Facoltµa di Ingegneriamocenni/dispense_07_08.pdfL’importante teorema di Hartman-Grobman...

Universita degli Studi di Siena

Facolta di Ingegneria

Modellistica e simulazione

Corso di Laurea Specialistica in Ingegneria Informatica

Corso di laurea Specialistica in Ingegneria Gestionale

Dispense del corso

Prof. Chiara MocenniIn collaborazione con Sergio Trimboli

A. A. 2007/2008

Indice

Introduzione 1

1 I sistemi lineari 2

2 Biforcazioni di Sistemi Dinamici 52.1 Biforcazione Nodo-Sella . . . . . . . . . . . . . . . . . . . . . . . . . 62.2 Biforcazione Transcritica . . . . . . . . . . . . . . . . . . . . . . . . . 72.3 Biforcazione Pitchfork Supercritica . . . . . . . . . . . . . . . . . . . 72.4 Biforcazione Pitchfork Subcritica . . . . . . . . . . . . . . . . . . . . 82.5 Biforcazione di Hopf . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

3 Mappe e Punti Fissi 123.1 La mappa logistica . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13

4 Mappe di Poincare 154.1 Stabilita lineare di orbite periodiche . . . . . . . . . . . . . . . . . . 16

5 Flussi sul cerchio e simulazioni 175.1 Pendolo semplice . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 175.2 Pendolo sovrasmorzato . . . . . . . . . . . . . . . . . . . . . . . . . . 205.3 Beat phenomenon . . . . . . . . . . . . . . . . . . . . . . . . . . . . 215.4 Lucciole . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

6 Il modello di Lorenz 31

7 Controllo del caos 407.1 Controllo di un processo caotico in continua evoluzione . . . . . . . . 417.2 Targeting . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 427.3 Controllo del sistema di Lorenz . . . . . . . . . . . . . . . . . . . . . 43

Bibliografia 48

i

Introduzione

[... Negli ultimi venti anni, si e parlato molto dei metodi e dei risultati matematiciche hanno portato alla definizione di caos deterministico. Questi risultati sono statiottenuti nell’ambito di quel settore della Matematica noto come Teoria qualitativadei sistemi dinamici e sono stati stimolati dall’esigenza di rappresentare, mediantemodelli matematici, i sistemi reali che evolvono nel tempo come il moto dei pianeti, leoscillazioni di un pendolo, il flusso delle correnti atmosferiche, lo scorrere piu o menoregolare dell’acqua in un fiume, il numero di insetti che anno dopo anno popolano unacerta regione, l’andamento giornaliero dei prezzi delle azioni nei mercati finanziarie cosı via. L’apparente contraddizione (o paradosso) contenuto nel termine caosdeterministico, ha molto incuriosito anche il pubblico dei non specialisti. I modellimatematici di tipo deterministico vengono in genere associati all’idea di fenomeniregolari, prevedibili, che si ripetono nel tempo, mentre il termine caotico viene riferitoa situazioni caratterizzate da assenza di regole e da imprevedibilita. La scopertadel caos deterministico spezza questa dicotomia, in quanto mostra come modellimatematici deterministici (cioe privi di ogni elemento aleatorio nelle equazioni cheli definiscono) sono in grado di generare andamenti estremamente complessi, sottomolti aspetti imprevedibili, tanto da risultare quasi indistinguibili da sequenze dieventi generati attraverso processi aleatori... ] Gian-Italo Bischi, Rosa Carini, LauraGardini e Paolo Tenti.

1

Capitolo 1

I sistemi lineari

Si consideri il sistema dinamico descritto dalla seguente equazione differenziale or-dinaria:

dx

dt=x = f(x), x ∈ Rn. (1.1)

Una soluzione di (1.1) e una funzione x(x0, t) a valori in Rn dipendente da t econdizione iniziale x(0) = x0. Possiamo interpretare la parte destra di (1.1) comeun flusso φt : U → Rn, U ⊆ Rn, dove φt(x) = φ(x, t) e una funzione smooth, chesoddisfa la (1.1), definita per ogni x in U e t in qualche intervallo I = (a, b) ⊆ R.

Se il sistema e lineare e autonomo, cioe non dipende esplicitamente dal tempo t,f(x) = Ax, dove A e una matrice a coefficienti costanti.

In questo caso la soluzione della (1.1) e

x(x0, t) = etAx0, (1.2)

dove etA e la matrice definita come segue:

etA = [I + tA +t2

2!A2 + · · ·+ tn

n!An + . . . ]. (1.3)

Una soluzione generale x(t) della (1.1) puo essere ottenuta comme combinazione din soluzioni indipendenti

x1(t), . . . , xn(t)

:

x(t) =n∑

j=1

cjxj(t),

dove le n costanti incognite possono essere calcolate dalle condizioni iniziali.Se A ha n autovettori linearmente indipendenti vj , j = 1, . . . , n, allora possiamo

utilizzare come base dello spazio delle soluzioni le funzioni

xj(t) = eλjtvj ,

2

3

dove λj e l’autovalore associato all’autovettore vj . Per gli autovalori complessi senzamolteplicita, λj = αj±iβj , con autovettori vR±ivI , possiamo prendere come coppiadi soluzioni reali indipendenti ad essi associate

xj = eαjt(vRcos(βt)− vIsin(βt)),

xj+1 = eαjt(vRcos(βt) + vIsin(βt)).

Denotiamo conX(t) = [x1(t), . . . , xn(t)]

la matrice fondamentale delle soluzioni che ha le n soluzioni nelle colonne. Le colonnexj(t), j = 1, . . . , n di X(t) formano una base dello spazio delle soluzioni della (1.1),e inoltre

etA = X(t)X(0)−1.

Nel caso in cui gli autovalori abbiano molteplicita, e possibile trovare una tra-sformazione invertibile che diagonalizza A o almeno permette di metterla in formanormale di Jordan. L’equazione (1.1) diventa

y = Jy,

dove J = T−1AT e x = Ty e le colonne sono gli autovettori generalizzati di A.L’esponenziale di A diventa allora

etA = TetJT−1,

dove gli esponenziali sono calcolati per le tre matrici 2 x 2 in forma di Jordan:

A =(

λ1 00 λ2

)etA =

(eλ1t 00 eλ2t

);

A =(

α −ββ α

)etA = eαt

(cos(βt) −sin(βt)sin(βt) cos(βt)

);

A =(

λ 01 λ

)etA = eλt

(1 0t 1

).

4

Consideriamo ora il caso in cui n = 2 e

A =(

a bc d

).

Il determinante della equazione caratteristica e

λ2 − τλ + ∆ = 0,

doveτ = traccia(A) = a + d;

∆ = det(A) = ad− bc.

Abbiamo

λ1,2 =τ ±√τ2 − 4∆

2, ∆ = λ1λ2, τ = λ1 + λ2.

Questo permette di classificare la stabilita delle soluzioni per ogni combinazionedegli autovalori e di visualizzare questa classificazione nel piano (∆, τ) come segue:

• ∆ < 0: gli autovalori hanno segno opposto. Il punto di equilibrio x = 0 e unpunto di sella;

• ∆ > 0: gli autovalori hanno lo stesso segno o sono complessi coniugati. Il puntodi equilibrio e un nodo se τ2−4∆ > 0 e una spirale se τ2−4∆ < 0. La parabolaτ2 − 4∆ = 0 e il caso limite in cui si trovano le stelle e i nodi degeneri. Lastabilita dei nodi e delle spirali e determinata da τ : se τ < 0 tutti e due gliautovalori hanno parte reale negativa, quindi il punto di equilibrio e stabile; seinvece τ > 0 si hanno nodi e spirali instabili. I centri si trovano nell’asse τ = 0;

• ∆ = 0: c’e almeno un autovalore nullo, cio significa che l’origine non e un puntodi equilibrio isolato: esiste almeno una linea di punti di equilibrio o addiritturaun piano di punti di equilibrio.

Esempio: Studiare il sistema x = x + y; y = 4x− 2y, (x0, y0) = (2,−3).

Capitolo 2

Biforcazioni di Sistemi Dinamici

Un punto di equilibrio di un sistema dinamico di ordine n e detto iperbolico se tuttigli autovalori della sua linearizzazione si trovano fuori dall’asse immaginario.L’importante teorema di Hartman-Grobman stabilisce che il phase portrait (cioe ildiagramma di tutte le traiettorie del sistema dinamico nello spazio delle fasi) vicinoad un punto iperbolico e topologicamente equivalente al phase portrait del sistemalinearizzato. Topologicamente equivalente significa che esiste un omeomorfismo (cioeuna deformazione continua con inversa continua) che trasforma un phase portraitlocale nell’altro, mappa traiettorie in traiettorie e mantiene la direzione temporale.In particolare, le proprieta di stabilita del punto di equilibrio sono completamentecaratterizzate dalla linearizzazione.Un phase portrait viene detto strutturalmente stabile se la sua topologia non puoessere modificata da una perturbazione arbitrariamente piccola del campo vettoriale.Se analizziamo le proprieta di un sistema dinamico in funzione dei parametri checompaiono nelle equazioni, puo accadere che vi siano dei cambiamenti nella strutturaqualitativa delle soluzioni per qualche valore dei parametri. Questi cambiamentivengono chiamati biforcazioni e i parametri che danno loro luogo prendono il nomedi valori di biforcazione.Definizione 1. Dato il sistema dinamico

x = fµ(x), x ∈ Rn, µ ∈ Rk, (2.1)

un valore µ0 di (2.1) per cui il flusso di (2.1) non e strutturalmente stabile e unvalore di biforcazione di µ.Definizione 2. Il luogo dei punti nel piano (µ, x) degli insiemi invarianti di (2.1) sichiama diagramma di biforcazione.Di seguito vengono elencati i tipi piu comuni di biforcazioni locali in termini dellaloro forma normale. La forma normale di una biforcazione e una rappresentazioneprototipale di tutte le biforcazioni di quel tipo. Per chiarezza esaminiamo sistemidi ordine 1. Nelle equazioni che seguono r e un parametro che puo essere positivo,negativo o nullo.

5

2.1. BIFORCAZIONE NODO-SELLA 6

2.1 Biforcazione Nodo-Sella

x = r + x2. (2.2)

−6 −4 −2 0 2 4 6−15

−10

−5

0

5

10

15

20

Figura 2.1: Biforcazione Nodo-Sella. r0 = 0.

Figura 2.2: Diagramma di biforcazione: biforcazione nodo-sella.

2.2. BIFORCAZIONE TRANSCRITICA 7

2.2 Biforcazione Transcritica

x = rx− x2. (2.3)

−3 −2 −1 0 1 2 3−3

−2.5

−2

−1.5

−1

−0.5

0

0.5

1

1.5

2

Figura 2.3: Biforcazione Transcritica. r0 = 0.

Figura 2.4: Diagramma di biforcazione: biforcazione transcritica.

2.3 Biforcazione Pitchfork Supercritica

x = rx− x3. (2.4)

2.4. BIFORCAZIONE PITCHFORK SUBCRITICA 8

−2 −1.5 −1 −0.5 0 0.5 1 1.5 2−3

−2

−1

0

1

2

3

Figura 2.5: Biforcazione Pitchfork Supercritica. r0 = 0, α = −1, ω = 1, a = −1, b =0.

2.4 Biforcazione Pitchfork Subcritica

x = rx + x3. (2.5)

−2 −1.5 −1 −0.5 0 0.5 1 1.5 2−3

−2

−1

0

1

2

3

Figura 2.6: Biforcazione Pitchfork Subcritica. r0 = 0.

2.5 Biforcazione di Hopf

In sistemi di ordine almeno 2 si incontra anche un altro tipo di biforcazione: laBiforcazione di Hopf. Supponiamo che il sistema linearizzato abbia una coppia

2.5. BIFORCAZIONE DI HOPF 9

Figura 2.7: Diagramma di biforcazione: biforcazioni pitchfork supercritica esubcritica.

di autovalori complessi coniugati per il valore critico µ0 = 0 del parametro µ:

λ = µ + iω, λ = µ− iω

La forma normale della biforcazione di Hopf e la seguente:

x = µx− ωy + (−x− by)(x2 + y2) + O(|x|5, |y|5),y = ωx + µy + (bx− y)(x2 + y2) + O(|x|5, |y|5). (2.6)

L’equazione (2.6) in coordinate polari (x = rcos(θ), y = rsin(θ)) diventa:

r = µr − r3 + O(|r|5),θ = ω + br2 + O(|r|4). (2.7)


x ’ = mu x − omega y + ( − x − b y) (x2 + y2)y ’ = omega x + mu y + (b x − y) (x2 + y2)

omega = 1b = 1

mu = − 1

−3 −2 −1 0 1 2 3

−3

−2

−1

0

1

2

3

x

y

Figura 2.8: Biforcazione di Hopf Supercritica. µ < 0.


omega = 1b = 1

mu = 0

−3 −2 −1 0 1 2 3

−3

−2

−1

0

1

2

3

x

y

Figura 2.9: Biforcazione di Hopf Supercritica. µ = 0.



omega = 1b = 1

mu = 1

−3 −2 −1 0 1 2 3

−3

−2

−1

0

1

2

3

x

y

Figura 2.10: Biforcazione di Hopf Supercritica. µ > 0.

Figura 2.11: Diagramma di biforcazione: biforcazione di Hopf supercritica.

Figura 2.12: Diagramma di biforcazione: biforcazione di Hopf subcritica.

Capitolo 3

Mappe e Punti Fissi

Una mappa e descritta da una equazione del tipo

xn+1 = f(xn), (3.1)

dove f e una funzione continua con derivata continua, f : R→ R. Vengono chiamatemappe le equazioni alle differenze, le iterazioni, e in generale le equazioni in cui iltempo e discreto.

Definizione. x∗ si dice punto fisso per la mappa (3.1) se:

f(x∗) = x∗. (3.2)

Per la (3.2), se per qualche n accade che xn = x∗, allora xn+1 = f(xn) = f(x∗) = x∗,cioe l’orbita di (3.1) rimane uguale a x∗ per tutte le iterazioni successive.

Per studiare la stabilita di x∗ si considera un’orbita vicina xn = x∗ + ηn e sianalizza il fatto che essa sia attratta o respinta da x∗ stessa, cioe si si va a vedere sela perturbazione ηn cresce o si riduce. Sostituendo nella (3.1) otteniamo:

x∗ + ηn+1 = xn+1 = f(x∗ + ηn) = f(x∗) + f ′(x∗)ηn + O(η2n). (3.3)

Poiche f(x∗) = x∗, l’equazione (3.3) diventa

ηn+1 = f ′(x∗)ηn + O(η2n). (3.4)

Se trascuriamo i termini in η2 otteniamo una equazione lineare in η con autovalore,detto anche in questo caso moltiplicatore, λ = f ′(x∗).

La soluzione della mappa lineare si trova iterando a partire da un valore inizialeη0:

ηn = λnη0. (3.5)

Si distinguono tre casi:

i) |λ| = |f ′(x∗) < 1. In questo caso ηn → 0 per n → ∞ e il punto fisso elinearmente stabile.

12

3.1. LA MAPPA LOGISTICA 13

ii) |λ| = |f ′(x∗) > 1. In questo caso ηn → ∞ per n → ∞ e il punto fisso elinearmente instabile.

iii) |λ| = |f ′(x∗) = 1. In questo caso il punto fisso e marginalmente stabile e la suastabilita e determinata dai termini non lineari.

3.1 La mappa logistica

Esempio. Determinare la stabilita dei punti fissi della mappa logistica:

xn+1 = rxn(1− xn). (3.6)

Assumiamo x ≥ 0. I punti fissi della mappa sono: x∗1 = 0 e x∗2 = 1− 1/r. x∗1 e unpunto fisso per ogni valore di r, mentre x∗2 lo e solo per r ≥ 1.

La stabilita dipende dal moltiplicatore f ′(x∗) = r − 2rx∗.

• x∗1 = 0. Poiche f ′(x∗1) = r, l’origine e stabile per r < 1 ed instabile per r > 1.

• x∗2. f ′(x∗2) = 2−r, quindi x∗2 e stabile per −1 < (2−r) < 1, cioe per 1 < r < 3,e instabile per r > 3.

Per r = 3 f ′(x∗2) = −1 cioe r raggiunge un valore critico superato il qualeil punto x∗2 perde la sua stabilita e le traiettorie del sistema cominciano adoscillare indefinitamente tra due valori.

Questo tipo di perdita di stabilita prende il nome di biforcazione flip. Quelloche succede al sistema quando avviene una biforcazione flip e che un ciclo di periodoT perde la sua stabilita dando luogo ad un ciclo di periodo 2T . Per capire questofatto osserviamo che un ciclo di periodo 2T e un ciclo di periodo T per la secondaiterata della mappa f2(x) = f(f(x)). Poiche in questo caso f(x) e un polinomioquadratico, f2(x) e un polinomio di grado 4 che ha quindi 4 intersezioni con labisettrice.

La biforcazione flip puo manifestarsi in ogni tipo di mappe.Il ciclo di periodo 2 rimane stabile per 3 < r < 3.449. Oltre questo valore

anch’esso perde la propria stabilita dando luogo ad un nuovo ciclo di periodo 4T .Per comprendere questo fenomeno basta analizzare la quarta iterata della mappaf(x) o la seconda iterata della mappa f2(x).

Questo fenomeno si ripete per valori sempre crescenti di r, fino a che per r > 3.57si osserva un comportamento aperiodico denominato caos deterministico.

E’ interessante osservare che la mappa logistica ha la seguente proprieta:

limn→∞rn − rn−1

rn+1 − rn= 4.669,

cioe esiste una costante universale δ = 4.669 chiamata costante di Feigenbaumche fissa il rapporto che lega valori consecutivi della costante di biforcazione.

3.1. LA MAPPA LOGISTICA 14

Figura 3.1: Diagramma di biforcazione della mappa logistica.

Il diagramma di biforcazione della mappa logistica e mostrato nella seguentefigura:

Siti web:- http://www.lboro.ac.uk/departments/ma/gallery/doubling/- http://math.la.asu.edu/ chaos/logistic.html

Capitolo 4

Mappe di Poincare

Sia S una superficie di dimensione n − 1 trasversale al flusso definito da x = f(x).Sia inoltre xk la k − esima intersezione del flusso con la superficie S. Si definisceMappa di Poincare la mappa seguente:

xk+1 = P (xk),

dove P : S → S rappresenta la successione delle intersezioni successive del flussocon S (vedi figura 4.1).

Figura 4.1: Mappa di Poincare.

Per definizione, se x∗ e un punto fisso della mappa di Poincare, allora P (x∗) = x∗.Cio significa che la traiettoria del sistema e una curva chiusa.

Le mappe di Poincare riducono il problema dell’esistenza di una curva chiusa perun sistema dinamico all’esistenza di un punto fisso per una mappa di dimensioneinferiore.

Si consideri ad esempio il sistema r = r(1 − r2); θ = 1. Se prendiamo S ugualeall’asse x positiva, il primo ritorno su S avviene dopo t = 2π. La mappa di Poincarer1 = P (r0) e definita da

∫ r1

r0

drr(1−r2)

=∫ 2π0 dt = 2π.

15

4.1. STABILITA LINEARE DI ORBITE PERIODICHE 16

4.1 Stabilita lineare di orbite periodiche

Sia x = f(x) un sistema dinamico con un’orbita chiusa. Ci poniamo la domandase questa orbita sia stabile, ovvero se il punto fisso x∗ della mappa di Poincare siastabile. Si applica una perturbazione infinitesima v0 ad x∗ tale che x∗ + v0 sia in S.

Dopo il primo ritorno si ha che:

x∗ + v1 = P (x∗ + v0) = P (x∗) + [DP (x∗)] v0 + O(||v0||2),dove DP (x∗) e la matrice della mappa di Poincare linearizzata intorno a x∗. Allora lastabilita del punto fisso x∗ si esprime in base agli autovalori λj di DP (x∗). L’orbitachiusa e stabile sse |λj | < 1, ∀j = 1, . . . , n − 1. Gli autovalori λj prendono anche ilnome di Moltiplicatori di Floquet.

Si consideri l’esempio precedente r = r(1 − r2); θ = 1. Sappiamo che esiste unciclo limite r∗ = 1. Sia r = 1 + η con η infinitesimo. Allora

r = η = (1 + η)(1− (1 + η)2).

Se trascuriamo i termini non lineari, otteniamo

η = −2η,

cioeη(t) = η0e

−2t.

Dopo 2π, η1 = e−4πη0 e e−4π e il moltiplicatore di Floquet. Poiche |e−4π| < 1 sideduce che il ciclo limite r∗ = 1 e stabile.

Capitolo 5

Flussi sul cerchio e simulazioni

In questo capitolo viene descritta la la dinamica di sistemi oscillanti, quali il pendolosemplice ed il pendolo sovrasmorzato. Si introduce inoltre il concetto di sincroniz-zazione di sistemi dinamici.

5.1 Pendolo semplice

Si consideri il sistema pendolo governato dalle seguenti equazioni dinamiche:

θ = ω

ω = −sin(θ)−Dω

Con θ si rappresenta la posizione angolare del pendolo e con ω si intende lavelocita angolare dello stesso (in realta e la velocita angolare del centro di massadel sistema fisico pendolo). Il parametro D e la costante di proporzionalita relati-va all’attrito viscoso del sistema, infatti la componente Dω dipende dalla velocitaangolare del pendolo. Si considerino due casi distinti da analizzare nel piano del-le fasi. Il primo in cui D = 0 e il secondo in cui D 6= 0 (ma positivo). Il pianodelle fasi e caratterizzato da un campo vettoriale periodico di periodo 2π. Questosignifica che si possono considerare due punti di equilibrio adiacenti per descrivereil comportamento del piano delle fasi.

Caso 1: D = 0. Si consideri solo una finestra del piano delle fasi compresa fra ipunti (θ,ω)=(−π,0) e (2π,0), riportata in Figura 5.1. In questo spazio vi sono quattropunti di equilibrio. Il primo, di coordinate x0=(θ,ω)=(−π,0), e un punto di equilibrioinstabile perche non e ne un attrattore ne Lyapunov stabile (Re(λ) > 0). Il secondopunto, di coordinate x1=(θ,ω)=(0,0), e un punto di equilirio marginalmente stabile(Re(λ) = 0), in quanto per condizioni iniziali nell’intorno di x1, le traiettorie tendonoa formare un centro in x1. Quindi il punto x1 e stabile ma non attrattivo (oscillatorearmonico). Il terzo punto, di coordinate x2=(θ,ω)=(π,0), e un punto di equilibrioinstabile (e un punto di sella) perche per una condizione iniziale nell’intorno della

17

5.1. PENDOLO SEMPLICE 18

posizione di equilibrio x2, la traiettoria associata diverge dal punto di equilibrio.Dato che il sistema non e smorzato (D = 0), una condizione iniziale nell’intorno dellaposizione di equilibrio instabile fa si che la traiettoria ad essa associata, per t →∞,tenda a raggiungere la posizione di equilibrio instabile adiacente, oppure a portareil sistema in oscillazione nell’intorno della posizione di equilibrio stabile pi vicina(dipende dalla velocita iniziale). Fisicamente, considerando un filo inestensibile a cuie collegata una massa, libero di ruotare attorno ad un fulcro, portando il pendoloin posizione di equilibrio instabile e perturbando la sua posizione iniziale (velocitanon nulla), esso tendera ad allontanarsi dal punto di equilibrio instabile e continueraad oscillare per t → ∞. Il quarto punto, di coordinate x3=(2π,0) e un punto diequilibrio instabile con caratteristiche uguali al secondo punto.

θ ’ = ω ω ’ = − sin(θ) − D ω

D = 0

−3 −2 −1 0 1 2 3 4 5 6

−4

−3

−2

−1

0

1

2

3

4

θ

ω

t tendente a −inf

CentroP.to L.−stabile (centro)

P.to instabile (sella)

Figura 5.1: Piano delle fasi del pendolo: D = 0.

Caso 2: D 6= 0. Si sceglie D = 0.25 per la simulazione. In riferimento allaFigura 5.1, il piano delle fasi risultante e diverso dal precedente per cio che riguardai punti di equilibrio stabili. Infatti considerando una condizione iniziale nell’intornodi un punto di equilibrio instabile, la traiettoria associata converge verso il puntodi equilibrio stabile adiacente formando una spirale (λ ∈ C). In altri termini, lavelocita del pendolo diminuisce gradualmente per t → ∞ e la posizione angolare si

5.1. PENDOLO SEMPLICE 19

assesta nell’intorno della posizione di equilibrio stabile, al contrario di quanto accadeper D = 0, caso in cui il pendolo oscilla sempre con la stessa ampiezza (oscillatorearmonico). I punti di equilibrio stabili sono anche attrattori perche le traiettorie concondizione iniziale nell’intorno dei punti di equilibrio stabili convergono ad essi pert → ∞. Aumentando il parametro D si nota che diminuisce il tempo con cui unatraiettoria giunge nel punto di equilibrio stabile, partendo da un intorno di esso, pert →∞ e con movimento a spirale.

θ ’ = ω ω ’ = − sin(θ) − D ω

D = 0.25

−3 −2 −1 0 1 2 3 4 5 6

−4

−3

−2

−1

0

1

2

3

4

θ

ω

t tendente a −infSpirale

P.to L.−stabile e attrattore (spirale)

P.to instabile (sella)

Figura 5.2: Piano delle fasi del pendolo: D = 0.25.

5.2. PENDOLO SOVRASMORZATO 20

5.2 Pendolo sovrasmorzato

In questa sezione si analizza un sistema oscillante al variare dei parametri caratte-ristici, mettendo in luce la perdita di uniformita dinamica del sistema stesso, cioe ilmoto del sistema non e piu simmetrico.

Descrizione del modello

Si consideri il modello del pendolo fisico descritto dalla seguente equazione differen-ziale del secondo ordine:

mL2θ + bθ + mgLsin(θ) = Γ (5.1)

in cui ogni parametro e positivo. In particolare m e la massa del pendolo, L ela lunghezza dello stesso, b e il coefficiente di attrito viscoso del sistema e Γ e unacoppia costante positiva. Per inserire il modello nell’interfaccia di pplane7 si puoriscrivere l’equazione differenziale del secondo ordine come sistema di equazioni delprimo ordine come segue:

θ = ω ω = −gsin(θ)L

− b

mL2ω +

F

mL2

A questo punto si considera il caso in cui il pendolo sia sovrasmorzato imponendob grande e il modello (5.1) diventa il seguente:

b

mgLθ =

ΓmgL

sin(θ)

doveτ =

mgL

bt; γ =

ΓmgL

Quindi si ottiene:θ′ = γ − sin(θ)

doveθ′ =

∂θ

∂τ

Simulazione

Gli angoli sono misurati in senso antiorario. La posizione di equilibrio instabile valeπ e quella stabile vale 0. Si analizzano i piani delle fasi risultanti al variare di γ cherappresenta il rapporto fra la coppia esterna applicata e la massima coppia gravita-zionale.

5.3. BEAT PHENOMENON 21

Caso 1: γ → 1+. Il piano delle fasi e mostrato in Figura 5.3.a. Sono visibilianche le nullclines, le quali non si intersecano mai, dati i parametri associati almodello (5.1), ma sono molto vicine ad ogni periodo della sinusoide in arancione(θ = 0). Non esistono punti di equilibrio (F non viene mai annichilita), pertantoper ogni condizione iniziale il pendolo ruota per t → ∞. Il moto del pendolo e piurapido in discesa, partendo da una condizione iniziale pari a π per raggiungere 2π,piuttosto che in salita, partendo da 0 verso π (ω >> 0). Infatti in discesa agisconosul pendolo sia la coppia F che la coppia gravitazionale, mentre in salita l’energiafornita dalla coppia F viene in parte persa con l’azione opponente della coppiagravitazionale e dell’attrito viscoso (supposto elevato). La Figura 5.3.b mostra unasoluzione particolare del piano delle fasi in funzione del tempo. In particolare conb = 5 e γ → 1+.

Caso 2: γ > 1. Piu e elevato il rapporto γ e piu le nullclines si distanziano.Non esistono punti di equilibrio. In questo caso il tempo di salita e il tempo didiscesa sono pressoche uguali, in quanto la coppia F e cosı elevata da rendere quasitrascurabili i termini relativi alla forza peso e all’attrito viscoso. Caso 3: γ = 1.Appare un punto di equilibrio in posizione θ = π/2. In Figura 5.4 e rappresentato ilpiano delle fasi con γ = 1. Le nullclines si intersecano ad ogni periodo della sinusoideed e evidenziato in rosso il punto di equilibrio (nodal skin) periodico.Caso 4: γ < 1. In questo caso i punti di equilibrio sono due e sono tanto piu distantiquanto γ → 0. Il primo e uno spiral skin e il secondo e una sella. Il piano delle fasie mostrato in Figura 5.5.Conclusioni. Si nota come al variare del rapporto γ si possa osservare un sistemaoscillante non uniforme. Partendo da γ = 0 (caso non contemplato in quanto ilsistema si riduce ad un oscillatore con due punti di equilibrio di cui il minore e stabile,vedi il Capitolo 1), per γ → 1.5 si osserva come i due punti di equilibrio si concentrinoin uno solo (per γ = 1) per poi svanire completamente. L’analisi qui presentata puoessere fatta anche considerando un campo vettoriale sul cerchio, considerando unsistema oscillante non uniforme semplice (di primo grado) e considerando il flussosul cerchio, piuttosto che su di una retta.

5.3 Beat phenomenon

In questo esempio si mostra un primo esempio di sincronizzazione fra sistemi.


Nei sistemi oscillanti con due variabili indipendenti e interessante capire quando ledue variabili assumono lo stesso valore. Si consideri l’esempio riportato nel libro ditesto, in cui vi sono due corridori che partono insieme su un circuito circolare. Ilprimo impiega T1 secondi per compiere un giro, mentre il secondo impiega T2 > T1

secondi (e piu lento). Bisogna capire quando il primo corridore sorpassa il secondo.Dato che la velocita dei due corridori e costante, i sistemi che li rappresentano sono

5.4. LUCCIOLE 22

descritti dalle equazioni:θi = ωi

dove i = 1, 2 sono i due corridori e ωi = 2π/Ti e la pulsazione naturale dei duesistemi. Definendo la differenza di fase φ = θ1 − θ2, si vuole capire quando essa vale2π, cioe quando il primo corridore, al tempo t, ha compiuto un giro in piu rispettoal secondo corridore. Quindi

φ = θ1 − θ1 = ω1 − ω2 = 2π

In questo caso θ1 e in fase con θ2 ogni

T =2π

ω1 − ω2=

T1T2

T1 − T2

Simulazione

Per la simulazione di questo sistema si utilizza un codice Matlab creato apposita-mente. Il sistema viene riscritto utilizzando il sistema di riferimento in coordinatepolari: (

xy

)=

(cos(ωt)sin(ωt)

)

Si associa ad ogni corridore un sistema in coordinate polari. Per la simulazione sisceglie la seguente configurazione:

• ω1 = 2;

• ω2 = 1;

Cio significa che il primo corridore doppiera il secondo ogni

Tlap =2π

ω1 − ω2= 2π

La condizione iniziale e (x1,2, y1,2) = (0, 0). Si considera uno spazio temporale diquattro periodi e si tracciano le curve φx = x1 − x2 e φy = y1 − y2 dei due sistemi.Come mostrato in Figura 5.6, il primo corridore doppiera il secondo ogni volta che(φx, φy) = (0, 0), cioe ogni volta che i due grafici si intersecano contemporaneamentesulla retta y = 0.

5.4 Lucciole

In questa sezione si analizza un secondo esempio di sincronizzazione tra sistemi.

5.4. LUCCIOLE 23


In correlazione con l’esempio precedente, vi e il fenomeno naturale della luce pe-riodica emessa dalle lucciole. La differenza sostanziale e che in questo caso, oltre alfatto che la luce emessa dalla lucciola e fortemente collegata alla presenza di unostimolo esterno, la differenza di fase tra l’emissione di luce e lo stimolo non aumentalinearmente. Di seguito si fa una traduzione di alcune parti del libro di testo per unamigliore comprensione .Modello. Supponendo che θ(t) sia la fase con cui la lucciola si illumina, dove θ = 0corrisponde al momento in cui viene emessa la luce, e che in assenza di stimoli esternila lucciola emette luce a frequenza ω, si puo scrivere:

Θ = Ω

dove Θ e un generico stimolo periodico. La risposta della lucciola puo essere mo-dellata come segue: se lo stimolo e avanti nel ciclo, allora si assume che la lucciolaaumenti la frequenza di emissione per sincronizzarsi; se la lucciola emette luce a unafrequenza maggiore dello stimolo, allora tentera di diminuire la sua frequenza. Unmodello semplice e il seguente:

θ = ω + Asin(Θ− θ)

in cui A > 0. Per esempio, se Θ e avanti rispetto a θ (0 < Θ − θ < π) allora lalucciola aumenta la frequenza di emissione di luce (θ > ω). Il parametro A e unamisura dell’abilita della lucciola di modificare la sua frequenza istantanea (resettingstrength).Analisi. Si consideri la differenza di fase φ:

φ = Θ− θ = Ω− ω −Asin(φ)

In questa forma e un oscillatore non uniforme. Come per il pendolo, si introduconodue parametri che rendono il sistema adimensionale:

τ = At; µ =Ω− ω

A

Alloraφ′ = µ− sin(φ)

dove φ′ = ∂φ/∂τ .

Simulazione

A questo punto si simula il modello con pplane7 al variare del parametro µ.Caso 1: µ = 0. Come mostra la Figura 5.7, le traiettorie tendono a zero per t →∞(in verde). Infatti φ = 0 e un punto di equilibrio in cui si concentrano tutte letraiettorie con condizione iniziale nell’intorno di zero. Invece φ = π/2 e un punto

5.4. LUCCIOLE 24

di equilibrio instabile, in quanto le traiettorie tendono a divergere per t → ∞ (inrosso).Caso 2: 0 < µ < 1. I punti di equilibrio si avvicinano, ma mantengono le stessecaratteristiche di stabilita, come mostra la Figura 5.8. Le nullclines del sistemalucciola-stimolo si allontanano fra loro. In realta il fatto che il punto di equilibriostabile φ∗ > 0 significa che la lucciola e lo stimolo esterno hanno la stessa frequenzama sono sfasati (phase-locked).Caso 3: µ = 1. In Figura 5.9 si nota come le traiettorie convergano all’unico punto diequilibrio se hanno condizione iniziale minore di esso e divergano se hanno condizioneiniziale maggiore dello stesso. Il punto di equilibrio e mezzo-stabile.Caso 4: µ > 1. In questo caso le traiettorie divergono per t → ∞ per qualunquecondizione iniziale (Figura5.10. Non esistono punti di equilibrio e le nullclines delsistema lucciola-stimolo non si intersecano mai. L’emissione di luce della lucciola nonsara mai in fase con l’emissione di luce esterna. Questo e dovuto al fatto che lo stimoloesterno ha una frequenza elevata e il sistema lucciola non riesce ad agganciarsi infase ad esso.

5.4. LUCCIOLE 25

θ ’ = ω ω ’ = − g/L sin(θ) − b/(m L2) ω + F/(m L2)

g = 9.81L = 1

m = 1F = 10

b = 5

−5 −4 −3 −2 −1 0 1 2 3 4 5−4

−3

−2

−1

0

1

2

3

4

θ

ω

(a) Piano delle fasi

θ

ω


g = 9.81L = 1

m = 1F = 10

b = 5

0 2 4 6 8 10 12 14 16 18 20

0

2

4

6

8

10

12

14

16

18

t

θ an

d ω

X: 15.07Y: 3.177

X: 16.15Y: 6.337

Salita (lenta)

Discesa (rapida)

(b) Soluzione particolare in funzione del tempo

Figura 5.3: Piano delle fasi con γ → 1+

5.4. LUCCIOLE 26


g = 9.81L = 1

m = 1F = 9.81

b = 5

−10 −8 −6 −4 −2 0 2 4 6 8 10

−4

−3

−2

−1

0

1

2

3

4

θ

ω

Figura 5.4: Piano delle fasi con γ = 1

5.4. LUCCIOLE 27


g = 9.81L = 1

m = 1F = 5

b = 5

−5 −4 −3 −2 −1 0 1 2 3 4 5

−4

−3

−2

−1

0

1

2

3

4

θ

ω

Figura 5.5: Piano delle fasi con γ < 1

5.4. LUCCIOLE 28

0 5 10 15

−1.5

−1

−0.5

0

0.5

1

1.5

2

X: 6.293Y: 0.009814

t

φ x , φ y

(φx,φ

y)=(0,0) (φ

x,φ

y)=(0,0)

C.I.=(0,0)

Figura 5.6: Beat phenomenon

−2 −1 0 1 2 3 4−2

−1

0

1

2

3

4

t

phi

φ ’ = O − o − A sin(φ) O = 1 A = 1 o = 1

Instabile

Stabile

Figura 5.7: Lucciole: µ = 0

5.4. LUCCIOLE 29

−2 0 2 4 6 8 10

−4

−3

−2

−1

0

1

2

3

4

t

phi

φ ’ = O − o − A sin(φ) O = 1 A = 1 o = 0.5

Instabile

Stabile

Figura 5.8: Lucciole: 0 < µ < 1

−2 0 2 4 6 8 10

−4

−3

−2

−1

0

1

2

3

4

t

phi

φ ’ = O − o − A sin(φ) O = 2 A = 1 o = 1

Instabile

Stabile

Figura 5.9: Lucciole: µ = 0

5.4. LUCCIOLE 30

−2 0 2 4 6 8 10

−4

−3

−2

−1

0

1

2

3

4

t

phi

φ ’ = O − o − A sin(φ) O = 3 A = 1 o = 1

Figura 5.10: Lucciole: µ > 1

Capitolo 6

Il modello di Lorenz

In un articolo del 1963, Deterministic Nonperiodic Flow, Lorenz cita le equazioni diSaltzman (1962), derivate dall’analisi di un sistema composto da un fluido di profon-dita (altezza) uniforme, con differenti temperature tra la parte superiore e la parteinferiore, in particolare con una variazione di temperatura lineare. SuccessivamenteRayleigh scoprı un punto critico per il quale le equazioni di Saltzman mostrano unmoto convettivo, basato sul numero di Rayleigh. In questo capitolo si analizza ilsistema di Lorenz a partire dal modello matematico.

Modello matematico

Lorenz definisce un set di tre equazioni dinamiche. In particolare x e proporzionaleall’intensita del moto convettivo, y e proporzionale alla differenza di temperatura trale correnti ascensionali e discensionali e z e proporzionale alla distorsione del profiloverticale di temperatura (dalla linearita). Il sistema tridimensionale e il seguente:

x = −σx + σy

y = −xz + rx− y

z = xy − bz

dove, in questo caso, σ = 10 (numero di Prandtl), r e il numero di Rayleigh eb = 8/3 e la largezza approssimativa della regione di lavoro. Di seguito si analizzanoalcune peculiarita del sistema di Lorenz, con particolare attenzione alle biforcazionidel sistema al variare del numero di Rayleigh 1.

1Il modello della convezione del fluido e realistico per valori di r compresi fra 0 e 1. Per valorimaggiori di r si modellano sistemi come spiking irregolare nei lasers, convezione in regioni toroidali,etc.

31

32

Analisi con matcont per 0 ≤ r < 30

In questa sezione si analizza empiricamente il sistema di Lorenz al variare del para-metro r, con σ e b fissati.

Caso 1: 0 ≤ r < 1. Il punto di equilibrio stabile, di coordinate (x, y, z) =(0, 0, 0), e anche un attrattore globale. In Figura 6.1 si riportano tre simulazioni perr = 0; 0.5; 1 − ε. Per comodita viene rappresentata la variabile x(t) in funzione deltempo.

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5−10

−8

−6

−4

−2

0

2

4

6

8

10

t

x

0 < r <1

Figura 6.1: Caso 1: 0 ≤ r < 1

Caso 2: 1 < r < 24. Per r = 1 si verifica la prima biforcazione del sistema.In Figura 6.2 si puo osservare come le traiettorie, per t → ∞, generate a partireda condizioni iniziali opposte in segno, si allontanino definitivamente dalla prece-dente soluzione di equilibrio stabile ((0, 0, 0)) e si distanzino fra loro al variare dir. L’incremento del parametro r comporta la perdita di stabilita della soluzione diequilibrio (0, 0, 0) e la generazione di due nuove soluzioni di equilibrio stabili e sim-metriche nel diagramma di biforcazione. Risulta chiaro che la biforcazione e di tipopitchfork subcritica. Nel caso riportato in Figura 6.3, in cui r = 22, le oscillazioni delsistema sono molto sostenute, ma comunque smorzate (spirale). Inoltre, osservandoil transitorio di x(t), si nota la presenza di uno spike. Di seguito si comprende comemai vi sia questo fenomeno.

Caso 3: 24 < r < 30. Per r = 24.7 si verifica la seconda biforcazione del sistema.

33

0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5−10

−8

−6

−4

−2

0

2

4

6

8

10

t

x

r è circa 1 e le traiettorie sono stabili verso (0,0,0)

r incrementa

traiettoria sovrasmorzata

Figura 6.2: Caso 2: 1 < r < 24

In realta, per la simmetria della biforcazione pitchfork, le biforcazioni sono due e sonosimmetriche rispetto all’asse x = 0. Si nota che i due autovalori complessi e coniugatiλ2, λ3 hanno parte reale strettamente minore di zero per r = 24.7 − δ e assumonoparte reale strettamente maggiore di zero per r = 24.7+δ. In questo caso l’equilibrioe instabile e la biforcazione e di Hopf subcritica. Le traiettorie tridimensionali, comeenuncia la teoria classica, dovrebbero allontanarsi definitivamente dagli autovettoriinstabili. In realta nessuna traiettoria diverge, ma si nota che esse entrano in dueregioni diverse senza un criterio intuitivo, con moto a spirale, senza mai stabilizzarsinell’origine della spirale stessa. In riferimento alla Figura 6.5, si osserva il motoa spirale caotica in tridimensione. Il segnale e assolutamente aperiodico, rendendoimpredicibile il comportamento del sistema per valori futuri del tempo, rispetto adun tempo t0 generico. Da un sistema deterministico e totalmente predicibile si passaad un sistema caotico impredicibile e, come mostra la Figura 6.4, molto sensibile allecondizioni iniziali, per un r fissato. Infatti per una piccolissima variazione (dell’ordinedi 10−10) della condizione iniziale relativa alla coordinata z si ottengono due rispostetemporali differenti.In Figura 6.6 e rappresentato il diagramma di biforcazione del modello di Lorenzche riassume i tre casi analizzati.

34

Altri casi notevoli

Si consideri il caso in cui 144 < r < 166. In Figura 6.7 e rappresentata la risposta delsistema per r = 145. Si nota che le traiettorie tendono a stabilizzarsi su un’orbitaperiodica (per r crescente), ma tornano ad essere aperiodiche dopo un determinatointervallo di tempo. Questo fenomeno e chiamato caos intermittente e accade peruna finestra di valori di r, tale per cui si verifica il period doubling. Se r e dentroquesto range, il caos intermittente puo eventualmente cessare, ma se r assume unvalore prossimo all’ upper bound (r = 166.3), allora il fenomeno si verifica per ognivalore di t. Infatti per r crescente da 145 a 166 si nota che il periodo in cui domina ilcaos intermittente tende ad aumentare fino a dominare l’intera evoluzione temporale.

Nel caso in cui r = 212 si osserva un secondo fenomeno di period doubling, macon periodo raddoppiato rispetto alla finestra considerata in precedenza. La rispostadinamica del sistema segue una traiettoria stabile, ma in questo caso il caos non edi tipo intermittente, in quanto, per r crescente, la periodicita della soluzione x(t) edi esattamente un periodo per ogni intervallo temporale.In Figura 6.8 e rappresentata la realizzazione del sistema per r = 212.

Conclusioni

Dall’analisi proposta emergono alcune delle caratteristiche principali del sistema diLorenz.Aperiodicita. Per un determinato valore del numero di Rayleigh si osserva il fe-nomeno dell’aperiodicita. E’ del tutto contro-intuitivo pensare che da un sistemalineare si riesca a generare un segnale aperiodico, non predicibile e apparentementecasuale (in un certo senso e anche non causale).Sensitivita. Il sistema risulta sensibile alle condizioni iniziali per quanto riguardala risposta aperiodica. Un fattore di differenza molto piccolo genera due segnali to-talmente differenti dopo un tempo t. Dall’abstract dell’articolo di Lorenz del 1963:[... soluzioni nonperiodiche sono instabili rispetto a piccole perturbazioni, cosicchecondizioni iniziali di poco differenti possono evolvere in stati considerevolmente di-versi...].Biforcazioni. Al variare di r vi sono tre biforcazioni di cui due sono simmetriche,ma le traiettorie, per qualsiasi valore di r, non divergono verso l’infinito. Un sistemacon questo comportamento si chiama dissipativo.Simmetria. Il sistema e invariante se si inverte il segno delle variabili di stato (tran-ne z). L’invarianza dell’asse z implica che tutte le traiettorie relative ad esso restanosullo stesso e si avvicinano all’origine. Un fenomeno di simmetria si osserva ancheper le biforcazioni di Hopf, dal momento che si generano a partire dai due ramistabili della biforcazione pitchfork.

35

0 5 10 15 20 25 30−20

−15

−10

−5

0

5

10

15

20

t

x

(a) x(t)

0

5

10

15

20−20 −15 −10 −5 0 5 10 15 20

0

5

10

15

20

25

30

35

40

yx

z

(b) Plot tridimensionale

Figura 6.3: Caso 2: r=22

36

0 5 10 15 20 25 30 35 40t

si "gira" nella seconda spirale

si "gira" nella prima spirale

scambio

differenziazione traiettorie

Figura 6.4: Risposte temporali per r = 27 e condizioni iniziali differenti per 10−10

0

5

10

15

20

−20−15−10−505101520

x

y

Figura 6.5: Moto caotico a spirale

37

−10 −8 −6 −4 −2 0 2 4 6 8 10

0

5

10

15

20

25

30

x

r

BP

H

BP

H

Figura 6.6: Biforcazioni al variare del numero di Rayleigh

38

0 2 4 6 8 10 12−40

−30

−20

−10

0

10

20

30

40

t

x

(a) x(t)

−1000

100200

300

−100−50050100150

0

50

100

150

200

250

yx

z

(b) Plot tridimensionale

Figura 6.7: Risposta del sistema per r = 145

39

−100

0

100

200

300

400

−100−50

050

100150

0

50

100

150

200

250

x

y

z

(a) Plot tridimensionale

0 2 4 6 8 10 12 14 16 18 20−50

−40

−30

−20

−10

0

10

20

30

40

t

x

(b) x(t)

Figura 6.8: Risposta dinamica per r = 212

Capitolo 7

Controllo del caos

Per completare la panoramica sul caos deterministico, si analizza il problema delcontrollo del caos. Si nota che al giorno d’oggi si utilizzano tecniche di controllo chederivano senza alterazioni dalla teoria classica del controllo.

Introduzione

Un sistema caotico e caratterizzato dalla dipendenza dalle condizioni iniziali. Duetraiettorie che partono da punti molto vicini nello spazio divergono, ad esempio,esponenzialmente. Infatti, in un processo di identificazione di un sistema (black-box ), spesso non si riesce a distinguere esattamente se si identifica un sistema affettoda rumore oppure un sistema caotico. Un quantificatore della divergenza delle duetraiettorie e l’esponente di Lyapunov, in particolare il suo massimo. Un secondoproblema e la possibilita di sincronizzazione fra due sistemi caotici identici. Recentistudi dimostrano che effettivamente puo avvenire la sincronizzazione, con importantiapplicazioni come telecomunicazioni digitali non lineari (Pecora e Carrol, 1990).

Il problema messo in luce in questo lavoro e la possibilita di stabilizzare unsistema caotico (di Lorenz) con un semplice regolatore PI (proporzionale-integrale).

Teoria del controllo

Con la parola controllo si intende controllo in retroazione. Si consideri un sistemanon lineare. Si hanno alcune variabili dipendenti dal tempo e si puo pensare dicontrollarle, ossia di modificarne lo stato e la sua evoluzione, conoscendo la storiadel sistema e il suo stato attuale. Il controllo puo essere inteso in tre modi diversi,elencati di seguito:

1. Dato un sistema caotico in costante evoluzione, aumentare la sua performancemedia nel tempo, dove performance e un termine generico per definire unacaratteristica particolare che si vuole ottenere dal sistema;

40

7.1. CONTROLLO DI UN PROCESSO CAOTICO IN CONTINUAEVOLUZIONE 41

2. Dato un sistema caotico con una condizione iniziale associata, trascinare rapi-damente la traiettoria rapidamente verso lo stato desiderato;

3. Dato un processo caotico in costante evoluzione, prevenire che diverga versouna regione non buona dello spazio delle soluzioni.

L’esistenza di un’orbita complessa e molto rilevante per ottenere il target asse-gnato nel primo dei tre problemi. Per il secondo e il terzo target e piu importantela sensibilita esponenziale alle condizioni iniziali. Per tutti i targets e da notare cheil controllo puo essere ottenuto soltanto con piccoli controlli, cioe con piccole va-riazioni dello stato attuale. L’interesse per il terzo target risiede nella possibilita dicontrollare e prevenire eventi catastrofici che si verificano nel momento in cui lo statoconverge in una precisa zona dello spazio delle soluzioni.

Ovviamente esiste un approccio diverso per ognuno dei tre targets. Nelle sezioniseguenti si tratta del primo e del secondo tipo di approccio rispettivamente per ilprimo e per il secondo target.

7.1 Controllo di un processo caotico in continua evolu-zione

Il punto chiave per questo problema e la presenza di un’ orbita complessa nel caos.Con questo si intende che esistono molte orbite topologicamente distinte all’internodi un set caotico invariante. Quindi vi sono tipicamente molte orbite periodicheinstabili (UPO - Unstable periodic orbit) nascoste in un attrattore caotico. Infatti,per un attrattore iperbolico, il numero di UPOs con periodo minore o uguale di Pincrementa esponenzialmente con P ,

#(P ) ∼ exp(hT P ) (7.1)

dove hT e l’entropia topologica (Katok, 1980). Per risolvere il problema esiste l’ap-proccio seguente (Ott et al. 1990a,b):

1. Esaminare le dinamiche non controllare del sistema e determinare alcune UPOsdi periodo piccolo nascoste nell’attrattore;

2. Esaminare queste UPOs e trovare la caratteristica risultante dalla dinamicatemporale di ogni UPO. Confrontare le caratteristiche delle UPOs e sceglierela migliore.

3. Formulare un algoritmo di controllo che stabilizzi la UPO selezionata in unintorno di essa.

In sostanza si cerca di far convergere una traiettoria verso un’orbita propria delsistema e mantenere questo comportamento per tutta la durata del processo. Nelcaso in cui una traiettoria tenti di divergere da quest’orbita, si cerca di intervenire

7.2. TARGETING 42

tempestivamente per riportarla sull’orbita. La tempestivita e importante perche ilsistema puo essere controllato in maniera semplice solo in un intorno molto ristret-to della UPO selezionata. Si pensi all’attrattore di Lorenz. Sostanzialmente esisteun’orbita in cui le traiettorie convergono e successivamente divergono e si spostanoverso una seconda orbita, tentando di convergere, ma, per la naturale evoluzione delsistema non controllato di Lorenz, divergono anche da essa e si sposteranno versola prima orbita, reiniziando il ciclo. Se si identificasse una regione Ω con l’approcciodi Ott, si potrebbe pensare di stabilizzare la traiettoria nell’intorno di una delledue orbite caotiche. Data la natura del sistema non si puo pensare di stabilizzare latraiettoria nell’origine, in quanto non si saprebbe come governare il sistema versoun punto che non e piu di equilibrio (nel caso caotico). A questo punto diventerebbeun problema di inseguimento e non piu di stabilizzazione, ma in realta un sistemacaotico non ha ingressi esterni (e autonomo di ed ha caos autonomo).

Per identificare una UPO e le sue caratteristiche si considerano delle mappericavate sezionando con un piano un’orbita dell’attrattore e analizzando le serietemporali derivate. In questo modo si trasferisce il problema tridimensionale in unproblema bidimensionale. Ovviamente la riduzione del problema ad un caso piusemplice porta ad avere un’analisi semplificata. Nella fase in cui si analizzano leserie temporali vi e in gioco il fattore rumore. Come precedentemente illustrato eun fattore di non facile identificazione, soprattutto quando non si conosce il sistemada cui proviene l’attrattore strano. Per stabilizzare una UPO si utilizza la teoria delcontrollo lineare. Si considerano quindi approcci come la tecnica del pole placementoppure un regolatore generico PI, come descritto nell’esempio di questa relazione.

7.2 Targeting

Per la seconda tipologia di controllo del caos si analizza il problema di direzionarel’orbita di un sistema caotico in una zona desiderata del piano delle fasi. Alla basedel targeting vi e l’idea che un sistema caotico e esponenzialmente sensibile a piccoleperturbazioni. Una piccola perturbazione ha un largo effetto in un tempo relati-vamente piccolo. Se le perturbazioni sono rumori ignoti, allora il sistema diventaimpredicibile. Ma se il rumore presente e relativamente piccolo, si possono tran-quillamente inserire piccole perturbazioni per controllare il sistema. Una possibilestrategia di controllo e il metodo forward-backward (Shinbrot et al., 1990), descrittodi seguito. Si consideri una mappa bidimensionale con una direzione di contrazione(un esponente di Lyapunov negetivo) e una direzione di espansione (un esponentedi Lyapunov positivo). Considerando uno scalare di controllo p, si puo perturbareil punto di partenza verso un set di regioni nello spazio delle soluzioni, cioe si cercadi indirizzare il punto di partenza verso la regione desiderata. L’azione di control-lo e manifestata nel sistema come perturbazione iniziale, cioe il controllo e attivosolo inizialmente (per t = 0). A questo punto si utilizza la conoscenza a priori delsistema, quindi la conoscenza del set di traiettorie che il punto potrebbe percorrere

7.3. CONTROLLO DEL SISTEMA DI LORENZ 43

evolvendo in un tempo T. Simultaneamente si fa evolvere la regione a cui si puntaindietro nel tempo T . Incrementando T , l’immagine del segmento di controllo inizia-le (in avanti nel tempo) e l’immagine della regione puntata (indietro nel tempo) siespandono esponenzialmente in lunghezza. Nel momento in cui la preimmagine dellaregione puntata e l’immagine del segmento di intersecano, si identifica un valore delsegmento iniziale p. Idealmente la traiettoria del punto iniziale tendera a colliderecon la regione puntata in un tempo 2T . In un sistema reale, affetto da rumore, sicommette un errore anche molto piccolo quando si imposta il segmento iniziale: l’er-rore aumenta esponenzialmente nel tempo ed e e probabile che il controllo adottatonon sia robusto. Questo sistema funzione per una dimensione del sistema non troppoelevata. Per sistemi di dimensione alta esiste il metodo di Kostelich (Kostelich etal., 1993). La prima applicazione del controllo formalizzato da Kostelich si ha nel1985 da parte della NASA per il targeting di una sonda spaziale con la cometa diGiacobini-Zinner.

7.3 Controllo del sistema di Lorenz

In questa sezione si analizza la possibilita di controllare alcuni punti instabili delsistema di Lorenz e stabilizzarli con la teoria del controllo lineare, in particolarecon un regolatore proporzionale e integrale, utilizzando l’approccio di Ott (primometodo). Si consideri il sistema seguente:

x = Ax + g(x) + u (7.2)

dove x ∈ Rn e il vettore di stato, u ∈ Rn e il comando da progettare, A ∈ Rn×n euna matrice costante e g(x) e una funzione non lineare continuamente differenziabile.Considerando xs un punto di equilibrio instabile del sistema, tale che

Axs + g(xs) = 0 (7.3)

si puo progettare un controllore u che stabilizzi la traiettoria verso xs. Il controlloreu e un PI della forma

u = λ[B(Kx + k

∫ t

0(y − ys)dτ)] (7.4)

dove B ∈ R1×n e un vettore costante di guadagni, K ∈ R1 timesn e il vettoredi guadagni di retroazione proporzionale allo stato, k ∈ R e il guadagno integrale,y = Cx e l’uscita con una matrice costante C ∈ R1×n, ys = Cxs e l’osservazionedel punto di equilibrio target xs, e λ ∈ 0, 1 e una variabile binaria di attivazioneche vale 1 se il vettore di stato ad un certo tempo t e dentro la regione di controlloΩxs di xs. Quindi il controllore u e attivo se e solo se la traiettoria del punto x(t) edentro la regione ΩxS . Si riscrive il sistema come segue:

x = Ax + g(x) + λ

[B

(Kx + k

∫ t

0(Cx− Cxs)dτ

)](7.5)


Si consideri il sistema di Lorenz. Imponendo le condizioni per i punti di equilibrio(derivate prime uguali a zero), si identificano tre punti di equilibrio dipendenti daiparametri: (−

√β(ρ− 1),−

√β(ρ− 1), (ρ−1)), (0, 0, 0) e

√β(ρ− 1),

√β(ρ− 1), (ρ−

1)). Riscrivendo il sistema nella forma 7.5, segue che l’errore di inseguimento edescritto dalla seguente formula:

g(x)− g(xs) =

0 0 0−zs 0 −xs

ys xs 0

ex

ey

ez

+

0−exez

exey

= g′(xs)e + g(e) (7.6)

La dinamica dell’ errore incrementale e descritto dalle seguenti equazioni di stato:

e1 = (A + λBK)e1 + g(e1) (7.7)

dove

e1 =[

eq

](7.8)

g(e1) =[

g(e)0

]=

0−exez

exey

0

(7.9)

A = A + g′(xs) =

−σ σ 0ρ− zs −1 −xs

ys xs −β

(7.10)

A =[

A 0C 0

]=

−σ σ 0 0ρ− zs −1 −xs 0

ys xs −β 0c1 c2 c3 0

(7.11)

B =[

B0

]=

b1

b2

b3

0

(7.12)

K =[

k1 k2 k3 k]

(7.13)

Se la matrice A+BK e stabile in xs , allora si puo affermare che il controllore stabiliz-za l’evoluzione del sistema verso il target (la 7.7 e asintoticamente stabile). Se la stes-sa matrice e stabile nel punto di equilibrio P1 =(−

√β(ρ− 1),−

√β(ρ− 1), (ρ−1))),

allora la 7.7 e asintoticamente stabile e il sistema di Lorenz puo essere stabilizzatoverso il punto di equilibrio P1. Le stesse osservazioni valgono per P2 = (0, 0, 0) eP3 = (

√β(ρ− 1),

√β(ρ− 1), (ρ − 1)). Ad esempio si puo pensare di stabilizzare il

sistema di Lorenz per σ = 16, ρ = 45.6 e β = 4, per i quali il sistema ha un com-portamento caotico. Il sistema e stabilizzabile con il controllore PI verso il punto di


equilibrio P3 perche la matrice A + BK e stabile in xs ha autovalori a parte realeminore di zero, quindi e stabile (λ1 = −2.8, λ2 = −2.7, λ3 = −2.3, λ4 = −2.5). Inquesto caso le matrici del sistema di controllo sono le seguenti

B =

100

(7.14)

C =[ −1 1 0

](7.15)

K =[

10.70 −19.18 12.91], k = 0.12 (7.16)

In Figura 7.1 e riportato lo schema Simulink realizzato ad hoc per la simula-zione. Le traiettorie delle tre variabili di stato tendono asintoticamente al punto diequilibrio P3 = (13.35, 13.35, 44.60) pertanto l’obiettivo proposto e stato raggiunto.Nelle Figure 7.2, 7.3 e 7.4 si mostrano le traiettorie temporali rispettivamente perx(t), y(t) e z(t).

Figura 7.1: Modello Simulink del sistema di Lorenz controllato con un PI


50 100 150 200 250 300 350 400 450 500t

Figura 7.2: x(t)

0 50 100 150 200 250 300 350 400 450 500−10

−5

0

5

10

15

20

y

x

Figura 7.3: y(t)


0 50 100 150 200 250 300 350 400 450 50030

35

40

45

50

55

60

t

z

Figura 7.4: z(t)

Bibliografia

[1] Edward Ott Chaos in dynamical systems Cambridge University Press, 1993.

[2] G. Jiang, G. Chen, W.K. Tang Stabilizing unstable equilibrium points of aclass of chaotic systems using a state PI regulator (IEEE Trans. on Circuits andSystems, VOL. 49, NO. 12, 2002)

[3] Steven H. Strogatz Nonlinear Dynamics and Chaos Westview, 1994.

48

Facoltµa di Ingegneriamocenni/dispense_07_08.pdfL’importante teorema di Hartman-Grobman...

Documents

Transcript of Facoltµa di Ingegneriamocenni/dispense_07_08.pdfL’importante teorema di Hartman-Grobman...