Metody numeryczne w10.pdf

(148 KB) Pobierz
Microsoft Word - Metody numerycznew10.doc
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Równania rózniczkowe zwyczajne
Zagadnienie początkowe
dy
=
f
(
y
(
x
),
x
),
y
n
+ 1
y
n
=
f
(
y
,
x
)
n
n
dx
h
rozwiązanie analityczne y(x)
rozwiązanie numeryczne: na przedziale [a,b]
punkty: a=x 0 , x 1 , x 2 , ....., x n =b, x i -x i-1 =h i
SCHEMAT
RÓŻNICOWY
wartości przybliżone: y 0 , y 1 , y 2 , ....., y n
ANALIZA BŁĘDU
?
wartości dokładne: y(x 0 ), y(x 1 ), .......,y(x n )
Błędy
Wybrane schematy różnicowe
W10-1
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Metoda EULERA
dy
=
f
(
y
(
x
),
x
),
y
n
+
1
y
n
f
(
y
,
x
)
y
n
y
n
+
1
=
f
(
y
,
x
)
dx
h
n
n
h
n
+
1
n
+
1
y
n
+
1
=
y
n
+
f
(
y
n
,
x
n
)
y
n
+
1
=
y
n
+
f
(
y
n
+
1
,
x
n
+
1
)
otwarta (jawna) zamknięta (niejawna)
Zmodyfikowana metoda Eulera
pół kroku w kierunku pochodnej - poprawienie pochodnej
cały krok w poprawionym kierunku
h
h
f
=
f
(
y
,
x
)
y
=
y
+
f
f
=
f
(
y
,
x
+
)
n
n
n
n
+
1
n
n
2
n
+
1
n
+
1
n
2
2
2
2
y
=
y
+
f
h
n
+
1
n
1
n
+
2
Schemat jednokrokowy
y
n
+
1
=
y
n
+
Φ
f
(
y
n
+
1
,
y
n
,
x
n
,
h
)
- niejawny (zamknięty)
y
n
+ 1
=
y
n
+
Φ
f
(
y
n
,
x
n
,
h
)
- jawny (otwarty)
W10-2
SCHEMAT
RÓŻNICOWY
115603068.007.png
 
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Schematy Rungego-Kutty r -poziomowy (etapowy) schemat R-K
=
r
Φ
f
(
y
n
,
x
n
,
h
)
=
c
i
K
i
(
y
n
,
x
n
,
h
)
i
1
r
r
K
i
(
y
,
x
,
h
)
=
f
(
y
+
h
b
ij
K
j
,
x
+
h
b
ij
),
i
=
1
2
,...,
r
j
=
1
j
=
1
otwarty (jawny):
b ij
= 0
j
i
i
1
1
K
=
f
(
y
,
x
),
K
(
y
,
x
,
h
)
=
f
(
y
+
h
b
K
,
x
+
h
b
),
i
=
2
,...,
r
1
i
ij
j
ij
j
=
1
j
=
1
K
1
=
f
(
y
n
,
x
n
),
K
=
f
(
y
+
1
K
,
x
+
1
h
),
2
n
2
1
n
2
RK4:
K
=
f
(
y
+
1
K
,
x
+
1
h
),
3
n
2
2
n
2
K
4
=
f
(
y
n
+
K
3
,
x
n
+
h
),
y
=
y
+
1
h
K
+
2
K
+
2
K
+
K
)
n
+
1
n
6
1
2
3
4
W10-3
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Metody wielokrokowe
Metody Adamsa: Adamsa-Bashfortha (jawne), Adamsa –Moultona
(niejawne)
x
n
+
1
y
(
x
n
+
1
)
=
y
(
x
n
)
+
f
(
y
(
x
),
x
)
dx
x
n
przybliżamy funkcję podcałkową wielomianem interpolacyjnym
Lagrange’a w węzłach
(
x
,
f
(
y
,
x
))
i
=
1
,
0
,
1
,...
k
metoda
niejawna
n
i
n
i
n
i
i
=
0
,
1
,...
k
metoda
jawna
całkujemy wielomian i otrzymujemy wzory postaci
k
~
k
~
y
=
α
y
+
h
β
f
(
y
,
x
)
dla metody jawnej
n
+
1
j
n
j
j
n
j
n
j
j
=
0
j
=
0
k
k
y
=
α
y
+
h
β
f
(
y
,
x
)
+
h
β
f
(
y
,
x
)
dla
n
+
1
j
n
j
j
n
j
n
j
1
n
+
1
n
+
1
j
=
0
j
=
0
metody niejawnej
W10-4
,
i
115603068.008.png
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Metoda wstecznego różniczkowania
przybliżamy rozwiązanie y(x) wielomianem interpolacyjnym W(x)
zbudowanym na węzłach
x
,
y
)
i
=
1
,
0
,
1
,...
k
metoda
niejawna
n
i
n
i
i
=
0
,
1
,...
k
metoda
jawna
obliczamy pochodną tego wielomianu W’(x) , która przybliża
pochodną rozwiązania,
z równości
y
(
x
n
)
=
f
(
y
n
,
x
n
)
W
'
(
x
n
)
otrzymujemy wzór
+
k
~
~
y
=
α +
y
h
β
f
(
y
,
x
)
dla metody jawnej
n
1
j
n
j
0
n
n
j
=
0
a z
y
(
x
n
+
1
)
=
f
(
y
n
+
1
,
x
n
+
1
)
W
'
(
x
n
+
1
)
=
k
y
α
y
+
h
β
f
(
y
,
x
)
dla metody niejawnej
n
+
1
j
n
j
1
n
+
1
n
+
1
j
=
0
W10-5
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Stosowanie metod niejawnych
k
k
y
=
α
y
+
h
β
f
(
y
,
x
)
+
h
β
f
(
y
,
x
)
n
+
1
j
n
j
j
n
j
n
j
1
n
+
1
n
+
1
j
=
0
j
=
0
y i jest
rozwiązywane metodą iteracyjną (najczęściej iteracji prostej):
n
+
1
y
[
n
i
]
=
k
α
y
+
h
k
β
f
(
y
,
x
)
+
h
β
f
(
y
[
n
+
1
1
]
,
x
)
, i=1,2,...
+
1
j
n
j
j
n
j
n
j
1
n
+
1
j
=
0
j
=
0
Wymaga ona punktu startowego, np. z metody Eulera
)
[
n
0
]
=
y
+
hf
(
y
,
x
. Uzyskanie dokładnego rozwiązania wymaga
+
1
n
n
n
wtedy wielu iteracji.
Metody typu PREDYKTOR-KOREKTOR wyznaczają przybliżenie
początkowe jawną metodą wielokrokową o takiej samej liczbie
kroków, a następnie stosują kilka iteracji rozwiązujących równanie
nieliniowe.
W10-6
(
jest nieliniowym równaniem algebraicznym względem
i
y
115603068.009.png
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Zbieżność metod jednokrokowych i sterowanie długością kroku
)
y
n
+ 1
=
y
n
+
Φ
f
(
y
n
,
x
n
,
h
y
n
+ 1
=
y
(
x
n
)
+
Φ
f
(
y
(
x
n
),
x
n
,
h
)
r
n
+ 1
=
y
(
x
n
+
h
)
[
y
(
x
n
)
+
Φ
f
(
y
(
x
n
),
x
n
,
h
)
]
- błąd lokalny
Metoda rządu p :
r
(
h
)
(
y
(
x
),
x
)
p
+
1
+
O
(
h
p
+
2
)
n
+
1
n
n
metoda
p
Eulera
1
zmodyfikowana Eulera
2
RK2,3,4
2,3,4
RKm m=5,6,..
p<m
x
n
,
y
(
x
n
)
h
y
n
+
1
x
n
,
y
(
x
n
)
h
h
y
1
1
n
+
+
2
2
2
2
h
p
+
1
y
(
x
+
h
)
y
=
ϕ
h
p
+
1
+
...
y
(
x
+
h
)
y
=
2 ϕ
+
...
n
n
+
1
n
1
1
2
n
+
+
2
2
W10-7
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
2
p
ϕ
h
p
+
1
(
y
y
)
2
p
1
n
+
1
+
1
n
+
1
2
2
16
dla p=4 : BŁĄD=
(
y
y
)
15
n
+
1
+
1
n
+
1
2
2
16
y
y
n
+
1
+
1
n
+
1
ERR < y max RELREER+ABSERR
y
:
=
2
2
n
+
1
15
ERR > y max RELREER+ABSERR zmniejszyć h
inaczej
ERR <
y
n
+
y
n
+
1
RELREER+ABSERR lub
2
ERR
ETOL =
<
1
y
+
y
n
n
+
1
RELERR
+
ABSERR
2
α
szukamy α , by
ETOL α :
(
h
)
1
ETOL
(
α ≈
h
)
ETOL
(
h
)
α
1
0.9 ,
h nowe
=
α
h
5
ETOL
(
h
)
W10-8
115603068.001.png 115603068.002.png 115603068.003.png 115603068.004.png
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Jeżeli
ETOL
(
h
)
<
1
, to
α
=
min
0
.
9
,
5
5
ETOL
(
h
)
Metoda wymaga dla RK4 4+7=11 obliczeń prawej strony.
Algorytm Rungego-Kutty-Fehlberga stosuje dwa schematy
Rungego-Kutty m+1 i m etapowy z odpowiednio dobranymi
współczynnikami. Schemat m-etapowy jest rzędu p a schemat m+1
etapowy jest rzędu p+1, a współczynniki
K są jednakowe.
n +
=
m
1
y
=
y
+
Φ
(
y
,
x
,
h
)
=
y
+
c
K
(
y
,
x
,
h
)
n
+ 1
n
f
n
n
i
i
n
n
i
1
~
n =
m
~
~
n
+ 1
=
y
n
+
Φ
f
(
y
n
,
x
n
,
h
)
=
y
+
i
K
i
(
y
n
,
x
n
,
h
)
m
+
1 =
0
i
1
= +
=
m
1
~
Wtedy
ERR
(
h
)
h
(
c
i
i
)
i
1
Algorytm GEAR’A
W10-9
Instytut Automatyki Politechniki Łódzkiej - Metody Numeryczne wykład 10
Stabilność i sztywność
Rozważmy układ n równań
dy
=
Ay
(
x
),
y
(
0
)
=
y
dx
0
n
=
rozwiązanie:
y
(
x
)
=
e
Ax
y
=
c
e
λ
i
x
dąży do 0 dla
Re λ
<
0
0
i
i
i
1
Metodą Eulera
y
n
+ 1
=
y
n
+
hAy
n
=
(
I
+
hA
)
y
n
:
y +
1
=
I
hA
)
y
0
,
y
=
(
I
+
hA
)
y
=
(
I
+
hA
)
2
y
2
1
0
y
=
(
I
+
hA
)
y
=
(
I
+
hA
)
3
y
3
2
0
.............................................
n
lim
y
n
=
0
1
+
λ
i
h
<
1
i
=
1
,...,
Obszar stabilności schematu różnicowego: zbiór wartości
zespolonych
λ , dla których wszystkie rozwiązania zadania
testowego są ograniczone dla
h
n . Jeżeli obszar stabilności
zawiera punkt 0, to metodę nazywamy stabilną.
Metoda A-stabilna A(a)-stabilna , A( α ) stabilna
W10-10
n
115603068.005.png 115603068.006.png
Zgłoś jeśli naruszono regulamin