1. Sposoby sformułowania zagadnień brzegowych (lokalne, globalne) i koncepcje
budowy rozwiązań przybliżonych.
Zagadnienie (problem) brzegowe: dany jest obszar
Ω
, w którym poszukiwane jest
rozwiązanie, układ równań różniczkowych cząstkowych oraz warunki początkowo –
brzegowe nałożone na zbiór punktów należących do brzegu
∂Ω obszaru.
Ω
∂Ω
( )
P x
W rozważanym obszarze poszukiwana jest funkcja
( )
u x w każdym punkcie ( )
P x . Można
stosować następujące sformułowania zagadnień brzegowych:
• Sformułowanie lokalne (mocne, silne): szukane jest rozwiązanie układu równań
różniczkowych w każdym z punktów obszaru osobno:
u
f
dla P
u g dla P
=
∈
=
∈
L
B
Ω
∂Ω
∈∂Ω
gdzie
są operatorami różniczkowymi odpowiednio w obszarze i na jego
brzegu. Równanie
nosi nazwę warunków brzegowych. Jeżeli
są one nałożone na funkcję (tzn.
i
L B
u g dla P
=
B
1
≡
B
), noszą nazwę podstawowych warunków
brzegowych Dirichleta, natomiast dowolna kombinacja warunków brzegowych
złożona z pochodnych nosi nazwę nosi nazwę naturalnych warunków brzegowych
Neumanna.
• Sformułowanie globalne : może być formułowanie jako problem optymalizacji
funkcjonału lub jako zasada wariacyjna.
o
Minimalizacja funkcjonału:
1
( )
( , )
( )
2
I u
u u
=
−
B
L u
1
W funkcjonałach energetycznych pierwszy składnik prezentuje energię
wewnętrzną układu, podczas gdy drugi jest równy pracy wykonanej przez siły
zewnętrzne. Nieznana funkcja
może przedstawiać sobą przemieszczenia
, odkształcenia
( )
u P
u
ε , naprężenia σ lub wszystkie z nich. Funkcja
u
realizująca
ekstremum (minimum, punkt stacjonarny) funkcjonału
( )
min ( )
u
I u
jest szukana.
Można rozważać problem optymalizacji funkcjonału bez ograniczeń (w całej
przestrzeni rozwiązań dopuszczalnych) lub z ograniczeniami (ekstremum jest
szukane w podprzestrzeni narzuconych ograniczeń).
o
Zasada wariacyjna
( ,
)
( )
u u
u
dla
u V
∂ =
∂
∂ ∈
B
L
W mechanice powyższe równanie może mieć sens np. zasady prac
wirtualnych. Sformułowanie wariacyjne (tzw. słabe) ma podstawowe
znaczenie przy konstruowaniu rozwiązań przybliżonych. Można go uzyskać ze
sformułowania mocnego w czterech krokach:
Przemnożenie równania różniczkowego przez dowolną funkcję (tzw.
funkcja testująca),
Przecałkowanie wyniku po rozważanym obszarze Ω ,
Całkowanie przez części z wykorzystaniem twierdzenia Greena w celu
zredukowania pochodnych do minimalnego rzędu,
Wprowadzenie do funkcjonału warunków brzegowych Neumanna.
Sformułowania globalne wymagają dodatkowego całkowania po obszarze.
Sformułowanie wariacyjne jest ogólniejsze, gdyż możliwe jest w przypadku
wszystkich zagadnień brzegowych, podczas gdy ułożenie funkcjonału możliwe
jest tylko dla niektórych zadań mechaniki, np. dla zadań liniowej sprężystości
(funkcjonał Lagrange’a, Hamiltona, Reissnera, Castigliano, itp.).
• Możliwe są również podejścia mieszane, polegające np. na podziale obszaru Ω na
podobszary, gdzie stosuje się różne sformułowania wraz z odpowiednimi warunkami
ograniczającymi.
Budowa rozwiązania przybliżonego problemu brzegowego zależy przede wszystkim od
wybranej metody dyskretnej. Można wyróżnić dwie główne koncepcje:
• Rozwiązanie dyskretne w postaci kombinacji liniowej współczynników liczbowych
oraz funkcji bazowych:
1 1
2 2
1
( )
( )
( ) ...
( )
( )
n
n n
i i
i
p x
a
x
a
x
a
x
a
x
ϕ
ϕ
ϕ
ϕ
=
=
+
+ +
=
∑
.
Funkcje bazowe (najczęściej: wielomiany, funkcje trygonometryczne, funkcje
specjalne) muszą być liniowo niezależne, odpowiednio ciągłe oraz muszą spełniać
jednorodne warunki brzegowe rozważanego problemu (jednorodne warunki to takie, w
których po prawej stronie stoi 0, (np.
0
0
( ) 0,
( ) 0
I
u x
u x
=
= ). Przy takim zapisie
postaci rozwiązania przybliżonego można szukać budując odpowiednie residua, (czyli
wyrażenia świadczące o spełnieniu przez rozwiązanie przybliżone wyjściowych
równań różniczkowych) odpowiednio w obszarze i na brzegu:
2
( )
( )
,
( )
( )
d
b
x
p x
f
x
p x
g
ε
ε
=
−
=
L
B
−
Funkcjonał ważący powyższe wyrażenia ma postać:
( )
d
d
b
b
I p
w d
w d
ε
ε
Ω
∂Ω
=
Ω +
∫
∫
∂Ω
b
.
Wagi
świadczą o odejściu p(x) od wyniku ścisłego odpowiednio w obszarze i
na jego brzegi. Dla metod residuów ważonych (metoda Bubnowa - Galerkina, metoda
najmniejszych kwadratów, metoda kolokacji) i metod energetycznych (metoda
Rayleigha – Ritza) zakłada się błąd na brzegu
d
w i w
0
b
ε
= (ścisłe spełnienie warunków
brzegowych) i rozważa jedynie
( )
d
d
I p
w d
ε
Ω
=
Ω
∫
. Odmienną koncepcję prezentują
tzw. metody Trefza, w których zakłada się ścisłe spełnienie równania wewnątrz
obszaru a rozwiązań przybliżonych poszukuje na jego brzegu.
• Rozwiązanie dyskretne w wybranych punktach obszaru (lub/i jego brzegu) zwanych
węzłami. W tej koncepcji niezbędna jest dyskretyzacja obszaru (na węzły, elementy
itp.), gdzie zastępuje się wielkości ciągłe wielkościami dyskretnymi. Numeryczne
wyniki dyskretne można aproksymować funkcją ciągłą w ramach tzw. postprocesingu.
Do tych metod należą: metoda różnic skończonych (MRS, zamiana operatorów
różniczkowych na różnicowe, poszukiwanie wartości węzłowych funkcji szukanej,
aproksymacja metodami najmniejszych kwadratów), metoda elementów skończonych
(MES, podział na elementy i aproksymacja funkcjami kształtu) oraz metod elementów
brzegowych (MEB, podział brzegu na odcinki, obliczanie całek brzegowych).
Przykład
Belka swobodnie podparta obciążona obciążeniem ciągłym równomiernie rozłożonym.
Sformułowanie lokalne:
2
2
( )
( )
( )
( )
( )
0
2
d
M x
x
w x
f x
f x
x
l
dx
EJ
=
=
= −
≤
L
≤
1
( )
(2
)
(0) 0
(2 ) 0
2
M x
qx l x
w
w l
=
−
=
=
3
Sformułowanie globalne:
• W postaci funkcjonału:
2
2
0
1
( )
min ( )
( )
[ (
)
] ,
(0)
(2 ) 0
2
l
w
dw
M x
I w
I w
w dx
w
w l
dx
EJ
⇒
=
−
=
∫
=
• W postaci zasady wariacyjnej:
2
2
2
0
( )
[
] ( )
0,
(0)
(2 ) 0
( ) funkcja próbna, odpowiednio ciągla, spelnia warunki brzegowe: v(0)
(2 ) 0
l
dw
M x
v x dx
w
w l
dx
EJ
v x
v l
+
=
=
=
−
=
∫
=
lub po przecałkowaniu przez części (sformułowanie słabe):
2
0
( )
[
( )]
0,
(0)
(2 ) 0,
(0)
(2 ) 0
l
dw dv M x
v x dx
w
w l
v
v l
dx dx
EJ
−
=
=
=
=
∫
=
Rozwiązanie przybliżone dla metod residualnych:
• Funkcje bazowe:
2
1
2
( )
(
2 ),
( )
(
2 )
x
x x
l
x
x x
l
ϕ
ϕ
=
−
=
−
,
• Rozwiązanie próbne:
2
1
2
( )
( )
( )
(
2 )
(
2 )
p x
a
x
b
x
a x x
l
b x x
l
ϕ
ϕ
=
+
=
−
+
−
,
• Residuum w obszarze:
2
2
( )
( )
1
( )
2
(6
4 )
(2
)
2
d p x
M x
x
a b x
l
qx l x
dx
EJ
ε
=
+
=
+
−
−
−
• Dla metody Bubnowa - Galerkina:
2
2
1
0
0
2
2
2
2
0
0
1
( )
( )
0
[2
(6
4 )
(2
)] (
2 )
0
2
,
...
1
( )
( )
0
[2
(6
4 )
(2
)] (
2 )
0
2
l
l
l
l
x
x dx
a b x
l
qx l x x x
l dx
a b
x
x dx
a b x
l
qx l x x x
l dx
ε
ϕ
ε
ϕ
⎧
⎧
⋅
=
+
−
−
−
−
=
⎪
⎪
⎪
⎪
⇒
⇒
⎨
⎨
⎪
⎪
⋅
=
+
−
−
−
−
=
⎪
⎪
⎩
⎩
∫
∫
∫
∫
=
• Dla metody najmniejszych kwadratów:
2
( , )
0
2
2
0
( , )
( ) ( )
min ( , )
( , ) 0
1
( , )
[2
(6
4 )
(2
)]
,
...
2
( , ) 0
l
a b
l
I a b
x
x dx
I a b
I a b
a
I a b
a b x
l
qx l x dx
a b
I a b
b
ε
ε
=
⋅
⇒
∂
⎧
=
⎪⎪∂
=
+
−
−
−
⇒
⇒
=
⎨ ∂
⎪
=
⎪∂
⎩
∫
∫
• Dla metody kolokacji (punkty kolokacji:
1
2
2
,
3
3
l
l
x
x
=
=
):
4
2
1
1
1
1
0
1
2
2
1
2
2
2
0
1
( ) (
)
0
2
(6
4 )
(2
) 0
( ) 0
2
,
...
( ) 0
1
2
(6
4 )
(2
) 0
( ) (
)
0
2
l
l
x
x x dx
a b x
l
qx
l x
x
a b
x
a b x
l
qx
l x
x
x x dx
ε
δ
ε
ε
ε
δ
⎧
⎧
⋅
−
=
⎪
+
−
−
−
=
⎪
=
⎧
⎪
⎪
⇒
⇒
⇒
⎨
⎨
⎨
=
⎩
⎪
⎪
+
−
−
−
=
⋅
−
=
⎪
⎪
⎩
⎩
∫
∫
=
W metodzie różnic skończonych MRS wprowadzono w ramach dyskretyzacji obszaru 3
węzły (patrz: rysunek). Z trzech wartości węzłowych dwie z nich stanowią warunki
brzegowe:
, pozostaje do obliczenia wartość
. Przy sformułowaniu lokalnym
zamianie na operator różnicowy ulega operator różniczkowy na drugą pochodną:
0
2
0
w
w
=
=
1
w
2
0
1
1
1
2
2
II
x l
w
w
w
d w
w
Lw
dx
l
=
−
+
=
≈
=
2
2
. Równania różnicowe generuje się metodą kolokacji
(ścisłe spełnienie równania w węzłach obszaru):
}
}
0
0
2
0
1
2
1
1
1
2
2
1
...
2
w
w
w
ql
Lw
f
w
l
EJ
−
+
=
⇒
=
⇒
=
W sformułowaniu globalnym można ułożyć funkcjonał energii potencjalnej układu. Po jego
dyskretyzacji (całkowanie kwadraturą Newtona-Cotesa między węzłami) otrzymuje się:
2
2
1
0
2
1
0
1
2
0
0
1 1
1 1
2
2
1
( , ,
)
[(
)
(
)
(
)
(
)
]
2
2
w
w
w
w
l
l
I w w w
l
l
M w
M w
M w
M w
l
l
EJ
EJ
−
−
=
+
−
+
−
+
2
Niewiadomą
(oczywiście
) otrzymuje się minimalizując powyższy funkcjonał
względem
:
1
w
0
2
0
w
w
=
=
1
w
1
1
1
( ) 0
...
d
I w
w
dw
=
⇒
=
Przy sformułowaniu wariacyjnym słabym (funkcja testowa: (0)
(2 ) 0
v
v l
=
= ) od razu
otrzymuje się gotowe równanie różnicowe:
1
0
1
0
2
1
2
1
0
0
1 1
1 1
2
2
(
)
(
)
2
2
w
w v
v
w
w v
v
l
l
l
l
M w
M w
M w
M w
l
l
l
l
EJ
EJ
−
−
−
−
+
−
+
−
+
0
= .
Podstawiając
i przyrównując wyrażenie stojące przy
dowolnym do zera otrzymuje się wartość
.
0
2
0
2
0 oraz
0
w
w
v
v
=
=
=
=
1
v
1
w
2. Wybrane metody numeryczne, które mają zastosowanie w komputerowych
modelach mechaniki.
Podczas modelowania komputerowego ważnym etapami są: stworzenie modelu
numerycznego badanego obiektu rzeczywistego, co wiąże się z zamianą wielkości ciągłych na
wielkości dyskretne (punktowe), a następnie z implementacją komputerową algorytmów
działających na tych dyskretnych wielkościach. Niezbędne są tu metody numeryczne, które w
sposób zrozumiały dla maszyny cyfrowej pozwalają na rozwiązywanie podstawowych
problemów algebry i analizy funkcjonalnej.
5
• Rozwiązywanie równań algebraicznych
Równania nieliniowe (z jedną zmienną) pojawiają się w przypadku materiału
nieliniowego (np. w plastyczności) lub mogą mieć pochodzenie geometryczne (np.
uwzględnienie dużych ugięć). Rozwiązuje się je metodami iteracyjnymi, w których na
podstawie znajomości tzw. schematu iteracyjnego
1
( )
n
n
x
f x
+
=
, punktu startowego
0
x
oraz żądanej dokładności rozwiązania otrzymuje się ciąg przybliżeń miejsca zerowego
równania ( ) 0
F x
= . Istnieje kilka metod iteracyjnych, np.:
o
Metoda Newtona (metoda stycznych):
1
0
( )
,
...
'( )
n
n
n
n
F x
x
x
x
F x
+
=
−
= ,
o
Metoda siecznych:
1
1
1
( )
(
),
...,
...
( )
(
)
n
n
n
n
n
n
n
F x
x
x
x
x
x
x
F x
F x
+
−
−
=
−
−
=
=
−
0
1
o
Metoda relaksacji (optymalnie szybka metoda wokół punktu startowego):
0
1
0
0
0
'( )
1
( )
,
...,
1
'( )
1
'( )
n
n
n
f x
x
f x
x
f x
f x
+
=
−
−
−
x
=
Każda z metod ma swoją interpretację geometryczną.
• Rozwiązywanie układów równań algebraicznych.
Konieczność taka występuje właściwie w każdej metodzie dyskretnej, jest zazwyczaj
ostatnim zabiegiem algebraicznym, po rozwiązaniu którego otrzymuje się niewiadome
pierwotne, np. w statyce budowli – przemieszczenia węzłowe układu konstrukcyjnego
(MRS, MES, Metoda sił, Metoda przemieszczeń).
W przypadku układu równań nieliniowych (nieliniowość materiałowa lub
geometryczna) stosuje się najczęściej iteracyjną metodę Newtona - Rhapsona. Dla
rozwiązania układu równań:
( ) 0,
1, 2,...
i
F
i
=
=
x
n
przyjmuje się punkt startowy
obliczeń
(0)
(0)
(0)
0
1
2
(
,
,...,
)
n
x
x
x
=
x
buduje się następujący schemat iteracyjny:
1
( )
( )
( )
n
n
n
n
n
+
⋅
=
⋅
−
J x
x
J x
x
F x
, gdzie:
1
1
1
1
( ) ...
( )
( )
...
...
...
( ) ...
( )
n
n
n
n
d
d
F
F
dx
dx
d
d
F
F
dx
dx
⎡
⎤
⎢
⎥
⎢
⎥
= ⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎣
⎦
x
x
J x
x
x
- macierz Jacobiego.
W przypadku układu równań liniowych:
(
) ( 1)
( 1)
(det
0)
n n n
n
×
×
×
=
≠
A x
b
A
metody jego
rozwiązywania można podzielić na eliminacyjne (dają wyniki ścisłe, wymagają
dużego nakładu pracy), iteracyjne (dają wyniki obarczone błędem, łatwe w
implementacji), kombinowane i specjalne (metody macierzy rzadkich, metody analizy
frontalnej). Do metod eliminacyjnych należą m.in. popularna metoda eliminacji
Gaussa (polegająca na doprowadzeniu macierzy A do postaci diagonalnej) i metoda
eliminacji Choleskiego (dla macierzy symetrycznych, dodatnio określonych, polega na
6
rozkładzie macierzy A na macierze samosprzężone: dolno- i górno-trójkątną). Do
metod iteracyjnych zaliczyć można metody iteracji Jacobiego i Gaussa – Seidle’a
będące uogólnieniem metody iteracji prostej dla równań algebraicznych. Powyższe
metody mogą służyć również do odwracania macierzy.
• Rozwiązywanie problemu własnego.
Konieczność taka pojawia się np. w zadaniach stateczności konstrukcji (poszukiwanie
siły krytycznej) lub w zadaniach dynamiki budowli (poszukiwanie częstości drgań
własnych). Algebraicznie problem własny można zapisać następująco:
(
) ( 1)
( 1)
(
) 0,
n n
n
n
λ
λ
λ
×
×
×
⋅
= ⋅
⇒
−
=
∈
A
x
x
x A
I
ℜ
.
Aby powyższy układ miał niezerowe rozwiązanie, (ale niejednoznaczne!), musi być
spełniony warunek:
det(
) 0
λ
−
=
A
I
. Otrzymuje się w ten sposób równanie
wielomianową na niewiadomą
λ , zwane równaniem charakterystycznym. Po jego
rozwiązaniu dostaje się wartości własne ,
1, 2,...,
i
i
n
λ
=
. Wstawiając je kolejno do
wyjściowego równania otrzymuje się n wektorów własnych
,
1, 2,...,
i
i
n
=
x
,
jednoznacznych, co do długości bądź kierunku. W przypadku, gdy wszystkie wartości
własne są równe 0, to zbiorem wektorów własnych jest cała hiperpłaszczyzna.
Problemu własnego nie rozwiązuje się korzystając z definicji, jako iż jest to
czasochłonne. Ponieważ większość problemów własnych z dziedzin mechaniki
dotyczy macierzy symetrycznych, tak więc metody numeryczne są właśnie dla nich
opracowane. Metodą numeryczną można rozwiązywać problem własny kilkoma
wariantami: można szukać wybranej wartości własnej (i odpowiadającego jej wektora
własnego), wartości własnych z wybranego przedziału lub wszystkich wartości
własnych i wektorów własnych (np. metoda Jacobiego). Do pierwszej grupy metod
należy iteracyjna metoda potęgowa i jej modyfikacje. Wykorzystuje ona własności
tzw. ilorazu Rayleigha. Czysta metoda potęgowa prowadzi do znalezienia wartości
własnej największej co do modułu oraz odpowiadającego jej wektora własnego. Jej
algorytm można zapisać następująco:
0
0
0
0
1
1
1
1
1
max
( )
( )
1
1
...
,
0,1, 2,...
,
,
k
k
k
k
k
T
k
k
k
k
k
k
k
k
k
k
v
v
k
λ
λ
λ
λ
λ
ε
ε
λ
+
+
+
+
+
−
−
=
=
⎧
=
=
⎪
⎪
⎪
=
⇒
=
⇒
⎨
⎪
−
⎪
=
=
−
⎪⎩
max
v
x
x
x
x
x
Av
x
v x
v
v
v
Aby znaleźć kolejne wartości własne, można zastosować powyższy algorytm z
odpowiednimi modyfikacjami, np. rozwiązując problem własny metodą potęgową, ale
dla macierzy odwrotnej
-1
A znajdzie się wartość własną najbliższą zeru a przesuwając
dodatkowo widmo (zbiór wartości własnych) macierzy o zadaną liczbę znajdzie się
wartość własną najbliższą danemu przesunięciu .
d%
d%
7
W zagadnieniach mechaniki teorii elementów skończonych ma się też często do
czynienia z tzw. uogólnionym problemem własnym:
(
)
(
)
( 1)
( 1)
,
n n
n n
n
n
λ
λ
×
×
×
×
⋅
= ⋅
⋅
∈
A
x
B
x
ℜ
,
różniącym się od postaci zwykłej obecnością macierzy (z założenia symetrycznej).
B
Najczęściej problem uogólniony sprowadza się do postaci standardowej:
1
T
λ
−
−
⋅ = ⋅
=
⋅ ⋅
y,
L
A
C y
C
L ,
gdzie macierze dolno- i górnotrójkątne (
T
L , L
) pochodzą z rozkładu macierzy
i rozwiązuje dalej jedną z poznanych wcześniej metod Istnieją też
odpowiednie modyfikacje metody potęgowej dla uogólnionego problemu, bez
potrzeby jego liniowej transformacji do postaci zwykłej.
T
=
⋅
B
L L
• Aproksymacja funkcji
Zagadnienie aproksymacji funkcji f(x) (ciągłej – danej wzorem analitycznym bądź
dyskretnej – danej zbiorem punktów) polega na jej zastąpieniu przez inną funkcję p(x)
w danym przedziale
[ ]
,
x
a b
∈
. Poszukiwana funkcja aproksymująca powinna mieć
łatwy wzór, aby jej dalsza analiza była stosunkowo prosta. Zagadnienie aproksymacji
zapisuje się następująco.
0 0
1 1
0
( )
( )
( )
( ) ...
( )
( )
m
m m
i i
i
f x
p x
a
x
a
x
a
x
a
x
ϕ
ϕ
ϕ
ϕ
=
≈
=
+
+ +
=
∑
m oznacza stopień aproksymacji – jest to liczba związana z ilością liniowo
niezależnych funkcji bazowych ( ),
0,1,...,
i
x
i
m
ϕ
=
przyjmowanych z góry
(najczęściej wielomiany, mogą to być także funkcje trygonometryczne lub funkcje
specjalne np. wielomiany Lagrange’a, Czebyszewa, Legendre’a, Bessela, itp.).
Zadanie aproksymacji sprowadza się do znalezienia liczbowych współczynników
, a od strony algebry do rozwiązania układu równań liniowych.
Najczęściej spotykana jest aproksymacja dyskretna, bazująca na danych punktach
,
0,1,...,
i
a
i
=
m
( , ),
1, 2,...,
i
i
x f
i
=
n , do których funkcja p(x) ma się dopasować. Współrzędne
i
x zwane są węzłami, a
i
f - wartościami węzłowymi. Aproksymacja, w której
sumaryczny błąd odejścia funkcji p(x) od zadanych punktów ( , )
i
i
x f podlega
minimalizacji, zwana jest najlepszą aproksymacją. Jeżeli używa się do tego normy
średniokwadratowej (
2
1
n
i
i
x
x
=
=
∑
), to najlepsza aproksymacja zwana jest metodą
najmniejszych kwadratów, co dla jednomianowych funkcji bazowych sprowadza się
do minimalizacji następującego funkcjonału błędu:
0
1
0
1
1
2
0
1
0
1
( , ,...,
)
1
0
...
...
( , ,...,
)
(
)
min
( , ,...,
)
...
...
m
n
m
j
m
j i
i
m
a a
a
i
j
m
a
a
B a a
a
a x
f
B a a
a
a
−
=
=
=
=
=
−
⇒
⇒
=
∑ ∑
8
Można też stosować metodę ważoną polegającą na wprowadzeniu wag
do funkcjonału
,
1, 2,...,
i
w
i
=
n
i
i
1
2
0
1
1
0
( , ,...,
)
(
)
n
m
j
m
j i
i
j
B a a
a
a x
f w
−
=
=
=
−
∑ ∑
regulujących
odejście wielomianu aproksymującego
1
0
( )
m
j
j
i
p x
a x
−
=
=
∑
od oryginalnych wartości
węzłowych
i
f . Im większa waga, tym bliżej będzie przechodzić krzywa danego
punktu.
Szczególnym przypadkiem aproksymacji jest interpolacja, gdzie w ogólności od
krzywej p(x) wymaga się jedynie ścisłego przejścia przez wszystkie dane punkty
( , ),
1, 2,...,
i
i
x f
i
=
n , co można zapisać: ( )
,
1, 2,...,
i
i
p x
f
i
n
=
=
i co implikuje od
razu stopień interpolacji m=n-1. Krzywa interpolacyjna będzie miała wzór:
1
( )
( )
n
i i
i
p x
a
ϕ
=
=
∑
x . Współczynniki ,
1, 2,...,
i
a
i
n
=
znajduje się z następującego
układu równań:
⋅ = F
Φ a
, gdzie:
( ),
,
,
1, 2,...,
ij
j
i
i
i
x
F
f
i j
n
ϕ
Φ =
=
=
. Jeżeli funkcje
bazowe przyjmie się tak, aby spełniały następujący warunek
0,
( )
( ),
( )
1,
i
i
i
j
ij
dla i
j
x
L x
L x
dla i
j
ϕ
δ
≠
⎧
=
=
= ⎨
=
⎩
, to wtedy prawdziwy jest następujący
wzór (tzw. wzór interpolacyjny Lagrange’a):
1
( )
( )
n
i
i
i
p x
f L
=
=
⋅
∑
x
n
n
. Rozwiązywanie
układu równań nie jest wtedy potrzebne. Uogólnieniem interpolacji Lagrange’a jest
tzw. interpolacja l’Hermitte’a, gdzie w węzłach oprócz wartości funkcji mogą być
również dane wartości pochodnych.
Interpolacja wielomianami coraz wyższych stopni (przy coraz większej liczbie danych
punktów) nie jest korzysta ze względu na zachowanie międzywęzłowe funkcji
interpolujących. Dlatego też przy zachowaniu własności interpolacyjnych można
budować tzw. funkcje sklejane (ang. spline) z wielomianów niskich rzędów na
odcinkach wyznaczanych przez kolejne pary sąsiednich węzłów, stosując odpowiednie
kryteria ciągłości dla tych funkcji.
Problem aproksymacji (w tym także interpolacji) funkcji pojawia się zarówno przy
obróbce danych eksperymentalnych jak i numerycznych; często też dokonuje się
aproksymacji funkcji nieznanej np. przy budowie funkcji kształtu w MES lub przy
generacji wzorów różnicowych.
• Numeryczne różniczkowanie i całkowanie
Numeryczne różniczkowanie i całkowanie najczęściej wykonywane przez algorytm
komputerowy numeryczne operacje przy rozwiązywaniu problemów brzegowych.
Numeryczne różniczkowanie polega na znalezieniu na podstawie wartości
dyskretnych funkcji wartości pochodnej w wybranych węzłach. Najprostszym
rozwiązaniem problemu jest więc aproksymacja funkcji i jej zróżniczkowanie, a
następnie obliczenie wartości pochodnej we wskazanym punkcie, ale ze względu na
złożoność takiego rozwiązania stosuje się inne metody generacji tzw. wzorów
różnicowych. Pochodną numeryczną rzędu k przedstawia się jako kombinację liniową
wartości funkcji
i współczynników liczbowych
:
,
1, 2,...,
i
w
i
=
,
1, 2,...,
i
a
i
=
9
( )
1
,
1, 2,...,
(1, 2,..., )
n
k
j
i
i
i
w
a w
k
j
=
=
=
∈
∑
n .
Współczynniki
można znaleźć rozwijając wartości węzłowe w szereg Taylora
wokół węzła centralnego (tego, w którym liczona jest pochodna, czyli (j)), a następnie
przemnożeniu rozwinięć przez odpowiedni współczynnik, dodaniu stronami i
porównaniu z oryginalnym operatorem różniczkowym
i
a
( )
( )
k
i
k
d
w
dx
. Taki sposób noszący
nazwę metody współczynników nieoznaczonych pozwala też na oszacowanie
dokładności wzoru różnicowego poprzez zebranie pierwszych odrzuconych w
rozwinięciach w szereg wyrazów wyższych rzędów. W bezsiatkowej wersji MRS
generacja wzorów różnicowych odbywa się za pomocą techniki metody najmniejszych
ważonych kroczących kwadratów (ang. Moving Weighted Least Squares – MWSL).
Należy również zaznaczyć, iż wzory różnicowe mogą być budowane nie tylko na
samych wartościach funkcji w węzłach, ale też na wartościach pochodnych dowolnych
rzędów. Takie schematy różnicowe noszą nazwę uogólnionych wzorów różnicowych.
Numeryczne całkowanie funkcji może odbywać się na dwa sposoby. Pierwszy,
noszący nazwę wzorów (tzw. kwadratur) Newtona – Cotesa, polega na zastępowaniu
funkcji podcałkowej wielomianami coraz wyższych rzędów. Do tej grupy należą
popularne wzory prostokątów, trapezów i Simpsona:
( )
calka oryginalna
( ) ,
metoda prostokątów
1
[ ( )
( )] ,
metoda trapezów
2
1
[ ( ) 4
(
)
( )] ,
metoda Simpsona
3
2
2
b
a
p
t
S
I
f x dx
I
f a h h b a
I
f a
f b
h h b a
a b
b a
I
f a
f
f b
h h
=
−
⎧
⎪
=
⋅
= −
−
⎪
⎪
=
+
⋅
= −
−
⎨
⎪
+
−
⎪ =
+ ⋅
+
⋅
=
−
⎪⎩
∫
Inną koncepcję prezentują kwadratury Gaussa, które zakładają wartość całki jako
kombinację liniową wartości funkcji podcałkowej f(x) w punktach zwanych punktami
całkowania ,
1, 2,...,
i
x
i
=
N oraz wag ,
1, 2,...,
i
w
i
N
=
:
1
1
2
2
1
( )
( ) ...
( )
( )
N
G
N
N
i
i
i
I
w f x
w f x
w
f x
w f x
=
=
⋅
+
⋅
+ +
⋅
=
⋅
∑
.
Wagi dobierane są według zasady, by wzór całkowania przybliżony był wzorem
ścisłym dla wielomianu możliwe wysokiego stopnia ortogonalnego (z wagą) w
przedziale
[
]
1,1
−
. Dobór odpowiednich wielomianów implikuje rodzaj kwadratur
Gaussa. (np. kwadratury Gaussa - Legendre’a, Gaussa – Laguerra itp.). Miejsca
zerowe i wagi są najczęściej tablicowane. Ilość miejsc zerowych (wag) N świadczy o
punktowości wzoru Gaussa (np. N=2 – wzór 2-punktowy). Przed zastosowaniem
kwadratury Gaussa należy wyjściowy przedział całkowania
[ ]
,
a b
przetransformować
do przedziału
[
]
1,1
−
.
10
W celu podniesienia dokładności wyniku numerycznego całkowania dzieli się
wyjściowy przedział
[ ]
,
a b
na podprzedziały i do każdego z osobna stosuje się
odpowiednio wybrany wzór niskiego rzędu. Następnie wyniki z poszczególnych
podprzedziałów dodaje się do siebie.
Generacja wzorów różnicowych jest podstawą metod różnicowych (np. MRS),
natomiast całkowanie numeryczne występuje niemal w każdej metodzie dyskretnej,
zwłaszcza przy sformułowaniu globalnym, w którym występuje konieczność
całkowania po zadanym obszarze.
• Rozwiązywanie równań różniczkowych
Problemy, w których występują równania różniczkowe (oprócz nich dany jest obszar,
brzeg i odpowiednie warunki) można podzielić na dwie grupy: problemy początkowe
(warunki – początkowe- zadane są w jednym punkcie obszaru na funkcje niewiadomą
i jej kolejne pochodne – klasyczne zagadnienie Cauchy’ego) oraz problemy brzegowe
(warunki – brzegowe – zadane są w kilku punktach obszaru). Rozwiązywanie
problemów brzegowych zostało omówione wyczerpująco w pytaniu 1, tu omówione
zostanie rozwiązywanie numeryczne problemów początkowych. Poniższe metody
opracowane są dla problemów brzegowych rzędu 1-szego, wszystkie wyższe rzędy
można bowiem do niego sprowadzić. Oryginalne zadanie ma postać:
0
0
( , ),
( )
dy
f x y
y x
y
dx
=
=
.
Ogólny zapis poszukiwania wartości dyskretnych funkcji y(x) ma postać wspólną dla
wszystkich metod:
1
1
( , )
n
n
x
n
n
n
x
y
y
f x y dx
y
+
+
=
+
=
+ Δ
∫
n
.
Poszczególne funkcje różnią się między sobą sposobem obliczania przyrostu
n
Δ na
odcinku
[
]
1
,
n
n
x x
+
. I tak: metody jednokrokowe (metoda Eulera, metoda Rungego -
Kutty) bazują przy obliczaniu
n
Δ na jednym punkcie wstecz, podczas gdy metody
wielokrokowe (metoda Adamsa – Bashfortha, metoda Adama – Moultona) bazują na
znajomości kilku punktów wstecz. Najprostsza z nich wszystkich, metoda Eulera,
zakłada stałą postać funkcji podcałkowej na danych odcinku
[
]
1
,
n
n
x x
+
(
):
1
n
n
h x
x
+
=
−
1
( , )
n
n
n
n
y
y
h f x y
+
=
+ ⋅
.
Następne metody aproksymują funkcję podcałkową wielomianami coraz wyższych
rzędów, np. metoda Rungego – Kutty rzędu II:
1
2
1
1
1
( , )
( ,
)
1
(
)
2
n
n
n
n
n
n
K
h f x y
K
h f x y
K
2
y
y
K
K
+
= ⋅
= ⋅
+
=
+
+
lub metoda Adamsa – Bashfortha rzędu II:
1
1
[23
( ) 16
(
) 5
(
)]
12
n
n
n
n
n
h
y
y
f x
f x
f x
+
−
=
+
⋅
− ⋅
+ ⋅
2
−
.
11
Metoda jednokrokowa wysokiego rzędu oraz para metod wielokrokowych (otwarta i
zamknięta, – czyli bazująca na informacjach odpowiednio znanych i nieznanych
a’priori) stanowią aparat do rozwiązywania problemów początkowych zwanych
metodą predyktor - korektor.
3. Idea tworzenia modeli dyskretnych zgodnie z koncepcją MES na przykładzie 1- i
2- wymiarowego układu ciągłego.
Każda metoda dyskretna (w tym także MES) dokonuje dyskretyzacji rozważanego modelu
ciągłego (matematycznego) danego obiektu rzeczywistego. Dyskretyzacji (zamianie problemu
ciągłego na dyskretny – punktowy) podlegają: obszar (pręt, płyta, tarza itp.), funkcja
niewiadoma (poprzez wprowadzenie tzw. wartości węzłowych) oraz warunki brzegowe.
Podstawową cechą MES jest podział obszaru na elementy – elementy skończone. Są to figury
geometryczne, w których skład wchodzą tzw. funkcje kształtu oraz stopnie swobody (w
węzłach – wierzchołkach elementu). Stopnie swobody to najczęściej wartości węzłowe
funkcji niewiadomej, ale mogą one być również dowolną kombinacją jej pochodnych.
Funkcje kształtu to zestaw funkcji liniowo niezależnych przyporządkowanych elementowi w
taki sposób, że ich wartości w węzłach (wierzchołkach) elementu są równe wartościom
węzłowym. Ogólnie funkcje kształtu (zazwyczaj wielomiany Lagrange lub l’Hermitte’a – ale
też np. funkcje trygonometryczne lub funkcje spline) spełniają warunki reproduktywności
(
( )
,
,
1, 2,..., ,
liczba węzlów
i
j
ij
N x
i j
N
N
δ
=
=
=
) oraz konsystentności (
1
( ) 1
N
i
i
N x
=
=
∑
).
Elementy nie mogą na siebie zachodzić, węzły jednego elementu nie mogą znajdować się
wewnątrz drugiego ani na żadnym z jego boków – węzły mogą jedynie pokrywać się w
wierzchołkach. Suma elementów powinna jak najlepiej pokrywać się ze kształtem obszaru.
Topologię układu m MES opisuje się za pomocą tzw. incydencji, czyli zależności między
elementem, jego węzłami i numerami stopni swobody w węzłach. Każdemu elementowi
przypisuje się wiele wielkości wektorowych i macierzowych. Jedną z nich jest macierz
transformacji, informująca o położeniu elementu w stosunku do globalnego układu
współrzędnych.
Podstawowe znaczenie dla dalszej analizy dyskretnej MES ma wybór sformułowania
problemu oraz model numeryczny. Najczęściej stosuje się sformułowanie słabe,
rozpoczynając rozważania od zasady prac wirtualnych lub globalne, gdzie podstawą analizy
jest odpowiedni funkcjonał energetyczny. Model MES może być np. w przypadku zadań
statyki przemieszczeniowy, naprężeniowy lub mieszany - hybrydowy.
W MES postępowanie mające na celu obliczenie globalnych wielkości i uzyskanie
niewiadomych pierwotnych (np. przemieszczeń) i wtórnych (np. naprężeń i odkształceń)
nazywa się syntezą natomiast przeniesienie tych informacji do układów lokalnych nazywa się
analizą.
W metodzie elementów skończonych można wyróżnić następujące etapy analizy konstrukcji:
analiza na poziomie punktu P, przekroju S, ciała C oraz układu U. Synteza zwie się w MES
powrotem do elementu.
Przystępując do rozwiązania zadania za pomocą MES dokonujemy kolejno:
• Dyskretyzacji układu: wprowadzenie węzłów, podział na elementy, wprowadzenie
stopni swobody w węzłach, opis topologii układu (incydencja!),
• Obliczenia wielkości macierzowych i wektorowych na poziomie elementu (w
układach lokalnych): macierzy sztywności, macierzy transformacji (przejścia od
układu lokalnego do globalnego), ew. macierzy mas (dla zadań dynamiki) i macierzy
wstępnych naprężeń (dla zadań stateczności), wektora zastępników elementowych
oraz wektora wstępnych reakcji.
12
• Transformacja wielkości elementowych do ew. układów węzłowych i układu
globalnego,
• Agregacja elementów w układ (złożenie wielkości elementowych zapisanych w
układzie globalnym z pominięciem działania na wielkościach zerowych),
• Budowa układu równań
: K – globalna macierz sztywności, q – wektor
niewiadomych pierwotnych, P – wektor obciążeń węzłowych, R – wektor reakcji
węzłów (niewiadome wtórne),
K q P R
⋅ = +
• Uwzględnienie podstawowych i naturalnych warunków brzegowych,
• Rozwiązanie układu równań, otrzymanie niewiadomych pierwotnych i wtórnych,
• Powrót do elementu, obliczenie wielkości elementowych (reakcji, sił przekrojowych)
na podstawie rozwiązań globalnych (przemieszczeń, reakcji węzłów),
• Ew. proces adaptacyjny związany z obliczeniem poziomu błędu i przyjęcia nowej
siatki elementów.
Przykład 1D – kratownica
Dyskretyzacja (numeracja węzłów, elementów, globalnych stopni swobody, relacje
przylegania elementów do węzłów):
Y,V
X,U
P = 8N
3
1
2
Q
6
Q
5
3
Q
4
Q
3
2
Q
2
Q
1
1
Incydencja:
• element nr 1 - węzły nr 1 i 2 – stopnie swobody 1,2,3,4,
• element nr 2 - węzły nr 2 i 3 – stopnie swobody 3,4,5,6,
• element nr 3 - węzły nr 1 i 3 – stopnie swobody 1,2,5,6.
Dodatkowo określa się położenie elementów w stosunku do globalnego układu
współrzędnych poprzez podanie ich dostaw kierunkowych (sinus i cosinus). Warunki
brzegowe
Przyjęto elementy prętowe - kratowe (2 węzły, 2 stopnie swobody w elemencie, po jednym
na każdy węzeł – związane z poziomym ruchem sztywnym węzłów).
Uwzględnienie siły skupionej w węźle nr 3 odbywa się na poziomie globalnego układu
równań w wektorze obciążeń węzłowych. Uwzględnienie warunków brzegowych implikuje
zerowe przemieszczenia:
i przebudowę układu równań (skreślenie trzech
kolumn i wierszy odpowiadających tym niewiadomym).
1
2
4
0
Q
Q
Q
=
=
=
13
Przykład 2D – tarcza
Tarczę dyskretyzowano elementami tarczowymi (w płaskim stanie naprężenia) CST
(Constant Strain Triangle – element stałego odkształcenia) trójwęzłowymi. Stopniami
swobody w węzłach są ich przemieszczenia poziome i pionowe. Każdy węzeł na 2 stopnie
swobody. Element ma ich 6. Układ ma
2 N
⋅
stopni swobody, gdzie
oznacza liczbę
węzłów.
N
Ponieważ do tarczy przyłożone jest obciążenie ciągłe, będzie ono musiało zostać zamienione
na tzw. zastępniki obciążenia elementowego – siły węzłowe wywołujące ten sam efekt
kinematyczny, co oryginalne obciążenie. Inna nazwa: równoważniki (zastępniki
przeskalowane przez ‘-‘).
q [kN/m]
i
1
2
3
Y,V
X,U
Uwzględnienie warunków brzegowych przebiega podobnie jak w przypadku 1.
Na podstawie obliczonych przemieszczeń węzłów można aproksymować pole przemieszczeń
dla całej tarczy lub poszczególnych elementów za pomocą wprowadzonej wcześniej bazy
wielomianowej funkcji kształtu (funkcje liniowe – opisują przemieszczenia, ich pochodne to
funkcje stałe a mają charakter odkształceń – stąd nazwa: elementy stałego odkształcenia).
Można też stosować elementy wyższego rzędu (np. elementy trójkątne LST – 6-węzłowe
elementy liniowego odkształcenia) lub elementy tarczowe innego kształtu geometrycznego
(najczęściej obok trójkątnych: prostokątne).
4. Zasadnicze podobieństwa i różnice w koncepcjach MES i MRS.
Różnice i podobieństwa zostaną przedstawione na podstawie cech charakterystycznych
obydwu metod dyskretnych. MRS – metoda różnic skończonych jest najstarszą metodą
14
dyskretną do rozwiązywania zadań brzegowych, pochodzi jeszcze sprzed ery komputerów –
była to wersja na siatkach regularnych. Wraz z rozwojem komputeryzacji powstała jej wersja
na siatki dowolnie nieregularne (tzw. bezsiatkowa MRS) należąca do szerokiej grupy metod
bezsiatkowych. MES natomiast jest najbardziej rozpowszechnioną obecnie metodą
komputerową, mającą bardzo szerokie zastosowanie w problemach technicznych.
Sformułowanie
Obydwie metody służą do rozwiązywania problemów brzegowych mechaniki. Mogą one być
różnorako sformułowane. MRS może być równorzędnie użyta w każdym z nich (lokalne,
globalne, słabe), ale najlepiej znana jest jej wersja w sformułowaniu lokalnym (układ równań
różnicowych). Natomiast MES najlepiej pracuje w sformułowaniach wymagających
całkowania po obszarze – przy zasadzie wariacyjnej lub w zagadnieniach minimalizacji
funkcjonału.
Dyskretyzacja
W MRS wszystkie wielkości dyskretne sprowadza się do węzłów. W wyniku działania MRS
otrzymuje się wyniki węzłowe – dopiero ich późniejsza aproksymacja, niezależna od
wcześniejszych zabiegów – daje rozwiązanie w całym obszarze. Węzły w wersji klasycznej
MRS musiały być rozłożone regularnie, w równych odstępach – tworzyły tzw. siatkę
różnicową, natomiast w wersji bezsiatkowej – węzły mogą być rozłożone dowolnie – bez
żadnych dodatkowych warunków typu regularność czy element. Taka cecha pozwala na łatwe
dokładanie węzłów, przesuwanie oraz usuwanie. Topologia siatki w MRS opisywana jest za
pomocą sąsiadów węzłów: mocnych i słabych. Określenia te pochodzą z podziału obszaru na
podobszary przypisywane węzłom (tzw. wielokąty Voronoi). Na tej podstawie dokonywana
też jest triangularyzacja obszaru. W sformułowaniach globalnych potrzebny jest dodatkowy
podział obszaru na podobszary do całkowania, co jest robione na różne sposoby.
Topologię w MES stanowią węzły i elementy razem budujące siatkę MES. To wybór rodzaju
elementu i aproksymacji między jego węzłami decyduje o końcowym wyniku węzłowym.
Elementy skończone stanowią wielką bibliotekę w zależności od potrzeb można wybierać
elementy pod względem kształtu, liczby węzłów, liczby i rodzaju stopni swobody, rzędu i
jakości funkcji kształtu i innych cech. Do generacji siatek w MES stosuje się specjalne
programy zwane generatorami siatek.
Aproksymacja funkcji niewiadomej
W MRS podstawą metody jest zamiana operatorów różniczkowych na operatory różnicowe,
czyli obliczające wartość danej pochodnej na podstawie informacji z otaczających dany węzeł
innych węzłów. Razem ta konfiguracja węzłów wraz z węzłem centralnym nosi nazwę
gwiazdy różnicowej. Są różne kryteria doboru węzłów do gwiazd, np. kryterium
odległościowe – należy przy tym pamiętać, że zły dobór węzłów spowoduje, iż schemat
różnicowy może okazać się osobliwy lub źle określony. W klasycznej MRS (na siatkach
regularnych) schematy różnicowe generuje się najczęściej metodą współczynników
nieoznaczonych. W bezsiatkowej MRS do generacji całego kompletu wzorów różnicowych
od pierwszej pochodnej aż do ostatniego wymaganego rzędu służy metoda zwana Metodą
Najmniejszych Ważonych Kroczących Kwadratów (ang. Moving Weighted Least Squares –
MWSL). Posługuje się ona lokalną aproksymacją funkcji, dla której różnice między
wartościami węzłowymi są regulowane przez funkcje wagowe. Wpływ danego węzła zanika
wraz z jego odległością od węzła centralnego. Wagi osobliwe wprowadzają interpolację.
Ogólnie aproksymację MWLS w gwieździe z węzłem centralnym
0
x opisuje zależność:
15
1
0
0
0
1
( , )
( )
( ) ( )
m
i
t
i
i
u x x
a x
x
a x p x
−
=
=
⋅
=
∑
Współczynniki ,
1, 2,...,
i
a
i
m
=
(m – rząd lokalnej aproksymacji) znajduje się z układu
równań (n – liczba węzłów w gwieździe) :
, gdzie:
2
t
t
P W Pa P W f
=
2
P – macierz wartości interpelantów (jednomiany) w każdym z węzłów,
W – diagonalna macierz wagowa,
f – wektor wartości węzłowych,
a – wektor niewiadomych współczynników,
2
1
1
1
2
2
2
2
2
(
)
3
3
3
2
1
1
1
...
2
!
1
1
1
...
2
!
1
1
1
...
2
!
... ...
...
...
...
1
1
1
...
2
!
m
m
m
n m
m
n
n
n
h
h
h
m
h
h
h
m
P
h
h
h
m
h
h
h
m
×
⎡
⎤
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
= ⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎣
⎦
,
1
1
2
2
0
(
)
( 1)
0
0
0
0
0
0
,
,
0
0
...
0
...
0
0
0
i
i
n n
n
n
n
w
f
w
f
W
f
w
f
×
×
⎡
⎤
⎡ ⎤
⎢
⎥
⎢ ⎥
⎢
⎥
⎢ ⎥ h x x
=
=
=
⎢
⎥
⎢ ⎥
⎢
⎥
⎢ ⎥
⎣
⎦
⎣ ⎦
−
,
Wtedy współczynniki można wyrazić:
, natomiast aproksymację w
gwieździe:
. Funkcje
noszą nazwę pseudo funkcji
kształtu, dlatego iż nie odtwarzają przy dowolnych wagach wartości węzłowych (
2
1
2
(
)
t
t
a
P W P
P W f
−
=
( )
( )
t
u x
a p N x f
=
=
⋅
%
( ),
1, 2,...,
i
N x
i
m
=
%
( )
i
i
N x
f
≠
%
).
Ich dobór jest podstawą klasyfikacji wszystkich metod bezsiatkowych.
W MES funkcje kształtu zapewniają ciągłą aproksymację w elemencie, w jego wierzchołkach
odtwarzają wartości węzłowe funkcji niewiadomej. Za pomocą funkcji kształtu i ich
pochodnych można opisać dowolne pole fizyczne w elemencie: ( )
( )
u x
N x f
=
⋅ . Dzięki temu
wszystkie rozważania sprowadza się do węzłów, a między nimi (wewnątrz elementu)
wszystko zależy od przyjętej aproksymacji.
Układ równań algebraicznych
W MRS rzadki lub pasmowy, może być symetryczny. Kłopotliwe jest uwzględnianie
warunków brzegowych ze względu na postać funkcji bazowych (ilorazy wielomianów).
Układ równań jest otrzymywany zależnie od sformułowania zagadnienia metodą kolokacji,
Bubnowa – Galerkina albo Rayleigha – Ritza. Dokładne różniczkowanie i całkowanie funkcji
bazowych jest w tej metodzie czasochłonne.
W MES układ jest na ogół pasmowy i symetryczny, uzyskiwany metodą Bubnowa –
Galerkina albo Rayleigha – Ritza. Stosunkowo łatwo oblicza się współczynniki układu
charakterystyczną dla MES metodą agregacji, często połączoną z eliminacją Gaussa (metoda
frontalna). Łatwo wymusza się warunki brzegowe Dirichleta.
Opracowanie wyników (tzw. postprocessing)
W MRS sprowadza się ono do aproksymacji (uciąglenia) danych dyskretnych w węzłach (co
czyni się najczęściej metodami najmniejszych kwadratów – w tym MWLS, gdzie przyjęcie
nieosobliwych wag powoduje dodatkowe wygładzenie wyniku) lub zróżniczkowania
numerycznego wyniku (np. w celu obliczenia odkształceń należy zróżniczkować
16
przemieszczenia węzłowe – tu obserwowana jest tzw. nadzbieżność – lepsza dokładność
pochodnej numerycznej niż samej funkcji wyjściowej).
W MRS gotowy aparat aproksymacji elementowej pozwala na szybkie znajdywanie wartości
funkcji w dowolnym punkcie obszaru.
W razie, gdy konieczna jest dalsza analiza zadania, np. w wyniku adaptacji (polepszenia
dyskretyzacji w miejscu występowania największego błędu) charakterystyczna dla MRS jest
adaptacja typu h, czyli zagęszczenie siatki węzłów, podczas gdy w MES najczęściej podnosi
się stopień lokalnej (elementowej) aproksymacji (adaptacja typu p), gdyż dokładanie nowych
węzłów i elementów wiąże się z kłopotliwą przebudową wyjściowej siatki.
W MES są dobrze rozwinięte metody szacowania błędu. W MRS te techniki są w ciągłym
rozwoju.
5. Klasyfikacja źródeł błędów w metodach komputerowych
Każda metoda komputerowa zakłada budowę ciągu modeli poprzez przyjmowanie różnych
założeń upraszczających pozwalających na zamianę rzeczywistego obiektu na model
komputerowy. Poszczególne etapy modelowania wiążą się więc z budową kolejnych modeli.
W ciągu tego procesu powstają błędy charakterystyczne dla każdego modelu z osobona,
związane z coraz większym upraszczaniem sytuacji wyjściowej dla danego modelu.
Poszczególne etapy modelowania:
• Punktem wyjścia jest OBIEKT RZECZYWISTY ( istniejący bądź projektowany )
będący przedmiotem analizy statycznej lub dynamicznej albo też pewien PROCES
FIZYCZNY odbywający się w określonym środowisku.,
• MODEL MECHANICZNY będący idealizacją rzeczywistego obiektu poprzez
przyjęcie uzasadnionych założeń upraszczających ( dot. materiału, geometrii,
obciążenia, przyszłych warunków eksploatacji), podczas jego tworzenia popełnia się
błąd modelu zwany nieuniknionym; może się wiązać ze złym przyjęciem parametrów
do modelu,
• MODEL MATEMATYCZNY budowany na podstawie modelu mechanicznego,
ujmujący powyższe idealizowane cechy w opisie matematycznym, na ogół w postaci
układu różniczkowych równań cząstkowych lub też funkcjonału albo zasady
wariacyjnej,
• MODEL DYSKRETNY, dyktowany możliwościami obliczeniowymi komputera (
proste operacje algebraiczne i logiczne), polega na przetworzeniu modelu
matematycznego w układ równań algebraicznych; przy jego tworzeniu niezwykle
ważny jest wybór właściwej metody dyskretnej do analizy problemów brzegowych;
• MODEL NUMERYCZNY, wiążący się z właściwym wyborem nowoczesnych
technik numerycznych, służących rozwiązywaniu dużych układów równań
algebraicznych liniowych i nieliniowych, uogólnionych zagadnień własnych, a także
zagadnień programowania liniowego i nieliniowego. Błędy numeryczne to największa
zmora dzisiejszych modeli komputerowych. Najczęściej wiążą się one z tzw. błędami
obcięcia (obcięcie określonej liczby wyrazów szeregu przy rozwinięciu danej funkcji)
i błędami zaokrągleń. Można też mówić o błędzie residuum (kontrola spełnienia przez
dane rozwiązanie wyjściowego równania – zapobiega przed zbieżnością wyniku do
17
liczby niebędącej rozwiązaniem) lub błędzie zbieżności (kontrola zbieżności pozwala
uniknąć rozbieżności procesu – wynik rozbiega się do nieskończoności lub
niestabilności procesu).
• MODEL KOMPUTEROWY, związany z wyborem odpowiedniego sprzętu
komputerowego do obliczeń oraz napisaniu programu w jednym ze znanych języków (
tzw. modelowanie otwarte) w celu kontrolowania implementowanych procedur
numerycznych lub wyborze dostępnego na rynku programu komercyjnego,
traktowanego na ogół przez modelującego jako „czarna skrzynka” ( tzw. modelowanie
zamknięte).
• PRZYGOTOWANIE DANYCH WEJŚCIOWYCH I WYKONANIE OBLICZEŃ
• WERYFIKACJA POSZCZEGÓLNYCH MODELI NA PODSTAWIE
OTRZYMANYCH WYNIKÓW, bardzo ważny etap, który każe patrzeć krytycznie na
otrzymane wyniki z komputera, które mogą być całkowicie błędne ( „śmieci na
wejściu – śmieci na wyjściu”) lub niekiedy nie spełniają wcześniejszych założeń w
poszczególnych modelach. Analiza błędu zakłada istnienie dwóch rodzajów błędów:
a’priori („przed faktem”) i a’posteriori („po fakcie”). Błąd a’priori liczony na
podstawie teorii danej metody pozwala na oszacowanie dokładności wyniku przed
rozwiązaniem zadania (np. w MES – można określić dla zadań liniowych błąd na
podstawie rozmiaru siatki oraz stopnia aproksymacji). Błąd a’posteriori liczony po
rozwiązaniu danego zadania może dotyczyć samego rozwiązania (kontrolna
zbieżności rozwiązania w węzłach wspólnych dla dwóch siatek) lub tzw. residuum,
czyli dokładności spełnienia wyjściowego równania różniczkowego w wybranych
punktach obszaru. Duży poziom błędu residualnego jest najczęściej podstawą do
ponownego rozwiązania tego samego zadania, ale dla innej, lepszej dyskretyzacji.
Należy też zaznaczyć, iż w miarę możliwości obliczenia, zwłaszcza w przypadku
zadań skomplikowanych, należy zawsze weryfikować przez porównanie wyników np.
z różnych metod dyskretnych lub pochodzących z różnych programów
komputerowych.
Przykład : zamodelować ściskany wspornik i przy pomocy MRS obliczyć siłę krytyczną
Eulera.
0
1
f
w
f
w
1
w
0
18
• Etap I : Model mechaniczny wspornika : belka ( ściskanie + zginanie ) zamocowana
jednostronnie, założenia upraszczające :
1. Geometria : pręt o długości L ( ciało 1D ) – rozważania sprowadzone do osi
obojętnej,
2. Materiał : ograniczenie do zakresu sprężystego ( moduł Younga E ),
3. Przekrój : pręt lity, znane wymiary ( charakterystyki : A,I ),
4. Zastosowano liniową teorię sprężystości ( liniowe związki fizyczne i
kinematyczne ), stąd założenie małych przemieszczeń,
5. Obciążenie : rezygnacja z zasady zesztywnienia aby uwzględnić efekt
wyboczenia : warunki równowagi rozpatruje się w konfiguracji odkształconej,
obciążenie jest zachowawcze ( nie zmienia kierunku podczas wyboczenia )
• Etap II : Model matematyczny : zdecydowano się na model matematyczny w
sformułowaniu lokalnym, w postaci równania różniczkowego razem z warunkami
brzegowymi (układ współrzędnych kartezjańskich zaczepiono we wsporniku );
równanie uwzględnia zasadę małych przekrojów i małych przemieszczeń z modelu
mechanicznego :
sformułowanie problemu (stateczność zlinearyzowana) :
l
x
EJ
w
x
M
w
x
f
w
x
f
x
w
dx
d
x
≤
≤
−
=
=
=
0
)
,
(
)
,
(
,
)
,
(
)
(
)
(
2
2
L
gdzie (wzór na moment zginający uwzględnia odkształcenie pręta) :
:
0
)
0
(
'
0
)
0
(
)
(
)
,
(
=
=
−
−
=
w
w
w
f
P
w
x
M
f oznacza strzałkę ugięcia pręta ( dla x = L ).
• Etap III : Model dyskretny : zdecydowano się użyć Metodę Różnic Skończonych,
zdyskretyzowano pręt dwoma węzłami („1” – x = 0, „2” – x = L ), warunki brzegowe i
nieznaną funkcję przemieszczeń punktów pręta zapisano za pomocą przemieszczeń
węzłowych ( w
1
= 0 i w
1
’ = 0 – warunek na pochodną zapisano za pomocą fikcyjnego
węzła „0” – x = -L : w
0
= w
2
, niewiadoma : w
2
= ? ), równanie różniczkowe
zamieniono na różnicowe ( operator na drugą pochodną zastąpiono operatorem
różnicowym centralnym )
w’’ = (1/h^2) * ( w(i-1) – 2*w(i) + w(i+1) )
• Etap IV : Obliczenia : kollokacja w węźle „1” :
( w0 – 2*w1 + w2 ) / (1/L^2) + ( P / EJ )*w1 = P*f / EJ
uwzględnienie warunków podparcia oraz w2 = f prowadzi do równania na siłę
krytyczną :
2*f / (L^2) = P*f / EJ
19
a
stąd :
P = 2*EJ / (L^2)
• Etap V : Weryfikacja wyników :
Analiza błędu : przez porównanie ze znanym wynikiem analitycznym :
Pe= (1/4) * (
π^2) * EJ / (L^2) = 2.467*EJ / (L^2)
E = | Pe – P | / | Pe | = 0.189
Uwaga! Wynik analityczny jest ścisły tylko dla przyjętego wyżej modelu
mechanicznego!
Źródło błędu w modelu mechanicznym : model jest wystarczający dla potrzeb
inżynierskich, wynik dokładniejszy byłby po przedsięwzięciu środków zaradczych :
uwzględnienie dużych ugięć, rozważanie wyboczenia plastycznego ( hipoteza
Karmana ), uwzględnienie rzeczywistej geometrii pręta ( wyboczenie pręta
cienkościennego – hipotezy Własowa ), przeprowadzenie doświadczenia : zbadanie
siły krytycznej.
Źródło błędu w modelu matematycznym : nieuwzględnienie wpływu ściskania
(podwyższenie siły krytycznej) oraz ścinania w okolicach podpory ( hipoteza
Timoszenki ), środki zaradcze : uwzględnienie w równaniu ugięcia również skrócenia
od ściskania ( 1 / ( E*J*L ).
Źródło błędu w modelu dyskretnym : zbyt mała liczba węzłów, niski rząd operatora
różnicowego ( w obszarze i na brzegu ), środki zaradcze : zwiększyć liczbę węzłów (
przy trzech węzłach na belce rozwiązanie wynosi P = 0.95*Pe ), przeprowadzić
analizę adaptacyjną ( typu h, p i mieszane ), zwiększyć rząd operatorów ( np. poprzez
technikę aproksymacji wyższego stopnia ), zmienić metodę ( np. MES – stateczność
zlinearyzowana - porównać wyniki ), odwołać się do danych z eksperymentu.
6. Zachowanie się układów prętowych przy obciążeniach termicznych i
geometrycznych.
Deformację konstrukcji prętowej mogą spowodować różnego rodzaju czynniki, które noszą
nazwę obciążeń zewnętrznych. Wyróżnia się ich trzy zasadnicze rodzaje:
• Obciążenia statyczne (indeks p) czyli wszelkiego rodzaju siły uogólnione,
• Obciążenia termiczne (indeks t) w postaci równomiernego i nierównomiernego
ogrzania elementów konstrukcyjnych,
• Obciążenia geometryczne (indeks g) pochodzące od wszelkiego rodzaju
wymuszeń kinematycznych (np. osiadania podpór, niedokładności montażu
itp.).
Wpływ obciążeń statycznych na pracę konstrukcji zostanie pominięty, gdyż pytanie nie
dotyczy tego problemu bezpośrednio.
20
Obciążenia termiczne
Obciążenia termiczne ( j t
= ) powodują odkształcenia elementarnego wycinka pręta o
długości
, spowodowane zmianą temperatury. Wpływ temperatury może się uwydatnić w
zmianie wartości trzech sił przekrojowych na odcinku
pręta: momentu gnącego M, siły
podłużnej N oraz siły poprzecznej Q. Przyrosty powyższych wielkości wynoszą odpowiednio:
ds
ds
0
,
,
M
N
s
w
t
t
t
t
ds
t ds
h
α
α
Δ − Δ
Δ =
Δ = ⋅ Δ ⋅
Δ = 0.
Q
t
Zmianie mogą ulec, co najwyżej wartości przyrostów momentu zginającego i siły podłużnej
natomiast temperatura nie ma wpływu na wartość zmiany siły poprzecznej – a raczej ten
wpływ jest pomijalny w przyjętym poprzednio modelu prętowym (ciało 1D). Przyrost
momentu zginającego związany jest z nierównomiernym ogrzaniem dwóch stron pręta – inna
temperatura jest po stronie tzw. „spodów” (indeks „s”) a inna po stronie „wierzchów” (indeks
„w”). Rozróżnienie „wierzchów” i „spodów” jest zarówno w mechanice budowli jak i
wytrzymałości materiałów założeniem umownym i związane jest z przyjęciem globalnego
układu współrzędnych – pozwala na rysowanie dodatnich wartości momentu zginającego po
stronie włókien rozciąganych – „spodów”. Podsumowując: nierównomierne ogrzanie
„spodów” i „wierzchów” powoduje powstanie dodatkowego niezerowego momentu na belce i
faktyczne jej ugięcie – w kierunku temperatury wyższej. Natomiast ogrzanie osi obojętnej
pręta powoduje jego wydłużenie (temperatura dodatnia) lub skrócenie (temperatura ujemna).
Efekty te pojawiają się niezależnie od stopnia statycznej niewyznaczalności układu.
Zakłada się, że przyrosty temperatur
t
Δ
zmieniają się liniowo wzdłuż wysokości przekroju
poprzecznego. Ich wartości należy obliczać ze wzorów:
h
0
0
,
,
s
s
m
m
w
m
m
t
t
t
t
t
t
t
t
t
Δ = −
Δ = −
Δ = − ,
w których przez oznaczono tzw. temperaturę montażu konstrukcji, a przez
-
aktualne temperatury włókien odpowiednio dolnych (spody), górnych (wierzchy) i
położonych w osi obojętnej przekroju poprzecznego. Współczynnik
m
t
0
, i
s
w
t t
t
α jest współczynnikiem
rozszerzalności termicznej (ma wymiar
0
1
C
).
Drugim efektem statycznym oprócz powstania naprężeń w przekroju pręta jest wywołanie
przemieszczeń punktów pręta pochodzenia termicznego – nie mylić z odkształceniami.
Przekrój może się odkształcić tylko wtedy, gdy się przemieści, – ale może się przemieścić bez
powstania odkształceń – taka sytuacja ma miejsce w mechanizmach i była omawiana na
przedmiocie „mechanika ogólna” – w postaci ruchu sztywnego. W przypadku działania
obciążeń termicznych tzw. dopełniająca praca wirtualna (zasada prac wirtualnych w
mechanice budowli jest podstawą do wyprowadzenia wszelkiego rodzaju wzorów na
przemieszczenia) uogólnionych sił przekrojowych (
,
,
i
i
i
)
M N Q wynosi dla całego układu:
0
( )
[
( )
(
t
s
w
w
i
u
t
t
L
M s
t
h
δ
α
Δ − Δ
= −
+ Δ
∫
)] .
i
N s ds
Obciążenia geometryczne
Przy działaniu obciążeń geometrycznych ( j g
=
)
w układach statycznie wyznaczalnych nie
powstają żadne odkształcenia, co bynajmniej nie wyklucza występowania przemieszczeń
spowodowanych przez te obciążenia. Wynika to z faktu, że układy takie nie są
przesztywnione i wszelkie wymuszenia kinematyczne sprawiają, że zachowują się one tak jak
21
mechanizmy ruchome, tzn. przemieszczają się bez występowania jakichkolwiek odkształceń
ich ogniw (elementów). Można więc w rozważanym przypadku (dla USW – układów
statycznie wyznaczalnych) zapisać:
0
0,
M
N
Q
g
g
g
Δ =
Δ =
Δ = 0.
Stąd również praca wirtualna odpowiadająca temu typowi obciążeń wynosi:
0.
g
w
L
δ
=
Natomiast w układach statycznie niewyznaczalnych (USN) każda pojedyncza osiadająca
podpora (osiadanie odbywa się po kierunku nieruchomego przemieszczenia uogólnionego)
zwalnia jeden więz kinematyczny, co nie musi koniecznie oznaczać przejścia konstrukcji w
mechanizm. Stąd ogólnie dla USN obciążenia geometryczne wywołują niezerowe
przemieszczenia, naprężenia i odkształcenia.
Praktyczne uwzględnienie obciążeń geometrycznych (dla USW i USN) odbywa się poprzez
końcowy efekt zastosowania zasady prac wirtualnych dla całego układu będący wzorem
pozwalającym obliczyć przemieszczenie dowolnego węzła (tzw. wzór Maxwella – Mohra):
0
1
( )
( )
( )
( )
( )
( )
( )
[
(
n
i
p
i
p
i
p
s
w
i
i
i
k
u
M s M s
N s N s
k Q s Q s
t
t
p
M
EI
EA
GA
h
α
α
=
⋅
⋅
⋅
⋅
Δ − Δ
=
+
+
+
+ Δ
−
∑
∫
)
( )]
ik
k
s
t N s ds
R
Δ
Wzór ten określa wartość poszukiwanego rzeczywistego
i
p o właściwym dla niego wymiarze
fizycznym tylko wtedy, gdy w i-tym stanie pomocniczym (wirtualnym) posłużymy się
odpowiednim wirtualnym, jednostkowym i bezwymiarowym obciążeniem uogólnionym
1
i
P
= . Wzór uwzględnia wpływy od wszystkich trzech rodzajów obciążeń. Ostatni składnik
wzoru – już poza znakiem całki – jest liczbowym uwzględnieniem wielkości osiadań
kinematycznych.
Przykłady
• Rama portalowa z obciążeniem termicznym
A
B
C
D
1
w
t
t
=
2
s
t
t
=
2
t
2
t
Rama z obciążeniem termicznym
(tzw. stan j=t)
A = A’
B
C
D
Rysunek odkształconej ramy
C’
B’
D’
22
Rama, pokazana na rys., poddana jest obciążeniu termicznemu scharakteryzowanemu przez
łókien skrajnych:
C ,
- współczynnik rozszerzalności termicznej:
następujące wartości:
- temperatura montażu
0
5
m
t
C
= −
,
- aktualna temperatura w
0
0
1
2
5
,
15
t
C
t
= +
= +
5
0
1
10
C
α
−
=
.
Wysokość przekroju prostokątnego słupów i rygla ram
0.2
h
m
=
y wynosi
.
amy:
m równym
i przyrostem
.
- pręt BC jest
ie ogrzany – wygięty i w dłużony z następującymi
Rozważane będą obciążenia termiczne każdego z trzech prętów r
- pręty AB oraz DC są równomiernie ogrzane z wydłużeniem termiczny
0
N
t
t
ε
α
= ⋅ Δ
0
0
2
15 ( 5) 20
m
t
t
t
C
Δ = − =
− − =
nierównomiern
y
odkształceniami:
0
oraz
M
N
s
w
t
t
t
t
t
h
κ
α
ε
α
Δ − Δ
= ⋅
= ⋅ Δ ,
związanymi z wartościami przyrostów:
C
= − = − =
− − =
Δ = − = − = − − =
Δ = Δ − Δ
=
−
+
−
=
− −
+
− −
=
+
=
Na rysunku przedstawiono też szkic ramy odkształconej i przemieszczonej.
• Rama portalowa z obciążeniem geometrycznym
y przemieszczonej wydać charakterystyczną cechę statycznie wyznaczalnych
0
2
0
1
0
0
1
2
15 ( 5) 20
5 ( 5) 10
(
) / 2 [(
) (
)]/ 2 {[5 ( 5)] [15 ( 5)]}/ 2 (10 20) / 2 15
s
s
m
m
w
w
m
m
w
s
m
m
t
t
t
t
t
C
t
t
t
t
t
C
t
t
t
t
t
t
t
Δ
Na rysunku ram
układów poddanych obciążeniu geometrycznemu, a mianowicie proste, niezakrzywione i
niewydłużone poszczególne pręty ramy.
A
B
C
D
A
B
C
D
C’
B’
D’
Rama z obciążeniem geometrycznym
(tzw. stan j=g)
Rysunek przemieszczonej ramy
(1 )
A
A
v
Δ
=
(2 )
A
A
ϕ
Δ
=
(3 )
D
D
v
Δ
=
A’
23
7. Drgania własne i wymuszone, rezonans i tłumienie + Klasyfikacja wpływów
zagadnieniach statyki (bądź też stateczności) budowli zakłada się, że obciążenia są
rakter losowy. W
ń dynamicznych należy wymienić przede wszystkim: takie, które
) okresowe,
ednoznaczne opisanie ruchu konstrukcji wymaga znajomości nie tylko jej własności
dynamicznych, charakterystyki dynamiczne konstrukcji, uwzględnienie
wpływów dynamicznych zagadnieniach konstrukcjach budowlanych
W
przykładane do konstrukcji w sposób statyczny, – czyli energia kinetyczna obciążanego
układu jest pomijalnie mała w porównaniu z jego energią potencjalną. W zadaniach, gdzie jej
wpływ nie może być pominięty, należy analizować konstrukcję metodami dynamiki budowli.
Różnorodność obciążeń zmieniających się w czasie, jakie działają na konstrukcje budowlane
w trakcie ich eksploatacji, jest ogromna. Zmianom mogą ulegać wartości tych obciążeń, ich
zwroty kierunki, także punkty przyłożenia. Każda zmiana, jeżeli tylko jest dostatecznie duża i
szybka, wywołuje tzw. efekty dynamiczne związane z ruchem konstrukcji.
Działające na konstrukcję obciążenia mają w większości przypadków cha
dalszych rozważaniach przedstawione zostaną obciążenia deterministyczne o ściśle
określonych parametrach.
Z najważniejszych obciąże
występują w procesach technologicznych i pochodzą od pracy różnych maszyn i urządzeń,
działanie wiatru na budowle wysokie, ruch pojazdów po mostach i dźwigów po torach
podsuwnicowych oraz trzęsienia Ziemi, efekty parasejsmiczne pochodzące od ruchu
drogowego szybkich i ciężkich pojazdów kołowych, odstrzały w kamieniołomach i tzw.
tąpnięcia górotworu na skutek eksploatacji górniczej prowadzonej w kopalniach.
Na wykresie poniżej pokazano typowe obciążenia dynamiczne: a) harmoniczne, b
c) losowe, d) skokowe i e) impulsowe:
t
J
sprężystych, ale również rozkładu jej masy. Konieczna jest także znajomość efektywnego
tłumienia drgań analizowanego układu. O jego wartości decydują zarówno własności tłumiące
t
t
a) b)
c)
t
d)
t
e)
24
samego materiału konstrukcyjnego, jak również parametry charakteryzujące interakcję układu
drgającego z otaczającym go ośrodkiem fizycznym (powietrzem, wodą). Efektywne tłumienie
drgań jest wielkością fizyczną niewyznaczalną na drodze analitycznej. Dlatego też musi ono
być określane eksperymentalnie, za pomocą odpowiedniej aparatury pomiarowej.
Podstawowym problemem dynamiki budowli jest obliczenie tzw. odpowiedzi konstrukcji Q(t)
ajprostszymi modelami obliczeniowymi konstrukcji drgających są modele jednomasowe o
na działanie danego obciążenia zewnętrznego P(t). Przy zadaniach statyki pomijało się efekty
związane ze zmianą energii kinetycznej przy przykładaniu obciążenia. W dynamice to
założenie odrzuca się jednocześnie rozróżniając drgania własne konstrukcji (wywołane
ciężarem konstrukcji) oraz wymuszone (wywołane przez czynniki zewnętrzne – np. drgające
silniki spoczywające na konstrukcji). Elementy powyższych wektorów interpretuje się
odpowiednio jako obciążenie zewnętrzne i przemieszczenie związane z i-tym stopniem
swobody. Macierze M,C,K określają własności modelu odpowiednio: masowe
(bezwładnościowe), tłumiące i sprężyste.
N
prętach nieważkich (tzn. pozbawionych masy). Stanowią one punkt wyjścia do tworzenia
innych modeli, w tym modeli wielomasowych – modeli statycznych uzupełnionych o pewną
liczbę punktowych mas skupionych. Cechą takiego modelu jest tzw. liczba stopni swobody
drgań (LSSD) – określająca w sposób jednoznaczny położenia wszystkich mas układu
drgającego w dowolnej chwili czasu. Dla układu z rysunku powyżej przy założeniu EA
= ∞
jest LSSD = n. Uwzględnienie rzeczywistych własności mechanicznych pręta ( EA
≠ ∞ )
powoduje zmianę modelu na bardziej realny o LSSD = 2n.
Jedną z podstawowych cech układu drgającego jest jego tzw. widmo częstości drgań
ł W. Kolousek. Zaproponował on przyjęcie dla
Q
(t)
P
(t)
własnych. Niskie częstości drgań własnych związane są z dominującym efektem zginania
prętów. Wraz ze wzrostem efektu rozciągania bądź ściskania, częstości rosną i są tym wyższe,
im większa jest sztywność EA. W większości przypadków, z jakimi spotykamy się w praktyce
obliczeniowej najważniejszą rolę odgrywają najniższe częstości drgań własnych. One to,
bowiem najczęściej decydują o występowaniu takich, niekorzystnych dla konstrukcji,
zjawisk, jak rezonans objawiający się nadmiernym wzrostem jej przemieszczeń, a co za tym
idzie – także przyspieszeń i towarzyszących im sił bezwładności powodujących czasem nawet
bardzo znaczne przeciążenie konstrukcji.
Autorem modeli bardziej złożonych by
każdego pręta LSSS =
∞ . Oznaczało to uwzględnienie ciągłego rozkładu masy wzdłuż osi
pręta. W konsekwencji otrzymał dla pojedynczego pręta jego równanie równowagi w postaci
równania różniczkowego cząstkowego. Przejście od analizy pojedynczych prętów do analizy
układu złożonego odbywa się według algorytmu metody przemieszczeń. Natomiast rewolucję
w zakresie budowania modeli dynamicznych wywołała metoda elementów skończonych
MES, gdzie zarówno modele o masach skupionych jak i ciągłe analizowane są na podstawie
równań różniczkowych zwyczajnych.
M,C,K
1
( )
Q t
( )
Q t
( )
Q t
1
m
n
2
EA
= ∞
2
m
m
( )
P t
( )
P t
n
1
2
( )
P t
n
25
Dla dowolnego układu drgającego warunek „równowagi dynamicznej” układu można zapisać
dzie:
•
korzystając z zasady d’Alemberta. Dla układu o jednym stopniu swobody jest ono w postaci:
( )
( )
( )
( ) 0
B t
T t
S t
P t
+
+
+
= ,
g
( )
( )
B t
M Q
••
= −
t
- siła bezwładności d’Alemberta drgającej masy,
•
- si
- siła wym
tąd otrzymuje się następujące równanie różniczkowe ruchu układu:
ła tłumienia występującego w układzie,
( )
( )
T t
C Q t
•
= −
•
- siła reakcji sprężystej układu,
( )
( )
S t
K Q t
= −
•
uszenia (obciążenia) zewnętrznego.
( )
P t
S
M
( )
( )
( )
( )
Q t
C Q t
K Q t
P t
••
•
+
+
=
.
dpowiedź Q(t), wywołana przez dane wymuszenie P(t), zależy, jak wiadomo, od warunków
O
początkowych
( )
Q t
Q
0
0
=
,
( )
Q t
Q
•
•
=
, określających w sposób jednoznaczny stan układu w
chwili początko
0
0
wej
0
t t
= .
Równanie można też
d
prze stawić w postaci:
2
1
( ) 2
( )
( )
( )
Q t
Q t
Q t
P t
M
γ
ω
••
•
+
+
=
,
po wprowadzeniu oznaczeń
2
2
,
C
K
M
M
γ
ω
=
=
a
i korzystając z zadady superpozycji
(t) jako sumę trzech odpowiedzi:
Q a t
Q b t
Q t
•
=
+
+
a trzy niezależne efekty wywołujące drgania:
• Wychylenie początkowe (dla
przedstawić można rozwiązanie Q
0
0
( )
( )
( )
( )
p
Q t
n
0
t t
= ):
0
0
( )
0
Q t
Q
=
≠ ,
• Prędkość początkowa (dla
):
•
0
t t
=
0
0
( )
0
Q t
Q
•
≠
,
=
• Obciążenie zewnętrzne (dla
0
t t
= ): P(t).
yniki można zapisać w postaci:
W
2
2
( )
[cos(
)
sin(
)],
t
d
d
d
d
a t
e
t
t
γ
γ
ω
ω
ω
ω
ω
−
=
+
=
γ
−
,
1
( )
sin(
),
d
d
b t
e
t
γ
ω
ω
−
=
26
0
1
( )
( ) (
)
t
p
Q t
P
b t
d
M
τ
τ τ
=
−
∫
- tzw. całka Duhamela.
Przypadkiem szczególnym jest brak tłumienia (
0
γ
= ). Wtedy otrzymuje się:
0
1
1
( ) cos( ),
( )
sin( ),
( )
( )sin[ (
)]
t
p
a t
t
b t
t
Q t
P
t
d
M
ω
ω
τ
ω
ω
ω
=
=
=
−
∫
τ τ
]
.
W omawianych powyżej wzorach i zależnościach występowały następujące wielkości
fizyczne:
[
/
rad s
ω
- kołowa częstość drgań (własnych),
2
2
1
d
2
ω
ω
γ
ω
ξ
=
−
=
−
- kołowa częstość drgań tłumionych (
/
ξ γ ω
=
),
pozwalające wyznaczyć wielkości z nimi związane:
2
[ ]
T
s
π
ω
=
- okres drgań własnych,
1
[
lub
/ ]
f
Hz
cykl s
T
=
- częstość drgań (własnych).
Powyższe wzory pozwalają obliczyć odpowiedź układu o LSSD = 1 dla dowolnych
warunków początkowych ruchu i dowolnego obciążenia zewnętrznego. Całkę Duhamela
oblicza się zazwyczaj numerycznie, chyba że rozwiązanie udaje się przedstawić w postaci
analitycznej.
W przypadku, gdy do układu znajdującego się w spoczynku (
0
0
( ) 0 , ( ) 0
Q t
Q t
•
=
=
) zostaje, w
chwili , przyłożone obciążenie o stałej wartości
0
p
t
>
P const
=
, otrzymuje się odpowiedź w
postaci:
(
)
( )
sin[
(
)]
p
t
t
d
d t
P
Q t
e
w t
d
M
γ
τ
τ τ
ω
−
−
=
−
∫
,
a dla braku tłumienia (
0
γ
= ):
( )
{1 cos[ (
)]}
st
Q t
Q
t
ω
τ
=
−
−
,
przy czym
2
/
(
st
Q
P K K
M
)
ω
=
=
jest przemieszczeniem analizowanego układu, wywołanym
przez statycznie działanie siły P. Miarą przewyższenia odpowiedzi statycznej przez
dynamiczną jest tzw. współczynnik dynamiczny
μ
wynoszący w rozważanym przypadku:
2
( )
max
2.0
st
t
st
st
Q
Q t
Q
Q
μ
=
=
=
.
Przy występowaniu tłumienia (
0
γ
≠ ) charakterystyczne są dwa stany występujące podczas
odpowiedzi układu w czasie. Pierwszy to zjawisko dążenia odpowiedzi układu do osiągnięcia
pewnego stanu ustalonego. Stanem ustalonym może być stan spoczynkowej równowagi
27
statycznej układu oraz stan równowagi statycznej. Przed osiągnięciem stanu spoczynkowego
odpowiedź Q(t) ma charakter nieustalony, przejściowy (tzw. stan nieustalony). Przechodzenie
od stanu nieustalonego do ustalonego odbywa się dzięki uwzględnianiu rozpraszania energii
drgań do otoczenia na skutek tłumienia obecnego w każdym układzie fizycznym.
Bardzo często występującym wymuszeniem jest wymuszenie harmoniczne ( )
sin( )
P t
P
t
θ
=
.
Wtedy odpowiedzią jest:
(
)
0
1
( )
sin( )
sin[
(
)]
t
t
d
d
Q t
P
e
w t
d
M
γ
τ
θτ
τ
ω
−
−
=
−
∫
τ
,
a dla braku tłumienia (
0
γ
= ):
( )
[sin( )
sin( )]
st
Q t
Q
t
t
μ
θ
ν
ω
=
−
, przy czym:
2
2
1/(1
),
/
/
/(
)
st
m
Q
P K
P M
μ
ν
ν θ ω
ω
=
−
=
=
=
,
Szczególnym przypadkiem jest stan, kiedy różnica częstości drgań własnych i harmonicznych
θ ω
− = Δ
ω maleje, zbliżając się do przypadku granicznego (rezonansowego) θ ω
=
. Taki
charakterystyczny rodzaj drgań nazywa się dudnieniami, które powstają zawsze wtedy, gdy
mamy do czynienia z superpozycją (nakładaniem się) drgań harmonicznych o mało
różniących się częstościach. Wraz z przybliżaniem się do rezonansu mamy
,
1
θ ω ν
≈
≈ , co
oznacza się zbliżanie się amplitudy drgań (maksymalnej odpowiedzi) do sinusoidy o bardzo
długim okresie T
.
Δ
≈ ∞
W przypadku szczególnym
θ ω
=
mamy do czynienia z wymuszeniem ( )
sin( )
P t
P
t
ω
=
,
wywołującym w układzie drgającym tzw. zjawisko rezonansu na jego częstości własnej
ω .
Wtedy cechą charakterystyczną odpowiedzi jest liniowo-asymptotyczny, względem czasu,
wzrost amplitudy do nieskończoności. Nieograniczone narastanie amplitudy nie występuje
jednak nigdy w układach rzeczywistych, a to ze względu n fakt występowania tłumienia.
Drgania własne układu o wielu stopniach swobody
Postać ogólna równania ruchu układu drgającego o wielu (
n) stopniach swobody:
( )
( )
( )
( )
t
t
t
••
•
+
+
=
M Q
C Q
KQ
P t .
Macierze M,C,K są macierzami kwadratowymi o wymiarach (LSSD
× LSSD), natomiast P i
Q
to macierze jednokolumnowe (wektory) o wymiarach (LSSD
×1). Elementy macierzy
sztywności K i bezwładności M oblicza się na podstawie ich następujących definicji
(
):
,
1, 2,...,
i j
LSSD
=
Element
macierzy K jest reakcją i-tego więzu kinematycznego (siłą uogólnioną)
równoważącą działanie obciążenia (wymuszenia) kinematycznego pochodzącego od
jednostkowego przemieszczenia
j-tego więzu;
ij
k
1
j
Δ =
28
Element
macierzy M jest reakcją i-tego więzu kinematycznego (siłą uogólnioną)
równoważącą działanie sił bezwładności d’Alemberta obciążających układ na skuter
jednostkowego przyspieszenia
, wymuszonego za pomocą j-tego więzu kinematycznego.
ij
m
1
j
Q
••
=
Elementy macierzy tłumienia C oblicza się najczęściej ze wzoru
α
β
=
+
C
M
K .
Wyznaczenie stałych ,
α β
odbywa się na drodze eksperymentalnej.
Dla drgań własnych układu, nietłumionych jest:
0, ( ) 0
t
=
≡
C
P
, co prowadzi do równania
ruchu:
( )
( ) 0
t
t
••
+
=
M Q
KQ
.
Podstawiając możliwą odpowiedź układu w postaci ( )
sin(
)
t
t
ω ϕ
=
+
Q
A
otrzymujemy
równanie:
2
( -
)
0
ω
=
K
M A
reprezentujące problem własny dynamiki budowli. Z algebry wiadomo, iż
ma on niezerowe rozwiązanie, gdy
2
Det( -
) 0
ω
=
K
M
. Po rozwinięciu wyznacznika
otrzymuje się równanie wielomianowe mające dokładnie n dodatnich i rzeczywistych
pierwiastków zwanych częstościami (kołowymi) drgań własnych układu. Zbiór tych
częstości, uporządkowany niemalejąco:
1
2
0
...
n
ω ω
ω
<
≤
≤ ≤
nazywa się widmem częstości drgań własnych układu. Z równania problemu własnego można
otrzymać, po podstawieniu
i
ω ω
=
, odpowiednie niezerowe rozwiązanie
i
A , określające i-tą
postać drgań własnych, która jest i-tym wektorem własnym problemu. Pełne rozwiązanie
problemu własnego tworzy zbiór n tzw. par własnych ( , )
i
i
ω
A . Po wprowadzeniu dwóch
dodatkowych macierzy: spektralnej (widmowej)
Ω i modalnej (własnej)
V
:
1
2
n
ω
ω
ω
⎡
⎤
⎢
⎥
⎢
⎥
=
⎢
⎥
⎢
⎥
⎣
⎦
O
Ω
,
,
1
2
[ ,
,...,
]
n
=
V
A A
A
można sformułować tzw. warunek ortogonalności:
T
V MV = M , w którym macierz M jest
macierzą diagonalną o elementach
i
m , nazywaną modalną macierzą mas.
W zapisie problemu własnego zamiast macierzy sztywności K można stosować macierz
podatności D, o następującej własności: KD = DK = I. Wtedy problem własny przybiera
postać:
2
2
( -
)
0 lub (
-
)
0,
1/
ω
λ
λ
=
=
I
DM A
DM
I A
ω
=
.
29
Rozwiązanie problemu własnego stanowi klucz do utworzenia ważnego algorytmu
dekompozycji modalnej przy obliczaniu odpowiedzi Q(t). Algorytm ten nazywa się analizą
modalną układów drgających.
Równania drgań wymuszonych harmonicznie
Zakłada się, na konstrukcję działa obciążenie harmoniczne
0
( )
sin( )
t
t
θ
=
P
P
wywołujące
zmienne w czasie przemieszczenia, a także siły bezwładności. Równanie macierzowe
opisujące drgania wymuszone harmonicznie:
2
0
( -
)
θ
=
K
M Q
P
0
.
Dla danego wektora amplitud obciążenia wymuszającego
0
P należy obliczyć wektor
amplitud przemieszczeń
.
0
Q
Chcąc wyznaczyć amplitudy przemieszczeń
w czasie drgań wymuszonych harmonicznie,
przy danej macierzy podatności D i macierzy mas M, należy rozwiązać równanie
macierzowe:
0
Q
2
0
( -
)
θ
=
0
I
DM Q
DP .
W analizie drgań wymuszonych używane jest także nastające równanie macierzowe:
-1
0
0
2
1
( -
)
0
θ
=
D
M
B + DP
,
zawierające macierz podatności D, macierz mas M oraz – jako poszukiwane wielkości –
amplitudy sił bezwładności
. Równanie powyższe ma także skróconą postać:
0
B
0
0
0
p
+
=
*
D B
Δ
dzięki następującym podstawieniom:
-1
0
0
2
1
,
p
θ
= −
=
*
D
D
M
D
Δ
P .
8. Założenia i podstawy analizy statycznej prętów cienkościennych
Pręt cienkościenny to taki pręt, którego jeden z wymiarów określających przekrój poprzeczny
(grubość) jest nieporównywalnie mały w stosunku do drugiego. Pręty cienkościenne
wymagają innej analizy statycznej niż pręty lite, gdyż w ich przypadku nie obowiązują pewne
założenia upraszczające dla prętów litych: hipoteza Bernoulliego (przekrój płaski przed
przyłożeniem obciążenia pozostaje płaski po obciążeniu i prostopadły do ugiętej osi) i zasada
de Saint Venanta (przyłożone obciążenie można zastąpić innym statycznie równoważnym –
naprężenia, przemieszczenia i odkształcenia nie zmienią się poza niewielkim obszarem
zawierającym obciążone pole). Pręty cienkościenne można podzielić na pręty o profilu
otwartym i zamkniętym.
Analiza statyczna pręta cienkościennego zakłada wprowadzenie tzw. wycinkowych wielkości
geometrycznych na jego przekroju poprzecznym na bazie tzw. teorii grafu. Na profilu pręta
cienkościennego wybiera się tzw. korzeń, od którego odchodzą dendryty (drzewa) do każdego
innego węzła. Definiuje się lokalny krzywoliniowy układ współrzędnych ortogonalnych,
30
którego współrzędną skierowaną wzdłuż krawędzi profilu jest tzw. współrzędna łukowa.
Licząc całkę od bieguna (dowolny punkt płaszczyzny przekroju) aż do dowolnego punktu na
drodze skierowanej otrzyma się tzw. współrzędną wycinkową, stanowiącą podstawę
wyprowadzeń wszystkich innych wielkości geometrycznych profilu (tzw. charakterystyki
wycinkowe i geometryczne): momentu statycznego, momentów odśrodkowych, momentu
bezwładności i innych. Ważnymi pojęciami są także środek zginania (takie położenie
bieguna, dla którego obydwa wycinkowe momenty odśrodkowe zerują się) oraz główny punkt
zerowej współrzędnej wynikowej (główny korzeń).
Dla prętów cienkościennych otwartych poddanych prostemu skręcaniu wprowadzono
następujące założenia upraszczające:
• Przekrój poprzeczny aproksymuje się skończoną ilością prostokątów,
• Prostokąty pracują niezależnie – tak jednak, że jednostkowy kąt skręcenia każdego
nich jest taki sam.
Innymi pojęciami charakterystycznymi dla teorii prętów cienkościennych są bimoment i
bipara. Na ich podstawie powstaje odpowiednik zasady de Saint Venanta w prętach
cienkościennych (o kinematycznej równoważności układów).
U podstaw teorii pręta cienkościennego, oprócz założeń poprzednich, leżą następujące
hipotezy sformułowane przez Własowa:
• Powierzchnia środkowa (powierzchnia równoodległa od powierzchni górnej i dolnej)
deformuje się tak, jakby w płaszczyźnie każdego przekroju rozpostarta była na linii
środkowej sztywna tarcza, idealnie jednak wiotka dla deformacji w kierunku
prostopadłym do tego przekroju,
• Odkształcenia kątowe w punktach powierzchni środkowej są małe (
0
xs
γ
≈ ),
• Naprężenia normalne
n
σ
do powierzchni równoległych do powierzchni środkowej są
pomijalnie małe w stosunku do dwóch pozostałych naprężeń normalnych
,
x
s
σ σ
.
Do wyprowadzenia równania różniczkowego rządzącego (np. równania ruchu) stosuje się
najczęściej w teorii pręta cienkościennego podejście kinematyczne (przypuszczenie funkcji
przemieszczeń).
9. Siły przekrojowe w ustrojach prętowych
Rozważany będzie dowolny pręt określony w układzie współrzędnych globalnych
(
1
2
,
,
3
X X X ). Do pręta przyłożony jest dowolny układ sił zewnętrznych, również określony w
tym układzie poprzez współrzędne sił
1
2
3
( , , )
i
i
i
P
P P P
=
i
)
oraz punkty ich przyłożenia
1
2
3
(
,
,
i
i
i
i
A
X X X
=
1, 2,...
i
=
,
. W dowolnym punkcie Q dokonuje się przekroju poprzecznego
pręta: pręt rozpada się na dwie części: lokalny układ współrzędnych (
1
2
3
, ,
x x x ) powiązany jest
z II-gą częścią. Układ sił z części II zredukowany do bieguna Q składa się z dwóch wektorów
zaczepionych w Q:
• Sumy sił układu S :
i
i II
S
P
∈
=
∑
,
• Momentu układu względem bieguna Q:
,
i
i
i
i II
i
M
r P
r
QA
∈
=
×
=
∑
.
31
1
X
2
X
3
X
Przekrój poprzeczny
pręta płaszczyzną
Globalny układ
współrzędnych
3
x
1
x
2
x
II
I
Lokalny układ
współrzędnych
1
P
1
A
2
P
2
A
Q
S
M
Składowe S i M nazywa się siłami przekrojowymi (siła i moment). Słuszna jest następująca
definicja sił przekrojowych:
Siłami przekrojowymi przyłożonymi do przekroju poprzecznego, o wersorze normalnej
zewnętrznej
x
v
e
= , nazywać będzie się składowe zredukowanego, do środka ciężkości
przekroju poprzecznego, układu sił zewnętrznych przyłożonych do drugiej części przeciętej
bryły, odniesione do lokalnego układu współrzędnych.
Po rozkładzie składowych S i M w układzie lokalnym dostaje się następujący zestaw sił
przekrojowych:
I
;
z
z
F M
;
x
x
F M
;
y
y
F M
;
z
z
F M
;
x
x
F M
•
,
x
x
F M - siła normalna, moment skręcający,
•
- siły styczne (poprzeczne), momenty zginające.
, ,
,
y
z
y
F F M M
n
Postać kanoniczna:
,
,
x
y
z
x
y
z
x
y
z
x
y
z
S
F
F
F
F i
F j F k
M
M
M
M
M i
M j M k
=
+
+
=
+
+
=
+
+
=
+
+
32
Współrzędne:
( ,
, ),
(
,
,
)
x
y
z
x
y
z
S F F F
M M M M
,
Zapis tradycyjny:
(
,
,
),
(
,
,
)
x
y
z
x
y
z
S N Q Q
M M M M
.
Znakowanie sił przekrojowych:
Wartościom siły przekrojowej przypisuje się znak „+” wówczas, gdy:
• Normalna zewnętrzna przekroju ma zwrot zgodny (niezgodny) ze zwrotem osi, do
której jest równoległa oraz:
• Siła przekrojowa ma zwrot zgodny (niezgodny) ze zwrotem osi, do której jest
równoległa.
W przypadku, gdy nie zachodzi ta podwójna zgodność (lub podwójna niezgodność) sile
przekrojowej przypisujemy znak „-”.
Poprzez siły przekrojowe wyraża się rozkład sił wewnętrznych w przekroju pręta. W tym celu
buduje się funkcje sił przekrojowych w zależności od kolejnych położeń przekroju
poprzecznego (w układach płaskich - zmienna niezależna
x):
( ),
( ),
( )
( ),
( ),
( )
x
x
y
y
z
z
x
x
y
y
z
z
F
F x
F
F x
F
F x
M
M x
M
M x
M
M x
=
=
=
=
=
=
Opisów powyższych funkcji dokonuje się w tzw. przedziałach charakterystycznych,
wyznaczonych przez kolejne punkty charakterystyczne (początek i koniec pręta, punkty
przyłożenia obciążeń skupionych, początek i koniec obciążenia ciągłego, punkty nieciągłości
funkcji lub jej kolejnych pochodnych – np. przegub). Przy budowaniu funkcji sił
przekrojowych należy pamiętać, iż:
• W punktach przyłożenia skupionej siły zewnętrznej o niezerowej składowej
prostopadłej do osi pręta funkcja siły podłużnej jest nieciągła,
• W punktach przyłożenia skupionej siły zewnętrznej o niezerowej składowej
równoległej do osi funkcja siły podłużnej jest nieciągła,
• W punktach przyłożenia skupionego momentu zewnętrznego funkcja momentu jest
nieciągła.
W punktach nieciągłości nie definiuje się żadnej wartości tych sił. Z punktu widzenia
matematyki wartość wartości kolejnych pochodnych w punktach nieciągłości są opisywane
przez pseudo-funkcję Diraca (impulsowe przyłożenie obciążenia).
10. Analiza stanu naprężenia i odkształcenia w punkcie
Stan naprężenia
Stanem naprężenia w określonym punkcie ciała nazywamy zbiór naprężeń na wszystkich
dowolnie zorientowanych elementach przekrojów zawierających ten punkt.
Tensor stanu naprężenia
33
Z wnętrza bryły pozostającej w równowadze przy obciążeniu układem sił zewnętrznych
wycinamy element objętościowy w kształcie czworościanu, którego trzy krawędzie są
odpowiednio równoległe do osi układu współrzędnych. Wektor v normalnej zewnętrznej do
ściany ukośnej będzie wersorem o współrzędnych
α
ni
czyli:
n [
α
n1
;
α
n2
;
α
n3
]
Ponieważ długość wersora jest równa jedności, przeto współrzędne jego są cosinusami kątów
między wersorem n a osiami współrzędnych i oczywiście zachodzi równość:
α
n1
2
+
α
n2
2
+
α
n3
2
= 1
Ponieważ cosinusy kąta między dwiema ścianami jest równy cosinusowi kąta miedzy
normalnymi do tych ścian, więc możemy zapisać:
cos(n, x
i
) =
α
ni
Na każdy punkt powierzchni wyciętego czworościanu działają siły wewnętrzne, pochodzące
od punktów materialnych bryły. Zgodnie z poznanym twierdzeniem o równoważności
układów sił wewnętrznych i zewnętrznych, układ działający na ścianki czworościanu jest
równoważny układowi zerowemu.
Wynika stąd, że suma sił jest równa zeru, oraz moment względem dowolnego punktu jest
również równy zeru. Z równań równowagi sił otrzymujemy:
σ
n1
=
σ
11
α
n1
+
σ
21
α
n2
+
σ
31
α
n3
σ
n2
=
σ
12
α
n1
+
σ
22
α
n2
+
σ
32
α
n3
σ
n3
=
σ
13
α
n1
+
σ
23
α
n2
+
σ
33
α
n3
Zapisując drugi z warunków równowagi czworościanu, mianowicie przyrównując do zera
współrzędne M, otrzymujemy:
σ
12
=
σ
21
σ
13
=
σ
31
34
σ
23
=
σ
32
uwzględniając symetrie możemy zapisać:
σ
n1
=
σ
11
α
n1
+
σ
12
α
n2
+
σ
13
α
n3
σ
n2
=
σ
21
α
n1
+
σ
22
α
n2
+
σ
23
α
n3
σ
n3
=
σ
31
α
n1
+
σ
32
α
n2
+
σ
33
α
n3
Wykorzystując umowę sumacyjną, powyższy układ równań możemy zapisać w następującej,
skondensowanej formie:
σ
ni
=
σ
ij
α
nj
Tak więc, jeśli dana jest macierz naprężeń w punkcie oraz dane są współrzędne wersora
normalnego dowolnej płaszczyzny, którą przecinamy bryłę, przechodzącą przez ten punkt, to
za pomocą powyższego układu równań obliczymy współrzędne wektora naprężenia,
odpowiadającego temu przecięciu.
Do określenia stanu naprężenia wystarcza znajomość macierzy naprężeń w każdym punkcie,
a więc znajomość sześciu funkcji:
σ
ij
=
σ
ij
(x
1
; x
2
; x
3
)
Powyższy układ równań możemy zapisać w postaci macierzowej:
p
n
= T
σ
⋅ v
Z powyższego zapisu odczytujemy, że macierz naprężeń jest tensorem drugiego rzędu.
Zgodnie bowiem z definicją, tensorem drugiego rzędu nazywamy macierz, która pomnożona
przez wektor daje w wyniku wektor.
Zauważmy, że macierz ta w wyniku twierdzenia o równowartości naprężeń stycznych, jest
macierzą symetryczną względem przekątnej.
Transformacja macierzy naprężeń do innego układu.
Macierz naprężeń jest zbiorem współrzędnych trzech funkcji wektorowych p
1
; p
2
; p
3
które
otrzymaliśmy przecinając bryłę płaszczyznami prostopadłymi kolejno do osi x
1
; x
2
; x
3
.
Gdybyśmy bryłę przecięli trzema innymi płaszczyznami odpowiednio prostopadłymi do osi
ortogonalnego układu (x’
1
; x’
2
; x’
3
), wówczas każdej z tych nowych płaszczyzn będą
przyporządkowane wektory naprężenia odpowiednio p’
1
; p’
2
; p’
3
których współrzędne w
układzie (x’
i
) zapiszemy w postaci macierzy:
σ’
11
σ’
12
σ’
13
T’
σ
=
σ’
21
σ’
22
σ’
23
σ’
31
σ’
32
σ’
33
Macierz T
σ
jest tensorem, wobec tego elementy jego transformują się do nowego układu
według prawa tensorowego.
Podczas transformacji naprężeń z jednego układu kierunków (i, j) na drugi (k, l) znajduje
zastosowanie prawo transformacji tensorowej. Jest to prawo, zgodnie z którym
transformują się współrzędne tensora przy obrocie układu współrzędnych. Jeśli obiekt
transformuje się zgodnie z tym prawem, to jest (z definicji) tensorem, to dla tensora
drugiego rzędu prawo transformacji ma postać:
35
σ
kl
=
α
ki
α
lj
σ
ij
gdzie: i, j = x, y, z – współrzędne tensora w starym układzie
k, l = ξ,
η, ζ – współrzędne tensora w nowym układzie
α
ki
i
α
lj
- elementy macierzy przejścia (kosinusy kierunkowe transformacji)
Pod pojęciem macierzy przejścia rozumiemy macierz współczynników dostaw
kierunkowych (kosinusów kierunkowych) nowego układu współrzędnych w starym
układzie (albo na odwrót). Kosinusy kierunkowe są to kosinusy kątów pomiędzy osiami
układów współrzędnych: pierwotnego i obróconego, są to współrzędne wersorów jednego
układu w drugim.
α
11
α
12
α
13
α
21
α
22
α
23
α
31
α
32
α
33
Kolejne wiersze macierzy przejścia są współrzędnymi wersorów (e’
i
) nowych osi,
odniesionymi do układu starego (x
i
), kolumny zaś przedstawiają współrzędne wersorów
starych osi e
i
określonych w układzie nowym (x’
i
).
Tensor stanu naprężenia jest tensorem symetrycznym tzn. tensorem, którego współrzędne
spełniają warunek:
σ
ij
=
σ
ji
.
2.3
Warunki równowagi wewnętrznej
Postulat zachowania równowagi zarówno przez ciało jako całość, jak i przez jego części
myślowo odcięte, formułuje się w klasycznej wytrzymałości materiałów w postaci
całkowej. W teorii sprężystości badamy równowagę „nieskończenie małych” elementów
ciała (czyli o wymiarach bliskich zeru), wynikiem czego są warunki równowagi
wewnętrznej w formie równań różniczkowych. Postać badanych elementów zależy od
przyjętego układu współrzędnych. I tak np. w układzie kartezjańskim element taki będzie
prostopadłościanem, w układzie współrzędnych walcowych będzie odcinkiem „klina” itd.
W dowolnym punkcie A ciała, które pozostaje w równowadze pod wpływem działających
na nie sił, przyjmijmy początek układu współrzędnych kartezjańskich x
1
, x
2
, x
3
.
Płaszczyzny odniesienia uważajmy za płaszczyzny trzech parami prostopadłych
przekrojów, przeprowadzonych przez rozpatrywany punkt. Nieskończenie blisko tego
poprowadźmy równolegle drugą trójkę parami prostopadłych przekrojów, odpowiednio w
odległościach: dx
1
, dx
2
, dx
3
. Odcinki te będą krawędziami utworzonego w ten sposób
elementarnego prostopadłościanu. Wydzielmy go całkowicie z wnętrza ciała, odrzucając
myślowo materię przylegająca do jego ścianek. Aby nasz element pozostał nadal w
równowadze, musimy teraz zastąpić oddziaływanie nań materii otaczającej siłami
powierzchniowymi (naprężeniami) i objętościowymi (masowymi). Wypadkowe sił
powierzchniowych zaczepiamy w środkach ciężkości poszczególnych ścianek biorąc pod
uwagę średnią wartość naprężeń w obrębie każdej ścianki. Wypadkową sił objętościowych
36
– odniesioną do jednostki objętości – zaczepiamy w środku ciężkości całego elementu.
Rzuty tej wypadkowej na osie współrzędnych oznaczamy odpowiednio przez X, Y, Z.
Przyjmując niejednorodny stan naprężenia musimy uwzględnić zmianę wartości
naprężenia (a więc i jego składowych) w kierunkach osi współrzędnych. Zmianę tę
zapiszemy jako różniczki wyrażające przyrost danej składowej, wywołany przyrostem
odpowiedniej współrzędnej, przy przejściu w kierunku prostopadłym do ścianki
przeciwległej. I tak jeżeli w przekrojach przecinających się w punkcie A występują
składowe (i, j, k -wersory osi x
1
, x
2
, x
3
):
a) na ściance o normalnej zewnętrznej – i:
σ
11
σ
12
σ
13
b) na ściance o normalnej zewnętrznej – j:
σ
21
σ
22
σ
23
c) na ściance o normalnej zewnętrznej – k:
σ
31
σ
32
σ
33
to na ścianki przeciwległe elementu działać będą:
a) na ściankę o normalnej zewnętrznej – i:
σ`
11
σ`
12
σ`
13
b) na ściankę o normalnej zewnętrznej – j:
σ`
21
σ`
22
σ`
23
c) na ściankę o normalnej zewnętrznej – k:
σ`
31
σ`
32
σ`
33
Przy tym:
σ`
11
=
σ
11
+
1
11
x
∂
∂
σ
dx
1
σ`
12
=
σ
12
+
1
12
x
∂
∂
σ
dx
1
σ`
13
=
σ
13
+
1
13
x
∂
∂
σ
dx
1
σ`
21
=
σ
21
+
2
21
x
∂
∂
σ
dx
2
σ`
22
=
σ
22
+
2
22
x
∂
∂
σ
dx
2
σ`
23
=
σ
23
+
2
23
x
∂
∂
σ
dx
2
σ`
31
=
σ
31
+
3
31
x
∂
∂
σ
dx
3
σ`
32
=
σ
32
+
3
32
x
∂
∂
σ
dx
3
σ`
33
=
σ
33
+
3
33
x
∂
∂
σ
dx
3
Przestrzenny układ sił, obciążający nasz element, musi spełniać sześć warunków
równowagi: trzy warunki momentów i trzy warunki rzutów na osie.
37
Warunki momentów:
Warunki momentów formułujemy względem osi ξ, η, ζ, przeprowadzonych przez środek
ciężkości prostopadłościanu równolegle do osi x
1
, x
2
, x
3
układu współrzędnych. W ten
sposób wyeliminujemy z rachunku siłę masową, siły powierzchniowe normalne i siły
powierzchniowe styczne, zarówno równoległe do danych osi, jak i przechodzące przez
odpowiednią oś.
Warunek momentów względem osi ξ przybierze postać:
σ
23
dx
3
dx
1
2
2
dx
+
σ`
23
dx
3
dx
1
2
2
dx
-
σ
32
dx
1
dx
2
2
3
dx
-
σ`
32
dx
1
dx
2
2
3
dx
= 0
σ
23
dx
3
dx
1
2
2
dx
+ (
σ
23
+
2
23
x
∂
∂
σ
dx
2
)dx
3
dx
1
2
2
dx
-
σ
32
dx
1
dx
2
2
3
dx
- (
σ
32
+
3
32
x
∂
∂
σ
dx
3
)dx
1
dx
2
2
3
dx
= 0
Dzieląc całe równanie przez połowę objętości elementu 0,5 dx
1
dx
2
dx
3
otrzymamy
równanie, zawierające wielkości skończone
σ
23
i
σ
32
oraz wielkości nieskończenie małe,
określające przyrosty odnośnych naprężeń w obrębie badanego elementu. Przy zmierzaniu
do zera z wymiarami tego elementu, różnice naprężeń dążyć będą również do zera i w
granicy otrzymamy równość
σ
23
=
σ
32
.
Analogiczne równości
σ
13
=
σ
31
i
σ
12
=
σ
21
otrzymamy z następujących dwóch warunków
momentów względem osi η i ζ. Te trzy równania noszą nazwę warunków Cauchy’ego i
stanowią twierdzenie o równoważności naprężeń stycznych w przekrojach wzajemnie
prostopadłych.
Warunki rzutów sił:
Warunek rzutów na oś x
1
:
-
σ
11
dx
2
dx
3
+
σ`
11
dx
2
dx
3
-
σ
21
dx
1
dx
3
+
σ`
21
dx
1
dx
3
-
σ
31
dx
1
dx
2
+
σ`
31
dx
1
dx
2
+
+ Xdx
1
dx
2
dx
3
= 0
-
σ
11
dx
2
dx
3
+ (
σ
11
+
1
11
x
∂
∂
σ
dx
1
)dx
2
dx
3
-
σ
21
dx
1
dx
3
+(
σ
21
+
2
21
x
∂
∂
σ
dx
2
) dx
1
dx
3
-
σ
31
dx
1
dx
2
+
+ (
σ
31
+
3
31
x
∂
∂
σ
dx
3
)dx
1
dx
2
+ Xdx
1
dx
2
dx
3
= 0
Po otwarciu nawiasów, zredukowaniu i podzieleniu pozostałych wyrazów przez objętość
elementu dx
1
, dx
2
, dx
3
otrzymamy:
1
11
x
∂
∂
σ
+
2
21
x
∂
∂
σ
+
3
31
x
∂
∂
σ
+ X = 0
W analogiczny sposób sformułowane warunki rzutów na oś x
2
i oś x
3
dadzą dwa następne
równania. Równania te, wraz z poprzednim, tworzą układ równań, zwanych ogólnie
warunkami równowagi wewnętrznej lub warunkami Naviera w postaci różniczkowej:
1
11
x
∂
∂
σ
+
2
21
x
∂
∂
σ
+
3
31
x
∂
∂
σ
+ X = 0
38
1
12
x
∂
∂
σ
+
2
22
x
∂
∂
σ
+
3
32
x
∂
∂
σ
+ Y = 0
1
13
x
∂
∂
σ
+
2
23
x
∂
∂
σ
+
3
33
x
∂
∂
σ
+ Z = 0
Warunki równowagi wewnętrznej zostały wyprowadzone przy założeniu, że rozpatrywane
ciało znajduje się w spoczynku (w równowadze statycznej).
Ze względu na warunek Cauchy`ego równania możemy zapisać:
1
11
x
∂
∂
σ
+
2
12
x
∂
∂
σ
+
3
13
x
∂
∂
σ
+ X = 0
1
21
x
∂
∂
σ
+
2
22
x
∂
∂
σ
+
3
23
x
∂
∂
σ
+ Y = 0
1
31
x
∂
∂
σ
+
2
32
x
∂
∂
σ
+
3
33
x
∂
∂
σ
+ Z = 0
W przypadku szczególnym kiedy mamy stan płaskich naprężeń np. w
płaszczyźnie x
1
, x
2
, wszystkie składowe w kierunku osi x
3
będą się zerować, czyli
σ
11
=
σ
31
=
σ
32
= Z= 0.
Wówczas warunki równowagi wewnętrznej Nawiera sprowadzą się do dwóch równań:
1
11
x
∂
∂
σ
+
2
12
x
∂
∂
σ
+ X = 0;
1
21
x
∂
∂
σ
+
2
22
x
∂
∂
σ
+ Y = 0
Jeżeli nie ma udziału sił masowych, równania te upraszczają się do postaci:
1
11
x
∂
∂
σ
+
2
12
x
∂
∂
σ
= 0;
1
21
x
∂
∂
σ
+
2
22
x
∂
∂
σ
= 0
2.4
Warunki brzegowe
Warunki brzegowe będziemy chcieli sformułować w ten sposób, aby wyrażały zależności
pomiędzy naprężeniami i siłami zewnętrznymi obciążającymi ciało – w elementach
ograniczonych częściowo powierzchnią zewnętrzną ciała. W tym celu rozpatrywane ciało
należy podzielić na elementy płaszczyznami równoległymi do płaszczyzn odniesienia.
Wówczas elementy położone na skraju bryły będą ograniczone od zewnątrz trzema
wzajemnie prostopadłymi płaszczyznami przekroju, a od zewnątrz częścią powierzchni
ciała. Przy dostatecznie małych wymiarach elementu powierzchnię tę możemy
aproksymować płaszczyzną, otrzymując w rezultacie elementarny czworościan wraz z
przyłożonymi do jego ścianek naprężeniami. Należy przy tym zaznaczyć, że siła
powierzchniowa p
n
(p
n
x
1
, p
n
x
2
, p
n
x
3
) na czwartej ukośnej ściance nie jest, w ścisłym tego
słowa znaczeniu naprężeniem, tylko jednostkową siłą zewnętrzną o wymiarze naprężenia.
Należy również wspomnieć, że pominięcie sił masowych jest uzasadnione tym, że będąc
proporcjonalne do objętości ciała wartości tych sił dążą do zera, wraz z wymiarami
elementu, szybciej od naprężeń.
39
Orientacje zewnętrznej ukośnej ścianki elementu określają kosinusy kierunkowe normalnej
zewnętrznej n, kosinusy te oznaczamy dla skrócenia zapisu przez:
cos (n, x
1
) = n
1
cos(n, x
2
) = n
2
cos(n, x
3
) = n
3
Pole ścianki ukośnej oznaczymy dF
n
, a pola ścianek prostopadłych o normalnych x
1
, x
2
, x
3
oznaczymy odpowiednio dFx
1,
dFx
2,
dFx
3
. Pola te możemy uważać za rzuty pola dF
n
w
odpowiednich kierunkach, skąd wynikają następujące związki:
dFx
1
= n
1
·dF
n,
dFx
2
= n
2
·dF
n,
dFx
3
= n
3
·dF
n
Za pomocą powyższych związków i oznaczeń warunki rzutów na osie x
1
, x
2
, x
3
przybiorą
postać:
σ
11
· n
1
·dF
n
+
σ
21
· n
2
·dF
n
+
σ
31
·n
3
·dF
n
- p
n
x
1
dF
n
= 0
σ
22
· n
2
·dF
n
+
σ
12
·n
1
·dF
n
+
σ
32
·n
3
·dF
n
- p
n
x
2
dF
n
= 0
σ
33
·n
3
·dF
n
+
σ
13
·n
1
·dF
n
+
σ
23
·n
2
dF
n
- p
n
x
3
dF
n
= 0
Po podzieleniu tych równań przez wspólny czynnik dF
n
, otrzymamy ostatecznie:
p
n
x
1
=
σ
11
·n
1
+
σ
21
·n
2
+
σ
31
·n
3
p
n
x
2
=
σ
12
·n
1
+
σ
22
·n
2
+
σ
32
·n
3
p
n
x
3
=
σ
13
·n
1
+
σ
23
·n
2
+
σ
33
·n
3
Wykorzystując warunki Cauchy’ego równania możemy zapisać następująco:
p
n
x
1
=
σ
11
·n
1
+
σ
12
·n
2
+
σ
13
·n
3
p
n
x
2
=
σ
21
·n
1
+
σ
22
·n
2
+
σ
23
·n
3
p
n
x
3
=
σ
31
·n
1
+
σ
32
·n
2
+
σ
33
·n
3
gdzie: cos (n, x
1
) = n
1
, cos(n, x
2
) = n
2
, cos(n, x
3
) = n
3
Dwa powyższe układy równań formułują poszukiwane warunki brzegowe dla naprężeń, w
przyjętym układzie współrzędnych kartezjańskich. Z uwagi na kształt badanego elementu
ciała, nazywane bywają one również warunkami czworościanu.
Przy użyciu konwencji sumacyjnej możemy równania zastąpić jednym równaniem:
p
ni
=
σ
ji
n
j
gdzie: i, j = 1, 2, 3
i – wskaźnik kierunków rzutów.
j – wskaźnik osi układu jako normalnej zewnętrznej przekroju.
40
Główne kierunki naprężeń.
Określmy naprężenie normalne na ściance ukośnej. Ponieważ normalna zewnętrzna
zachowuje oznaczenie n, a kosinusy kierunkowe (n,x
1
), (n,x
2
), (n,x
3
) nie posiadają
wskaźników, dlatego naprężenie to nazwiemy
σ
n
i do zapisania go możemy użyć wprost
jedne go ze wzorów transformacyjnych. Po uwzględnieniu powyższych uwag wzór ten
przyjmuje następującą postać:
1
3
31
3
2
23
2
1
12
2
3
33
2
2
22
2
1
11
2
2
2
n
n
n
n
n
n
n
n
n
n
α
α
σ
α
α
σ
α
α
σ
α
σ
α
σ
α
σ
σ
+
+
+
+
+
=
Naprężenie
σ
n
jest ciągłą i różniczkowalną funkcją kosinusów kierunkowych
α
n1
,
α
n2
,
α
n3
normalnej n do przekroju
(
)
3
2
1
,
,
n
n
n
n
f
α
α
α
σ
=
Ekstrema takiej funkcji są jednocześnie ekstremami naprężenia normalnego w ustalonym
punkcie ciała ze względu na kierunek normalnej n do badanego przekroju. Dodatkowym
warunkiem geometrycznym, którym związane są 3 zmienne jest funkcja:
(
)
0
1
,
,
2
3
2
2
2
1
3
2
1
=
−
+
+
=
n
n
n
n
n
n
g
α
α
α
α
α
α
W myśl ogólnych reguł badania ekstremów warunkowych funkcji wielu zmiennych
budujemy funkcję pomocniczą
(
)
(
)
(
)
3
2
1
3
2
1
3
2
1
,
,
,
,
,
,
n
n
n
n
n
n
n
n
n
g
f
F
α
α
α
λ
α
α
α
α
α
α
+
=
gdzie
λ jest czynnikiem Lagrange’a, który jest poszukiwaną. Warunki konieczne istnienia
ekstremum
0
1
=
∂
∂
n
F
α
0
2
=
∂
∂
n
F
α
0
3
=
∂
∂
n
F
α
Po wykonaniu działań otrzymujemy układ równań warunkowych
41
⎪
⎩
⎪
⎨
⎧
−
=
+
+
−
=
+
+
−
=
+
+
3
3
33
2
32
1
31
2
3
23
2
22
1
21
1
3
13
2
12
1
11
n
n
n
n
n
n
n
n
n
n
n
n
λα
α
σ
α
σ
α
σ
λα
α
σ
α
σ
α
σ
λα
α
σ
α
σ
α
σ
co można także zapisać
p
nx
= -
λα
n1
,
p
ny
= -
λα
n2
, p
nz
= -
λα
n3
,
1
3
31
3
2
23
2
1
12
2
3
33
2
2
22
2
1
11
2
2
2
n
n
n
n
n
n
n
n
n
α
α
σ
α
α
σ
α
α
σ
α
σ
α
σ
α
σ
λ
+
+
+
+
+
=
−
po przemnożeniu równań: pierwszego przez
α
n1
, drugiego przez
α
n2
, a trzeciego przez
α
n3
i
dodaniu, otrzymujemy
lub po prostu
λ
σ
−
=
nextr
Wynika z tego wprost wartość czynnika Lagrange’a
σ
σ
λ
−
=
−
=
nextr
Symbol
σ wprowadzono dla uproszczenia oznaczenia ekstremalnego naprężenia
normalnego do przekroju o normalnej zewnętrznej n. Po przekształceniu otrzymamy
p
nx
=
σα
n1
, p
ny
=
σα
n2
,
p
nz
=
σα
n3
,
Podniesienie równania do kwadratu i dodanie daje
2
2
2
2
σ
=
+
+
nz
ny
nx
p
p
p
Równocześnie wiemy, że
2
2
2
2
n
nz
ny
nx
p
p
p
p
=
+
+
Ostatnie dwa równania dają podstawy do wysnucia wniosku
n
p
=
σ
czyli:
naprężenie normalne o wartości ekstremalnej jest zarazem naprężeniem całkowitym
(wypadkowym) w odpowiadającym mu przekroju, czyli że nie towarzyszy mu składowa
styczna.
Naprężenia ekstremalne z uwagi na położenie przekroju nazywamy naprężeniami
głównymi. Przepiszmy równania warunkowe w zmodyfikowanej postaci – wykorzystamy
w nich zależności
42
(
)
(
)
(
)
⎪
⎩
⎪
⎨
⎧
=
−
+
+
=
+
−
+
=
+
+
−
0
0
0
3
33
2
32
1
31
3
23
2
22
1
21
3
13
2
12
1
11
n
n
n
n
n
n
n
n
n
α
σ
σ
α
σ
α
σ
α
σ
α
σ
σ
α
σ
α
σ
α
σ
α
σ
σ
Układ równań jest
• Jednorodny ze względu na niewiadome a, b, c;
• Ze względu na istniejący warunek geometryczny nie może posiadać rozwiązań
zerowych;
Z teorii równań wynika, że układ równań jednorodnych posiada rozwiązanie niezerowe
wtedy, gdy wyznacznik główny tego układu jest równy 0
(
)
(
)
(
)
0
33
32
31
23
22
21
13
12
11
=
−
−
−
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
Równanie tego typu nazywa się charakterystycznym, sekularnym lub wiekowym.
Rozwinięcie wyznacznika we wzorze daje następujące równanie ze względu na
σ
(kubiczne)
(
)
(
)
[
]
(
)
[
]
0
2
2
12
33
2
31
22
2
23
11
31
23
12
33
22
11
2
31
2
23
2
12
11
33
33
22
22
11
2
33
22
11
3
=
+
+
−
+
−
−
+
+
−
+
+
+
+
+
−
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
Równanie powyższe jest równoważne równaniu
(
)(
)(
)
0
3
2
1
=
−
−
−
σ
σ
σ
σ
σ
σ
w którym
σ
1
,
σ
2
,
σ
3
są pierwiastkami i mają fizyczny wymiar naprężeń głównych.
Po rozwinięciu
(
)
(
)
0
3
2
1
1
3
3
2
2
1
2
3
2
1
3
=
−
+
+
+
+
+
−
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
co jest równoważne
0
3
2
2
1
3
=
−
+
−
J
J
J
σ
σ
σ
Niezmienniki tensora stanu naprężenia
Z porównania równań o tych samych pierwiastkach widać, że mają one ten sam czynnik 1
przy niewiadomej przy najwyższej potędze. Wynika stąd równość przy niższych potęgach
– 2, 1 i 0 . W takim razie współczynniki w równaniu nie mogą zależeć od orientacji układu,
czyli muszą być niezmiennikami obrotu. Nazywamy je niezmiennikami stanu naprężenia
i oznaczamy
43
3
2
1
33
32
31
23
22
21
13
12
11
3
0
0
0
0
0
0
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
=
=
J
1
3
3
2
2
1
11
13
31
33
33
32
23
22
22
21
12
11
2
0
0
0
0
0
0
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
+
+
=
+
+
=
J
m
J
σ
3
1
=
zaś
(
)
(
)
3
2
1
33
22
11
3
1
3
1
σ
σ
σ
σ
σ
σ
σ
+
+
=
+
+
=
m
σ
m
- naprężenie średnie
Rozkład tensora stanu naprężenia
Rozkład stanu naprężenia na stany podstawowe
• Stan hydrostatyczny
3
2
1
σ
σ
σ
=
=
• Czyste ścinanie
0
3
2
1
=
+
+
σ
σ
σ
Celem jest wykazanie, że każdy (dowolny) stan naprężenia można rozłożyć na dwa stany
podstawowe (jw)
Dowód:
Przyjmijmy, że stan dowolny naprężenia N określony jest naprężeniami głównymi
σ
1
,
σ
2,
σ
3
, albo – inaczej mówiąc – tensorem naprężeń
3
2
1
0
0
0
0
0
0
σ
σ
σ
τ
σ
=
Wprowadzając
σ
m
jako naprężenie średnie
)
(
3
1
3
2
1
σ
σ
σ
σ
+
+
=
m
możemy przedstawić naprężenia główne jako tożsamości
(
)
,,
1
,
1
1
1
σ
σ
σ
σ
σ
σ
+
=
−
+
=
m
m
(
)
,,
2
,
2
2
2
σ
σ
σ
σ
σ
σ
+
=
−
+
=
m
m
44
(
)
,,
3
,
3
3
3
σ
σ
σ
σ
σ
σ
+
=
−
+
=
m
m
przy czym
m
σ
σ
=
,
1
m
σ
σ
=
,
2
m
σ
σ
=
,
3
m
σ
σ
σ
−
=
1
,,
1
m
σ
σ
σ
−
=
2
,,
2
m
σ
σ
σ
−
=
3
,,
3
Rozkładamy stan naprężenia N na dwa stany A i B.
'
3
'
2
'
1
'
3
'
2
'
1
σ
σ
σ
σ
σ
σ
=
=
⇒
⎪
⎭
⎪
⎬
⎫
A
STAN
0
''
3
''
2
''
1
''
3
''
2
''
1
=
+
+
⇒
⎪
⎭
⎪
⎬
⎫
σ
σ
σ
σ
σ
σ
B
STAN
Stan A nazywamy hydrostatycznym, stan B nazywamy czystym ścinaniem.
Przechodząc do rachunku tensorowego
3
''
2
''
1
''
3
'
2
'
1
'
3
2
1
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
σ
σ
σ
σ
σ
σ
σ
σ
σ
+
=
po uwzględnieniu
)
(
0
0
0
)
(
0
0
0
)
(
0
0
0
0
0
0
0
0
0
0
0
0
3
2
1
3
2
1
m
m
m
m
m
m
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
−
−
−
+
=
• Tensor w którym elementy położone na przekątnej głównej są sobie równe a pozostałe
elementy równe są zeru nosi nazwę tensora kulistego lub aksjatora
• Tensor, w którym suma elementów położonych na przekątnej głównej równa się zeru –
nosi nazwę dewiatora
W odniesieniu do stanu naprężenia używać będziemy nazw i symboli:
• Aksjator naprężeń
A
σ
,
• Dewiator naprężeń
D
σ
.
Symbolicznie rozkład taki przedstawiamy jako
σ
σ
σ
D
A
T
+
=
czyli każdy tensor można rozłożyć na aksjator i dewiator naprężeń
45
33
32
31
23
22
21
13
12
11
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
=
T
Odpowiednio
''
33
32
31
23
''
22
21
13
12
''
11
'
33
'
22
'
11
33
32
31
23
22
21
13
12
11
0
0
0
0
0
0
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
+
=
gdzie
m
m
m
m
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
−
=
−
=
−
=
+
+
=
=
=
=
33
''
33
22
''
22
11
''
11
33
22
11
'
33
'
22
'
11
)
(
3
1
Ostatecznie
(
)
(
)
(
)
4
4
4
4
4
4
3
4
4
4
4
4
4
2
1
4
4 3
4
4 2
1
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
D
m
m
m
A
m
m
m
−
−
−
+
=
33
32
31
23
22
21
13
12
11
33
32
31
23
22
21
13
12
11
0
0
0
0
0
0
Analiza odkształcenia
Rozważania dotyczące odkształcenia przestrzeni zajętej przez ciało o dowolnym kształcie
można prowadzić znając charakter zmian geometrycznych jakim ulegają poszczególne części
tego ciała. Zmiany te spowodowane są przemieszczeniem punktów odkształcalnego ciała pod
wpływem obciążeń. Przemieszczenie jest definiowane jako odcinek skierowany (wektor)
łączący dwa kolejne położenia punktu – przed i po deformacji.
Odkształceniem liniowym
włókna, przechodzącego przez dany punkt, będziemy nazywać
jego względną zmianę długości na wskutek przyłożenia sił.
46
Składowe liniowe odkształcenia:
3
3
3
2
2
2
1
1
1
x
u
x
u
x
u
x
x
x
∂
∂
=
∂
∂
=
∂
∂
=
ε
ε
ε
Odkształceniem kątowym
nazywamy połowę kąta, o jaki zmieni się kąt prosty między
dwoma włóknami przechodzącymi przez wspólny punkt i wzajemnie prostopadłymi przed
przyłożeniem obciążenia.
Składowe kątowe odkształcenia:
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
3
1
1
3
31
2
3
3
2
23
1
2
2
1
12
2
1
2
1
2
1
x
u
x
u
x
u
x
u
x
u
x
u
ε
ε
ε
47
Równania Cauchy’ego.
Równania opisujące odkształcenia liniowe i kątowe tworzą tzw.
równania Cauchy’ego, które mają następującą postać:
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
3
1
1
3
13
3
3
33
2
3
3
2
23
2
2
22
1
2
2
1
12
1
1
11
2
1
'
2
1
'
2
1
'
x
u
x
u
x
u
x
u
x
u
x
u
x
u
x
u
x
u
ε
ε
ε
ε
ε
ε
a w zapisie sumacyjnym:
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
i
j
j
i
ij
i
i
ii
x
u
x
u
x
u
1
2
1
'
ε
ε
Tensor odkształceń
Odkształcenia kątowe i liniowe stanowią zbiór funkcji, transformują się wg tego samego
prawa co składowe stanu naprężenia, który można przedstawić w postaci macierzy,
umieszczając na przekątnej odkształcenia liniowe, a poza nią odkształcenia kątowe.
Ten wniosek i twierdzenie o równowartości posunięć pozwala uważać tę macierz za tensor
symetryczny II rzędu zwany tensorem odkształceń
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
=
33
32
31
23
22
21
13
12
11
,
,
,
,
,
,
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
T
Z definicji odkształcenia kątowego wynika, że
ji
ij
ε
ε
=
, zatem macierz odkształceń jest
symetryczna.
Biorąc pod uwagę dwa dowolne punkty A i B bryły, które po deformacji przyjmują położenie
A' i B', współrzędne tych punktów, wektora A, B oraz funkcje wektorową przemieszczeń u(r)
możemy określić w układzie x
i.
( ) (
)
( )
i
i
i
i
i
x
u
u
x
AB
x
x
B
x
A
=
Δ
+
),
(
,
,
48
To różnica kwadratów odległości punktów A i B wyrażać się będzie następująco
j
i
ij
x
x
e
AB
B
A
Δ
⋅
Δ
⋅
⋅
=
−
~
2
'
'
2
2
Jeśli współrzędne punktów A i B wektora AB oraz u określimy w układzie x
i
'
( ) (
)
( )
i
i
i
i
i
x
u
u
x
AB
x
x
B
x
A
'
),
'
(
,
'
'
,
'
=
Δ
+
to różnice kwadratów długości zapiszemy
'
'
'
~
2
'
'
2
2
j
i
ij
x
x
e
AB
B
A
Δ
⋅
Δ
⋅
⋅
=
−
Ponieważ różnica kwadratów długości tych samych dwu wektorów nie zależy od przyjętego
układu dlatego
'
'
'
~
2
~
2
j
i
ij
j
i
ij
x
x
e
x
x
e
Δ
⋅
Δ
⋅
⋅
=
Δ
⋅
Δ
⋅
⋅
Występujące w tym wzorze
Δx'
i
Δx
i
są współrzędnymi wektora w układach x
i
i x'
i
. Jeśli układ
x'
i
określany jest w układzie x
i
macierzą o elementach
α
ij
to współrzędne wektora AB
związane są wzorem transformacyjnym
k
ki
i
x
x
Δ
⋅
=
Δ
α
podstawiając, otrzymujemy:
r
rj
k
ki
ij
j
i
ij
x
x
e
x
x
e
'
'
~
'
~
Δ
⋅
⋅
Δ
⋅
⋅
=
Δ
⋅
Δ
⋅
α
α
Zmieniając po prawej stronie wskaźniki sumacyjne i na k i j na r i odwrotnie otrzymujemy
jr
ik
kr
ij
j
jr
i
ik
kr
j
i
ij
e
e
x
x
e
x
x
e
α
α
α
α
⋅
⋅
=
Δ
⋅
⋅
Δ
⋅
⋅
=
Δ
⋅
Δ
⋅
~
'
~
'
'
~
'
~
Ostanie równanie jest prawem tensorowym, a wiec macierz jest tensorem. Biorąc obustronną
49
granicę relacji przy B
→A otrzymamy
jr
ik
kr
ij
e
e
α
α
⋅
⋅
=
'
Przy założeniu małych pochodnych e
ij
=
ε
ij
, a zatem mamy:
jr
ik
kr
ij
α
α
ε
ε
⋅
⋅
=
'
Transformacja macierzy T
ε na nowy układ, odbywa się po powyższym wzorze.
Tensor odkształceń można rozłożyć na aksjator i dewiator odkształceń
ε
ε
ε
D
A
T
+
=
po rozpisaniu
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
+
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
−
−
−
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
+
+
=
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
−
=
−
=
−
=
=
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
=
=
=
=
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
+
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
m
m
m
m
m
m
m
m
m
m
m
m
m
m
m
m
def
D
A
def
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
0
0
0
0
0
0
3
"
"
"
0
0
0
0
0
0
'
'
,
'
'
0
0
0
'
0
0
0
'
"
"
"
33
32
31
23
22
21
13
12
11
33
32
31
23
22
21
13
12
11
33
22
11
33
33
32
31
23
22
22
21
13
12
11
11
33
22
11
33
22
11
33
32
31
23
22
21
13
12
11
33
32
31
23
22
21
13
12
11
Niezmienniki stanu odkształcenia
W każdym punkcie ciała istnieją trzy wzajemnie prostopadłe kierunki, w których
odkształcenia liniowe osiągają wartości ekstremalne, zwane odkształceniami głównymi a ich
kierunki – kierunkami głównymi odkształceń. W kierunkach tych odkształcenia znikają.
Znaczy to, że element prostopadłościenny doznaje jedynie wydłużenia boków
33
22
11
ε
ε
ε
≥
≥
50
Odkształceniami głównymi nazywamy wartości główne tensora odkształcenia
ε
ij.
Obliczamy
je jako pierwiastki równania algebraicznego stopnia trzeciego.
0
3
2
2
1
3
=
−
⋅
+
⋅
−
ε
ε
ε
ε
ε
ε
I
I
I
i
i
i
gdzie
są niezmiennikami stanu odkształcenia
ε
ε
ε
3
2
1
,
,
I
I
I
(
)
(
)
)
(
2
2
1
2
12
33
2
13
22
2
23
11
31
23
12
33
22
11
3
2
13
2
23
2
12
11
33
33
22
22
11
2
33
22
11
1
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
⋅
+
⋅
+
⋅
−
⋅
⋅
+
⋅
⋅
=
+
+
−
⋅
+
⋅
+
⋅
=
−
=
+
+
=
=
I
I
I
ji
ij
jj
ii
ii
Kinematyczne warunki brzegowe.
Równania geometryczne są równaniami różniczkowymi. Obliczenie przemieszczeń przy
danych funkcjach odkształceń daje całą rodzinę rozwiązań. Z rodziny tej wybieramy te
rozwiązania, które spełniają zadane więzy kinematyczne bryły. Warunki jakie narzucają
więzy na przemieszczenia, nazywamy kinematycznymi warunkami brzegowymi. Z uwagi na
różnorodność typów więzów, kinematycznych warunków brzegowych nie da się ująć
generalnym zapisem, jak to było w przypadku statycznych warunków brzegowych.
Warunki nierozdzielności.
Musimy zwrócić uwagę, że rozwiązywanie równań geometrycznych przy danych funkcjach
odkształceń nie zawsze ma rozwiązanie. Wynika to z tego, że równań geometrycznych jest
sześć a szukanych funkcji tylko trzy. Dlatego jeżeli mamy sześć równań, a tylko trzy
niewiadome, to rozwiązanie możliwe jest tylko wtedy, jeśli trzy spośród równań są
kombinacją liniową trzech pozostałych, co prowadzi między innymi do pewnych zależności
między wyrazami wolnymi.
W przypadku równań geometrycznych również muszą istnieć zależności między wyrazami
wolnymi, czyli nie dla każdej macierzy odkształceń istnieje rozwiązanie dla funkcji
przemieszczeń. Z geometrycznej interpretacji macierzy odkształceń wynika, że w każdym
punkcie określa ona deformację jednostkowego sześcianu.
Gdyby obliczyć i narysować deformację każdego sześcianu tak, aby otrzymać kształt
zdeformowanej bryły, zauważymy, że deformacje poszczególnych nie mogą być niezależne,
jeśli zachowana ma być ciągłość bryły.
Związki te nazywamy warunkami nierozdzielności. Warunki nierozdzielności otrzymamy
drogą następujących operacji różniczkowych:
51
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
∂
∂
+
∂
∂
=
i
j
j
i
ij
x
u
x
u
2
1
ε
zróżniczkujemy dwukrotnie stronami względem x
r
oraz x
k
.
Ponadto odkształcenie ε
kr
wyrażone odpowiednio przez przemieszczenia zróżniczkujemy
dwukrotnie względem x
i
i x
j
, następnie ε
ik
względem x
j
i x
r
oraz ε
rj
względem x
i
i x
k
.
Otrzymamy : u
i
i u
j
– składowe przemieszczeń w punkcie
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
∂
∂
∂
∂
+
∂
∂
∂
∂
=
∂
∂
∂
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
∂
∂
∂
∂
+
∂
∂
∂
∂
=
∂
∂
∂
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
∂
∂
∂
∂
+
∂
∂
∂
∂
=
∂
∂
∂
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
∂
∂
∂
∂
+
∂
∂
∂
∂
=
∂
∂
∂
r
k
i
r
r
k
i
j
k
i
jr
r
k
j
k
r
k
j
i
r
j
ik
r
j
i
r
r
j
i
k
j
i
kr
r
k
j
j
r
k
j
i
r
k
ij
x
x
x
u
x
x
x
u
x
x
x
x
x
u
x
x
x
u
x
x
x
x
x
u
x
x
x
u
x
x
x
x
x
u
x
x
x
u
x
x
3
3
2
3
3
2
3
3
2
3
3
2
2
1
2
1
2
1
'
2
1
ε
ε
ε
ε
Jeżeli pomnożymy dwa ostatnie związki stronami przez –1 i dodamy stronami dwoma
pierwszymi, to otrzymamy:
0
2
2
2
2
=
∂
∂
∂
−
∂
∂
∂
−
∂
∂
∂
−
∂
∂
∂
k
i
jr
r
j
ik
j
i
kr
r
k
ij
x
x
x
x
x
x
x
x
ε
ε
ε
ε
związki te dla i,j,k,r=1,2,3 są szukanymi warunkami nierozdzielności. Jak widać warunków
tych jest 3
4
=81. Równań niezależnych jest natomiast sześć.
Mają one postać:
2
1
33
2
2
13
1
23
3
12
3
3
1
22
2
3
12
1
23
2
13
2
3
2
11
2
3
12
2
13
1
23
1
3
2
23
2
2
2
33
2
2
3
22
2
3
1
13
2
2
1
33
2
2
3
11
2
2
1
12
2
2
1
22
2
2
2
11
2
2
2
2
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
x
∂
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
+
∂
∂
−
∂
∂
∂
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
+
∂
∂
−
∂
∂
∂
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
+
∂
∂
−
∂
∂
∂
∂
∂
=
∂
∂
+
∂
∂
∂
∂
∂
=
∂
∂
+
∂
∂
∂
∂
∂
=
∂
∂
+
∂
∂
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
Równanie Hooke'a.
52
Stan naprężenia i odkształcenia można rozpatrywać zupełnie oddzielnie. Dopiero chcąc
ustalić związki pomiędzy składowymi obu stanów dochodzimy do właściwych podstaw teorii
sprężystości. Uważając składowe tensora naprężeń za przyczyny a składowe tensora
odkształceń za skutki, i stosując zasadę superpozycji, otrzymujemy sześć równań
określających każdą składową stanu naprężenia w formie funkcji liniowej wszystkich
składowych stanu naprężenia.
23
66
13
65
12
64
33
63
22
62
11
61
23
23
56
13
55
12
54
33
53
22
52
11
51
13
23
46
13
45
12
44
33
43
22
42
11
41
12
23
36
13
35
12
34
33
33
22
32
11
31
33
23
26
13
25
12
24
33
23
22
22
11
21
22
23
16
13
15
12
14
33
13
22
12
11
11
11
σ
σ
σ
σ
σ
σ
ε
σ
σ
σ
σ
σ
σ
ε
σ
σ
σ
σ
σ
σ
ε
σ
σ
σ
σ
σ
σ
ε
σ
σ
σ
σ
σ
σ
ε
σ
σ
σ
σ
σ
σ
ε
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
a
Współczynniki związków są funkcjami trzech zmiennych
)
,
,
(
3
2
1
x
x
x
a
a
ij
ij
=
. Odkształcenia
mają miejsce przy pojawieniu się naprężenia i znikają całkowicie przy zniknięciu naprężeń.
Rozwiązując powyższe równania ze względu na naprężenia , otrzymujemy:
23
66
13
65
12
64
33
63
22
62
11
61
23
23
56
13
55
12
54
33
53
22
52
11
51
13
23
46
13
45
12
44
33
43
22
42
11
41
12
23
36
13
35
12
34
33
33
22
32
11
31
33
23
26
13
25
12
24
33
23
22
22
11
21
22
23
16
13
15
12
14
33
13
22
12
11
11
11
ε
ε
ε
ε
ε
ε
σ
ε
ε
ε
ε
ε
ε
σ
ε
ε
ε
ε
ε
ε
σ
ε
ε
ε
ε
ε
ε
σ
ε
ε
ε
ε
ε
ε
σ
ε
ε
ε
ε
ε
ε
σ
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
⋅
+
⋅
+
⋅
+
+
⋅
+
⋅
=
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
b
Współczynniki b
ij
nazywamy współczynnikami sprężystości (zachodzi równość b
ij
=b
ji
)
Materiałem izotropowym
nazywamy materiał, który w dowolnym punkcie ma jednakowe
własności we wszystkich kierunkach. Związki pomiędzy naprężeniami i odkształceniami w
danym punkcie są jednakowe
)
,
,
(
)
,
,
(
3
2
1
3
2
1
'
x
x
x
b
x
x
x
b
ij
ij
=
′
′
′
Można również wykazać , że w miarę ograniczania ogólnej anizotropii pewnymi
uprzywilejowanymi kierunkami i płaszczyznami, liczba ta ulegnie dalszej stopniowej
redukcji. Dla ciała izotropowego otrzymujemy w końcu trzy różne stałe:
53
31
23
12
66
55
44
33
22
11
b
b
b
b
b
b
b
b
b
=
=
=
=
=
=
przy czym stałe są powiązane zależnością :
44
12
11
2
b
b
b
+
=
współczynniki b
12
i b
44
nazywane są stałymi
Lamego i oznaczane symbolami:
b
12
=
λ , b
44
=
μ
Po podstawieniu stałych Lamego otrzymujemy ogólne prawo Hooke’a dla ciała
izotropowego :
(
)
(
)
(
)
(
(
)
(
31
31
23
23
12
12
11
22
33
33
33
11
22
22
33
22
11
11
2
2
2
με
σ
με
σ
με
σ
ε
ε
λ
ε
μ
λ
σ
ε
ε
λ
ε
μ
λ
σ
)
)
ε
ε
λ
ε
μ
λ
σ
=
=
=
+
+
+
=
+
+
+
=
+
+
+
=
a po rozwiązaniu ze względu na odkształcenia :
(
)
(
)(
)
(
)
(
)(
)
(
)
(
)(
)
μ
σ
ε
μ
σ
ε
μ
σ
ε
σ
σ
μ
λ
λ
σ
μ
λ
μ
μ
λ
ε
σ
σ
μ
λ
λ
σ
μ
λ
μ
μ
λ
ε
σ
σ
μ
λ
λ
σ
μ
λ
μ
μ
λ
ε
31
31
23
23
12
12
11
22
33
33
33
11
22
22
33
22
11
11
2
2
3
2
2
3
2
2
3
=
=
=
⎥
⎦
⎤
⎢
⎣
⎡
+
+
−
+
+
=
⎥
⎦
⎤
⎢
⎣
⎡
+
+
−
+
+
=
⎥
⎦
⎤
⎢
⎣
⎡
+
+
−
+
+
=
ν
wprowadzając stałe:
(
)
(
)
(
)
G
E
E
=
+
=
=
+
=
+
+
ν
μ
ν
μ
λ
λ
μ
λ
μ
λ
μ
1
2
2
2
3
gdzie: E-moduł Younga,
ν- stała Poissona, G- moduł Kirchhoffa.
przy użyciu tych stałych prawo Hooke’a można przedstawić:
54
[
]
[
]
[
]
31
31
31
23
23
23
12
12
12
22
11
33
33
33
11
22
22
33
22
11
11
1
1
2
1
1
2
1
1
2
)
(
1
)
(
1
)
(
1
σ
σ
ν
ε
σ
σ
ν
ε
σ
σ
ν
ε
σ
σ
ν
σ
ε
σ
σ
ν
σ
ε
σ
σ
ν
σ
ε
⋅
=
⋅
+
=
⋅
=
⋅
+
=
⋅
=
⋅
+
=
+
⋅
−
⋅
=
+
⋅
−
⋅
=
+
⋅
−
⋅
=
G
E
G
E
G
E
E
E
E
Trzy pierwsze równania wyrażają odkształcenie liniowe w funkcji wszystkich trzech
naprężeń normalnych. Pozostałe równania określają związki między odpowiadającymi sobie
odkształceniami postaciowymi i naprężeniami stycznymi.
Pierwszą grupę równań zapisuje się często w innej formie.
(
)
[
]
(
)
[
]
(
)
[
]
3
:
3
1
1
3
1
1
3
1
1
33
22
11
33
33
22
22
11
11
σ
σ
σ
σ
νσ
σ
ν
ε
νσ
σ
ν
ε
νσ
σ
ν
ε
+
+
=
−
+
⋅
=
−
+
⋅
=
−
+
⋅
=
m
m
m
m
gdzie
E
E
E
1.9. Prawo zmiany objętości
Dodając powyższe równanie stronami otrzymujemy:
(
)
(
)
m
m
m
E
E
σ
ν
θ
ε
ε
ε
θ
σ
ν
ε
2
1
3
2
1
3
3
33
22
11
−
=
+
+
=
−
=
Wzór ten wyraża związek pomiędzy dylatacją i średnim naprężeniem normalnym , wzór ten
daje odpowiedź na pytanie , kiedy objętość nie ulega zmianie:
a. w przypadku gdy E=
∞, wówczas mamy do czynienia z ciałem doskonale sztywnym:
b. w przypadku gdy
ν=0.5, takie ciało nazywamy doskonale nieściśliwym (teoria
plastyczności):
c. w przypadku gdy
σ
m
=0, jest to typowy stan naprężenia zwany czystym ścinaniem, w
którym odkształcenie polega wyłącznie na zmianie postaci bez zmiany objętości.
55
Wprowadzając nowe stałe materiałowe:
-moduł Helmholtza
(
)
ν
2
1
2
−
=
E
H
-moduł ściśliwości
(
)
ν
2
1
3
−
=
E
K
otrzymujemy:
θ
σ
K
m
=
lub
m
m
H
ε
σ
2
=
Z matematycznego punktu widzenia prawo zmiany objętości mówi, że składowe aksjatora
naprężeń są proporcjonalne do składowych aksjatora odkształceń
ε
σ
HA
A
2
=
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
=
m
m
m
A
σ
σ
σ
σ
0
0
0
0
0
0
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
=
m
m
m
A
ε
ε
ε
ε
0
0
0
0
0
0
Prawo zmiany postaci
Odejmując w równaniach prawa Hooke’a po lewej stronie odkształcenie
ε
m
a po prawej jego
wartość określoną wzorem
(
)
m
m
E
σ
ν
ε
2
1
−
=
i dokonując następujących przekształceń :
(
)
[
]
(
)
[
]
(
)
(
)
[
]
(
)
(
)
m
m
m
m
m
m
m
m
m
G
E
E
E
E
E
σ
σ
σ
σ
ν
σ
ν
σ
ν
νσ
σ
νσ
σ
ν
σ
ν
νσ
σ
ν
ε
ε
−
=
−
+
=
=
+
−
+
=
+
−
−
+
=
=
−
−
−
+
=
−
11
11
11
11
11
11
2
1
1
1
1
1
2
3
1
1
2
1
3
1
1
a więc podobnie :
(
)
(
)
(
)
m
m
m
m
m
m
G
G
G
σ
σ
ε
ε
σ
σ
ε
ε
σ
σ
ε
ε
−
=
−
−
=
−
−
=
−
33
33
22
22
11
11
2
1
2
1
2
1
Ostatnie trzy równania określają, jak widać zależności liniowe pomiędzy położonymi na
przekątnej głównej elementami dewiatorów D
ε
i D
σ
Z zależności tych i ze związków:
31
31
31
31
23
23
23
23
12
12
12
12
1
1
1
ε
σ
σ
ε
ε
σ
σ
ε
ε
σ
σ
ε
⋅
=
⇒
⋅
=
⋅
=
⇒
⋅
=
⋅
=
⇒
⋅
=
G
G
G
G
G
G
56
otrzymamy zależności liniowe pomiedzy wszystkimi elementami tych dewiatorów, w formie
następujących rownań:
(
)
(
)
(
)
2
2
2
2
2
2
2
2
2
31
31
23
23
12
12
33
33
22
22
11
11
ε
σ
ε
σ
ε
σ
ε
ε
σ
σ
ε
ε
σ
σ
ε
ε
σ
σ
⋅
=
⋅
=
⋅
=
−
⋅
=
−
−
⋅
=
−
−
⋅
=
−
G
G
G
G
G
G
m
m
m
m
m
m
równania wyrażają prawo zmiany postaci
ε
σ
GD
D
2
=
po rozpisaniu mamy:
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
−
−
−
=
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
−
−
−
m
m
m
m
m
m
G
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
33
32
31
23
22
21
13
12
11
33
32
31
23
22
21
13
12
11
2
Wnioski:
- całą zmianę objętości opisuje aksjator tensora odkształcenia, nie opisuje on natomiast
zmiany postaci;
- zmianę postaci opisuje dewiator tensora odkształcenia, nie opisuje on natomiast
zmiany objętości
Równania liniowej teorii sprężystości.
Poniższe równania stanowią formalny, matematyczny zapis służący do rozwiązywania zadań.
Równania statyczne tworzą trzy warunki równowagi wewnętrznej Naviera. Równaniom tym
przyporządkowane są warunki brzegowe dla naprężeń.
57
⎪
⎪
⎪
⎩
⎪
⎪
⎪
⎨
⎧
=
+
∂
∂
+
∂
∂
+
∂
∂
=
+
∂
∂
+
∂
∂
+
∂
∂
=
+
∂
∂
+
∂
∂
+
∂
∂
0
0
0
3
3
33
2
32
1
31
2
3
23
2
22
1
21
1
3
13
2
12
1
11
P
x
x
x
P
x
x
x
P
x
x
x
σ
σ
σ
σ
σ
σ
σ
σ
σ
⎪
⎪
⎪
⎭
⎪
⎪
⎪
⎬
⎫
0
=
+
∂
∂
i
j
ij
P
x
σ
σ
ij
=
σ
ji
statyczne warunki brzegowe:
⎪
⎩
⎪
⎨
⎧
+
+
=
+
+
=
+
+
=
3
33
2
32
1
31
3
3
23
2
22
1
21
2
3
13
2
12
1
11
1
v
v
v
v
v
v
v
v
v
v
v
v
p
p
p
α
σ
α
σ
α
σ
α
σ
α
σ
α
σ
α
σ
α
σ
α
σ
⎪
⎭
⎪
⎬
⎫
vj
ij
vi
p
α
σ
=
Równania geometryczne tworzą sześć związków Cauchy’ego
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
⎟⎟
⎠
⎞
⎜⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
3
2
2
3
32
3
3
33
3
1
1
3
13
1
1
22
2
1
1
2
12
1
1
11
2
1
,
2
1
,
2
1
,
x
u
x
u
x
u
x
u
x
u
x
u
x
u
x
u
x
u
ε
ε
ε
ε
ε
ε
W zapisie sumacyjnym:
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
∂
∂
+
∂
∂
=
∂
∂
=
i
j
j
i
ij
i
i
ii
x
u
x
u
x
u
1
2
1
,
ε
ε
Plus kinematyczne warunki brzegowe (z uwagi na różnorodność typów więzów,
kinematycznych warunków brzegowych nie da się ująć generalnym zapisem).
Równania fizykalne tworzą sześć równań wyrażających prawo sprężystości Hooke’a
Ogólne prawo Hooke’a dla ciała izotropowego:
(
)
(
(
)
(
(
)
(
11
22
33
33
33
11
22
22
33
22
11
11
2
2
2
ε
ε
λ
ε
μ
λ
σ
ε
ε
λ
ε
μ
λ
σ
)
)
)
ε
ε
λ
ε
μ
λ
σ
+
+
+
=
+
+
+
=
+
+
+
=
31
31
23
23
12
12
με
σ
με
σ
με
σ
=
=
=
58
Wprowadzając stałe:
(
)
(
)
(
)
G
E
=
=
ν
E
+
=
=
+
+
+
ν
μ
μ
λ
λ
μ
λ
μ
λ
μ
1
2
2
2
3
gdzie: E-moduł Younga,
ν- stała Poissona,
G- moduł Kirchhoffa.
przy użyciu tych stałych prawo Hooke’a można przedstawić:
[
]
[
]
[
]
31
31
31
23
23
23
12
12
12
22
11
33
33
33
11
22
22
33
22
11
11
1
1
2
1
1
2
1
1
2
)
(
1
)
(
1
)
(
1
σ
σ
ν
ε
σ
σ
ν
ε
σ
σ
ν
ε
σ
σ
ν
σ
ε
σ
σ
ν
σ
ε
σ
σ
ν
σ
ε
⋅
=
⋅
+
=
⋅
=
⋅
+
=
⋅
=
⋅
+
=
+
⋅
−
⋅
=
+
⋅
−
⋅
=
+
⋅
−
⋅
=
G
E
G
E
G
E
E
E
E
Trzy pierwsze równania wyrażają odkształcenie liniowe w funkcji wszystkich trzech naprężeń
normalnych. Pozostałe równania określają związki między odpowiadającymi sobie
odkształceniami postaciowymi i naprężeniami stycznymi.
Związek pomiędzy dylatacją i średnim naprężeniem normalnym:
(
)
m
E
σ
ν
θ
2
1
3
−
=
Prawo zmiany objętości:
ε
σ
HA
A
2
=
-moduł Helmholtza
(
)
ν
2
1
2
−
=
E
H
⎟
⎟
⎟
⎠
⎞
m
σ
0
0
⎜
⎜
⎜
⎝
⎛
=
m
m
A
σ
σ
σ
0
0
0
0
⎜
⎜
⎜
⎝
⎛
=
m
m
A
ε
ε
ε
0
0
0
0
⎟
⎟
⎟
⎠
⎞
m
ε
0
0
Prawo zmiany objętości:
ε
σ
GD
D
2
=
59
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
−
−
−
=
⎟
⎟
⎟
⎠
⎞
⎜
⎜
⎜
⎝
⎛
−
−
−
m
m
m
m
m
m
G
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
ε
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
σ
33
32
31
23
22
21
13
12
11
33
32
31
23
22
21
13
12
11
2
Płaski stan naprężenia i odkształcenia
Jeżeli jedno z naprężeń (odkształceń) głównych jest równe zero, to taki stan naprężenia
(odkształcenia) nazywa się płaskim. Wówczas niezależnie od od obranego układu osi
niezmiennik I
3
=0. Jest to wystarczającym warunkiem dla stwierdzenia płaskości stanu
naprężenia (odkształcenia).
11. Stateczność prętów prostych, postacie utraty stateczności, metody określania sił
krytycznych
jl
60
Document Outline