Definicja modelu numerycznego

Powszechnie wiadomo, że analizy numeryczne są bardzo wrażliwe na rodzaj zastosowanego elementu kontaktowego oraz gęstość zastosowanej siatki. Do dys-kretyzacji analizowanych belek użyto trójwymiarowych elementów skończonych typu solid.

W celu wybrania odpowiedniego elementu skończonego, który został użyty w dalszych analizach do weryfikacji modelu, przeprowadzono najpierw analizę zbieżności wyników dla różnych elementów skończonych. Analiza ta obejmowała nie tylko gęstość siatki, rodzaj użytego elementu, ale także liczbę węzłów elemen-tów oraz czas obliczeniowy (cpu time). Elementy, które zostały poddane analizie oraz ich charakterystyki, zostały pokazane w tabeli 4.2.

Jako rozwiązanie porównawcze zastosowano analityczne rozwiązanie belki obustronnie utwierdzonej, którą obciążono w środku rozpiętości siłą skupioną P (Rysunek 4.1).

Rysunek 4.1 – Konfiguracja początkowa belki w analizie zbieżności

Do obliczeń przyjęto belkę o przekroju kołowym, o promieniu przekroju r = 100 i długości 2000. Przyjęto następujące dane materiałowe: moduł Younga E = 2500 oraz współczynnik Poissona v = 0.4. Belka została obciążona siłą skupio-ną o wartości P = 1000. Charakterystyka geometryczna i materiałowa belki spełnia założenia teorii belki Eulera – Bernoullego. Wyniki rozwiązania analitycznego przedstawione zostały w tabeli 4.1.

Założono trzy kryteria weryfikacji: przemieszczenie w punkcie B, naprężenie normalne w punkcie B oraz czas obliczeń procesora (cpu), w zależności od gęstości siatki wykorzystanej w modelu. Wyniki przeprowadzonych analiz zostały pokaza-ne na wykresach (rysupokaza-nek 4.2 – 4.4).

Tabela 4.1 Rozwiązanie analityczne Przemieszczenie

u_B

Naprężenie σB

Rozwiązanie

analityczne 3.3953 2.5465

Tabela 4.2 Elementy skończone poddane analizie zbieżności Nr Nazwa

elementu

Typ elementu

Funkcja kształtu

Liczba węzlów

Całkowanie

zredukowane Wielkość siatki

1 C3D8

brick

liniowa 8 Nie 25 / 20 / 10 / 5

2 C3D8R liniowa 8 Tak 25 / 20 / 10 / 5

3 C3D20 kwadratowa 20 Nie 25 / 20 / 10

4 C3D20R kwadratowa 20 Tak 25 / 20 / 10

5 C3D4

tetrahedron

liniowa 4 Nie 25 / 20 / 10 / 5

6 C3D10 kwadratowa 10 Nie 20 / 20 / 10

Rysunek 4.2 Przemieszczenie punktu środkowego B w analizie zbieżności

Rysunek 4.3 Naprężenie normalne w punkcie środkowym B w analizie zbieżności 1.6

Rysunek 4.4 Czas obliczeniowy w analizie zbieżności

Na wykresie, na rysunku 4.2 pokazane zostały przemieszczenia punktu B, ja-kie uzyskano przy użyciu różnych elementów skończonych w zależności od za-gęszczenia siatki elementów skończonych. Wykres na rysunku 4.3 pokazuje war-tość naprężeń normalnych w punkcie B uzyskaną w analizach zaś, wykres na ry-sunku 4.4 – czas potrzebny do jej wykonania. Można zauważyć, że wszystkie za-stosowane elementy, jeżeli są wystarczająco małe, dają bardzo dobrą zgodność wyniku z rozwiązaniem analitycznym, które na wykresie oznaczono linią przery-waną.

Przy wyborze elementu ważne było także, aby wyniki analizy były jak naj-mniej wrażliwe na rozmiar przyjętej siatki. Analizując wyniki można zauważyć, że elementy z liniową funkcją kształtu: czworościenny element C3D4 oraz elementy C3D8 i C3D8R, nie spełniają powyższego wymogu. Dwa pierwsze wykresy (rysu-nek 4.2 i 4.3) pokazują, że elementy skończone o kwadratowych funkcjach kształtu (C3D10 i C3D20, C3D20R) nie są zbyt wrażliwe na rozmiar siatki i są one z tego punktu widzenia najkorzystniejsze.

Ostatni wykres przedstawia czas, jaki potrzebny jest do wykonania analizy, czyli tzw. koszt obliczeń. Kryterium te okazało się bardzo istotne, gdyż wiele prze-prowadzonych symulacji komputerowych okazało się bardzo czasochłonnych.

Biorąc pod uwagę te trzy kryteria do przeprowadzenia analiz wybrano jako najkorzystniejsze dwa skończone elementy: C3D10 oraz C3D20R.

0

4.3 Przykłady numeryczne

Ponieważ użycie w analizach elementu C3D10, który jest elementem ogólne-go przeznaczenia, daje bardzo dobre wyniki i dokładność jest wystarczająca w przypadku wartości sił normalnych w miejscu kontaktu, został on użyty w pierwszym etapie weryfikacji, który miał na celu sprawdzić prawidłowość zasto-sowanych modeli (Przykład 1 – Przykład 3). W kolejnych analizach (Przykład 4) skupiono się na sprawdzeniu dokładności obliczeń deformacji przekroju i w tych analizach posłużono się elementem dwudziestowęzłowym C3D20R, który umożli-wia otrzymanie bardziej dokładnych wyników niż C3D10 przy tej samej liczbie stopni swobody i wydaje się najlepszy do liniowych obliczeń sprężystych.

Przykład 1

W tym przykładzie rozpatrywany jest kontakt między dwiema obustronnie utwierdzonymi belkami, których konfiguracja początkowa została pokazana na rysunku 4.5. Belki ułożone są w taki sposób, że osie belek są do siebie prostopadłe, a odległość początkowa między belkami wynosi 0.001. Końce belki 2 podlegają przemieszczeniom pionowym o wartości ∆ = 0.1, które przyłożono w 40 równych przyrostach. Belki mają przekrój kołowy o równych promieniach r1 = r2 = 0.1, a ich długość wynosi 6.0. Do obliczeń przyjęto jednakowe dane materiałowe dla obu belek: moduł Younga E = 210· 10⁹, współczynnik Poissona v = 0.3. W przy-padku analiz przeprowadzonych przy pomocy elementów belkowych zapropono-wanych w pracy każdą z belek podzielono na 20 jednakowych elementów. W ana-lizach przy użyciu programu Abaqus, każdą z belek podzielono na 91006 elemen-tów trójwymiarowych C3D10 otrzymując 8 elemenelemen-tów przypadających na wyso-kość przekroju.

Analizę przeprowadzono dla dziewięciu różnych konfiguracji położenia belki drugiej. Dla pierwszej konfiguracji belki umieszczono symetrycznie względem siebie, następnie zmieniano położenie belki drugiej przesuwając ją za każdym ra-zem o wartość δ = 0.3 w lewa stronę, jak pokazano na rysunku 4.5. Wyniki obli-czeń zostały przedstawione w tabelach 4.3 – 4.5. Wartości sił normalnych dla dziewięciu różnych konfiguracji położenia belki 2 zostały przedstawione w tabeli 4.3. W kolumnie drugiej podano końcowe wartości siły normalnej otrzymane w analizie z użyciem skończonych elementów belkowych uwzględniajacych deformację przekroju: Model 1, Model 2 oraz modelu, który jej nie uwzględnia, a w kolumnie trzeciej wyniki otrzymane z obliczeń numerycznych wykonanych przy użyciu programu Abaqus. Ostatnia kolumna zawiera różnice procentowe między siłami otrzymanymi przy uzyciu modeli uwzględniających deformację i siłami otrzymanymi w programie Abaqus. W tabeli 4.4 pokazano otrzymane war-tości deformacji przekroju dla symetrycznego położenia belek, tzn. gdy δ = 0

Rysunek 4.5 Przykład 1 – Konfiguracja początkowa osi belek

Rysunek 4.6 Przykład 1 – Konfiguracja odkształcona osi belek dla δ = 1.5

Na rysunku 4.6 przedstawiono konfigurację odkształconą osi belek dla przy-padku, gdy przesunięcie belki 2 względem pierwszej wynosi δ = 1.5.

W przypadku siły normalnej otrzymano bardzo dobrą zgodność wyników dla wszystkich dziewięciu położeń belki 2. Różnice między wartościami sił normal-nych otrzymanormal-nych przy użyciu Modelu 1 i Modelu 2 oraz wartościami uzyskanymi podczas analizy z użyciem programu Abaqus nie przekracza 2.5 %. W przypadku wartości deformacji d dla symetrycznego położenia belek, czyli δ = 0 zmierzone deformacje są wartościami tego samego rzędu (patrz tabela 4.4). Natomiast w przypadku niesymetrycznych układów tzn., gdy δ > 0 różnice w zmierzonych deformacjach są tak duże, że są one nieporównywalne. Wynika to ze sposobu defi-niowania punktów kontrolnych w przypadku analiz przeprowadzonych za pomocą Abaqus/CAE. Punkty kontrolne A, B przedstawione na rysunku 4.7 zdefiniowane zostały w konfiguracji początkowej belek, czyli przed odkształceniem.

W przypadku symetrycznego położenia belek punkty kontrolne pokrywają się z punktami kontaktowymi w konfiguracji odkształconej. Dla niesymetrycznych układów punkty te nie pokrywają się z rzeczywistymi punktami styku, co uniemoż-liwia wyznaczenie rzeczywistej deformacji przekroju, która wyznaczona został na podstawie przemieszczeń punktów kontrolnych A i B. Deformacja d obliczana jest, jako różnica miedzy przemieszczeniami punktów kontrolnych dla nieodkształco-nych i odkształconieodkształco-nych konfiguracji belek. Konfiguracja belek przed i po odkształ-ceniu wraz z zaznaczonymi punktami kontrolnymi została przedstawiona na rysun-ku 4.7.

Możemy zapisać, że:

B A

AB u r gap r u

u′ + =2⋅ + +2⋅ + (4.1)

Deformacja przekroju rozumiana, jako skrócenie promieni belek jest oblicza-na z zależności:

uAB

r

d = 4⋅ − ′ (4.2)

Rysunek 4.7 Przykład 1 – Położenie punktów kontrolnych dla odkształconych i nieod-kształconych belek

Tabela 4.3 Przykład 1 – Wartość siły normalnej Fn

Położenie belki 2

δ

Beam-to-beam

Abaqus

Różnica Model

bez deformcji Model 1 Model 2 Model 1 /Abaqus

Model 2 /Abaqus 0.0

0.3

778000.73 790019.50

772599.83 784474.00

773266.24 785167.67

760881.0 772169.0

1.54 % 1.59 %

1.63 % 1.68 % 0.6 827423.28 821050.40 821809.95 818256.0 0.34 % 0.43 % 0.9 891905.63 885091.58 892470.60 881160.0 0.45 % 1.28 % 1.2 989118.01 981131.56 982076.59 975741.0 0.55 % 0.65 % 1.5 1123948.79 1113632.64 1114874.23 1106910.0 0.61 % 0.72 % 1.8 1297335.08 1283809.65 1285520.15 1276190.0 0.60 % 0.73 % 2.1 1496822.50 1478721.28 1481229.91 1471200.0 0.51 % 0.68 % 2.4 1680115.57 1657221.07 1660364.75 1651970.0 0.32 % 0.51 %

Tabela 4.4 Przykład 1 – Wartość deformacji d dla δ = 0

Dla symetrycznego położenia belek (δ = 0) wykonano także dodatkowe obli-czenia przyjmując inne wartości modułu Younga oraz współczynnika Poissona. W zadaniu przyjęto dane: E =250· 10⁷ oraz ν = 0.4. Wartość siły normalnej oraz de-formacji zostały przedstawione w tabeli 4.5. Jak poprzednio otrzymano dobrą zgodność wyników dla sił normalnych oraz deformacji.

Tabela 4.5 Przykład 1 – Wyniki analizy dla E =250·10⁷ i ν = 0.4

W tym przykładzie analizowany jest kontakt między dwiema belkami wspor-nikowymi, których konfiguracja początkowa została przedstawiona na rysunku 4.8.

Belki mają przekrój kołowy o promieniach: r1 = r2 = 0.1, długość belek wynosi 6.0, odległość początkowa między belkami 0.001. Do obliczeń przyjęto: wartość modu-łu Younga E = 250·10⁵, współczynnik Poissona v = 0.3. Każdą z belek podzielono na 20 jednakowych elementów a w przypadku analiz w programie Abaqus, doko-nano podziału na 91006 elementów trójwymiarowych C3D10 otrzymując 8 ele-mentów przypadających na wysokość przekroju. Wolny koniec belki drugiej pod-dany jest przemieszczeniom pionowym ∆, przykładanym w 60 równych przyro-stach o rożnej wartości dla kolejnych konfiguracji położenia belek. Wartości tych przemieszczeń wraz z przyjętymi współczynnikami kary w przypadku analizy z użyciem Modelu 1 zostały przedstawione w tabeli 4.6. Przeprowadzono oblicze-nia dla trzynastu różnych konfiguracji położeoblicze-nia osi belek. Zmieoblicze-niano położenie obu belek przesuwając je za każdym razem o wartość δ = 0.3 od konfiguracji po-czątkowej, w kierunku pokazanym na rysunku 4.8.

Wyniki analizy zostały przedstawione w tabelach 4.6 – 4.8. W tabeli 4.6 po-kazano końcowe wartości siły normalnej w miejscu kontaktu, jakie uzyskano w analizach przy użyciu Modelu 1, Modelu 2 oraz programu Abaqus. Jak w przy-padku poprzednim uzyskano bardzo dobrą zgodność wyników. Różnica pomiędzy ich wartościami nie przekracza 3.1%.

Rysunek 4.8 Przykład 2 – Konfiguracja początkowa osi belek

Rysunek 4.9 Przykład 2 – Konfiguracja odkształcona osi belek dla δ = 3.0

Tabela 4.6 Przykład 2 – Wartość siły normalnej Fn

Położenie belki 2

x

Przemie- szczenie

∆

ε

^N Beam-to-beam

Abaqus/CEA

Różnica

Model 1 Model 2 Model 1

/Abaqus

Model 2 /Abaqus 1.8 0.09 1·10⁷ 5.80803 5.79701 5.63627 3.05% 2.85%

2.1 0.09 1·10⁷ 5.25181 5.25211 5.14890 2.00% 2.00%

2.4 0.09 1·10⁷ 4.79779 4.79810 4.72851 1.47% 1.47%

2.7 0.1 1·10⁷ 4.93046 4.93012 4.88361 0.96% 0.95%

3 0.1 1·10⁷ 4.56138 4.55994 4.52692 0.76% 0.73%

3.3 0.1 1·10⁷ 4.23675 4.23686 4.22213 0.35% 0.35%

3.6 0.1 1·10⁷ 3.94825 3.94914 3.93403 0.36% 0.38%

3.9 0.09 1·10⁴ 3.28588 3.31114 3.29694 0.34% 0.43%

4.2 0.09 1·10⁴ 3.07698 3.09639 3.08498 0.26% 0.37%

4.5 0.08 1·10⁴ 2.55464 2.56805 2.55779 0.12% 0.40%

4.8 0.07 1·10⁴ 2.08615 2.09535 2.08630 0.01% 0.43%

5.1 0.06 1·10⁴ 1.66632 1.67256 1.66470 0.10% 0.47%

5.4 0.05 1·10⁴ 1.29119 1.29533 1.28851 0.21% 0.53%

Dodatkowo w tym przykładzie dla trzech wybranych konfiguracji pokazano przemieszczenia wolnych końców belek, jakie otrzymano w analizach. Można zauważyć, że obliczone przemieszczenia pozostają w dobrej zgodności miedzy porównywanymi metodami obliczeniowymi. Różnica nie przekracza kilku

punk-tów procentowych i wynosi maksymalnie 4.9%. Ze względu na niesymetryczne położenie belek wielkość deformacji przekroju d nie była porównywana. Na ry-sunku 4.9 przedstawiono konfigurację odkształconą osi belek dla δ = 3.0, na któ-rym zaznaczono także elementy kontaktowe.

Tabela 4.7 Przykład 2 – Wartości przemieszczenia wolnych końca belek dla Modelu 1 Położenie

belki 2 x

Model 1 Abaqus/CEA Różnica

Belka 1 Belka 2 Belka 1 Belka 2 Belka 1 Belka 2 1.8 0.03088 0.09 0.029766 0.09 3.73% 0 3.0 0.06220 0.10 0.616966 0.10 0.82% 0 4.5 0.07053 0.08 0.0704139 0.08 0.16% 0

Tabela 4.8 Przykład 2 – Wartości przemieszczenia wolnych końca belek dla Modelu 2 Położenie

belki 2 x

Model 2 Abaqus/CEA Różnica

Belka 1 Belka 2 Belka 1 Belka 2 Belka 1 Belka 2 1.8 0.03087 0.09 0.029766 0.09 3.72% 0 3.0 0.06220 0.10 0.616966 0.10 0.82% 0 4.5 0.07089 0.08 0.0704139 0.08 0.67% 0

Przykład 3

W tym przykładzie rozpatrywany jest kontakt miedzy czterema identycznymi belkami. Belki maja przekrój kołowy o promieniu r = 0.1, a ich długość wynosi 8.0. Do obliczeń przyjęto moduł Younga E = 250·10⁷ oraz współczynnik Poissona v = 0.3. Wolne końce belek poddano przemieszczeniom pionowym równym

∆ = 0.5 przyłożonym w 30 identycznych przyrostach.

W każdej belce punkty środkowe belki zostały częściowo utwierdzono umoż-liwiając jedynie obrót względem osi prostopadłej do osi belki, odpowiednio x lub y.

Konfiguracja początkowa osi belek została pokazana na rysunku 4.10. Współczyn-nik kary w przypadku analizy z użyciem modelu pierwszego

ε

^{n = 2·106}. W przy-padku analiz przeprowadzonych przy pomocy elementów belkowych każdą z belek podzielono na 10 jednakowych elementów. W analizach przy użyciu programu Abaqus, dokonano dyskretyzacji belek elementami C3D10, w taki sposób, aby otrzymać 8 elementów na wysokość poprzeczną belki.

Rysunek 4.10 Przykład 3 – Konfiguracje początkowa osi belek

Jak w poprzednich przykładach otrzymano zgodne wartości siły normalnej, które przedstawiono w tabeli 4.9. Także i w tym zadaniu ze względu przemiesz-czanie się punktów kontrolnych, nie można było określić deformacji przekroju na podstawie ich przemieszczeń i nie zostały one porównane. Konfiguracja odkształ-conych osi belek została pokazana na rysunku 4.11, zaś na rysunku 4.12 pokazany został widok odkształconych belek uzyskany podczas analizy w programie Abaqus.

Rysunek 4.11 Przykład 3 – Konfiguracja odkształcona osi belek

Rysunek 4.12 Przykład 3 – Konfiguracja odkształcona belek

Tabela 4.9 Przykład 3 – Wartość siły normalnej Fn Beam-to-beam

Abaqus/CEA

Różnica

Model 1 Model 2 Model 1

/Abaqus

Model 2 /Abaqus Siła

no-malna Fn 25858.12 26933.00 26369.00 1.98% 2.14%

Przykład 4

W tym zadaniu skupiono się jak sprawdzeniu poprawności obliczanej według Modelu 1 i Modelu 2 deformacji przekroju w strefie kontaktu. W celu zweryfiko-wania otrzymywanych wyników rozwiązano zadanie, w którym rozpatrywane były dwie, symetrycznie ułożone belki. Obie belki były obustronnie utwierdzone i miały jednakową długość równą 6.0.

Konfiguracja początkowa osi belek jest identyczna jak dla belek omawianych w przykładzie pierwszym niniejszego rozdziału (patrz Rysunek 4.5). Belki ułożone są w taki sposób, że osie belek są do siebie prostopadłe a odległość początkowa między belkami wynosi 0.001. Końce belki 2 podlegają przemieszczeniom piono-wym o wartości ∆ = 0.1, które w przypadku zaproponowanych modeli przykładano w 40 równych przyrostach. Każda z belek została podzielona na 20 równych ele-mentów. W przypadku analizy przeprowadzanej przy użyciu programu Abaqus zastosowano dwudziestowęzłowy element C3D20R, który jest elementem o dużo większej dokładności.

Tabela 4.10 Przykład 4 – Zestawienie danych materiałowych, geometrycznych oraz koń-cowej wartości siły normalnej Fn

Dane Wartość siły normalnej Fn

Nr E v r₁ r₂

W ramach tego przykładu przygotowano dziesięć zadań, które różniły się da-nymi materiałowymi oraz przyjętymi promieniami belek. Zmianie ulegał moduł Younga oraz współczynnik Poissona. Zadania te zostały odpowiednio oznaczone, jako: A lub B. Zestawienie zadań wraz z ich numeracją przedstawiono w tabeli 4.10. W ramach każdego zadania sprawdzano wartość deformacji przekroju dla różnej gęstości siatki w programie Abaqus, sprawdzając w ten sposób zbieżność wyników z wartościami deformacji uzyskanymi przy zastosowaniu zaproponowa-nych modeli elementów belkowych.

W tabeli 4.11 dla wybranych gęstości siatek zestawiono liczbę elementów oraz liczbę węzłów uzyskaną w poszczególnych zadaniach w przypadku analiz wykonywanych przy pomocy programu Abaqus oraz przy zastosowaniu skończo-nych elementów kontaktowych zaproponowaskończo-nych w pracy. W pierwszej kolumnie podano rozmiar siatki (tzw. wielkość oczek) zastosowanej do podziału belek na elementy, w kolejnych ilość: elementów w zadaniu, węzłów oraz elementów przy-padająca na przekrój poprzeczny belki. Choć bardziej istotna jest ilość elementów, na jaką podzielono belki oraz ilość elementów przypadająca na przekrój poprzecz-ny belki, w celu uproszczenia zapisu w dalszej części rozważań posługiwano się tzw. rozmiarem siatki.

Tabela 4.11 Przykład 4 – Liczba elementów oraz węzłów w przykładowych zadaniach, w zależności od gęstości nałożonej siatki

Rozmiar

siatki Liczba

Beam-to beam Abaqus

A0/B0 - A6/B6 A0/B0 - A3/B3 A4/B4 A5/B5

25

elementów 40 28800 42240 19200

węzłów 42 133850 192110 90998

elem. / wysokość 1 8 12 i 8 8 i 4

20

elementów 40 105840 94800 38400

węzłów 42 471154 413822 176834

elem. / wysokość 1 12 16 i 10 10 i 6

15

elementów 40 144000 206400 89600

węzłów 42 626762 881638 396734

elem. / wysokość 1 14 20 i 14 14 i 8

10

elementów 40 460800 672000 276000

węzłów 42 1963234 2809090 1187070

elem. / wysokość 1 20 30 i 20 20 i 10

Przeglądając tabele 4.11 można zauważyć, że w przypadku analizy z wyko-rzystaniem zaproponowanych w pracy kontaktowych elementów belkowych roz-miar zadania jest nieporównywalnie mniejszy, niż zadania obliczanego przy użyciu elementów typu solid w programie Abaqus, co stanowi główną zaletę proponowa-nych modeli. Przykładowo dla zadań, w których analizowano belki o jednakowych promieniach (zadania A0 – A3, B0 – B3) podzielenie belek na 20 równych elemen-tów w analizie Beam-to-beam daje zgodne wyniki z wartościami, jakie uzyskano w analizie z użyciem pełnych elementów C3D20R, w której każdą z belek należało podzielić aż na 14400 elementów przy rozmiarze siatki 25 (patrz tabela 4.11).

Na rysunku 4.13 – 4.15 pokazano przykładową dyskretyzacje belek oraz czę-ściowe wyniki, jakie uzyskano w Abaqusie, w zadaniu A3. Na rysunku 4.13 przed-stawiono podział belek na elementy skończone C3D20R, jakie otrzymano w pro-gramie Abaqus/CAE przy zastosowaniu siatki o rozmiarze 20 oraz 10. Na rysunku tym zaznaczono także punkty kontrolne A i B, które służyły do wyznaczenia warto-ści deformacji przekroju w strefie kontaktu.

Rozkład końcowych przemieszczeń w belkach pokazano na rysunku 4.14 zaś rozkład naprężeń – na rysunku 4.15. W celu pokazania koncentracji naprężeń, któ-re występują w miejscu kontaktu, rozkład ten został pokazany na przekroju po-przecznym belek.

Rysunek 4.13 Przykład 4 – Podział belek na elementy skończone dla siatki 20 i 10 w zada-niach A3 dla konfiguracji odkształconej

Rysunek 4.14 Przykład 4 – Przemieszczenie belek w zadaniu A3 przy siatce 20

Rysunek 4.15 Przykład 4 – Rozkład naprężeń w belkach w zadaniu A3 przy siatce 20 – przekrój poprzeczny

Wyniki przeprowadzonych analiz zostały przedstawione na wykresach oraz w tabelach 4.10 – 4.13. W tabeli 4.10 zestawiono końcowe wartości siły normalnej w miejscu kontaktu uzyskane w analizach przy użyciu modelu pierwszego, drugie-go oraz programu Abaqus przy zastosowaniu siatki o rozmiarze 25. Można zauwa-żyć, że otrzymano bardzo dobrą zgodność wyników. Różnice miedzy końcowymi wartościami sił nie przekraczają 2%. Przebieg przyrostu siły normalnej w czasie trwania procesu dla Modelu 2 oraz dla analizy z wykorzystaniem elementów C3D20R pokazano na wykresach (Rysunek 4.13 i 4.14). Ze względu na bardzo zbliżone wartości sił dla obu modeli (Modelu 1 i Modelu 2), dla modelu pierwsze-go przebiegu nie pokazano.

Rysunek 4.16 Przykład 4 – Zmiany siły normalnej dla zadań A0 – A5

Rysunek 4.17 Przykład 4 – Zmiany siły normalnej dla zadań B0 – B5 0

200000 400000 600000 800000 1000000 1200000 1400000

0 20 40 60 80 100

Siła normalna Fn

Procent końcowego przemieszczenia ∆ [%]

A0-A3 Model 2 A4 Model 2 A5 Model 2 A0-A3 Abaqus A4 Abaqus A5 Abaqus

0 2000 4000 6000 8000 10000 12000 14000 16000

0 20 40 60 80 100

Siła normalna Fn

Procent końcowego przemieszczenia ∆[%]

B0-B3 Model 2 B4 Model 2 B5 Model 2 B0-B3 Abaqus B4 Abaqus B5 Abaqus

Rysunek 4.18 Przykład 4 – Zmiany deformacji przekroju dla zadań A0 – A5

Rysunek 4.19 Przykład 4 – Zmiany deformacji przekroju dla zadań B0 – B5 0

Procent końcowego przemieszczenia ∆[%]

A0 Model 2

Procent końcowego przemieszczenia ∆[%]

B0 Model 2

Wartości deformacji przekroju w strefie kontaktu, jakie otrzymano w anali-zach przy użyciu obu modeli zostały przedstawione w tabeli 4.12. Zostały one ze-stawione wraz z wartościami deformacji, jakie uzyskano z programu Abaqus przy dyskretyzacji belek siatką o rozmiarze 25. W kolumnie ostatniej pokazano różnice procentowe pomiędzy otrzymanymi wartościami. W przypadku zastosowania siatki o gęstości 25 różnice te są znaczne i w skrajnych przypadkach osiągają wartość 35%.

W związku z tym, że dokładność obliczeń w metodzie elementów skończo-nych w znacznej mierze zależy od liczby elementów w dyskretyzacji, przeprowa-dzono kolejne symulacje numeryczne w Abaqusie, w których zwiększano gęstość nakładanej na belki siatki. Zagęszczanie siatki powoduje otrzymywanie dokład-niejszych wyników jednak znacznie zwiększa rozmiar zadania, a co za tym idzie 1 – czas obliczeniowy. W przypadku zadań, w których użyto Model 1 i Model 2 wykorzystywano własne, autorskie programy obliczeniowe. W przypadku użycia proponowanych skończonych elementów belkowych czas obliczeniowy nie prze-kraczał 1 minuty, natomiast w przypadku analiz przy użyciu pełnych elementów dwudziestowęzłowych w programie Abaqus, czas ten wynosił od kilkudziesięciu minut dla siatki o rozmiarze 20, do kilkudziesięciu godzin dla siatek o rozmiarze mniejszym niż 5.

Tabela 4.12 Przykład 4 – Deformacje przekroju w miejscu kontaktu

Nr

Beam-to-beam

Abaqus (siatka 25)

Różnica

Model 1 Model 2 Model 1

/Abaqus

Model 2 /Abaqus A0 0.00065015 0.00066851 0.00053100 22.44% 25.90%

A1 0.00062265 0.00065500 0.00053200 17.04% 23.12%

A2 0.00062240 0.00063216 0.00052700 18.10% 19.95%

A3 0.00058095 0.00059946 0.00051700 12.37% 15.95%

A4 0.00084662 0.00086501 0.00077400 9.38% 11.76%

A5 0.00022505 0.00023171 0.00033700 -33.22% -31.24%

B0 0.00065604 0.00066851 0.00053100 23.55% 25.90%

B1 0.00064279 0.00065500 0.00053200 20.82% 23.12%

B2 0.00062037 0.00063216 0.00052700 17.72% 19.95%

B3 0.00058828 0.00059946 0.00051700 13.79% 15.95%

B4 0.00080379 0.00082036 0.00075900 5.90% 8.08%

B5 0.00021350 0.00021970 0.00033900 -37.02% -35.19%

W tabeli 4.13 przedstawiono różnice procentowe pomiędzy wartościami de-formacji w zależności od przyjętego rozmiaru siatki, a co za tym idzie – liczbie elementów przypadających na wysokość belki.

Ponieważ w obu modelach belkowych otrzymano podobne wartości deforma-cji przekroju w strefie kontaktu, w tabeli 4.13 pokazano różnice procentowe, jakie otrzymano w analizach dla Modelu 2 i otrzymane w programie Abaqus/CAE.

W kolumnie pierwszej pokazano deformację przekroju, jaką uzyskano używa-jąc Modelu 2, a w kolejnych kolumnach różnicę miedzy tą wartością, a wartościa-mi uzyskanywartościa-mi przy użyciu programu Abaqus dla różnych rozwartościa-miarów siatek. Po-kazane niżej wyniki stanowią tylko część wyników symulacji, jakie przeprowadzo-no. Dla analizowanego zadania przeprowadzono także symulacje numeryczne przy użyciu gęstszych siatek o rozmiarze: 7, 5 z 40 elementami na wysokości jednak ze względu na nakładanie się błędów numerycznych otrzymywano wyniki znacznie różniące się od wyników otrzymywanych przy użyciu siatek o rozmiarze 10 i rzad-szych. W związku z tym zostały one odrzucone, jako błędne i nie były brane pod uwagę w dalszych analizach.

Przebieg zmiany wartości deformacji w zależności od siatki pokazano na ry-sunku 4.17. Można zauważyć, ze wraz ze wzrostem liczby elementów różnice te maleją i są bliskie rozwiązaniom otrzymywanym przy użyciu Modelu 2. Najmniej Tabela 4.13 Przykład 4 – Różnica procentowa między wartościami deformacji w zależności od gęstości siatek

Beam-to-beam Abaqus

(rozmiar siatki)

Model 2 25 20 17 16 15 10

A0 0.00066851 25.90% 17.08% 17.08% 1.91% -1.55% 5.61%

A1 0.00065500 23.12% 14.71% 2.83% 0.31% -3.11% 4.80%

A2 0.00063216 19.95% 12.28% 5.71% -1.38% -4.65% 2.46%

A3 0.00059946 15.95% 9.59% 3.71% -3.00% -6.19% 0.58%

A4 0.00086501 11.76% 0.58% 0.58% -7.88% -3.89% -

A5 0.00023171 -31.24% 7.27% 7.27% 7.27% 7.27% 5.32%

B0 0.00066851 25.90% 17.08% 4.62% 1.91% -1.55% 4.78%

B1 0.00065500 23.12% 14.71% 2.83% 0.31% -3.11% 3.80%

B2 0.00063216 19.95% 10.71% 5.71% -1.38% - -

B3 0.00059946 15.95% 9.59% 3.71% -3.00% -6.19% 0.58%

B4 0.00082036 8.08% -1.16% -10.05% -10.64% -9.35% 1.15%

B5 0.00021970 -35.19% -11.77% 1.24% 1.71% 1.71% 1.71%

Rysunek 4.20 Przykład 4 – Przebieg zmiany różnic procentowych między wartościami deformacji w zależności od gęstości siatek

korzystne wyniki otrzymano w zadaniu A5 i B6, gdzie różnice przy siatce o roz-miarze 25 wynoszą ponad 30 %. Błąd ten wynika z nałożenia na obie belki jedna-kowej siatki w wyniku, czego dla belki o mniejszym promieniu r2= 0.05 otrzymano podział belki na wysokości na zaledwie 4 elementy. Dla tego zadania przyjęcie siatki gęstszej (20 i mniej) daje bardziej korzystne wyniki, które zostały zaprezen-towane także na wykresie. Na wykresie nie przedstawiono skrajnych wyników, które zostały odrzucone ze względu na występujące w analizach niestabilności numeryczne. Wartości te zostały w tabeli 4.13 przekreślone.

Przedstawione wyżej wyniki pokazują, iż zaproponowane modele skończo-nych elementów belkowych dają poprawne wyniki. Różnice w przypadku wartości sił normalnych w strefie kontaktu otrzymywane przy wykorzystaniu obu modeli oraz sił obliczonych w programie Abaqus są mniejsze niż 2 %. Także w przypadku wartości deformacji przekroju otrzymywane wyniki są poprawne. Przy dobrze dobranym typie elementu oraz rozmiarze siatki maleją one nawet do 0.58 %. Poja-wiające się większe różnice w otrzymywanych wartościach wynikają albo z błę-dów numerycznych, albo z niestabilności numerycznej. Można, więc stwierdzić, iż proponowane modele belkowe dobrze odzwierciedlają zachowanie belek będących w kontakcie.

-15%

-5%

5%

15%

25%

10 15 20 25

Róznica wartości deformacji

Róznica wartości deformacji

W dokumencie Olga Kawa (Stron 79-102)