NOTE SU EQUAZIONI DIFFERENZIALI LINEARI Armando Bazzani Meccanica Analitica e Modelli

57 Pages·2017·0.35 MB·Italian

NOTE SU EQUAZIONI DIFFERENZIALI LINEARI
Armando Bazzani
Meccanica Analitica e Modelli Numerici per la Fisica
Dipartimento di Fisica e Astronomia - Bologna
Soluzione delle equazioni differenziali lineare: il flusso di fase
Le equazioni differenziali lineari costituiscono un esempio importante per lo studio dei
sistemi dinamici in quanto da una parte `e possibile uno studio esaustivo delle propriet`a
delle soluzioni, e dall’altra sono uno strumento per l’approssimazione locale di un sistema
dinamico generico. Consideriamo il sistema di equazioni differenziali lineari a coefficienti
costanti
x˙ = Ax x ∈ Rn (1.0)
dove A ´e una matrice data n×n. La soluzione all’equazione (1.0) ´e determinata una volta
fissata la condizione iniziale x(0) = x . Infatti possiamo trasformare il sistema (1.0) in
0
un’equazione integrale
(cid:90) t
x(t) = x +A x(s)ds (1.1)
0
0
definendo la soluzione in forma implicita. L’esistenza della soluzione si pu´o dimostrare
utilizzando un principio di punto fisso per la mappa
(cid:90) t
y(cid:48)(t) = x +A y(s)ds
0
0
nello spazio delle funzioni y(t) regolari nell’intervallo [0,t]. Possiamo anche scrivere la
soluzione andando ad approssimare l’integrale secondo la definizione
x(∆t) = (I +A∆t+O(∆t)2)x
0
e iterando la procedura fino ad ottenere
x(t) = (I +A∆t+O(∆t)2)t/∆tx
0
Abbiamo quindi il limite per le matrici analogo a quanto succede per i numeri reali
x(t) = lim (I +A∆t)t/∆tx = exp(At)x (1.2)
0 0
∆t→0
1
Tale relazione risulta interessante in quanto consente di scrivere degli integratori numerici
per il sistema (1.0). Per il calcolo della soluzione possiamo anche sviluppare x(t) in serie
di Taylor a partire da t = 0
(cid:88) tk dkx (cid:88) tk
x(t) = x + (0) = x + Akx (1.3)
0 k! dtk 0 k! 0
k≥1 k≥1
dove abbiamo utilizzato l’eguaglianza
dkx
(0) = Akx
dtk 0
Riconosciamo nella (1.3) la funzione esponenziale della matrice At definita dalla serie
(cid:88)∞ tk
exp(At) = Ak
k!
k=0
`
E possibile dimostrare la convergenza della serie introducendo la norma nello spazio delle
matrici
(cid:107)Ax(cid:107)
(cid:107)A(cid:107) = sup (cid:107)Ax(cid:107) = sup
(cid:107)x(cid:107)=1 x(cid:54)=0 (cid:107)x(cid:107)
(dove (cid:107)x(cid:107) indica l’usuale norma Euclidea) e utilizzando il fatto che lo spazio delle matrici
´
´e completo (vedi Appendice). E facile dimostrare che se la dimensione dello spazio `e finita
anche la norma di una matrice `e finita.
Dalla definizione della norma, vale la relazione
(cid:88)∞ tk
(cid:107)exp(At)(cid:107) ≤ (cid:107)A(cid:107)k = exp((cid:107)A(cid:107)t)
k!
k=0
e l’esponenziale converge su tutto l’asse reale. Allora il vettore x(t) = exp(At)x `e ben
0
definito ed `e la soluzione dell’equazione differenziale (1.0). La definizione dell’operatore
esponenziale expAt definisce in modo naturale un gruppo abeliano di operatori che prende
il nome di flusso di fase associato all’equazione (1.0). Abbiamo infatti le propriet`a
exp(O) = I exp(At)exp(As) = exp(A(t+s))
La matrice A prende il nome di generatore del gruppo. Vale anche il viceversa: sia Φtx
un gruppo abeliano di operatori lineari dipendenti dal parametro reale, possiamo definire
il generatore del gruppo come
Φ∆t −I
lim = A
∆t→0 ∆t
e rappresentare il gruppo Φt = exp(At) in quanto la funzione x(t) = Φtx , x(t) risolve
0
l’equazione differenziale
d
x˙ = Φtx = Ax
0
dt
2
e si pu`o scrivere nelle forma exp(At)x e vale un teorema unicit`a per le soluzioni dell’equa-
0
zione differenziale (1.0). Abbiamo pertanto una corrispondenza tra gruppi di operatori
lineari abeliani e equazioni differenziali. Notiamo che la matrice S(t) = exp(At) che
caratterizza il flusso di fase `e soluzione dell’equazione matriciale
dS
= AS(t)
dt
Notiamo infine che la linearit´a dell’operatore (1.2) implica che ogni combinazione lineare
(a coefficienti costanti) di soluzioni dell’equazione differenziale (1.0) ´e ancora soluzione con
condizione iniziale data dalla combinazione delle condizioni iniziali con gli stessi coeffcienti.
Quindi esiste un isomorfismo tra lo spazio delle condizioni iniziali e lo spazio delle soluzioni.
In particolare quest’ultimo risulta uno spazio lineare di dimensione n ed esiste una base per
lo spazio delle soluzioni che corrisponde ad una base dello spazio delle condizioni inziali.
Per un calcolo esplicito dell’esponenziale di una matrice usiamo il seguente Lemma:
Se λ `e un autovalore della matrice A corrispondente all’autovettore v allora vale
exp(At)v = eλtv
ovvero eλt `e autovalore della matrice exp(At).
Se la matrice A `e diagonalizzabile mediante una trasformazione lineare T
A = TΛT−1
allora possiamo scrivere
exp(At) = TeΛtT−1
Tale risultato segue in quanto vale l’eguaglianza
TΛnT−1 = (TΛT−1)n = An
Abbiamo quindi un metodo del tutto generale per risolvere le equazioni lineari a coefficienti
costanti. Una volta risolta l’equazione caratteristica
det(λI −A) = 0
cerchiamo le soluzione dell’equazione (1.0) nella forma
x(t) = T exp(Λt)c (1.4)
dove T `e la matrice che diagonalizza A (ovvero la matrice formata dagli autovettori) e c `e
un vettore costante che viene definito dalle condizioni condizioni iniziali. Sostituendo nel
sistema (1.0) si ottiene
TΛexp(Λt)c = AT exp(Λt)c
3
che risulta soddisfatta qualunque sia la scelta di c se vale
(TΛ−AT) = 0 ⇒ A = TΛT−1
Infine imponendo le condizioni iniziali otteniamo il sistema
x(0) = Tc
che ammette una sola soluzione. L’equazione
(TΛ−AT) = 0 (1.5)
si definisce equazione di coniugazione per la matrice A. Tale equazione risolve il problema
ditrovareunatrasformazioneT cheriducaAallaformapiu´semplicepossible(inparticolare
alla forma diagonale), per cui il calcolo delle soluzioni dell’equazione (1.0) diventa banale.
La soluzione di tale problema richiede l’utilizzo di risultati algebrici circa le forme normali
delle matrici: diremo infatti che l’equazione di coniugazione consente di ridurre in forma
normalelamatriceA. RicordiamocheΛeAdevonoaverelostessopolinomiocaratteristico
P(ζ): infatti dalla relazione TΛT−1 = A si ricava
P(ζ) = det(Iζ −A) = det(Iζ −TΛT−1) = det(T)det(Iζ −Λ)det(T−1) = det(Iζ −Λ)
Pertanto gli autovalori della matrice A nel sistema (1.0) definiscono delle propriet`a in-
trinseche della dinamica che non dipendono dalle variabili scelte e consentono una clas-
sificazione a seconda del valore degli autovalori. Come `e noto se il polinomio ha tutte
radici distinte allora Λ `e una diagonale formata dagli autovalori e T `e la matrice le cui
colonne sono le coordinate degli autovettori corrispondenti (nel caso di autovalori multipli
la formula `e piu` complicata e avremo le forme normali di Jordan).
In generale, possiamo risolvere il sistema (1.0) cercando soluzioni nella forma esponenziale
x(t) = (cid:126)veλt con (cid:126)v vettore. La linearit´a del problema implica che date delle soluzioni
(cid:80)
x (t) ogni combinazione lineare c x (t) con c coefficienti reali o complessi, ´e ancora
k k k k k
soluzione. Il sistam di soluzioni x (t) si dice linearmente indipendente se l’equazione
k
(cid:88)
c x (t) = 0
k k
k
ha un’unica soluzione c = 0 per ogni t. Si dimostra che se le soluzioni sono lineramente
k
indipendenti per t = 0 lo sono anche per ogni tempo: possiamo avere una prova diretta
dal fatto che il flusso di fase ´e invertibile e quindi non altera la condizione di indipendenza
lineare se questa ´e valida al tempo t = 0. Si ottiene quindi
λ(cid:126)veλt = A(cid:126)veλt
che equivale all’equazione di autovalori ed autovettori per la matrice A. Se troviamo una
base di autovettori per la matrice A, avremo automaticamente un sistema di soluzioni
4
linearmente indipedenti e il calcolo dello costanti c consente di esprimere qualunque con-
k
dizione iniziale. Lo spettro degli autovalori della matrice A consente quindi di ottenere
informazione sul comportamento delle soluzioni. Per i sistemi meccanici (ovvero per cui
l’equazione di Newton) si dimostra che la matrice A si fattorizza A = JS con J matrice
antisimmetrica con J2 = −I ed S matrice simmetrica. Per i sistemi meccanici J ´e
(cid:18) (cid:19)
0 I
J =
−I 0
In tal caso Sx diventa il gradiente di una forma quadratica H = x·Sx/2.
Si dimostra il Lemma: se λ `e radice dell’equazione caratteristica Det[λI −JS] = 0 allora
anche −λ lo `e.
Consideriamo la seguente relazione
0 = Det[λI −JS] = Det[λI +SJ] = Det[J]2Det[−λI −JS] = Det[−λI −JS]
Il primo passaggio si ottiene prendendo la trasposta e quindi si usano le propriet´a di J.
Data quindi una matrice reale nella forma JS gli autovalori si presentano a coppie reali
o a coppie immaginari puri o come quadruple di opposti e complessi coniugati. Passando
alla matrice esponenziale expJS le propriet´a discusse si riflettono nella relazione
[exp(JS)]TJ exp(JS) = J (1.6)
che caratterizza la cosidette matrici simplettiche. La dimostrazione segue da
[exp(JS)]TJ = JJ−1exp(−SJ)J = J exp(−JS)
ricordando che T−1exp(A)T = exp(T−1AT).
Le matrici simplettiche sono importanti per lo studio dei flussi di fase dei sistemi meccanici
Hamiltoniani e hanno la propriet´a che se λ ´e autovalore anche λ−1 lo ´e.
Alcuni esempi e applicazioni
Si consideri la successione di numeri interi (successione di Fibonacci)
a = a +a n ≥ 1
n+1 n n−1
con a = 0 ed a = 1. Rappresentare in modo matriciale la ricorrenza e determinare
0 1
ln(a )
n
lim
n→∞ n
Possiamo scrivere la ricorrenza in forma matriciale
(cid:18) (cid:19) (cid:18) (cid:19)(cid:18) (cid:19)
a 1 1 a
n+1 = n
a 1 0 a
n n−1
5
per cui
(cid:18) (cid:19) (cid:18) (cid:19)n(cid:18) (cid:19)
a 1 1 1
n+1 =
a 0 1 0
n
Calcoliamo gli autovalori della matrice dal polinomio caratteristico
√
1± 5
λ2 −λ−1 = 0 λ =
±
2
per cui |λ | > 1 e |λ | < 1. Allora per n → ∞ la componente lungo l’autovettore
+ −
corrispondente a λ tende a zero, e solo la componente lungo l’autovettore associato a λ
− +
sopravvive. Per valutare il limite dobbiamo scomporre nella base degli autovettori
√ √
(cid:18) (cid:19) (cid:18) (cid:19)
(1+ 5)/2 (1− 5)/2
v = v =
+ 1 − 1
´
che risultano ortogonali. E facile vedere che
(cid:18) (cid:19)
1 1
= √ (v −v )
0 5 + −
per cui
1
a = √ (λn −λn)
n + −
5
e vale il limite √
lna 1+ 5
n
lim = lnλ = ln
+
n→∞ n 2
λ `e anche noto come il numero d’oro.
+
Un caso di notevole interesse sono i sistemi dinamici lineari associati alle matrici
stocastiche π : ovvero matrici con la propriet´a che 0 < π < 1 e
jk jk
(cid:88)
π = 1 (a)
jk
j
In base alle ipotesi fatti gli elementi π si possono interpetare come probabilit´a di tran-
ij
sizionetraunostatoj diuncertosistemaallostatoi. Talimatricisonolegateall’evoluzione
dei processi di Markov secondo l’equazione
(cid:88)
ρ (t+1) = π ρ (t) (b)
j jk k
k
Vale il Lemma: in un caso non denegere la matrice Π ha un autovalore unitario mentre
tutti gli altri sono in modulo strettamente minori di 1.
Dallapropriet´a(a)seguecheilvettore(1,1,...,1)´eautovettoredellamatricetrasposta
con autovalore unitario, allora anche la matrice Π ha un autovalore eguale a uno.
6
Inoltre abbiamo un integrale primo del moto lineare dato da
(cid:88)
ρ (t) = cost.
j
j
Da questo segue che: se Π ammette un autovettore ρ∗ con (cid:80) ρ∗ (cid:54)= 0 allora il corrispon-
j j
dente autovalore deve essere 1. In effetti dall’integrale primo del moto segue
(cid:88) (cid:88) (cid:88)
ρ∗ = π ρ∗ = λ ρ∗
k kj j k
k k,j k
che pu´o essere vera solo se λ = 1. Dal momento che tutti gli elementi di Π sono positivi
il primo quadrante risulta invariante per la dinamica (b) e il determinante deve essere
strettamente minore di 1. Per dimostrare tale affermazione ricordiamo che (ε ´e il simbolo
di Levi-Civita in n dimensioni)
(cid:88)
DetΠ = ε π ...π
j ,...,j 1,j n,j
1 n 1 n
j
che risulta un sottoinsieme in senso stretto di tutti i prodotti positivi che si ottegono
dall’espressione
(cid:32) (cid:33)
(cid:89) (cid:88)
π = 1
kj
j k
Il valore 1 si ottiene solo se Π ´e una matrice diagonale, il che vuol dire che coincide con
l’identit´a. Pertanto la dinamica (b) ´e una contrazione nell’insieme ottenuto intersecando
(cid:80)
il primo quadrante con sottospazio ρ = 1 e quindi ammette un punto fisso che cor-
k k
risponde all’autovettore di autovalore unitario. L’autovettore ρ∗ pu´o essere interpretato
come distribuzione di probabilit´a associata ad uno spazio discreto di eventi (in corrispon-
denza con le componenti). Seguono quindi le propriet´a:
(cid:80)
1) il sottospazio lineare definito dal piano v = 0 ´e un sottospazio invariante per la
k k
dinamica (b) ed ´e possibile dividere lo spazio vettoriale come somma diretta di tale
sottospazio ed lo spazio associato all’autovettore ρ∗;
(cid:80)
2) gli autovettori nel sottospazio v = 0 hanno autovalori strettamente minori di 1
k k
nei casi non degeneri.
Per provare questa affermazione dato un generico vettore v consideriamo separata-
mente le sue componenti positive v e negative −v (tali insiemi devono essere entrambi
k+ k−
non vuoti) e introduciamo la norma
(cid:88) (cid:88) (cid:88)
|v| = |v | = v + v
k k+ k−
k k+ k−
Vale che
(cid:12) (cid:12) (cid:32) (cid:33)
(cid:88)(cid:12)(cid:88) (cid:12) (cid:88) (cid:88) (cid:88)
|Πv| = (cid:12) π v (cid:12) < π v + v
(cid:12) jk k(cid:12) jk+ k+ k−
(cid:12) (cid:12)
j k j k+ k−
7
Scambiando le sommatorie e utilizzando la propriet´a (a) si ottiene
(cid:88) (cid:88)
|Πv| < v + v = |v|
k+ k−
k+ k−
dove, salvo matrici degeneri che hanno autovalore 1 con molteplicit´a maggiore di uno, il
segno di minore ´e stretto. Ne segue che dato un autovettore v∗
|Πv∗| = |λv∗| = |λ||v∗| < |v∗|
e quindi il corrispondente autovalore deve essere strettamente minore di 1 in modulo. Dal
principio di contrazione segue che iterando la dinamica (b) a partire da un qualunque
vettore, che si pu´o scrivere nella forma
u = cρ∗ +v
(cid:80) (cid:80)
con c = u e v appartenente al sottospazio invariante v = 0, abbiamo
j j k k
lim Πnu = cρ∗ + lim Πnv → cρ∗
n→∞ n→∞
La successione converge verso l’autovettore di autovalore unitario e la velocit´a di conver-
genza ´e proporzionale a |λ |n dove λ ´e il piu´ grande degli autovalori in modulo
max max
esclusa l’unit´a.
Associate alle matrici stocastiche vi sono le cosidette marici Laplaciane L con la propriet´a
(cid:88)
l = 0
jk
j
e |l | ≤ 1. Le matrici Laplaciane discretizzano l’operatore di Laplace e si richiede che gli
jk
autovalori siano tutti maggiori o uguali di zero: sicuramente 0 ´e autovalore dal momento
che vi ´e una relazione lineare tra le colonne. Le propriet´a delle matrici Laplaciane si
riconducono a quelle delle matrici stocastiche dimostrando che la matrice esponenziale
exp(tL) ´e una matrice stocastica. Si vede facilmente la propriet´a
(cid:88)
[exp(tL)] = 1
jk
j
Occorre dimostrare che gli elementi della matrice sono positivi; consideriamo il caso di
∆t (cid:28) 1, si ha
[exp(∆tL)] = δ +∆tl +O(∆t2)
jk jk jk
per cui l’asserto segue se l > 0 per j (cid:54)= k (per cui possono essere interpretati come tassi
jk
di transizione per unit´a di tempo). Da momento che vale una propriet´a di gruppo
n
exp(n∆tL) = (exp(∆tL)
8
allora il risultato precedente vale anche per un tempo finito. Nelle ipotesi precedenti
possiamo associare una matrice stocastica ad una matrice Laplaciana con
π = δ +l ∆t
jk jk jk
purch´e 0 < l ∆t < 1. In tal caso abbiamo che una matrice stocastica con gli stessi
jk
autovettori di L e autovalori collegati dalla relazione
λ = 1+λ
Π L
Dal momento che λ < 1 ad eccezzione dell’autovalore unitario, ne segue che tutti gli
Π
autovalori di L hanno parte reale negativa (salvo in casi degeneri) tranne il primo che
´e nullo. L’autovettore associato all’autovalore nullo ´e diretto nel primo quadrante e pu´o
essere interpretato come distribuzione di probabilit´a. Le matrici Laplaciane entrano nella
definizione della Master Equation
(cid:88)
ρ˙ = l ρ (d)
j jk k
k
utilizzata sia in Meccanica Statistica che in Meccanica Quantistica, in cui ρ (t) ha il
k
sinificato di probbalit´a di essere nello stato k al tempo t. La soluzione formale formale si
scrive
ρ(t) = exp(Lt)ρ(0)
(cid:80)
e vale che ρ ´e un integrale primo del moto per la dinamica (d) (conservazione della
j j
probabilit´a totale). Da quanto detto l’equazione (d) ammette uno stato stazionario che
pu´o essere interpretato come distribuzione di probabilit´a e risulta attrattivo
lim ρ(t) = ρ
s
t→∞
Notiamo infine che la composizione di due matrici stocastiche Π(1) e Π(2) ´e ancora una
matrice stocastica
(cid:88)(cid:88) (cid:88)
(1) (2) (2)
π π = π = 1
ij jk jk
i j j
Quindi le soluzioni della Master Equation si possono approssimare utilizzando la compo-
sizione di matrici stocastiche. In particolare possiamo approssimare
(cid:18) ∆t2 ∆tn(cid:19)
ρ(t+∆t) = I +L∆t+L2 +.....+Ln ρ(t)+O(∆tn+1)
2 n!
e comporre in successioni le matrici stocastiche al membro destro.
Matrici SimpletticheUnsecondocasodiinteressesonoisistemidiequazionidifferenziali
lienari associate alle matrici simplettiche caratterizzate dalla relazione
MTJM = J
9
con J matrice antisimmetrica con J2 = 1 di dimensione pari. Le matrici simplettiche
hanno la propriet´a che se λ ´e autovalore anche λ−1 lo ´e:
Dimostrazione
Det(λI −MT) = Det(λI +JM−1J) = Det(λI −M−1)
Quindi λ ´e autovalore di M e della sua inversa; ci´o ´e possibile solo se gli autovalori si
presentano a coppie (λ,λ−1).
´
E facile dimostrare che se poniamo
M = exp(JS) (b)
con S matrice simmetrica, allora M ´e simplettica. Vale infatti l’eguaglianza
exp(−SJ)J exp(JS) = J (cid:0)J−1exp(−SJ)J(cid:1)exp(JS) = J exp(−JS)exp(JS) = J
Ne segue che le matrici simplettiche rappresentano le soluzioni del sistema lineare
dM
x˙ = JSx ovvero = JSM(t)
dt
con S matrice simmetrica. Dalle definizioni `e facile riconoscere che il sistema precedente
ha la forma canonica di un sistema Hamiltoniano con H = (xSx)/2. Le propriet´a di
gruppo delle matrici simplettiche consentono di mantenere la forma canonica dell equazioni
del moto. Se infatti operiamo un cambio di variabili x = My con M trasformazione
simplettica, il sistema cambia mediante la relazione
y˙ = M−1JSMy
e dato che dalla condizione di simpletticit´a segue che M−1 = −JMTJ, abbiamo
(−JMTJ)JSM = JMTSM = JS(cid:48)
dove S(cid:48) = MTSM `e la nuova matrice simmetrica associata al Hamiltoniano del sistema.
Si pu´o mostrare che formalmente vale anche il viceversa: data una matrice M ´e sempre
possibile una rappresentazione esponenziale del tipo (b). In effetti ´e formalmente possibile
invertire l’equazione (b), nella forma
log(M) = log(I +(M −I)) = JS ⇒ −J log(M) = S
Dal momento che S ´e simmetrica deve accadere che
log(MT)J = −J log(M) ⇒ J log(MT)J = log(M)
Ora dalla condzione di simpletticit´a abbiamo −(JMTJ)M = I che implica M−1 =
−(JMTJ). Si verifica quindi la relazione
JMnJ = −(−JMJ)n
10