Plot ellipsoid from matrix

Question

I would like to plot an ellipsoid from a (3x3) matrix in matplotlib

I have this code:

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.patches import FancyArrowPatch
from mpl_toolkits.mplot3d import proj3d

fontSize = 14
lineWidth = 2.0
plt.rc('text', usetex=True)
plt.rc('font', family='serif')
plt.rcParams["figure.figsize"] = (10,8)
plt.rcParams['axes.linewidth'] = lineWidth 

class Arrow3D(FancyArrowPatch):
    def __init__(self, xs, ys, zs, *args, **kwargs):
        super().__init__((0,0), (0,0), *args, **kwargs)
        self._verts3d = xs, ys, zs

    def do_3d_projection(self, renderer=None):
        xs3d, ys3d, zs3d = self._verts3d
        xs, ys, zs = proj3d.proj_transform(xs3d, ys3d, zs3d, self.axes.M)
        self.set_positions((xs[0],ys[0]),(xs[1],ys[1]))

        return np.min(zs)

def main_axis(ax, s=1.0):
    hat_1 = np.array([s*1, 0, 0])
    hat_2 = np.array([0, s*1, 0])
    hat_3 = np.array([0, 0, s*1])
    arrow_prop_dict = dict(mutation_scale=10, arrowstyle='-|>', color='k', linewidth=1)

    ax.add_artist(Arrow3D([0, s*1],[0, 0],[0, 0], **arrow_prop_dict))
    ax.text(*(s*hat_1),'x',fontsize=14) 

    ax.add_artist(Arrow3D([0, 0],[0, s*1],[0, 0], **arrow_prop_dict))
    ax.text(*(s*hat_2),'y',fontsize=14) 
    
    ax.add_artist(Arrow3D([0, 0],[0, 0],[0, s*1], **arrow_prop_dict))
    ax.text(*(s*hat_3),'z',fontsize=14) 
    

def plot(mat):      

    eigVal, eigVec = np.linalg.eig(mat)
    
    # Create figure
    fig = plt.figure()
    ax = fig.add_subplot(projection = '3d')
    
    main_axis(ax, s=1.15)
        
    nDiv = 30 
    theta = np.linspace(0, np.pi, nDiv).reshape(nDiv, 1)   
    phi = np.linspace(0, 2*np.pi, nDiv).reshape(1,1,nDiv) 
    theta, phi = np.meshgrid(theta, phi)

    x = eigVal[0]*(np.sin(theta) * np.cos(phi))
    y = eigVal[1]*(np.sin(theta) * np.sin(phi))
    z = eigVal[2]*(np.cos(theta))
   
    ax.plot_surface(x, y, z,  rstride=4, cstride=4, color='b', alpha=0.7)
    
    for vec in eigVec.T:
        drawvec = Arrow3D([0, vec[0]], [0, vec[1]], [0, vec[2]],
                          mutation_scale=20, lw=2, arrowstyle="-|>", color='r')
        ax.add_artist(drawvec)

mat = np.array([ [0.5, 0.0, 0.25],
                [0.0, 0.4, 0.0],
                [0.25, 0.0, 0.1]                
                ])

plot(mat)

How can I make the ellipsoid be described by the eigenvalues and the eigenvectors?

score 0 · Answer 1 · answered Jun 21 '23 at 20:24

I would to it that way

    u0=eigVec[:,0][None,None,:]
    u1=eigVec[:,1][None,None,:]
    u2=eigVec[:,2][None,None,:]

    pts = eigVal[0]*(np.sin(theta) * np.cos(phi))[:,:,None]*u0 + eigVal[1]*(np.sin(theta) * np.sin(phi))[:,:,None]*u1 + eigVal[2]*(np.cos(theta))[:,:,None]*u2
    x=pts[:,:,0]
    y=pts[:,:,1]
    z=pts[:,:,2]

    ax.plot_surface(x, y, z,  rstride=4, cstride=4, color='b', alpha=0.7)

Note: the fact that your ellipsoid is very flat, just comes from the fact that one of the eigenvalues is very small. I obviously did not solve that (that comes from the data); I just multiplied the said eigenvalue by 10 to have a more visible dimension (I could and should have used another input matrix whose eigenvalues are all big enough. But I was just to lazy for that :D).

So, what my code solve is the fact that your ellipsoid has the size of your eigenvalues, but have the canonical x,y,z axes. You can easily see that your can never make use of eigen vectors at all (for the ellipsoid, I mean. You use them for the arrows). So, no chance that the ellipsoid axes are aligned with the arrows.

The points of an ellipsoid are, as you almost coded,
λ₀.sin(θ).cos(φ)×u₀ + λ₁.sin(Θ).sin(φ))×u₁ + λ₂.cos(Θ))×u₂

Where λᵢ are the axes sizes, and uᵢ are the axes direction. Here, obviously, as you can guess both by my choice of λ and u, and by the fact that this is just your code, transformed into a math formula, with vector addition, axes sizes are the eigenvalues, and axes direction the eigenvectors.

What you did by plotting with just x,y,z, is the same as what my code would do, by replacing u₀ by [1,0,0], u₁ by [0,1,0], and u₂ by [0,0,1], aka i, j, k, aka the canonical axes.

Note that I had to play with some broadcasting/axis to have the computation done on all combination of θ and φ, and on all 3 canonical dimensions i,j,k (x,y,z).

Plot ellipsoid from matrix

1 Answers1