1. WST
P. METODA EULERA
1
����Ł�
1. WST
P. METODA EULERA
Wprowadzenie
Motywacja poznawania metod numerycznych:
1. Rozwiązywanie bardzo dużych konstrukcji o złożonej geometrii, wielu stopniach swobody (powyżej
1miliona), przy różnorodnym zachowaniu materiałów.
2. Świadome wykorzystywanie oprogramowania (wymaganiem jest posiadanie podstawowej wiedzy
na temat metod numerycznych).
3. Metody numeryczne pozwalają na  naukę matematyki wyższej.
Metody numeryczne pozwalają na formułowanie matematycznych problemów i rozwiązywanie ich w
efektywny sposób. Cechą charakterystyczną metod numerycznych jest zastąpienie zwykle pracochłonnych
obliczeń innymi metodami, wymagającymi zastosowania arytmetycznych działań.
Przed rozwojem przemysłu komputerowego najwięcej czasu poświęcano na rozwiązywanie
problemów (rozwiązywanie równań algebraicznych, różniczkowych), a znacznie mniej czasu poświęcano na
samo formułowanie problemu i interpretację wyników. Było to efektem trudności w otrzymywaniu
rozwiązań.
Obecnie komputery i metody numeryczne stanowią alternatywę dla metod tradycyjnych.
Wykorzystując komputery nie jesteśmy zmuszeni do robienia wielu założeń upraszczających co zwiększa
możliwość rozwiązywania złożonych problemów. W końcowym efekcie więcej czasu można poświęcić na
formułowanie problemu i interpretację wyników.
Równanie różniczkowe nazywamy zwyczajnym, gdy mamy do czynienia z pochodnymi tylko jednej
zmiennej niezależnej. Będziemy rozważać równania różniczkowe pierwszego rzędu:
dy
= f śą x , yźą
dx
lub układ równań
dy1
=� f1(x, y1K� yn )
dx
dy2
=� f2(x, y1K� yn )
dy
(1.1)
M�
dyn
=� fn(x, y1K� yn )
dx
z warunkiem początkowym:
yi (x0) =� yi,0 i =�1,K�n
(1.2)
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
1. WST
P. METODA EULERA
2
����Ł�
Zagadnienia powyższe możemy zapisać w postaci wektorowej:
dy
=� f(x,y)
dx (1.3)
y(x0 ) =� y
0
y1
�� ł�
f1 y1,0
�� ł�
�� ł�
ę�y ś�
ę� ś�
ę� ś�
2
f2 y2,0
ę� ś�
y =� ę� ś�
ę� ś�
f =� y =�
(1.4)
0
ę� M� ś�
ę� ś�
ę� M�ś� M�
ę�y ś�
ę� ś�
ę� ś�
n
�� �� fn yn,0
�� �� �� ��
W przypadku gdy n=1 zagadnienie przyjmuje postać jednego równania skalarnego:
dy
=� f(x,y)
dx (1.5)
y(x0 ) =� y0
W tym rozdziale będziemy rozważać metody całkowania równań różniczkowych z grupy metod
Rungego-Kutty, dla których wzór ogólny jest postaci:
yi+�1 =� yi +�j�h
(1.6)
h =� xi+�1 -� xi
(1.7)
gdzie h jest krokiem całkowania funkcji a Ć  kierunkiem poszukiwanego kolejnego punktu rozwiązania.
Metody Rungego  Kutty charakteryzują się tym, że poszukiwane rozwiązanie otrzymujemy za
pomocą wcześniej obliczonych wartości.
y
wielkość kroku całkowania - h
h=xśąi �ą1źą-xśąi źą
h
h
y
i �ą1
yśąi�ą1źą=yi�ąf śąxi , yiźąh
yi
Ć
xi xi�ą1
x
Rys.1.1 Interpretacja graficzna metody stycznych
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
1. WST
P. METODA EULERA
3
����Ł�
Metoda Eulera
Jedną z najprostszych metod rozwiązywania zagadnienia (1.5) jest metoda Eulera, według której
 kierunek jest wyrażony przez pierwszą pochodną funkcji:
dy
j� =� =� f (xi, yi ) (1.8)
dx
W takim razie mamy:
yi+�1 =� yi +� f (xi , yi )h
(1.9)
h =� xi+�1 -� xi
(1.10)
Tak więc pierwsza pochodna funkcji, to znaczy prawa strona równania różniczkowego, wyznacza
kierunek położenia nowego punktu rozwiązania. Jest to tzw. predyktor.
y
punkt przybliżony
y
i �ą1
predyktor
błąd kroku
y
i
punkt szukany
xi xi �ą1
x
Rys.1.2 Interpretacja graficzna metody Eulera
Odległość między predyktorem a rozwiązaniem dokładnym stanowi błąd metody. Ze względu na
swoją prostotę metoda Eulera jest łatwa do zastosowania. Dobrze oddaje charakter rozwiązania ale może
być obarczona dużym błędem.
Przykład 1.1:
Rozpatrzmy równanie różniczkowe postaci:
dy
=� -�2x3 +�12x2 -� 20x +� 8,5
dx
Znajdz
my rozwiązanie dla przedziału od x=0 do x=4 i warunku początkowego y(0)=1 i kroku
całkowania h=0,5.
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
hj
1. WST
P. METODA EULERA
4
����Ł�
Rozwiązanie dokładne:
3
��dy =� ��(-�2x +�12x2 -� 20x +� 8,5)dx
1
4 2
y =� -� x +� 4x3 -�10x +� 8,5x +� C
2
dla warunku początkowego:
y(0) =� 1
C =� 1
więc rozwiązanie dokładne ma postać:
1
y =� -� x4 +� 4x3 -�10x2 +� 8,5x +�1
2
Rozwiązanie metodą Eulera dla x = 0,5
y(0,5) =� y(0) +� f (0;1)��0,5 na podst. wzoru (1.9)
f (0;1) =� -�2(0)3 +�12(0)2 -� 20(0) +� 8,5 =� 8,5
y(0,5) =� 1+� 8,5��0,5 =� 5,25 rozwiazanie numeryczne
1
y(0,5) =� -� (0,5)4 +� 4(0,5)3 -�10(0,5)2 +� 8,5(0,5) +�1 =� 3,21875 rozw.dokladne
2
Błąd całkowity i procentowy w punkcie x = 0,5
Et =� ydokladne -� ynumeryczne
Et =� 3,21875 -� 5,25 =� -�2,03125
Et
Ek =� ��100%
ydokladne
Ek =� -�63,1%
Rozwiązanie metodą Eulera dla x = 1
y(1) =� y(0,5) +� f (0,5;5,25) �� 0,5 na podst. wzoru (1.9)
f (0,5;5,25) =� -�2(0,5)3 +� 12(0,5)2 -� 20(0,5) +� 8,5 =� 1,25
y(1) =� 5,25 +� 1,25 �� 0,5 =� 5,875 rozw.numeryczne
1
y(1) =� -� �� (1)4 +� 4 �� (1)3 -�10 �� (1)2 +� 8,5 ��1 +� 1 =� 3 roz.dokl.
2
Błąd całkowity i procentowy w punkcie x = 1
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
1. WST
P. METODA EULERA
5
����Ł�
Et =� 3 -� 5,875 =� -�2,875
Ek =� -�95,8%
Uzyskane rozwiązanie dla całego przedziału przedstawiono na (Rys 1.3) gdzie dla porównania
zamieszczono również rozwiązanie dokładne.
Metoda Eulera
7
6
5
4
3
2
1
0
0 0,5 1 1,5 2 2,5 3 3,5
Rozw. dokładne h=0,5
Rys.1.3 Rozwiązanie obliczone dla kroku h=0,5
Poprawienie dokładności można osiągnąć np. przez zmniejszenie kroku h. Na Rys1.4 przedstawiono
rozwiązania dla kroku całkowania h=0,1 i h=0,5. Widać wyraz
nie poprawienie dokładności dla mniejszego
kroku całkowania.
Metoda Eulera
7
6
5
4
3
2
1
0
0 0,5 1 1,5 2 2,5 3 3,5
Rozw. dokładne h=0,1 h=0,5
Rys.1.4 Rozwiązania dla różnych wielkości kroku
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
1. WST
P. METODA EULERA
6
����Ł�
liczba kroków
5 500
50
5000
100
10
1
0,1
1 0,1 0,01 0,001
wielkość kroku
Rys.1.5 Zależność między wielkością kroku, ilością kroków a dokładnością rozwiązania.
Oszacowanie błędu w metodzie Eulera
W analizie numerycznej wyróżniamy zasadniczo dwa typy błędów:
�� Błędy zaokragleń wynikające z wykonywania działań arytmetycznych
�� Błąd metody (dyskretyzacji, obcięcia). Określimy ten błąd dla metody Eulera.
Dowolną funkcję można przedstawić za pomocą rozwinięcia w szereg Taylora. Gdy oznaczymy
funkcję przez y mamy:
yi '' yi (n)
(1.11)
yi+�1 =� yi +� yi 'h +� �� h2 +�K� +� �� hn +� Rn
2! n!
gdzie:
yi(n +�1)
Rn =� (x� ) ��h(n+�1) (1.12)
(n +� 1)!
h =� xi+�1 -� xi
xi Ł� x� Ł� xi+�1
dy
Ponieważ nasze równanie różniczkowe ma postać =� f (x, y) zatem możemy napisać:
dx
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
Błąd względny [%]
1. WST
P. METODA EULERA
7
����Ł�
��
f '(xi , yi )�� h2
yi+�1 =� yi +� f (xi , yi ) �� h +� +�K� +�
ę�
2!
��
(1.13)
( n-�1)
ł�
f (xi , yi ) �� hn
n+�1
+� +� 0 �� (h )ś�
n!
��
Widać z równania (1.13), że błąd metody Eulera wynosi:
f '(xi , yi )�� h2
(1.14)
Et =� +�K� +� 0 �� (hn+�1)
2!
a błąd przybliżony ograniczamy do pierwszego wyrazu błędu obcięcia:
f ' (xi , yi )�� h2
(1.15)
Ea =�
2!
Jest to błąd lokalny na jednym kroku proporcjonalny do kwadratu kroku. Można wykazać, że błąd globalny
(całkowity) jest proporcjonalny do pierwszej potęgi kroku.
Przykład 1.2:
Wyznaczyć błąd metody Eulera dla wielomianu z poprzedniego zadania.
(3)
f '(xi , yi ) �� h2 f ''(xi , yi )�� h3 f (xi , yi ) �� h4
Et =� +� +� ze wzoru (1.14)
2! 3! 4!
f '(x, y) =� -�6x2 +� 24x -� 20
f ''(x, y) =� -�12x +� 24
(3)
f (x, y) =� -�12
Obliczamy poszczególne błędy metody dla 1,2 i 3 pochodnej w punkcie x= 0
-� 6(0)2 +� 24(0) -� 20
Et ,2 =� (0,5)2 =� -�2,5
2
-� 12(0) +� 24
Et ,3 =� (0,5)3 =� 0,5
6
-�12
Et ,4 =� (0,5)4 =� -�0,03125
24
Całkowity błąd równy jest sumie błędów dla poszczególnych pochodnych:
Et =� -�2,5 +� 0,5 -� 0,03125 =� -�2,03125
W ten sposób oszacowaliśmy całkowity błąd w punkcie x=0 nie znając dokładnego rozwiązania.
Zwróćmy uwagę, że błąd ten jest równy dokładnie błędowi otrzymanemu podczas rozwiązywania zadania w
przykładzie wyżej.
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
1. WST
P. METODA EULERA
8
����Ł�
Et =� ydok -� ynum =� -�2,03125
Widać, że wraz ze wzrostem rzędu pochodnej (kolejne wyrazy szeregu Taylora) wartość błędu maleje,
dlatego ograniczamy błąd do pochodnej drugiego rzędu. Dla metody Eulera błąd lokalny jest więc
proporcjonalny do kwadratu kroku a błąd globalny o jeden rząd niżej (liniowy).
W celu zmniejszenia błędu w metodzie Eulera można uwzględnić wyrazy wyższego rzędu.
Przypatrzmy się, jak to będzie w przypadku wyrazu II rzędu. I tak mamy:
2
f '(x, y) �� h
(1.16)
yi+�1 =� yi +� f (x, y)h +�
2
f ''(x, y)�� h3
(1.17)
Ea �
3!
Pierwszą pochodną musimy wyznaczyć jak dla funkcji uwikłanej. Jest ona następująca:
ś�f (x, y) ś�f (x, y) dy
f '(x, y) =� +�
(1.18)
ś�x ś�y dx
Druga pochodna będzie miała postać:
2 2
ć� ��
ś�2 f ś� f dy ś�2 f ś�2 f dy dy ś�f d y
f ''(x, y) =� +� +� �� +� �� +�
(1.19)
2 2 2
�� ��
ś�x ś�xś�y dx ś�yś�x ś�y dx dx ś�y dx
Ł� ł�
Widać, że obliczenie wyższych pochodnych wymaga dużego nakładu pracy i jest skomplikowane, dlatego
unika się takiego podejścia w praktyce.
Zastąpmy jeszcze we wzorze (1.16) wyrażenie na dpochodną wyrażeniem różnicowym:
yi+�1 -� yi
'
f �
h
Otrzymamy wówczas:
h2 yi+�1 -� yi
yi+�1 =� yi +� f (xi, yi )h +� ( )
2 h
i ostatecznie:
yi+�1 -� yi
yi+�1 =� yi +� f (xi , yi )h +� h
2
Powyższe równanie możemy zapisać również jako:
f (xi+�1, yi+�1) -� f (xi, yi )
yi+�1 =� yi +� f (xi, yi )h +� h
2
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa
1. WST
P. METODA EULERA
9
����Ł�
f (xi, yi ) +� f (xi+�1, yi+�1)
yi+�1 =� yi +� h
2
Zwróćmy uwagę, że powyższa formuła całkowania ma dokładność O(h3). Jest to tzw. zmodyfikowana
metoda Eulera. Powrócimy do tej metody w kolejnym rozdziale.
Metody Komputerowe - Dominika Mejbaum, Anna Snela, Marek Komosa