充液贮箱液固耦合振动模态与阻尼比研究

沈　骥， 张美艳， 唐国安

(复旦大学 力学与工程科学系，上海　200433)



充液贮箱液固耦合振动模态与阻尼比研究

沈骥， 张美艳， 唐国安

(复旦大学 力学与工程科学系，上海200433)

摘要：对于液体运载火箭纵向振动的稳定性分析，在建立贮箱内液体的振动模型时，不仅需要考虑结构的弹性特性，还需要考虑液体晃动时的阻尼特性。由于模态阻尼比与一个振动周期内阻尼力做功、平均机械能及模态固有频率相关。因此可以先通过结构与液体的耦合振动分析确定其模态，进而分别计算液体与结构接触面、自由面以及液体内部阻尼力做功，同时计算平均机械能，从而得到模态阻尼比。在计算模态阻尼比时，对耦合振动的模态分析、液体内部流场计算均采用等效比拟的方法，为充分利用通用程序计算创造了条件。此方法不仅能够提高效率，同时也便于工程上的应用。

关键词：液固耦合振动；充液贮箱；模态；阻尼比

充液贮箱在许多工程领域都有所利用，例如航天、石化、船舶、水利水电等领域。特别在航天领域，随着液体火箭任务目的、航程、性能和材料的发展，其贮液箱内液体所占比重越来越大，液体和结构耦合作用对箭体的影响越发明显，因此对耦合作用的研究也愈发重要。在进行液固耦合问题的动力学分析时，以往学者关注于模态频率和振型的分析。但是在分析稳定性与响应时需要对耦合阻尼进行更深入的研究[1]。

在流体力学中，因为N-S方程的非线性性具有一定的复杂度。通常通过简化流场区域、给定边界条件和一些近似假设而获得近似结果。在模态分析中，Fox等[2]讨论了不同形状的刚性二维贮箱问题，并对求解方法进行了归类。朱琳等[3]建立了基于有限元方法的液固耦合模态分析方法，并对三维圆柱形贮箱建模。朱昶帆等[4]建立了用通用软件实现液固耦合模态分析的方法。而阻尼分析中，Case等[5]计算了不可压缩流体的表面波阻尼。Henderson等[6]采用边界层假设获得三维圆柱形刚性贮箱的阻尼比。Miles等[7]在Martel等[8]的基础上用内部耗散对阻尼比进行修正并获得更符合实验的结果，王为等[9-10]引入自由面洁净系数并与几种液体的实验结果相比较。夏恒新等[11]计算了带半球形底的刚体圆柱并与实验结果对比。

本文在以上方法的基础上，对三维圆柱形贮箱液固耦合问题进行分析。通过小幅晃动和理想流体假设，用标准的有限元方法获得非对称形式液固耦合方程组。将此方程组经过自由面缩聚和主坐标变换后获得对称形式高度降阶的主坐标下的液固耦合方程组。将经过特征值计算获得的自由面主坐标模态通过热流比拟的方法还原成全流场模态，进而采用边界层理论处理粘性耗散得到模态阻尼比。

1充液贮箱液固耦合振动的模态阻尼比

1.1弹性结构与理想流体的耦合振动模态分析

贮箱液固耦合问题可以近似为无黏、无旋、不可压缩的理想流体在弹性贮箱结构内的小幅线性晃动问题。

贮箱内液体满足Laplace方程，液固接触面上满足运动边界条件，自由面上满足几何边界条件：

▽2p=0 液体所占区域内部

(1)

(2)

(3)

式中：p是液体的压强，n是液体的外法线方向，ρf为液体的密度，un为结构沿外法向的位移。

贮箱结构运动方程：

(4)

式中：Ms,Ks分别为结构的质量、刚度矩阵，qsf为液体作用在结构上的等效节点力向量。

根据液体的控制方程和边界条件及结构运动方程，用有限元加权参数方法获得液固耦合非对称形式的方程：

(5)

式中：Mf,Kf分别是液体的质量、刚度矩阵，L为液体对结构的作用矩阵。

经过自由面缩聚和主坐标变换，得到对称形式的液体与结构耦合振动方程：

(6)

式中，ξE为自由面主坐标，u为结构的位移向量。

为获得模态频率和振型，朱昶帆等[4]将自由面主坐标与结构位移向量时空分离：

(7)

将式(7)代入式(6)获得特征方程：

(8)

求解特征方程(8)获得各阶模态固有频率及对应的自由面主坐标和结构位移的模态振型。

为了求模态阻尼比需要全流场的压力模态振型。

假定压力函数也是时空分离的：

p=Peiωt

(9)

(10)

根据式(7)获得结构模态加速度：

(11)

进行热流比拟。因为控制方程式(1)与温度的Laplace方程形式一致，所以将节点压力模态P比拟为节点温度模态T。自由面上采用固定温度的边界条件，将式(10)中自由面节点压力比拟为自由面节点温度。接触面采用热流输入的边界条件，因为边界条件式(2)与热流输入边界条件形式一致，将式(11)中结构位移加速度比拟为热流输入。通过Nastran内置的热流模块计算出全流场节点温度T，即全流场节点压力P。

1.2黏性耗散及其模态阻尼比计算

根据文献[6]，液体的黏性耗散在不考虑毛细效应时，体现为液固接触面耗散、自由面耗散和内部耗散。根据Stokes边界层假设，在小Re数假设下，液体在边界层内流动遵从扩散方程：

(12)

在流体内部满足理想流体假设，存在速度势函数，自由晃动条件下，将势函数时空分离，同时因为速度势和速度的关系，u0也可以表示为时空分离的形式：

(13)

在液固接触面上一个晃动周期内的平均能量耗散率为：

(14)

根据文献[6]，在表面张力及自由面洁净情况等影响下自由面耗散和接触面耗散是同量级的，可以采用同样的边界层假设，文献[10]中引入系数α，通过的α变化表征自由面耗散情况，当α=0时表征自由面完全无耗散。

在自由面上一个晃动周期内的平均能量耗散率为：

(15)

式中，Γf为液体自由面。

则边界层内总耗散为液固接触面耗散与自由面耗散之和：

D=Dw+Df

(16)

液体的内部耗散：

根据文献[7]，液体的内部耗散和表面耗散相比是可比较的。忽略高阶项时可采用近似公式：

(17)

式中，∂n为外法线方向偏导数。

液体的机械能：

采用奥-高公式，将体积分转为面积分，则液体的一个周期内的平均机械能为：

液体晃动一个周期的阻尼比：

(18)

由于式(18)的模态阻尼比是关于势函数的，而上节推导出的流场模态是关于压力函数，所以需要得到模态阻尼比关于压力函数的阻尼比表达。

拉格朗日积分得到液体势函数与压力函数关系：

(19)

由于压力函数和势函数都是时空分离的，得到

(20)

将式(20)代入公式(18)，得到

(21)

因此在计算阻尼比的时候可以采用压力函数替代势函数。运用有限元方法对理想流体问题得到的流域节点压力进行插值，获得所需的节点一阶及二阶导数，进而得到阻尼比。

综上所述，在计算弹性贮箱液固耦合振动阻尼比时，可采用先将问题简化为理想流体与弹性结构的液固耦合，获得理想流体的流场节点压力，再引入边界层耗散用节点压力及其导数获得阻尼比。

2计算实例

2.1无限长二维弹性底板矩形槽模型

结构为如图1所示的无限长矩形截面槽，截面两侧壁面为刚性，底面为弹性薄板，弹性模量E=200 GPa，密度ρs=7 800 kg/m3，槽内充液的液体密度为ρf=1 000 kg/m3，运动黏度系数ν=1×10-6m2/s。假定自由面为洁净的，洁净系数α=0。在二维情况下，底边无限长弹性薄板的平面应变运动方程可以等效为梁的运动方程。通过截断级数的半解析方法获得流场，并计算阻尼比，与相同几何形状的底面也为刚体模型的解析结果比较。

图1　无限长二维弹性底板矩形槽简图Fig.1 Infinite long 2D elastic bottom rectangular canal

表1、2为弹性底板矩形槽模型与相同几何形状的全刚性壁面矩形槽模型的结果对比。分别对比浅水和深水(充液深度0.2 m和0.8 m)的前三阶情况，对应的模态振型为图2。其中，R1为边界层阻尼比，R2为边界层与内部之和的总阻尼比，i/s为内部耗散Di与边界层耗散D的比值。

表1　浅水时(0.2 m)弹性底板矩形槽与全刚矩形槽比较

表2　深水时(0.8 m)弹性底板矩形槽与全刚矩形槽比较

在充液深度为0.2 m时，第2阶模态，底部弹性模型和刚性模型的频率有较大差别，边界层阻尼比及总阻尼比也有较大差别，从模态振型(图2)可以看出，和第1、第3阶模态相比，第2阶模态液固耦合比较明显，属于耦合振动模态。在充液深度为0.8 m时，虽然弹性底部模型和刚性模型第3阶模态的频率基本相同，但是模态振型图中该阶模态液固耦合明显，在边界层阻尼比及总阻尼比上也具有较大差别。在浅水和深水情况下，在液固耦合不明显的各阶模态中(0.2 m的第1第3阶，0.8 m的第1第2阶)，弹性底部模型和刚性模型阻尼比基本一致，说明在液固耦合不明显的模态上，使用刚性模型进行近似也是合理和比较准确的。但是在液固耦合明显的模态上，使用弹性模型会有较大改进。

图2　弹性底板矩形槽模态振型图Fig.2 Elastic bottom canal modal shapes

2.2三维圆柱形贮箱模型

此模型模态频率和自由面振型取自文献[3]。

三维圆柱形贮箱模型几何尺寸如图3，材料弹性模量E=68 GPa，泊松比σ=0.33，前后短柱壳厚度6.92 mm，前后底扁球壳厚度2.0 mm，中间柱壳厚度4.95 mm。液体自由面高度为8 000 mm，密度ρ=1 000 kg/m3，运动黏度系数ν=1×10-6m2/s，洁净系数α=1。

图3　三维圆柱形贮箱模型简图Fig.3 3D circular cylinder tank model sketch

运动黏度为ν=1×10-6m2/s时，三维圆柱形贮箱模型中，各阶模态边界层阻尼比(R1)均为10的-4次方量级，内部耗散与边界层耗散比值为10的-2次方量级，见表3。在液固耦合不明显的模态上(第5、8、14、19阶)，内部耗散与边界层耗散的比值大约2%～3%，与文献[6]中，刚性壁面假设时，两者比值相近。

表3　三维圆柱形贮箱一些模态阻尼比

2.3三维圆柱形贮箱模型2

模型采用和2.2中相同模型贮箱，在模型设定相同的情况下，将底部材料弹性模量增加10倍，材料弹性模量E=6.8×1011Pa。将刚性底部情况与2.2中的标准情况相比较。

在液固耦合不明显的模态上，刚底情况与标准情况模态和阻尼比基本一致，见表4。说明在液固耦合不明显的模态上，采用刚性假设也能获得较为精确的模态频率和阻尼比结果。

表4　刚底情况与标准情况比较

而在液固耦合振动模态上，通过模态振型比较，得到刚底情况第64、67、91，106阶与标准情况第64、67、91、104阶模态分别对应，见表5。说明在可以比较的模态上，两种情况的模态频率和阻尼比计算结果差别不大。但是，底部刚度对贮箱液固耦合部分模态振型影响较大，会出现部分模态增加和减少的现象。因此在计算三维圆柱形贮箱模态阻尼比时需要考虑弹性贮箱液固耦合的影响。

表5　刚底情况与标准情况比较2

3结论

本文在弹性贮箱液固耦合理论和刚性容器液体小幅晃动阻尼理论基础上，建立了液体在弹性容器内的小幅晃动阻尼计算方法。用热流比拟的方法，处理通过自由面缩聚和主坐标表示的液固耦合动力学方程的模态，得到流场。并用小幅晃动理论获得模态阻尼比。

本文通过比较二维半解析模型与刚体模型，发现在液固耦合模态下的弹性模型与刚体模型的阻尼比有一定的差别，验证了在某些频率上采用液固耦合模型在阻尼比精度上比刚体假设有较大改善。而三维圆柱形贮箱模型的有限元计算表明，应用本文的方法，能够高效、即时获取一定液高下液固耦合模态、流场和模态阻尼比，为工程上更高精度建模进行动力学分析，设计防晃装置以降低液体火箭贮箱晃动等实践提供了分析基础。

参 考 文 献

[ 1 ] 马驰骋，张希农，罗亚军，等. 变质量矩形贮箱流固耦合系统动力学特性分析[J]. 西安交通大学学报，2014，48(7):109-116.

MA Chi-cheng, ZHANG Xi-nong, LUO Ya-jun, et al. Dynamic analysis on a fluid-structure coupling system with variable mass in a rectangular tank[J]. Journal of Xi’an Jiaotong University,2014,48(7):109-116.

[ 2 ] Fox D W, Kuttler J R. Sloshing frequencies[J]. Journal of Mechanics and Applied Mathematics and Physics,1983,34(5):668-696.

[ 3 ] 朱琳，唐国安，柳征勇，等. 推进剂与贮箱液固耦合振动的动力学分析[J]. 振动与冲击，2012，31(5):139-144.

ZHU Lin,TANG Guo-an, LIU Zheng-yong, et al. Dynamic analysis for fluid-structure coupled vibration of propellants and fuel tank [J]. Journal of Vibration and Shock, 2012,31(5):139-144.

[ 4 ] 朱昶帆，唐国安，张美艳. 考虑推进剂晃动的火箭液固耦合分析的比拟算法[J]. 动力学与控制学报，2014，12(3):239-242.

ZHU Chang-fan, TANG Guo-an, ZHANG Mei-yan. Theanalogy method to coupling analysis of rockets comsidering propellant sloshing[J]. Journal of Dynamics and Control, 2014,12(3):239-242.

[ 5 ] Case K M, Parkinson W C. Damping of surface waves in an incompressible liquid[J].J Fluid Mech,1957,2:172-184.

[ 6 ] Henderson D M, Miles J W. Surface-wave damping in a circular cylinder with a fixed contact line[J]. J Fluid Mech,1994,275:285-299.

[ 7 ] Miles J W, Henderson D M. A note on interior vs. boundary-layer damping of surface waves in a circular cylinder[J].J Fluid Mech,1998,364:319-323.

[ 8 ] Martel C, Nicolas J A, Vega J M. Surface-wave damping in a brimful circular cylinder[J].J Fluid Mech,1998,360:213-228.

[ 9 ] 王为，李俊峰，王天舒. 航天器贮箱内液体晃动阻尼研究(一)[J]. 宇航学报，2005，26(6):687-693.

WANG Wei, LI Jun-feng, WANG Tian-shu. Study on sloshing damping in a contanier of spacecraft(Part I): theory[J]. Journal of Astronautics, 2005, 26(6):687-693.

[10] 王为，李俊峰，王天舒. 航天器贮箱内液体晃动阻尼研究(二)[J]. 宇航学报，2005，27(2):177-181.

WANG Wei, LI Jun-feng, WANG Tian-shu. Study on sloshing damping in a contanier of spacecraft(Part Ⅱ): theory[J]. Journal of Astronautics, 2005, 27(2):177-181.

[11] 夏恒新，王为，宝音贺西,等. 带球形底的圆柱容器中液体小幅晃动的实验研究和有限元分析[J]. 动力学与控制学报，2009，7(1):66-70.

XIA Heng-xin, WANG Wei, BAOYIN He-xi, et al. Experiment and FEM analysis of liquid sloshing with small amplitude in a circular cylindrical tank with a spherical bottom[J]. Journal of Dynamics and Control, 2009, 7(1):66-70.

Fluid-structure coupled vibration modes and damping ratios for a liquid-filled contained

SHENJi,ZHANGMei-yan,TANGGuo-an

(Department of Mechanics & Engineering Science, Fudan University, Shanghai 200433, China)

Abstract:When conducting the stability analysis for longitudinal vibration of a liquid-carried rocket and establishing the vibration model of liquid inside the tank, the elastic characteristics of the structure should be considered, the damping character of liquid shaking should also be considered. Due to the fact that a modal damping ratio is correlated to the work of damping force per vibration period, mean mechanical energy, and the modal natural frequency, the fluid-structure coupled vibration modes were determined via analyzing the coupled vibrations of structure and liquid, and then the works done by the damping force on the contact surface between liquid and structure that on the free surface and the inner damping force were calculated, the mean mechanical energy was also computed, so as to achieve the modal damping ratio. The heat flux analogy method was adopted in the modal analysis of coupled vibration and the calculation of liquid inner fluid-field to take full advantages of a general program. This method not only improved the efficiency, but also facilitated engineering applications.

Key words:liquid-structure coupled vibration; tank; mode; damping ratio

中图分类号:O327；O353

文献标志码：A

DOI:10.13465/j.cnki.jvs.2016.04.015

通信作者唐国安 男，教授，博士生导师，1962年10月生

收稿日期：2014-12-19修改稿收到日期：2015-03-07

基金项目：国家自然科学基金(11202052)

第一作者 沈骥 男，硕士生，1989年12月生