曹 源，金先龍，李 政

（上海交通大學機械系統與振動國家重點實驗室，上海 200240）

儲液容器廣泛應用于儲運、航空等領域，隨著材料性能的提高，在許多抗沖擊場合，如空投、爆炸中，柔性儲液容器開始被越來越多的應用。許多學者已對容器內液面晃動、容器沖擊響應以及流固耦合等問題進行了長期深入的研究［1－3］。與早期研究相比，柔性容器沖擊響應問題的難點是大變形和非線性，使得其在沖擊過程中的動態(tài)響應規(guī)律更加復雜。對于這類問題的研究，目前還是以實驗為主。隨著計算方法的改進和計算機性能的提高，數值模擬手段已有可能發(fā)揮更大的作用。

對于大變形的流固耦合問題，很難用單一的Lagrange方法或Euler方法進行研究。所以ALE方法開始應用于流固耦合分析［4－5］。ZHANG Ai－nian等［6］分別使用附加質量法、Lagrange方法和 ALE方法分析了船舶碰撞時流固耦合作用，并認為與ALE方法相比，附加質量方法和Lagrange方法低估了液體與結構之間相互作用的影響。M.Anghileri等［7］研究了儲液箱體與地面碰撞時，箱內液面大幅晃動及箱體變形的問題，并將ALE方法與其他數值模擬方法進行了比較，并且與實驗結果吻合很好。

本文中針對跌落或空投中的柔性儲液罐動態(tài)響應進行實驗研究和數值模擬。采用基于Lagrange方法描述下的Belytschko－Tsay（B.T）膜單元建立柔性儲液容器模型，流體部分則是基于ALE方法建模，并利用罰函數法實現流固耦合。在實驗基礎上，使用數值模擬的方法深入分析了空投高度、裝水量等因素對容器動態(tài)響應的影響規(guī)律。

1 流固耦合方法與原理

利用ALE方法建模與罰函數法相結合對沖擊過程中的流固耦合現象進行數值模擬。罰函數法的優(yōu)點在于保證了耦合接觸過程中的能量守恒。首先，對于每個結構節(jié)點，搜索包含該節(jié)點的ALE單元。然后，通過計算穿透速度與時間來計算流體節(jié)點的穿透深度。最后，通過正比穿透深度計算基于罰函數的耦合力。耦合力可以看作是節(jié)點總力中的一個外部力。在此基礎上，通過計算每一個時間步長的總節(jié)點力，可以得到流體和結構體耦合界面上每一點的速度、位移和應力等物理量，從而實現流固耦合界面上的相互作用，即耦合計算。ALE方法則綜合了Lagrange方法與Euler方法的優(yōu)點，在材料域與空間域外引入了參考域，并在參考域網格上求解，既解決了Lagrange方法中材料嚴重扭曲的問題，又解決了Euler方法中移動邊界引起的復雜性問題。

基于ALE方法的連續(xù)性方程和動量方程為

式中：χ為ALE坐標；vi為水的流動速度；xi和xj為空間坐標；ci和cj均為對流速度；ρ為流體密度；bi為流體體力，σij為應力張量。

對儲液容器變形的描述使用彈性體的連續(xù)方程

式中：X是Lagrange坐標，ρs為容器密度，fi為體力，u為固體結構位移。

對于水這種弱可壓縮流體，引入線性牛頓流體本構方程

式中：p為水的靜態(tài)壓力，μ＝μ（sij）為動力粘性系數，sij為應變率張量，δij為克羅內克常數。

流固耦合需在液體和結構界面上滿足幾何相容條件和力的平衡條件，即

式中：vi和vs分別為流固耦合界面上水的流動速度和固體材料速度，Ff和Fs分別是流體和固體結構作用在流固耦合界面上的力，這些力可由罰函數方法計算獲得

2 實驗簡介

圖1 儲液罐結構Fig.1 Stucture of container

為驗證本文中方法和模型的正確性，首先對裝水量η＝80%的柔性儲液罐進行100m高度的空投沖擊實驗。所使用的儲液罐為?0.5m×0.8m的圓柱體，如圖1所示。儲液罐壁厚為10mm。罐體采用簾布材料，該簾布采用尼龍織物，基體采用天然橡膠。其頂端由鋼蓋密封。在100m高度，以零初始速度垂直拋落，整個拋落過程由高速照相機記錄。

實驗中，罐體肩部（位置1）破裂，儲液罐內的水濺出。所以，罐體肩部被確定為脆弱位置之一。此外，罐體中部（位置2）和底部（位置3）也是易受沖擊和易損部位，這些部位也是數值模擬時關注的重點。

3 數值模擬與分析

3.1 數值模型

以實驗為基礎，建立三維有限元數值模型。該模型包括：地面、罐狀柔性儲液容器、水和容器內空氣。相對于柔性容器，沖擊過程中混凝土地面變形較小，模型中將地面材料作為彈塑性體處理。容器的鋼蓋處理為彈性體。各項參數如表1所示，其中，ρ為密度，E為彈性模量，ν為泊松比，σy為屈服應力，G為剪切模量，β為硬化系數。

容器主體復合材料結構的力學性能是數值模擬的關鍵，本文中采用 XIA Zhi－h(huán)ui等［8］提出的周期性邊界加載方法計算獲得簾布材料的力學性能參數，如表2所示。其中E1、E2和E3分別為材料3個方向的彈性模量；G12、G23和G13分別為材料3個方向的剪切模量；ν12、ν23和ν13分別為材料3個方向的泊松比。由表2可知材料呈橫觀各向異性。由于容器壁較薄，所以采用基于Lagrange方法描述的B.T膜單元，沿膜厚度方向取3個積分點。

表1 材料參數Table1 Material porperties and parameters

表2 容器材料力學參數Table2 Mechanical properties of the composite

為了節(jié)省計算時間，忽略罐體墜落過程中的空氣阻力，將拋落高度換算為罐體沖擊地面時的速度，作為初始條件。對于內部流體采用基于ALE方法描述的六面體單元。為了確保結果的準確性，所有單元尺寸和形狀均仔細控制。整個模型單元數為228 033，節(jié)點數為238 484。

數值模擬采用基于顯示計算方法的ANSYS／LS－DYNA有限元軟件，為確保計算的穩(wěn)定性，時間步長由模型的最小單元尺寸控制，為約2.0μs。在流固耦合的模擬中，過多的耦合點將導致數值計算不穩(wěn)定，而耦合點較少則易產生泄漏現象。由于柔性容器為曲面結構，存在較小的夾角，故耦合點數取為4。此外，耦合最小體積參數決定了流體物質達到某個單元體積時，流體與結構發(fā)生耦合作用。由于空投物體在觸地時屬于高速碰撞沖擊問題，極易產生泄漏現象，所以計算中將耦合最小體積參數值設為0.1，可以有效地減小流體的泄漏現象。

3.2 數值模擬與實驗結果比較

對相同條件下的數值模擬和實驗結果進行比較，以驗證方法和模型的有效性。定義罐體直徑變化量Δd來評估柔性罐體的變形情況。所謂直徑變化量，是指沖擊過程中不同時刻罐體中間直徑d與盛有相應水量的罐體靜置于地面時罐體中間直徑d0之差

圖2 變形量比較Fig.2 The comparison of the deformation

對罐體徑向變形量的實驗結果和柔性罐狀容器的數值模擬結果相比較發(fā)現，數值模擬結果相對于物理實驗結果要小一些，但仍處于可接受的范圍，如圖2所示。

不同時刻下柔性容器應力分布和內部流體的變化，如圖3所示。從圖中可以看出，流體形態(tài)與柔性容器變形一致，其內部的空氣先是被壓縮，隨著柔性容器的反彈，空氣完全被水包圍。當沖擊開始時，柔性容器表面最大應力出現在最先與地面接觸的底部；反彈時，由于鋼蓋慣性較大，在頂部形成“彈坑”形狀，此時，容器肩部變形最大，應力也最大。

圖3 柔性容器表面應力分布及流體形態(tài)Fig.3 Stress distribution of the container and configuration of the fluid

柔性容器上3個位置主應力σp變化曲線，如圖4所示。從曲線中可以看出，3個位置的應力依次達到最大值，與這些位置受沖擊的先后順序一致。其中，肩部的最大應力幾乎是其他2個位置的2倍。這也與實驗中罐體反彈過程中，肩部破裂的結果十分吻合。

3.3 響應規(guī)律分析及討論

從拋落高度和裝水量等2個方面對容器的沖擊響應進行研究。對裝水量η＝80%時，100、200和400m等3個拋落高度進行數值模擬。圖5顯示了這3種高度下，3個位置上最大主應力σm的比較。

圖4 關鍵位置處主應力曲線Fig.4 Principal stresses at key positions

圖5 不同高度下最大應力Fig.5 Principal stresses at different heights

柔性容器肩部的應力一直大于另外2個位置的，而隨著高度的增加，容器中部的應力增加最快。當高度達到400m時，該位置的最大應力幾乎與肩部相同。所以，隨著高度的增加，容器中部也成為危險部位。另一方面，容器底部的應力隨高度的增加，升高緩慢。特別是當高度低于200m時，容器底部應力變化不大。

圖6給出了3個拋落高度下柔性容器直徑的變化曲線。容器直徑變化量和最大幅度都隨高度的增加而變大。其中容器直徑變化量的最大值從30.6mm增加到64mm。圖中還可以看出，盡管拋落高度不同，但振蕩收斂的時間沒有大的改變。

圖7顯示的是空拋落高度100m時，裝水量η＝30%，50%，80%和100%等4種情況下容器上3個位置的最大應力。在裝水量不同的情況下仍然是容器肩部應力最大，且與裝水量呈“V”型關系。當η＝80%時，肩部應力是4種情況下最小的。原因可能是：當裝水量少時，鋼蓋可以上下大幅運動，導致容器肩部產生大變形。隨著裝水量的增加，肩部變形減小，應力也降低。但是當容器100%充滿水時，沖擊產生的力和能量可以通過水直接作用在肩部。由于水的可壓縮性很小，剛度又大于柔性材料，減震效果小于容器本身，所以導致肩部的應力急劇增加。與肩部應力的變化不同，容器底部和容器中部的應力與裝水量的變化的關系不大。還發(fā)現容器的變形量隨著裝水量增加而呈線性增加。

圖6 不同高度下柔性容器直徑變化Fig.6 Container deformations at different heights

圖7 不同裝水量下最大應力Fig.7 Principal stresses under different amount of water

4 結 論

在對柔性儲液罐進行空投實驗的基礎上，采用基于ALE方法描述的流固耦合狀態(tài)和三維有限元模型對柔性儲液容器的沖擊響應進行研究。實驗和數值模擬結果都顯示，柔性儲液容器的肩部是其薄弱部位。當達到一定的拋落高度時（如400m），容器中部也會變成另外一個薄弱部位。容器的裝水量對其抗沖擊性有很大影響。其中肩部的應力與裝水量呈“V”型關系。而柔性容器的變形則基本與裝水量呈線性關系。在其他條件不變的情況下，裝水量為80%時，柔性儲液容器的抗沖擊性最好。

［1］Moiseev N N.Introduction to the theory of oscillations of liquid containing bodies［J］.Advances in Applied Mechanics，1964，8（5）：233－289.

［2］Ma D C，Su T C.Fluid－structure Vibration and Liquid Sloshing［C］∥Proceedings of Pressure Vessels and Piping Conference.New York：American Society of Mechanical Engineers，1987.

［3］Shimizu N.Advances and trends in seismic design of cylindrical liquid storage tanks［J］.Japan Society of Mechanical Engineers International Journal，1990，33（2）：111－124.

［4］Nitikitpaiboon C，Bathe K J.An arbitrary Lagrangian－Eulerian velocity potential formulation for fluid－structure interaction［J］.Computers ＆Structures，1993，47（4／5）：871－891.

［5］Nomura T.ALE finite element computations of fluid－structure interaction problems［J］.Computer Methods in Applied Mechanics and Enginneering，1994，112（1／2／3／4）：291－308.

［6］ZHANG Ai－nian，Suzuki K.A comparative study of numerical simulations for fluid－structure interaction of liquidfilled tank during ship collision［J］.Ocean Engineering，2007，34（5／6）：645－652.

［7］Anghileri M，Luigi M L，Castelletti，Tirelli M.Fluid－structure interaction of water filled tanks during the impact with the ground［J］.International Journal of Impact Engineering，2005，31（3）：235－254.

［8］XIA Zhi－h(huán)ui，ZHANG Yun－fang，Ellyin F.A unified periodical boundary conditions for representative volume elements of composites and applications［J］.International Journal of Solids and Structure，2003，40（8）：1907－1921.