Is there a visualization tool that shows matrices transforming (multiplying) ellipses?

Question

Matrices send ellipses to ellipses. Is there any online visualization tool where I can enter a matrix $M$ and see a plot of ellipse $e$ being sent to $Me$? This would be very helpful to visualize classes of matrices such as symmetric, skew symmetric, orthogonal, etc.

Ideally, the tool should show ellipses in 2D and ellipsoids in 3D.

It doesn't have to be a dedicated tool: If there's a way to do this with a more general tool (e.g. Wolfram Alpha), that will be helpful as well.

Answer 1

Here's a rudimentary Python script to get a 2-D visualization.

import numpy as np
import matplotlib.pyplot as plt
from numpy import sin, cos, pi

A = np.array([        # matrix to visualize
    [1, 2],
    [3, 4]
])
angle = 30 * (pi/180) # ellipse rotation angle
axis_1 = 5            # length of axis 1
axis_2 = 3            # length of axis 2

th_range = np.linspace(0,2*pi,61)

# get points in a circle
x,y = cos(th_range), sin(th_range)
pts = np.vstack([x,y])

# stretch and rotate to form elipse
pts = pts * np.array([[axis_1],[axis_2]])
rot = np.array([
    [cos(angle), -sin(angle)],
    [sin(angle), cos(angle)]
])
pts = rot @ pts

# get points after transformation is applied

pts_after = A @ pts


plt.plot(*pts)
plt.plot(*pts_after)

plt.legend(['Before Transformation', 'After Transformation'])
plt.show()

Sample result:

---

More sophisticated version:

def fix_aspect(height=5, margin=0.1):
    ax = plt.gca()
    fig = plt.gcf() # added this
    ax.set_aspect(1)
    box = ax.get_tightbbox(fig.canvas.get_renderer()) # added single argument
    box_width, box_height = box.x1 - box.x0, box.y1 - box.y0
    aspect = box_height / box_width
    width = height / aspect
#     fig = plt.gcf()
    fig.set_figwidth(width)
    fig.set_figheight(height)
    fig.tight_layout(pad=margin)

def draw_trans(A, angle = 0, axis_1 = 1, axis_2 = 1):
    
    th_range = np.linspace(0,2*pi,61)

    # get points in a circle
    x,y = cos(th_range), sin(th_range)
    pts = np.vstack([x,y])

    # stretch and rotate to form elipse
    pts = pts * np.array([[axis_1],[axis_2]])
    rot = np.array([
        [cos(angle), -sin(angle)],
        [sin(angle), cos(angle)]
    ])
    pts = rot @ pts

    # get points after transformation is applied
    pts_after = A @ pts
    
    # get arrows
    arrow_1 = np.array([cos(angle),sin(angle)])*axis_1
    arrow_2 = np.array([cos(angle + pi/2),sin(angle + pi/2)])*axis_2
    arrow_3 = A@arrow_1
    arrow_4 = A@arrow_2

    # plots
    plt.plot(*pts)
    plt.plot(*pts_after)


    plt.annotate("", xy=arrow_1, xytext=(0, 0), arrowprops=dict(
            arrowstyle="-|>", linewidth = 2, linestyle = '-', color = 'blue'))
    plt.annotate("", xy=arrow_2, xytext=(0, 0), arrowprops=dict(
            arrowstyle="-|>", linewidth = 2, linestyle = ':', color = 'blue'))
    plt.annotate("", xy=arrow_3, xytext=(0, 0), arrowprops=dict(
            arrowstyle="-|>", linewidth = 2, linestyle = '-', color = 'green'))
    plt.annotate("", xy=arrow_4, xytext=(0, 0), arrowprops=dict(
            arrowstyle="-|>", linewidth = 2, linestyle = ':', color = 'green'))
    
    plt.legend(['Before Transformation', 'After Transformation'])
    fix_aspect()
    
    plt.show()
    
A = np.array([        # matrix to visualize
    [1, 2],
    [3, 4]
])
draw_trans(A, axis_2 = 2)

Result:

Answer 2

Here's a 3D plotter, albeit with limited functionality. It only can use the unit sphere as its starting ellipsoid, and it shows what happens to the transformation's principal axes or to the $i,j,k$ axes (but not to a user-specified set).

def draw_principle_3d(A, mode = 'principal'):
    lincount = 20
    u, v = np.mgrid[0:2*pi:lincount*2j, 0:pi:lincount*1j] 
    # u, v = np.mgrid[0:pi:20j, 0:pi:10j]
    X = cos(u)*sin(v)
    Y = sin(u)*sin(v)
    Z = cos(v)

    pts = np.array([x_vals,y_vals,z_vals]).transpose()
    pts_new = np.dot(pts, A.T).transpose()

    X_new,Y_new,Z_new = pts_new

    U,s,VT = la.svd(A)

    fig,ax = plt.subplots(figsize = (10,10),subplot_kw = {"projection":"3d"})
    ax.plot_wireframe(X, Y, Z, alpha = .1, color = 'tab:blue')
    ax.plot_wireframe(X_new, Y_new, Z_new, alpha = .1, color = 'tab:orange')

    org = np.zeros(3)
    
    if mode == 'principal':
        arrows_before = VT
        arrows_after = (s*U).T
    else:
        arrows_before = np.eye(3)
        arrows_after = A.T
    
    
    for row in arrows_before:
        a = Arrow3D(*zip(org,row), mutation_scale=20, 
                    lw=2, arrowstyle="-|>", color='blue')
        ax.add_artist(a)

    for row in arrows_after:
        a = Arrow3D(*zip(org,row), mutation_scale=20, 
                    lw=2, arrowstyle="-|>", color='green')
        ax.add_artist(a)
    
    plt.show()

Running this in the default mode, as in the code below, produces the picture below.

A = 0.7 * np.diag([1,2,3])
Q1,_ = la.qr(np.random.randn(3,3)) # random rotation
Q2,_ = la.qr(np.random.randn(3,3)) # random rotation

A = Q1@A@Q2

draw_principle_3d(A)

Result:

Running this in the "ijk" mode, as in the code below, produces the picture below.

A = 0.7 * np.diag([1,2,3])
Q1,_ = la.qr(np.random.randn(3,3)) # random rotation
Q2,_ = la.qr(np.random.randn(3,3)) # random rotation

A = Q1@A@Q2

draw_principle_3d(A, mode = 'ijk')