Modul 2 Formiddag

Modul 2 Formiddag#

Introduktion til SVD#

Singular Value Decomposition (SVD) af reelle matricer baserer sig på følgende sætning:

Enhver reel $m \times n$ matrix $A$ kan skrives på SVD-formen

\[A = U\,\Sigma\, V^T\]

hvor $V$ er en ortogonal $n\times n$ matrix, $\Sigma$ er en $m\times n$ diagonal matrix, og $U$ er en ortogonal $m\times m$ matrix. Elementerne i $\Sigma$ kaldes singulærværdierne, de er alle $\geq 0$ og står i ikke-stigende rækkefølge.

I første del af Modul 2 (”formiddag”) undersøger vi grundigt $2\times 2$ matricer med fuld rang. Bemærk af der ved SVD af en $2\times 2$ matrix kun indgår $2\times 2$ matricer. Her er et eksempel på SVD af en $2\times 2$ matrix:

\[\begin{split} \begin{bmatrix} 1 & -2 \\ -2 & -2 \end{bmatrix} = \begin{bmatrix} -\frac{4\sqrt{5}}{15} & \frac{3\sqrt{5}}{10} \\ \frac{2\sqrt{5}}{15} & \frac{3\sqrt{5}}{5} \end{bmatrix} \begin{bmatrix} 3 & 0 \\ 0 & 2 \end{bmatrix} \begin{bmatrix} -\frac{2\sqrt{5}}{5} & \frac{\sqrt{5}}{5} \\ -\frac{\sqrt{5}}{5} & -\frac{2\sqrt{5}}{5} \end{bmatrix} \end{split}\]

Opsætning#

Vi kalder de python pakker og procedurer vi får brug for nedenfor:

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.patches import Circle
from ipywidgets import interact, FloatSlider, Layout
from IPython.display import clear_output
from scipy.optimize import brentq

Eksistens af løsning.#

Lad $A$ være en reel 2x2 matrix med fuld rang. Vi skal vise et vigtigt resultat på vejen mod SVD: Der findes en ortonormal basis $\,(\boldsymbol{v_1},\boldsymbol{v_2})\,$ for $\Bbb R^2$ som opfylder at basisvektorernes billeder $\,A\boldsymbol{v_1}\,$ og $\,A\boldsymbol{v_2}\,$ er ortogonale.

For vilkårligt $\,t\in\Bbb R\,$ betragter vi vektorerne

\[\begin{split} \boldsymbol{x_1(t)} = \begin{bmatrix} \cos(t) \\ \sin(t) \end{bmatrix} \quad \text{og} \quad \boldsymbol{x_2(t)} = \begin{bmatrix} -\sin(t) \\ \cos(t) \end{bmatrix}\,. \end{split}\]

Gør rede for at $\,\boldsymbol{x_1(t)}\,$ og $\,\boldsymbol{x_2(t)}\,$ er ortogonale for alle $t.$

Vis at der findes et $\,t_0\in \left[0,\frac {\pi}2\right]$ der gør billedvektorerne $\, \boldsymbol{y_1(t_0)} = A\,\boldsymbol{x_1(t_0)}\,$ og $\,\boldsymbol{y_2(t_0)} = A \,\boldsymbol{x_2(t_0)}\,$ ortogonale. Vink: Gør rede for at funktionen $\,s(t)=\boldsymbol{y_1(t)}\cdot \boldsymbol{y_2(t)}\,$ er kontinuert, overvej dens værdi i intervallets endepunkter og brug Bolzanos sætning.

Med tallet $t_0$ og funktionen $\,s(t)\,$ fra forrige spørgsmål: Vis at $\,s(t_0+ \frac{\pi}{2})=0\,.$
Vink: Udnyt at der for alle $t$ gælder $\,\cos(t+\frac{\pi}{2})=-\sin(t)\,$ og $\,\sin(t+\frac{\pi}{2})=\cos(t)\,.$

Forklar at vi hermed har løst opgaven.

Geometrisk løsning af SVD for 2x2 matrix $\,A\,$ med fuld rang#

Vi opstiller her en procedure “enhedscirkel_ellipse” som visualiserer transformationen af enhedscirklen ved en given 2x2 matrix $\,A\,:$

def enhedscirkel_ellipse(A, slider_width='80%', figsize=(14, 7)):
    """
    Opretter en interaktiv visualisering af enhedscirklen og dens billede under matrix A.
    Parametre:
    - A: 2x2 numpy array, transformationsmatrix
    - slider_width: bredde på slider (f.eks. '80%')
    - figsize: tuple med figurstørrelse (bredde, højde)
    """
    # Basis- og billedfunktioner
    def x1(t): return np.array([np.cos(t), np.sin(t)])
    def x2(t): return np.array([-np.sin(t), np.cos(t)])
    def y1(t): return A @ x1(t)
    def y2(t): return A @ x2(t)

    # Forudberegn ellipsepunkter
    ts = np.linspace(0, 2*np.pi, 361)
    ellipse_pts = np.column_stack([A @ x1(tt) for tt in ts]).T

    # Plot-funktion
    def plot_vectors(t):
        clear_output(wait=True)
        fig, (ax1, ax2) = plt.subplots(1, 2, figsize=figsize)
        for ax in (ax1, ax2):
            ax.set_aspect('equal')
            ax.grid(True, linestyle='--', linewidth=0.4, alpha=0.5)
        # Fælles aksebegrænsninger
        all_pts = np.vstack([ellipse_pts, x1(t), x2(t), y1(t), y2(t)])
        lim = 1.15 * np.max(np.abs(all_pts))
        for ax in (ax1, ax2):
            ax.set_xlim(-lim, lim)
            ax.set_ylim(-lim, lim)

        # Venstre: enhedscirkel
        circle = Circle((0, 0), radius=1, edgecolor='lightgrey', facecolor='none', linewidth=1.5, zorder=0)
        ax1.add_patch(circle)
        ax1.plot(np.cos(ts), np.sin(ts), color='lightgrey', linewidth=1.5, zorder=0)
        # Højre: ellipse
        ax2.plot(ellipse_pts[:, 0], ellipse_pts[:, 1], color='lightgrey', linewidth=1.5, zorder=0)

        # Funktion til pile og labels
        def draw_arrow(ax, v, color, label, z):
            ax.quiver(0, 0, v[0], v[1], angles='xy', scale_units='xy', scale=1,
                      color=color, width=0.012, zorder=z)
            ax.text(v[0], v[1], f'  {label}{tuple(np.round(v, 2))}',
                    va='bottom', ha='left', fontsize=9, color=color)

        # Tegn vektorer
        draw_arrow(ax1, x1(t), 'tab:blue', 'x₁', 3)
        draw_arrow(ax1, x2(t), 'tab:blue', 'x₂', 3)
        draw_arrow(ax2, y1(t), 'tab:red',  'y₁', 4)
        draw_arrow(ax2, y2(t), 'tab:red',  'y₂', 4)

        ax1.set_title('Enhedscirklen S')
        ax2.set_title('Ellipsen A.S')
        fig.suptitle(rf"$t={t:.2f}$   (træk slideren)", fontsize=14)
        plt.show()

    # Slider
    slider = FloatSlider(
        value=0,
        min=0,
        max=2*np.pi,
        step=0.01,
        description='t',
        continuous_update=True,
        readout_format='.2f',
        layout=Layout(width=slider_width)
    )
    interact(plot_vectors, t=slider)

I den følgende celle skal du indskrive den matrix du vil afprøve.
Prøv $\,A=\left[\begin{array}{c}1 & 2\\-1 & 1\end{array}\right]\,$ for grafisk løsning med slider.

Kør cellen, så proceduren enhedscirkel_ellipse aktiveres.

A = np.array([[1,2],
              [-1,1]])              
enhedscirkel_ellipse(A)

Vi tager udgangspunkt i situation og betegnelserne som i øvelse 1. Gå således frem:

Find med slideren et $\,t_0,\,$ således at $\,\boldsymbol{y_1(t_0)}\,$ og $\,\boldsymbol{y_2(t_0)}\,$ er ortogonale.

Sæt $\,\boldsymbol{v_1}=\boldsymbol{x_1(t_0)}\,$ og $\,\boldsymbol{v_2}=\boldsymbol{x_2(t_0)}\,$

Lad $\sigma_1$ og $\sigma_2$ være længderne $\,\sigma_1=\left\|\boldsymbol{y_1(t_0)}\right\|\,$ og $\,\sigma_2=\left\|\boldsymbol{y_2(t_0)}\right\|\,$

Sæt $\,\boldsymbol{u_1}=\frac{1}{\sigma_1}\,\boldsymbol{y_1(t_0)}\,$ og $\,\boldsymbol{u_2}=\frac{1}{\sigma_2}\,\boldsymbol{y_2(t_0)}\,$

Gør rede for at $\,V=\left[\boldsymbol{v_1}\,\,\boldsymbol{v_2}\right]\,$ og $\,U=\left[\boldsymbol{u_1}\,\,\boldsymbol{u_2}\right]\,$ er ortogonale matricer.

Lad $\Sigma$ være diagonalmatricen $\,\Sigma=\begin{bmatrix}\sigma_1&0\\0&\sigma_2\end{bmatrix},$ og gør rede for at $\,A\,V=U\,\Sigma\,$

Gør rede for og verificér i det givne eksempel, at $\,A=U\,\Sigma\,V^T\,$

Konvention: Det er fast konvention i SVD, at man ønsker singulærværdierne i faldende orden, det vil her sige, at $\,\sigma_1 \geq \sigma_2\,.$
Meget bekvemt fandt vi i første eksempel et $\,t_0\in \left[0,\frac {\pi}2\right],\,$ hvor $\,\left\|\boldsymbol{y_1(t_0)}\right\|>\left\|\boldsymbol{y_2(t_0)}\right\|\,.$
Hvis det modsatte gælder, trækker man blot slideren videre til næste ortogonalitetspunkt. Det opstår ved $\,t$-vinklen $\,t_0+\frac{\pi}2\,$ hvor længderne byttes om.

Prøv for eksempel $\,A=\left[\begin{array}{c}1 & -1\\0 & 2\end{array}\right]\,$ eller $\,A=\left[\begin{array}{c}1 & -1\\0 & -2\end{array}\right]\,.$

Håndregningsversion hvor du har mulighed for at regne med eksakte tal:

Udfør SVD på $\,A=\left[\begin{array}{c}5 & 3\\0 & 4\end{array}\right]\,.$

Vi opstiller nu en procedure SVD_2x2_matrix som samler punkterne 1-7, og som du kan bruge til at teste dine 2x2 matricer:
Bemærk specielt at vi sørger for at $\,\|\boldsymbol{y_1(t_0)}\|\,$ er større end eller lig med $\,\|\boldsymbol{y_2(t_0)}\|\,.$

def SVD_2x2_matrix(A, print_output = False):

    # Basis vektorer som funktion af t
    def x1(t): return np.array([np.cos(t), np.sin(t)])
    def x2(t): return np.array([-np.sin(t), np.cos(t)])
    def y1(t): return A.dot(x1(t))
    def y2(t): return A.dot(x2(t))

    # Dotprodukt-funktion
    f = lambda t: np.dot(y1(t), y2(t))

    # Find rod i (0, π/2)
    t0 = brentq(f, 1e-6, np.pi/2 - 1e-6)

    # Sørg for ||y1|| >= ||y2||
    if np.linalg.norm(y1(t0)) < np.linalg.norm(y2(t0)):
        t0 += np.pi/2

    # Så kan vi sætte værdierne direkte:
    sigma1 = np.linalg.norm(y1(t0))
    sigma2 = np.linalg.norm(y2(t0))
    Sigma = np.diag([sigma1, sigma2])
    V = np.column_stack([x1(t0), x2(t0)])
    U = np.column_stack([y1(t0)/sigma1, y2(t0)/sigma2])

    # printer outputtet, hvis tilvalgt
    if print_output:
        print("t0 =", t0)
        print("Σ =\n", Sigma)
        print("V =\n", V)
        print("U =\n", U)

    return  t0, U, Sigma, V

I den følgende celle afprøves proceduren SVD_2x2_matrix. Som output får du i denne rækkefølge: $\,t_0,\, U, \,\Sigma\,$ og $\,V\,.$

Indtast din test matrix $A$ og tjek i en følgende celle, at vi opnår den ønskede dekomposition $\,A = U\, \Sigma\,V^T\,.$

A = np.array([[1,2],
              [-1,1]])
t0, U, Sigma, V = SVD_2x2_matrix(A, print_output=True) 

Afprøvning af særlige matricer#

Find på simple test-eksempler på særlige 2x2 matricer:

Symmetrisk matrix, dvs. typen $\,\left[\begin{array}{c}a & b\\b & c\end{array}\right]$
Skævsymmetrisk matrix af typen $\,\left[\begin{array}{c}a & b\\-b & a\end{array}\right]$
Ortogonal matrix (du kan f.eks. bruge matrix $U$ lagret fra foregående øvelse)

Test matricerne skrives i følgende celle. Når du kører cellen aktiverer du begge procedurer enhedscirkel_ellipse og SVD_2x2_matrix.

Eksperimenter via slideren og studér det algebraiske output på matricer af type 1-3. Diskutér hvad der er det særlige!

A = np.array([[1,2],
              [-1,1]])
enhedscirkel_ellipse(A)
_ = SVD_2x2_matrix(A, print_output=True)

Gilbert Strangs visualisering af SVD.#

Verdens mest kendte underviser i lineær algebra, Gilbert Strang fra MIT, bliver i et 5 minutters interview på youtube spurgt om, hvad der er smukkest i lineær algebra, og hans svar er SVD. Spurgt ind til, hvad der er smukt i SVD, svarer han, at enhver matrix kan opfattes som en rotation efterfulgt af strækning i ortogonale retninger og afsluttet med en ny rotation. Eftersom rotationer ikke ændrer formen af et geometrisk objekt, er det altså midterdelen, strækningerne, der er det afgørende. Og de er givet ved singulærværdierne. Det skal vi nu studere nærmere via en figur fra en artikel af Gilbert Strang:

Gilbert

I stedet for at betragte dens transformation af enhedscirklen i ét hug, som ved at gå direkte fra venstre side til højre side i figuren, fanger vi den i tre trin bestemt ved højresiden af SVD-en $\,A=U\,\Sigma\,V^T\,.$ Inden vi går igang, minder vi om at en ellipse kan opfattes som en cirkel der er blevet strakt i to ortogonale retning. På standard-parameterform er ellipsen $e$ derfor givet ved $$\,e(t)=\big (a\cos(t),b\sin(t)\big )\,,\,\,t\in [0,2\pi]\,$$ hvor $a$ og $b$ er strækningsfaktorerne også kaldet halvakserne.

Her er et eksempel du kan arbejde med:

2x2SVD

Transformationen af enhedscirklen $C$ er naturligvis uafhængig af den parameterfremstilling der vælges for den. Vi vælger derfor, inspireret af Gilbert Strang-figuren, standard-paramaterfremstillingen for $c_v$ i det nye koordinatsystem bestemt ved basisvektorerne $\boldsymbol{y_1(t_0)}$ og $\boldsymbol{y_2(t_0)}$ (dvs. søjlerne i $V$).

Opgaverne i dette afsnit er tænkeopgaver. Til opgaverne er givet køreklare implementeringer, der plotter de relevante parametriseringer.

Hvilken parameterfremstilling $c_e$ har $C$ så i det sædvanlige koordinatsystem bestemt ved basisvektorerne $\boldsymbol{v_1}$ og $\boldsymbol{x_1(t_0)}\,?$
Vink: $V$ er basisskiftematricen, der skifter $v$-koordinater til sædvanlige $e$-koordinater.

def plot_parametric(func):
    """
    Plotter en parametrisk kurve defineret af func(t) for t i [0, 2π],
    med en farvegradient langs stien.

    Parameters:
    - func: lambda function of t, returning (x, y)
    """
    N = 500
    t = np.linspace(0, 2*np.pi, N)
    x, y = func(t)

    # Normalize t for colormap
    colors = plt.get_cmap("viridis")(np.linspace(0, 1, N))

    for i in range(N-1):
        plt.plot(x[i:i+2], y[i:i+2], color=colors[i])

    plt.axis("equal")
    plt.show()

V = np.array([[1,-1],
              [1, 1]]) * np.sqrt(2)/2
cv = lambda t: np.array([np.cos(t), np.sin(t)])  # funktion for standardparametrisering af enhedscirkel i v basis
ce = lambda t: V @ cv(t)  # enhedscirkel i standard basis

plot_parametric(ce)

Udfør den første transformation $V^T\,c_e\,.$ Kommentér resultatet.

# første transformation
T1 = lambda t: V.T @ ce(t)
plot_parametric(T1)

Udfør den næste transformation $\,\Sigma\,(V^T\,c_e)\,.$ Hvilken geometrisk figur er der tale om?

Sigma = np.diag([2,1])*np.sqrt(10)

# anden transformation
T2 = lambda t: Sigma @ T1(t)
plot_parametric(T2)

Kan du allerede nu se hvorfor $\,A\, c_e=U\,(\Sigma\,(V^T\,c_e))\,$ er en ellipse? Hvad er dens halvakser?

I den følgende celle vises det endelige billede $A\,C$.

A = np.array([[5,3],
              [0,4]])
plot_parametric(lambda t: A @ ce(t))

Geometrisk løsning ved 3x3 matricer: “Apollonios’ metode”#

Transformation af enhedskugle og tangentplan.

En 3x3 matrix $A$ med fuld rang er givet ved

\[\begin{split} A=\begin{bmatrix} 1 & 2 & 1 \\ 1 & 1 & -1 \\ -1 & 0 & 0 \end{bmatrix} \end{split}\]

Målet med opgaven er at finde en SVD for $A$. Til det formål indfører vi en lineær afbildning $f: \mathbb{R}^3\to \mathbb{R}^3$ givet ved $f(\boldsymbol{v})=A\, \boldsymbol{v}$ for alle $\boldsymbol{v}$ i $\mathbb{R}^3$.

Vi ønsker at finde SVD af $A$ ved at undersøge hvordan $f$ transformerer enhedskuglen $K = \{ (x,y,z) \in \mathbb{R}^3 \;\;|\;\; x^2 + y^2 + z^2 = 1 \}$. Billedet $f(K)$ betegner vi med $E$.

Til rådighed for opgaven haves følgende viden:

En lineær afbildning givet ved en $3\times 3$ matrix med fuld rang er tangentplan-bevarende. Det betyder, at tangentplanen for $K$ i et punkt $P$ på $K$ afbildes i tangentplanen for $E$ i $f(P)$.
Der findes 6 punkter på $E$, hvor stedvektoren for hvert punkt er ortogonal på tangentplanen for $E$ i punktet.
De 6 punkter optræder i 3 par af modsatte punkter.
Hvis man vælger et punkt fra hvert af de tre par, vil de tre punkters stedvektorer være indbyrdes ortogonale, og med dem kan vi opbygge SVD-en.

Spørgsmål 1#

Du skal finde tre punkter på $E$, hvis stedvektorer er indbyrdes ortogonale. Lad $P=(x_0,y_0,z_0)$ være et vilkårligt punkt på $K$. Ved følgende fremgangsmåde opstiller vi et ulineært ligningssystem i $x_0,y_0,z_0$, som vi vil løse med scipy’s fsolve, der kan finde rødderne til en funktion numerisk. Det opstilles senere i opgaven i funktionen equations_LS(P).

from scipy.optimize import fsolve
A = np.array([[1,2,1],
              [1,1,-1],
              [-1,0,0]])

Gør rede for, at tangentplanen $TK_P$ for $K$ i $P$ er givet ved

\[ x_0\cdot x+y_0\cdot y+z_0\cdot z=1\;. \]

Antag nu, at $x_0$ er foskellig fra 0 og fortsæt.

Gør rede for, at vektorerne

\[\begin{split} \boldsymbol{r}_1=\begin{bmatrix} -\frac{y_0}{x_0}\\1\\0\end{bmatrix}\quad \text{ og }\quad\boldsymbol{r}_2=\begin{bmatrix} -\frac{z_0}{x_0}\\0\\1\end{bmatrix} \end{split}\]

er lineært uafhængige retningsvektorer for $TK_P$.
Vink: Bestem en parameterfremstilling for $TK_P$.

I følgende celle defineres en funktion, der giver $\boldsymbol{r}_1$ og $\boldsymbol{r}_2$ givet $P$. Dette er en hjælpefunktion, der skal bruges senere.

def TK_vectors(P):
    "Returnerer to lineært uafhængige retningsvektorer for tangentplanen til K i punktet P"
    x0, y0, z0 = P
    r1 = np.array([-y0/x0, 1, 0])
    r2 = np.array([-z0/x0, 0, 1])
    return r1, r2

Lad nu $Q$ være billedpunktet $Q=f(P)$ (Q=A@P), og lad $TE_Q$ være tangentplanen for $E$ i $Q$.

Angiv to lineært uafhængige retningsvektorer $\boldsymbol{s}_1$ og $\boldsymbol{s}_2$ for $TE_Q$. Vink: Overvej igen parameterfremstillingen for $TK_P$.

I følgende celle defineres endnu en hjælpefunktion, der giver $\boldsymbol{s}_1$ og $\boldsymbol{s}_2$ givet punktet $P$.

Færdiggør funktionen. Vink: Gør evt. brug af funktionen TK_vectors fra før.

def TE_vectors(P):
    "Returnerer to lineært uafhængige retningsvektorer for tangentplanen til E i punktet Q=f(P)"
    x0, y0, z0 = P
    r1, r2 = TK_vectors(P)
    # EGEN KODE: Definer s1 og s2

    return s1, s2

Kan du begrunde, at det faktisk gælder, at $\boldsymbol{s}_1$ og $\boldsymbol{s}_1$ er lineært uafhængige?

Med dette fås, at $TE_Q$ kan parametriseres som $Q+t_1\boldsymbol{s}_1+t_2\boldsymbol{s}_2$.

Gør rede for at stedvektoren $OQ$ for $Q$ er ortogonal på $TE_Q$, netop når de nedenstående 3 ligninger er opfyldte.

$\quad OQ \cdot \boldsymbol{s}_1 = 0$
$\quad OQ \cdot \boldsymbol{s}_2 = 0$
$\quad x_0^2 + y_0^2 + z_0^2 - 1 = 0$

$OQ$ bestemmes af $P$ gennem afbildningen $f$, og $\boldsymbol{s}_1$ og $\boldsymbol{s}_2$ bestemmes tilsvarende af $P$ som vist i TE_vectors. Dermed afhænger alle tre ligninger i lingingssystemet alene af den ubekendte $P=(x_0,y_0,z_0)$.

I den følgende celle defineres en funktion, der beregner venstresiderne af de 3 ovenstående ligninger. Det er dermed den funktion, som vi ønsker at finde roden af, da alle højresiderne er 0.

Færdiggør funktionen. Vink: Brug evt. Numpys funktioner np.dot(v,w), np.linalg.norm(v). I Python skrives $x^a$ som x**a.

def equations_LS(P):
    "Returnerer værdierne af venstresiderne af de tre ligninger i det ulineære ligningssystem"
    x0, y0, z0 = P
    Q = A @ P
    s1, s2 = TE_vectors(P)
    # definer værdierne af venstresiderne i ligningssystemet
    eq1_LS =
    eq2_LS =
    eq3_LS =
    return (eq1_LS, eq2_LS, eq3_LS)

I den følgende celle finder vi rødderne til ligningssytemet med fsolve. fsolve tager en funktion og et startgæt som input. Derfor defineres også tilfældige startgæt på $K$. Løsningerne der findes er således også punkter $(x_0,y_0,z_0)$ på $K$.

Kør cellen.

solutions = [] # liste til løsninger

# forskellige startgæt til fsolve
guesses = [np.random.normal(size=3) for _ in range(50)] # 50 tilfældige gæt
guesses = [g/np.linalg.norm(g) for g in guesses] # normaliserer gættene, så de ligger på K

# finder løsningerne for hvert startgæt
for guess in guesses:
    # undgå x0 tæt på 0 (antagelsen!) - startgættet skippes
    if abs(guess[0]) <= 0.05:  
        continue

    # finder løsning med fsolve
    sol = fsolve(equations_LS, guess)

    # tjekker om løsningen er ny (ikke allerede fundet)
    if not any(np.allclose(sol, s) for s in solutions):
        # gemmer løsningen
        solutions.append(sol)

# printer løsningerne
i = 0
for solution in solutions:
    print(f'solutions[{i}] =', solution)
    i += 1

Verificer, at der er fundet 6 løsninger.

Hvis nej, så kan du have været uheldig med de tilfældige startgæt eller vi kan have taget fejl om antagelsen om at $x_0$ er forskellig fra $0$. Genoptag undersøgelsen for tilfældet $x_0=0$ og find de resterende.

Er der fundet 6 løsninger, gør så rede for, at spørgsmål 1 er besvaret.

Spørgsmål 2#

Bestem en SVD for $A$. Vi skal finde de tre $3\times 3$ matricer $\Sigma$, $V$ og $U$, hvor $\Sigma$ er en diagonalmatrix, og $V$ og $U$ er ortogonale.
Gå frem efter følgende opskrift.

Vælg vha. dine fundne løsninger til equations_LS 3 punkter på $E$, hvis stedvektorer er indbyrdes ortogonale, og bestem de tre vektorers længder. Vink: Brug evt. np.linalg.norm(v).

e1 = A @ solutions[]
e2 = A @ solutions[]
e3 = A @ solutions[]
print("e1 =", e1, "with length", np.linalg.norm(e1))
print("e2 =", e2, "with length", np.linalg.norm(e2))
print("e3 =", e3, "with length", np.linalg.norm(e3))

Opstil de tre længder i ikke-stigende rækkefølge og kald dem $\sigma_1$, $\sigma_2$ og $\sigma_3$. Dan diagonalmaticen $\Sigma$. Vink: Brug evt. np.diag([a,b,c])

sigma1 = np.linalg.norm()
sigma2 = np.linalg.norm()
sigma3 = np.linalg.norm()
Sigma = np.diag([sigma1, sigma2, sigma3])
print("Sigma =\n", Sigma)

Find de tre stedvektorer til punkter på $K$ hvis billeder har længderne $\sigma_1$, $\sigma_2$ og $\sigma_3$, og kald dem hhv. $\boldsymbol{v}_1$, $\boldsymbol{v}_2$ og $\boldsymbol{v_3}$. Dan matricen $V$.
Vink: Brug evt. np.columnstack([u,v,w]).

v1 = solutions[]
v2 = solutions[]
v3 = solutions[]
V = np.column_stack([v1, v2, v3])
print("V =\n", V)

Bestem vektorerne $\boldsymbol{u}_1=\frac{1}{\sigma_1}f(\boldsymbol{v}_1)$, $\boldsymbol{u}_1=\frac{1}{\sigma_2}f(\boldsymbol{v}_2)$ og $\boldsymbol{u}_3=\frac{1}{\sigma_3}f(\boldsymbol{v}_3)$. Dan matricen $U$.

u1 = (1/sigma1) * (A @ v1)
u2 = (1/sigma2) * (A @ v2)
u3 = (1/sigma3) * (A @ v3)
U = np.column_stack([u1, u2, u3])
print("U =\n", U)

Verificer, at $A=U\,\Sigma\,V^T$.

print("A =\n",A)
print("U@Sigma@V.T =\n",U @ Sigma @ V.T)