Metody Numeryczne Wykład 9

Metody Numeryczne Wykład 9

Kwadratury Na pocz¸ atku zdefiniujemy poj¸ecie kwadratury K . Kwadratur¸ a K opart¸ a na w¸ezłach x _j nazywamy

K =

n

X

j=0

a _j f (x _j ) ≈ Z b

a

f (x)dx

Współczynniki a j nazywamy wagami kwadratury.

Podstawowe kwadratury tzw. kwadratury Newtona-Cotesa oparte s¸ a na w¸ezłach wielo- mianu interpolacyjnego Lagrange

p _n (x) =

n

X

j=0

f (x _j )L _j (x)

gdzie funkcje ( porównaj wykłady 4-5 )

L _j (x) = (x − x ₀ )(x − x ₁ ) . . . (x − x _j−1 )(x − x _j+1 ) . . . (x − x _n ) (x _j − x ₀ )(x _j − x ₁ ) . . . (x _j − x _j−1 )(x _j − x _j+1 ) . . . (x _j − x _n ) . tworz¸ a tzw. baz¸ e Lagrange

L _j (x _i ) =  0 i 6= j, 1 i = j.

St¸ ad

Z b a

f (x)dx ≈ Z b

a

p _n (x)dx = Z b

a n

X

j=0

f (x _j )L _j (x)dx

=

n

X

j=0

f (x j ) Z b

a

L j (x)dx.

Przykład 1

Kład¸ ac n = 1 i interpoluj¸ ac funkcj¸e f w w¸ezłach x ₀ = a, x ₁ = b, otrzymujemy L ₀ (x) = x − b

a − b , L ₁ (x) = x − a

b − a .

a ₀ = Z b

a

x − b

a − b dx = b − a 2 , a ₁ =

Z b a

x − a

b − a dx = b − a 2 .

prost¸ a kwadratur¸ e trapezów:

Z b a

f (x)dx ≈ b − a

2 [f (a) + f (b)] .

Przykład 2

Kład¸ ac n = 2, x 0 = a, x 1 = ^b+a ₂ , x 2 = b, otrzymujemy prost¸ a kwadratur¸ e Simpsona (parabol)

Z b a

f (x)dx ≈ b − a 6



f (a) + 4f ( b + a

2 ) + f (b)

 .

Kwadratury złożone

Napiszemy w OCTAVE skrypt funkcyjny o nazwie trapez.m obliczaj¸ acy całki danych funkcji na przykład f (x) = e ^x , g(x) = sin2πx, h(x) = x ^1/3 , na przedziale [a, b] = [0, 1]

złożon¸ a kwadratur¸ a trapezów Z x

x

f (x)dx ≈ h (f (x _i−1 ) + f (x _i )

2 ,

gdzie

h = (b − a)/N jest długości¸ a kroku ,

N - ilości¸ a równych podziałów przedziału całkowania a b.

Sumuj¸ ac pola trapezów w każdym w przedziale [x i−1 , x _i ] i wykorzystuj¸ ac własność ad- dytywności całki oznaczonej - wzór na złożon¸ a kwadratur¸e trapezów można zapisać w postaci

Z b a

f (x)dx ≈ h f (a) + f (b)

2 +

N −1

X

j=1

f (a + jh)

!

Sumuj¸ ac zaś proste kwadratury Simpsona w każdym z przedziałów [x i , x _i+1 ], otrzy- mujemy złożon¸ a kwadratur¸ e Simpsona

Z b a

f (x)dx ≈ h

3 (f ₀ + 4f ₁ + 2f ₂ + 4f ₃ + 2f ₄ + ... + 2f _n−2 + 4f _{N −1} + f _N ) gdzie

f i = f (x i ) = f (a + ih), h = ^b−a _N .

Napiszemy prosty skypt funkcyjny w OCTAVE o nazwie it simpson.m obliczaj¸ acy całk¸e oznaczon¸ a złożon¸ a kwadratur¸ a Simpsona.

(Program simpson.m zał¸acznik 2 ) Kwadratury Gaussa

Kwadraturami Gaussa nazywamy kwadratury postaci Z 1

−1

f (t)dt =

n

X

i=1

w i f (t i ) gdzie

wagi kwadratury w _i i w¸ezły t _i s¸ a tak dobrane,że kwadratura jest dokładna dla wszystkich wielomianów f (t) stopnia 2n-1 i istnieje wielomian stopnia 2n+2 dla którego kwadratura nie jest dokładna.

Przykład

Jeśli f(t) stanowi baz¸e pot¸egow¸ a {1, t, t ² , t ³ } ( patrz wykłady 4-5 ), to R 1

−1 1dt = 1 − (−1) = 2 = w ₁ + w ₂ ; R 1

−1 tdt = 1/2 − 1/2 = 0 = w ₁ t ₁ + w ₂ t ₂ ; R 1

−1 t ² dt = 1/3 − (−1/3) = 2/3 = w ₁ t ² ₁ + w ₂ t ² ₂ ; R 1

−1 t ³ dt = 1/4 − 1/4 = 0 = w ₁ t ³ ₁ + w ₂ t ³ ₂ ;

Rozwi¸ azaniem układu s¸ a w 1 = w 2 = 1 i t 1 = − ^{√ 1} ₃ , t 2 = ^{√ 1} ₃ . St¸ ad otrzymujemy dwupunktow¸ a kwadratur¸e Gaussa

Z 1

−1

f (t)dt = f



− 1

√ 3

 + f

 1

√ 3

 .

Do obliczenia całki w dowolnym przedziale [a, b] - dokonujemy liniowego ptrzekształcenia

x = mt + c przedziału [a, b] na przedział [-1, 1]

x =  b − a 2



t +  b + a 2

 . Wówczas

Z b a

f (x)dx = b − a 2

Z 1

−1

f  b − a

2 t + b + a 2

 dt

Napiszemy w OCTAVE skrypt o nazwie gaussint.m obliczaj¸ acy dowoln¸ a całk¸e ozna- czon¸ a kwadratur¸ a Gaussa.

Zał¸ acznik 3 do wykładu.

Program OCTAVE posiada instrukcj¸e wewn¸etrzn¸ a o nazwie quad.m - numerycznego cał- kowania funkcji jednej zmiennej.

Sposób wywołania tej instrukcji jest nast¸epuj¸ acy:

[v, ier, nf un, err] = quad(f, a, b, tol, sing)

gdzie

f - jest nazw¸ a funkcji najcześciej zapisywana w odzielnym m-file;

a, b - granic¸ a doln¸ a i górn¸ a całkowania;

tol - opcjonalny argument dokładności obliczania całki - jest to wektor dwuelementowy, którego pierwszy element oznacza bład bezwzgl¸edny , drugi - bł¸ ad wzgł¸edny obliczeń sing - opcjonalny argument - wektor przedziału, w którym obliczana cal ka jest zbieżna (singular) ( gł’ownie ma zastosowanie dla dla całek niewłaściwych;

Otrzymujemy po wywołaniu instrukcji ENTER v - wartość całki;

ier - wskaźnik bł¸edu obliczeń, na przykład liczba 0 - obliczenie całki dokładne ; nfun - wskaźnik ilości obliczeń numerycznych funkcji, do otrzymania wyniku;

err - bł¸ ad obliczenia całki Na przykład

[v, ier, nf un, err] = quad( ⁰ exp ⁰ , −1, 1, [1.e − 3 1.e − 3], [1 0]) + EN T ER v = 2.3504

ier = 0 nfun = 63 err = 2.6095e-14

»

Na tym kończymy przegl¸ ad najważniejszych kwadratur numerycznych, choć autor wy-

kładu zdaje sobie spraw¸e, że nie przedstawił takich zagadnień jak dokładność i zbieżność

kwadratur numerycznych a także kwadratur Romberga. Wyszedł bowiem z założenia, że

w praktyce inżynierskiej przedstawione powyżej kwadratury s¸ a wystarczaj¸ ace.