pobieranie: Modelowanie i Symulacja, Rozwiązywanie równań różniczkowych zwyczajnych metodą klasyczną
Pobierz dokument pdfPrzeglądaj wersję html pliku:
Laboratorium Modelowania i symulacji
2008 r.
Wydział Elektryczny Zespół Automatyki (ZTMAiPC) ZERiA
LABORATORIUM MODELOWANIA I SYMULACJI Ćwiczenie 2
Rozwiązywanie równań róŜniczkowych zwyczajnych metodą klasyczną. Metoda Eulera. 1. Cel ćwiczenia
Celem ćwiczenia jest zapoznanie się z metodami rozwiązywania liniowych równań róŜniczkowych zwyczajnych (w tym metodami numerycznymi) w języku Matlab.
2. Wprowadzenie.
Równaniem róŜniczkowym liniowym nazywamy równanie postaci:
& & a n y ( n ) + a n −1 y ( n −1) + ... + a 2 && + a1 y + a 0 y = b0 x + b1 x + b2 && + bm −1 x ( m −1) + bm x ( m ) y x
w równaniu tym: x – wymuszenie y – odpowiedź
(1)
Rozwiązaniem równania (1) jest całka będąca sumą całki stanowiącej rozwiązanie równania róŜniczkowego jednorodnego oraz jednej z całek szczególnych będącej rozwiązaniem równania róŜniczkowego niejednorodnego.:
y (t ) = y p (t ) + y u (t )
(2) (3)
yp(t) – składowa przejściowa (swobodna) jest rozwiązaniem równania (3):
& a n y ( n ) + a n −1 y ( n −1) + ... + a 2 && + a1 y + a 0 y = 0 y
yu(t) – składowa wymuszona Rozwiązanie równania jednorodnego (II rzędu): Będziemy rozpatrywać równanie postaci:
a
podstawiamy:
d2y dy + b + cy = 0 2 dx dx
(4)
y = e rx ⇒ y ' = re rx , y ' ' = r 2 e rx
mamy:
(5) (6)
ar 2 e rx + bre rx + ce rx = 0
i dalej:
(7) Równanie (7) zwane jest równaniem charakterystycznym. Istnieją następujące warianty rozwiązania: 1. Dwa róŜne pierwiastki (rzeczywiste) r1 i r2
ar 2 + br + c = 0
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-1-
Laboratorium Modelowania i symulacji
2008 r.
y = C1e r1x + C 2 e r 2 x
2. Jeden pierwiastek podwójny r1,2
(8) (9) (10)
y = (C1 x + C 2 )e 1, 2
r x
3. Dwa pierwiastki zespolone, sprzęŜone r1 = α+iβ, r1 = α-iβ Stałe C1, C2 wyznacza się dla zadanych warunków początkowych.
y = eαx (C1 cos βx + C 2 sin βx)
Rozwiązanie równania niejednorodnego (I rzędu) – metoda uzmienniania stałej. Sposób rozwiązania równania niejednorodnego pierwszego rzędu zostanie przedstawiony na przykładzie. Rozpatrzmy równanie:
x
dy − y = 2 ⋅ x3 dx x dy −y=0 dx
(11)
Krok pierwszy: Rozwiązanie równania jednorodnego (metoda rozdzielenia zmiennych). (12)
1 dy 1 = y dx x 1 1 ∫ y dy = ∫ x dx ln( y ) = ln( x ) + C
Podstawienie: C = ln( C1 ); C1 ≠ 0 , daje następującą całkę ogólną równania (11):
(13)
y = C1 ⋅ x
Krok drugi: Uzmiennienie stałej C1: C1 = u(x)
(14)
y = u ( x) ⋅ x
dy du ( x) dy du ( x) = ⋅ x + u ( x) ⇒ = u ( x) + ⋅x dx dx dx dx
Krok trzeci: Podstawienie (15) do (11)
(15)
du ( x) du ( x) 2 x ⋅ u ( x) + ⋅ x − u ( x) ⋅ x = 2 ⋅ x 3 ⇒ ⋅ x = 2 ⋅ x3 dx dx du ( x) = 2⋅ x dx u ( x) = x 2 + C
(16)
Krok czwarty: Podstawienie (16) do (15):
dy dy = (x2 + C) + 2 ⋅ x 2 ⇒ = 3⋅ x2 + C dx dx
Krok Piąty: Rozwiązanie równania (11)
(17)
y = ∫ 3 ⋅ x 2 + C dx = 3∫ x 2 dx + C ∫ dx = x 3 + C ⋅ x
(
)
(18)
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-2-
Laboratorium Modelowania i symulacji
2008 r.
15 C= 3 C= 1 C = -3
10
5
Y
0
-5
-10
-15 -2
-1.5
-1
-0.5
0 X
0.5
1
1.5
2
Wykres 1. Ilustracja rozwiązania równania (11) dla trzech róŜnych wartości współczynnika C
3. Metoda Eulera Rozpatrywać będziemy zagadnienia początkowe, tzn. będziemy chcieli znaleźć rozwiązanie równania róŜniczkowego dla zadanej wartości u0(t0) w punkcie początkowym t0:
du (t ) = f (u (t ), t ) dt u (t 0 ) = u 0
Wprowadźmy oznaczenia: ti = t0+i∆t, ui=u(ti), fi=f(ui,ti). Ogólną metodą na rozwiązanie równania róŜniczkowego, jest zapisanie go w postaci całkowej a następnie odcałkowanie i zastosowanie przybliŜenia na całkę występującą po lewej stronie. W metodzie Eulera pochodną zastępuje się ilorazem róŜnicowym w przód opartym na węzłach tn i tn+1. Całka po lewej stronie przybliŜana jest w związku z tym iloczynem wartości funkcji w początku przedziału i jego długości. Dla metody Eulera mamy:
u n +1 − u n = f (u n , t n ) ⇒ u n+1 = u n + ∆tf (u n , t n ) ∆t
(19)
Wzór przedstawia zaleŜność pomiędzy wartością następną a poprzednią wyznaczanego rozwiązania. ∆t oznacza krok całkowania (dyskretyzacji). Ulepszona metoda Eulera (MIDPOINT) Z uwagi na wolną zbieŜność metody Eulera, aby zachować duŜą dokładność obliczeń tą metodą trzeba stosować bardzo mały krok całkowania. Zwiększa to ilość wykonywanych operacji a w następstwie wydłuŜa czas potrzebny na uzyskanie rozwiązania. Zwiększeniu ulega teŜ wymagana ilość pamięci operacyjnej niezbędna do wykonania całkowania. PowyŜsze niedogodności powodują, Ŝe częściej stosuje się ulepszoną metodę Eulera (MIDPOINT). Polega ona na wprowadzeniu dodatkowego punktu (środek przedziału). Odpowiednie wzory przedstawiono poniŜej:
u
n+
1 2
= un +
∆t f (u n , t n ) 2
n+ 1 2
u n +1 = u n + ∆t ⋅ f (u
(20)
n+ 1 2
,t
)
W metodzie tej wartość funkcji dla tn+1/2 oblicza się z dwa razy mniejszym krokiem. Są równieŜ moŜliwe inne modyfikacje metody Eulera.
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych -3-
Laboratorium Modelowania i symulacji
2008 r.
4. Program ćwiczenia: W ćwiczeniu rozwiązywane będą dwa równania róŜniczkowe: & a) y + 2 y = 10 + 10 cos 2t przy warunku pocz. y (0) = 2
(21) (22)
& & b) && + 2 y + 5 y = t y przy warunku pocz. y (0) = 0, y (0) = 2 Zostaną zaprezentowane trzy sposoby rozwiązania tych równań: - metodą analityczną
-
-
metodą Eulera przy wykorzystaniu istniejących w systemie Matlab funkcji wspierających rozwiązywanie równań róŜniczkowych (ODE – Ordinary Differential Equations)
I. Rozwiązywanie równania róŜniczkowego (21): 1. Rozwiązanie równania o współczynnikach stałych met. klasyczną: a) Rozwiązanie: rozwiązanie ogólne (równanie jednorodne): równanie charakterystyczne:
& y + 2 y = 10 + 10 cos 2t , y (0) = 2
r + 2 = 0 ⇒ r = −2
całka ogólna (rozwiązanie równania jednorodnego):
y p (t ) = C1e −2⋅t
yp - składowa swobodna. przewidujemy całkę szczególną postaci:
A + B cos(2 ⋅ t ) + C ⋅ sin(2t )
Po podstawieniu całki szczególnej do postaci ogólnej rozwiązywanego równania otrzymano: A = 5, B = 5/2, C = 5/2 Z warunku początkowego:, C1= -11/2, więc:
y (t ) = −
11 − 2t 5 e + 5 + (cos 2t + sin 2t ) 2 2
b) W oknie poleceń Matlaba wprowadzić następujące polecenia: >> t1 = 0:0.1:10 ; % wektor wartości zmiennej niezaleŜnej >> y1 = (-11/2)*exp(-2.*t1)+5+(5/2)*(cos(2.*t1)+ ... sin(2.*t1)); % wektor rozwiązań >> plot(t1, y1); % wykres Opisać i skopiować wykres. 2. Rozwiązanie równania róŜniczkowego metodą Eulera: a) W systemie Matlab utworzyć nowy m-plik (File | New | M-file). Wprowadzić następujący skrypt: deltat = 0.1; t_euler = 0:deltat:10; y_euler = zeros(size(t_euler)); y_euler(1,1) = 2 ; for i = 1:length(t_euler)-1 y_euler(1,i+1) = y_euler(1,i) + deltat.*(10+10.* ... cos(2.*t_euler(1,i)) - 2.*y_euler(1,i)); end
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-4-
Laboratorium Modelowania i symulacji
2008 r.
Zapisać skrypt w pliku euler_y1.m, w katalogu podanym przez prowadzącego. Z linii poleceń Matlaba wpisać: >> euler_y1 ; >> plot(t_euler, y_euler); Opisać i skopiować wykres. Dokonać zmian w skrypcie zgodnie ze wskazówkami prowadzącego (zmiana kroku całkowania, ew. zmiana całkowanej funkcji, warunku początkowego, zmiana metody („midpoint”)) Korzystając ze znanych poleceń umieścić na jednym rysunku i odpowiednio opisać wykresy będące rozwiązaniami równania (21) metodą klasyczną i metodą Eulera. 3. Rozwiązanie równania róŜniczkowego metodą ODE: Zadanie 1 Dla równania (21) napisać program rozwiązujący je metodą ODE
II. Rozwiązywanie równania róŜniczkowego (22): 1. Rozwiązanie równania o współczynnikach stałych met. klasyczną: c) Rozwiązanie: Rozwiązanie ogólne (równanie jednorodne): Równanie charakterystyczne:
&& + 2 y + 5 y = t , y (0) = 0, y (0) = 2 & & y
r2 + 2⋅r + 5 = 0
WyróŜnik równania charakterystycznego jest mniejszy od zera, a więc równanie posiada dwa pierwiastki zespolone, sprzęŜone:
r1 = −1 + 2 j , r2 = −1 − 2 j
całka ogólna (rozwiązanie równania jednorodnego):
y p (t ) = e − t ⋅ (C1 cos(2 ⋅ t ) + C 2 ⋅ sin(2 ⋅ t ) )
A⋅t + B
Z uwagi na to, Ŝe po prawej stronie rozwiązywanego równania występuje liniowa funkcja t przewidujemy całkę szczególną postaci: Podstawiając całkę szczególną do postaci ogólnej, otrzymano następujące stałe:
1 2 A = ,B = − . 5 25
Z warunków początkowych: C1=2/25, C2=47/50, więc:
yu (t ) =
i ostatecznie:
1 2 ⋅t − , 5 25
y (t ) = y p (t ) + y u (t ) =
1 2 ⋅t − + e −t 5 25
47 2 ⋅ cos(2 ⋅ t ) + sin(2 ⋅ t ) 50 25
d) W oknie poleceń Matlaba wprowadzić następujące polecenia: >> t2 = 0:0.1:10; >> y2 = (1/5).*t2-(2/25)+exp(-t2).*((2/25).*cos(2.*t2)+ ... (47/50).*sin(2.*t2)); >> plot(t2,y2); 2. Rozwiązanie równania róŜniczkowego metodą Eulera:
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-5-
Laboratorium Modelowania i symulacji
2008 r.
Zadanie 2 Napisać algorytm rozwiązania równania róŜniczkowego (22) metodą Eulera. Wskazówka Drugą pochodną, przedstawić w postaci ilorazu róŜnicowego II rzędu. Korzystając ze znanych poleceń umieścić na jednym rysunku i odpowiednio opisać wykresy będące rozwiązaniami równania (22) metodą klasyczną i metodą Eulera. 3. Rozwiązanie równania róŜniczkowego metodą ODE: Aby wykorzystać ODE naleŜy równanie n-tego rzędu zamienić na układ n równań I rzędu. Dla równania II rzędu:
&& + 2 y + 5 y = t , y (0) = 0, y (0) = 2 & & y
dokonuje się zamiany na układ 2 równań I rzędu wykorzystując następujące podstawienia y = y1, & oraz y 2 = y1 . Wtedy powyŜsze równanie moŜna zapisać jako:
& y1 = y 2 & y 2 = −5 y1 − 2 y 2 + t
Wykorzystując edytor systemu Matlab (File | New | M-file) zapisać powyŜszy układ równań w pliku rownanie2.m w następujący sposób: function dy = rownanie2(t,y) dy=[y(2);(-5)*y(1)-2*y(2)+t]; y(n) – oznacza tu pochodną zmiennej y rzędu (n-1), y(1) – oznacza szukaną funkcję y(t). W celu rozwiązania układu równań róŜniczkowych 1-go rzędu naleŜy wpisać z linii poleceń Matlaba: >> [t_ode,y_ode]=ode45('rownanie2',[0 10],[0;2]); >> plot(t_ode, y_ode(:,1), t2, y2) >> legend('ODE', 'Rozwiązanie analityczne'); Pierwszy parametr wywołania funkcji ode45 – nazwa funkcji opisującej rozwiązywany układ równań róŜniczkowych, drugi parametr – zakres zmian zmiennej niezaleŜnej (t), trzeci parametr – kolumnowy wektor warunków początkowych. W wyniku obliczeń zwrócony zostanie: t – wektor wartości zmiennej niezaleŜnej, dla których wyznaczono rozwiązanie, y – macierz wartości funkcji y(t) oraz y’(t). Korzystając z funkcji pomocy zapoznaj się z innymi metodami rozwiązywania równań róŜniczkowych zaimplementowanymi w Matlabie w postaci funkcji. Opracowanie sprawozdania W sprawozdaniu naleŜy umieścić wykresy rozwiązań omawianych równań. Porównać rozwiązania analityczne z numerycznymi dla róŜnych parametrów (krok całkowania, metoda – Eulera, ODE). Omówić jak zmiana parametrów rozwiązań numerycznych wpływa na dokładność rozwiązania. Zadanie 3 Napisać skrypt w języku Matlab rozwiązujący równanie róŜniczkowe (11) – wstęp teoretyczny (dowolna omawiana metoda). Literatura 1. B. Mrozek, Z. Mrozek: MATLAB i Simulink: poradnik uŜytkownika. Helion, Gliwice, 2004. 2. A. Zalewski, R. Cegieła: Matlab - obliczenia numeryczne i ich zastosowania.Wydawnictwo Nakom, Poznań, 2000. 3. J. Brzózka, L. Dorobczyński: Programowanie w Matlab. Wydawnictwo Mikom,Warszawa, 1998.
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-6-
Modelowanie i Symulacja, Rozwiązywanie równań różniczkowych zwyczajnych metodą klasyczną
Laboratorium Modelowania i symulacji
2008 r.
Wydział Elektryczny Zespół Automatyki (ZTMAiPC) ZERiA
LABORATORIUM MODELOWANIA I SYMULACJI Ćwiczenie 2
Rozwiązywanie równań róŜniczkowych zwyczajnych metodą klasyczną. Metoda Eulera. 1. Cel ćwiczenia
Celem ćwiczenia jest zapoznanie się z metodami rozwiązywania liniowych równań róŜniczkowych zwyczajnych (w tym metodami numerycznymi) w języku Matlab.
2. Wprowadzenie.
Równaniem róŜniczkowym liniowym nazywamy równanie postaci:
& & a n y ( n ) + a n −1 y ( n −1) + ... + a 2 && + a1 y + a 0 y = b0 x + b1 x + b2 && + bm −1 x ( m −1) + bm x ( m ) y x
w równaniu tym: x – wymuszenie y – odpowiedź
(1)
Rozwiązaniem równania (1) jest całka będąca sumą całki stanowiącej rozwiązanie równania róŜniczkowego jednorodnego oraz jednej z całek szczególnych będącej rozwiązaniem równania róŜniczkowego niejednorodnego.:
y (t ) = y p (t ) + y u (t )
(2) (3)
yp(t) – składowa przejściowa (swobodna) jest rozwiązaniem równania (3):
& a n y ( n ) + a n −1 y ( n −1) + ... + a 2 && + a1 y + a 0 y = 0 y
yu(t) – składowa wymuszona Rozwiązanie równania jednorodnego (II rzędu): Będziemy rozpatrywać równanie postaci:
a
podstawiamy:
d2y dy + b + cy = 0 2 dx dx
(4)
y = e rx ⇒ y ' = re rx , y ' ' = r 2 e rx
mamy:
(5) (6)
ar 2 e rx + bre rx + ce rx = 0
i dalej:
(7) Równanie (7) zwane jest równaniem charakterystycznym. Istnieją następujące warianty rozwiązania: 1. Dwa róŜne pierwiastki (rzeczywiste) r1 i r2
ar 2 + br + c = 0
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-1-
Laboratorium Modelowania i symulacji
2008 r.
y = C1e r1x + C 2 e r 2 x
2. Jeden pierwiastek podwójny r1,2
(8) (9) (10)
y = (C1 x + C 2 )e 1, 2
r x
3. Dwa pierwiastki zespolone, sprzęŜone r1 = α+iβ, r1 = α-iβ Stałe C1, C2 wyznacza się dla zadanych warunków początkowych.
y = eαx (C1 cos βx + C 2 sin βx)
Rozwiązanie równania niejednorodnego (I rzędu) – metoda uzmienniania stałej. Sposób rozwiązania równania niejednorodnego pierwszego rzędu zostanie przedstawiony na przykładzie. Rozpatrzmy równanie:
x
dy − y = 2 ⋅ x3 dx x dy −y=0 dx
(11)
Krok pierwszy: Rozwiązanie równania jednorodnego (metoda rozdzielenia zmiennych). (12)
1 dy 1 = y dx x 1 1 ∫ y dy = ∫ x dx ln( y ) = ln( x ) + C
Podstawienie: C = ln( C1 ); C1 ≠ 0 , daje następującą całkę ogólną równania (11):
(13)
y = C1 ⋅ x
Krok drugi: Uzmiennienie stałej C1: C1 = u(x)
(14)
y = u ( x) ⋅ x
dy du ( x) dy du ( x) = ⋅ x + u ( x) ⇒ = u ( x) + ⋅x dx dx dx dx
Krok trzeci: Podstawienie (15) do (11)
(15)
du ( x) du ( x) 2 x ⋅ u ( x) + ⋅ x − u ( x) ⋅ x = 2 ⋅ x 3 ⇒ ⋅ x = 2 ⋅ x3 dx dx du ( x) = 2⋅ x dx u ( x) = x 2 + C
(16)
Krok czwarty: Podstawienie (16) do (15):
dy dy = (x2 + C) + 2 ⋅ x 2 ⇒ = 3⋅ x2 + C dx dx
Krok Piąty: Rozwiązanie równania (11)
(17)
y = ∫ 3 ⋅ x 2 + C dx = 3∫ x 2 dx + C ∫ dx = x 3 + C ⋅ x
(
)
(18)
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-2-
Laboratorium Modelowania i symulacji
2008 r.
15 C= 3 C= 1 C = -3
10
5
Y
0
-5
-10
-15 -2
-1.5
-1
-0.5
0 X
0.5
1
1.5
2
Wykres 1. Ilustracja rozwiązania równania (11) dla trzech róŜnych wartości współczynnika C
3. Metoda Eulera Rozpatrywać będziemy zagadnienia początkowe, tzn. będziemy chcieli znaleźć rozwiązanie równania róŜniczkowego dla zadanej wartości u0(t0) w punkcie początkowym t0:
du (t ) = f (u (t ), t ) dt u (t 0 ) = u 0
Wprowadźmy oznaczenia: ti = t0+i∆t, ui=u(ti), fi=f(ui,ti). Ogólną metodą na rozwiązanie równania róŜniczkowego, jest zapisanie go w postaci całkowej a następnie odcałkowanie i zastosowanie przybliŜenia na całkę występującą po lewej stronie. W metodzie Eulera pochodną zastępuje się ilorazem róŜnicowym w przód opartym na węzłach tn i tn+1. Całka po lewej stronie przybliŜana jest w związku z tym iloczynem wartości funkcji w początku przedziału i jego długości. Dla metody Eulera mamy:
u n +1 − u n = f (u n , t n ) ⇒ u n+1 = u n + ∆tf (u n , t n ) ∆t
(19)
Wzór przedstawia zaleŜność pomiędzy wartością następną a poprzednią wyznaczanego rozwiązania. ∆t oznacza krok całkowania (dyskretyzacji). Ulepszona metoda Eulera (MIDPOINT) Z uwagi na wolną zbieŜność metody Eulera, aby zachować duŜą dokładność obliczeń tą metodą trzeba stosować bardzo mały krok całkowania. Zwiększa to ilość wykonywanych operacji a w następstwie wydłuŜa czas potrzebny na uzyskanie rozwiązania. Zwiększeniu ulega teŜ wymagana ilość pamięci operacyjnej niezbędna do wykonania całkowania. PowyŜsze niedogodności powodują, Ŝe częściej stosuje się ulepszoną metodę Eulera (MIDPOINT). Polega ona na wprowadzeniu dodatkowego punktu (środek przedziału). Odpowiednie wzory przedstawiono poniŜej:
u
n+
1 2
= un +
∆t f (u n , t n ) 2
n+ 1 2
u n +1 = u n + ∆t ⋅ f (u
(20)
n+ 1 2
,t
)
W metodzie tej wartość funkcji dla tn+1/2 oblicza się z dwa razy mniejszym krokiem. Są równieŜ moŜliwe inne modyfikacje metody Eulera.
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych -3-
Laboratorium Modelowania i symulacji
2008 r.
4. Program ćwiczenia: W ćwiczeniu rozwiązywane będą dwa równania róŜniczkowe: & a) y + 2 y = 10 + 10 cos 2t przy warunku pocz. y (0) = 2
(21) (22)
& & b) && + 2 y + 5 y = t y przy warunku pocz. y (0) = 0, y (0) = 2 Zostaną zaprezentowane trzy sposoby rozwiązania tych równań: - metodą analityczną
-
-
metodą Eulera przy wykorzystaniu istniejących w systemie Matlab funkcji wspierających rozwiązywanie równań róŜniczkowych (ODE – Ordinary Differential Equations)
I. Rozwiązywanie równania róŜniczkowego (21): 1. Rozwiązanie równania o współczynnikach stałych met. klasyczną: a) Rozwiązanie: rozwiązanie ogólne (równanie jednorodne): równanie charakterystyczne:
& y + 2 y = 10 + 10 cos 2t , y (0) = 2
r + 2 = 0 ⇒ r = −2
całka ogólna (rozwiązanie równania jednorodnego):
y p (t ) = C1e −2⋅t
yp - składowa swobodna. przewidujemy całkę szczególną postaci:
A + B cos(2 ⋅ t ) + C ⋅ sin(2t )
Po podstawieniu całki szczególnej do postaci ogólnej rozwiązywanego równania otrzymano: A = 5, B = 5/2, C = 5/2 Z warunku początkowego:, C1= -11/2, więc:
y (t ) = −
11 − 2t 5 e + 5 + (cos 2t + sin 2t ) 2 2
b) W oknie poleceń Matlaba wprowadzić następujące polecenia: >> t1 = 0:0.1:10 ; % wektor wartości zmiennej niezaleŜnej >> y1 = (-11/2)*exp(-2.*t1)+5+(5/2)*(cos(2.*t1)+ ... sin(2.*t1)); % wektor rozwiązań >> plot(t1, y1); % wykres Opisać i skopiować wykres. 2. Rozwiązanie równania róŜniczkowego metodą Eulera: a) W systemie Matlab utworzyć nowy m-plik (File | New | M-file). Wprowadzić następujący skrypt: deltat = 0.1; t_euler = 0:deltat:10; y_euler = zeros(size(t_euler)); y_euler(1,1) = 2 ; for i = 1:length(t_euler)-1 y_euler(1,i+1) = y_euler(1,i) + deltat.*(10+10.* ... cos(2.*t_euler(1,i)) - 2.*y_euler(1,i)); end
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-4-
Laboratorium Modelowania i symulacji
2008 r.
Zapisać skrypt w pliku euler_y1.m, w katalogu podanym przez prowadzącego. Z linii poleceń Matlaba wpisać: >> euler_y1 ; >> plot(t_euler, y_euler); Opisać i skopiować wykres. Dokonać zmian w skrypcie zgodnie ze wskazówkami prowadzącego (zmiana kroku całkowania, ew. zmiana całkowanej funkcji, warunku początkowego, zmiana metody („midpoint”)) Korzystając ze znanych poleceń umieścić na jednym rysunku i odpowiednio opisać wykresy będące rozwiązaniami równania (21) metodą klasyczną i metodą Eulera. 3. Rozwiązanie równania róŜniczkowego metodą ODE: Zadanie 1 Dla równania (21) napisać program rozwiązujący je metodą ODE
II. Rozwiązywanie równania róŜniczkowego (22): 1. Rozwiązanie równania o współczynnikach stałych met. klasyczną: c) Rozwiązanie: Rozwiązanie ogólne (równanie jednorodne): Równanie charakterystyczne:
&& + 2 y + 5 y = t , y (0) = 0, y (0) = 2 & & y
r2 + 2⋅r + 5 = 0
WyróŜnik równania charakterystycznego jest mniejszy od zera, a więc równanie posiada dwa pierwiastki zespolone, sprzęŜone:
r1 = −1 + 2 j , r2 = −1 − 2 j
całka ogólna (rozwiązanie równania jednorodnego):
y p (t ) = e − t ⋅ (C1 cos(2 ⋅ t ) + C 2 ⋅ sin(2 ⋅ t ) )
A⋅t + B
Z uwagi na to, Ŝe po prawej stronie rozwiązywanego równania występuje liniowa funkcja t przewidujemy całkę szczególną postaci: Podstawiając całkę szczególną do postaci ogólnej, otrzymano następujące stałe:
1 2 A = ,B = − . 5 25
Z warunków początkowych: C1=2/25, C2=47/50, więc:
yu (t ) =
i ostatecznie:
1 2 ⋅t − , 5 25
y (t ) = y p (t ) + y u (t ) =
1 2 ⋅t − + e −t 5 25
47 2 ⋅ cos(2 ⋅ t ) + sin(2 ⋅ t ) 50 25
d) W oknie poleceń Matlaba wprowadzić następujące polecenia: >> t2 = 0:0.1:10; >> y2 = (1/5).*t2-(2/25)+exp(-t2).*((2/25).*cos(2.*t2)+ ... (47/50).*sin(2.*t2)); >> plot(t2,y2); 2. Rozwiązanie równania róŜniczkowego metodą Eulera:
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-5-
Laboratorium Modelowania i symulacji
2008 r.
Zadanie 2 Napisać algorytm rozwiązania równania róŜniczkowego (22) metodą Eulera. Wskazówka Drugą pochodną, przedstawić w postaci ilorazu róŜnicowego II rzędu. Korzystając ze znanych poleceń umieścić na jednym rysunku i odpowiednio opisać wykresy będące rozwiązaniami równania (22) metodą klasyczną i metodą Eulera. 3. Rozwiązanie równania róŜniczkowego metodą ODE: Aby wykorzystać ODE naleŜy równanie n-tego rzędu zamienić na układ n równań I rzędu. Dla równania II rzędu:
&& + 2 y + 5 y = t , y (0) = 0, y (0) = 2 & & y
dokonuje się zamiany na układ 2 równań I rzędu wykorzystując następujące podstawienia y = y1, & oraz y 2 = y1 . Wtedy powyŜsze równanie moŜna zapisać jako:
& y1 = y 2 & y 2 = −5 y1 − 2 y 2 + t
Wykorzystując edytor systemu Matlab (File | New | M-file) zapisać powyŜszy układ równań w pliku rownanie2.m w następujący sposób: function dy = rownanie2(t,y) dy=[y(2);(-5)*y(1)-2*y(2)+t]; y(n) – oznacza tu pochodną zmiennej y rzędu (n-1), y(1) – oznacza szukaną funkcję y(t). W celu rozwiązania układu równań róŜniczkowych 1-go rzędu naleŜy wpisać z linii poleceń Matlaba: >> [t_ode,y_ode]=ode45('rownanie2',[0 10],[0;2]); >> plot(t_ode, y_ode(:,1), t2, y2) >> legend('ODE', 'Rozwiązanie analityczne'); Pierwszy parametr wywołania funkcji ode45 – nazwa funkcji opisującej rozwiązywany układ równań róŜniczkowych, drugi parametr – zakres zmian zmiennej niezaleŜnej (t), trzeci parametr – kolumnowy wektor warunków początkowych. W wyniku obliczeń zwrócony zostanie: t – wektor wartości zmiennej niezaleŜnej, dla których wyznaczono rozwiązanie, y – macierz wartości funkcji y(t) oraz y’(t). Korzystając z funkcji pomocy zapoznaj się z innymi metodami rozwiązywania równań róŜniczkowych zaimplementowanymi w Matlabie w postaci funkcji. Opracowanie sprawozdania W sprawozdaniu naleŜy umieścić wykresy rozwiązań omawianych równań. Porównać rozwiązania analityczne z numerycznymi dla róŜnych parametrów (krok całkowania, metoda – Eulera, ODE). Omówić jak zmiana parametrów rozwiązań numerycznych wpływa na dokładność rozwiązania. Zadanie 3 Napisać skrypt w języku Matlab rozwiązujący równanie róŜniczkowe (11) – wstęp teoretyczny (dowolna omawiana metoda). Literatura 1. B. Mrozek, Z. Mrozek: MATLAB i Simulink: poradnik uŜytkownika. Helion, Gliwice, 2004. 2. A. Zalewski, R. Cegieła: Matlab - obliczenia numeryczne i ich zastosowania.Wydawnictwo Nakom, Poznań, 2000. 3. J. Brzózka, L. Dorobczyński: Programowanie w Matlab. Wydawnictwo Mikom,Warszawa, 1998.
Ćwiczenie 2 – Rozwiązywanie równań róŜniczkowych
-6-
