Esercizio - Cholesky della matrice di rigidezza di una catena di tre masse

Teoria: Fattorizzazione LU, pivoting e CholeskyRisolvere $A\mathbf x=\mathbf b$ con Gauss equivale a fattorizzare $A=LU$ ($L$ triangolare inferiore con $1$ sulla diagonale e i moltiplicatori $l_{ik}=\frac{a_{ik}^{(k)}}{a_{kk}^{(k)}}$ sotto la diagonale, $U$ la matrice a scala finale) e poi risolvere $L\mathbf y=\mathbf b$ (in avanti) e $U\mathbf x=\mathbf y$ (all'indietro). Esiste (ed è unica) se i minori principali di testa $\det A_k$, $k<n$, sono $\ne0$. Costo $\frac23n^3$ flop per la fattorizzazione, $2n^2$ per ogni termine noto; $\det A=\prod u_{kk}$. Con un pivot piccolo l'algoritmo è instabile: il pivoting parziale (scambio di righe che porta in posizione di pivot l'elemento di modulo massimo della colonna) dà $PA=LU$ con $|l_{ik}|\le1$ ed esiste per ogni $A$ invertibile. Se $A$ è simmetrica definita positiva esiste $A=LL^T$ (Cholesky), con $l_{jj}=\sqrt{a_{jj}-\sum_{k<j}l_{jk}^2}$ e $l_{ij}=\frac{a_{ij}-\sum_{k<j}l_{ik}l_{jk}}{l_{jj}}$; costo $\frac13n^3$, senza pivoting.Fattorizzazione LU, pivoting e Cholesky →. Il sistema massa-molla è anche in Sistemi lineari n per n - matrice fondamentale, esponenziale e variazione delle costantiPer y' = A(t)y + g(t) di dimensione n una matrice fondamentale Φ(t) ha per colonne n soluzioni indipendenti (D: lo sono se e solo se det Φ(t₀) ≠ 0); l'integrale generale dell'omogeneo è y = Φ(t)c. Con A costante: se Ah = λh, e^(λt)h è soluzione; autovalori complessi λ = μ + iω danno parte reale e immaginaria di e^(λt)h; autovalore non diagonalizzabile: e^(λt)(h₂ + t h₁); l'esponenziale e^(tA) = Σ tᵏAᵏ/k! risolve y(t₀) = y₀ con y = e^((t−t₀)A)y₀. Variazione delle costanti (D): y = Φ(t)Φ(t₀)⁻¹y₀ + Φ(t)∫Φ(s)⁻¹g(s)ds, per A costante y = e^((t−t₀)A)y₀ + ∫e^((t−s)A)g(s)ds. Sistemi massa-molla x'' = −Kx: modi normali, pulsazioni ω_j² = autovalori di K, battimenti. Stabilità: Re λ < 0 per ogni autovalore dà y → 0.Sistemi lineari n per n - matrice fondamentale, esponenziale e variazione delle costanti →.

Testo. Tre carrelli sono in fila su una guida liscia, collegati tra loro e alle due pareti da quattro molle uguali di rigidezza $k=2000$ N/m (parete - $m_1$ - $m_2$ - $m_3$ - parete). Sul terzo carrello agisce una forza costante $F=10$ N verso destra. Siano $x_1,x_2,x_3$ gli spostamenti statici dei carrelli dall'equilibrio.

1. Scrivere le equazioni di equilibrio $K\mathbf x=\mathbf F$ e mostrare che $K$ è simmetrica definita positiva.
2. Calcolare la fattorizzazione di Cholesky di $K$.
3. Calcolare gli spostamenti e l'energia elastica immagazzinata.
4. Cosa accade se si tolgono le due molle di parete?

Soluzione

1. Equilibrio. Il carrello $i$ riceve da ciascuna molla adiacente una forza proporzionale all'allungamento relativo. Per $m_1$: $-kx_1+k(x_2-x_1)=0$; per $m_2$: $-k(x_2-x_1)+k(x_3-x_2)=0$; per $m_3$: $-k(x_3-x_2)-kx_3+F=0$. In forma matriciale $$K\mathbf x=\mathbf F,\qquad K=k\begin{pmatrix}2&-1&0\\-1&2&-1\\0&-1&2\end{pmatrix},\qquad\mathbf F=\begin{pmatrix}0\\0\\10\end{pmatrix}\text{ N}.$$ $K$ è simmetrica. Definita positiva: minori di testa (della matrice $K_0=K/k$) $2$, $2\cdot2-1=3$, $\det K_0=2(4-1)-(-1)(-2)=6-2=4$, tutti positivi (positivi anche moltiplicati per $k>0$). Fisicamente, $\mathbf x^TK\mathbf x=2U>0$ per ogni spostamento non nullo, perché è il doppio dell'energia elastica.

2. Cholesky di $K_0$ (la scala $k$ si toglie e si rimette alla fine: $K=kK_0=(\sqrt kL)(\sqrt kL)^T$).

• Colonna 1: $l_{11}=\sqrt2=1{,}414214$; $l_{21}=\frac{-1}{\sqrt2}=-0{,}707107$; $l_{31}=\frac0{\sqrt2}=0$.
• Colonna 2: $l_{22}=\sqrt{2-l_{21}^2}=\sqrt{2-\frac12}=\sqrt{1{,}5}=1{,}224745$; $l_{32}=\frac{-1-l_{31}l_{21}}{l_{22}}=\frac{-1}{1{,}224745}=-0{,}816497$.
• Colonna 3: $l_{33}=\sqrt{2-l_{31}^2-l_{32}^2}=\sqrt{2-\frac23}=\sqrt{\frac43}=1{,}154701$.

$$L=\begin{pmatrix}1{,}414214&0&0\\-0{,}707107&1{,}224745&0\\0&-0{,}816497&1{,}154701\end{pmatrix},\qquad\det K_0=(1{,}414214\cdot1{,}224745\cdot1{,}154701)^2=2^2=4 .$$ Coincide con $\det K_0=4$ calcolato sopra ✓. Il fatto che i radicandi siano $2,\ \frac32,\ \frac43$ (rapporti $\frac{\det A_j}{\det A_{j-1}}=\frac21,\frac32,\frac43$) conferma la definita positività. La matrice è tridiagonale e $L$ è bidiagonale: la struttura a banda si conserva.

3. Spostamenti. Si risolve $K_0\mathbf x=\frac{\mathbf F}k=(0,0,0{,}005)^T$ m, cioè con $\frac{10}{2000}=0{,}005$.

• Avanti $L\mathbf y=(0,0,0{,}005)^T$: $y_1=0$; $y_2=\frac{0-l_{21}y_1}{l_{22}}=0$; $y_3=\frac{0{,}005}{1{,}154701}=0{,}004330$.
• Indietro $L^T\mathbf x=\mathbf y$: $x_3=\frac{y_3}{l_{33}}=\frac{0{,}004330}{1{,}154701}=0{,}00375$; $x_2=\frac{y_2-l_{32}x_3}{l_{22}}=\frac{0+0{,}816497\cdot0{,}00375}{1{,}224745}=0{,}0025$; $x_1=\frac{y_1-l_{21}x_2}{l_{11}}=\frac{0{,}707107\cdot0{,}0025}{1{,}414214}=0{,}00125$.

$$\mathbf x=(1{,}25;\ 2{,}5;\ 3{,}75)\ \text{mm}.$$ Controllo con la prima equazione: $2000\,(2\cdot0{,}00125-0{,}0025)=0$ ✓; la terza: $2000\,(-0{,}0025+2\cdot0{,}00375)=2000\cdot0{,}005=10$ N ✓. Energia elastica: $U=\frac12\mathbf F^T\mathbf x=\frac12\cdot10\cdot0{,}00375=0{,}01875$ J (18,75 mJ).

4. Senza le molle di parete la matrice diventa $K'=k\begin{pmatrix}1&-1&0\\-1&2&-1\\0&-1&1\end{pmatrix}$, con $\det K'=0$: Cholesky darebbe $l_{11}=1$, $l_{21}=-1$, $l_{22}=\sqrt{2-1}=1$, $l_{32}=-1$ e un terzo radicando $l_{33}^2=1-l_{32}^2=0$, e l'algoritmo si ferma. È la matrice di una struttura labile: tutti i carrelli possono traslare insieme senza deformare le molle (autovettore $(1,1,1)$ con autovalore $0$), e il sistema $K'\mathbf x=\mathbf F$ non ha soluzione (la forza netta non è nulla) o ne ha infinite. Il tentativo di Cholesky è quindi anche un test fisico di vincolo.

Risultati: $\mathbf x=(1{,}25;\ 2{,}5;\ 3{,}75)$ mm, $U=18{,}75$ mJ, $\det K=4k^3$; radicandi di Cholesky $2,\frac32,\frac43$.

Errori tipici. Dimenticare la scala $k$ (spostamenti 2000 volte troppo grandi); scambiare $L$ e $L^T$ nella seconda sostituzione; saltare il controllo di definitezza prima di applicare Cholesky.