Metoda Elementów Skończonych Studium magisterskie
SFORMUA
OWANIE WARIACYJNE PB
aproksymacja skończenie wymiarowa
metoda Ritza
WYKA
AD 4
Wersja elektroniczna, http://www.okno.pg.gda.pl.
Literatura
KLEIBER M.: Wprowadzenie do metody elementów skończonych. PAN IPPT, PWN Warszawa  Poznań 1989,
Rozdz. 3  4.2, str. 20  42.
CHRÓŚCIELEWSKI J., MAKOWSKI J., PIETRASZKIEWICZ W.: Statyka i Dynamika Powłok wielogałęziowych. Nieliniowa
teoria i metoda elementów skończonych. Wydawnictwo IPPT PAN, Warszawa2004,
Rozdz. 4.1, str. 269  282.
WILiŚ Politechnika Gdańska
J.Chróścielewski, M.Miśkiewicz, W.Witkowski: Katedra Mechaniki Budowli i Mostów
MES sformułowanie słabe
Charakterystyka sformułowania silnego (S)
statyka liniowa:
DT� + f = 0 na B statyczne równania równowagi (3r), RR;
� = Du na B *""B równania geometryczne (6r), RG;
� = E� na B *""B równania konstytutywne (6r), RK;
u = u* na "Bd kinematyczne warunki brzegowe, KB;
nT� = t* na "Bf statyczne (dynamiczne) warunki brzegowe SB.
Układ zamknięty, liczba niewiadomych pól (15, u 3, � 6, � 6)
równa się liczbie równań (15, RR,RG,RK)
i po uwzględnieniu warunków brzegowych (KB,SB) może być rozwiązany.
Poszukiwanie rozwiązań w formie silnej (S) problemu brzegowego MOC jest bardzo uciążliwe.
Istnieje niewielka liczba zadań o rozwiązaniach zamkniętych (analitycznych, ścisłych).
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/4
MES sformułowanie słabe
Charakterystyka sformułowania słabego
Ułatwieniem i podstawą wielu rozwiązań numerycznych jest tzw.
sformułowanie wariacyjne, nazywane także
słabym lub ogólnym.
Istota sformułowania słabego polega na:
skonstruowanie pewnego funkcjonału,
którego warunkami stacjonarności są
odpowiednie równania pól (RR, RG, RK)
i warunki brzegowe (KB, SB) uzyskane jako
równania Eulera Lagrange a tegoż funkcjonału.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/6
MES sformułowanie słabe
Charakterystyka sformułowania słabego
W ogólności problem brzegowy w słabej formie może być wyrażony jako
albo
tożsamość całkowa, np.: zasada wariacyjna, np.:
zasada pracy wirtualnej zasada stacjonarności energii potencjalnej,
zasada wirtualnych przemieszczeń, zasada Hellingera Meissnera
zasada wirtualnej pracy komplementarnej. zasada Hu-Washizu,
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/9
MES sformułowanie słabe
Zalety sformułowania słabego:
1) całkowa reprezentacja problemu zamiast różniczkowej,
wymaga spełnienia słabszych warunków ograniczeń na regularność zarówno danych jak i rozwiązań,
2) w zwarty i jasny sposób koncentruje w jednym wyrażeniu wszystkie równania problemu,
3) postać niezmiennicza  nie zależy od układu, funkcjonał ma zawsze określony sens fizyczny,
4) naturalna baza do rozwiązań aproksymacyjnych, umożliwia badanie błędów aproksymacji.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/13
MES sformułowanie słabe
Pojęcie wariacji, pochodna kierunkowa
DEFINICJA Dowolne, dostatecznie regularna pole �u"TuC(B) spełniające jednorodne warunki �u = 0
na brzegu "Bu ciała jest nazywane kinematycznie dopuszczalną wariacją.
Zwyczajowo wariacje �u nazywane są wirtualnymi przemieszczeniami (oznaczaną też przez w a" �u)
Odnotować warto, że jeśli w przestrzeni rozwiązań C(B) wprowadzimy definicję krzywej u�= u(�) = u + ��u,
przez punkt u�(0) = u to operacja wyznaczenia pochodnej kierunkowej w punkcie u w kierunku �u daje:
d
�[u;�u] = (u + ��u) |�=0= �u a" w
d�
czyli wariację możemy rozumieć jako pochodną kierunkową funkcji.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/15
MES sformułowanie słabe
Zasada pracy wirtualnej (zarys wyprowadzenia)
Odkształcenia wirtualne ��.
a) Niech pole u"C(B) opisuje konfigurację ciała B.
b) Niech �u a" w "TuC(B) jako element przestrzeni stycznej do C(B) w punkcie u
oznacza pole wirtualnych przemieszczeń.
Wirtualne odkształcenia ��, na podstawie pochodnej kierunkowej odkształceń zapisanych dla krzywej
u�= u(�) = u + ��u a" u + �w , u�(0) = u
na przestrzeni konfiguracyjnej C(B) parametryzowanej przez �" , mają formę
d
w naszym przypadku
�� = ��[u;�u] = �(u + ��u) |�=0a" Lin(�u) Ż#Ż#Ż#Ż#Ż#Ż# �� = D�u a" Dw .

d�
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/17
MES sformułowanie słabe
Zasada pracy wirtualnej (zarys wyprowadzenia)
Zasada pracy wirtualnej stwierdza, że obowiązuje
,
��T�dV = �uTfdV + �uTt*dA
+" +"+"
BB "Bf



wirtualna praca
wirtualna praca obciążeń zewnętrznych
wewnętrzna (naprężeń)
wtedy i tylko wtedy gdy:
a) spełnione są lokalne równania równowagi (RR),
b) spełnione są statyczne warunki brzegowe (SB),
c) przemieszczenia wirtualne są kinematycznie dopuszczalne to między innymi znaczy, że �u = 0 na "Bd .
Zasada pracy wirtualnej jest niezależna od definicji materiału.
Uzupełniając ją zależnością odkształcenie przemieszczenie (RG) (niekoniecznie liniową) oraz
relacjami konstytutywnymi (RK) (niekoniecznie liniowosprężystymi)
otrzymuje się słabe sformułowanie problemu obowiązujące także w zagadnieniach nieliniowych.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/19
MES sformułowanie słabe
Zasada wirtualnych przemieszczeń (zarys wyprowadzenia)
Wychodząc z zasady pracy wirtualnej, wyrażając
wirtualną pracę wewnętrzną ��T�dV ,
+"
B
poprzez przemieszczenia tj. naprężenia jako funkcję przemieszczeń,
wykorzystując równania konstytutywne (RK) i zależność geometryczną (RG):
w naszym przypadku
� = � �(u),x Ż#Ż#Ż#Ż#Ż#Ż# �(u) = E� = EDu), przy czym ET = E,
(
()
otrzymuje się funkcjonał zasady wirtualnych przemieszczeń
G[u;�u] = Gi[u;�u] - Ge[u;�u] = 0,
gdzie
Gi[u;�u] = (D�u)TE(Du)dV .
+"
B
Ge[u;�u] = �uTfdV + �uTt*dA.
+"+"
B "Bf
Człon opisujący pracę wirtualną sił zewnętrznych jest tu taki sam jak w zasadzie pracy wirtualnej.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/22
MES sformułowanie słabe
Przykład
Sformułuj zasadę wirtualnych przemieszczeń dla liniowo sprężystego pręta prostego jak na rysunku
jedna współrzędna opisu dziedziny (s a" x), jeden parametr teorii (u), jeden węzłowy stopień swobody (ua).
(S)  silne, klasyczne sformułowanie problemu brzegowego (WM)
2 2
N + f = 0, � = u , N = EA�, + warunki brzegowe: u|x = 0 = u* , N|x = L = t* .
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/23
MES sformułowanie słabe
Rozwiązanie:
2
a) Mnożymy równanie równowagi N + f = 0 przez �u i całkujemy wzdłuż długości pręta
(N,x + f )�uds = 0 �! N,x �uds + f �uds = 0
+" +" +"
L L L
b) Obliczamy pochodną iloczynu (N�u),x = N,x �u + N�u,x �! N,x �u = (N�u),x -N�u,x
c) Formalnie stosujemy twierdzenie o dywergencji do członu (N�u),x ds = (�uN)|x = L
+"
L
d) Wracamy do punktu a) i otrzymujemy
N,x �uds + f �uds = 0 �! - N�u,x ds + f �uds + (�uN)|x = L = 0 �! N�u,x ds = f �uds + (�uN)|x = L
+"+" +" +" +" +"
L L L L L L
e) Ponieważ � = u,x �! �� = �u,x oraz N = EA� to otrzymujemy poszukiwany wynik w postaci
(W)  wariacyjne, słabe sformułowanie problemu brzegowego,
��Nds = �uf ds + (�uN)|x = L �! ��EA�ds = w f ds + wt* , wx) a" �u(x),
(
+" +" +" +"
L L L L
Tym samym, w toku ścisłych matematycznych przekształceń dowiedliśmy równoważności
(S) �! (W )


ciągły
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04A/28
MES aproksymacja skończenie wymiarowa
Metody numeryczne (MN) rozwiązywanie zadań matematycznych za pomocą działań arytmetycznych.
Zachodzi potrzeba przybliżania wielkości niearytmetycznych wielkościami
arytmetycznymi i badania błędów wywołanych takim postępowaniem.
Aproksymacja, uwagi ogólne
1. Podstawowym celem aproksymacji jest osiągnięcie pewnej żądanej dokładności obliczenia
przybliżanej wielkości, a także możliwość jej oszacowania.
Pytanie: Jaki jest dopuszczalna błąd wyniku?
Pytanie ogólne: Co rozumiemy przez błąd rozwiązania?
(np. błąd względy, bezwzględny, przedziałowy, błąd średni - wieloznaczność pojęcia).
2. Określenie ilości obliczeń potrzebnych do uzyskania przybliżenia.
Pytanie: Jak szybko można otrzymać rozwiązanie używając danej metody?
Zakładamy, że
a) żądaną dokładność można osiągnąć,
b) wszystkie rozpatrywane metody ,,działają (tzn. prowadzą do sensownego wyniku, np. zbieżność iteracji).
Klasy funkcji aproksymacyjnych
Wiele zadań MN polega na przybliżaniu funkcji f (x) kombinacją (najczęściej liniową) funkcji należących do
pewnej szczególnej klasy, np.:
- N pierwszych wyrazów szeregu Taylora,
- wielomianami stopnia n (np. wzór trapezów, MES),
- klasa funkcji trygonometrycznych - dla funkcji okresowych,
- funkcje wykładnicze, itp.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04B/5
MES aproksymacja skończenie wymiarowa
Problem ciągły (C), sformułowanie
Zastępowanie układu typu kontynualnego (o nieskończonej liczbie stopni swobody)
np. funkcjonału prac wirtualnych G :U �TU 
G[u;�u] = Gi[u;�u] - Ge[u;�u] = (D�u)TE(Du)dV - ( �uTfdV + �uTt*dA),
+" +"+"
BB "Bf
układem dyskretnym (o skończonej liczbie stopni swobody).
Problem C: Znalez
ć kinematycznie dopuszczalne pole przemieszczeń u"U �"C(B) takie, że
funkcjonał
G[u;�u] = 0,
dla wszystkich kinematycznie dopuszczalnych wirtualnych przemieszczeń
�u"TuU �" TC(B), to znaczy spełniających jednorodne przemieszczeniowe warunki
brzegowe �u = 0 na "Bd .
Tutaj C określa problem Ciągły.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04B/7
MES aproksymacja skończenie wymiarowa
Problem dyskretny (D), aproksymacja skończenie wymiarowa
a) Niech Ch(Bh) jest skończenie wymiarową podprzestrzenią C(B).
b) Niech TCh(Bh) oznacza odpowiednią wiązkę styczną do niej.
c) Oznaczamy skończenie wymiarowe odpowiedniki pól xh"Bh dla x ,
uh"Ch(Bh) dla u,
�uh"Tu Ch(Bh) dla �u.
h
h  parametr charakterystyczny aproksymacji taki, że
Ch(Bh) C(B) jeśli h " (alternatywnie h 0) w sensie definicji (przy Bh B lub Bha" B ).
Aproksymację skończenie wymiarową (problem dyskretny) do Problemu C formułuje się jako:
Problem D: Znalez
ć kinematycznie dopuszczalne pole przemieszczeń uh"Uh �"Ch(Bh) takie, że
funkcjonał
G[uh;�uh] = 0
dla wszystkich kinematycznie dopuszczalnych wirtualnych przemieszczeń
�uh "Tu Uh �" TCh(Bh), to znaczy spełniających jednorodne przemieszczeniowe
h
warunki brzegowe �uh = 0 na ("Bd )h .
Mówi się, że rozwiązanie uh problemu D jest aproksymacją zbieżną do rozwiązania u
problemu C jeśli uh u przy h ".
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04B/13
MES metoda Ritza
Metoda Rayleigha-Ritza (R-R), koncepcja
Aproksymacja poszukiwanego rozwiązania u(x)"U na całym rozpatrywanym obszarze x " B za pomocą
kombinacji liniowej układu h liniowo niezależnych znanych funkcji Ćn (n =1,2,...,h) (tworzących bazą
skończenie wymiarowej przestrzeni rozwiązań dopuszczalnych Uh ):
N
uh(x) = ha (x)aa = h(x)a,
"
a=1
x " B , uh "Uh �" Ch(B).
Jednoczesna i jednolita aproksymacja Ća
wszystkich trzech składowych przemieszczeń:
Ća (x) 0 0 aa
Ą# ń# ż# �#
�#b �#
ó# Ą#
ha (x) = 0 Ća (x) 0 , aa = ,
�# Ź#
a
ó# Ą#
�#c �#
ó# 0 0 Ća (x)Ą#
Ł# Ś# �# a �#
a1
ż# �#
�#a �#
�# �#
2
h(x) = h1(x) h2(x) ... hN (x) , a = .
[ ]
�# Ź#
...
�# �#
�#aN �#
�# �#
Funkcje h, (Ćn )  funkcje testowe,
parametry a  współrzędne uogólnione.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04C/3
MES metoda Ritza
Metoda Rayleigha-Ritza (R-R), tok postępowania
a) Poszukiwana funkcja uh musi spełniać kinematyczne warunki brzegowe,
tj. tworzące ją funkcje dobierane są tak, aby aproksymacja uh= u* na brzegu "Bd
niezależnie od wartości parametrów a.
b) Dziedzina problemu B i jej brzeg "Bd nie podlega aproksymacji.
N
uh(x) = ha(x)aa = h(x)a, x " B i uh(x) = u*, x ""Bd , uh "Uh �" Ch(B).
"
a=1
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04C/6
MES metoda Ritza
Metoda Rayleigha-Ritza (R-R), tok postępowania
c) Zakłada się (w R R) jednorodne warunki brzegowe (kinematyczne) tj. u = uh= u*= �u a" 0 wzdłuż "Bd .
Pozwala to na przyjęcie tej samej bazy aproksymacyjnej dla przemieszczeń wirtualnych
N
�uh(x) = ha(x)�aa = h(x)�a , x " B i �uh(x) = 0, x ""Bd , �uh"Tu Uh.
"
a=1 h
d) Podstawienie postulowanej postaci rozwiązania do funkcjonału G :U �TuU 
G[uh;�uh] a" G[a;�a] = (Dh�a)TE(Dha)dV - ( (h�a)T fdV + (h�a)T t*dA)
+" +"+"
BB "Bf
.
= �aT[[ hTDTEDhdV ]a - ( hTfdV + hTt*dA)]= �aT[Ka - r ]
+"+" +"
BB "Bf




K
r
e) Warunek G[a;�a] = 0 wobec dowolności �a, a więc dla każdego �a `" 0, prowadzi do
układu h = 3N równań algebraicznych względem parametrów a
G[a;�a] = �aT Ka - r = 0 �! Ka = r .
[ ]
f) Wyznaczenie poszukiwanego rozwiązania przybliżonego
a = K-1r �! uh(x) = h(x)a = h(x)K-1r .
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04C/10
MES metoda Ritza
Problemy związane ze stosowaniem metody Rayleigha Ritza:
a) Dobór funkcji aproksymujących rozwiązanie w całym obszarze x " B . Funkcje te musza nie tylko spełniać
warunki brzegowe, ale także umożliwiać opis wszystkich charakterystyk zadania np.: geometrii układu,
stałych materiałowych itp.
Efektywne generowanie Ćn w zagadnieniach praktycznych jest bardzo uciążliwe.
Na przykład, jakie funkcje przyjąć do opisu geometrii następującego zadania?
b) Całkowanie wielomianów aproksymujących wysokich rzędów jest uciążliwe.
c) Brak przesłanek, jakie funkcje aproksymujące należy dodać, aby poprawić dokładność rozwiązania.
d) Trudności z interpretacją współrzędnych uogólnionych, które nie mają czytelnej interpretacji fizycznej.
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04C/15
MES metoda Ritza
Przykład
Pręt wspornikowy, o stałym przekroju poprzecznym obciążony stycznym obciążeniem o rozkładzie parabolicznym.
Dane: EA =1000, L = 5, p0 = 1.
Specyfikacja macierzy zadania:
u(x) = {u(x)}, x" B ; f(x) = { f (x)} = p0{x2}, x" B ;
u*(x) = {u *(x)}, x ""Bd �! u*|x=0 = 0; t*(x) = {t *(x)}, x ""Bf �! t*|x=5 = 0;
"
Ą# ń#, E = EA, K = EA L(Dh)T DhdV , r = LhTf dA + LhTt*dA
D = DT =
+" +" +"
ó#"x Ą#
0 00
Ł# Ś#
Rozwiązań: poszukujemy przyjmując (arbitralnie) funkcje bazowe w postaci xn .
Rozpatrujemy dwa przypadki macierzy h(x):
a1
ż# �#
a1
ż# �#
�#a �#.
h(x) = [x1 x2], a =
�#a Ź# oraz h(x) = [x1 x2 x3], a = �# Ź#
2
�# 2 �#
�#a �#
�# 3 �#
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04C/18
MES metoda Ritza
Rozwiązanie z h(x) = [x1 x2]
Rozwiązanie analityczne
0.15
0.125
0.1
0.075
0.05
0.025
1 2 3 4 5
Długość pręta
uaproks (x =1) = 0.04625, uanalit (x =1) = 0.0415833
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04C/19
MES metoda Ritza
Rozwiązanie z h(x) = [x1 x2 x3]
Rozwiązanie analityczne
0.15
0.125
0.1
0.075
0.05
0.025
1 2 3 4 5
Długość pręta
uaproks (x =1) = 0.04125, uanalit (x =1) = 0.0415833
WILiŚ Politechnika Gdańska
J. Chróścielewski, M.Miśkiewicz, W. Witkowski: KMBiM W04C/20
Metoda Elementów Skończonych Studium magisterskie
Dziękuję za uwagę
cdn.
WILiŚ Politechnika Gdańska
J.Chróścielewski, M.Miśkiewicz, W.Witkowski, Katedra Mechaniki Budowli i Mostów