容器匀加速运动状态下液体大幅晃动的解析计算

张海涛　 孙蓓蓓

(1东南大学机械工程学院, 南京 211189)(2南昌大学机电工程学院, 南昌 330031)

液体晃动是一种常见的流体运动现象.在油气储运、航空航天等工业领域,液体晃动可能会对结构安全性和系统运动稳定性产生较大的影响,因此液体晃动特性分析及其抑制措施已成为工程实际中的重要问题.小幅线性晃动的研究起步较早,现已形成非常成熟的理论[1].而液体大幅晃动通常具有较强的非线性,运动形式趋于多样化,研究难度较大.目前关于液体晃动的研究主要集中于大幅非线性晃动.

罐式车辆在启动或制动过程中,惯性作用可能会引起液体产生大幅晃动,并对罐壁造成较大冲击.Shahravi等[2]采用ALE有限元法和基于等效力学模型的多体动力学方法计算容器内的液体晃动,并对结果进行了比较，同时模拟了罐车在转弯和加减速等工况下的晃动现象,分析了一些重要晃动参数的影响作用.Yu等[3]研究了罐车在急刹和急转弯等危险工况下的行驶稳定性,讨论了晃动力与罐体形状、充液率和晃动固有频率等因素之间的关联性.尚春雨等[4]用Fluent软件计算了容器做匀加速运动时的液面波动曲线.Nicolsen等[5]建立了一种可与多体系统相融合的新型拉格朗日连续介质晃动模型,其中液体位移和形状改变的核心计算方法是有限元绝对坐标法.他们使用该模型分析了罐车在转弯和减速等工况下的动力学行为和车辆稳定性.

容器匀变速运动引发的晃动现象与简谐激励或随机激励的受迫晃动有很大的不同.由于晃动的复杂性,学者们大多采用数值模拟方法计算该情况下的受迫晃动.数值方法计算结果较为准确,但计算耗时较长,且不易揭示出晃动非线性的内在机理.解析方法是一种传统的非线性理论分析方法,通过寻求非线性微分方程组的解析近似解来探讨系统的动力学特性.对于矩形容器内液体自由晃动和简谐激励下的受迫晃动,学者们已获得了解析近似解[6].本文采用解析方法计算矩形容器在匀变速运动状态下的液体大幅晃动,对自由液面波动和流体压力的非线性项进行理论分析,并探讨外部相关参数(加速度大小、充液率)变化对晃动非线性的影响.

1　数学模型

设定充液容器为矩形箱体,长为b,内部充有深度为h的液体(密度设为ρ=1 000 kg/m3).设定液体为无黏不可压缩的理想流体,忽略表面张力的影响,采用势流模型描述.假设在晃动过程中,自由液面的运动是连续的,即无翻转、破碎等现象.设定容器在绝对坐标系下向右做匀加速直线运动,其加速度为β.建立充液容器的相对坐标系xoz,坐标原点设定为左壁面和静止水平液面的交点,向右和向上为正方向.整体晃动模型如图1所示.在相对坐标系下考虑液体晃动,得到如下方程组[7]:

(1)

(2)

(3)

(4)

(5)

式中,φ为流体速度势;ξ为自由液面波高函数,反映了各自由液面质点偏离静平面的程度;g为重力加速度.式(1)为理想流体的连续性方程,式(2)、(3)为固壁边界条件,式(4)、(5)为自由液面边界条件.假设初始时刻液体和容器均静止,自由液面为未扰动水平面,由此可得如下初始条件:

φ(x,z,0)=0,ξ(x,0)=0

(6)

对上述微分方程进行求解后,内部流体压力p可根据如下表达式[8]进行计算:

(7)

图1　充液容器大幅晃动示意图

2　解析方法

首先对非线性方程进行无量纲化.对于液体晃动方程,一般选位移和频率作为量纲参数,进行后续处理.然而匀加速激励不同于简谐激励,没有实际的位移量，所以定义如下虚拟振幅:

(8)

式中,a为虚拟振幅,为无量纲化而设定的,没有实际意义;Ω1为液体线性晃动第一阶固有频率.由此可选取如下无量纲量:

(9)

并设定如下无量纲小参数:

(10)

(11)

采用基于变量分离的小参数法对无量纲晃动方程求解,设定近似解形式为

(12)

将式(12)代入无量纲晃动方程,令小参数ε同次幂系数相等,可将非线性方程转化为一系列线性方程组.

首先考虑ε=0的情况,该零次方程组正是线性晃动方程,可采用分离变量法求解.根据线性晃动理论,其解的形式可设定为[7]

(13)

式(13)可满足连续性方程和固壁边界条件.同时,做如下傅里叶级数展开[7]:

(14)

根据线性晃动理论[8]可知,晃动第一阶固有频率对应的第一阶反对称模态对受迫晃动有显著影响,后续奇数阶反对称模态的影响效果越来越弱,而偶数阶对称模态对受迫晃动没有影响.在式(14)中,除直流分量A0外,其他偶数项系数均为零,进而可算出偶数阶模态方程的解也为零,这与线性晃动理论完全相符.此处求解直流分量方程(解得速度势函数和波高函数分别为-bt/2和0),以及第一模态和第三模态方程,根据式(13),将零次方程组的解近似表示为以上离散解的叠加，即

(15)

式中

(16)

其次,考虑ε=1对应的一次方程组,自由液面边界条件分别为

(17)

(18)

(19)

式(18)表达式也与之类似.为了分离变量,可设定如下形式的解:

(20)

(21)

根据上述求解过程以及式(12),可将液面运动和流体速度势解析近似解表示为

(22)

将解出的流体速度势代入无量纲压力表达式(11),即可求得晃动过程中的液体内部压力.

3　计算结果与讨论

作为一类重要的非线性分析方法,解析近似法能够快速获得计算结果,不仅能确定非线性系统的运动随时间变化的规律,而且能够得到运动特性与系统参数之间的依赖关系[9].根据已获得的解析近似解,可进一步讨论晃动非线性因素,并研究外部参数(加速度和充液率)变化对晃动非线性效应的影响.

3.1　基于解析解的一般性分析

为验证该解析近似解的准确性,在如下算例中将数值解与解析近似解做对比,其中数值解采用文献[10]的方法计算.设定b=2 m,h=1 m,β=3 m/s2,代入晃动模型.图2(a)显示了2种方法获得的左壁处自由液面波高函数的时间历程,图2(b)显示了2种方法获得的左壁底部的液体压力变化时间历程.通过比较可看出,数值解和解析解符合程度较好,说明该解析解是比较准确的.随着时间的增加,两者之间存在着一定的相位偏差.Frandsen[7]对液体自由晃动的解析解和数值解做比较时,也发现了类似的现象.他认为,如果采用更高次数的解析近似解,相位偏差将会减小.

(a) 左壁液面波高函数时间历程比较

(b) 左壁底部流体压力时间历程比较

匀加速运动下的受迫晃动与简谐激励作用下的受迫晃动有很大的区别.由图2(a)可知,容器做匀加速运动时,左壁处自由液面质点在振动过程中总是位于静止水平面上;而在简谐激励作用下,液面质点的位移可能会高于或低于静止水平面.此外,图2(a)显示振动周期约为1.7 s,对应圆频率为3.7 rad/s,与线性晃动第一阶固有频率(Ω1=3.76 rad/s)非常接近.可见对于容器匀加速运动引发的受迫晃动而言,其液面质点的振动形式比较单一;而在简谐激励受迫晃动下,液面质点的振动与外激励频率极为相关,可能会引发拍振和共振等现象[10].

表1　波高函数各非线性项的最大振幅

从表中可以看出,仅由第三模态解生成的第2项非常微小,几乎可以忽略不计;第3项和第4项由第一模态解和第三模态解相互耦合作用生成,它们的波动幅度较为接近.因此,左壁自由液面波高函数的非线性部分主要由第1项、第3项和第4项构成,其中第1项占主导地位.

人类身体的运动属性是与生俱来的能力，体育行为是后天生存活动中的经验与总结。体育既是来源于人类潜在的本能需求，又是人类对身体健康的内在需求，体育搭建起主体寻求强蛮体魄与自由存在的桥梁，在人类本性抒发中完成自身肉体的解放。体育发挥主、客体关联的中介作用来彰显主体的姿态。“在人类的自然需求中，体育运动是一种寻求生存平衡的身体本能”[9]。生命的平衡要求人类内在向度达到精神满足与肉体需求的统一，外在尺度保持实践行为与自然开发的相互和谐。体育运动在人类的顽强精神中彰显躯体的健硕体态，在人类的本能行为中呈现出自然的动物属性，在顺应人类肉身本质的同时获取精神的自由。

图3　左壁底部流体压力非线性项时间历程

(23a)

(23b)

3.2　加速度对晃动非线性的影响

(24)

可见,波高函数ξ的非线性项与加速度平方成正比.设定参数b=2 m,h=1 m,加速度分别为β=3.0, 4.5 m/s2进行计算.计算结果显示,随着加速度的增大,左壁面波高函数非线性项的函数形状未发生变化,但其数值扩大了2.25倍,这与理论分析的结果一致.

对流体压力做类似分析可知,非线性项同样与加速度平方成正比.综上所述,加速度的改变将对晃动非线性效应产生显著影响.

3.3　充液率对晃动非线性的影响

(25)

分析可知,当h/b>0.7时,其函数曲线的斜率接近于零,即充液率的增加仅使得固有频率Ω1发生微小的变化.因固有频率几乎不变,故解析近似解也基本上保持不变.当0.5

设定充液率h/b分别为0.4,0.5,0.75和1.0,计算相应左壁液面波高函数的非线性项,计算结果如图4所示.图4显示,h/b=0.75和h/b=1.0的函数图像非常接近;尽管充液率仅减小了0.1,h/b=0.4与h/b=0.5的函数图像有明显区别,流体的减少加强了晃动非线性效应.事实上,充液率过低，更偏向于浅水晃动,其晃动动力学特性愈加复杂,解析解不一定适用.

(a) 较高充液率

(b) 较低充液率

4　结语

采用解析方法计算了匀加速运动状态下二维矩形容器内的液体大幅晃动问题,研究了自由液面波动和底部流体压力的非线性特征,并讨论了加速度和充液率变化对晃动非线性效应的影响.通过与数值计算结果的对比,可以看出解析近似解的精度较高.后续分析结果显示,左壁自由液面波高函数非线性部分的第1项(即仅由第一模态解生成的非线性项)为其主要影响成分;对于左壁底部流体压力非线性部分,其时间项为主要影响成分,“双峰”现象由流体压力线性部分和非线性部分共同作用形成.自由液面波高函数与流体压力的非线性项与加速度的平方成正比.当充液率较大时,充液率的变化对晃动非线性影响较小,而当充液率较小时(h/b<0.5)时,充液率的减小使得晃动非线性效应增强.

参考文献(References)

[1] Abramson H N. The dynamic behavior of liquids in moving containers[R]. Washington DC: NASA SP-106, 1966.

[2] Shahravi M, Sajjadi M R, Feizi M M. Investigation of sloshing coefficient by Arbitrary Lagrange-Euler methods in partially filled tankers[J].AdvancesinRailwayEngineering,anInternationalJournal, 2014,2(1): 1-11.

[3] Yu D, Li X, Liu H, et al. Theory and experiments on driving stability of tank trucks under dangerous working conditions[J].JournalofVibroengineering, 2015,17(5): 2521-2534.

[4] 尚春雨,赵金城.用FLUENT分析刚性容器内液面晃动问题[J].上海交通大学学报,2008,42(6):953-956.

Shang Chunyu, Zhao Jincheng. Studies on liquid sloshing in rigid containers using FLUENT code[J].JournalofShanghaiJiaotongUniversity, 2008,42(6): 953-956. (in Chinese)

[5] Nicolsen B, Wang L, Shabana A. Nonlinear finite element analysis of liquid sloshing in complex vehicle motion scenarios[J].JournalofSoundandVibration, 2017,405: 208-233. DOI:10.1016/j.jsv.2017.05.021.

[6] Ibrahim R A.Liquidsloshingdynamics:Theoryandapplications[M]. Cambridge,UK: Cambridge University Press, 2005:253-263.

[7] Frandsen J B. Sloshing motions in excited tanks[J].JournalofComputationalPhysics, 2004,196(1): 53-87. DOI:10.1016/j.jcp.2003.10.031.

[8] Dodge F T.Thenewdynamicbehaviorofliquidsinmovingcontainers[M]. San Antonio, Texas,USA: Southwest Research Institute, 2000:1-10.

[9] 刘延柱,陈立群.非线性振动[M].北京:高等教育出版社,2001:2-4.

[10] 张海涛,孙蓓蓓,陈建栋.基于自由液面预测的非线性液体晃动问题的数值模拟[J].东南大学学报(自然科学版),2014,44(2):277-282.DOI: 10.3969/j.ssn.1001-0505.2014.02.010.

Zhang Haitao, Sun Beibei, Chen Jiandong. Numerical simulation of nonlinear liquid sloshing problems based on forecast of free surface[J].JournalofSoutheastUniversity(NaturalScienceEdition), 2014,44(2): 277-282. DOI: 10.3969/j.ssn.1001-0505.2014.02.010.(in Chinese)

[11] 周宏.液体晃动数值模拟及刚-液耦合动力学研究[D].北京:清华大学航天航空学院,2008.

[12] Wu G X, Ma Q W, Eatock Taylor R. Numerical simulation of sloshing waves in a 3D tank based on a finite element method[J].AppliedOceanResearch, 1998,20(6): 337-355. DOI:10.1016/s0141-1187(98)00030-3.