Metoda Gaussa

Jak dobrać punkty podziału [a,b] przedziału na N części

aby uzyskać dokładność 2(N+1)?





b

a

dx

)

x

(

f

)

x

(

p

I

Rozwiązanie problemu można znaleźć za pomocą

wielomianów ortogonalnych.

Ciąg wielomianów:

)

x

(

P

),...,

x

(

P

),

x

(

P

N

1

0

jest ortogonalny w przedziale [a,b] z wagą p(x) jeżeli zachodzi:















k

=

i

dla

A

k

i

dla

0

dx

)

x

(

P

)

x

(

P

)

x

(

p

))

x

(

P

),

x

(

P

(

k

k

i

b

a

)

x

(

p

k

i

Waga p(x)=1, przedział [-1,1].

Wielomiany Legendre’a:

n

2

n

n

n

n

)

1

x

(

dx

d

!

n

2

1

)

x

(

P





 

 
 





i.t.d.

1

x

3

4

1

x

P

x

x

P

1

x

P

2

2

1

0









Interpolując funkcję podcałkową f(x) za pomocą

wielomianu Legendre’a otrzymujemy:











N

0

k

k

k

1

1

)

x

(

f

A

dx

)

x

(

f

gdzie

)

x

(

P

)

x

(

P

)

2

N

(

2

A

k

'

1

N

k

2

N

k











a x

k

jest k-tym pierwiastkiem równania:

0

)

x

(

P

1

N





i błąd metody jest:

]

1

,

1

[

);

(

f

)

)!

2

N

2

)((

3

N

2

(

)

)!

1

N

((

2

e

)

2

N

2

(

3

4

3

N

2

N





















Pierwiastki i współczynniki kwadratury Gaussa z wagą p(x)=1:

N

k

1

0

0

2

2

0
1

-

0.5773502691896

26

0.5773502691896

26

1
1

3

0
1
2

-

0.7745966692414

83

0

0.7745966692414

83

5/9
8/9
5/9

4

0
1
2
3

-0.861136311594053
-0.339981043584856

0.339981043584856
0.861136311594053

0.347854845137454
0.652145154862546
0.652145154862546
0.347854845137454

5

0
1
2
3

4 

-0.906179845938664
-0.538469310105683

0

0.538469310105683
0.906179845938664

0.236926885056189
0.478628670499366
0.568888888888889
0.478628670499366
0.236926885056189

węzły x

k

współczynniki A

k

W przypadku całkowania w przedziale [a,b] dokonujemy
jego transformacji na przedział [-1,1]:

  



1

,

1

b

,

a









dt

2

a

b

dx

1

t

2

a

b

a

x













Przykład

4

dx

x

1

I

1

0

2











n

1

2

3

4

5

I

n

0.86602

54

0.79611

301

0.78725

969

0.78705

74

0.78629

544

błąd
%

10.3

1.36

0.237

0.211

0.114

Ta sama całka, ale liczona w ten sposób, że przedział całkowania

jest podzielony na 5 równych części i w każdym przedziale

zastosowano aproksymację N=1

I

1

=0.792996956 z błędem 0.97%

Całki typu:



 

  











b

a

dx

x

f

a

x

x

b

Transformacja przedziału całkowania:

  



1

,

1

b

,

a









dt

2

a

b

dx

1

t

2

a

b

a

x













Pozwala zapisać całkę w postaci:



 

  













1

1

0

dt

t

f

t

1

t

1

A

gdzie



















 



1

0

2

a

b

A

Wystarczy więc rozważyć całkę:



 

  















1

1

dt

t

f

t

1

t

1

I

Mamy obliczyć całkę funkcji f(t) z wagą



















t

1

t

1

t

w

Wielomiany ortogonalne na przedziale [-1,1] z wagą w(t)

są nazywane wielomianami Jacobiego i mają postać:





    

 















n

n

n

n

n

n

,

n

t

1

t

1

dt

d

t

1

t

1

!

n

2

1

t

P

































i całkę I można obliczyć z zależności:



 

  

 



















N

0

k

k

k

1

1

t

f

A

dt

t

f

t

1

t

1

gdzie



 

 





 

 







 





 

k

,

2

N

k

,

1

N

k

t

P

t

P

2

N

2

N

!

2

N

2

N

2

N

4

N

2

2

A









































































a t

k

jest pierwiastkiem równania:





 

0

t

P

k

,

1

N









Przypadki szczególne:

i

 

 













0

1

x

dt

t

exp

t

x





!

n

1

n







dla n=0,1,2,...





























 



n

2

1

n

2

5

3

1

2

1

n



dla n=1,2,3,...

 





 n

dla n=0,-1,-2,...

Typowe funkcje wagowe



















t

1

t

1

t

w

.

1 ==-0.5

czyli

 

2

t

1

1

t

w





a całka ma postać:

 









1

1

2

dt

t

1

t

f

I

W tym przypadku wielomiany Jacobiego są powiązane

z wielomianami Czebyszewa T

n

(t) związkiem:





 

 

t

T

C

t

P

n

n

5

.

0

,

5

.

0

n







gdzie C

n

– stała , a

 





t

arccos

n

cos

t

T

n





Łatwo więc znaleźć pierwiastki t

k

z równania:

 





0

t

arccos

N

cos

t

T

k

k

N







i mamy:



















N

2

1

k

2

cos

t

k

dla k=1,2,3,...,N

Korzystając z odpowiednich tożsamości dla funkcji Czebyszewa
i Jacobiego można wykazać, że współczynniki wagi są:





 

2

2

N

2

N

k

N

C

N

!

N

5

.

0

N

2

A









a ponieważ nie są zależne od k, więc dla danego N powinny być
stałe i równe A. Stałą A łatwo wyznaczymy jeżeli przyjmiemy, że
f(t)=1 i mamy:

NA

A

dt

t

1

N

1

k

k

1

1

2



















i stąd

N

A 



mamy ostatecznie:

 

 

f

R

N

2

1

k

2

cos

f

N

dt

t

1

t

f

N

N

1

k

1

1

2









































gdzie R

N

(f) jest błędem aproksymacji wynoszącym:

 

 

 

 

!

N

2

f

2

f

R

N

2

1

N

2

N











i –1<<1

Przykład













0

2

2

cos

k

1

d

I

gdzie |k|<1

Podstawiamy: cos=t i

2

t

1

dt

d







czyli:

 











1

1

2

2

t

1

kt

1

dt

I

w tym przypadku:

 











1

1

2

2

t

1

kt

1

dt

I

funkcja f(t) ma postać:

 

 

2

kt

1

1

t

f





i na mocy zależności:

 

 

f

R

N

2

1

k

2

cos

f

N

dt

t

1

t

f

N

N

1

k

1

1

2









































mamy:







































N

1

i

2

2

0

2

2

N

2

1

i

2

cos

k

1

1

N

cos

k

1

d

I

0

0.2

0.4

0.6

0.8

1

2

4

6

8

I k

( )

Ip k 2



(

)

Ip k 4



(

)

Ip k 6



(

)

k

0

0.2

0.4

0.6

0.8

1

20

6.67

33.33

60

eps k 2



(

)

eps k 4



(

)

eps k 6



(

)

eps k 8



(

)

k

%

100

)

k

(

I

)

n

,

k

(

Ip

)

k

(

I

)

n

,

k

(

eps





2. ==0.5

Funkcja wagowa:

 

2

t

1

t

w





 









1

1

2

dt

t

f

t

1

I

a całka ma postać:

W tym przypadku funkcje Jacobiego można wyrazić przez
funkcje Czebyszewa II – go rodzaju U

n

(t):





 









 

t

U

!

1

n

!

n

2

!

1

n

2

t

P

n

n

2

5

.

0

,

5

.

0

n







gdzie

 





 





2

n

t

1

t

arccos

1

n

sin

t

U







i wzór dla kwadratur ma postać:

 

 

f

R

1

N

k

cos

f

1

N

k

sin

1

N

dt

t

f

t

1

N

N

1

k

2

1

1

2



























































gdzie błąd określa wzór:

 

 

 

 

!

N

2

f

2

f

R

N

2

N

2

N







-1<<1

Przykład:















1

0

2

2

dx

4

x

4

x

x

x

2

I

transformacja przedziału [0,1] do przedziału [-1,1]:1

t

x





i całka przyjmuje postać:

















1

1

2

2

dt

1

t

2

t

t

1

I

Z oceny błędu widać, że wystarczy wybrać N=3, a wynik
będzie dokładny i mamy: I=1.96349541

Jeszcze jeden przypadek szczególny

3. =0.5 =-0.5

Funkcja wagowa jest:

 

t

1

t

1

t

w







i liczymy całkę:

 











1

1

dt

t

f

t

1

t

1

I

Również w tym przypadku można wyrazić wielomiany
Jacobiego przez wielomiany Czebyszewa:





 

 

 

 

 

1

t

t

T

t

T

!

n

2

!

n

2

t

P

n

1

n

2

n

2

5

.

0

,

5

.

0

n











i kwadraturę można zapisać:

 

 

f

R

1

N

2

k

2

cos

f

1

N

2

k

sin

1

N

2

4

dt

t

f

t

1

t

1

N

N

1

k

2

1

1





























































gdzie

 

 

 

 







N

2

N

2

N

f

!

N

2

2

f

R

gdzie -1<<1

Całki na przedziale nieograniczonym

Można wybrać specjalne reguły całkowania za pomocą

wielomianów ortogonalnych Laguerre’a lub Hermite’a

w przedziale [0,] lub [-,+] odpowiednio.

 

 





















dx

x

f

I

dx

x

f

I

2

a

1

Całki postaci:

 









0

x

dx

x

f

e

x

gdzie >-1

wagą jest

 

x

e

x

x

w







Wielomianami ortogonalnymi na półosi [0,) z wagą w(x)

są wielomiany Laguerre’a:

 

   





x

n

n

n

x

n

n

e

x

dx

d

e

x

1

x

L

















Dla =0 mamy:

 

 













N

1

k

k

k

0

x

x

f

A

dx

x

f

e

a x

k

– zera wielomianu Laguerre’a i współczynniki A

k

podaje

tablica:

N=
1

x

1

=1.0000000

0

A

1

=1.00000000

N=
2

x

1

=0.5857864

37627
x

2

=3.4142135

62373

A

1

=0.85355339059

3
A

2

=0.14644660940

7

N=
3

x

1

=0.2157745

56783
x

2

=2.2942803

60279
x

3

=6.2899450

82937

A

1

=0.71109300992

9
A

2

=0.27851773356

9
A

3

=0.01038925650

16

N=

4

x

1

=0.3225476

89619
x

2

=1.7457611

01158
x

3

=4.5366202

96921
x

4

=9.3950709

12301

A

1

=0.60315410434

2
A

2

=0.35741869243

8
A

3

=0.03888790851

5
A

4

=0.53929470556

1·10

-3

N=
5

x

1

=0.26356031

9718
x

2

=1.41340305

9107
x

3

=3.59642577

1041
x

4

=7.08581000

5859
x

5

=12.6408008

44276

A

1

=0.52175561058

3
A

2

=0.39866681108

3
A

3

=0.07594244968

17
A

4

=0.00361175867

992
A

5

=233699723858·

10

-4

N=
6

x

1

=0.22284660

4179
x

2

=1.18893210

1673
x

3

=2.99273632

6059
x

4

=5.77514356

9105
x

5

=9.83746741

8383
x

6

=15.9828739

80602

A

1

=0.45896467395

0
A

2

=0.41700083077

2
A

3

=0.11337338207

4
A

4

=0.01039919745

31
A

5

=0.26101720281

5·10

-3

A

6

=0.89854790643·

10

-6

Przykład:
Dana jest całka:













0

x

2

x

dx

e

1

xe

I

Dokładna wartość całki jest:

8

I

2





Obliczenia metodą trapezów przy obcięciu górnej granicy:

 











x

0

t

2

t

dt

e

1

te

x

In

Błąd obliczeń liczymy z zależności:

 

100

I

I

x

In

)

x

(

err





Wyniki przedstawia tabela:

x=1

-61.6%

x=5

-3.28%

x=10

-0.04%

x=100

-2.28·10

-

10

%

x=1000

2.26·10

-

9

%

dlaczego?

Obliczenia korzystając z wielomianu Laguerre’a

N=1 -6.26%
N=2 -0.7%
N=3 0.068%
N=4 0.048%
N=5 7.82·10

-

3

%

N=6 -2.52·10

-

3

%

Można też zastosować następujące postępowanie dla całki:

 

0

a

;

dx

x

f

I

a

1











podstawiamy:

t

dt

dx

e

t

x









i mamy:

 









a

e

0

t

1

t

dt

e

f

I

Warunkiem istnienia całki I

1

jest:

 

0

t

e

f

t

0

t

lim

















a więc stosując np. metodę trapezów trzeba w punkcie t=0

uwzględnić wartość graniczną.

Całkę

 













dx

x

f

I

2

najwygodniej rozbić na dwie:

 

 

 





























0

0

2

dx

x

f

dx

x

f

dx

x

f

I

które liczymy jak poprzednio.

Funkcja podcałkowa posiada osobliwość wewnątrz

przedziału całkowania

Np.:







2

0

x

1

dx

I

ogólnie:

 





b

a

dx

x

f

I

i w punkcie x

0

[a,b] występuje osobliwość.

1. Jeżeli istnieje

   

0

x

x

x

f

x

f

lim

0





to wydzielamy osobliwość:

 

   







  

0

b

a

0

b

a

x

f

b

a

dx

x

f

x

f

dx

x

f

I















Przykład:















2

2

2

2

dx

5

x

6

x

4

x

3

x

I

wewnątrz przedziału całkowania mianownik ma miejsce
zerowe x

0

=1, a granica funkcji w tym punkcie wynosi:

25

.

1

5

x

6

x

4

x

3

x

lim

2

2

1

x















czyli

4

25

.

1

dx

25

.

1

5

x

6

x

4

x

3

x

dx

5

x

6

x

4

x

3

x

I

2

2

2

2

2

2

2

2











































5

dx

25

.

1

5

x

6

x

4

x

3

x

I

2

2

2

2



























w obliczeniach numerycznych mamy funkcję podcałkową:

 

























1

x

dla

0

1

x

dla

25

.

1

5

x

6

x

4

x

3

x

x

f

2

2

z której całkę w przedziale [-2,2] liczymy dowolną metodą.

2. Funkcja f(x) posiada słabą osobliwość w punkcie x

0

[a,b]:

 







0

x

x

x

g

)

x

(

f

gdzie 0<<1 oraz g(x

0

)0.

i mamy:

 

   

 



























b

a

b

a

0

0

0

0

b

a

0

x

x

dx

x

g

dx

x

x

x

g

x

g

dx

x

x

x

g

I

Ostatnią całkę można obliczyć analitycznie:

 

 



























1

x

a

x

b

x

g

x

x

dx

x

g

1

0

1

0

0

b

a

0

0

Natomiast funkcja podcałkowa:

 

   









0

0

x

x

x

g

x

g

x

f

ma w punkcie x

0

granicę:

 

0

x

f

lim

0

x

x





jeżeli







0

x

x

dx

dg

Równanie różniczkowe rzędu n-go

 





0

x

,...,

x

,

x

,

t

f

n





przykład:

























t

sin

E

x

x

x

1

a

x

2

2





Układ równań różniczkowych I-go rzędu w postaci
normalnej:













n

2

1

n

n

n

2

1

2

2

n

2

1

1

1

x

,...,

x

,

x

,

t

f

x

........

..........

..........

x

,...,

x

,

x

,

t

f

x

x

,...,

x

,

x

,

t

f

x













Numeryczne metody całkowania równań różniczkowych

Sprowadzenie równania rzędu n do układu równań na
przykładzie:

























t

sin

E

x

x

x

1

a

x

2

2





oznaczamy:

 

 

2

2

2

1

1

x

x

dt

d

x

dt

d

x

x

x

x

x

x























i ostatecznie:









1

2

2

2

1

2

2

1

x

x

x

1

a

t

sin

E

x

x

x























Równanie różniczkowe rzędu n może praktycznie zawsze

być zapisane w formie układu n równań I-go rzędu

mówimy, że mamy układ równań w postaci normalnej.

Zagadnienie początkowe dla równania różniczkowego:

Przykład:

 

t

f

x

x









całka ogólna:

0

x

x













t

exp

A

x







Rozwiązanie:

Całka szczególna metodą uzmienniania stałej:



  

t

f

t

exp

A









czyli

 

 













t

0

d

exp

f

C

A

Rozwiązanie równania:

 

t

f

x

x









ma postać:





 































t

0

d

t

exp

f

t

exp

C

x

Stałą C wyznaczamy z warunków początkowych

Dla równania różniczkowego I-go rzędu potrzebujemy
jeden warunek początkowy:





0

x

0

t

x





i dla wyznaczenia stałej C mamy równanie:

C

x

0



i rozwiązanie ma postać:

 





 































t

0

0

d

t

exp

f

t

exp

x

t

x

Przykład równania różniczkowego II rzędu:

R

0

E

C

R

L

t=0

i(t)

u

c

c

u

C

i





0

u

u

RC

u

LC

c

c

c











Warunki początkowe:

 

E

R

R

R

0

u

0

c





0

R

R

E

)

0

(

i







 

E

R

R

R

0

u

0

c





0

R

R

E

)

0

(

i







Biorąc pod uwagę, że

c

u

C

i



zapisujemy warunki w postaci:

 





C

R

R

E

u

E

R

R

R

0

u

0

0

t

c

0

c















Warunków początkowych należy postawić tyle i nie więcej

ile wynosi rząd równania różniczkowego

0

u

u

RC

u

LC

c

c

c











 





C

R

R

E

u

E

R

R

R

0

u

0

0

t

c

0

c















Rozwiązanie:

Niech równanie charakterystyczne:

0

1

RC

LC

2











ma dwa różne pierwiastki rzeczywiste

2

1

,



wtedy całka ogólna ma postać:

 

t

2

t

1

c

2

1

e

A

e

A

t

u









Stałe A

1

i A

2

wyznaczamy z warunków początkowych:

 





C

R

R

E

u

E

R

R

R

0

u

0

0

t

c

0

c















czyli





C

R

R

E

A

A

E

R

R

R

A

A

0

2

2

1

1

0

2

1



















Rozwiązując powyższy układ równań wyznaczamy stałe
a następnie znajdujemy napięcie na kondensatorze.

Podsumowanie:

1. Mamy kłopoty z rozwiązywaniem równań o stałych
współczynnikach, jeżeli rząd wyższy od dwóch
równanie charakterystyczne jest wielomianem
rzędu n i nie znamy wzorów na pierwiastki.
2. Równań o zmiennych w czasie współczynnikach
i nieliniowych nie potrafimy rozwiązać w ogólnym
przypadku często już dla rzędu I-go

Konieczność użycia metod numerycznych

Ze względu na łatwość zbudowania ogólnego algorytmu

obliczenia prowadzimy dla układu równań różniczkowych

I-go rzędu w postaci normalnej

który po wprowadzeniu formalnego zapisu w postaci wektorów:

 

 

 

 

























t

x

t

x

t

x

t

X

n

2

1



 

 

 

 

























t

x

t

x

t

x

t

X

n

2

1























n

2

1

n

n

n

2

1

2

2

n

2

1

1

1

x

,...,

x

,

x

,

t

f

x

........

..........

..........

x

,...,

x

,

x

,

t

f

x

x

,...,

x

,

x

,

t

f

x





















































n

2

1

n

n

2

1

2

n

2

1

1

x

,...,

x

,

x

,

t

f

x

,...,

x

,

x

,

t

f

x

,...,

x

,

x

,

t

f

X

,

t

F



i układ zapisujemy:





X

,

t

F

X 







X

,

t

F

X 



Wniosek – Jeżeli opanujemy metodę numerycznego
rozwiązywania jednego równania różniczkowego
pierwszego rzędu, czyli

 

x

,

t

f

x 



to wyniki łatwo uogólnimy na układ n równań
I-go rzędu w postaci normalnej

Rozpoczynamy od następującego zadania:

 

x

,

t

f

x 



z warunkiem początkowym:





0

x

0

t

x





Należy znaleźć rozwiązanie dla t[0,T]

Metoda aproksymacji wielomianowej

metody wielokrokowe

 

t

,

x

f

x 



Przyjmujemy algorytm w postaci























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

t

k

=kh

Algorytm nazywamy

jawnym

, jeżeli sumowanie rozpoczyna

się od 0 w przeciwnym przypadku mówimy, że algorytm jest

niejawny

Wyznaczenie współczynników a

i

, b

i

.

Metoda jest dokładna jeżeli rozwiązaniem równania:

 

t

,

x

f

x 



jest wielomian stopnia zerowego, czyli równanie ma postać:

0

x

.



którego rozwiązaniem jest:

 

0

c

t

x









p

1

i

0

i

0

c

a

c

podstawiając do























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

mamy

ponieważ f(x,t)=0. Dzieląc przez c

0

Wyznaczenie współczynników a

i

, b

i

.

Metoda jest dokładna dla wielomianu stopnia zerowego
jeżeli rozwiązaniem równania:

 

t

,

x

f

x 



jest wielomian stopnia zerowego, czyli równanie ma postać:

0

x

.



którego rozwiązaniem jest:

 

0

c

t

x









p

1

i

0

i

0

c

a

c

podstawiając do























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

mamy

ponieważ f(x,t)=0. Dzieląc przez c

0

otrzymujemy warunek:

1

a

p

0

i

i







Jeżeli współczynniki a

i

dobierzemy tak, aby spełnić

powyższy warunek, to algorytm jest dokładny dla
wielomianów stopnia zerowego.

Żądamy, jeżeli rozwiązaniem równania

0

1

c

t

c

x





 

t

,

x

f

x 



jest wielomian pierwszego stopnia

to algorytm























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

jest dokładny.

Równanie spełniane przez wielomian

0

1

c

t

c

x





ma postać:

1

.

c

x 

czyli

c

)

t

,

x

(

f

1



a więc algorytm:























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

uwzględniając, że t

i

=ih przyjmuje postać:























p

1

i

i

i

p

0

i

0

1

i

0

1

c

b

h

]

c

h

)

i

n

(

c

[

a

c

h

)

1

n

(

c

Biorąc pod uwagę, że dla wielomianu stopnia zerowego
mamy warunek:

1

a

p

0

i

i







Z równania:























p

1

i

i

i

p

0

i

0

1

i

0

1

c

b

h

]

c

h

)

i

n

(

c

[

a

c

h

)

1

n

(

c

otrzymujemy:

1

b

ia

p

1

i

i

p

0

i

i

















Powtórzmy jeszcze jako ćwiczenie rozumowanie dla
wielomianu II-go stopnia:

0

1

2

2

c

t

c

t

c

x







który spełnia równanie

1

2

.

c

t

c

2

x





czyli

1

2

c

t

c

2

)

t

,

x

(

f





Przyjmując dla skrócenia zapisu t

n

=0 algorytm:























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

zapisujemy:

i uwzględniając:

1

2

i

n

i

n

0

1

2

2

i

n

0

1

2

2

1

n

c

)

ih

(

c

2

)

t

,

x

(

f

c

)

ih

(

c

)

ih

(

c

x

c

h

c

h

c

x

























































p

1

i

1

2

i

p

0

i

0

1

2

2

i

0

1

2

2

]

c

ih

c

2

[

b

h

]

c

ih

c

)

ih

(

c

[

a

c

h

c

h

c

a biorąc pod uwagę dwa poprzednie warunki:

1

a

p

0

i

i







1

b

ia

p

1

i

i

p

0

i

i

















otrzymujemy:













1

i

b

2

i

a

p

1

i

i

2

i

Dla wielomianu stopnia k:

0

1

1

k

1

k

k

k

c

t

c

...

t

c

t

c

x















równie różniczkowe ma postać

1

2

k

1

k

1

k

k

.

c

...

t

c

)

1

k

(

t

kc

x

















Przyjmując podobnie jak dla drugiego stopniaih

t

i



































k

1

m

1

m

m

i

n

i

n

k

0

m

m

m

i

n

k

0

m

m

m

1

n

)

ih

(

mc

)

t

,

x

(

f

)

ih

(

c

x

h

c

x

i podstawiając do:























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

otrzymujemy:

 

 

























p

1

i

k

1

m

1

m

m

i

p

0

i

k

0

m

m

m

i

k

0

m

m

m

)

ih

(

mc

b

h

)

ih

(

c

a

h

c

Porównując współczynniki przy tych samych potęgach
h i poprzednie k-1 równań znajdujemy:

1

)

i

(

b

k

)

i

(

a

p

1

i

1

k

i

p

0

i

k

i





















Algorytm wielokrokowy:























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

jest dokładny dla wielomianu stopnia k, jeżeli współczynniki
a

i

, b

i

spełniają następujący układ k równań:

1

a

p

0

i

i







1

b

ia

p

1

i

i

p

0

i

i





























1

i

b

2

i

a

p

1

i

i

2

i

...............................

1

)

i

(

b

k

)

i

(

a

p

1

i

1

k

i

p

0

i

k

i





















Każdy algorytm spełniający warunki dla k>1 nazywamy

zwartym. Jeżeli metoda jest dokładna dla wielomianu stopnia

k, to mówimy o metodzie rzędu k-go.

Liczba niewiadomych do wyznaczenia w algorytmie:























p

0

i

p

1

i

i

n

i

n

i

i

n

i

1

n

)

t

,

x

(

f

b

h

x

a

x

Współczynników a

i

wynosi: p+1

Współczynników b

i

wynosi: p+2

Całkowita liczba niewiadomych: 2p+3

1

)

i

(

b

k

)

i

(

a

p

1

i

1

k

i

p

0

i

k

i





















Układ równań:

k=0,1,...,n

można spełnić pod warunkiem, że liczba równań n+1  2p+3

W metodzie Adamsa - Bashfortha przyjmujemy p=n-1

Liczba niewiadomych wynosi: 2(n-1)+3=2n+1, więc przy

spełnieniu n+1 równań możemy dowolnie wybrać

n współczynników

Przyjmujemy n-1 współczynników a

1

=a

2

=....=a

n-1

=0

Pierwsze

z równań przy spełnieniu powyższych warunków

przyjmuje postać:

1

a

p

0

i

i







1

a

0



Jako n-ty dowolnie wybrany współczynnik przyjmujemy b

-1

=0

A więc schemat Adamsa - Bashfortha jest

schematem jawnym

Dla n=1 mamy p=0 i równanie wyznaczające b

0

1

)

i

(

b

k

)

i

(

a

p

1

i

1

k

i

p

0

i

k

i





















jest

1

b

0



Przy n=2 będzie p=1 i mamy dwie niewiadome b

0

i b

1

spełniające układ równań:

2

1

b

1

b

b

1

1

0









czyli

2

1

b

;

2

3

b

1

0



































































































































n

1

.

3

1

2

1

1

b

.

b

b

b

))

1

n

(

(

...

)

2

(

)

1

(

0

.

...

.

.

.

)

1

n

(

...

4

1

0

)

1

n

(

...

2

1

0

1

...

1

1

1

1

n

2

1

0

1

n

n

n

2

Ogólnie jeżeli mamy n niewiadomych współczynników b

i

,

to wyznaczamy je z równań

Rozwiązanie powyższego układu równań przedstawia

tabela:

Współczynniki b

k

metoda Adamsa -

Bashfortha

n

b

0

b

1

b

2

b

3

b

4

b

5

1

1

2

3/2

-1/2

3

23/12

-16/12

5/12

4

55/24

-59/24

37/24

-9/24

5

1901/72

0

-

2774/72

0

2616/72

0

-

1274/720

251/720

6

4277/14

40

-

7923/14

40

9982/14

40

-

7298/144

0

2877/14

40

-

475/144

0

Algorytm Adamsa - Bashfortha rzędu pierwszego jest

nazywany algorytmem Eulera i ma postać:





i

i

i

1

i

t

,

x

hf

x

x







Przykład:

 

t

30

sin

10

x

3

x







z warunkiem początkowym x(0)=0

Zapisujemy równanie w postaci normalnej:

 

x

3

t

30

sin

10

x







czyli w tym przypadku funkcja f(x,t) jest

 

 

x

3

t

30

sin

10

t

,

x

f





Wybór kroku całkowania h:

Wybór kroku całkowania h:

Równanie jednorodne ma postać:

0

x

3

x







Całka ogólna tego równania ma postać:





t

3

exp

A

x





Przebieg rozwiązania jest scharakteryzowany przez

wielkość tłumienia, które w tym przypadku wynosi

a=3 i charakterystyczny czas wynosi 1/a=1/3 s.

Krok czasowy h należy wybrać co najwyżej

h<1/(10a), co w tym przypadku pozwala wybrać

h=0.03s.

Należy również uważnie rozważyć funkcję wymuszającą,
która w rozpatrywanym przypadku ma postać:

 

t

30

sin

10

Okres T zmian analizowanej funkcji wynosi:

s

21

.

0

30

2

T







Rysunek funkcji sinus jest „gładki”, jeżeli podzielimy

okres na co najmniej 20 kroków, czyli krok h powinien

wynosić h<T/20, a więc h=0.01s.

t

0

h

2h

3h

4h

Z obu ograniczeń h=0.03s i h=0.01s wybieramy mniejszy
czyli w tym przypadku przyjmujemy h=0.01s.
Algorytm przy przyjętym kroku h ma postać:









n

n

1

n

x

3

n

01

.

0

30

sin

10

01

.

0

x

x













x

0

=0

i mamy:









02955

.

0

0

3

01

.

0

30

sin

10

01

.

0

x

01

.

0

t

1

1















następny:









08513

.

0

02955

.

0

3

02

.

30

sin

10

01

.

02955

.

0

x

02

.

0

t

2

2















itd..

Równanie drugiego rzędu na przykładzie wahadła
matematycznego o długości l:

0

sin

l

g

dt

d

2

2











Wprowadzamy nowe zmienne
celem zastąpienia układem
równań I-go rzędu:















sin

l

g

dt

d

dt

d

Przyjmujemy wahadło o długości l=1m i g=10m/s

2

Warunki początkowe:

 

 

0

dt

d

0

3

0

0

t

















Jeżeli zlinearyzujemy równanie:

0

sin

l

g

dt

d

2

2









to przyjmie on postać:

0

l

g

dt

d

2

2









którego całka ogólna ma postać:































t

l

g

cos

A

t

l

g

sin

A

2

1

Okres T drgań zlinearyzowanego wahadła wynosi:

s

2

g

l

2

T







czyli krok należy przyjąć

1

.

0

h

20

T

h





ponieważ mamy nieliniowe równanie, więc dla bezpieczeństwa
przyjmujemy h=0.05s.

Jawny schemat Eulera dla układu równań















sin

l

g

dt

d

dt

d

przyjmuje postać:

n

n

1

n

n

n

1

n

sin

l

g

h

h

























z warunkiem startowym:

0

3

0

0









