《红外诱饵弹干扰特性有限元建模.doc》由会员分享，可在线阅读，更多相关《红外诱饵弹干扰特性有限元建模.doc（14页珍藏版）》请在三一办公上搜索。
1、精品论文红外诱饵弹干扰特性有限元建模朱敏，杨春玲，张国亮（哈尔滨工业大学电气工程及自动化学院）5摘要：针对国内诱饵弹半实物仿真技术和红外制导技术研究的迫切需要，提出针对诱饵弹干 扰特性的有限元建模方法。通过分析诱饵弹的工作机理和运行状态，借助有限元方法建立诱 饵弹的温度流场分布和空气变阻力运动特性模型，并通过仿真实验验证所提出模型的有效性 和正确性，该方法不仅能为红外诱饵弹的设计和红外制导算法的研究提供依据，也对红外半 实物仿真系统的搭建有重要指导意义。10关键词：红外诱饵弹；干扰特性；有限元建模；仿真中图分类号：TN216The infrared decoys interference ch
2、aracteristics of finite element modeling15zhumin, yangchunling, zhangguoliang(School of Electrical Engineering and Automation, Harbin Institute of Technology)Abstract: An urgent need for domestic decoys loop simulation technology and infrared guidancetechnology research, a finite element modeling me
3、thod interference characteristics for decoys is proposed in this paper. The mechanism and the operating status of the decoys are analyzed, the20flow field of the decoys are solved accurately based on finite element theory. Then the changes of air resistance in the process of decoys sports movement a
4、nd the distribution of temperature fieldare studied in the paper, and the decoys movement characteristics and radiative properties of the finite element model are setted up, which provides a theoretical basis for the hardware-in-loop simulation. And simulation verify the validity and correctness of
5、the proposed model. The method25can not only provide the basis for the infrared decoys design and infrared-guided algorithmresearch, but also have important guiding significance of infrared decoy HIL simulation system.(10 Points, Times New Roman)Key words: the infrared decoys; interference character
6、istics; finite element modeling; simulation300引言目前，国外相关机构已经在红外诱饵弹特性模型和仿真系统的研制方面展开工作。如美 国 Naval Air Warfare Center 建立的红外成像仿真实验室，美国 CSA 实验中心（Computer Scienceand Application Inc）设计开发的红外景像模拟系统 IRSP 等，都已经得到实际应用1,2。 AIISPES 系统是美国陆军在 20 世纪末建立的武器性能评价系统，可以对常规飞行武器、红35外对抗武器、机载舰载诱饵武器等多种红外制导反导武器进行打靶实验和仿真测试3。英国 C
7、hemring 防务公司为便携式防空武器开发了专业的目标特征模拟系统和配套的仿真软件 COUNTERSIM4。该系统可以真实的再现仿真目标的工作状态。20 世纪 80 年代，以色列研 制并生产出一款红外模拟器，该模拟器专门用于动态目标的红外图像仿真，可以对不同背景 下的红外目标进行实时模拟5。40虽然国外已经对红外目标半实物仿真系统进行了大量的研究，但针对仿真目标，特别是基金项目：高等学校博士学科点专项科研基金（编号 20092302110026）作者简介：朱敏（1981-），男，讲师，主要研究红外图像的仿真，模式识别等通信联系人：杨春玲（1965-），女，教授，主要研究智能图像处理，模式识别
8、等. E-mail: yangcl1- 14 -红外诱饵弹的自身特性研究比较少，现有的仿真系统多是使用简化仿真模型对其运功和辐射特性进行仿真，使得仿真系统的可信度有所降低6-9。如美国空军研究实验室（Air Force Research Laboratory）对红外诱饵弹特性进行了研究，提出了诱饵弹运动的简化模型和单波 段红外辐射模型10。弗罗茨瓦夫理工大学（Wroclaw University of Technology）和波兰空军45材料研究院（Air Force Material Command），对直升机自我防卫系统进行了研究，建立了直 升机和机载诱饵等防卫武器在天空背景下的红外辐射模
9、型11。国内对红外诱饵弹的特性也有一定的研究。在诱饵弹特性建模方面，空军工程大学、西 北工业大学和东北电子技术研究所等高等院校和科研机构，对诱饵弹的运动规律及辐射特性 进行了研究和仿真12-17。空军工程大学的刘加从等人，对诱饵弹运动简化模型进行了研究，50着重讨论了诱饵弹运动过程中的受力情况18。东北电子技术研究室的研究人员给出了面源 式诱饵弹的辐射特性及其仿真模型，并对该类型的诱饵弹干扰性能进行了评价19。空军工 程大学的张文华等人，以匀直流点源模型为基础，建立了诱饵弹的流体力学模型，并推导出 此模型下诱饵弹的相关气动参数，同时提出了等效分子假设，给出了诱饵弹的辐射特性模型 20。55总之
10、，国内各单位在诱饵弹半实物仿真技术的研究上还处于起步阶段，特别是在诱饵弹 精确模型的建立和新型红外模拟器的研制方面，与国外先进水平有很大的差距。因此，有必 要对诱饵弹半实物仿真技术进行更深入的研究，以适应制导武器的发展需要。本文对诱饵弹 的流场分布进行研究，建立诱饵弹空气动力特性的有限元建模，研究诱饵弹的运动规律和辐 射特性，并建立相应的干扰特性模型。601诱饵弹气动特性有限元分析1.1诱饵弹工作流场的分析与简化诱饵弹工作时，内部燃料迅速燃烧，产生大量的高温燃气，喷射到弹体周围。这一特性 使得诱饵弹周围流场的分布变得十分复杂，不同于一般弹丸运行时周围流场的匀直流情况。 根据空气动力学基础理论，
11、诱饵弹运动时，周围空气流场的分布满足流体力学基本方程组： + v=0t65de +p d ( 1 )=q(1)dtdt R dv = F -p+2 v dt式中 为气体密度，v 是气体流速，e 是气体内能， p 是压力张量，qR 是气体辐射热，F 是 气体所受外力， 是气体粘性系数。自然流区和燃气扰动区的流场有相同性质：首先，流场的粘性可忽略，整个流场为无旋场。其次，流场的参数不随时间变化，是定常流场。最后，流场中气体的外力很小（主要是70重力），可忽略不计。因此，可以对流场的基本方程进行化简：忽略连续性方程中的外力和 时间导数，式(1a)可表示为： ( v ) + ( v ) + ( v )
12、 = 0xx yy zz在直角坐标系中，无粘性气体的欧拉运动方程如式(3)所示：(2)v vx +vvx +vvx =- 1 p x xy yz z x75vvy +vvy +vvy =- 1 p(3) x xy yz z y v vz +vvz +vvz =- 1 p该方程可以缩写成(4)的形式：x xy yz z z根据声速公式，可推导出：(v )v= - 1 p(4)1 (v ) = v 1 = 1 v p= - 1 v (v )v(5)c2 c280将式(5)代入式(2)，可以得出诱饵弹无粘性流动区域中流速所满足的微分方程：v2 vv2 vv2 v(1- x ) x +(1- y )y
13、 +(1- z ) z c2xc2yc2 zvx vy vvy-(x +)-vy vz ( vz + vy )(6)c2yxc2 yz- vz vx ( vx + vz ) = 0c2 zx由前面的讨论可知诱饵弹流场是无旋场，根据场论的有关知识，在无旋场中可以引入速度势函数 ，使其满足以下条件： vx = x = x，vy = y = y，vz = z =z85将式(7)代入式(6)，整理可得用速度势函数表示的流场运动方程：(7) 2(1- x )c2xx+(1- 2y )c2yy 2+(1- z )c2zz(8)- 2 x y c2xy - 2 y z c2yz - 2z x c2zx =
14、0由于诱饵弹的轴向对称性，可以将式(8)中 Z 轴方向的分量消去，从而得到诱饵弹流场的二 维运动方程： 2 2 x y x y y(1-) +(1-) - 2c2 xx c2 yy c2 xy += 0y(9)90式中的 y / y 项是由于取 X 轴为对称轴而引入的。至此，已推导出诱饵弹流场关于势函数的 运动方程，只要求出势函数 ，就能利用式(7)的约束条件求出流场中的流速，然后再利用理想气体方程求出压强、温度等其他气动参数。1.2诱饵弹流场有限元模型的建立用速度势函数表示的流场运动方程本质上依然是二阶非线性微分方程，求解起来有很大95难度。为了求解流场运动方程，本文引入了有限元方法21-2
15、5。1.2.1诱饵弹流场的有限单元剖分本文求解区域是以诱饵弹表面为内边界，以弹体直径 5 倍长度为半径的空心圆柱区域， 子单元为任意四边形。对诱饵弹流场解算区域进行有限单元剖分，剖分后的求解区域共有1715 个单元，3514 个结点，剖分结果如图 1 所示。100105110图 1 求解区域有限元剖分结果.Fig.1 Solving regional finite element split results.1.2.2单元插值基函数的确定 借助标准等边四边形的基函数选取方法，再通过变换得到任意四边形单元的基函数表达式。变换过程如图 2 所示。(a) 实际坐标系中的任意四边形(b) 标准坐标系中
16、的标准四边形 图 2 任意四边形的标准化.Fig.2 Standardization of arbitrary quadrilateral.根据图 2，在实际坐标系 XOY 中，取任意一个四边形 ABCD，其结点坐标分别为A( x1 , y1 ) ， B( x2 , y2 ) ，C (x3 , y3 ) ， D(x4 , y4 ) 。建立一个与之对应的标准坐标系O ，将四边形 ABCD 映射成一个标准的等边四边形 ABCD，各结点坐标为 A(-1,1)，B(-1,-1)，C(1,-1)，D(1,1)。那么这两个坐标系之间的变换关系如式(2-10)所示： x = 1 (a+ b + c + d
17、)41 1 1 1(10) y = 1 (a+ b + c + d )115其中各系数满足如下关系：42 2 2 2 a1 = x1 + x2 + x3 + x4 b1 = -x1 + x2 + x3 - x4 c1 = -x1 - x2 + x3 + x4 d1 =x1 - x2 + x3 - x4(11) a2 = y1 + y2 + y3 + y4 b2 = - y1 + y2 + y3 - y4 c2 = - y1 - y2 + y3 + y4d = 2 y1 - y2 + y3 - y4在标准坐标系中，等边四边形的插值基函数有多种取法，只要是相互正交的函数即可。这里 取最为常用的线性
18、插值基函数，其形式如式(12)所示： = 1 (1 - )(1 - ) 14 = 1 (1 + )(1 - )120 24 = 1 (1 + )(1 + )(12) 34 = 1 (1 - )(1 + ) 44根据式(10)和(12)，可以得到标准基函数对 x，y 的导数： i =1 ( A+ B + C ) x8 J xi xizi(i = 1, 2, 3, 4)(13) i =1 ( A+ B + C ) y8 J其中各系数满足以下关系：yi yi yi Axi y2 -y4y3 -y1y4 -y2y1 -y3 Bxi y4 -y3Cxi y3 -y2=y3 -y4y1 -y4y1 -y2
19、y4 -y1y2 -y1 y2 -y3 (14) Ayi x4 -x2 Byi x3 -x4x1 -x3x4 -x3x2 -x4x2 -x1x3 -x1 x1 -x2 Cyix2 -x3x4 -x1x1 -x4x3 -x2125将标准坐标系下的基函数及其导数的表达式带入到有限元计算方程中，便可以求得标准坐标系下的有限元解计算结果，然后坐标变换，即可得到实际坐标系中的计算结果。1.2.3有限单元分析 将式(9)进行积分，可以得出诱饵弹流场势函数的单元积分表达式： 2 2 x y x y y(1- c2 ) xx +(1- c2 ) yy - 2 c2xy + 2 yd = 0y(15)130其中
20、 表示子单元的求解域， 表示子单元的近似解。插值函数i 和结点值i 的线性组合 表示子单元的解 如式(16)所示： =iii=1,2,3,4(16)将式(16)代入式(15)，并把实际坐标系中的积分变量(x，y)替换成标准坐标系中的( ， )，可以得到带有系数矩阵的近似解方程：1351 12H(， ) - 1 G (， ) yd d = 0(17)-1 -1ijc2 ij i其中 Hij (， ) 和 Gij (， ) 的表达式如如(18)所示：H (， ) = 1 i + i j + i j ijy dyi x x y y (18) G (， ) = ( i )3j + ( i )3j +
21、2( i )3 j ij xx记yyxyF (， ) = H(， ) - 1 G (， )ijijij1 1c2 ij(19)140则式(17)可以表示为：Ae = 2-1 -1Fij (， )yd de 0Aiji =(20)其中 Ae 称为有限单元积分方程的系数矩阵，上标 e 表示 Ae 是第 e 个单元的系数矩阵，e=1，ij ij1452，N，N 为结点总数。将坐标变换关系式(10)代入式(19)，可以得到有限单元系数矩阵 的最终形式。1.2.4总体有限元方程的建立ij在求出每个单元的系数矩阵 Ae 之后，要通过总体合成，求出整个求解区域所对应的总 体系数矩阵 Aij ，然后，以所有结
22、点值i 为未知量，建立求解区域的总体有限元方程：Aiji = 0(21)150假设某一子单元 e，其结点 i、j 对应的总体结点序号为 m，n，在总体合成时按如下规 则进行：A= A + Ae nmnmij(22)在求解域中任取两个相邻单元 e1 、e2 ，其总体结点序号为 n1 至 n6 ，单元 e1 的结点i1 、i2 、i3 、i4 分别对应总体结点 n1 、n2 、n3 、n4 ，单元 e2 的结点i1 、i2 、i3 、i4 分别对应总体结点 n4 、n3 、 n6 、 n5 ，具体结构如图 3 所示。155图 3 两节点合成示意图Fig.3 Enter figures captio
23、n herenm根据式(22)，可以得到两个单元进行总体合成以后的系数矩阵 A ，如式(23)所示： A11A12A13A14A15A16 A21 AAA=31A22A32A23A33A24A34A25A35A26 A36 nm AAAAAA 41 42 43 44 45 46 51A52A53A54A55A56 A61A62A63A64A65A66(23) A(1)A(1)A(1)A(1)00 11 12 1314 A(1)A(1)A(1)A(1)00 21 22 23 24 A(1)A(1)A(1) +A(2)A(1) +A(2)A(2)A(2) = 31 323322 34 21 24
24、26 A(1)A(1)A(1) +A(2)A(1) +A(2)A(2)A(2) 41 42 43 12 44 11 14 16 00A(2)A(2)A(2)A(2) 42 41 44 43 00A(2)A(2)A(2)A(2) 32 31 34 33 160165最终的总体有限元方程是一个含有 3514 个未知数的线性方程组。对于诱饵弹的工作环境来讲，影响其流场分布的主要边界条件有两个：一是诱饵弹的运动速度，一是燃气出口状态。 在不同边界条件下，通过编程求解总体方程，可以得到流场速度势函数在每个结点的值，根 据式(7)对势函数求导，便可以求出流场的速度分布 v(x，y)。表 1 是取诱饵弹运动
25、速度 300m/s， 燃气出口压强为 3 个标准大气压，温度 1000K 时，总体有限元方程计算结果的一部分，图 4 是用后处理软件绘制的诱饵弹流场速度分布。表 1 有限元方程求解部分结果Tab.1 The finite element equation solver partial results结点序号X 坐标(m)Y 坐标(m)vx(m/s)vy(m/s)10.050.55350.000.0020.050.53350.510.0030.050.51350.700.1040.050.49350.820.121000.310.55369.284.591010.290.55367.534.13
26、1020.270.55365.503.431030.250.55363.352.6035110.630.17372.8581.8135120.640.15411.1438.9735130.660.17382.5488.5135140.680.19428.4737.37170175180图 4 诱饵弹流场速度分布Fig.4 The decoys flow velocity distribution从图中可以明显的看出，诱饵弹流场速度分布式不均匀的，弹体前后表面有很大的速度 差，这是诱饵弹释放的燃气与自然流空气相互扰动而引起的。诱饵弹的流场分布情况直接影响其运动和辐射特性。2诱饵弹运动模型的建立及
27、特性仿真2.1诱饵弹空气阻力的离线求解为了提高仿真系统的实时性，本文采取“离线求解，在线拟合”的方法，将不同边界条 件下的空气阻力在仿真之前，离线地求解出来，把结果编制成若干数据表，然后通过数据拟 合得到诱饵弹阻力与速度的函数关系式，由仿真程序调用。这种方法将用时较长的求解过程 和实时性要求很高的仿真过程分离开来，将复杂的求解转换成了简单的函数关系式，缩短了 求解时间，提高了系统的实时性。表 2 给出了不同速度条件下，诱饵弹空气阻力的计算结果。根据数据表给出的数据，利用最小二乘拟合，可以得到阻力与速度的函数关系，如式(24)所示：185f (v) = 3.2e9v5 5.5e6v4 +(24)
28、d0.0033 e8 v3 0.88v2 + 1.1e2 v 4.3e3表 2 空气阻力与速度关系.Tab.2 Air resistance and velocity relations.速度80130180230280330380430480530(m/s)空气阻力41896713771905273635335332585463516524(N)190图 5 则给出了拟合结果与原始数据的对比，由图可以看出，拟合结果与原始数据非常接近。800070006000原始数据 拟合曲线空气阻力（N）500040003000200010000100 150 200 250 300 350 400 450
29、 500速度（m/s）图 5 阻力数据拟合结果对比.Fig.5 Resistance data fitting results compared1952.2诱饵弹运动方程的建立及仿真在求出诱饵弹所受空气阻力 fd 以后，可以根据牛顿运动定律建立诱饵弹的运动方程组。 首先建立诱饵弹受力分析的笛卡尔坐标系：以地面上某定点为参考点，令 X 轴与地面平行， 指向观测方向，Y 轴竖直向上，Z 轴与 X、Y 轴满足右手定则。诱饵弹在该坐标系内受力分 析如图 6 所示：YvfdO mg X200205Z图 6 诱饵弹受力分析Fig.6 Decoys force analysis其中 fd 为空气阻力，方向与
30、速度 v 相反，mg 为重力，方向竖直向下。不妨设诱饵弹与地面参考系的俯仰角为 ，偏航角为 。将阻力 fd 向坐标系的三个坐标轴投影，根据牛顿运动定律，可以得到诱饵弹加速度的分量表示形式： dvx = f cos( )cos( ) / m dt dvy dtd= fd sin( ) / m g(25) dvz = f dtdcos( ) sin( ) / m将微分方程转化为差分方程，便可根据运动学定律得到诱饵弹运动方程组： ax = - fd cos( ) cos ( ) / m ay = - fd sin( ) / m - g az = fd cos( ) sin( ) / m vt = v
31、t -1 + a t x x x vt= vt -1 + a t(26) y y y vt = vt -1 + a t z z zx X t = X t -1 + vt -1 t Y t = Y t -1 + vt -1 ty Z t = Z t -1 + vt -1 tz210取合适的仿真步长 ，通过迭代计算可以得到诱饵弹每一时刻的速度和位置，图 7 给出了不同初始条件下，诱饵弹运动速度的仿真结果。图 8 是位置的仿真结果。图 7 不同初始条件下速度仿真结果Fig.7 Different initial conditions. speed simulation results215图 8 不
32、同初始条件下位置仿真结果Fig.8 The location of the simulation results under different initial conditions3诱饵弹干扰辐射特性模型建立与仿真T诱饵弹的辐射特性遵循普朗克定律26,27，该定律描述了辐射度与波长和温度 T 的函数 关系式，以波长表示的普朗克公式为：220M ( ,T ) = c1 1(27)式中 h普朗克常数，b 5 ec2 / T 1h = ( 6.6260755 0.0000040) 1034 J s；k 玻尔兹曼常数，k = (1.380658 0.000012) 1023 J K-1；c1 第一辐
33、射常数，1c = 2 hc2 = (3.741774 0.0000022) 1016 W m2；c2 第二辐射常数，c2 = hc / k = (1.4387869 0.00000012) 102 m K。一般的红外探测器只对某些波段的红外辐射比较敏感（通常在 2.43.1m 或 4.34.55m 范围内），因此仿真时只需计算相应波段的辐射度即可，此时，目标的总辐射度可表示 为： c1c T 4 (c T )3 c 225M= 1 d = 1 2 d 2 (28)0 0 5 ec2 / T 14 c / T c2 / T c2 e 12设 x = c2 / T ，上式可化简为：=M 1 21c
34、 T 44c2 / 1Tc / T3xxdx(29)2 2c2 e 1230其中 1 和 2 分别表示波段的上下限。目标表面的温度分布 T(x)与目标燃料的类型和燃烧情况有关，根据干扰弹的具体型号来确定。由于干扰弹的目标较小，为了便于仿真计算，可以用平均温度Ts 代替各点实际温度，代入式(29)并化简，便可求得目标总辐射强度 I0 与平均温 度Ts 的关系：1I = G(T ) = 02 c11d S cos(30)0s 5 ec2 / Ts 10 0235其中 0 为弹体表面的发射率，与干扰弹表面材料的特性有关； S0 为干扰弹的截面积；0 为探测方向与红外辐射面法向量的夹角。仿照诱饵弹空气
35、阻力离线求解的过程，可以将诱饵弹辐射强度与温度的积分关系拟合成简单的多项式函数：13 57 4 3 2 3I0 = 6.2et + 1.1e t 0.00018t+ 0.11t 30t + 3e(31)图 9 给出了拟合函数与原始数据的对比，由图可以看出，拟合结果与原始数据非常接近。16000原始数据14000 拟合结果12000辐射强度（W/Sr）1000080006000400020000-20003004005006007008009001000温度（K）240图 9 辐射强度数据拟合结果对比Fig.9 Radiation intensity fitting results contrast目标产生的红外辐射经过大气衰减，最终被红外探测器所检测到。大气对不同波段红外辐射的衰减可以用衰减系数 ( ) 表示，红外光谱透过率 ( ) 和衰减系数 ( ) 之间的关系可用布盖尔-朗伯定律表示： a ( ) = exp(a( ) R)(32)245式中 R 为目标与红外系统之间的距离， 为波长。在仿真波段内，目标的平均透过率和平均衰减系数为：12 a = 1
