Ilustracja metody
MONTE CARLO
obliczania
całek podwójnych
Często jest tak, iż wiemy, że istnieje całka oznaczona z funkcji f
jednak nie potrafimy jej analitycznie policzyć. Konieczne jest
wtedy zastosowanie jakiejÅ› metody numerycznej.
Metody całkowania oparte o użycie liczb losowych pojawiły się
w latach czterdziestych XX wieku. W Los Alamos w ramach
"Projektu Manhattan" (nad budowÄ… bomby jÄ…drowej pracowali
wtedy m.in. John von Neumann, Stanisław Ulam, Richard
Feynman) potrzebne było obliczenie całek dotyczących
rozpraszania i absorpcji neutronów. Wykorzystana została
do tego metoda Monte Carlo.
Nazwa Monte Carlo, kasyna w Monaco, nawiÄ…zuje
do losowości gier hazardowych.
Projekt prezentuje zastosowanie metody Monte Carlo do
obliczania całek podwójnych. Zgodnie z Mocnym Prawem
Wielkich Liczb Kołmogorowa, całka dana jest wzorem:
b
n
b - a
n
f (X ) çÅ‚ çÅ‚çÅ‚"
"
k
+"f (x)dx
n
k = 1
a
Przybliżenie całki dla całek podwójnych wygląda następująco:
(b - a)2 n
f (xk , yk )
"
+"+"f (x, y)dxdy H"
n
k = 1
D
Gdzie xk i yk to odpowiednie losowe liczby
yk = a + (b - a)vk
xk = a + (b - a)uk
Przykłady
całek
x3 + y3 dxdy
+"+"
D
x3 + y3 dxdy
+"+"
D
D = {(x, y) : 2 e" x e" - 2,0 d" y d" 4 - x2 }
x3 + y3 dxdy
+"+"
D
Obliczenia analityczne
4- x2
2 4- x2 2
1
ëÅ‚ öÅ‚
(x3 + y3)dy = dx =
+"dx +" +"ìÅ‚ 4 y4 ÷Å‚
íÅ‚ Å‚Å‚
- 2 0 - 2
0
2 2
1 1
+"(4 - x2)2 dx = 4 +"(16 - 8x2 + x4)dx =
4
- 2 - 2
2
1 8 1 128
ëÅ‚
16x - x3 + x5 öÅ‚ =
ìÅ‚ ÷Å‚
4 3 5 15
íÅ‚ Å‚Å‚
- 2
Metoda Monte Carlo
D'= [- 2;2]× [- 2;2]
~obszar D zawarty jest w kwadracie
~n liczb ui i vi losujemy z rozkładu jednostajnego U (0;1)
I przekształcamy według wzorów.
~wykorzystując wylosowane liczby możemy przybliżyć całkę:
42 n
3
( (xk 3 + yk 3)
"
+"+"x + y3)dxdy H" n
k = 1
D
Stworzony został wykres przedstawiający rozkład wyników dla wykonanych 50 prób dla każdego n,
wyskalowany został półlogarytmicznie, pominięto ukazanie rozrzutu wartości dla n=1 i n=10 ze
względu na ich duże różnice
Wyniki prób Monte Carlo Wartość całki
10
9
y = 8,53333
8
7
6
10 100 1000 10000 100000 1000000
Ilosc prób Monte Carlo n
Wartość obliczona metodami standardowymi to 8,5(3) . Metodą Monte Carlo, dla 50 prób
przy n=1000000 uzyskujemy wartość średnią 8,535893, co stanowi bardzo dobre
przybliżenie rzeczywistej wartości całki.
Wartosci uzyskane
sin( y)
dxdy
+"+"
y
D
sin(y)
dxdy
+"+"
y
D
D = {(x, y) : 0 d" x d" Ä„ , x d" y d" Ä„ }
sin(y)
dxdy
+"+"
y
D
Obliczenia analityczne
x= y
y
Ä„ Ä„
ëÅ‚ öÅ‚
sin(y) sin(y) sin(y)
dxdy = dx = dy =
+"+" +"dy+" +"ìÅ‚ x y ÷Å‚
ìÅ‚ ÷Å‚
y y
íÅ‚ Å‚Å‚
D 0 0 0
x= 0
Ä„
=
+"sin( y)dy = - [cos(Ä„ ) - cos(0)] = 2
0
Metoda Monte Carlo
D'= [0;Ä„ ]× [0;Ä„ ]
~obszar D zawarty jest w kwadracie
~n liczb ui i vi losujemy z rozkładu jednostajnego
U (0;1)
I przekształcamy według wzorów na xi i yi
~wykorzystując wylosowane liczby możemy przybliżyć całkę:
2
n
sin(y) 4Ä„ sin( yk )
dxdy H"
"
+"+"
y n yk
k = 1
D
Wyniki prób Monte Carlo Wartość całki
4
3,5
3
2,5
y = 2
2
1,5
1
10 100 1000 10000 100000 1000000
Ilosc prób Monte Carlo n
Wykres przedstawiający rozkład wyników dla wykonanych
100 prób dla każdego n. Wartość obliczona metodami
standardowymi to 2 . Metodą Monte Carlo, dla 100 prób przy
n=1000000 uzyskujemy wartość średnią 1,99978.
Wartosci uzyskane
2
(
+"+"x + y2)dxdy
D
2
(
+"+"x + y2)dxdy
D
D = {(x, y) : y e" 0, y d" x2 + y2 d" x}
2
(
+"+"x + y2)dxdy
D
Obliczenia analityczne
Aby obliczyć całkę analitycznie,
konieczne jest przejście na współrzędne biegunowe.
Obszar, po którym całkujemy, jest fragmentem koła
y d" x2 + y2 d" x
2
Ô! Á sinÕ d" Á d" Á cosÕ
Ô! sinÕ d" Á d" cosÕ
we współrzędnych biegunowych odpowiada on:
Ä„
Å„Å‚ üÅ‚
" = (Á ,Õ ) : 0 d" Õ d" ,sinÕ d" Á d" cosÕ
òÅ‚ żł
4
ół þÅ‚
Ä„ Ä„
Á = cosÕ
cosÕ
4 4
1 öÅ‚
ëÅ‚
2 2 2 4
( Á Á dÁ dÕ =
÷Å‚
+"+"x + y2)dxdy = +"+" +"dÕ +"Á Á dÁ = +"dÕ ìÅ‚ 4 Á ÷Å‚
íÅ‚
Å‚Å‚
D " 0 sinÕ 0
Á = sinÕ
Ä„
Ä„
Õ =
4
4
1 1 1
= (cos4 Õ - sin4 Õ )dÕ = (cosÕ sinÕ ) =
+"
4 4 8
Õ = 0
0
w obliczenia wykorzystane zostaly calki :
1
4
+"cos Õ dÕ = 32 (12Õ + 8sin(2Õ ) + sin(4Õ ))
1
4
+"sin Õ dÕ = 32 (12Õ - 8sin(2Õ ) + sin(4Õ ))
1
(cos4 Õ - sin4 Õ )Õ dÕ = *16sin(2Õ ) = sinÕ cosÕ
+"
32
Metoda Monte Carlo
~obszar D zawarty jest w kwadracie
D'= [0;1]× [0;1]
~n liczb xi i yi losujemy bezpośrednio
z rozkładu jednostajnego
U (0;1)
~wykorzystując wylosowane liczby możemy przybliżyć całkę:
n
1
2
( (xk 2 + yk 2)
"
+"+"x + y2)dxdy H" n
k = 1
D
Wyniki prób Monte Carlo Wartość całki
0,25
0,2
0,15
y = 0,125
0,1
0,05
0
10 100 1000 10000 100000 1000000
Ilosc prób Monte Carlo n
Wartosci uzyskane
e- ( x2 + y2 )dxdy
+"+"
D
e- ( x2 + y2 )dxdy
+"+"
D
D = {(x, y) : x e" 0, y d" x, x2 + y2 d" ln3}
e- ( x2 + y2 )dxdy
+"+"
D
Obliczenia analityczne
e- ( x2 + y2 )dxdy
D = {(x, y) : x e" 0, y d" x, x2 + y2 d" ln3}
+"+"
D
D we współrzędnych biegunowych odpowiada
Ä„ Ä„
" = {( f , p) : - d" f d" ,0 d" p d" ln 3}
2 4
Ä„
ln 3
4
- (x2 + y2 ) - p2 - p2
dxdy = df
+"+"e +"+"e dfdp = +" +"pe dp = *
Ä„
D " 0
-
2
Ä„ Ä„
îÅ‚ - p2
Å‚Å‚
z =
Ä„
4 4
ln 3
1 1 Ä„
4
ïÅ‚ śł
- p2
= (- e- p2 0 df = f = * z'= - 2 p =
+"pe dp = ïÅ‚ śł
+" +"
Ä„
-
2 3 4
2
Ä„ Ä„
ïÅ‚
dz = - 2 pdpśł
- -
ðÅ‚ ûÅ‚
2 2
1 1
z p2
-
+"e dz = - 2 e-
2
Metoda Monte Carlo
~obszar D zawarty jest w kwadracie D'= [- ln 3;ln 3]× [- ln 3;ln 3]
U (0;1)
~n liczb ui i vi losujemy z rozkładu jednostajnego
i przekształcamy według wzorów.
~wykorzystując wylosowane liczby możemy przybliżyć całkę:
n
2 2
4ln2 3
e- ( x2 + y2 )dxdy H" e- ( xk + yk )
"
+"+"
n
k = 1
D
Wyniki prób Monte Carlo Przewidywana wartość całki
1,75
1,25
y = Ä„/4
0,75
0,25
1 10 100 1000 10000 100000 1000000
Ilosc prób Monte Carlo n
Wykres przedstawiający rozkład wyników dla wykonanych 80
prób dla każdego n
Wartosci uzyskane
10xxsin( y)dxdy
+"+"
D
10xxsin( y)dxdy
+"+"
D
D = {(x, y) : 0 d" x d" 1; 0 d" y d" x}
10xxsin( y)dxdy
+"+"
D
Obliczenia analityczne
Całki tej nie da się wyrazić za pomocą funkcji
elementarnych. Tutaj z pomocÄ… przychodzi metoda Monte
Carlo, dzięki której możemy poznać jej przybliżoną wartość.
Metoda Monte Carlo
~obszar D zawarty jest w kwadracie
D'= [0;1]× [0;1]
~n liczb xi i yi losujemy bezpośrednio
z rozkładu jednostajnego
U (0;1)
~wykorzystując wylosowane liczby możemy przybliżyć całkę:
n
10
10xxsin( y)dxdy H" xk xk sin( yk )
"
+"+"
n
k = 1
D
Wyniki prób Monte Carlo Przewidywana wartość całki
4
3
2
y = 1.31947
1
0
1 10 100 1000 10000 100000 1000000
Ilosc prób Monte Carlo n
Wartosci uzyskane
Autorzy projektu:
" Dorota Moskal 171695
" Mateusz Sawicki 171620
" Mariusz Orda 171665
" Michał Piórek 171677
" Mateusz Fuławka 171623
" Miłosz Karolonek 171595