陈晓丽，牛 禄，潘科玮，王志新，陈浩然
(上海航天动力技术研究所，上海 201109)
固体姿轨控发动机利用燃气直接横向喷射,获得作用于飞行器的横向推力,从而产生机动控制力或力矩,可应用于导弹[1]、逃逸系统[2]等。喉栓式可调喷管是固体姿轨控发动机实现推力调节的重要部件,通过调节喉栓位置,可控制工质流量,实现推力的调节。
20世纪60年代,美国开始研究喉栓式推力可调固体发动机,并将喉栓式可调喷管应用于标准-3导弹的推力连续调节固体姿轨控动力系统[3], 2013年完成了鉴定试验,2014年标准-3 Block IB导弹进入海军正式部署阶段。国内在喉栓式可调喷管的理论分析、数值模拟和试验研究等方面也开展了较多研究[4-6],但工程应用尚有差距。喷管推力是喉栓式可调喷管的重要性能,也是固体姿轨控发动机推力调控系统设计的基础。喉栓式可调喷管的推力常采用数值模拟的方法获得[7-9],该方法在获得喷管推力等性能的同时,还能对喷管内流场进行分析,但固体姿轨控动力系统存在诸多非稳态,喷管将面临多种工况,增加了数值模拟的工作量及工作难度,数值模拟的大量处理及计算工作也不利于喷管的设计优化。此外,现有喉栓式可调喷管的推力测量手段较少,测量难度较大,也给喉栓式可调喷管的推力性能研究带来困难。
针对喉栓式可调喷管,建立了几何喉部面积计算模型和理论推力模型,研究几何喉部面积和喷管理论推力随喉栓位置的变化情况,并建立了三维仿真模型,对喉栓在不同位置下的喷管内流场进行了数值模拟,分析喉栓在不同位置时的喷管内流场分布,研究喷管质量流率、推力随喉栓位置的变化规律,与理论推力进行了对比。建立了喉栓式可调喷管气动等效喉部面积的计算方法,对等效喉部面积计算模型和推力模型进行了优化,并进行了冷气试验验证。研究结果可为后续喉栓式可调喷管的设计及优化、固体姿轨控发动机推力调控系统的设计提供参考。
1.1 喉栓式可调喷管工作原理
喉栓式可调喷管的结构如图1所示,进气段与喷管轴线成一定角度,喉栓与喷管内壁形成环形气体通道,即等效喉部,喉栓式可调喷管通过改变喉栓的位置改变等效喉部面积,从而控制喷管输出推力。喉栓处于初始位置时,喷管开度最大,此时喷管的等效喉部面积最大,喷管的推力达到最大值,喉栓向外移动,喷管开度减小,通过的工质流量减小,喷管推力减小,当喉栓达到最大位移时,喷管完全关闭,喷管推力下降为0。
图1 非同轴喉栓式可调喷管示意图Fig.1 Non-coaxial pintle flow regulation nozzle schematic diagram
在固体姿轨控发动机推力调控系统的设计及仿真研究中[10-13],等效喉部面积计算模型和推力模型是整个系统的重点,等效喉部面积取喉栓和喷管型面构成的几何喉部面积,即喷管流道面积最小处的面积,推力取基于等熵流动假设的理论推力。
1.2 几何喉部面积计算模型
当喷管收敛段和喉栓头部型面均为圆弧时,设喷管轴线为x轴,方向如图2所示,O1,O2分别为喷管收敛段和喉栓头部型面圆弧中心,其半径分别为R1,R2,圆心初始坐标为(xo1,yo1),(xo2,yo2),设喉栓向外运动了x(0 ≤x≤xmax,xmax为喉栓伸入喉部的最大位移,此时喉栓型面与喷管型面相切,喷管处于关闭状态),则O2坐标变为(xo2+x,yo2),A,B分别圆O1和圆O2上的一点,坐标为(xA,yA),(xB,yB),RA,RB分别为A,B到喷管轴线的距离,则其环道面积S为[10]:
(1)
S的最小值即喷管的几何喉部面积At,At的求解是一个有界的非线性规划问题,可以采用内点法进行求解。At是x的函数,且随着x近似线性变化[9]。
图2 几何喉部示意图Fig.2 Geometric throat schematic diagram
1.3 理论推力模型
喷管的推力是喷管内、外表面所受气体压力的合力,根据动量守恒定律,可得推力公式[14]:
(2)
基于一维等熵流动假设,喷管质量流率和喷气速度分别为[14]:
(3)
(4)
(5)
式中:R,k,Tf,pc分别为气体常数、比热比、推进剂绝热燃烧温度和喷管入口截面处的压强。代入推力公式得喷管的理论推力F为[14]:
F=CFpcAt
(6)
(7)
入口压强pc、喉部面积At、出口面积Ae已知,求得喷管出口压强即可得到推力大小。膨胀压强比与喷管面积比关系为[14]:
(8)
根据压强比和面积比的关系,可由喷管面积比迭代求得对应的膨胀压强比,进而得到喉栓在不同位置时的喷管推力输出。
设喷管入口压强为9 MPa,出口压强为大气压强,工质比热比为1.25,喉栓从位置0运动到xmax处,喷管几何喉部面积和理论推力等性能的变化如图3所示。其中,F0和At0分别为喉栓处于初始位置时喷管的理论推力和几何喉部面积。喷管理论推力和几何喉部面积在喉栓处于初始位置时达到最大值,两者随喉栓向外运动而近似线性减小,几何喉部面积在喉栓达到最大位移处减小为0,而理论推力在喉栓达到最大位移前就已减小为0,此后出现负值,这是由于上述推力模型仅适用于喷管的完全膨胀、欠膨胀和轻微过膨胀状态,在喉栓运动后期,喷管扩张比剧增,喷管处于过膨胀状态,公式已不再适用。
图3 喷管性能随喉栓位置的变化Fig.3 The variation of nozzle performance with the position of pintle
为了研究喉栓位置对喷管性能的影响,建立更为准确的推力模型,对喉栓在不同位置下的喷管内流场进行了数值模拟,分析喉栓在不同位置时的喷管内流场分布,研究喷管质量流率、推力等随喉栓位置的变化规律。
2.1 数值计算方法
采用FLUENT进行内流场数值模拟,RANS方程表达式为[15]:
1)连续方程
(9)
2)动量方程
(10)
3)能量方程
(11)
采用RNGk-ε模型,可以更好地处理高应变率及流线弯曲程度较大的流动,其湍动能k和耗散率ε方程分别为[15]:
(12)
(13)
式中:Gk是由平均速度梯度引起的湍动能;Gb是由浮力引起的湍动能;YM代表可压缩湍流中脉动膨胀的贡献;C1ε,C2ε,C3ε为经验常数,αk,αε分别为湍动能和耗散率对应的普朗特数;Sk和Sε为源项。
2.2 网格划分及无关性验证
几何模型由喷管内型面、喉栓头部型面等组成,根据结构的对称性,计算采用二分之一三维模型。对计算区域分块进行网格划分,在进气段等部位建立结构网格,对几何形状较为复杂的部位采用非结构网格,计算区域通过交界面进行参数传递。流体流经喉栓头部和喷管收敛段时流速会发生较大变化,因此需要对网格进行局部加密,局部网格示意图见图4。
图4 喷管内流场局部网格示意图Fig.4 Local grid diagram of flow field in nozzle
将边界条件设置为:
1)入口边界条件,采用压力入口边界条件,设置入口压强为9 000 000 Pa,总温为1 800 K;
2)出口边界条件,采用压力出口边界条件,设置出口压强为101 325 Pa,总温为300 K;
3)固体壁面条件,采用标准壁面函数,绝热无滑移条件。
4)假设流体介质为单一理想气体,不考虑多相流动,不考虑重力影响。
为排除各个工况因网格数量不同而造成的差异,需对已建立的网格进行无关性验证。将网格数分别为758 046和1 846 890的网格模型进行比较,验证喷管质量流率和推力大小,结果见表1。
表1 两种网格模型的计算结果Table 1 Calculation results of two grid models
两种网格模型计算得出的质量流率差异小于0.1%,推力值差异小于0.01%,可认为模型1已达到网格无关。
2.3 计算结果及讨论
通过计算得到了喷管内流场、喷管推力等随着喉栓伸入喉部的变化,图5给出了喉栓向外运动时喷管对称面上的压强和马赫数分布变化。
喉栓距离喷管喉部较远时对流场产生的扰动很小,如图5(a)所示,在喷管喉部及扩张段前半部分,压强逐渐减小,由于存在轻微的过膨胀,在扩张段接近出口处压强基本不再变化,马赫数沿着喷管轴线方向均匀增加,在喷管喉部达到Ma1,此时喷管每秒质量流率和推力最大。
喉栓向外运动到一定位置后,燃气在喉栓表面与喷管收敛段内壁之间加速膨胀,在喉栓头部附近形成等效喉部,由于非同轴进气,等效喉部先在进气一侧形成,如图5(b)所示,此时等效喉部面积与几何喉部面积相差较大。
喉栓继续伸入喉部,在喉栓与喷管内壁围成的环形通道处形成等效喉部,在喉栓头部形成了涡流区,且喉栓头部附近流场分布非对称性明显,如图5(c)所示。随着喷管开度减小,形成等效喉部的位置向上游移动且等效喉部面积减小,气流在等效喉部后的膨胀加剧,喉栓头部的涡流区缩小。喉栓处于全关位置附近时,喷管内流场分布愈加紊乱,喷管处于过膨胀状态,出现气流分离现象,如图5(d)所示。
图5 喉栓在不同位置时喷管对称面上的压强和马赫数分布图Fig.5 The pressure and Mach number distribution on the symmetry plane of the nozzle with different positions
喷管质量流率和推力随喉栓位置的变化如图6、图7所示,随着喉栓向外运动,喷管每秒质量流率和几何喉部面积近似线性减小。喷管推力在喉栓运动初期近似线性减小,而随着喉栓位移增大,对流场的扰动增加,且喉部面积减小,喷管扩张比增大,喷管过膨胀,推力随喉栓位置的变化率呈减小的趋势。
图6 每秒质量流率随喉栓位置的变化Fig.6 The change of mass flow rate per second with the position of pintle
理论推力随喉栓位置的变化率大于CFD仿真结果,在喉栓运动初期,理论推力大于CFD结果,且随着喉栓向外运动,两者差距逐渐减小,当喉栓运动至某一位置时,理论推力与CFD计算结果相吻合,而喉栓继续伸入喉部时,理论推力小于CFD计算结果且两者差距随喉栓伸入喉部而增大。
图7 理论推力和CFD仿真结果对比Fig.7 Comparison of theoretical thrust and CFD simulation result
经分析,等效喉部在接近喉栓头部的位置形成,燃气经过喉栓头部后流动方向发生较大变化,气体的惯性作用较为明显,存在流动损失,且由于非同轴的影响,喉栓上下型面形成等效喉部的位置不对称,使得几何喉部面积大于实际等效喉部面积,理论推力计算结果大于CFD仿真结果。随着喉栓伸入喉部,喉部出现狭缝式喉径,导致喷管的扩张比剧增,喷管处于过膨胀状态,气流分离加剧,不符合等熵流动假设,1.3节中的推力公式已不适用。因此,需要对推力模型进行优化。
3.1 气动等效喉部面积
等效喉部面积计算模型直接影响推力模型,是推力模型优化的基础。定义喷管的气动等效喉部面积计算公式为:
(14)
代入CFD仿真结果,可得喷管的气动等效喉部面积。改变入口边界条件和工质参数,出口边界条件不变,得到喷管在不同位置时的等效扩张比Ae/At,q如表2所示。
表2 不同喉栓位置和工况下的等效扩张比Table 2 The equivalent expansion ratio under different pintle positions and working conditions
表3 各工况条件Table 3 Working conditions
由结果可知,喉栓在初始位置时,当入口压强由9 MPa降为3 MPa,气动喉部面积略减小,喷管等效扩张比增加了0.046 2(0.36%,工质为燃气)和0.067 2(0.53%,工质为冷气);当入口压强不变,工质由燃气变为冷气时,气动喉部面积略有增加,喷管等效扩张比减小了0.093 8(0.74%,入口压强为9 MPa)和0.072 8(0.57%,入口压强3 MPa),各工况下的喷管等效扩张比标准差为0.056。
喉栓向外运动1个单位后,入口压强和工质发生变化时,对应的喷管等效扩张比变化率分别为0.28%,0.31%,0.47%,0.44%,各工况下的喷管等效扩张比标准差为0.058,较喉栓在初始位置时各工况下等效扩张比的标准差略有增加。
3.2 等效喉部面积计算模型
由于喷管内存在流动损失,喷管几何喉部面积大于实际等效喉部面积,给喷管推力计算带来误差,为优化喷管推力计算模型,采用喷管气动等效喉部代替推力模型中的几何喉部。
由CFD仿真结果可知,喉栓在同一位置而工况不同时,喷管等效扩张比的差异很小,因此,忽略工况对气动等效喉部面积的影响,近似认为气动等效喉部面积随喉栓位置的变化规律只与喉栓和喷管型面相关。
由图3可知,几何喉部面积与喉栓位置x近似成一次方函数关系,而根据CFD仿真结果,气动等效喉部面积与喉栓位置x近似成三次方函数关系,因此建立以下修正公式:
(15)
应用非线性最小二乘法,求得系数a1=0.92,a2=-0.75/(yo1-yo2),a3=-1.25(R1+R2)/(yo1-yo2)。
几何喉部面积、修正后的等效喉部面积及工况1下气动等效喉部面积对比如图8所示,气动等效喉部面积小于几何喉部面积,结合图5,随着喉栓向外运动,等效喉部形成的位置向上游移动,与流速方向接近垂直,与几何喉部面积差距减小。修正后的等效喉部面积计算公式能较好拟合气动喉部面积随喉栓位置的变化规律。
图8 等效喉部面积计算结果对比Fig.8 Comparison of equivalent throat area calculation results
3.3 推力模型优化
将修正后的等效喉部面积代入理论推力计算公式,计算得到工况1下喷管推力随喉栓位置变化如图9所示。
图9 初步修正后的推力计算结果Fig.9 The thrust calculation result after preliminary correction
由图可知,在喉栓运动初期,推力得到较好修正,与CFD仿真结果相吻合,当喉栓运动至某一位置后,推力计算结果与CFD仿真结果的偏差逐渐增大,此时喷管内流场不符合等熵流动假设,需要对推力计算模型再次进行修正,最终推力计算模型为:
(16)
(17)
(18)
(19)
式中:b1=46;b2=0.015;b3=-22.5/k;b4=-0.09。
上述推力模型对动推力和静推力两部分皆进行了修正,利用优化后的推力计算模型,计算得到工况1下喷管推力随喉栓位置变化如图10所示。
图10 优化后的推力模型计算结果Fig.10 The optimized thrust model calculation results
表4 不同喉栓位置下的推力计算结果对比Table 4 Comparison of thrust calculation results under different pintle positions
将理论推力、优化模型计算结果与CFD仿真结果进行对比,结果如表4所示,均方根误差由0.063减小至0.002,优化后推力模型计算结果与CFD仿真结果吻合较好。
4.1 试验系统
试验系统由对称螺旋进气管、压强传感器、电磁阀、静态标定组件、喷管、板簧、螺杆和力传感器等组成,如图11所示。进气管路在结构上采用了螺旋方式进气且在两侧对称分布,高压气瓶由减压阀、高压软管与螺旋管路连接,通过减压阀可初步调整气体流出的压强,进气端通过连接板支架固定在底板上,由于气流流经对称螺旋管路后同时流出,能对气流起缓冲和稳定的作用,最后经由三通管将稳定的气流送入主管道。
图11 试验系统示意图Fig.11 Diagram of test system
高压气体通过两端对称螺旋管路汇入主管路,经过压强传感器,待气流稳定后由压强传感器记录腔体压强值,再由电磁阀控制主管路气流的开闭,保障主管路气流的稳定性,同时可以检验管路的密闭性是否良好。喉栓作动时,喷管朝上喷气,通过板簧和螺杆将力传递给力传感器。对整个试验系统进行标定后,可进行喷管冷流试验,测得不同腔体压强下的喷管推力。
4.2 试验结果及分析
测量喉栓在不同位置时喷管的稳定输出推力,测试结果见表5。
腔体压强分别在9 MPa和3 MPa左右时,推力模型计算结果和试验结果的最大绝对百分比误差分别为10.72%和 11.07%,均发生在x/xmax=0.9处,即喷管即将关闭时。经分析,推力模型计算结果与试验结果存在偏差是由于推力模型的建立过程中存在简化,比如忽略了工况对气动喉部面积的影响、喷管外流场的影响等,且随着喉栓伸入喉部,对流场的扰动增加,喷管内部流动复杂,给推力计算带来了困难。此外,试验系统存在进气扰动、管路和板簧等带来的推力传递损失等,也给推力测量带来了误差。
腔体压强分别在9 MPa和3 MPa左右时,不同喉栓位置的喷管推力计算结果与试验结果的平均绝对百分比误差分别为3.10%和4.99%,修正后的推力模型计算结果与试验结果基本吻合。
表5 试验结果与推力模型计算结果Table 5 Test results and thrust model calculation results
1)建立了三维仿真模型,可准确地模拟喉栓对喷管内流场的影响,随着喉栓深入喉部,喷管推力和出口质量流率减小,喉栓对流场的扰动加剧,喷管扩张段出现气流分离现象。
2)建立了气动等效喉部面积的计算方法,修正了喉部面积计算模型,气动等效喉部面积与入口边界条件、工质参数、喉栓和喷管型面相关,入口边界条件与工质参数对等效喉部面积的影响在一定范围内可以忽略。
3)优化了推力模型,该模型可预示喷管推力随喉栓位置的变化规律,喉栓向外运动时,喷管推力先近似线性减小,后喉栓伸入喉部,喷管推力不满足线性变化,随着喉栓向外运动逐渐降为0。
4)进行了喷管冷气推力测试,获得了不同喉栓位置和腔体压强下的喷管推力,试验结果与推力模型计算结果基本吻合。
猜你喜欢 流率喉部气动 中寰气动执行机构天然气与石油(2022年4期)2022-09-21轴排凝汽器喉部设计今日自动化(2022年11期)2022-02-22CDN 缓存空间调整优化方法广东通信技术(2021年9期)2021-10-12基于NACA0030的波纹状翼型气动特性探索北京航空航天大学学报(2021年6期)2021-07-20带交换和跳跃的一维双向自驱动系统的仿真研究安庆师范大学学报(自然科学版)(2021年1期)2021-03-14电子喉镜联合窄带成像技术对喉部早期恶性病变的诊断价值研讨健康大视野(2020年19期)2020-10-26多束光在皮肤组织中传输的蒙特卡罗模拟泰山学院学报(2018年6期)2018-12-18凝汽器喉部流场数值模拟发电设备(2018年6期)2018-11-29基于反馈线性化的RLV气动控制一体化设计测控技术(2018年9期)2018-11-25聚合物气辅共挤中熔体体积流率对挤出胀大的影响上饶师范学院学报(2016年3期)2016-08-02