Zagadnienia

4. Równania z małym parametrem

Mały parametr w równaniu różniczkowym może pojawiać się zasadniczo na dwa sposoby: albo z prawej strony albo z lewej strony. W pierwszym przypadku mamy do czynienia z małym zaburzeniem układu, o którym na ogół sporo wiemy a w drugim przypadku pojawiają się tzw. drgania relaksacyjne. Oba przypadki są omawiane kolejno w następnych rozdziałach.

4.1. Uśrednianie

Przykładem układu z pierwszej grupy jest znany układ van der Pola

x´=y,  ⁢y˙=-x+ε⁢1-x2⁢y,

gdzie ε jest naszym małym parametrem. Jest to szczególny przypadek zaburzenia układu hamiltonowskiego o jednym stopniu swobody

x˙=Hy′+ε⁢P⁢x,y,   ⁢y˙=-Hx′+ε⁢Q⁢x,y. (4.1)

W zastosowaniach często pojawiają się układy hamiltonowskie z wieloma stopniami swobody postaci

q˙i=Hpi′,   ⁢p˙i=-Hqi′,    ⁢i=1,…⁢n, (4.2)

gdzie H⁢q1,…,qn,p1,…,pn jest funkcją Hamiltona, lub hamiltonianem (Zadania 4.11 i 4.12). Na ogół układ (4.2) nie daje się rozwiązać. Jednak istnieje klasa układów hamiltonowkich w pełni rozwiązywalnych.

Definicja 4.1. Układ (4.2) nazywa się zupełnie całkowalnym jeśli istnieje układ funkcjonalnie niezależnych całek pierwszych  F1=H, F2,…, Fn taki, że każda funkcja Fj jest całką pierwszą dla innych układów hamiltonowskich generowanych przez inne funkcje Fi. Mówi się też, że funkcje Fj są w inwolucji.

Przykładami układów zupełnie całkowalnych jest zagadnienie Keplera i potok geodezyjny na powierzchni elipsoidy (patrz [4]); oba mają dwa stopnie swobody.

Dla układów spełniających warunek z Definicji 4.1 zachodzi następujące twierdzenie, które przytaczamy bez dowodu (patrz [4]).

Twierdzenie 4.2 (Liouville–Arnold). Jeśli wspólne poziomice {F1=c1,…, Fn=cn} zupełnie całkowalnego układu hamiltonowskiego są zwarte i gładkie, to są one torusami Tn.

Ponadto  w otoczeniu danego takiego torusa istnieje układ współrzędnych I1,…,In,φ1,…,φn, tzw. zmienne działanie–kąt, w których układ (4.2) przyjmuje następującą postać hamiltonowską

I˙j=0,   ⁢φ˙j=ωj⁢I=∂⁡H0/∂⁡Ij,     ⁢j=1,…,n, (4.3)

gdzie H⁢q,p=H0⁢I1,…,In jest hamiltonianem po zamianie. W szczgólności ruch na torusach I1=d1,…,In=dn, które są parametryzowane przez kąty φj mod 2⁢π, jest okresowy lub prawie-okresowy (patrz Rysunek 4.1):

φj⁢t=φj⁢0+ωj⁢I⁢t.

Przykład 4.3. Dla układu van der Pola z ε=0 i H=12⁢x2+y2 zmienne działanie–kąt są następujące: I=H i φ=arg⁡x+i⁢y.

Rozważmy teraz następujące zaburzenie układu (4.3)

I˙=ε⁢g⁢I,φ,   ⁢φ˙=ω⁢I+ε⁢f⁢I,φ, (4.4)

gdzie I=I1,…,In, φ=φ1,…,φn and ω=ω1,…,ωn. Naturalne jest spodziewać się, że rozwiązanie układu (4.4) po czasie rzędu O⁢1 różni sę od rozwiązania układu (4.3) z tymi samymi warunkami początkowymi o wielkość rzędu O⁢ε. Tymczasem poniższe Twierdzenie 4.4 mówi, że taką samą wielkość O⁢ε można uzyskać po czasie dążącym do nieskończoności przy ε→0. Tego rodzaju zjawisko ma miejsce dzięki tzw. uśrednieniu.

Rys. 4.1. Dynamika prawie okresowa na torusie.

Idea uśrednienia wiąże się z faktem, że na większości torusów Tn={I=d} trajektori układu niezaburzonego jest gęsta (jak na Rysunku 4.1). Zatem średnie odchylenie działania I⁢t można wyliczyć (w przybliżeniu) poprzez uśrednienie po torusie.

Definijemy układ uśredniony

J˙=ε⁢G⁢J, (4.5)

gdzie

G⁢J=12⁢πn⁢∫02⁢π…⁢∫02⁢πg⁢J,φ⁢d⁢φ1⁢…⁢d⁢φn

jest uśrednioną po Tn wielkością prędkości zmian działania.

Twierdzenie 4.4 (O uśrednianiu). Niech n=1 i funkcje ω,f,g będą klasy C1 i ω⁢I>0 na otwartym podzbiorze R1×T1. Jeśli I⁢t,φ⁢t i J⁢t,ψ⁢t są rozwiązaniami układów (4.4) i (4.5) takimi, że I⁢0=J⁢0, to dla

0<t<1/ε

amy

I⁢t-J⁢t<C⋅ε,

gdzie stała C zależy tylko odω,f,g.

Dowód. Dokonajmy zamiany

K=I+ε⁢k⁢I,φ (4.6)

tak, aby zachodziło K˙=O⁢ε2. Wyliczenie k⁢J,φ przebiega następująco

K˙=I˙+ε⁢∂⁡k∂⁡I⁢I˙+ε⁢∂⁡k∂⁡φ⁢φ˙=ε⁢g+∂⁡k∂⁡φ⁢ω+O⁢ε2=ε⁢g⁢K,φ+∂⁡k∂⁡φ⁢K,φ⁢ω⁢K+O⁢ε2.

Zatem chcemy rozwiązać równanie

∂⁡k∂⁡φ⁢K,φ⁢ω⁢K=-g⁢K,φ,

z oczywistym rozwiązaniem g⁢K,φ=-1ω⁢K⁢∫0φg⁢K,ψ⁢d⁢ψ. Niestety, na ogół to rozwiązanie nie jest jednoznaczną (czyli okresową) funkcją od φ. Przeszkodą jest wielkość ∫02⁢πg⁢h,ψ⁢d⁢ψ, która może być niezerowa.

Ale, zapisując

g⁢K,φ=G⁢K+g~⁢K,φ

tak, że ∫02⁢πg~⁢h,ψ⁢d⁢ψ=0, możemy zdefiniować jednoznaczną funkcję

g⁢K,φ=-1ω⁢K⁢∫0φg~⁢K,ψ⁢d⁢ψ.

Dostajemy równanie

K˙=ε⁢G⁢K+O⁢ε2.

Widać, że po czasie O⁢1/ε różnica pomiędzy J⁢t i K⁢t jest rzędu O⁢ε. Z drugiej strony, różnica pomiędzy K⁢t i I⁢t jest rzędu O⁢ε, dzięki zamianie (4.6). ∎

Dla zaburzeń typu (4.4) zupełnie całkowalnych układów hamiltonoskich z wieloma stopniami swobody oszacowania są słabsze niż w tezie Twierdzenia 4.4. Okazuje się, że po czasie czędu O⁢1/εa, dla warunków początkowych spoza zbioru o mierze Lebesque's O⁢εb, odchylenie J⁢t od I⁢t nie przekracza O⁢εc, gdzie a,b,c>0 są wykładnikami zależnymi od ω,f,g. Po więcej informacji odsyłam czytelnika do [5].

Przykład 4.5 (Całki abelowe). Rozważmy następujące zaburzenie dwuwymiarowego układu hamiltonowskiego

x˙=Hx′+ε⁢P⁢x,y,    ⁢y˙=-Hx′+ε⁢Q⁢x,y.

Dla ε=0 krzywe fazowe leżą w poziomicach funkcji Hamiktona H⁢x,y. W pewnym obszarze przestrzeni fazowej te krzywe są zamknięte. Jak już robiliśmy to kilkakrotnie, badanie cykli granicznych układu zaburzonego polega na analizie przekształcenia powrotu Poincarégo z⁢S (cięcie transweralne do krzywych fazowych) do S. Parametryzując S za pomocą HS warunek cyklu granicznego to Δ⁢H=H⁢B-H⁢A=0 (patrz Rysunek 4.2). Mamy

Rys. 4.2. Przekształcenie powrotu dla zaburzenia układu hamiltonowskiego.
Δ⁢H=∫0Td⁢Hd⁢tdt=∫0T{Hx′(Hy′+εP)+Hy′(-Hx′+εQ}dt=ε⁢∫P⁢H′⁢x+Q⁢Hy′⁢d⁢t=ε⁢∫P⁢Hx′-ε⁢Q+Q⁢Hy′+ε⁢P=ε⁢∫Γ⁢hQ⁢d⁢x-P⁢d⁢y=ε⁢∮H=hQ⁢d⁢x-P⁢d⁢y+O⁢ε2,

gdzie T jest czasem powrotu do S a Γ⁢h jest krzywą fazową układu zaburzonego startującą z A∈S takiego, że H⁢A=h.

Wyrażenie

I⁢h=∮H=hQ⁢d⁢x-P⁢d⁢y (4.7)

jest tzw. całką abelową.15Pojęcie całki abelowej wywodzi się z zespolonej geometrii algebraicznej. Są to całki z 1-form meromorficznych wzdłuż pewnych zamkniętych krzywych na zespolonych krzywych algebraicznych (powierzchniach Riemanna). Gdy H⁢x,y, P⁢x,y i Q(x,y są wielomianami, to powierzchnia Riemanna jest zespoloną krzywą {H(x,y)=h}⊂C2 a 1-forma to ω=Q⁢d⁢x-P⁢d⁢yH=h. Z Twierdzenia o funkcjach uwikłanych wynika, że jeśli I⁢h0=0 i I′⁢h0≠0, to dla ε≠0 i małego istnieje cykl graniczny γε, który dąży do krzywej H=h0 przy ε→0. To podejście do problemu cykli granicznych jest szeroko stosowane w Jakościowej Teorii.

Nietrudno zauważyć, że funkcja I⁢h jest odpowiednikiem całki uśrednienia G⁢J, która występuje we wzorze (4.5).

4.2. Teoria KAM

Rozważmy układ hamiltonowski

q˙=Hp′,    ⁢p˙=-Hq′,

p=p1,…,pn,q=q1,…,qn, z hamiltonianem postaci

H⁢p,q=H0⁢p,q+ε⁢H1⁢p,q,

gdzie H0 jest hamiltonianem układu zupełnie całkowalnego, czyli układu typu (4.3) w zmiennych działanie–kąt. Z tą sytuacją wiąże się jedno z najważniejszych twierdzeń matematycznych drugiej połowy XX wieku. Przed jego sformułowaniem musimy wprowadzić jeszcze dwa założenia dotyczące niezdegenerowania niezaburzonego hamiltonianu H0:

det⁡∂⁡ωi∂⁡Ij=det⁡∂2⁡H0∂⁡Ii⁢∂⁡Ij≠0, (4.8)
det⁡∂2⁡H0∂⁡Ii⁢∂⁡Ij∂⁡H0∂⁡Ii∂⁡H0∂⁡Ij0. (4.9)

Warunek (4.8) oznacza, że częstości ωi⁢I znieniają się niezależnie i dosyć szybko wraz ze zmianą działań Ij, natomiast warunek (4.9) oznacza, że te częstości zmianiają się szybko i w miarę niezależnie po ograniczeniu do poziomic H0=c⁢o⁢n⁢s⁢t.

Twierdzenie 4.6 (Kołmogorow–Arnold–Moser). Jeśli są spełnione warunki niezdegenerowania (4.8) i (4.9) dla H0, to dla małego zaburzenia H=H0+ε⁢H1 większość torusów niezmienniczych I=c⁢o⁢n⁢s⁢t nie znika, ale tylko lekko deformuje się i ruch na nich jest dalej prawie okresowy.

To twierdzenie zostało sformułowane w 1954 roku na Międzynarodowym Kongresie Matematyków w Amsterdamie, ale na ścisły dowód musiało czekać do początku lat sześdziesiątych. Podali go niezależnie V. Arnold (w przypadku analitycznym) i J. Moser (w przypadku gładkim klasy C333). Później klasa gladkości została obniżona do C3. Oczywiście nie jestem w stanie przedstawić tego dowodu tutaj.

Przykład 4.7. (Płaski ograniczony problem trzech ciał) Płaskie ograniczone zagadnienie trzech ciał jest to układ w którym dwa ciała (oddziałujące na siebie siłą grawitacji) obracają się za stała prędkością kątową wokół ich środka masy (w początku układu współrzędnych) a trzecie ciało porusza się w płaszczyźnie obrotu dwu ciał i ma masę tak małą, że nie zakłóca ich ruchu. Na Rysunku 4.3 mamy taki układ, w którym S oznacza Słońce, J -Jowisz, A zaś jest Asteroidem. Jednostki czasu, długości i masy można dobrać tak, aby prędkość kątowa, suma mas S i J oraz stała grawitacyjna były równe 1. Wtedy też odległość między S i J też równa się 1. Jedynym parametrem charakteryzującym układ jest masa Jowisza μ.

Równania ruchu Asteroidu są hamiltonowskie z hamiltonianem

12⁢p12+p22-1-μρ1-μρ2,

gdzie ρ1 i ρ2 są odległościami A od S i J odpowiednio. Zauważmy, że położenia S i J zmieniają się z czasem: J=1-μ⁢cos⁡t,sin⁡t, S=-μ⁢cos⁡t,sin⁡t; zatem hamiltoniam zależy bezpośredno od czasu.

Rys. 4.3. Problem trzech ciał i niezmiennicze torusy.

Aby pozbyć się tej zależności od czasu, dokonujemy następującej zamiany (jednoczesny obrót współrzędnych i pędów)

q′=M⁢t⁢q,   ⁢p′=M⁢t⁢p,   ⁢M⁢t=cos⁡tsin⁡t-sin⁡tcos⁡t.

Okazuje się, że w nowych zmiennych układ nadal jest hamiltonowski z nowym hamiltonianem

H=12⁢p1′+q2′2+12⁢p2′-q1′2-V⁢q1′,q2′, (4.10)
V=q1′A2+q2′A22+1-μρ1+μρ2,

ρ12=q1′+μ2+q2′A2, ρ22=q1′+μ-12+q2′A2 (Zadanie 4.13). W nowych zmiennych q1′,q2′ ciała⁢S i J spoczywają.

Punkty równowagi układu hamiltonowsiego to punkty krytyczne funkcji hamiltona (Zadanie 4.14). W przypadku hamiltonianu (4.10) te punkty, które nazywamy względnymi położeniami równowagi, zadane są przez

p1′=-q2′,   ⁢p2′=q1′,   ⁢∂⁡V/∂⁡q1′=∂⁡V/∂⁡q2′=0.

Mamy

∂⁡V∂⁡q2′=q2′⁢1-1-μρ13-μρ23=q2′⁢f,
∂⁡V∂⁡q2′=q1′⁢f-μ⁢1-μ⁢1ρ13-1ρ23.

Mamy dwie możliwości:

1. q2′=0; tutaj znajdujemy trzy punkty tzw. współliniowe punkty libracji L1, L2, L3 (Zadanie 4.15), które okazują się niestabilne.

2. f=0 i ρ1=ρ2=1; tutaj mamy dwa tzw. trójkątne punkty libracji L4 i L5, które leżą w wierzchołkach dwóch trójkątów równobocznych o podstawie S⁢J¯.

Wyliczenia, których nie przeprowadzamy, pokazują, że dla 27⁢μ⁢1-μ>1 punkty L4,5 są niestabilne natomiast w przeciwnym przypadku, tj. dla μ<μ1=12⁢1-23/27≈0.03852, wartości własne części liniowej układu hamiltonowskiego są postaci ±i⁢ω1, ±i⁢ω2, gdzie ω1<0<ω2≠ω1. Jesteśmy na granicy obszaru stabilności.

Ponadto część kwadratowa H w punkcie L4 przyjmuje postać

H0=12⁢ω1⁢p~12+q~12+12⁢ω2⁢p~22+q~22

w odpowiednim układzie wspólrzędnych w otoczeniu L4 (patrz [19]). Jest to Hamiltonian układu zupełnie całkowalnego ze zmiannymi działanie–kąt I1=12⁢p~12+q~12, I2=12⁢p~22+q~22, φ1=arg⁡q~1+i⁢p~1, φ2=arg(q~2+ip~2) i z H0=ω1⁢I1+ω2⁢I2 (Zadanie 4.16).

Mamy sytuację jak w Twierdzeniu KAM: H=H0+H1, gdzie H0 jest zupełnie całkowalny a H1 zawiera wyrazy rzędu >2 ze względu na Ij (które są małe). Niestety, to nie wystarcza, ponieważ częstości ωj=∂⁡H0/∂⁡Ij są stałe, a z warunku niezdegenerowania (4.8) powinny się zmianiać wraz z Ij. Należy więc uwzględnić jeszcze dalsze wyrazy rozwinięcia H w otoczeniu L4.

Dokładniej, dokonujemy uproszczenia wyrazów rzędu trzeciego i czwartego w hamiltonianie H. To uproszczenie jest analogiem formy normalnej Poincarégo–Dulaca i zostało udowodnione przez G. Birkhoffa w Twierdzeniu 4.9 poniżej. Ta forma normalna  Birkhoffa w naszym przypadku ma następującą postać

H=H0+H1,  ⁢H0=ω1⁢I1+ω2⁢I2+∑ωi⁢j⁢Ii⁢Ij,  ⁢Ij=12⁢Pj2+Qj2, (4.11)

gdzie Pj=p~j+…, Qj=q~j+… są nowymi zmiennymi a H1 zawiera wyrazy rzędu ≥5 (oraz H0 i H1 są inne niż powyżej). W założeniu twierdzenia Birkhoffa pojawia się warunek braku relacji rezonansowych rzędu 4 i 3. Okazuje się, że takie relacje zachodzą dla wartości μ2=12⁢1-1833/45≈0.02429 i μ3=12⁢1-213/15≈0.01352; zatem te wartości parametru μ należy wykluczyć.

Hamiltonian H0=H0⁢I1,I2 jest hamiltonianem zupełnie całkowalnym i ma szansę na spełnienie warunków niezdegenerowania (4.8) i (4.9). Okazuje się, że tylko warunek (4.9) jest istotny. A. Leontowicz pokazał, że może on zostać naruszony tylko dla dyskretnego zbioru wartości parametru μ.W [19] można dowiedzieć się, że warunek (4.9) zostaje naruszony dla dokłanie jednej konkretnej wartości parametru μ. Ta wartość wynikała ze wzoru na wyznacznik w równaniu (4.9) podanego przez francuskich astronomów A. Deprit i A. Deprit-Bartholomé (i cytowanego w bardzo poważnych monografiach). Ostatnio z moją magistrantką W. Barwicz odkryliśmy, że ten wzór jest nieprawdziwy, a nawet sprzeczny z wyliczeniami Leontowicza. W istocie, tenże wyznacznik jest bardzo skomplikowaną funkcją algebraiczną od μ, króra nie jest tożsamościowo równa zeru. Załóżmy zatem, że μ spełnia wszystkie warunki wypisane powyżej, czyli jest prawdziwa teza twierdzenia KAM.

Jak z twierdzenia KAM wynika stabilność? Otóż znajdujemy się w przestrzeni 4-wymiarowej w otoczeniu punktu równowagi. Ponieważ układ jest hamiltonowski z hamiltonianem niezależnym od czasu, więc ruch odbywa się po powierzchniach H=const. Są one trójwymiarowe. Z twierdzenia KAM wynika, że każda taka powierzchnia jest prawie zapełniona torusami niezmienniczymi T2, których jest tym więcej im bliżej jesteśmy torusa I1=I2=0. Każdy torus niezmienniczy rozbija powierzchnię H=const na dwie części, swoje wnętrze i zewnętrze. Zaden punkt z wnętrza nie wychodzi zeń w trakcie ewolucji. Ponieważ w przestrzeni zmiennych P,Q torusy mogą być dowolnie małe, to wynika stąd stabilność w sensie Lapunowa.

Uzupełnimy powyższy przykład. Załóżmy, że mamy hamiltonian w postaci

H=∑ωj⋅12⁢pj2+qj2+….

Definicja 4.8. Mówimy, że `częstości' ωj spełniają relację resonansową rzędu d, jeśli isnieją liczby całkowite k1,…,kn z ∑kj=d takie, że

∑kj⁢ωj=0.

Twierdzenie 4.9 (Birkhoff). Jeśli częstości ωj nie spełniają żadnej relacji rezonansowej rzędu ≤2⁢m, to istnieje kanoniczna zamiana zmiennych p,q⟼P,Q=p+…,q+… prowadząca do hamiltonianu

H=∑l≤mal⁢Il+O⁢p,q2⁢m+1,

gdzie Ij=12⁢Pj2+Qj2 i sumowanie przebiega po wielowskaźmikach l1,…,ln z l=l1+…+ln iIl=I1l1⁢…⁢Inln.

Uwaga 4.10. Zamiana p,q⟼P,Q, występująca w powyższym twierdzeniu jest kanoniczna jeśli

∑d⁢pj∧d⁢qj=∑d⁢Pj∧d⁢Qj.

Okazuje się, że po kanonicznej zamianie zmiennych układ hamiltonowski przechodzi w układ hamiltonowski (patrz [4]).

ZADANIA

Zadanie 4.11. Pokazać, że jeśli funkcja hamiltona H nie zależy bezpośrednio od czasu, to jest całką pierwszą dla układu (4.2).

Zadanie 4.12. Pokazać, że pole wektorowe zadane wzorem (4.2) ma zerową dywergencję. Wywnioskować stąd, że odpowiedni potok fazowy zachowuje objętość.

Zadanie 4.13. Udowodnić wzór (4.10).

Zadanie 4.14. Pokazać, że jeśli H nie zależy bezpośrednio od czasu, to punkty równowagi układu (4.2) są dokładnie punktami krytycznymi funkcji H.

Zadanie 4.15. Pokazać, że istnieją dokładnie trzy współliniowe punkty libracji.

Zadanie 4.16. Pokazać, że hamiltonian postaci H0=ω1⁢I1+ω2⁢I2 (lub jak we wzorze (4.11)) jest hamiltonianem układu zupełnie całkowalnego.

Zadanie 4.17. Zastosować metodę całek abelowych (Przykład 4.5) do pokazania, że układ van der Pola x˙=y, y˙=-x-a⁢x2-1⁢y dla małego parametru a>0 posiada dokładnie jeden cykl graniczny.

4.3. Drgania relaksacyjne

Zacznijmy od znanego przykładu.

Przykład 4.18 (Układ van der Pola).

x˙=y-x3+x,   ⁢y˙=-ε⁢x.

(Gdy ε=1 i położyć y1=y-x3+x, to dostaje się x˙=y1, y˙1=-x-3⁢x2-1⁢y1; z dokładnością do przeskalowania jest to układ z Przykładu 2.34.)

Rys. 4.4. Układ van der Pola typu `wolny–szybki'.

Widać, że x zmienia się szybko w porównaniu z y; mówimy, że x jest szybką zmienną a y wolną. Dla ε=0 mamy y=const i w istocie mamy równanie na x zależne od parametru y (teoria bifurkacji się kłania, patrz Rysunek 4.4). Gdy ε≠0 (ale małe), to fizycy powiedzieliby, że parametr y `płynie'. Oczekuje się istnienia cyklu granicznego γε (w istocie γε jest stabilny) dążącego do kawałkmi gładkiej krzywej γ0 przedstawionej na Rysunku 4.5. Cykl γ0 składa się z:

— kawałków ruchu powolnego wzdłuż krzywej y=x3-x (gdzie x˙=0),

— odcinków skoku wzdłuż prostych y=const.

Taki ruch jest przykładem dragań relaksacyjnych (jak bicie serca).

Rys. 4.5. Drgania relaksacyjne.

Rozważmy teraz ogólną sytuację. Mamy układ niezaburzony

x˙=f⁢x,y,    ⁢y˙=0,

(x∈Rk, y∈Rl); tutaj x to szybkie współrzędne a y to wolne współrzędne. Mamy też układ zaburzony

x˙=F⁢x,y;ε,   ⁢y˙=ε⁢G⁢x,y;ε,    ⁢F⁢x,y;0=f⁢x,y.

Definicja 4.19. Powierzchnia S={f(x,y)=0} nazywa się powolną powierzchnią.

Powolna powierzchnia dzieli się na obszary stabilności i niestabilności układu niezaburzonego; odpowiadają one sytuacjom gdy Re⁢λj⁢A<0, j=1,…,k, A=∂⁡f∂⁡x, i gdy istnieje Re⁢λj⁢A>0.

Na powolnej powierzchni mamy pole wektorowe definiowane następująco. Bierzemy pole

∂∂⁡ε⁢F⁢∂x+ε⁢G⁢∂yε=0=f1⁢x,y⁢∂x+g⁢x,y⁢∂y

w punkcie x,y∈S i rzutujemy je na Tx,y⁢S wzdłuż zmiennych y. Jest to pole ruchu powolnego.

Przypomnę, że na początku tego rozdziału mówiliśmy, że drgania relaksacyjne charakteryzują się własnością, ze mały parametr występuje po lewej stronie. Aby się o tym przekonać wprowadzamy czas powolny τ=ε⁢t. Wtedy dostajemy układ

ε⁢d⁢xd⁢τ=f⁢x,y+O⁢ε,   ⁢d⁢yd⁢τ=g⁢x,y+O⁢ε.

Teraz równanie ruchu powolnego na S (lokalnie parametryzowanej przez y) jest postaci

d⁢yd⁢τ=h⁢y+O⁢ε

(z odpowiednią funkcją h).

Przeanalizujmy ruch typowego punktu x0,y0. Składa się on z kawałków trzech rodzajów: dochodzenie do powierzchni powolnej, ruch wzdłuż powierzchnii powolnej i ruch w obszarze przejściowym.

Rys. 4.6. Dochodzenie do powierzchni powolnej.

4.20. Dochodzenie do powierzchni powolnej. Niech punkt x0,y0 spoza S rzutuje się (wzdłuż współrzędnych y) na punkt x∗,y0, x∗=x∗⁢y0 , na S w obszarze stabilności (patrz Rysunek 4.6). To znaczy, że punkt x0 leży w basenie przyciągania punktu x∗ dla równania x˙=f⁢x,y0 (y0 stałe). Rozważmy obszar U={|x-x∗(y0)|<δ, y0∈V}, gdzie V jest pewnym obszarem odpowiadającum podzbiorowi obszaru stabilności w S. Okazuje się, że powolny czas dochodzenia rozwiązania z warunkiem początkowym x0,y0 do U jest rzędu τ1∼C1⁢ε⁢ln⁡ε, co odpowiada rzeczywistemu czasowi

t1∼C1⁢ln⁡ε

(stała C1 zależy od U i od F,G).

4.21. Ruch powolny. W obszarze U mamy ruch powolny, opisywany równaniem d⁢y/d⁢τ=h⁢y+O⁢ε. Trwa on do momentu τ2=T=O⁢1, co odpowiada długiemu czasowi rzeczywistemu t2=T/ε.

4.22. Ruch w obszarze przejściowym. Obszar przejściowy leży blisko granicy pomiędzy obszarami stabilności i niestabilności w S. Mamy dwie typowe możliwości (jak w teorii bifurkacji):

A.λ1⁢A=0 (gdzie A=∂⁡f∂⁡x|f=0);

B.Re⁢λ1,2=0.

A. Zryw. Ten przypadek (który odpowiada bifurkacji siodło–węzeł) zanalizujemy dla sytuacji gdy x∈R i y∈R (można do tego wszystko zredukować). Po odpowiednich przeskalowaniach mamy następujący układ

Rys. 4.7. Zjawisko `zrywu'.
x˙=x2-y+…,   ⁢y˙=-ε+…

Dokonujemy normalizacji

ε=μ3,   ⁢x=μ⁢X,   ⁢y=μ2⁢Y.

łatwo sprawdzić, że prowadzi to do pola

X˙=μ⁢X2-Y+O⁢μ,   ⁢Y˙=μ⁢-1+O⁢μ

orbitalnie równoważnemu polu X2-Y⁢∂X-∂Y. Jego portret fazowy jest zadany równaniem Riccatiego

d⁢X/d⁢Y=X2-Y (4.12)

i jest przedstawiony na Rysunku 4.717Równanie (4.12) jest chyba najprostszym przykładem równania różniczkowego, którego nie można rozwiązać w tzw. kwadraturach.. Zjawisko, które tutaj obserwujemy nosi nazwę zrywu.

B. Opóźnienie utraty stabilności. W tym przypadku, który odpowiada bifurkacji Andronowa–Hopfa, problem redukuje się do następującego modelowego układu

z˙=y+i⁢ω⁢z+c⁢z⁢z2,   ⁢y˙=ε, (4.13)

z⁢x1+i⁢x2∈C≃R2,y∈R. Oczywiście y=ε⁢t jest `płynącym' parametrem. Załóżmy jeszcze, że

c=-1;

przypadek c>0 jest mniej ciekawy. Dla amplitudy r=z dostajemy równanie Bernoulliego

r˙=r⁢ε⁢t-r2.

Połóżmy warunek początkowy

y⁢t0=-μ,    ⁢r⁢t0=r0,   ⁢t0=-μ/ε,

gdzie μ>0 jest ustaloną (nie za dużą i nie za małą) stałą. To zagadnienie początkowe ma następujące rozwiązanie

r⁢t=r0⁢eε⁢t02-t2+2⁢r02⁢∫t0teε⁢s2-t2⁢d⁢s-1/2 (4.14)

(Zadanie 4.24). Zbadamy asymptotyczne zachowanie się tego rozwiązania przy ε→0 dzieląc zakres czasu t na cztery obszary:

(a) 0<t-t0<O⁢1, czyli 0<y+μ<O⁢ε.

Niech u=t-t0. Wtedy ε⁢t02-t2=ε⁢t0+t⁢u≈2⁢μ⁢u i ε⁢s2-t2≈2⁢μ⁢u-v, gdzie v=s-t0. Zatem

∫t0teε⁢s2-t2⁢d⁢s≈∫0ue2⁢μ⁢u-v⁢d⁢v=12⁢μ⁢e2⁢μ⁢u-1

oraz

r(t)≈r0{e2⁢μ⁢u+r02(e2⁢μ⁢u-1)/μ}-1/2

jest malejącą funkcją od u.

(b) y=ε⁢t jest ustalone tak, że -μ<y<μ.

Tutaj eε⁢t02-t2≈eμ2-y2/ε→∞. Zatem

r⁢t<C1⁢e-C2/ε→0,

przy czym jest to bardzo szybkie dążenie do zera.

(c) 0<t0-t<O⁢1, czyli 0<μ-y<O⁢ε.

Wprowadźmy zmienną w=t0-t. Jak w punkcie (a) mamy eε⁢t02-t2≈e2⁢μ⁢w.

Obszar całkowania dla całki we wzorze (4.14) podzielimy na trzy odcinki: od t0 do t0/2<0, od t0/2 do t0/2 i od t0/2 do t. Przez I1, I2 i I3 oznaczymy odpowiednie całki. Podobnie jak w punkcie (a) pokazuje się, że I1=O⁢1 i I3=O⁢1. Z rachunków w punkcie (b) wynika, że I2→0 bardzo szybko. Zatem

r⁢t=O⁢1.

(d) t0<t, czyli y>μ i jest ustalone. Teraz exp⁡ε⁢t02-t2≈exp⁡-y2-μ2/ε→0. Następnie ε⁢s2-t2≈s-t⋅2⁢y dla s bliskich t, tj. dla tych s, dla których wkład do całki jest dominujący. Dostajemy ∫te2⁢y⁢s-t⁢d⁢s≈12⁢y. Stąd

r⁢t≈r0⁢r02/2⁢y-1/2=y.

Możemy podsumować powyższe obliczenia.

Rys. 4.8. Zjawisko opóźnienia utraty stabilności.

Twierdzenie 4.23.W przypadku B opisywanego układem (4.13) z c<0 zachodzi zjawisko opóźnienia utraty stabilności. Polega ono na tym, przy zmianie zmiennej y (która jest współczynnikiem stabilności ruchu niezaburzonego) od wartości ujemnej y⁢t0=-μ do wartości dodatniej μ układ (względem z) jest cały czas stabilny, a zmiana stabilności rozwiązania następuje dla parametru y=μ, przy czym dalej amblituda oscylacji rośnie jak w zwykłej bifurkacji Andronowa–Hopfa.

Zjawisko opóźnienia utraty stabilności można objaśnić fizycznie.Zmienna y jest ujemna przez bardzo długi czas, rzędu 1/ε. Wtedy układ fizyczny zdąży podejść bardzo blisko położenia równowagi; na tyle blisko, że potrzeba potem tyle samo czasu, aby od położenia równowago odejść (patrz Rysunek 4.8).

ZADANIA

Zadanie 4.24. Udowodnić wzór (4.14).

Treść automatycznie generowana z plików źródłowych LaTeXa za pomocą oprogramowania wykorzystującego LaTeXML.

Projekt współfinansowany przez Unię Europejską w ramach Europejskiego Funduszu Społecznego.

Projekt współfinansowany przez Ministerstwo Nauki i Szkolnictwa Wyższego i przez Uniwersytet Warszawski.