Katedra stavební mechanikyFakulta stavební, VŠB - Technická univerzita Ostrava
Spolehlivost a bezpečnost staveb, 4.ročník bakalářského studia
Téma 2Simulační metody typu Monte Carlo
• Princip simulačních metod typu Monte Carlo• Metoda Simulation Based Reliability Assessment (SBRA)• Ukázky posudku spolehlivosti na vybraných příkladech
2 / 42Princip simulačních metod typu Monte Carlo
Historie metody Monte Carlo
Jedním z nejstarších popsaných případů využití metody je problém tzv. Buffonovy jehly, nazvaný po francouzském
matematikovi Georges-Louis Leclerc Comte de Buffonovi, který se pokoušel odhadnout hodnotu Ludolfova čísla π
náhodným vrháním jehly na linkovaný papír.
Georges-Louis Leclerc, Comte de Buffon(1707-1788)
π2
=p
Pravděpodobnost jevu, kdy jehla stejné délky, jako je vzdálenost mezi linkami, po dopadu na papír zůstane ležet na papíře tak, že protíná některou z linek, je rovna:
3 / 42Princip simulačních metod typu Monte Carlo
Historie metody Monte Carlo
Pravděpodobně první systematické využití metody Monte Carlo s reálnými výsledky je datováno až k roku 1930, kdy Nobelovou cenou oceněný italský fyzik Enrico Fermi tento přístup využíval ke generování náhodných čísel k výpočtu vlastností v té době nově objevené částice – neutronu.
Enrico Fermi(1901-1954)
Dodnes jsou tak počátky rozvoje metody Monte Carlo spojovány se jmény Stanislaw Marcin Ulam, John von Neumann,
Nicholas Metropolis nebo již zmíněný Enrico Fermi.Pojmenována podle známého kasina v Monaku (Ulamův
strýc zde sázel). Dříve se používala pod označením „statistical sampling“ – statistický výběr.
Stanislaw Ulam(1909-1984)
4 / 42Princip simulačních metod typu Monte Carlo
Historie metody Monte Carlo
Náhodnost jevů a opakování jejich výskytu jsou identické k činnostem prováděných v kasinech (ruleta je jednoduchý generátor náhodných čísel).
Metoda Monte Carlo hrála klíčovou roli při simulacích, kterými se odhadovala štěpná reakce při vývoji atomové bomby v rámci projektu Manhattan (krycí název pro utajený americký vývoj
atomové bomby za 2. světové války).
Metoda je využívána pro výpočet integrálů hustot pravděpodobnostíspojitých náhodných veličin, zejména vícerozměrných, kde běžnémetody nejsou efektivní.
Je založena na provádění náhodných experimentů s modelem systému a jejich vyhodnocení. Výsledkem provedení velkého množství experimentů je obvykle pravděpodobnost určitého jevu.
5 / 42Princip simulačních metod typu Monte Carlo
Metoda Monte Carlo
kde N je počet náhodných experimentů
Výhodou je jednoduchá implementace, nevýhodou relativně malápřesnost.
Nerr 1
=
Pro zvýšení přesnosti výsledku o jeden řád je tedy nutno zvýšit počet simulací alespoň o dva řády. Pro zisk výsledku s přesností na 6 desetinných míst, což odpovídá přesnosti jiných metod, je tedy potřeba zpracovat 1012 historií.
6 / 42Princip simulačních metod typu Monte Carlo
Zákon velkých čísel
Při velkém počtu nezávislých pokusů je možné téměř jistě očekávat, že relativní četnost se bude blížit teoretické hodnotěpravděpodobnosti.
( )NN XXN
X ++= ...11
Lze popsat s pomocí střední hodnoty náhodné veličiny:
kde X1, X2, ..., XN představuje nekonečnou posloupnost vzájemněnezávislých náhodných čísel s konečnou střední hodnotou . Se zvyšujícím se počtem historií bude střední hodnota vygenerované posloupnosti konvergovat ke střední hodnotě , což lze demonstrovat na jednoduchém příkladu s hrací kostkou.
∞<μ∞→N
μ→nX
7 / 42Princip simulačních metod typu Monte Carlo
Zákon velkých čísel
V případě hrací kostky o šesti stranách je aritmetický průměr součtu čísel na jednotlivých stranách roven:
5,3621
6654321
==+++++
=μ
Střední hodnota vržených čísel
2,0
2,5
3,0
3,5
4,0
4,5
5,0
0 2000 4000 6000 8000 10000 12000 14000 16000 18000 20000
Počet hodů
Stře
dní h
odno
ta
Vývoj vypočtenéstřední hodnoty 20000 vržených čísel
8 / 42Princip simulačních metod typu Monte Carlo
Zákon velkých čísel
Počty zastoupení jednotlivých čísel v 50000 hodech kostkou
8206 8385 8223 8383 8458 8345
0
2000
4000
6000
8000
10000
1 2 3 4 5 6
Počty zastoupení vržených čísel v 50000 hodech kostkou
9 / 42Princip simulačních metod typu Monte Carlo
Zákon velkých čísel
Procentuální zastoupení vržených čísel(celkové maximum počtu hodů 65528 je limitováno kapacitními možnostmi tabulkového procesoru Excel)
12
3
4
5
6
10100
100010000
5000065528
0,00%
5,00%
10,00%
15,00%
20,00%
25,00%
30,00%
Proc
entu
ální
zas
toup
ení
Číslo
Celkový počet hodů
Procentuální zastoupení jednotlivých čísel
11 / 42Princip simulačních metod typu Monte Carlo
Kongruenční generátory náhodných čísel
Generování pseudonáhodných čísel s pomocí rekurentního vztahu:
( ) MCUAU nn mod1 +×=+
kde konstanty A, C, a M určují statistickou kvalitu generátoru(žádoucí nesoudělnost A a M).
Nejpoužívanější generátory náhodných čísel, poprvézavedené americkým matematikem Lehmerem v roce 1948. Slouží pro generování posloupností náhodných veličin s rovnoměrným rozdělením.
Derrick Henry Lehmer(1905-1991)
12 / 42
0
1
2
3
4
5
6
7
0 4 8 12 16 20 24 28 32 36 40 44 48 52 56 60 64 68 72 76 80Princip simulačních metod typu Monte Carlo
Vliv vstupních konstant na vygenerovaná čísla
7M
3C
1A
1x0
( ) MCUAU nn mod1 +×=+
Vstupní údaje
13 / 42Princip simulačních metod typu Monte Carlo
Vliv vstupních konstant na vygenerovaná čísla
23M
3C
7A
1x0
( ) MCUAU nn mod1 +×=+
Vstupní údaje
0
5
10
15
20
25
0 4 8 12 16 20 24 28 32 36 40 44 48 52 56 60 64 68 72 76 80
14 / 42Princip simulačních metod typu Monte Carlo
Vliv vstupních konstant na vygenerovaná čísla
1011M
10C
7A
1x0
( ) MCUAU nn mod1 +×=+
Vstupní údaje
0
200
400
600
800
1000
1200
0 4 8 12 16 20 24 28 32 36 40 44 48 52 56 60 64 68 72 76 80
15 / 42
0,00
0,25
0,50
0,75
1,00
0 4 8 12 16 20 24 28 32 36 40 44 48 52 56 60 64 68 72 76 80
Princip simulačních metod typu Monte Carlo
Vliv vstupních konstant na vygenerovaná čísla
1M
0,333C
758A
0,5x0
( ) MCUAU nn mod1 +×=+
Vstupní údaje
0,00
0,25
0,50
0,75
1,00
0 4 8 12 16 20 24 28 32 36 40 44 48 52 56 60 64 68 72 76 80
Setříděno
16 / 42Princip simulačních metod typu Monte Carlo
Integrace metodou Monte Carlo
Metoda Monte Carlo se využívá nejčastěji k řešení vícerozměrných integrálů.
( ) ( )∫∫∫ ==V
y
y
x
x
yxyxfyxyxfIh
d
h
d
...dd,...,...dd,...,
Numerické integrování s využitím metody Monte Carlo spočívá ve stanovení hodnoty funkce f v N náhodných bodech, ležících v integrovanéoblasti V .
( ) ∑−
=≈N
iifN
VfVNfI1
..;
fkde představuje střední hodnotu funkce f, vypočtenou v N náhodných bodech.
17 / 42Princip simulačních metod typu Monte Carlo
Integrace metodou Monte Carlo 3D Graf vygenerovaných bodů a hodnot funkce f(x,y)
( ) 222, yxryxf −−=
Např. objem polokoule.
Ukázkový výpočet byl proveden pro 1000
vygenerovaných pseudonáhodných čísel. Poloměr polokoule r se
rovná 1 m
Výsledná hodnota odhadu integrálu
I=2,17707 (přesnáhodnota 2,09440).
Směrodatná odchylka odhadu integrálu
σ=0,04347.
18 / 42Princip simulačních metod typu Monte Carlo
Integrace metodou Monte Carlo
Náhodná proměnná 2D
0,00
0,20
0,40
0,60
0,80
1,00
1,20
0,00 0,20 0,40 0,60 0,80 1,00 1,20x
y
Tisíc vygenerovaných pseudonáhodných dvojic čísel pro výpočet objemu polokoule, zobrazených jako graf 2D náhodné proměnné
(vstupní parametry: x0=0,5; A=758; C=0,333 a M=1; y0=0,5; A=239; C=0,666 a M=1)
19 / 42
Princip simulačních metod typu Monte Carlo
Generování omezených rozdělení a transformace na požadovanérozdělení
Princip simulačních metod typu Monte Carlo
( ) MCUAU nn mod1 +×=+
20 / 42
• Např. Marek a kol.CRC Press, 1995.
• Vstupní proměnnécharakterizují useknuténeparametrické histogramy.
• Analýza funkce spolehlivosti metodou Monte Carlo.
• Spolehlivost je vyjádřena jako Pf < Pd, kde Pf je pravděpodobnost poruchy, a Pd je v normová návrhovápravděpodobnost poruchy.
Pf = Σ / Σ < Pd
Metody pro určení pravděpodobnosti poruchy
Účinek zatížení S
Odo
lnos
t R
R - S = 0
Posudek spolehlivosti metodou SBRA
21 / 42
Náhodné veličiny
Metoda Simulation Based Reliability Assessment (SBRA)
Proměnné hodnoty zatížení, variability průřezu a pevnostní charakteristiky
Dlouhodobé nahodilé hL1
Dlouhodobé nahodilé hL2Stálé hD Sníh hSn
Krátkodobé nahodilé hS Vítr hS
Napětí na mezi kluzu fy
Reprezentace náhodně proměnných veličin neparametrickými omezenými rozděleními
22 / 42
Výpočet metodou SBRA, program AntHill
Pf = 0,0000039 < Pd = 0,0000080táhlo vyhoví – úroveň spolehlivosti zvýšená
Metoda Simulation Based Reliability Assessment (SBRA)
23 / 42
Koncepty posudku spolehlivosti
Koncept „Design Pointu“ (PFD) Pravděpodobnostní alternativa
S
R
Rd > Sd
Rd
Sd
Pf = (modré)/(zelené) body
Metoda Simulation Based Reliability Assessment (SBRA)
24 / 42
Podstata metody, závěry
Metoda Simulation Based Reliability Assessment (SBRA)
• Vstupní náhodné veličiny jsou vyjádřeny useknutými histogramy (neparametrickými omezenými rozděleními).
• Pravděpodobnost poruchy Pf je získána analýzou funkce spolehlivosti RF (Reliability function, Safety function)s využitím simulační techniky Monte Carlo.
• Spolehlivost je posouzena na základě nerovnosti Pf < Pd , kde Pd je návrhová pravděpodobnost daná normou,např. ČSN 73 1401 Navrhování ocelových konstrukcí
• Výsledek se pokaždé liší, důležité zvolit dostatečný počet simulačních kroků, jejichž počet je závislý zejména na počtu náhodných veličin
• Metoda univerzální, pro složitější výpočty málo efektivní
27 / 42Ukázky posudku spolehlivosti na vybraných příkladech
Příklad 1
WIN=70kN
D=80kN
L=293,5kN
S=80kN
SN=40kN
Funkce spolehlivosti
Účinek zatížení
Odolnost konstrukce
S = NSd =80.DL + 293,5.LL ++ 80.SL + 70.WIN + 40.SN
R = NRd = Anom . Avar . fy
Statické schéma táhla
RF = ( R – S ) Ocelovýprofil IPE
28 / 42
Výpočet metodou SBRA, program AntHill
Pf = 0,0000039 < Pd = 0,0000080táhlo vyhoví – úroveň spolehlivosti zvýšená
Ukázky posudku spolehlivosti na vybraných příkladech
29 / 42Ukázky posudku spolehlivosti na vybraných příkladech
Příklad 2 – Mezní stav únosnosti
S =45kN2S =45kN1
L =75kN1
D=5kN/m'
Funkce spolehlivosti : RF = ( R – S )
S ... účinek zatížení (ohybový moment)R ... odolnost konstrukceS = D.DL.L2/8 + L1.LL.L/4 + S1.SL.L/3R = Wnom . Wvar . fy
30 / 42Ukázky posudku spolehlivosti na vybraných příkladech
Příklad 3 – Mezní stav použitelnosti
S =45kN2S =45kN1
L =75kN1
D=5kN/m'
S ... účinek zatížení (průhyb)R ... odolnost konstrukceS = (5.D.DL.L4/384 + L1.LL. L3/48 +
+ 0.0355.S1.SL.L3)/(210000. Inom.Ivar)
R = L / 350
Funkce spolehlivosti : RF = ( R – S )
31 / 42Ukázky posudku spolehlivosti na vybraných příkladech
Posouzení spolehlivosti konstrukce
AN
±=σ
AfNN yRdSd .=≤
32 / 42Ukázky posudku spolehlivosti na vybraných příkladech
Posouzení spolehlivosti konstrukce
⎟⎟⎠
⎞⎜⎜⎝
⎛±±±=
z
z
y
y
WM
WM
ANσ
35 / 42Ukázky posudku spolehlivosti na vybraných příkladech
Posouzení spolehlivosti konstrukce
var,var,var zz
z
yy
y
hWWM
hWWM
hAANS
∗±
∗±
∗=
N = - ( 1,35.D.hD + 1,5.( L1.hL1 + SN.hSn + WINx.hW + Sx.hS + L2x.hL2 )) [kN]
My = b . 1,5 . ( Sz.hS - L2z.hL2 + WINz.hW ) [kNm]
Mz = b . 1,5 . ( Sy.hS - L2y.hL2 ) [kNm]
Funkce spolehlivosti : RF = ( R - S )
S ... účinek zatížení
kde
36 / 42
S(t)
R(t)
čas t
Závislost R a S na čase t
Ukázky posudku spolehlivosti na vybraných příkladech
Posouzení životnosti konstrukce
38 / 42
RF = {MR(t) – MS(t)} kdeMS (t) = = 1/8×(DL×DLvar)×L2 + 1/4×(LL×LLvar)×L + 1/3×(SL×SLvar)×L
a MR(t) = (As×As,var) × fy × z [kN.m]:
Ukázky posudku spolehlivosti na vybraných příkladech
Křivky trvání účinků zatížení
39 / 42
MS(t) [kN.m] 0 až 50 let
Ukázky posudku spolehlivosti na vybraných příkladech
Účinek zatížení
40 / 42Ukázky posudku spolehlivosti na vybraných příkladech
Odolnost konstrukce
MR(t) [kN.m] 0 až 50 let
41 / 42
0,0
Bezpečné
Porucha
Ukázky posudku spolehlivosti na vybraných příkladech
MR(t) [kN.m] 0 až 50 let
Funkce spolehlivosti