基于LS_DYNA的弹体撞水过程流固耦合动力分析

引 言
1
结构物撞水的研究，对于解决鱼雷、反潜导弹、深水炸
弹、破障炮弹入水的运动，宇宙飞船返回舱的水上回收，水
上飞机的起降等工程实际问题具有重要意义
[1-2]
。当结构物
与流体发生撞击作用时，不但流体会对结构物产生强大的撞
击作用力，结构物的运动和形变也会对流场产生显著的影
响，因此两者之间存在着强烈的流固耦合作用。而高速弹体
博导，研究方向为飞行器导航、制导与控制；林鹏(1978-), 男, 辽宁本
溪人, 讲师, 博士后，研究方向为飞行器结构设计、导弹先进控制理论
与应用。
单元网格与结构是重合的，网格随着结构的变形而变形。但
是，对于流固耦合问题，材料流动将造成有限元网格的严重
畸变，引起数值计算困难，导致仿真计算无法完成。为此，
考虑到弹体的高速运动，这样得到的流场区域往往很大，为
了将求解时间控制在可以接受的范围内，需要对时间步长进
行严格控制。
对于 LS-DYNA 程序的显式时间积分算法，若给定了模
型的材料特性，则仿真的最小时间步长由最小单元尺寸控
制：对于给定的网格划分，
min
Δt 取决于材料的纵波速，它
是材料性质的函数（密度、弹性模量和泊松比）。由于仿真
我们采用 LS-DYNA 程序的 ALE (Arbitrary Lagrange-Euler )
列式来解决这个问题。ALE 列式先执行一个或几个 Lagrange
时步计算，此时单元网格随材料流动而产生变形，然后执行
ALE 时步计算：首先，保持变形后的物体边界条件，对内
部单元进行重分网格，网格的拓扑关系保持不变，称为
学者认为是一种分析流固耦合的有效方法
[5-8]
。
1 基于 ALE 网格的流固耦合方法
固体结构物的有限元计算通常采用 Lagrange 列式，其
收稿日期：2009-10-14 修回日期：2010-03-01
作者简介：孙琦(1982-), 男, 湖南益阳人, 博士生, 研究方向为飞行器结
构设计、导弹控制系统设计与仿真；周军(1966-), 男, 江苏常州人, 教授,
and the water-gas interface is reasonable. Pressure applied from fluid on projectile is analyzed, important references for the
design of projectile structure. Projectile attitude and position variety under different initial state is compared, which is
第 22 卷第 6 期
系 统 仿 真 学 报©
Vol. 22 No. 6
2010 年 6 月
Journal of System Simulation Jun., 2010
• 1498 •
基于LS-DYNA的弹体撞水过程流固耦合动力分析
孙 琦, 周 军，林 鹏
（西北工业大学精确制导与控制研究所，西安 710072）
• 1499 •
和
L
∂X ∂t 分别表示网格点和流体质点的速度。
基于 ALE 网格结构，在 LS-DYNA 程序中可以很方便
的实现欧拉—拉格朗日耦合算法，从而实现流体—结构的耦
合计算。采用流固耦合算法模拟问题时，往往要对拉格朗日
结构进行约束，将结构的相关参量传递给流体单元。我们采
用罚函数约束的方法来实现耦合。罚函数耦合系数追踪拉格
Abstract: Multi-material ALE method was used to fluid-structure interaction dynamic analysis on water impact of
three-dimensional elastic projectile, by the explicit formulation dynamic analysis finite element software LS-DYNA. The
朗日节点（弹体结构，即从物质）和欧拉流体（水和空气，
即主物质）位置间的相对位移 d ，检查每一个从节点对主物
质表面的贯穿，如果从节点不出现贯穿，就不进行任何操作；
如果发生从节点对主物质表面的贯穿，界面力 F 就会分布
到欧拉流体的节点上。界面力的大小与发生的贯穿数量成正
比，即
i
F = k ⋅ d(3)
式中
i
k 表示基于主从节点质量模型特性的刚度系数。
2 弹体撞水计算模型
在进行弹体撞水动力分析时，我们所主要关心的是弹体
结构受到的流体反作用的压力峰值，以便为提高弹体结构的
抗撞强度提供一个参考依据；以及在入水过程中，弹体位置
和姿态的变化情况，通过选择合适的入水初态，使整个撞水
过程中，弹体姿态满足一定的要求。
influences of key control parameters were discussed. Numerical simulation results indicate that, this model can solve the
numerical computation difficulty originated from fluid great shape change in fluid-structure interaction dynamic analysis,
—应变关系仍然为线性，但是斜率发生变化。卸载曲线与加
载段曲线斜率相同，这样当完全卸载（应力等于零）后，材
料中将保留永久的塑性变形
p
ε 。一维情况下弹塑性材料本
构的描述，首先判断结构的应力状态是否达到屈服应力，如
果没有达到，则按照线弹性材料本构进行处理。如果材料内
的应力已经超过屈服应力，则按照塑性变形本构计算结构中
Smooth Step；然后，将变形网格中的单元变量（密度、能量、
应力张量等）和节点速度矢量输运到重分后的新网格中，称
为 Advection Step。
对于 ALE 网格，可以理解为有两层网格重叠在一起，
一层是空间点网格，另一层是附着在材料上并随着材料在空
间网格中运动的网格，设流体质点 K ( a , b , c )和网格点
弹体姿态和位置变化情况，可以为弹体撞水初态的选择提供依据。
关键词：流固耦合；撞水响应；LS-DYNA；数值模拟
中图分类号：TP391.9 文献标识码：A 文章编号：1004-731X (2010) 06-1498-04
Dynamic Analysis of Fluid-Structure Interaction
其时间积分采用显式中心差分方法，基本格式如下：
1
1/2
( ) [ ( ) ( ) ( ) ( )]
n n n n n
x t M P t F t H t Cx t
的应力—应变。在三维条件下，判断材料是否进入塑性可以
使用 V.Mises 屈服准则，即
2
1 1
0
2 3
ij ij s
s s − σ= (4)
式中，
ij ij m ij
s = σ −σ δ
为斜应力张量；
11 22 33
1
( )
3
m
σ = σ + σ + σ为平均应力。
I ( i , j , k )的运动可用下式分别表示：
( )
,
L
X = X K t， ( )
,
A
X =X I t
取 f 为流场内某一个物理量，它是坐标与时间变量的连续函
数，则 f 对网格 I 与质点 K 的随体导数可分别表示为
d
d
A
A
X
f f f X
t t X t
I
公式为
2 3 2
0 1 2 3 4 5 6
p = c + c μ + c μ + c μ + ( c + c μ + c μ)E(6)
2.3 显式分析求解
在建立弹体撞水的有限元动力分析模型时，外流场需包
含整个弹体运动区域，我们将流场边界到弹体的最近距离取
为弹体相应尺寸的 3~5 倍，流场外界采用无反射边界条件，
helpful to design appropriate initial attitude in water impact of projectile.
Key words: fluid-structure interaction; water impact of projectile; LS-DYNA; numerical simulation
2.1 弹体材料模型
在材料本构模型中，主要有弹性、塑性、粘性以及三者
的组合。考虑到弹体高速入水时承受的巨大撞击力，有时需
要通过一定的结构来进行缓冲，存在结构的变形及失效问
题，因此，我们采用弹塑性材料来对弹体进行建模。
弹塑性材料在加载段应力与应变保持线性，当应力大于
屈服应力
s
σ 时，材料进入塑性，此后如果继续加载，应力
的冲击入水，还需要考虑水—气交界面的变化情况，对于这
样具有复杂边界条件的实际工程问题，很难进行解析求解，
借助于计算机的数值分析方法进行近似求解是比较有效的