33, s. 107-112, Gliwice 2007
WYBRANE PRZYKŁADY MODELOWANIA ZAGADNIENIA SPRZĘŻENIA W ŚRODOWISKU ALE
JERZY MAŁACHOWSKI, JACEK ŁAZOWSKI
Katedra Mechaniki i Informatyki Stosowanej, Wojskowa Akademia Techniczna,
e-mail: jmalachowski@wme.wat.edu.pl, jlazowski@wme.wat.edu.pl
Streszczenie. Uzyskanie wyników o odpowiedniej jakości w przypadku komputerowej analizy układów poddanych oddziaływaniu ciśnienia spalanych gazów, tj. w warunkach gwałtownego wzrostu obciążenia, wymaga modelowania procesu interakcji pomiędzy ciałami stałymi a ośrodkiem płynnym (w tym przypadku gazem). Prawidłowe wykonanie zadania obciążania wszystkich elementów modelu numerycznego MES wymaga zastosowania sprzężenia opisu Eulera (gaz) i Lagrange’a (ciała stałe). W pracy zaproponowano dwa przypadki takiej analizy: koncepcję numerycznego modelowania procesu spalania ładunku miotającego, w wyniku którego jest generowane ciśnienie oddziaływające na pocisk oraz drugi przypadek, w którym numerycznie realizuje się proces interakcji gazu z ciałem stałym na przykładzie elementu rury poddanego oddziaływaniu silnemu impulsowi dynamicznemu, który powstaje w wyniku detonacji materiału wybuchowego .
1. WSTĘP
Uzyskanie wyników o odpowiedniej jakości w przypadku komputerowej analizy układów poddanych oddziaływaniu ciśnienia spalanych gazów, tj. w warunkach gwałtownego wzrostu obciążenia, wymaga modelowania procesu interakcji pomiędzy ciałami stałymi a ośrodkiem płynnym (w tym przypadku gazem). Prawidłowe wypełnienie zadania obciążania wszystkich elementów modelu numerycznego MES wymaga zastosowania sprzężenia opisu Eulera (gaz) i Lagrange’a (ciała stałe). W przypadku opisu Eulera ciała przemieszczają się na tle siatki, którą został opisany dany obszar, zaś w opisie Lagrange’a siatka odkształca się wraz z deformowanym ciałem. Szczególnym zagadnieniem związanym z gwałtownym wzrostem obciążenia działającego na układ lufa-nabój jest modelowanie procesu spalania materiału wybuchowego miotającego w lufie [3]. Rozwiązanie tego problemu i zastosowanie w opracowywanym modelu numerycznym układu lufa-nabój w pakiecie LS-DYNA opartym na metodzie elementów skończonych, może spowodować w dużym stopniu wyeliminowanie długotrwałych i kosztownych badań poligonowych.
Ostatnia sytuacja zarówno w kraju jak też i na świecie pokazuje, iż działań o charakterze terrorystycznym praktycznie możemy się spodziewać w każdej chwili, z uwagi na łatwość dostępu do instalacji służących do transportu gazu i paliw płynnych. Wielokrotnie też dochodzi do przypadków samozapłonów lub też eksplozji na skutek oddziaływań zewnętrznych.
W przypadku symulacji rozważanych zagadnień problemem do rozwiązania jest dobór
odpowiednich modeli konstytutywnych dla materiałów, które są wykorzystywane w budowie rurociągów i elementów ich infrastruktury. W pracy autorzy prezentują tematykę procesu interakcji gazu z ciałem stałym na przykładzie elementu rury poddanego oddziaływaniu silnego impulsu dynamicznego, który powstaje w wyniku detonacji materiałów wybuchowych [2,5].
Na tle tych dwóch przykładów będzie omówiona procedura numerycznego sprzężenia pomiędzy układem opisanym we współrzędnych Lagrange’a, a falą ciśnienia generowaną we współrzędnych Eulera, czyli tzw. ALE (ang. Arbitrary Lagrangian-Eulearian) [2,4].
2. PODSTAWY TEORETYCZNE ZAGADNIENIA SPRZĘŻENIA
Procedura ALE składa się z następujących po sobie kolejno kroków: kroku odwzorowawczego i kroku adwekcyjnego. Krok adwekcyjny przeprowadzany jest przy założeniu, że zmiany położenia węzłów są niewielkie (bardzo małe) w porównaniu z wielkościami charakterystycznymi (długościami krawędzi) elementów otaczających te węzły.
Dodatkowym atutem wykorzystania tej procedury jest zapewnienie stałej topologii siatki MES.
Dokładność tę uzyskuje się dzięki algorytmowi użytemu w niej do odwzorowania rozwiązania z siatki zniekształconej do wygładzonej, który realizuje to z dokładnością do małych drugiego rzędu. W podejściu teoretycznym, procedura ALE zawiera w sobie jako podzbiór formuły eulerowskie. Formuły te pozwalają na określenie parametrów dla więcej niż jednego materiału w pojedynczym elemencie. Jednakże wzrost liczby materiałów pociąga za sobą wzrost liczby niezbędnych parametrów materiałowych (stałych, współczynników, itp.), co sprawia, że pracochłonność jednego kroku drastycznie rośnie. Zakres problemów, jakie możemy realizować przy użyciu procedury ALE, zależy wyłącznie od stopnia skomplikowania algorytmów odpowiedzialnych za wygładzanie siatki MES. W trakcie realizacji kroku eulerowskiego większość czasu jest poświęcana na obliczenia związane z przenoszeniem materiału pomiędzy sąsiadującymi elementami, a jedynie niewielki ułamek poświęcany jest na rozwiązanie problemu, gdzie i jak zmodyfikować siatkę MES. Obecnie stosowane algorytmy adwekcyjne są skomplikowane i czasochłonne, jednakże pozwoliły na wyeliminowanie błędów, jakie pojawiały się w pierwszych algorytmach I rzędu dokładności (fałszywe oscylacje w otrzymanych wynikach, brak stabilności, ograniczenie zakresu parametrów, itp.).
W ogólności w procedurze ALE możemy wyróżnić:
1. Przeprowadzenie klasycznego kroku lagrange’owskiego.
2. Przeprowadzenie kroku adwekcyjnego, z którym związane jest:
a) podjęcie decyzji o tym, które węzły przemieścić, b) przemieszczenie skrajnych węzłów,
c) przemieszczenie węzłów znajdujących się wewnątrz,
d) przeliczenie wszystkich zmiennych odniesionych do elementów, e) przeliczenie wartości pędu oraz uaktualnienie prędkości.
Każda wielkość zmienna musi być „przetransportowana”. Te wielkości to m.in.: prędkość, gęstość, energia wewnętrzna, sześć składowych tensora naprężenia i odkształceń plastycznych oraz wzmocnienie kinematyczne. Należy pamiętać, że prędkość musi być „przenoszona”
oddzielnie, gdyż jest ona odnoszona do węzłów, a nie jak pozostałe zmienne – do elementów.
Dodatkowo dla każdego elementu muszą być zachowane zasady zachowania masy, pędu i energii wyrażone w następujących równaniach [1]:
( ) ( ) ( )
V t S t
dM d
dV w v ndS
dt =dt
∫
ρ =! ∫
ρ % − ⋅% % (1)( ) ( ) ( ) - V(t) + V(t)
V t S t
dQ d
vdV v w v ndS pdV gdV
dt% =dt
∫
ρ% =! ∫
ρ% % − ⋅% %∫
∇%∫
ρ% (2)V ( t ) S( t ) S(t) V(t)
dE d
edV e( w v ) ndS - pv ndS+ g vdV
dt =dt
∫
ρ =! ∫
ρ % − ⋅% %∫
% %⋅∫
ρ%⋅% (3) gdzie: ρ jest gęstością ośrodka płynnego (w naszym przypadku gazu), p jest ciśnieniem gazu,ρ% prędkością rozchodzenia się gazu, g
% przyspieszeniem ziemskim i e energią właściwą.
Wielkości takie jak M (masa całkowita), Q (pęd układu) i E (energia układu) dotyczą objętości V(t) zamkniętej przez płaszczyznę S. W przestrzeni Eulera (u nas jest to gaz) prędkość jest definiowana przez w
% . W układzie współrzędnych Lagrange’a wektor prędkości jest oznaczony przez v
% . Wektor n
% jest definiowany jako wektor normalny do powierzchni S.
Odpowiednia realizacja procesu sprzężenia pomiędzy ośrodkiem płynnym (gazem) a ośrodkiem stałym wymaga doboru odpowiedniej liczby punktów całkowania na granicy tych dwóch obszarów. W innym razie może dochodzić do tzw. sztucznego wypływu gazu i przenikania przez ośrodek Lagrange’a (rys. 1). Sprzężenie między obszarami realizowane jest za pomocą metody opartej na funkcji kary [4].
Rys. 1. Realizacja procesu sprzężenia pomiędzy gazem a ciałem stałym
3. PRZYKŁADY ANALIZ
3.1. Spalanie ładunku prochowego
Na bazie klasycznych równań rozwiązania problemu głównego balistyki wewnętrznej, opracowano dla pakietu LS-Dyna koncepcję modelowania spalania się materiału wybuchowego miotającego na podstawie klasycznych równań modelu matematycznego z uwzględnieniem równania bilansu energii, równania bilansu masy i równania ruchu pocisku [6]. Obliczenia przy wykorzystaniu powyższego układu równań można prowadzić w oddzielnym pakiecie obliczeniowym, np. Matlab, przy założeniu całkowitego spalenia się prochu zanim pocisk opuści lufę. Proces obliczeniowy można wykonywać metodą kolejnych przybliżeń, przerywając iteracje, gdy względna różnica pomiędzy założoną i obliczoną wartością względnej masy spalonego prochu Ψ w chwili, gdy pocisk opuszcza przewód lufy, będzie mniejsza od założonego błędu. Na rys. 2 przedstawiono sposób obliczania zmiany ciśnienia gazów prochowych w komorze nabojowej na podstawie znajomości intensywności powstawania gazów materiału wybuchowego miotającego. Wzrost ciśnienia gazów prochowych następuje, gdy materiał wybuchowy miotający ulega spalaniu. Po spaleniu (Ψ=1) nad układem lufa – nabój wykonywana jest praca adiabatycznie rozprężających się gazów.
Charakterystykę intensywności powstawania gazów Γ (rys. 2), która opisuje własności materiału wybuchowego miotającego (prochu) w funkcji względnej masy spalonego ładunku,
można otrzymać z eksperymentu przy użyciu lufy balistycznej bądź z literatury [6]. Istotnym elementem koncepcji jest również trójwymiarowy model dyskretny, na którym została rozpięta siatka węzłów wraz z elementami skończonymi (rys. 3). Ruch deformujących się ciał został opisany współrzędnymi Lagrange’a. Dodatkowo model został uzupełniony o obszary w opisie we współrzędnych Eulera, w którym występuje ruch płynów (gazy prochowe, powietrze) (rys. 3). Zmianę ciśnienia w pakiecie LS-DYNA realizowano w obszarze, który modelował proch, za pomocą następującej relacji:
) (t p
pLSDYNA = obl (4)
gdzie: pLSDYNA - ciśnienie realizowane w pakiecie LS-Dyna, pobl- ciśnienie gazów prochowych obliczone z równań [6].
t
ps- ciśnienie końca palenia się prochu
Γ
= 0
Ψ Ψ = 1
Is- impuls ciśnienia do chwili spalenia prochu.
pmax
P
t
ps- ciśnienie końca palenia się prochu
Γ
= 0
Ψ Ψ = 1
Is- impuls ciśnienia do chwili spalenia prochu.
pmax
P
t
ps- ciśnienie końca palenia się prochu
Γ
= 0
Ψ Ψ = 1
Is- impuls ciśnienia do chwili spalenia prochu.
pmax
P
Rys. 2. Wykres zmiany ciśnienia i intensywności opisujących powstawanie gazów prochowych
Rys. 3. Wstępny model dyskretny układu lufa –nabój w opisie eulerowsko-lagrange’owskim
3.2. Symulacja procesu detonacji
Dokonano sprzężenia dwóch ośrodków, tj. ośrodka Lagrange’a, w którym była modelowana rura z ośrodkiem Eulera, który ją otaczał, a w którym realizowano symulację detonacji i rozchodzenie się fali ciśnienia jej towarzyszącej. Fala ta stanowiła obciążenie dla ścianek badanej rury. Głównym problemem w symulacji tych oddziaływań jest odpowiednie zdefiniowanie warunków początkowo-brzegowych oraz poprawne opisanie zjawiska detonacji, a następnie procesu rozchodzenia się fali ciśnienia w ośrodku oraz oddziaływania jej na elementy znajdujące się w rozważanym obszarze.
Do opisu procesu spalania materiału wybuchowego (HE – high explosive) i rozchodzenia się powstałej w jego wyniku fali detonacyjnej zostały zastosowane algorytmy wykorzystujące metodę ALE. Analizowane w układzie współrzędnych Eulera ciało (w badanym przypadku gaz) przemieszcza się w obszarze opisanym przez nieruchomą przestrzenną siatkę przestrzenną. W analizowanym przypadku we wspomnianym układzie zostały opisane dwa ciała: powietrze i materiał wybuchowy. Wielkości fizyczne, które opisują te ciała, zostały zdefiniowane w odpowiednich równaniach stanu. Do opisania procesu spalania materiału detonacyjnego wykorzystano równanie Jones-Wilkins-Lee (JWL) [2,4,5].
W realizowanych badaniach numerycznych w układzie współrzędnych lagrange’owskich zostało opisane ciało - rurociąg, na który oddziałuje fala ciśnienia generowana w wyniku rozwiązania równania JWL. Ośrodek ten charakteryzuje się tym, że siatka elementów pokrywających badany obiekt (w naszym przypadku rurociąg) przemieszcza i deformuje się wraz z tym ciałem. W wyniku symulacji procesu detonacji otrzymuje się zjawisko rozchodzenia się ciśnienia fali detonacyjnej przez oczka siatki eulerowskiej, które napotyka na ruchomą siatkę lagrange’owską (rurociąg) i stanowi dla niej obciążenie w postaci impulsu ciśnienia. Efekt wstępnego oddziaływania wygenerowanej fali ciśnienia na rurociąg przedstawiony jest na rys. 4. Do opisu procesu zniszczenia w modelu MES rurociągu zastosowano odkształceniowe kryterium zniszczenia z uwzględnieniem zmiany prędkości odkształceń [7].
Rys. 4. Model układu numerycznego sprzężenia eulerowsko-lagrange’owskiego na przykładzie rury poddanej działaniu fali detonacyjnej oraz obraz deformacji wstępnej uzyskanej w wyniku
przeprowadzonych symulacji.
4. WNIOSKI
W pracy omówiono koncepcję realizacji sprzężenia eulerowsko-lagrange’owskiego zaimplementowanego w pakiecie LS-DYNA na przykładzie dwóch badanych problemów. Oba te przypadki charakteryzują się tym, że w trakcie symulacji realizowane jest wysokoenergetyczne obciążenie. Przeprowadzone analizy numeryczne pozwoliły na realizację
szybkozmiennego obciążenia generowanego w środowisku współrzędnych Eulera (gaz), a następnie propagację tego gazu i interakcję z obiektami fizycznymi (pocisk, lufa, rura) opisanymi we współrzędnych Lagrange’a.
Przedstawiona koncepcja pozwoli efektywniej wykorzystać modele komputerowe i uzyskiwać rozwiązania złożonych problemów nieosiągalne innymi metodami (np.
analitycznymi). Proces weryfikacji otrzymanych wyników jest aktualnie prowadzony w warunkach poligonowych.
PODZIĘKOWANIE
Pracę zrealizowano dzięki wsparciu finansowemu otrzymanemu z MNiSzW w ramach projektów 4T12C00628 i 0T00B01429.
LITERATURA
1. Casadei F., Halleux JP. :An algorithm for permament fluid-structure interaction in explicit transient dynamics. “Computer Methods in Applied Mechanics and Engineering”, 1995, 128 s. 231-289.
2. Cichocki K.: Effects of underwater blast loading on structures with protective elements.
“Int. J. of Impact Engineering”, 1999, 22, s. 609-617.
3. Dębski A.: Badania metalograficzne luf kalibru 5,56 mm po próbach trwałościowych. W: VI Międzynarodowa Konferencja Uzbrojenia. Waplewo 2006.
4. Hallquist JO.: LS-Dyna. Theoretical manual. California Livermore Software Technology Corporation 1998.
5. Lu Y. Xu K., Lim HS.: Numerical simulation of concrete break-up under explosive loading.
Proceedings of the Design and Analysis of Protective Structures against Impact/Impulsive/Shock Loads, Tokyo, Japan 2003, s. 278-286.
6. Łazowski J., Małachowski J., Niezgoda T.: Preliminary concept of numerical modeling of combustion process for gunpowder using Ls-Dyna code. “Journal of KONES Powertrain and Transport”, 2007 Vol. 14, No. 1, s. 255-260.
7. Małachowski J., Niezgoda T., Wekezer J.: Analytical damage assessment of pipelines due to H.E. material detonation. Proceedings of the 2nd Int. Conference on Thermal Eng. Theory and Applications, Al Ain, United Arab Emirates, 2006 s. 26-29.
CHOSEN EXAMPLES OF COUPLING PROCESS MODELLING WITH ALE METHOD
Summary. Research communities around the world are using numerical methods to study systems and structures under impact due to blast. Dynamic response of some mechanical units subjected to the shock wave produced by the detonation or gunpowder burning is presented in this paper. Coupled Euler and Lagrange formulation are used in the finite element analysis of such problems to accurately represent the wave propagation phenomenon. In the first part of this paper the preliminary concept of numerical modelling of combustion process for gunpowder is presented. In the second part dynamic response of a pipe element subjected to the shock wave produced by the detonation of high explosive (HE) materials is presented. In both cases coupling effects between Euler and Lagrange formulations were accounted for by using the “CONSTRAINED LAGRANGE IN SOLID” feature of the LS-Dyna code.
