I'm trying to understand a function, and its name isn't really helping me (I've found resources similar to this and this but this is either R or not exactly "ELLIPLOT" and I'm not sure I fully understand the concept of an elliplot).
I know what a Cholesky decomposition is but I can't understand the sort of link it has with ellipses. I guess I'm looking for resources or explanation as to what pan of mathematics this is related to. I've encountered this function in a Linear Programming course but we haven't done much more than writing it and I'm a bit confused as to why we did.
def ellipplot(P,gamma=None,xc=None):
"""Matlab ELLIPPLOT Internal function for plotting ellipsoid"""
if xc is None:
xc = np.zeros(len(P))
if gamma is None:
gamma = 1
P=P/gamma
P=np.linalg.cholesky(P)
theta = np.linspace(-3.14, 3.14, 1000)
z = np.array([np.cos(theta), np.sin(theta)])
x = np.linalg.inv(P)@z
for n in range(x.shape[1]):
x[:,n]=x[:,n]+xc
xe = x[0,:]
ye = x[1,:]
return xe, ye
xe, ye = ellipplot(
np.array(
[[2, 1],
[1, 2]]
)
)
plt.plot(xe, ye)
plt.show()
This is the matlab version of the function, I believe.
The ellipse is defined by $$ p_{11}x_1^2+p_{22}x_2^2+p_{12}x_1x_2+p_{21}x_1x_2=1\,. $$ In matrix/vector notation this is $$\tag{1} x^\top Px=1 $$ where $$ P=\left(\begin{matrix}p_{11}&p_{12}\\p_{21}&p_{22}\end{matrix}\right)\,,~~x=\left(\begin{matrix}x_1\\x_2\end{matrix}\right)\,. $$ To plot the points that satisfy (1) the code takes the Cholesky decomposition (square root) of $P=LL^\top$ so that (1) can be written as $$\tag{2} y^\top y=1\,~~\text{ where }y=L^\top x\,. $$ Clearly, $$ y=\left(\begin{matrix}\cos\theta\\\sin\theta\end{matrix}\right)\,\quad\quad\theta\in[-\pi,\pi] $$ satisfy (2). Therefore, to find the vector $x$ we onyly have to invert $y=L^\top x$ and that's what the code does.
Strictly speaking the code inverts $L$ (not $L^\top$) but I checked that it plots an ellipse.