XLVIII KONFERENCJA NAUKOWA
KOMITETU INŻYNIERII L
DOWEJ I WODNEJ PAN
I KOMITETU NAUKI PZITB
Opole  Krynica 2002
Joanna KALISZUK1
Zenon WASZCZYSZYN2
PROBABILISTYCZNA OCENA NIEZAWODNOÅšCI KONSTRUKCJI
METODAMI MONTE CARLO Z WYKORZYSTANIEM SSN
1. Wstęp
Ocena niezawodności konstrukcji, w odniesieniu do jej różnych stanów czy też różnych
czynników, jest złożonym zagadnieniem, które wymaga wyjścia poza sformułowania
deterministyczne. Losowość zmiennych i parametrów analizowanych modeli implikuje
posługiwanie się probabilistycznymi metodami analizy (pomijamy inne metody, w tym
metody półprobabilistyczne, por. [1]). Metody probabilistyczne, stosowane do analizy
konstrukcji inżynierskich są omawiane w bogatej literaturze, że wymienimy tylko dwie
polskie monografie [1, 2].
Złożoność analizy problemów niezawodności konstrukcji powoduje, że można je
analizować głównie numerycznie, stosując różne metody symulacyjne, a wśród nich różne
modyfikacje metod Monte Carlo [3]. Zostały one oparte na solidnych podstawach matema-
tycznych, por. [3, 4] i stale sÄ… wzbogacane nowymi propozycjami symulacji komputerowych
[5]. Najprostsza metoda klasyczna MC polega na generowaniu rozwiązań próbnych (dalej
nazywamy je krótko "próbkami"), i sprawdzaniu czy spełniają one przyjęte kryteria nieza-
wodności. Bardziej złożoną jest metoda średniej ważonej (ang. importance sampling), która
wymaga mniejszej liczby próbek do obliczenia prawdopodobieństwa niezawodności [3, 5].
W teorii konstrukcji próbki dla metod MC liczymy zazwyczaj metodą elementów
skończonych (MES), por. [6]. Wymienione symulacje numeryczne wymagają na ogół bardzo
dużej liczby próbek co przy bardziej złożonych konstrukcjach czyni metody MC
numerycznie nieefektywnymi. W tym zakresie interesującą jest próba zastosowania
sztucznych sieci neuronowych do generowania próbek [7].
W referacie rozwijamy pomysł z pracy [7] i sprawdzamy go na przykładzie oceny
probabilistycznej nośności płaskiej ramy sprężysto-plastycznej z losową granicą
plastyczności materiału. Do obliczenia wzorców do uczenia i testowania sieci stosujemy
własny program MES, oparty na pracy [8]. Dzięki prostej sieci jednokierunkowej istotnie
skrócono czas generowania próbek co umożliwiło szybkie wyznaczanie probabilistycznej
granicy nośności dla różnych przypadków zmienności granicy plastyczności w prętach
rozpatrywanej ramy.
1
Mgr inż., Wydział Inżynierii Lądowej i Środowiska Uniwersytetu Zielonogórskiego
2
Prof. zw. dr hab. inż., Wydział Inżynierii Lądowej Politechniki Krakowskiej
78
2. Niezawodność i zawodność konstrukcji
W dalszym ciągu pomijamy czynnik czasu. Za miarę niezawodności konstrukcji przyjmuje
się prawdopodobieństwo dla ustalonego czasu, zdefiniowane w następujący sposób:
Q = Prob {G (R, S ) > 0} a" Prob {R > S } , (1)
gdzie: G (S, R )  funkcja stanów (zapas bezpieczeństwa), R = R (XR)  losowa nośność
konstrukcji, S (XS)  losowe obciążenie, X = [XR, XS]  wektor losowych stanów kon-
strukcji, por. [9].
Do dalszych rozważań wprowadzamy też losową zawodność (niesprawność) konstrukcji:
G (X) a" R - S d" 0 , (2)
oraz prawdopodobieństwo zawodności :
q = Prob { G (X) d" 0 } = +" f (X ) dX , (3)
G(X) d" 0
gdzie: f (X)  funkcja gęstości prawdopodobieństwa zawodności.
3. Numeryczne symulacje metodami Monte Carlo
Symulacje metodami Monte Carlo służą numerycznemu obliczaniu wartości całki (3).
Zakłada się, że wszystkie zmienne losowe funkcji G (X) są zmiennymi niezależnymi. Przy
takim założeniu jest słuszne prawo wielkich liczb, a więc przybliżoną wartość prawdopodo-
bieństwa zawodności konstrukcji można obliczać wzorem [4]:
N
1
q = I(Xi ) (4)
"
N
i=1
gdzie: N  całkowita liczba symulacji (wylosowanych próbek), I (Xi) następująco zdefinio-
wany wskaz
nik:
1 dla G(Xi )d" 0 ,
Å„Å‚
I(Xi ) = (5)
òÅ‚
)
ół0 dla G(Xi > 0 .
Powyższa metoda, zwana metodą klasyczną MC (KMC lub metoda  orzeł-reszka ,
por.[4]), jest łatwa w stosowaniu, jednak dla uzyskania oszacowań obarczonych małym
błędem, wymaga bardzo dużej liczby symulacji. Przybliżony błąd tej metody określony jest
wzorem (6), por. [4]:
´ = 1 N , (6)
gdzie: N  całkowita liczba symulacji.
Istnieje wiele metod o znacznie lepszej efektywności [3, 4]. Jedną z nich jest metoda
 średniej ważonej (w skrócie WMC, zwana też metodą losowania istotnego, por. [4]). W
metodzie tej, łączne prawdopodobieństwo zawodności konstrukcji przybliża wzór:
79
N
1
q = (Xi) , (7)
"K
N
i=1
gdzie K (Xi)  zdefiniowano poniżej:
Å„Å‚ f (Xi ) g (Xi ) dla G(Xi )d" 0 ,
X X
K(Xi )= (8)
òÅ‚
0 dla G(Xi )> 0.
ół
gdzie: fX (X)  funkcja gęstości rozkładu prawdopodobieństwa zmiennej losowej, gX (X) 
przyjęta funkcja wagi.
4. Zastosowanie SSN do generowania próbek dla metod MC
Do generowania próbek okreś lających stany konstrukcji zazwyczaj stosuje się metodę
elementów skończonych. Stosowanie programu MES do obliczania próbek do metod MC
jest jednak na ogół nieefektywne numerycznie. Wynika to z konieczności obliczenia bardzo
wielkiej liczby próbek (rzędu kilkunastu a nawet kilkuset tysięcy), co może wymagać bardzo
długich czasów obliczeń. Z tego powodu w pracy [7] zaproponowano zastosowanie SSN. W
takim podejściu program MES jest używany do obliczenia wzorców do uczenia i testowania
SSN. Po nauczeniu sieci jest ona stosowana do generowania próbek dla metod MC. W [7]
wykazano, że takie podejście daje skrócenie czasu rzędu 10-3 czasu potrzebnego do
generowania próbek za pomocą programu MES w porównaniu z generowaniem próbek
sieciÄ… neuronowÄ….
Wyżej opisane podejście stosujemy do obliczania próbek odpowiadającej nośności
płaskich ram sprężysto-plastycznych, którą można obliczyć jednym z wielu programów
MES. W naszym referacie zastosowano program ELPLAF-v1,opracowany na podstawie [8]
Program uwzględnia rozłożoną plastyczność i geometryczne nieliniowości rzędu II. Do
obliczeń przyjęto model materiału idealnie plastycznego. Nośność graniczna odpowiada
osiągnięciu granicznego punktu na ścieżce równowagi, a więc nośność jest obciążeniem
krytycznym ogólnej utraty stateczności ramy.
5. Probabilistyczna ocena nośności płaskiej ramy stalowej
Celem sprawdzenia działania opisanej metody oceniania bezpieczeństwa ram stalowych
przeprowadzono test numeryczny. Do analizy przyjęto wielokondygnacyjną ramę
analizowaną w [6]. Schemat ramy i konfigurację jednoparametrowego obciążenia pokazano
na rys.1. Przyjęto, że rama jest złożona z jednakowych prętów o geometrycznych
chakterystykach przekroju poprzecznego i module sprężystoś ci jak na rysunku. Przyjęto, że
wartość Å›rednia granicy plastycznoÅ›ci wynosi µ = 240 Mpa, a odchylenie standardowe
à = 24MPa.
Obliczenia wykonane programem ELPLAF-v1 dla granicy platycznoÅ›ci µ , jednakowej
dla wszystkich prętów dały wartość nośności ramy Pk = 1042.3 kN, a więc niższą niż
nośność o wartości 1143.7 kN, obliczona w [6] metodą przegubów plastycznych na
podstawie teorii geometrycznie liniowej.
Dalej przyjęto, że granice plastyczności dla wszystkich słupów i rygli są dwoma
zmiennymi losowymi, odpowiednio X1 = R s i X 2 = R r . Wszystkie obliczenia, zarówno
dotyczące tworzenia sieci neuronowej jak i generowania próbek metodami MC wykonano w
środowisku MATLAB 5, [10].
80
P P P
0.1P
E = 21000 kN/cm2
0.1P
Ã0 = 24 kN/cm2
A = 149 cm2
0.1P
J = 25170 cm4
S = 1868 cm3
0.1P
Z = 1638 cm3
1200
1200
Rys. 1. Schemat ramy i dane geometryczne i materiałowe
5.1. Neuronowa predykcja nośności ramy
Przyjęto sieć typu WPB (ze wsteczną propagacja błędów) o strukturze 2-3-1 , a więc z jedną
warstwą ukrytą o dwuelementowym wektorze wejścia, trzech sigmoidalnych neuronach w
warstwie ukrytej i jednoelementowym wyjściu. Problem przewidywania obciążenia
krytycznego sformułowano jako odwzorowanie wektora wejś cia x w skalarne wyjście y :
x(2×1) = { Rs , Rr } y = Pk , (9)
gdzie: Rs  granica plastyczności słupów, Rr  granica plastyczności rygli, Pk  wartość
obciążenia krytycznego.
Do uczenia sieci przyjęto 256 różnych wzorców, zaś do jej testowania 100 wzorców.
Obliczenia wykonano za pomocÄ… symulatora neuronowego [10] korzystajÄ…c z metody
uczenia Back-Propagation z członem Momentum. Do obliczeń przyjęto 2000 epok. Tabl. 1
pokazuje dokładność predykcji neuronowej w porównaniu z wartościami Pk obliczonymi za
pomocą MES. W tablicy zestawiono błędy:
V
avr ep = 1/V " ep , max ep =max{epćł p = 1,...,V } dla ep =ćł1- y(p) / t(p)ćłÅ" 100 %
V
2
1
(p)
RMS = (f - y(p)) , (10)
"
V
p =1
gdzie: V = L , T - liczebność wzorców odpowiednio dla zbioru uczącego i testującego, t( p) ,
y( p) - wartości wyjściowe znane i obliczane siecią dla wzorca p . W Tabl. 1 podano też
wartości współczynnika korelacji dla zbioru par { t( p) , y( p) } .
Na rys. 2 pokazano wyniki predykcji neuronowej na tle wyników otrzymanych za
pomocą MES dla wartości obciążenia krytycznego ramy jako funkcji zmiennych granic
plastyczności rygli i słupów Pk (Rs , Rr ) .
Na rys. 2a pokazano rozkład wartości obciążenia krytycznego Pk , obliczonych siecią
neuronową podczas uczenia na L = 256 regularnie rozłożonych wartościach granic
plastycznoÅ›ci ( 16 × 16 wartoÅ›ci), a na Rys.2b sÄ… widoczne wyniki testowania na T = 100
wzorcach odpowiadających wylosowanym wartościom granic plastyczności. Przedstawione
400
400
400
400
81
wyniki wskazują na bardzo dobre właś ciwości predykcyjne sieci o dokładnościach podanych
w tabl. 2.
Tablica 1. Błędy aproksymacji neuronowej
Proces avr ep [%] max ep [%] RMS r
uczenie 0.379 2.128 0.004822 0.9996
testowanie 0.387 1.118 0.004760 0.9980
a) b)
1.3
1.4
1.2
1.2
1.1
1.0
1.0
0.8
0.9
0.8
300
400 0.6
300 350 300
250
300
250
Rr [MPa] 200
Rr [MPa]
200
200
200
Rs [MPa]
Rs [MPa]
100
150
100 150
Rys. 2. Porównanie wyników predykcji sieci z wynikami analizy sprężysto-plastycznej
dla wzorców: a) uczących i b) testujących, gdzie: Ą% - MES, " - sieć neuronowa
Na rys. 2a pokazano rozkład wartości obciążenia krytycznego Pk , obliczonych siecią
neuronową podczas uczenia na L = 256 regularnie rozłożonych wartościach granic
plastycznoÅ›ci ( 16 × 16 wartoÅ›ci), a na Rys.2b sÄ… widoczne wyniki testowania na T = 100
wzorcach odpowiadających wylosowanym wartościom granic plastyczności. Przedstawione
wyniki wskazują na bardzo dobre właś ciwości predykcyjne sieci o dokładnościach podanych
w tabl. 1.
5.2. Symulacja Monte Carlo
Rozpatrzono trzy przypadki różniące się przyjmowanymi do analizy zmiennymi
losowymi Xi. Dla każdej zmiennej losowej określono funkcję rozkładu gęstoś ci
prawdopodobieństwa f (Xi) oraz przyjęto parametry charakteryzujące rozkłady, tj. wartość
oczekiwanÄ… µ i odchylenie standardowe à . Charakterystyki poszczególnych przypadków
zestawiono w tabl. 2.
W przypadku 1 założono, że granice plastyczności zmieniają się jednakowo we
wszystkich słupach i ryglach, tj. losowane są wartości zmiennych X2 = Rs i X3 = Rr .W
przypadku 2 jest ustalona wartość granicy plastyczności w ryglach Rr = 240 MPa , a wartość
granicy plastyczności w słupach jest losowana jako zmienna X2 = Rs.. W przypadku 3 przy
stałej wartości granicy plastyczności słupów Rs = 240 MPa losuje się wartość granicy
plastyczności w ryglach X3 = Rr .
3
3
k
k
e
e
Obciażnie krytyczne P *10 [kN]
Obcią żnie krytyczne P *10 [kN]
82
Tablica 2. Symulowane przypadki
Przypa- zmienna losowa rozkład funkcji wartość oczekiw. odchyl. stand.
dek Xi f(Xi)
µ Ã
1 obciążenie S [kN] log.-normalny 6.592 0.2
Rs słupów [MPa] normalny 240 24
Rr rygli [MPa] normalny 240 24
2 obciążenie S [kN] log.-normalny 6.592 0.2
Rr słupów [MPa] normalny 240 24
3 obciążenie S [kN] log.-normalny 6.592 0.2
Rr rygli [MPa] normalny 240 24
We wszystkich trzech opisanych przypadkach, ocenę bezpieczeństwa konstrukcji
przeprowadzono dwiema metodami symulacji MC, tzn. metodÄ… klasycznÄ… KMC i metodÄ…
średniej ważonej WMC. We wszystkich symulacjach metodą WMC jako funkcję wagi g(Xi),
dla obciążenia funkcjÄ™ rozkÅ‚adu normalnego o parametrach µ =1042.3 kN i à =165 kN.
W tabl. 3 przedstawiono wyniki symulacji dla wszystkich rozpatrywanych
przypadków. Analiza wyników umieszczonych w tabl. 3. potwierdza lepszą skuteczność
metody WMC w porównaniu z metodą KMC.
Tablica 3. Wartości prawdopodobieństwa zawodności ramy obliczne metodami
Monte Carlo: klasyczna metoda MC (KMC) i metoda średniej ważonej (WMC
Liczba Prawdopodobieństwo zawodności q
próbek
Przypadek 1 Przypadek 2 Przypadek 3
KMC WMC KMC WMC KMC WMC
0.0000 0.0333
50 0.02000 0.05023 0.02000 0.05534
0.0400 0.0371
100 0.04000 0.04926 0.04000 0.05012
0.0550 0.0420
1000 0.05000 0.05393 0.06700 0.05560
0.0420 0.0415
2000 0.05850 0.05731 0.05900 0.05330
0.0446 0.0425
5000 0.05520 0.05697 0.06220 0.05374
0.0420 0.0415
10000 0.05550 0.05583 0.05860 0.05411
0.0418 0.0417
50000 0.05778 0.05623 0.05462 0.05627
0.0416 0.0416
100000 0.05509 0.055825 0.05498 0.05505
Rys. 3 pokazuje jak zmienia się prawdopodobieństwo niezawodnoś ci ramy dla
różnych wartoś ci obciążeń S = P, przyjmowanych jako stałe w symulacjach przypadków
1-3. Widać , że krzywe losowej niezawodnoś ci dla przypadków 1 i 2 pokrywają się.
Inaczej jest w przypadku 3 gdy losowo zmienia siÄ™ tylko granica plastycznoÅ› ci rygli.
Ustalona wartość granicy plastycznoś ci słupów powoduje, że od wartoś ci obciążenia
równego wartoś ci krytycznej Pkr = 1042 kN o noś ności ramy decyduje uplastycznienie
utwierdzenia słupów [6].
Wykresy losowej niezawodności ramy pozwalają oszacować wartość obciążenia jakie
można przykładać do ramy o założonym prawdopodobieństwie niezawodności. Na Rys. 3
pokazano przykładowo przyjęte prawdopodobieństwa Q = 0.9. W przypadku 1 losowej
zmienności granic plastyczności słupów i rygli można przykładać obciążenie P(0.9) = 927
kN, a w przypadku 3 obciążenie bezpieczne wynosi P(0.9) = 1023kN.
83
1
0.9
0.8
Ę%
Rr- zm.los.; Rs=const
"
Rr=const; Rs- zm.los.
0.6

Rr- zm.los.; Rs- zm.los.
0.4
0.2
927 1023
0.0
600 800 1000 1200 1400
Obciążenie S = P [kN]
Rys. 3. Prawdopodobieństwa niezawodności Q dla ustalonych wartości obciążeń S = P
5.3. Analiza czasów obliczeń
Obliczenia wykonywano na komputerze z procesorem AMD A1000 pod systemem
Windows98. W trakcie obliczeń pomierzono czasy generowania próbek (obliczanie
wartoÅ› ci krytycznej Pkr ( Rs , Rr ) ) za pomocÄ… programu ELPLAF-v1 oparto na [8] oraz
za pomocą sieci neuronowej BP+M: 2-3-1 programem [10]. Okazało się, że dla
przypadku jednego losowania wartoś ci granic plastycznoś ci, Rr i Rr czas obliczeń
programem MES wynosił od 10.13 do 83.16 sec i ś rednio ok. 32 sec ( dla 67 losowań
czas wynosił 2077 sec). Oznacza to, że dla przygotowania 356 wzorców do uczenia i
testowania sieci używając program MES był potrzebny czas obliczeń ok. 3.1 godz. W
taki sam sposób. Jeśli przyjmiemy metodę WMC i liczbę wzorców 1000 to stosując
program MES uśredniony czas obliczenia dla jednego rozwiązania (punkt na rys. 3)
wynosi 32 × 1000 =3.2 ×104 sec.
Obliczanie jednego rozwiÄ…zania (punkt na Rys.3) metodami Monte Carlo z
generowaniem 100000 próbek za pomocą sieci neuronowe wymagało 1.20 sec dla KMC i
1.80 sec dla WMC. Taki czas jest rzÄ™du 0.5 ×10-4 mniejszy niż czas obliczeÅ„ dla 1000
próbek za pomocą programu MES.
W pracy [6], gdzie zastosowano program MES oparty na założeniu przegubów
plastycznych i teorii rzÄ™du I dla 1000 próbek otrzymano czas ok. 0.6 ×104 sec. W pracy [7]
analizowano większe ramy niż rama na rys.3 i czasy obliczeń metodami MC dla 100000
próbek generowanych za pomocą sieci BPNN: 2-6-1 były rzędu 6.0 sec.
6. Wnioski końcowe
1. Zastosowanie sieci neuronowych do generowania próbek w metodach Monte Carlo
może istotnie skrócić czas obliczania krzywych prawdopodobieństwa niezawodności
konstrukcji.
2. W pracy analizowaliśmy stosunkowo prosty przykład ramy prostokątnej, dwunawo-
wej i czteropiętrowej dla dwóch losowych zmiennych materiałowych (granice plastyczności
wszystkich słupów i rygli). Bez istotnych zmian algorytmu można stosować opisane
podejście dla innych zmiennych losowych (np. zmian charakterystyk przekrojów,
geometrycznych imperfekcji prętów i całej ramy, losowych obciążeń nieproporcjonalnych
lub wieloparametrowych) i innych typów konstrukcji dla których możemy obliczać za
pomocÄ… MES wzorce do uczenia sieci.
prawd. niezawodności ramy Q
84
3. Pełniejsze oszacowanie numerycznej efektywności dyskutowanego podejścia wymaga
znacznie szerszej analizy numerycznej różnych konstrukcji i zapewne będzie istotnie zależne
od liczby i typu zmiennych losowych.
Literatura
[1] MURZEWSKI J., Niezawodność konstrukcji inżynierskich, Warszawa, Arkady, 1989.
[2] BIEGUS A., Probabilistyczna analiza konstrukcji stalowych, Warszawa-Wrocław,
PWN, 1999.
[3] RUBINSTEIN R.Y., Simulation and the Monte Carlo Method, New York, John Wiley
& Sons, 1981.
[4] ZIELIC
SKI R., Metody Monte Carlo, Warszawa, WNT, 1970.
[5] MAREK P., GUSTAR M., ANAGNOS Th., Simulation-Based Reliability Assessment
for Structural Engineers,Boca Raton, Florida, CRC Press, 1996.
[6] PULIDO J. E., JACOBS T. L., PRATES DE LIMA E. C., Structural reliability using
Monte Carlo simulation with variance reduction techniques on elastic-plastic structures.
Computers & Structures. 1992, Vol. 43, No. 3, s. 419-430.
[7] PAPADRAKAKIS M., PAPADOPOULOS V., LAGAROS N. D., Structural reliability
analysis of elastic-plastic structures using neural networks and Monte Carlo simulation.
Computer methods in applied mechanics and engineering. Vol. ... 1996, s. 145-163.
[8] WASZCZYSZYN Z., PABISEK E., Elastoplastic analysis of plane steel frames by a
new superelement., Archiwum Inżynierii Lądowej, (praca przyjęta do druku).
[9] MACHOWSKI A., Zagadnienia stanów granicznych i niezawodności szkieletów
budynków wielokondygnacjowych, Politechnika Krakowska, Seria Inż. Ląd., Monografia
262, 1999.
[10] DEMUTHN H., BEALE M., Neural Network Toolbox for Use with MATLAB, User's
Guide, Version 3, The Math Works, Inc., 1998.
PROBABILISTIC ESTIMATION OF STRUCTURAL RELIABILITY
BY NEURAL NETWORK SUPPORTED MONTE CARLO METHODS
Summary
In the paper a forward neural network is used for generating samples in the Monte Carlo
methods. The patterns for network training and testing are computed by a FE program.
A high numerical efficiency of the proposed approach is demonstrated on example of a plane
multi-storey elasto-plastic frame.
Praca została wykonana w ramach projektu badawczego KBN Nr 8 T07E 002 20 pt.
"Zastosowanie sztucznych sieci neuronowych do analizy konstrukcji stalowych".