基于淺水波理論液體晃動(dòng)初值問(wèn)題的數(shù)值模擬

薛偉順 熊淵博

（四川大學(xué)建筑與環(huán)境學(xué)院，成都 610065）

基于淺水波理論液體晃動(dòng)初值問(wèn)題的數(shù)值模擬

薛偉順 熊淵博?

（四川大學(xué)建筑與環(huán)境學(xué)院，成都 610065）

針對(duì)日益受到關(guān)注的液體晃動(dòng)問(wèn)題，提出了一種基于淺水波理論的研究方案.該方案采用淺水波理論而非勢(shì)流理論導(dǎo)出系統(tǒng)控制方程，并通過(guò)哈密頓體系表達(dá)；利用中心有限差分法和Stormer-Verlet算法進(jìn)行空間和時(shí)間離散；模擬了不同初值條件下的液體晃動(dòng)情況并對(duì)比分析了影響系統(tǒng)非線性響應(yīng)的主要因素.結(jié)果表明，基于淺水波理論能有效解決液體晃動(dòng)問(wèn)題；與Euler格式對(duì)比，Stormer-Verlet算法精度較高；除共振外對(duì)于系統(tǒng)非線性響應(yīng)的影響容器初始位移比初始速度更顯著；非共振情況一定條件下，充液容器運(yùn)動(dòng)過(guò)程中液體晃動(dòng)能起到阻尼作用.

液體晃動(dòng)，淺水波理論，初值問(wèn)題，數(shù)值模擬，非線性

引言

液體晃動(dòng)問(wèn)題普遍存在于工業(yè)流水線、航天航空、能源動(dòng)力、原油儲(chǔ)運(yùn)等眾多工程領(lǐng)域，且日益得到工程界和學(xué)術(shù)界的關(guān)注.尤其近些年來(lái)，隨著航空、船舶等領(lǐng)域的快速發(fā)展，艙內(nèi)液體晃動(dòng)問(wèn)題引起了人們極大地關(guān)注［1-2］.當(dāng)容器部分充液時(shí)，液體在外部激勵(lì)下會(huì)晃動(dòng)起來(lái)，可能對(duì)結(jié)構(gòu)安全造成嚴(yán)重的危害.因此，對(duì)不同激勵(lì)形式下充液系統(tǒng)響應(yīng)進(jìn)行研究有著非常重要的工程應(yīng)用價(jià)值.

在液體晃動(dòng)的描述和液體晃動(dòng)的物理特性等方面，國(guó)內(nèi)外許多專家學(xué)者開(kāi)展了大量的研究.盧軍等［3］基于勢(shì)流理論研究了平放圓柱形貯箱內(nèi)任意充液比液體晃動(dòng)問(wèn)題；J.H.Jung等［4］采用VOF方法對(duì)液體晃動(dòng)問(wèn)題進(jìn)行了數(shù)值計(jì)算；李裕龍等［5］模擬計(jì)算了艙內(nèi)液體的非線性晃蕩并建立了波浪中載液船舶耦合運(yùn)動(dòng)方程；R.Sygulski［6］假設(shè)流體為無(wú)黏、不可壓縮、微幅晃動(dòng)的理想流體，采用邊界元方法研究了三維液體的晃動(dòng)問(wèn)題；劉富等［7］采用SPH方法對(duì)棱形液艙不同充液比所對(duì)應(yīng)的艙內(nèi)液體晃動(dòng)進(jìn)行了三維數(shù)值模擬；歐陽(yáng)林輝等［8］利用有限體積法研究液體晃動(dòng)行為，對(duì)二維和三維矩形儲(chǔ)液容器液體晃動(dòng)進(jìn)行了分析；孫麗娜等［9］對(duì)罐式集裝箱液體晃動(dòng)過(guò)程進(jìn)行了數(shù)值模擬研究；李文盛等［10］基于非線性波動(dòng)理論模型，研究了儲(chǔ)液容器內(nèi)液體晃動(dòng)的動(dòng)力學(xué)響應(yīng)問(wèn)題.然而在以往的液體晃動(dòng)問(wèn)題的研究中，更多的是以勢(shì)流理論為基礎(chǔ)，導(dǎo)出系統(tǒng)的控制方程，考慮外部激勵(lì)頻率與容器一階固有頻率相近或相等時(shí)的液體晃動(dòng)問(wèn)題［12］，而對(duì)諸如容器初始位移、容器初始速度等更一般的初值條件對(duì)液體非線性響應(yīng)的影響分析還不多.

為了提出有效的數(shù)值方案，并比較分析影響液體晃動(dòng)非線性響應(yīng)的不同初始條件，從而為以后進(jìn)一步研究液體晃動(dòng)控制問(wèn)題提供可靠參考.本文嘗試以淺水波理論為基礎(chǔ)，導(dǎo)出系統(tǒng)的控制方程，進(jìn)而得到完整的數(shù)值方案.在此基礎(chǔ)上，對(duì)系統(tǒng)在不同初始條件下的非線性響應(yīng)做了數(shù)值仿真對(duì)比試驗(yàn).通過(guò)將該方案與成熟的數(shù)值方案進(jìn)行對(duì)比，證明了該方案能夠精確穩(wěn)定的對(duì)淺水晃動(dòng)與容器平動(dòng)的動(dòng)態(tài)耦合問(wèn)題進(jìn)行計(jì)算，并且具有較高的精度.

1　系統(tǒng)控制方程

容器水平運(yùn)動(dòng)，令q（t）表示空間坐標(biāo)（X，Y）下容器的水平位移，k為彈性常數(shù)，h0為靜水位，L為容器水平長(zhǎng)度，（x，y）為隨體坐標(biāo)系，則運(yùn)動(dòng)系統(tǒng)如圖1所示.

圖1　部分充液容器的運(yùn)動(dòng)原理圖Fig.1 Schematic of a moving vehicle partially filled with fluid

容器水平運(yùn)動(dòng)控制方程［11］

式中，ρ為液體密度，h（x，t）為液體高度，U（x，t）為液體水平表面速度，mf為液體質(zhì)量，mv為容器質(zhì)量，˙q為容器水平速度.

淺水方程［13］如下

式中，g為重力加速度，¨q為容器水平加速度.

聯(lián)合（1）和（2）得到剛性容器水平運(yùn)動(dòng)導(dǎo)致液體晃動(dòng)的系統(tǒng)動(dòng)態(tài)耦合控制方程.則哈密頓體系下系統(tǒng)控制方程如下

2　液體晃動(dòng)問(wèn)題的數(shù)值分析

2.1空間離散

離散參考空間如下

令xi（t）=x（ai，t），wi（t）=w（ai，t），控制方程（3）中第一個(gè)方程等式右側(cè)采用變分離散［14］；利用梯形法則離散（3）第三個(gè)和第四個(gè)方程中的積分項(xiàng).

由此得到半離散化的哈密頓形式的控制方程組如下

式中，η為梯形法則，具有廣義動(dòng)量p的量綱.

2.2時(shí)間離散

時(shí)間離散采用如下方案

這個(gè)方案即Stormer-Verlet算法，它是顯式的，并且在長(zhǎng)時(shí)間跨度上能保持能量守恒［15］.

對(duì)（6）在每個(gè)時(shí)間步n?n+1中應(yīng)用Stormer-Verlet算法，得到離散控制方程如下

式中

耦合系統(tǒng)的能量方程

當(dāng)轉(zhuǎn)化為以正則變量（q，x，p，w）表示的形式時(shí)即為哈密頓方程（4）.

2.3系統(tǒng)的數(shù)值仿真試驗(yàn)

取液體密度ρ=1000kg/m3，重力加速度g= 9.81 m/s2，容器水平長(zhǎng)度L=0.525m，為滿足淺水波方程水深遠(yuǎn)小于波長(zhǎng)的適用條件，水深在0.05m～0.08m之間取值，空間步長(zhǎng)Δa=1.05×10-4m，時(shí)間步長(zhǎng)Δt=10-4s.

其中sj是滿足方程

的有限正解，其中R是與質(zhì)量有關(guān)的無(wú)量綱參數(shù)，G是與能量有關(guān)的無(wú)量綱參數(shù).圖3（a）和圖4（a）分別給出系統(tǒng)能量誤差和液面波動(dòng)情況.由圖3（a）可以看到采用該方案系統(tǒng)能量誤差波動(dòng)穩(wěn)定并保持有界，量級(jí)保持在10-6，足夠小.由圖4（a）液面波動(dòng)情況可以看到仿真結(jié)果出現(xiàn)多峰波，表現(xiàn)出明顯的非線性符合共振對(duì)液體波動(dòng)的實(shí)際影響，說(shuō)明該方案是合理的.

圖2　能量Ev隨時(shí)間t的波動(dòng)圖Fig.2 Energy-time curves

圖3　能量誤差Fig.3 Energy error

上述兩種初始條件下，能量波動(dòng)和液面波動(dòng)都表現(xiàn)出了一定的非線性，但q（0）=0.03m時(shí)非線性表現(xiàn)更為明顯，說(shuō)明在相同條件下，即使容器的初始位移很小，系統(tǒng)也表現(xiàn)出明顯的非線性.q（0）=0.03m時(shí)雖然液面出現(xiàn)明顯的非線性波動(dòng)，能量波動(dòng)劇烈，但依然很快達(dá)到收斂，并且保持在有界范圍內(nèi)，計(jì)算穩(wěn)定性較好.兩種初始條件下能量標(biāo)準(zhǔn)誤差均保持了足夠小量級(jí)，保證了較高的精確度，說(shuō)明數(shù)值結(jié)果是穩(wěn)定可靠的.由能量誤差波動(dòng)曲線同樣可以看出，計(jì)算很快達(dá)到穩(wěn)定并能保持足夠小量級(jí)，˙q（0）=0.12m/s時(shí)誤差雖然出現(xiàn)了“漂移”，但計(jì)算表明漂移量很小并沒(méi)有發(fā)生數(shù)量級(jí)的跳躍，計(jì)算結(jié)果是可靠的.表1能量標(biāo)準(zhǔn)誤差再次印證了以上得出的結(jié)論，對(duì)比表中數(shù)據(jù)可以看出，相同條件下本文所采用的算法在精度上有顯著提高.

圖4　隨機(jī)捕捉的自由液面Fig.4 Snapshots of the free surface profile

最后考慮給液面一個(gè)初始位移，液體靜止，忽略彈性力，來(lái)探究液體晃動(dòng)是如何影響容器運(yùn)動(dòng)的.取液面初始傾角θ0=-0.5℃，規(guī)定逆時(shí)針旋轉(zhuǎn)為正，自由液面的初始位置為

容器水平位移及其細(xì)節(jié)在圖5中給出，能量波動(dòng)及其誤差在圖6和圖7中給出，能量標(biāo)準(zhǔn)誤差在表1中列出.

表1　能量標(biāo)準(zhǔn)差表Table 1 Standard deviation of the energy

圖5　容器的位移q（t）隨時(shí)間t的變化情況及0s～5s容器的位移q（t）隨時(shí)間t的變化情況Fig.5 Computed vehicle position-time curves of the container together with the corresponding curves in sub-interval of 0s～5s

圖6　能量Ev隨時(shí)間t的波動(dòng)圖Fig.6 Energy-time curve

圖7　能量誤差Fig.7 Energy error

當(dāng)k=0kg/s2時(shí)，積分q（t）表達(dá)式并代入初始條件得到

等式右端第一項(xiàng)是均勻位移，其中˙q（0）= 0.01m/s，這從圖4中0s～100s的宏觀圖像中容器位移和時(shí)間基本呈線性關(guān)系得到表現(xiàn)，第二項(xiàng)令容器位移產(chǎn)生波動(dòng)，這是由于液體質(zhì)心的晃動(dòng)所導(dǎo)致的，取一段短時(shí)間的位移進(jìn)行觀察，在圖4中0s～5s容器位移的微觀圖像中可以看出曲線有一定的波動(dòng)，說(shuō)明液體晃動(dòng)在系統(tǒng)運(yùn)動(dòng)過(guò)程中起到阻尼作用.另外，結(jié)合圖5圖6與表1可以看出此時(shí)能量出現(xiàn)不穩(wěn)定波動(dòng)，但依然保持有界，誤差很快達(dá)到穩(wěn)定并在足夠小范圍內(nèi)波動(dòng)，并且能量標(biāo)準(zhǔn)誤差足夠小，精度保持在10-8量級(jí)，表明算例的計(jì)算結(jié)果依然是可靠的.

3　結(jié)論

針對(duì)液體晃動(dòng)問(wèn)題，本文提出了基于淺水波理論，結(jié)合Stormer-Verlet算法模擬的方案.該方案采用淺水波理論導(dǎo)出系統(tǒng)的控制方程，利用中心有限差分進(jìn)行空間離散結(jié)合Stormer-Verlet算法進(jìn)行時(shí)間離散，得到了精確穩(wěn)定的數(shù)值離散方案，通過(guò)與Euler格式對(duì)比，本文方法在計(jì)算效率相近的情況下其精度顯著提高，并能夠快速收斂和保持較好的穩(wěn)定性.模擬結(jié)果表明，基于淺水波理論能有效地解決液體晃動(dòng)問(wèn)題.

利用哈密頓體系表達(dá)系統(tǒng)控制方程，為一般的導(dǎo)出簡(jiǎn)明有效的數(shù)值方案奠定了基礎(chǔ).首先利用共振算例驗(yàn)證算法精度及方案有效性，然后對(duì)初始狀態(tài)分別取容器初始位移不為零和容器初始速度不為零的兩種情況進(jìn)行了對(duì)比討論，得出除共振頻率外，容器初始位移對(duì)液體非線性晃動(dòng)的影響更為顯著，在工程應(yīng)用中應(yīng)該引起重視.然后探究了液體晃動(dòng)對(duì)剛性容器運(yùn)動(dòng)的影響，發(fā)現(xiàn)非共振情況一定條件下由于液體質(zhì)心的晃動(dòng)，液體在系統(tǒng)運(yùn)動(dòng)過(guò)程中能起到阻尼的作用.

本文模擬方案為進(jìn)一步研究更為復(fù)雜工況條件下充液系統(tǒng)響應(yīng)以及液體晃動(dòng)控制問(wèn)題奠定了基礎(chǔ).

1 Xue M A，Lin P Z.Numerical study of ring baffle effects on reducing violent liquid sloshing.Computers&Fluids，2011，52：116～129

2 Hasheminejad S M，Aghabeigi M.Sloshing characteristics in half-full horizontal elliptical tanks with vertical baffles. Applied Mathematical Modeling，2012，36（1）：57～71

3 盧軍，楊翊仁.平放圓柱體內(nèi)液體晃動(dòng)特性及橫向受力分析.科學(xué)技術(shù)與工程，2009，9（14）：3945～3949（Lu J，Yang Y R.Characteristics of liquid sloshing and transverse force analysis of horizontal liquid-filled cylindrical tank.Science Technology and Engineering，2009，9（14）：3945～3949（in Chinese））

4 Jung J H，Yoon H S，Lee C Y，et al.Effect of the vertical baffle height on the liquid sloshing in a three-dimensional rectangular tank.Ocean Engineering，2012，44：79～89

5 李裕龍，朱仁傳，繆國(guó)平等.波浪中載液船舶運(yùn)動(dòng)的時(shí)域模擬.上海交通大學(xué)學(xué)報(bào)，2012，46（3）：379～384（Li Y L，Zhu R C，Miao G P，et al.Simulation of ship motions coupled with tank sloshing in time domain.Journal of Shanghai Jiaotong University，2012，46（3）：379～384（in Chinese））

6 Sygulski R.Boundary element analysis of liquid sloshing in baffle tanks.Engineering Analysis with Boundary Element，2011，35（8）：978～983

7 劉富，童明波，陳建平.基于SPH方法的三維液體晃動(dòng)數(shù)值模擬.南京航空航天大學(xué)學(xué)報(bào)，2010，42（1）：122～126（Liu F，Tong M B，Chen J P.Numerical simulation of three-dimensional liquid sloshing based on SPH method. Journal of Nanjing University of Aeronautics&Astronautics，2010，42（1）：122～126（in Chinese））

8 歐陽(yáng)林輝，鄧曉濤，霍世慧.基于有限體積法的液體晃動(dòng)數(shù)值模擬.航空制造技術(shù)，2013，19：75～79（Ouyang L H，Deng X T，Huo S H.Numerical simulation of sloshing waves based on volume of fluid method.Aeronautical Manufacturing Technology，2013，19：75～79（in Chinese））

9 孫麗娜，周?chē)?guó)發(fā).罐式集裝箱液體晃動(dòng)過(guò)程數(shù)值模擬研究.振動(dòng)與沖擊，2012，31（22）：147～150（Sun L N，Zhou G F.Numerical simulation for liquid sloshing process of a tank truck.Journal of Vibration and Shock，2012，31（22）：147～150（in Chinese））

10 李文盛，趙友清，賈善坡等.儲(chǔ)液容器內(nèi)液體晃蕩的非線性動(dòng)力學(xué)分析.爆炸與沖擊，2014，34（1）：87～92（Li W S，Zhao Y Q，Jia S P，et al.Numerical analysis on liquid sloshing in storage container by nonlinear dynamics method.Explosion and Shock Waves，2014，34（1）：87～92（in Chinese））

11 寧德志，宋偉華，滕斌等.容器固有頻率對(duì)液體晃蕩的影響.海洋科學(xué)進(jìn)展，2012，30（1）：45～53（Ning D Z，Song W H，Teng B，et al.Effect of natural frequencies of container on fluid sloshing.Advances in Marine Science，2012，30（1）：45～53（in Chinese））

12 Cooker M J.Water waves in a suspended container.Wave Motion，1994，20（4）：358～395

13 Alemi Ardakani H，Bridges T J.Dynamic coupling between shallow-water sloshing and horizontal vehicle motion.European Journal of Applied Mathematics，2010，21（6）：479～517

14 Marsden J E，West M.Discrete mechanics and variational integrators.Acta Numerica，2001，10：357～514

15 Hairer E，Lubich C，Wanner G.Geometric numerical integration illustrated by the Stormer-Verlet method.Acta Numerica，2003，12：399～450

NUMERICAL SIMULATION ON THE INITIAL VALUE PROBLEMS OF LIQUID SLOSHING BY SHALLOW WATER WAVE THEORY

Xue Weishun Xiong Yuanbo?
（College of Architecture and Environment，Sichuan University，Chengdu 610065，China）

Based on the shallow water wave theory，a research scheme was proposed for liquid sloshing problems in this paper.In this scheme，the system control equations was derived using the shallow water wave theory rather than potential flow theory，and Hamiltonian system was used to express the equations.The centered finite difference and Stormer-Verlet algorithm were then used for discrete space and time，respectively.Moreover，liquid sloshing under different initial conditions were simulated and compared to obtain the main factor that affecting the nonlinear response of the system.The simulating results indicate that shallow water wave theory can effectively solve the problem of liquid sloshing.Compared with the Eulerian scheme，the Stormer-Verlet scheme had higher accuracy，while compared with the initial velocity of the vehicle，the initial displacement had more remarkable effect on the nonlinearity response of the system.It is also found that under certain conditions，liquid sloshing acted as damping in the process of the vehicle movement.

liquid sloshing，shallow water wave theory，initial value problems，numerical simulation，nonlinearity

23 January 2015，revised 8 April 2015.

E-mail：xyuanbo@163.com

10.6052/1672-6553-2015-020

2015-01-23收到第1稿，2015-04-08收到修改稿.

E-mail：xyuanbo@163.com.