pobieranie: 06 Ansys Analiza dokładności obliczeń MES
Pobierz dokument pdfPrzeglądaj wersję html pliku:
POLITECHNIKA SZCZECIŃSKA KATEDRA MECHANIKI I PODSTAW KONSTRUKCJI MASZYN
Ćwiczenie nr 6
Instrukcja do ćwiczeń laboratoryjnych Numeryczne metody analizy konstrukcji
Analiza dokładności obliczeń metodą elementów skończonych
Szczecin 1999
1
Cel ćwiczenia W niniejszej instrukcji przedstawiono zadanie (pręt o zmiennym przekroju obcią ony statycznie), rozwiązano je metodą analityczną, a następnie metodą elementów skończonych. Istotą ćwiczenia jest zobrazowanie wpływu podziału pręta na odkształcalne elementy skończone na dokładność otrzymanych wyników (przemieszczenia i odkształcenia) w odniesieniu do rozwiązania dokładnego uzyskanego metodą analityczną. 1. Opis zadania Dany jest sto kowy pręt o prostokątnym kształcie przekroju w płaszczyźnie x - y. Lewy koniec pręta jest utwierdzony, a prawy jest obcią ony osiową siłą P=20 kN. Moduł Younga jest równy 3000 kN/cm². Długość pręta wynosi 1m. Powierzchnia przekroju lewego końca pręta wynosi 10cm², a prawego - 1cm² i zmienia się liniowo według wzoru A(x)=10-0.09x.. Taka zmienność przekroju jest świadomie wybrana w celu wymuszenia wyraźnych zmian w przebiegu przemieszczeń i naprę eń. Poszukiwane są: wydłu enia całkowite wzdłu osi x – u(x), naprę enia normalne σ(x), jak równie reakcja R na lewym końcu pręta.
Rysunek 1.1. Przykład obliczeniowy 2. Rozwiązanie analityczne Aby to zadanie mo na było oszacować wykorzystując metodę elementów skończonych, konieczne jest wcześniejsze opracowanie dokładnego rozwiązania analitycznego. Uproszczenie zagadnienia do postaci idealnej W naszym przykładzie długość jest przewa ającym wymiarem w porównaniu z pozostałymi. W związku z tym, e obcią enie jest zadane na tym wymiarze mo emy zastosować teorię prętów. Wyprowadzenie równania ró niczkowego Aby wyprowadzić równanie ró niczkowe, nale y wziąć pod uwagę trzy warunki: -statyczny:
2
S ( x) = P = 20kN = const. S ( x) σ ( x) = A( x) -materiałowy: σ ( x) ε ( x) = ⇒ σ ( x) = ε ( x) ⋅ E E -geometryczny: du ε ( x) = ⇒ du = ε ( x) ⋅ dx dx Równanie ró niczkowe:
du σ ( x) S ( x) P = = = dx E A( x) ⋅ E A( x) ⋅ E P du = ⋅ dx A( x) ⋅ E
(2.1)
(2.2)
(2.3)
(2.4)
(2.5) (2.6)
Rozwiązanie równania ró niczkowego
u( x ) = ∫ du = ∫
x x 0 0
P P dx 20 dx 20 1 ⋅ dx = ⋅ ∫ = ⋅ = ⋅− ⋅ [ln( 10 − 0 ,09 ⋅ x ) − ln10] E ⋅ A( x ) E 0 A( x ) 3000 ∫ 10 − 0.09 ⋅ x 3000 0.09 0
x x
(2.7) Maksymalne wydłu enie całkowite końca pręta oblicza się następująco: 20 1 u ( x = 100) = ⋅− ⋅ [ln 1 − ln 10] = 0,17056cm . 3000 0,09
Przebieg wydłu enia całkowitego jest pokazany na rysunku 2.1.
Rys. 2.1. Przemieszczenia u(x) – rozwiązanie analityczne
Znając funkcję przemieszczeń u(x) mo na obliczyć naprę enia σ(x).
σ ( x) = E ⋅
du 20 = dx (10 − 0,09 ⋅ x )
(2.8)
3
dla x = 0
σ ( x = 0) = 2
kN cm 2
kN cm 2
dla x = 100
σ ( x = 100) = 20
Wykres zmienności naprę eń pokazano na rysunku 2.3. Zbie ność geometrii pręta powoduje szybki wzrost naprę eń przy prawym jego końcu.
Rys.2.3. Naprę enia σ(x) – rozwiązanie analityczne
3. Rozwiązanie za pomocą metody elementów skończonych
W tym rozdziale do rozwiązania przykładu zastosowana została metoda elementów skończonych. Aby uwydatnić tylko istotne cechy zagadnienia, podjęto tu szereg uproszczeń: obcią enie zostało sprowadzone do węzłów; geometria została uproszczona w taki sposób, e elementy tworzące pręt, na płaszczyźnie x-y są prostokątami. Na początku strukturę podzielono tylko na dwa elementy, później na więcej, wraz z zagęszczeniem elementów w obszarze wy szego gradientu naprę enia (w celu osiągnięcia lepszych wyników). Ze względu na strukturę elementu wybieramy element dwuwymiarowy le ący na płaszczyźnie x-y. Zaletą tego wyboru jest to, e element taki mo e być opisany wzdłu osi x. Pozwoli to zmniejszyć liczbę niewiadomych.
3.1. Podział na elementy (dyskretyzacja)
Pręt dzielimy najpierw na dwie części. Dla ka dego elementu wybieramy liniową postać przybli onego rozwiązania równania ró niczkowego dla nieznanej funkcji przemieszczenia u(x). Funkcja przemieszczenia u(x) zostanie zastąpiona przez funkcję liniową opisującą element. Dzięki temu, e elementy ściśle się stykają, zapewniona jest równa wartość przemieszenia w miejscu ich połączenia - czyli ciągłość przemieszczeń (patrz rys. 3.1).
4
Aproksymacja przebiegu przemieszczeń w elemencie Wybierzemy najprostszą postać funkcji liniowej u = b+a⋅x, gdzie a,b są wartościami stałymi. (3.1)
Zamiast stałych a i b wprowadzimy przemieszczenia węzłów, jako niewiadome u1 i u2. Przemieszczenia węzłów oznaczają tu te stopnie swobody. Funkcję przemieszczeń oblicza się jako iloczyn przemieszczeń węzłów i funkcji kształtu.
Rys.3.1. Model zło ony z dwóch elementów Dla
i u ( x = 0) = b = u1
u ( x = L) = b + a ⋅ L = u 2 , L – długość elementu, u 2 − u1 . L (3.2) (3.3)
otrzymujemy b = u1 i a = Następnie
u − u1 u 2 u u ( x) = u1 + 2 ⋅ x = ⋅ x + u1 − 1 ⋅ x L L L
x x u ( x) = 1 − ⋅ u1 + ⋅ u 2 L L
u ( x) = N 1 ( x) ⋅ u1 + N 2 ( x) ⋅ u 2
(3.4) (3.5)
5
gdzie N 1 ( x) = 1 −
x x i N 2 ( x) = L L są funkcjami kształtu. W zapisie macierzowym funkcja przemieszczeń ma postać:
u u ( x) = [N 1 ; N 2 ] ⋅ 1 = [N ] ⋅ {u} u 2 Funkcje kształtu Ni mo na obliczyć podstawiając za u1 = 1 i u2 = 0 , względnie u2 = 1 i u1 = 0 .
(3.6)
(3.7)
Jak widać na rys. 3.2 sformułowanie funkcji u(x) wg wz. (3.7), odpowiada wz. (3.1).
Rys.3.2. W ogólnym przypadku wzór (3.7) mo na zapisać jako: ux (3.8) u y = [N ] ⋅ {u} u z Oprócz przemieszczeń całkowitych będziemy potrzebować te przemieszczenia jednostkowe i naprę enia, wyra one za pomocą wielkości przemieszczeń w węzłach. du ( x) ε ( x) = (3.9) dx d ε ( x) = [N ] ⋅ {u} dx d [N ] d = [N 1 ; N 2 ] = [B1 ; B2 ] = [B ] , Wykorzystując dx dx otrzymujemy ε ( x) = [B ] ⋅ {u} (3.10)
W naszym przykładzie:
6
d x x 1 1 1 − L ; L = [B1 ; B2 ] = − L ; L . dx W ogólnym przypadku przemieszczenia jednostkowe mogą występować na ró nych kierunkach. Zatem wzór (3.10) mo na zapisać w postaci {ε } = [B ] ⋅ {u} . Naprę enie jest powiązane z przemieszczeniem jednostkowym następującą zale nością: (3.11) σ = E ⋅ε , gdzie E jest modułem sprę ystości. W ogólnym ujęciu wielkość skalarna E jest macierzą [D]. Jako macierz nie odnosi się do wartości naprę enia, ale do składowych naprę eń zestawionych razem w wektor naprę eń. {σ } = [D ] ⋅ {ε } . Podstawiając zale ność (3.10) otrzymujemy: {σ } = [D] ⋅ [B] ⋅ {u} . (3.13) Stąd mo na wyznaczyć przemieszczenia całkowite (3.7) oraz wielkości pochodne, jak przemieszczenia jednostkowe {ε} (3.11) i naprę enia {σ} (3.12).
[B] =
3.2. Tworzenie macierzy sztywności
Pozostałe do określenia niewiadome przemieszczenia węzłowe wyznacza się wykorzystując zasadę minimum energii potencjalnej. Potencjał п składa się z potencjału sił wewnętrznych пi i zewnętrznych пa. π = ∑ (π i ) e + (π a ) k + ∑ (π a ) e = min . (3.14) e e Potencjał sił wewnętrznych пi odpowiada sumie potencjałów pojedynczych elementów пie. Potencjał sił zewnętrznych пa równy jest sumie prac tych sił, пak działających na węzły i пae działających na element. Na początku siły zewnętrzne w obrębie elementu nie są brane pod uwagę (odpada więc człon пae). Pracę odkształcenia пie w obrębie elementu określa się jako: v 1 T (3.15) π ie = ∫ {ε } ⋅ {σ }dv . 20 Uwaga: Wektor przemieszczeń jednostkowych {ε} musi być transponowany, w celu zbudowania iloczynu skalarnego. Wykorzystując (3.11), (3.13) oraz zale ność: {ε }T = [B]T ⋅ {u}T , (3.16) otrzymujemy: v 1 T T π ie = ∫ {u} ⋅ [B ] ⋅ [D ] ⋅ [B ] ⋅ {u}dv . (3.17) 20 Podstawiając do powy szego wzoru:
[K ]e = ∫ [B]T ⋅ [D] ⋅ [B]dv ,
v
(3.18)
0
gdzie [K]e oznacza macierz sztywności, otrzymujemy: 1 T π ie = ⋅ {u} ⋅ [K ]e ⋅ {u}. 2
(3.19)
7
Potencjał sił zewnętrznych w węzłach oblicza się jako iloczyn przemieszczeń węzłowych i obcią eń w węzłach F: T π ak = −{u} ⋅ {F } . Potencjał globalny jest równy: 1 T T π = ⋅ ∑ {u} ⋅ [K ]e ⋅ {u} − {u} ⋅ {F } . (3.20) 2 e Z zasady minimum energii potencjalnej otrzymujemy macierz sztywności dla układu globalnego, przy pomocy której wyznaczymy nieznane przemieszczenia. ∂π = 0 ⇔ ∑ [K ]e ⋅ {u} − {F } = 0 . (3.21) ∂{u} e
(3.22) ∑ [K ]e ⋅ {u} = {F } . e [K ] ⋅ {u} = {F } , (3.23) gdzie macierz sztywności układu globalnego [K] jest sumą macierzy sztywności poszczególnych elementów [K]e, {u} jest wektorem przemieszczeń węzłowych, a {F} – wektorem sił w węzłach. W naszym przykładzie macierze [K]e i [K] oraz wektor {F} są dwuelementowe.
Wyznaczenie macierzy [K]e:
v 0
[K ]e = ∫ [B ]T ⋅ [D] ⋅ [B]dv = ∫ ∫ [B]T ⋅ [D] ⋅ [B]dadx = Am ⋅ ∫ [B]T ⋅ [D] ⋅ [B]dx ,
L A L 0 0 0
A1 + A2 - wartość średnia przekroju. 2 Macierz materiałowa [D] jest wielkością skalarną równą modułowi sprę ystości E. Mo emy więc zapisać: Am =
gdzie :
[B ]T ⋅ [D] ⋅ [B ] = [B ]T ⋅ E ⋅ [B ] .
[B] = − 1 L
1 ; L
[B ]
T
1 − = L , 1 L
E L2 = E ⋅ 1 − 1 , E L2 − 1 1 L2
L L
E 2 [B] ⋅ E ⋅ [B] = LE − 2 L
T
−
[K ]e
E 1 − 1 E ⋅ Am x − x = Am ⋅ ∫ [B ] ⋅ E ⋅ [B ]dx = Am ⋅ ∫ 2 dx = , − 1 1 L2 − x x 0 0 0 L
L T
1 − 1 . (3.24) L − 1 1 Macierz sztywności dla ka dego z elementów oblicza się wykorzystując dane z tabeli.
[K ]e = E ⋅ Am
8
Element i 1 2
E 3000 3000
Li 50 50
Ami 7,75 3,25
[K ]ei = E ⋅ Ami
Li
− 1 , − 1 1
1
[K ]e1 =
465 − 465 , − 465 465
[K ]e 2 =
195 − 195 . − 195 195
Macierze [K]ei budują główną macierz sztywności [K]. Wyznaczenie macierzy [K]: Globalny potencjał wynosi: π = (π i1 + π i 2 ) + π ak , Przy czym пi1 jest potencjałem wewnętrznym elementu 1, пi2 – potencjałem wewnętrznym elementu 2, natomiast пak zewnętrznym potencjałem sił węzłowych. T u1 F1 T T u u 1 u 1 u π = ⋅ 1 ⋅ [K ]1 ⋅ 1 + ⋅ 2 ⋅ [K ]2 ⋅ 2 − u 2 ⋅ F2 . u 2 2 u 3 2 u 2 u 3 u F 3 3 Wyra enie powy sze mo na zapisać w postaci: T T u1 u1 u1 F1 1 [K ] π = u 2 ⋅ 1 ⋅ u − u ⋅ F . (3.25) [K ]2 2 2 2 2 u u F u3 3 3 3 Podstawiając wartości liczbowe otrzymamy: T T − 465 0 u1 u1 F1 u1 465 1 π = u 2 ⋅ − 465 (465 + 195) − 195 ⋅ u 2 − u 2 ⋅ F2 . 2 − 195 195 u 3 u 3 F3 u 3 0 ∂π Zgodnie z warunkiem: = 0, ∂{u} otrzymujemy: 0 u1 R 465 − 465 − 465 660 − 195 ⋅ u = 0 (3.26) 2 u 20 0 − 195 195 3
[K ] ⋅ {u} = {F }
9
3.4. Rozwiązanie równania
Z równania: [K ] ⋅ {u} = {F }, mo na wyznaczyć nieznane przemieszczenia węzłowe {u}: [K ]−1 ⋅ [K ] ⋅ {u} = [K ]−1 ⋅ {F } ,
(3.27)
(3.28) W formie skalarnej wygląda to następująco: K ⋅u = F 1 u = ⋅F K W programach wykorzystujących metodę elementów skończonych, równanie będzie rozwiązane zgodnie z zasadą eliminacji Gaussa, co będzie pokazane w naszym przykładzie. Równanie wygląda następująco: 0 u1 R 465 − 465 − 465 660 − 195 ⋅ u = 0 2 0 − 195 195 u 3 20
[I ] ⋅ {u} = [K ]−1 ⋅ {F } {u} = [K ]−1 ⋅ {F } .
,
I
} II
W takiej postaci równania nie mo na rozwiązać, gdy nieznana jest reakcja R. W związku z tym równanie podzielimy na część zawierającą warunki brzegowe (II), oraz na część zawierającą nieznane przemieszczenia (I). Po rozwiązaniu (I) będzie mo na rozwiązać (II). Podział: 660 − 195 u 2 0 − 465 ⋅ u1 I − 195 195 ⋅ u = 20 − 0 3
II
[465
u1 − 465 0] ⋅ u 2 = {R} u 3
Określenie przemieszczeń węzłowych za pomocą algorytmu Gaussa: Transformacja do macierzy trójkątnej: 660 − 195 u 2 0 − 195 195 ⋅ u = 20 3 Dzielenie pierwszego wiersza przez 660: − 0.295 u 2 0 1 ⋅ = − 195 195 u 3 20 Dodanie do drugiego wiersza, wiersza pierwszego pomno onego przez 195: 1 − 0.295 u 2 0 0 137.386 ⋅ u = 20 3 Dzielenie drugiego wiersza przez 137.386:
10
1 − 0.295 u 2 0 ⋅ = 0 1 u 3 0.1456 Podstawienie:
u 3 = 0.1456
u 2 = 0 − (−0.295 ⋅ 0.1456) = 0.0428 Wynik z (I): u1 0 Wektor przemieszczeń: {u} = u 2 = 0.0428 u 0.1456 3
Wyznaczenie reakcji z (II):
(465 ⋅ u1 ) − (465 ⋅ u 2 ) + (0 ⋅ u3 ) = R R = 465 ⋅ (u1 − u 2 ) = 465 ⋅ (0 − 0.0428) = −19.902 ≈ −20[kN ]
3.5. Obliczenie naprę eń
Przemieszczenia lub naprę enia mo emy wyznaczyć z wektora przemieszczeń {u}. Element 1:
u1 0 Wektor przemieszczeń węzłowych: {u}1 = = u 2 0.0428 Przebieg przemieszczeń: x x u u ( x )1 = [N ] ⋅ {u}1 = 1 − ; ⋅ 1 L L u 2
x x u ( x )1 = 1 − ⋅ 0 + ⋅ 0.0428 = 0 + 0.000856 ⋅ x 50 50 u ( x = 0 )1 = 0cm
u ( x = 50 )1 = 0.0428cm .
E E 1 1 u1 σ ( x )1 = [D ] ⋅ [B ] ⋅ {u}1 = E ⋅ − ⋅ = − ⋅ u1 + ⋅ u 2 L L u 2 L L 1 σ ( x )1 = 3000 ⋅ ⋅ 0.0428 = 2.568[kN 2 ] . cm 50 Element 2: u 2 0.0428 Wektor przemieszczeń węzłowych: {u}2 = = u 3 0.1456 Przebieg przemieszczeń: x x u 2 u ( x )2 = [N ] ⋅ {u}2 = 1 − ; ⋅ L L u 3
11
Przebieg naprę eń:
x x u ( x )2 = 1 − ⋅ 0.0428 + ⋅ 0.1456 = 0.0428 + 0.0029 ⋅ x 50 50 u ( x = 0 )2 = 0.0428cm u ( x = 50 )2 = 0.1456cm . Przebieg naprę eń:
E E 1 1 u 2 σ ( x )2 = [D ] ⋅ [B] ⋅ {u}2 = E ⋅ − ⋅ u = − L ⋅ u 2 + L ⋅ u 3 L L 3 1 σ ( x )2 = 3000 ⋅ ⋅ ( −0.0428 + 0.1456) = 6.168[ kN 2 ] . cm 50
3.6. Ocena i kontrola wyników
Na rysunku rys. 3.3 przedstawiono przebieg przemieszczeń.
Rys.3.3. Przebieg przemieszczeń u(x) uzyskany w MES i analitycznie. Po porównaniu wyników uzyskanych metodą analityczną z wynikami metody elementów skończonych, mo na przejść do oceny przebiegu naprę eń.
12
Rys.3.4. Przebieg naprę eń σ(x) uzyskany w MES i analitycznie. Jak widać na rys. 3.4, metoda elementów skończonych pozwala określić stałe wartości na całej długości danego elementu. Nie ma więc zachowanej ciągłości naprę eń na granicy elementów, chocia jest zachowana ciągłość przemieszczeń (rys. 3.3). Warto zaznaczyć, e wartości naprę eń w środkach elementów równają się wartościom uzyskanym analitycznie. Jasnym jest teraz, e zwiększenie liczby elementów pozwoli uzyskać dokładniejsze wyniki.
4. Zwiększenie dokładności wyników
Dotychczas uzyskane wyniki nie są wystarczająco zadowalające. Ich jakość zwiększymy przez: - zagęszczenie elementów przez podział na większą ich ilość (4.1), - zagęszczenie elementów w obszarze większego gradientu naprę eń przez przesuwanie węzłów (4.1), - podwy szenie stopnia przybli onego rozwiązania równania ró niczkowego (4.2).
13
4.1.
Większa liczba elementów - 4 elementy skończone
Rys.4.1. Model zło ony z czterech elementów
Rys.4.2. Przebieg przemieszczeń u(x) w MES i uzyskany analitycznie.
Rys.4.3. Przebieg naprę eń σ(x) w MES i uzyskany analitycznie.
14
- 8 elementów o równej długości
Rys.4.4. Model zło ony z ośmiu elementów o równej długości.
Rys.4.5. Przebieg przemieszczeń u(x),dla ośmiu elementów o równej długości.
Rys.4.6. Przebieg naprę eń σ(x) dla ośmiu elementów o równej długości.
15
- 8 elementów zagęszczonych
Rys.4.7. Model zło ony z ośmiu elementów zagęszczonych.
Rys.4.8. Przebieg przemieszczeń u(x) dla ośmiu elementów zagęszczonych.
Rys.4.9. Przebieg naprę eń σ(x) dla ośmiu elementów zagęszczonych.
16
4.2.
Podwy szenie stopnia przybli onego rozwiązania równania ró niczkowego
Polepszenie wyników nie musi się wiązać tylko ze zmianą liczby elementów, lecz równie z podwy szeniem stopnia przybli onego rozwiązania równania ró niczkowego. Dlatego te w tej części równanie to będzie miało postać kwadratową. Taki wybór narzuci nam konieczność dodania dodatkowego węzła w środku elementu. Przybli one rozwiązanie równania ró niczkowego dla przemieszczeń mo na zapisać jako:
u ( x ) = a 0 + a1 ⋅ x + a 2 ⋅ x 2 .
(4.1)
Postępując analogicznie do poprzednich obliczeń zapiszemy:
u ( x ) = N1 ⋅ u1 + N 2 ⋅ u 2 + N 3 ⋅ u 3
u ( x ) = [N 1
N2
u ( x ) = [N ] ⋅ {u}
u1 N 3 ] ⋅ u 2 u 3
Rys.4.10. Zastosowanie kwadratowego równania ró niczkowego.
Odkształcenia jednostkowe oblicza się zgodnie ze wz. (3.9):
{ε } = du ( x ) = d ([N ] ⋅ {u}) = dN 1
dx dx
dx
dN 2 dx
dN 3 ⋅ {u} . dx
17
Przyjmując:
dN i = Bi , i = 1,2,3 , dx otrzymamy: {ε } = [B1 B2 B3 ] ⋅ {u} = [B ] ⋅ {u}.
W naszym przykładzie mamy:
{ε } = − 3 + 4 x L L2
4 8x − L L2
u1 1 4x − + 2 ⋅ u 2 . L L u 3
Elementy macierzy sztywności [K]e oblicza się wykorzystując wz. (3.18):
[K ]e = ∫ [B]T ⋅ [D] ⋅ [B]dv .
v
0
Poniewa : dv = A( x )dx , [D] = E to:
L
[K ]e = E ⋅ ∫ [B]T ⋅ [B] ⋅ A( x )dx .
B ⋅ B1 ⋅ A( x ) L 1 [K ]e = E ⋅ ∫ B2 ⋅ B1 ⋅ A( x ) 0 B3 ⋅ B1 ⋅ A( x )
0
B1 ⋅ B2 ⋅ A( x ) B1 ⋅ B3 ⋅ A( x ) B2 ⋅ B2 ⋅ A( x ) B2 ⋅ B3 ⋅ A( x ) ⋅ dx, B3 ⋅ B2 ⋅ A( x ) B3 ⋅B 3 ⋅ A( x ) f 13 ( x ) f 2 3 ( x ) ⋅ dx, f 3 3 ( x )
(4.2)
[K ]e
f11 ( x ) = E ⋅ ∫ f 21 ( x ) 0 f 31 ( x )
L
f 2 2 (x ) f 3 2 (x )
f1 2 (x )
(4.3)
[K ]e
K 11 = E K 21 K 31
L
K 12 K 22 K 32
K 13 K 23 K 33
(4.4)
Ogólnie współczynniki macierzy sztywności zapiszemy jako:
K ij = E ∫ Bi ⋅ B j ⋅ A( x )dx .
0
W celu usprawnienia obliczeń zastosujemy jedno przekształceń Gaussa, które przedstawia się następująco: L L 2 f ( x )dx ≈ ⋅ ∑ λi ⋅ f ( x )i , (4.5) ∫ 2 i =2 0 gdzie:
18
L ⋅ (ξ i + 1) , 2 λi – współczynnik wagowy. xi =
W tym przykładzie będą przeprowadzone obliczenia po dyskretyzacji na dwa elementy.
Podwy szenie stopnia przybli onego rozwiązania równania ró niczkowego
Dwa elementy. Elementy macierzy sztywności:
3 4x k11 = E ⋅ ∫ f11 ( x )dx = E ⋅ ∫ B1 ( x ) ⋅ B1 ( x ) ⋅ A( x )dx = E ⋅ ∫ − + 2 (10 − 0.09 x )dx L L 0 0 0
L L L 2
k11 ≈ E ⋅
L 2 ∑ λi ⋅ f11 ( xi ) . 2 1
Podstawiając za: λ1 = λ2 = 1.0
oraz ξ1 = +0.57735027; ξ 2 = −0.57735027 względnie: L 50 x1 = (ξ1 + 1) = (0.57735207 + 1) = 39.4337567 2 2 x 2 = 10.5662432, otrzymamy: k11≈1265 .
Komentarz [p1]:
Rys.4.11. Model zło ony z dwóch elementów.
19
Pozostałe współczynniki macierzy [K]1 wyznaczone analitycznie zestawiono w tabeli poni ej.
i;j
fij(x1) 0.000062 -0.000922 0.000860 -0.000922 0.013762 -0.012840 0.000860 -0.012840 0.011980
fij(x2) 0.016805 -0.018011 0.001207 -0.018011 0.019305 -0.001293 0.001207 -0.001293 0.000087
∑λ
2
v
v =1
⋅ f ij ( xv )
kij 1265 -1420 155 -1420 2480 -1060 155 -1060 905
1;1 1;2 1;3 2;1 2;2 2;3 3;1 3;2 3;3
0.016867 -0.018933 0.002067 -0.018933 0.033067 -0.014133 0.002067 -0.014133 0.012067
1265 − 1420 155 [K ]1 = − 1420 2480 − 1060. 155 − 1060 905
Według tego samego schematu oblicza się macierz [K]2: 65 635 − 700 [K ]2 = − 700 1040 − 340. 65 − 340 275 Globalna macierz sztywności przyjmie postać: 0 0 u1 R 1265 − 1420 155 − 1420 2480 − 1060 0 0 u 2 0 155 − 1060 1540 − 700 65 ⋅ u 3 = 0 0 − 700 1040 − 340 u 4 0 0 0 0 0 − 340 275 u 5 20
Rozwiązując powy szy układ otrzymamy wektor przemieszczeń: 0 0,018915 {u} = 0,044254 0,084160 0,166320 oraz reakcję R: R = −1420 ⋅ u 2 + 155 ⋅ u 3
} II I
R = −1420 ⋅ 0,018915 + 155 ⋅ 0,044254 = −19,99993 ≈ −20[kN ].
Wyniki:
20
Elementarne wektory przemieszczeń: 0 0.044254 {u}1 = 0,018915 {u}2 = 0,084160 . 0.044254 0.166320 Przebieg przemieszczeń: u(x)=[N]·{u}. Dla elementu 1.
u ( x )1 = N1 ⋅ u1 + N 2 ⋅ u 2 + N 3 ⋅ u 3 u ( x = 0)1 = 0cm u ( x = 25)1 = 0.01892cm u ( x = 50)1 = 0.04426cm
Dla elementu 2.
u ( x )1 = 0 + 6.282 ⋅ 10 −4 ⋅ x + 5.139 ⋅ 10 −6 ⋅ x 2
u ( x )2 = 0.04426 + 7.508 ⋅ 10 − 4 ⋅ x + 3.381 ⋅ 10 −5 ⋅ x 2 u ( x = 0)2 = 0.04426cm u ( x = 25)2 = 0.08416cm u ( x = 50)2 = 0.16632cm
u ( x )2 = N1 ⋅ u1 + N 2 ⋅ u 2 + N 3 ⋅ u 3
Rys.4.12. Przebieg przemieszczeń u(x). Przemieszczenia są określone tylko w węzłach. Powstałe na wykresie punkty są połączone odcinkami prostymi, chocia przybli one rozwiązanie równania ró niczkowego ma postać kwadratową.
21
Przebieg naprę eń:
{σ } = [D ] ⋅ [B ] ⋅ {u} [D ] ≡ E
L
{σ } = E ⋅ − 3 + 4 x 2
L
4 8x − L L2
−
1 4x + ⋅ {u} L L2
Dla elementu 1.
u1 E {σ }1 = 2 ⋅ [− 3L + 4 x 4 L − 8 x − L + 4 x] ⋅ u 2 . L u 3 Dla elementu 2:
{σ }2
u 3 E = 2 ⋅ [− 3L + 4 x 4 L − 8 x − L + 4 x ] ⋅ u 4 . L u 5
Przebieg naprę eń w elemencie jest w tym przypadku linowy.
Rys.4.13. Przebieg naprę eń σ(x).
Naprę enia zostały wyznaczone za pomocą MES w punktach całkowania. Na podstawie tych punktów ekstrapoluje się liniowo wartości naprę eń w węzłach. Ró nica wartości naprę eń na styku elementów pozwala oszacować na ile rozwiązanie jest dobre. Skok wartości naprę eń powinien przyjąć mo liwie najmniejszą wartość. Przez zwiększenie liczby elementów podnosi się dokładność wyników. Na poni szych rysunkach przedstawione są wyniki dla czterech i ośmiu elementów, uzyskane dla kwadratowej postaci przybli onego rozwiązania równania ró niczkowego.
22
Rys.4.14. Przebieg przemieszczeń u(x).
Rys.4.15. Przebieg naprę eń σ(x).
23
Rys.4.16. Przebieg przemieszczeń u(x).
Rys.4.17. Przebieg naprę eń σ(x).
5. Zestawienie wyników
W rozdziale tym przedstawiono i porównano wyniki z ró nych rozwiązań MES w celu wyciągnięcia wniosków końcowych. Tabela poni ej zawiera maksymalne przemieszczenia u(x) w [cm] na prawym końcu pręta:
24
Elementy 1 2 3 4 5 6 7 8 8 o ró nej długości
Aproksymacja linio- Aproksymacja kwaDokładne obliczenia wa dratowa 0,121212 0,145581 0,155437 0,160463 0,163371 0,165199 0,166419 0,167272 0,169490 0,156028 0,166320 0,168859 0,169750 0,170128 0,170310 0,170407 0,170461 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562
Na rys. 5.1 widać jak rośnie dokładność wyników wraz ze wzrostem liczby elementów. Tabela poni ej zawiera maksymalne naprę enia σ(x) w [kN/cm2]na końcu pręta: Elementy 1 2 3 4 5 6 7 8 8 o ró nej długości Aproksymacja linio- Aproksymacja kwa- Dokładne obliczenia wa dratowa 3,636 8,511 20,000 6,154 8,000 9,412 10,526 11,429 12,174 12,800 19,110 12,394 14,545 15,878 16,766 17,391 17,849 18,194 20,000 20,000 20,000 20,000 20,000 20,000 20,000 20,000
25
Rys.5.1. Maksymalne przemieszczenia węzłowe.
Rys.5.2. Maksymalne naprę enia elementów.
26
Rozwiązanie zadania w programie ANSYS
■ PREPROCESOR
1. Ustawienia preferencji Main Menu: Preferences Structural (will show) 2. Rysowanie przedmiotu wg rysunku
A
B
A1
D
C
(jedna powierzchnia stworzona z linii i punktów bazowych o współrzędnych jak na rysunku, przy czym aby stworzyć linie AB i CD nale y wskazać punkt punkt B oraz D jako pierwszy
3. Definiowanie typu elementu i opcji Main Menu: Preprocessor → Element Type → Add/Edit/Delete→ Element Types → Add... Wybierz element Quad 4node (PLANE 42)
Main Menu: Preprocessor → Element Type → Add/Edit/Delete→ Element Types → Options... K3: Plane strs w/thk
4. Nadawanie elementowi grubości Main Menu: Preprocessor → Real Constants→ Add.. → OK THK: 1 5. Definiowanie stałych materiałowych Main Menu: Preprocessor → Material Props → -Constant- Isotropic EX: 3000 NUXY:0.29 6. Utwierdzanie przedmiotu Main Menu: Preprocessor → Solution → -Loads- Apply → -Structural- Displacement → On Keypoints Wybierz dwa punkty bazowe znajdujące się po lewej stronie przedmiotu i ustaw pełne utwierdzenie (All DOF)
27
7. Definiowanie obcią enia Main Menu: Preprocessor → Solution → -Loads- Apply → -StructuralForce/Moment → On Keypoints Wybierz dwa punkty bazowe znajdujące się po prawej stronie przedmiotu i ustaw wartość siły FX = 10 (wypadkowa obu sił daje siłę 20kN umieszczoną na środku prawego końca) 8. Ustawienia liczby elementów 8.1 Włączenie numerowania elementów Utility Menu: PlotCtrls → Numbering... Elem / Attrib numbering: Element numbers 8.2 Zapisanie pracy Toolbar: SAVE_DB 8.3 Podział na dwa elementy o równej długości Main Menu: Preprocessor → -Meshing- Shape & Size → -Global- Size Wpisz w polu SIZE wartość nieco większą od 1/2 długości przedmiotu np.: SIZE: 51 8.4 Podział na cztery elementy o równej długości Main Menu: Preprocessor → -Meshing- Shape & Size → -Global- Size Wpisz w polu SIZE wartość nieco większą od 1/4 długości przedmiotu np.: SIZE: 26 8.5 Podział na osiem elementów o równej długości Main Menu: Preprocessor → -Meshing- Shape & Size → -Global- Size Wpisz w polu SIZE wartość nieco większą od 1/8 długości przedmiotu np.: SIZE: 13 9. Dyskretyzacja Main Menu: Preprocessor → -Meshing- Mesh → -Areas- Free 10. Rozwiązanie zadania Main Menu: Solution → -Solve- Current LS 11. Wykres przemieszczeń Main Menu: General Postproc → Path Operations → Define Path + Wska jeden z punktów na prawym i lewym końcu przedmiotu Path Operations → Map ontho Path... Wpisz etykietę: Lab User laber for item: UX Wybierz przemieszczenie względem osi X: Item, Comp: DOF solution; Translation UX Wyłącz uśrednianie wartości na długości elementu: Average results across element: No Path Operations → Plot Path Items... Wybierz przemieszczenia: UX
28
12. Wykres naprę eń Path Operations → Map ontho Path... Wpisz etykietę: Lab User laber for item: naprez Wybierz naprę enia względem osi X: Item, Comp: Stress; X-direction SX Wyłącz uśrednianie wartości na długości elementu: Average results across element: No Path Operations → Plot Path Items... Wybierz naprę enia: naprez 13. Reakcje Utility Menu: List → Results... → Reaction Soultion... All items
By przeprowadzić obliczenia dla większej liczby elementów nale y powtórzyć kroki od punktu 7.4. po uprzednim wczytaniu zbioru: Toolbar: RESUME_DB i odświe eniu ekranu: Utility Menu: Plot → Multi-Plots Dla obliczeń z u yciem ośmiu zagęszczonych elementów nale y podać sposób podziału linii AB i CD: Main Menu: Preprocessor → -Meshing- Shape and Size → -Global- Size w polu SIZE wpisać wartość większą lub równą 10 Main Menu: Preprocessor → -Meshing- Shape and Size → -Lines- Picked Lines wybrać linie AB i CD SIZE: 0 NDIV: 8 SPACE: 8 Main Menu: Preprocessor → -Meshing- Mesh → -Areas- Free
W przypadku u ycia elementów o wy szym stopniu równania ró niczkowego nale y zamiast czerowęzłowego elementu PLANE 42 u yć ośmiowęzłowego elementu PLANE 82 (punkt 2).
29
Ocena wyników
1. Zwiększenie liczby elementów pozwala uzyskać lepsze wyniki. 2. Zagęszczenie elementów w obszarze du ego gradientu naprę enia znacznie poprawia dokładność określenia przemieszczeń, a szczególnie naprę eń 3. Kwadratowa postać przybli onego rozwiązania równania ró niczkowego dostarcza dokładniejsze wyniki, ni liniowa jego postać. 4. Dla kwadratowej postaci przybli onego rozwiązania równania ró niczkowego uzyskuje się większą zbie ność wyników. 5. Na granicach elementów przemieszczenia mają równe wartości, natomiast naprę enia wykazują ró ne wartości. 6. Wartości naprę eń w środkach elementów lub w punktach całkowania, lepiej ni w węzłach zgadzają się z rzeczywistymi wartościami naprę eń. Przemieszczenia mo na określić dokładniej ni naprę enia (naprę enie jest wielkością pochodną). Przez zwiększenie liczby elementów wzrosła dokładność wyników. Rysunek 4.9 pokazuje, e większą dokładność mo na uzyskać przez zagęszczanie elementów w obszarze wy szego gradientu naprę eń. Nale y jeszcze raz podkreślić, e wybrany przykład szczególnie krytycznie podchodzi do obliczeń metodą elementów skończonych. Został on wybrany, aby zaznaczyć aproksymacyjny charakter metody.
30
06 Ansys Analiza dokładności obliczeń MES
POLITECHNIKA SZCZECIŃSKA KATEDRA MECHANIKI I PODSTAW KONSTRUKCJI MASZYN
Ćwiczenie nr 6
Instrukcja do ćwiczeń laboratoryjnych Numeryczne metody analizy konstrukcji
Analiza dokładności obliczeń metodą elementów skończonych
Szczecin 1999
1
Cel ćwiczenia W niniejszej instrukcji przedstawiono zadanie (pręt o zmiennym przekroju obcią ony statycznie), rozwiązano je metodą analityczną, a następnie metodą elementów skończonych. Istotą ćwiczenia jest zobrazowanie wpływu podziału pręta na odkształcalne elementy skończone na dokładność otrzymanych wyników (przemieszczenia i odkształcenia) w odniesieniu do rozwiązania dokładnego uzyskanego metodą analityczną. 1. Opis zadania Dany jest sto kowy pręt o prostokątnym kształcie przekroju w płaszczyźnie x - y. Lewy koniec pręta jest utwierdzony, a prawy jest obcią ony osiową siłą P=20 kN. Moduł Younga jest równy 3000 kN/cm². Długość pręta wynosi 1m. Powierzchnia przekroju lewego końca pręta wynosi 10cm², a prawego - 1cm² i zmienia się liniowo według wzoru A(x)=10-0.09x.. Taka zmienność przekroju jest świadomie wybrana w celu wymuszenia wyraźnych zmian w przebiegu przemieszczeń i naprę eń. Poszukiwane są: wydłu enia całkowite wzdłu osi x – u(x), naprę enia normalne σ(x), jak równie reakcja R na lewym końcu pręta.
Rysunek 1.1. Przykład obliczeniowy 2. Rozwiązanie analityczne Aby to zadanie mo na było oszacować wykorzystując metodę elementów skończonych, konieczne jest wcześniejsze opracowanie dokładnego rozwiązania analitycznego. Uproszczenie zagadnienia do postaci idealnej W naszym przykładzie długość jest przewa ającym wymiarem w porównaniu z pozostałymi. W związku z tym, e obcią enie jest zadane na tym wymiarze mo emy zastosować teorię prętów. Wyprowadzenie równania ró niczkowego Aby wyprowadzić równanie ró niczkowe, nale y wziąć pod uwagę trzy warunki: -statyczny:
2
S ( x) = P = 20kN = const. S ( x) σ ( x) = A( x) -materiałowy: σ ( x) ε ( x) = ⇒ σ ( x) = ε ( x) ⋅ E E -geometryczny: du ε ( x) = ⇒ du = ε ( x) ⋅ dx dx Równanie ró niczkowe:
du σ ( x) S ( x) P = = = dx E A( x) ⋅ E A( x) ⋅ E P du = ⋅ dx A( x) ⋅ E
(2.1)
(2.2)
(2.3)
(2.4)
(2.5) (2.6)
Rozwiązanie równania ró niczkowego
u( x ) = ∫ du = ∫
x x 0 0
P P dx 20 dx 20 1 ⋅ dx = ⋅ ∫ = ⋅ = ⋅− ⋅ [ln( 10 − 0 ,09 ⋅ x ) − ln10] E ⋅ A( x ) E 0 A( x ) 3000 ∫ 10 − 0.09 ⋅ x 3000 0.09 0
x x
(2.7) Maksymalne wydłu enie całkowite końca pręta oblicza się następująco: 20 1 u ( x = 100) = ⋅− ⋅ [ln 1 − ln 10] = 0,17056cm . 3000 0,09
Przebieg wydłu enia całkowitego jest pokazany na rysunku 2.1.
Rys. 2.1. Przemieszczenia u(x) – rozwiązanie analityczne
Znając funkcję przemieszczeń u(x) mo na obliczyć naprę enia σ(x).
σ ( x) = E ⋅
du 20 = dx (10 − 0,09 ⋅ x )
(2.8)
3
dla x = 0
σ ( x = 0) = 2
kN cm 2
kN cm 2
dla x = 100
σ ( x = 100) = 20
Wykres zmienności naprę eń pokazano na rysunku 2.3. Zbie ność geometrii pręta powoduje szybki wzrost naprę eń przy prawym jego końcu.
Rys.2.3. Naprę enia σ(x) – rozwiązanie analityczne
3. Rozwiązanie za pomocą metody elementów skończonych
W tym rozdziale do rozwiązania przykładu zastosowana została metoda elementów skończonych. Aby uwydatnić tylko istotne cechy zagadnienia, podjęto tu szereg uproszczeń: obcią enie zostało sprowadzone do węzłów; geometria została uproszczona w taki sposób, e elementy tworzące pręt, na płaszczyźnie x-y są prostokątami. Na początku strukturę podzielono tylko na dwa elementy, później na więcej, wraz z zagęszczeniem elementów w obszarze wy szego gradientu naprę enia (w celu osiągnięcia lepszych wyników). Ze względu na strukturę elementu wybieramy element dwuwymiarowy le ący na płaszczyźnie x-y. Zaletą tego wyboru jest to, e element taki mo e być opisany wzdłu osi x. Pozwoli to zmniejszyć liczbę niewiadomych.
3.1. Podział na elementy (dyskretyzacja)
Pręt dzielimy najpierw na dwie części. Dla ka dego elementu wybieramy liniową postać przybli onego rozwiązania równania ró niczkowego dla nieznanej funkcji przemieszczenia u(x). Funkcja przemieszczenia u(x) zostanie zastąpiona przez funkcję liniową opisującą element. Dzięki temu, e elementy ściśle się stykają, zapewniona jest równa wartość przemieszenia w miejscu ich połączenia - czyli ciągłość przemieszczeń (patrz rys. 3.1).
4
Aproksymacja przebiegu przemieszczeń w elemencie Wybierzemy najprostszą postać funkcji liniowej u = b+a⋅x, gdzie a,b są wartościami stałymi. (3.1)
Zamiast stałych a i b wprowadzimy przemieszczenia węzłów, jako niewiadome u1 i u2. Przemieszczenia węzłów oznaczają tu te stopnie swobody. Funkcję przemieszczeń oblicza się jako iloczyn przemieszczeń węzłów i funkcji kształtu.
Rys.3.1. Model zło ony z dwóch elementów Dla
i u ( x = 0) = b = u1
u ( x = L) = b + a ⋅ L = u 2 , L – długość elementu, u 2 − u1 . L (3.2) (3.3)
otrzymujemy b = u1 i a = Następnie
u − u1 u 2 u u ( x) = u1 + 2 ⋅ x = ⋅ x + u1 − 1 ⋅ x L L L
x x u ( x) = 1 − ⋅ u1 + ⋅ u 2 L L
u ( x) = N 1 ( x) ⋅ u1 + N 2 ( x) ⋅ u 2
(3.4) (3.5)
5
gdzie N 1 ( x) = 1 −
x x i N 2 ( x) = L L są funkcjami kształtu. W zapisie macierzowym funkcja przemieszczeń ma postać:
u u ( x) = [N 1 ; N 2 ] ⋅ 1 = [N ] ⋅ {u} u 2 Funkcje kształtu Ni mo na obliczyć podstawiając za u1 = 1 i u2 = 0 , względnie u2 = 1 i u1 = 0 .
(3.6)
(3.7)
Jak widać na rys. 3.2 sformułowanie funkcji u(x) wg wz. (3.7), odpowiada wz. (3.1).
Rys.3.2. W ogólnym przypadku wzór (3.7) mo na zapisać jako: ux (3.8) u y = [N ] ⋅ {u} u z Oprócz przemieszczeń całkowitych będziemy potrzebować te przemieszczenia jednostkowe i naprę enia, wyra one za pomocą wielkości przemieszczeń w węzłach. du ( x) ε ( x) = (3.9) dx d ε ( x) = [N ] ⋅ {u} dx d [N ] d = [N 1 ; N 2 ] = [B1 ; B2 ] = [B ] , Wykorzystując dx dx otrzymujemy ε ( x) = [B ] ⋅ {u} (3.10)
W naszym przykładzie:
6
d x x 1 1 1 − L ; L = [B1 ; B2 ] = − L ; L . dx W ogólnym przypadku przemieszczenia jednostkowe mogą występować na ró nych kierunkach. Zatem wzór (3.10) mo na zapisać w postaci {ε } = [B ] ⋅ {u} . Naprę enie jest powiązane z przemieszczeniem jednostkowym następującą zale nością: (3.11) σ = E ⋅ε , gdzie E jest modułem sprę ystości. W ogólnym ujęciu wielkość skalarna E jest macierzą [D]. Jako macierz nie odnosi się do wartości naprę enia, ale do składowych naprę eń zestawionych razem w wektor naprę eń. {σ } = [D ] ⋅ {ε } . Podstawiając zale ność (3.10) otrzymujemy: {σ } = [D] ⋅ [B] ⋅ {u} . (3.13) Stąd mo na wyznaczyć przemieszczenia całkowite (3.7) oraz wielkości pochodne, jak przemieszczenia jednostkowe {ε} (3.11) i naprę enia {σ} (3.12).
[B] =
3.2. Tworzenie macierzy sztywności
Pozostałe do określenia niewiadome przemieszczenia węzłowe wyznacza się wykorzystując zasadę minimum energii potencjalnej. Potencjał п składa się z potencjału sił wewnętrznych пi i zewnętrznych пa. π = ∑ (π i ) e + (π a ) k + ∑ (π a ) e = min . (3.14) e e Potencjał sił wewnętrznych пi odpowiada sumie potencjałów pojedynczych elementów пie. Potencjał sił zewnętrznych пa równy jest sumie prac tych sił, пak działających na węzły i пae działających na element. Na początku siły zewnętrzne w obrębie elementu nie są brane pod uwagę (odpada więc człon пae). Pracę odkształcenia пie w obrębie elementu określa się jako: v 1 T (3.15) π ie = ∫ {ε } ⋅ {σ }dv . 20 Uwaga: Wektor przemieszczeń jednostkowych {ε} musi być transponowany, w celu zbudowania iloczynu skalarnego. Wykorzystując (3.11), (3.13) oraz zale ność: {ε }T = [B]T ⋅ {u}T , (3.16) otrzymujemy: v 1 T T π ie = ∫ {u} ⋅ [B ] ⋅ [D ] ⋅ [B ] ⋅ {u}dv . (3.17) 20 Podstawiając do powy szego wzoru:
[K ]e = ∫ [B]T ⋅ [D] ⋅ [B]dv ,
v
(3.18)
0
gdzie [K]e oznacza macierz sztywności, otrzymujemy: 1 T π ie = ⋅ {u} ⋅ [K ]e ⋅ {u}. 2
(3.19)
7
Potencjał sił zewnętrznych w węzłach oblicza się jako iloczyn przemieszczeń węzłowych i obcią eń w węzłach F: T π ak = −{u} ⋅ {F } . Potencjał globalny jest równy: 1 T T π = ⋅ ∑ {u} ⋅ [K ]e ⋅ {u} − {u} ⋅ {F } . (3.20) 2 e Z zasady minimum energii potencjalnej otrzymujemy macierz sztywności dla układu globalnego, przy pomocy której wyznaczymy nieznane przemieszczenia. ∂π = 0 ⇔ ∑ [K ]e ⋅ {u} − {F } = 0 . (3.21) ∂{u} e
(3.22) ∑ [K ]e ⋅ {u} = {F } . e [K ] ⋅ {u} = {F } , (3.23) gdzie macierz sztywności układu globalnego [K] jest sumą macierzy sztywności poszczególnych elementów [K]e, {u} jest wektorem przemieszczeń węzłowych, a {F} – wektorem sił w węzłach. W naszym przykładzie macierze [K]e i [K] oraz wektor {F} są dwuelementowe.
Wyznaczenie macierzy [K]e:
v 0
[K ]e = ∫ [B ]T ⋅ [D] ⋅ [B]dv = ∫ ∫ [B]T ⋅ [D] ⋅ [B]dadx = Am ⋅ ∫ [B]T ⋅ [D] ⋅ [B]dx ,
L A L 0 0 0
A1 + A2 - wartość średnia przekroju. 2 Macierz materiałowa [D] jest wielkością skalarną równą modułowi sprę ystości E. Mo emy więc zapisać: Am =
gdzie :
[B ]T ⋅ [D] ⋅ [B ] = [B ]T ⋅ E ⋅ [B ] .
[B] = − 1 L
1 ; L
[B ]
T
1 − = L , 1 L
E L2 = E ⋅ 1 − 1 , E L2 − 1 1 L2
L L
E 2 [B] ⋅ E ⋅ [B] = LE − 2 L
T
−
[K ]e
E 1 − 1 E ⋅ Am x − x = Am ⋅ ∫ [B ] ⋅ E ⋅ [B ]dx = Am ⋅ ∫ 2 dx = , − 1 1 L2 − x x 0 0 0 L
L T
1 − 1 . (3.24) L − 1 1 Macierz sztywności dla ka dego z elementów oblicza się wykorzystując dane z tabeli.
[K ]e = E ⋅ Am
8
Element i 1 2
E 3000 3000
Li 50 50
Ami 7,75 3,25
[K ]ei = E ⋅ Ami
Li
− 1 , − 1 1
1
[K ]e1 =
465 − 465 , − 465 465
[K ]e 2 =
195 − 195 . − 195 195
Macierze [K]ei budują główną macierz sztywności [K]. Wyznaczenie macierzy [K]: Globalny potencjał wynosi: π = (π i1 + π i 2 ) + π ak , Przy czym пi1 jest potencjałem wewnętrznym elementu 1, пi2 – potencjałem wewnętrznym elementu 2, natomiast пak zewnętrznym potencjałem sił węzłowych. T u1 F1 T T u u 1 u 1 u π = ⋅ 1 ⋅ [K ]1 ⋅ 1 + ⋅ 2 ⋅ [K ]2 ⋅ 2 − u 2 ⋅ F2 . u 2 2 u 3 2 u 2 u 3 u F 3 3 Wyra enie powy sze mo na zapisać w postaci: T T u1 u1 u1 F1 1 [K ] π = u 2 ⋅ 1 ⋅ u − u ⋅ F . (3.25) [K ]2 2 2 2 2 u u F u3 3 3 3 Podstawiając wartości liczbowe otrzymamy: T T − 465 0 u1 u1 F1 u1 465 1 π = u 2 ⋅ − 465 (465 + 195) − 195 ⋅ u 2 − u 2 ⋅ F2 . 2 − 195 195 u 3 u 3 F3 u 3 0 ∂π Zgodnie z warunkiem: = 0, ∂{u} otrzymujemy: 0 u1 R 465 − 465 − 465 660 − 195 ⋅ u = 0 (3.26) 2 u 20 0 − 195 195 3
[K ] ⋅ {u} = {F }
9
3.4. Rozwiązanie równania
Z równania: [K ] ⋅ {u} = {F }, mo na wyznaczyć nieznane przemieszczenia węzłowe {u}: [K ]−1 ⋅ [K ] ⋅ {u} = [K ]−1 ⋅ {F } ,
(3.27)
(3.28) W formie skalarnej wygląda to następująco: K ⋅u = F 1 u = ⋅F K W programach wykorzystujących metodę elementów skończonych, równanie będzie rozwiązane zgodnie z zasadą eliminacji Gaussa, co będzie pokazane w naszym przykładzie. Równanie wygląda następująco: 0 u1 R 465 − 465 − 465 660 − 195 ⋅ u = 0 2 0 − 195 195 u 3 20
[I ] ⋅ {u} = [K ]−1 ⋅ {F } {u} = [K ]−1 ⋅ {F } .
,
I
} II
W takiej postaci równania nie mo na rozwiązać, gdy nieznana jest reakcja R. W związku z tym równanie podzielimy na część zawierającą warunki brzegowe (II), oraz na część zawierającą nieznane przemieszczenia (I). Po rozwiązaniu (I) będzie mo na rozwiązać (II). Podział: 660 − 195 u 2 0 − 465 ⋅ u1 I − 195 195 ⋅ u = 20 − 0 3
II
[465
u1 − 465 0] ⋅ u 2 = {R} u 3
Określenie przemieszczeń węzłowych za pomocą algorytmu Gaussa: Transformacja do macierzy trójkątnej: 660 − 195 u 2 0 − 195 195 ⋅ u = 20 3 Dzielenie pierwszego wiersza przez 660: − 0.295 u 2 0 1 ⋅ = − 195 195 u 3 20 Dodanie do drugiego wiersza, wiersza pierwszego pomno onego przez 195: 1 − 0.295 u 2 0 0 137.386 ⋅ u = 20 3 Dzielenie drugiego wiersza przez 137.386:
10
1 − 0.295 u 2 0 ⋅ = 0 1 u 3 0.1456 Podstawienie:
u 3 = 0.1456
u 2 = 0 − (−0.295 ⋅ 0.1456) = 0.0428 Wynik z (I): u1 0 Wektor przemieszczeń: {u} = u 2 = 0.0428 u 0.1456 3
Wyznaczenie reakcji z (II):
(465 ⋅ u1 ) − (465 ⋅ u 2 ) + (0 ⋅ u3 ) = R R = 465 ⋅ (u1 − u 2 ) = 465 ⋅ (0 − 0.0428) = −19.902 ≈ −20[kN ]
3.5. Obliczenie naprę eń
Przemieszczenia lub naprę enia mo emy wyznaczyć z wektora przemieszczeń {u}. Element 1:
u1 0 Wektor przemieszczeń węzłowych: {u}1 = = u 2 0.0428 Przebieg przemieszczeń: x x u u ( x )1 = [N ] ⋅ {u}1 = 1 − ; ⋅ 1 L L u 2
x x u ( x )1 = 1 − ⋅ 0 + ⋅ 0.0428 = 0 + 0.000856 ⋅ x 50 50 u ( x = 0 )1 = 0cm
u ( x = 50 )1 = 0.0428cm .
E E 1 1 u1 σ ( x )1 = [D ] ⋅ [B ] ⋅ {u}1 = E ⋅ − ⋅ = − ⋅ u1 + ⋅ u 2 L L u 2 L L 1 σ ( x )1 = 3000 ⋅ ⋅ 0.0428 = 2.568[kN 2 ] . cm 50 Element 2: u 2 0.0428 Wektor przemieszczeń węzłowych: {u}2 = = u 3 0.1456 Przebieg przemieszczeń: x x u 2 u ( x )2 = [N ] ⋅ {u}2 = 1 − ; ⋅ L L u 3
11
Przebieg naprę eń:
x x u ( x )2 = 1 − ⋅ 0.0428 + ⋅ 0.1456 = 0.0428 + 0.0029 ⋅ x 50 50 u ( x = 0 )2 = 0.0428cm u ( x = 50 )2 = 0.1456cm . Przebieg naprę eń:
E E 1 1 u 2 σ ( x )2 = [D ] ⋅ [B] ⋅ {u}2 = E ⋅ − ⋅ u = − L ⋅ u 2 + L ⋅ u 3 L L 3 1 σ ( x )2 = 3000 ⋅ ⋅ ( −0.0428 + 0.1456) = 6.168[ kN 2 ] . cm 50
3.6. Ocena i kontrola wyników
Na rysunku rys. 3.3 przedstawiono przebieg przemieszczeń.
Rys.3.3. Przebieg przemieszczeń u(x) uzyskany w MES i analitycznie. Po porównaniu wyników uzyskanych metodą analityczną z wynikami metody elementów skończonych, mo na przejść do oceny przebiegu naprę eń.
12
Rys.3.4. Przebieg naprę eń σ(x) uzyskany w MES i analitycznie. Jak widać na rys. 3.4, metoda elementów skończonych pozwala określić stałe wartości na całej długości danego elementu. Nie ma więc zachowanej ciągłości naprę eń na granicy elementów, chocia jest zachowana ciągłość przemieszczeń (rys. 3.3). Warto zaznaczyć, e wartości naprę eń w środkach elementów równają się wartościom uzyskanym analitycznie. Jasnym jest teraz, e zwiększenie liczby elementów pozwoli uzyskać dokładniejsze wyniki.
4. Zwiększenie dokładności wyników
Dotychczas uzyskane wyniki nie są wystarczająco zadowalające. Ich jakość zwiększymy przez: - zagęszczenie elementów przez podział na większą ich ilość (4.1), - zagęszczenie elementów w obszarze większego gradientu naprę eń przez przesuwanie węzłów (4.1), - podwy szenie stopnia przybli onego rozwiązania równania ró niczkowego (4.2).
13
4.1.
Większa liczba elementów - 4 elementy skończone
Rys.4.1. Model zło ony z czterech elementów
Rys.4.2. Przebieg przemieszczeń u(x) w MES i uzyskany analitycznie.
Rys.4.3. Przebieg naprę eń σ(x) w MES i uzyskany analitycznie.
14
- 8 elementów o równej długości
Rys.4.4. Model zło ony z ośmiu elementów o równej długości.
Rys.4.5. Przebieg przemieszczeń u(x),dla ośmiu elementów o równej długości.
Rys.4.6. Przebieg naprę eń σ(x) dla ośmiu elementów o równej długości.
15
- 8 elementów zagęszczonych
Rys.4.7. Model zło ony z ośmiu elementów zagęszczonych.
Rys.4.8. Przebieg przemieszczeń u(x) dla ośmiu elementów zagęszczonych.
Rys.4.9. Przebieg naprę eń σ(x) dla ośmiu elementów zagęszczonych.
16
4.2.
Podwy szenie stopnia przybli onego rozwiązania równania ró niczkowego
Polepszenie wyników nie musi się wiązać tylko ze zmianą liczby elementów, lecz równie z podwy szeniem stopnia przybli onego rozwiązania równania ró niczkowego. Dlatego te w tej części równanie to będzie miało postać kwadratową. Taki wybór narzuci nam konieczność dodania dodatkowego węzła w środku elementu. Przybli one rozwiązanie równania ró niczkowego dla przemieszczeń mo na zapisać jako:
u ( x ) = a 0 + a1 ⋅ x + a 2 ⋅ x 2 .
(4.1)
Postępując analogicznie do poprzednich obliczeń zapiszemy:
u ( x ) = N1 ⋅ u1 + N 2 ⋅ u 2 + N 3 ⋅ u 3
u ( x ) = [N 1
N2
u ( x ) = [N ] ⋅ {u}
u1 N 3 ] ⋅ u 2 u 3
Rys.4.10. Zastosowanie kwadratowego równania ró niczkowego.
Odkształcenia jednostkowe oblicza się zgodnie ze wz. (3.9):
{ε } = du ( x ) = d ([N ] ⋅ {u}) = dN 1
dx dx
dx
dN 2 dx
dN 3 ⋅ {u} . dx
17
Przyjmując:
dN i = Bi , i = 1,2,3 , dx otrzymamy: {ε } = [B1 B2 B3 ] ⋅ {u} = [B ] ⋅ {u}.
W naszym przykładzie mamy:
{ε } = − 3 + 4 x L L2
4 8x − L L2
u1 1 4x − + 2 ⋅ u 2 . L L u 3
Elementy macierzy sztywności [K]e oblicza się wykorzystując wz. (3.18):
[K ]e = ∫ [B]T ⋅ [D] ⋅ [B]dv .
v
0
Poniewa : dv = A( x )dx , [D] = E to:
L
[K ]e = E ⋅ ∫ [B]T ⋅ [B] ⋅ A( x )dx .
B ⋅ B1 ⋅ A( x ) L 1 [K ]e = E ⋅ ∫ B2 ⋅ B1 ⋅ A( x ) 0 B3 ⋅ B1 ⋅ A( x )
0
B1 ⋅ B2 ⋅ A( x ) B1 ⋅ B3 ⋅ A( x ) B2 ⋅ B2 ⋅ A( x ) B2 ⋅ B3 ⋅ A( x ) ⋅ dx, B3 ⋅ B2 ⋅ A( x ) B3 ⋅B 3 ⋅ A( x ) f 13 ( x ) f 2 3 ( x ) ⋅ dx, f 3 3 ( x )
(4.2)
[K ]e
f11 ( x ) = E ⋅ ∫ f 21 ( x ) 0 f 31 ( x )
L
f 2 2 (x ) f 3 2 (x )
f1 2 (x )
(4.3)
[K ]e
K 11 = E K 21 K 31
L
K 12 K 22 K 32
K 13 K 23 K 33
(4.4)
Ogólnie współczynniki macierzy sztywności zapiszemy jako:
K ij = E ∫ Bi ⋅ B j ⋅ A( x )dx .
0
W celu usprawnienia obliczeń zastosujemy jedno przekształceń Gaussa, które przedstawia się następująco: L L 2 f ( x )dx ≈ ⋅ ∑ λi ⋅ f ( x )i , (4.5) ∫ 2 i =2 0 gdzie:
18
L ⋅ (ξ i + 1) , 2 λi – współczynnik wagowy. xi =
W tym przykładzie będą przeprowadzone obliczenia po dyskretyzacji na dwa elementy.
Podwy szenie stopnia przybli onego rozwiązania równania ró niczkowego
Dwa elementy. Elementy macierzy sztywności:
3 4x k11 = E ⋅ ∫ f11 ( x )dx = E ⋅ ∫ B1 ( x ) ⋅ B1 ( x ) ⋅ A( x )dx = E ⋅ ∫ − + 2 (10 − 0.09 x )dx L L 0 0 0
L L L 2
k11 ≈ E ⋅
L 2 ∑ λi ⋅ f11 ( xi ) . 2 1
Podstawiając za: λ1 = λ2 = 1.0
oraz ξ1 = +0.57735027; ξ 2 = −0.57735027 względnie: L 50 x1 = (ξ1 + 1) = (0.57735207 + 1) = 39.4337567 2 2 x 2 = 10.5662432, otrzymamy: k11≈1265 .
Komentarz [p1]:
Rys.4.11. Model zło ony z dwóch elementów.
19
Pozostałe współczynniki macierzy [K]1 wyznaczone analitycznie zestawiono w tabeli poni ej.
i;j
fij(x1) 0.000062 -0.000922 0.000860 -0.000922 0.013762 -0.012840 0.000860 -0.012840 0.011980
fij(x2) 0.016805 -0.018011 0.001207 -0.018011 0.019305 -0.001293 0.001207 -0.001293 0.000087
∑λ
2
v
v =1
⋅ f ij ( xv )
kij 1265 -1420 155 -1420 2480 -1060 155 -1060 905
1;1 1;2 1;3 2;1 2;2 2;3 3;1 3;2 3;3
0.016867 -0.018933 0.002067 -0.018933 0.033067 -0.014133 0.002067 -0.014133 0.012067
1265 − 1420 155 [K ]1 = − 1420 2480 − 1060. 155 − 1060 905
Według tego samego schematu oblicza się macierz [K]2: 65 635 − 700 [K ]2 = − 700 1040 − 340. 65 − 340 275 Globalna macierz sztywności przyjmie postać: 0 0 u1 R 1265 − 1420 155 − 1420 2480 − 1060 0 0 u 2 0 155 − 1060 1540 − 700 65 ⋅ u 3 = 0 0 − 700 1040 − 340 u 4 0 0 0 0 0 − 340 275 u 5 20
Rozwiązując powy szy układ otrzymamy wektor przemieszczeń: 0 0,018915 {u} = 0,044254 0,084160 0,166320 oraz reakcję R: R = −1420 ⋅ u 2 + 155 ⋅ u 3
} II I
R = −1420 ⋅ 0,018915 + 155 ⋅ 0,044254 = −19,99993 ≈ −20[kN ].
Wyniki:
20
Elementarne wektory przemieszczeń: 0 0.044254 {u}1 = 0,018915 {u}2 = 0,084160 . 0.044254 0.166320 Przebieg przemieszczeń: u(x)=[N]·{u}. Dla elementu 1.
u ( x )1 = N1 ⋅ u1 + N 2 ⋅ u 2 + N 3 ⋅ u 3 u ( x = 0)1 = 0cm u ( x = 25)1 = 0.01892cm u ( x = 50)1 = 0.04426cm
Dla elementu 2.
u ( x )1 = 0 + 6.282 ⋅ 10 −4 ⋅ x + 5.139 ⋅ 10 −6 ⋅ x 2
u ( x )2 = 0.04426 + 7.508 ⋅ 10 − 4 ⋅ x + 3.381 ⋅ 10 −5 ⋅ x 2 u ( x = 0)2 = 0.04426cm u ( x = 25)2 = 0.08416cm u ( x = 50)2 = 0.16632cm
u ( x )2 = N1 ⋅ u1 + N 2 ⋅ u 2 + N 3 ⋅ u 3
Rys.4.12. Przebieg przemieszczeń u(x). Przemieszczenia są określone tylko w węzłach. Powstałe na wykresie punkty są połączone odcinkami prostymi, chocia przybli one rozwiązanie równania ró niczkowego ma postać kwadratową.
21
Przebieg naprę eń:
{σ } = [D ] ⋅ [B ] ⋅ {u} [D ] ≡ E
L
{σ } = E ⋅ − 3 + 4 x 2
L
4 8x − L L2
−
1 4x + ⋅ {u} L L2
Dla elementu 1.
u1 E {σ }1 = 2 ⋅ [− 3L + 4 x 4 L − 8 x − L + 4 x] ⋅ u 2 . L u 3 Dla elementu 2:
{σ }2
u 3 E = 2 ⋅ [− 3L + 4 x 4 L − 8 x − L + 4 x ] ⋅ u 4 . L u 5
Przebieg naprę eń w elemencie jest w tym przypadku linowy.
Rys.4.13. Przebieg naprę eń σ(x).
Naprę enia zostały wyznaczone za pomocą MES w punktach całkowania. Na podstawie tych punktów ekstrapoluje się liniowo wartości naprę eń w węzłach. Ró nica wartości naprę eń na styku elementów pozwala oszacować na ile rozwiązanie jest dobre. Skok wartości naprę eń powinien przyjąć mo liwie najmniejszą wartość. Przez zwiększenie liczby elementów podnosi się dokładność wyników. Na poni szych rysunkach przedstawione są wyniki dla czterech i ośmiu elementów, uzyskane dla kwadratowej postaci przybli onego rozwiązania równania ró niczkowego.
22
Rys.4.14. Przebieg przemieszczeń u(x).
Rys.4.15. Przebieg naprę eń σ(x).
23
Rys.4.16. Przebieg przemieszczeń u(x).
Rys.4.17. Przebieg naprę eń σ(x).
5. Zestawienie wyników
W rozdziale tym przedstawiono i porównano wyniki z ró nych rozwiązań MES w celu wyciągnięcia wniosków końcowych. Tabela poni ej zawiera maksymalne przemieszczenia u(x) w [cm] na prawym końcu pręta:
24
Elementy 1 2 3 4 5 6 7 8 8 o ró nej długości
Aproksymacja linio- Aproksymacja kwaDokładne obliczenia wa dratowa 0,121212 0,145581 0,155437 0,160463 0,163371 0,165199 0,166419 0,167272 0,169490 0,156028 0,166320 0,168859 0,169750 0,170128 0,170310 0,170407 0,170461 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562 0,170562
Na rys. 5.1 widać jak rośnie dokładność wyników wraz ze wzrostem liczby elementów. Tabela poni ej zawiera maksymalne naprę enia σ(x) w [kN/cm2]na końcu pręta: Elementy 1 2 3 4 5 6 7 8 8 o ró nej długości Aproksymacja linio- Aproksymacja kwa- Dokładne obliczenia wa dratowa 3,636 8,511 20,000 6,154 8,000 9,412 10,526 11,429 12,174 12,800 19,110 12,394 14,545 15,878 16,766 17,391 17,849 18,194 20,000 20,000 20,000 20,000 20,000 20,000 20,000 20,000
25
Rys.5.1. Maksymalne przemieszczenia węzłowe.
Rys.5.2. Maksymalne naprę enia elementów.
26
Rozwiązanie zadania w programie ANSYS
■ PREPROCESOR
1. Ustawienia preferencji Main Menu: Preferences Structural (will show) 2. Rysowanie przedmiotu wg rysunku
A
B
A1
D
C
(jedna powierzchnia stworzona z linii i punktów bazowych o współrzędnych jak na rysunku, przy czym aby stworzyć linie AB i CD nale y wskazać punkt punkt B oraz D jako pierwszy
3. Definiowanie typu elementu i opcji Main Menu: Preprocessor → Element Type → Add/Edit/Delete→ Element Types → Add... Wybierz element Quad 4node (PLANE 42)
Main Menu: Preprocessor → Element Type → Add/Edit/Delete→ Element Types → Options... K3: Plane strs w/thk
4. Nadawanie elementowi grubości Main Menu: Preprocessor → Real Constants→ Add.. → OK THK: 1 5. Definiowanie stałych materiałowych Main Menu: Preprocessor → Material Props → -Constant- Isotropic EX: 3000 NUXY:0.29 6. Utwierdzanie przedmiotu Main Menu: Preprocessor → Solution → -Loads- Apply → -Structural- Displacement → On Keypoints Wybierz dwa punkty bazowe znajdujące się po lewej stronie przedmiotu i ustaw pełne utwierdzenie (All DOF)
27
7. Definiowanie obcią enia Main Menu: Preprocessor → Solution → -Loads- Apply → -StructuralForce/Moment → On Keypoints Wybierz dwa punkty bazowe znajdujące się po prawej stronie przedmiotu i ustaw wartość siły FX = 10 (wypadkowa obu sił daje siłę 20kN umieszczoną na środku prawego końca) 8. Ustawienia liczby elementów 8.1 Włączenie numerowania elementów Utility Menu: PlotCtrls → Numbering... Elem / Attrib numbering: Element numbers 8.2 Zapisanie pracy Toolbar: SAVE_DB 8.3 Podział na dwa elementy o równej długości Main Menu: Preprocessor → -Meshing- Shape & Size → -Global- Size Wpisz w polu SIZE wartość nieco większą od 1/2 długości przedmiotu np.: SIZE: 51 8.4 Podział na cztery elementy o równej długości Main Menu: Preprocessor → -Meshing- Shape & Size → -Global- Size Wpisz w polu SIZE wartość nieco większą od 1/4 długości przedmiotu np.: SIZE: 26 8.5 Podział na osiem elementów o równej długości Main Menu: Preprocessor → -Meshing- Shape & Size → -Global- Size Wpisz w polu SIZE wartość nieco większą od 1/8 długości przedmiotu np.: SIZE: 13 9. Dyskretyzacja Main Menu: Preprocessor → -Meshing- Mesh → -Areas- Free 10. Rozwiązanie zadania Main Menu: Solution → -Solve- Current LS 11. Wykres przemieszczeń Main Menu: General Postproc → Path Operations → Define Path + Wska jeden z punktów na prawym i lewym końcu przedmiotu Path Operations → Map ontho Path... Wpisz etykietę: Lab User laber for item: UX Wybierz przemieszczenie względem osi X: Item, Comp: DOF solution; Translation UX Wyłącz uśrednianie wartości na długości elementu: Average results across element: No Path Operations → Plot Path Items... Wybierz przemieszczenia: UX
28
12. Wykres naprę eń Path Operations → Map ontho Path... Wpisz etykietę: Lab User laber for item: naprez Wybierz naprę enia względem osi X: Item, Comp: Stress; X-direction SX Wyłącz uśrednianie wartości na długości elementu: Average results across element: No Path Operations → Plot Path Items... Wybierz naprę enia: naprez 13. Reakcje Utility Menu: List → Results... → Reaction Soultion... All items
By przeprowadzić obliczenia dla większej liczby elementów nale y powtórzyć kroki od punktu 7.4. po uprzednim wczytaniu zbioru: Toolbar: RESUME_DB i odświe eniu ekranu: Utility Menu: Plot → Multi-Plots Dla obliczeń z u yciem ośmiu zagęszczonych elementów nale y podać sposób podziału linii AB i CD: Main Menu: Preprocessor → -Meshing- Shape and Size → -Global- Size w polu SIZE wpisać wartość większą lub równą 10 Main Menu: Preprocessor → -Meshing- Shape and Size → -Lines- Picked Lines wybrać linie AB i CD SIZE: 0 NDIV: 8 SPACE: 8 Main Menu: Preprocessor → -Meshing- Mesh → -Areas- Free
W przypadku u ycia elementów o wy szym stopniu równania ró niczkowego nale y zamiast czerowęzłowego elementu PLANE 42 u yć ośmiowęzłowego elementu PLANE 82 (punkt 2).
29
Ocena wyników
1. Zwiększenie liczby elementów pozwala uzyskać lepsze wyniki. 2. Zagęszczenie elementów w obszarze du ego gradientu naprę enia znacznie poprawia dokładność określenia przemieszczeń, a szczególnie naprę eń 3. Kwadratowa postać przybli onego rozwiązania równania ró niczkowego dostarcza dokładniejsze wyniki, ni liniowa jego postać. 4. Dla kwadratowej postaci przybli onego rozwiązania równania ró niczkowego uzyskuje się większą zbie ność wyników. 5. Na granicach elementów przemieszczenia mają równe wartości, natomiast naprę enia wykazują ró ne wartości. 6. Wartości naprę eń w środkach elementów lub w punktach całkowania, lepiej ni w węzłach zgadzają się z rzeczywistymi wartościami naprę eń. Przemieszczenia mo na określić dokładniej ni naprę enia (naprę enie jest wielkością pochodną). Przez zwiększenie liczby elementów wzrosła dokładność wyników. Rysunek 4.9 pokazuje, e większą dokładność mo na uzyskać przez zagęszczanie elementów w obszarze wy szego gradientu naprę eń. Nale y jeszcze raz podkreślić, e wybrany przykład szczególnie krytycznie podchodzi do obliczeń metodą elementów skończonych. Został on wybrany, aby zaznaczyć aproksymacyjny charakter metody.
30
