本文介绍了一种用于模拟固体火箭发动机内部弹道和轴向声学波动的非稳态准一维流动求解器。由于控制方程的非守恒形式和横截面积的非平滑轴向变化,准一维控制方程的高阶数值解容易产生数值振荡。在中心差分中添加人工耗散被发现不足以抑制这种振荡,因此使用了一个简单的低耗散激波捕捉方案(名为SLAU2)来代替。该方案固有的数值耗散有助于正确捕捉陡峭的声波(在触发的不稳定性开始时形成的),而不会对声波产生不良的数值阻尼。使用新的求解器,提出了一个计算特征频率、相应的模态形状和阻尼率的程序,并对一个具有圆柱形装药几何形状的发动机进行了验证。准一维公式可以适应装药边界的轴向速度的滑移或非滑移边界条件,并相应地预测阻尼率。用这两种边界条件得到的衰减率的差异被证明是流动转向的贡献。

命名法

1.引言

在大多数固体火箭发动机中,正常工作的流动都是一维的、准一维模型可以给出可接受的弹道预测。尽管一些多维效应,如与潜入喷管、熔渣沉积和涡流脱落有关的效应,确实需要多维建模,但其他效应,如边界层和喷管烧蚀损失,可以在一维模型中用经验公式来处理[1]。根据固体火箭发动机的许多实验,Blomshield[2]指出,在固体发动机中,高能径向模式非常罕见,而切向模式在使用双基推进剂时可能会出现,但在使用金属化推进剂时几乎没有出现过。这并不奇怪,因为固体火箭发动机的燃烧室有很大的长径比。因此,一维程序对于固体火箭发动机设计的快速线性稳定性评估非常有用。

大多数早期的弹道学程序是基于打靶法[1,3](在空间上有累积误差,碰到流道面积突变的情况时,误差会非常大)。它涉及到迭代调整头端压力,并在空间上整合稳定的准一维控制方程,以便喷管入口条件对应于给定喉部面积的拥塞流动条件。这种方法不适合模拟瞬态过程,如点火和波传播,不仅因为它使用了伪稳态和短喷管近似,而且还因为它的一阶空间精度。燃烧不稳定性研究中的声动力学是使用通过这种射击方法得到的平均流场上的扰动的线性化方程来处理的。对线性化方程的特征分析导致通过扰动的模态分解进一步简化了问题。只有几个主要的大尺度模式可以用一组常微分方程来跟踪,而不是解决波运动的偏微分方程。完整的特征分析和确定所有的特征模式在计算上非常昂贵。Culick[4]介绍的使用格林函数的线性波方程的空间平均技术提供了一种计算主要特征模的有效方法,并且一直是过去稳定性程序的支柱。这种相对低代价的线性稳定性分析对于快速评估和排除固有的不稳定设计非常有用。

触发的不稳定性,主要是由非线性机制引起的,[2,5-7],不能用线性稳定性分析来预测。Flandro等人[9]指出,使用固体火箭稳定性预测(SSP)程序[8]预测为线性稳定的发动机可以被脉冲到不稳定[9]。SSP经过快速修改[10],包括了一些非线性因素,以改善稳定性预测。

确定大量的模态和相应的频率来处理压缩波陡峭化为激波波的问题在计算上非常昂贵[4]。使用高阶有限差分法或有限体积法代替模态表示法计算扰动总是可能的,但必须处理与非光滑面积变化和气体动态不连续带来的数值问题。如果必须处理这些问题,可以计算完整的解决方案本身,而不是在一个平均值上的扰动[11, 12]。

用于捕捉激波而不产生欧拉和NS方程解的数值振荡的人工粘性方法[13]也可用于解决准一维方程[14,15]。然而,由于准一维方程的非守恒形式,它们的效果并不理想。如果面积变化不平滑,就会遇到问题,如本研究中所示。这就是为什么Baum等人[12,16]以及Levine和Baum[17]使用三种不同方案的组合来模拟非线性不稳定性的一个原因。Loncaric等人[18]、Montesano等人[19]、Baczynski和Greatrix[20]、Greatrix[21]以及Montesano等人[22]最近的工作是基于非确定性的随机选择方法(RCM),主要是为捕捉双曲守恒定律的不连续数值振荡解而开发的。

基于黎曼求解器的确定性方法[23-28]使用本质上非振荡(ENO)[24,25,27]或不连续Galerkin(DG)[26]方案扩展到高阶,也被开发并在公开文献中报告。Pekkan和Ucer[23]开发的Riemann求解器被用于稳态和脉动模拟,但只用于圆柱形和端燃的装药。前者的轴向面积变化沿轴线平滑,后者的轴向面积变化为零。端燃装药对应于圆柱形腔体的流入边界条件。d,Agostino和Andrenucci[24]以及Ferretti[27]的ENO方案被证明可以处理急剧的横截面积变化而不产生虚假的振荡。后者的工作[27]侧重于在一个准一维框架内对涡度-声学相互作用进行建模。Cavallini等人[28]也专注于这方面的研究,使用单调的上流中心守恒定律方案(MUSCL)方法,而不是ENO方案来实现二阶空间精度。Willcox等人[25]的三阶通量分割方案被用来测试准稳态假设(通常在进行零或一维弹道预测时进行),并研究侵蚀性燃烧的影响。据报道,用这种方法在面积突变附近的位置进行了非物理性的预测。除了Shimada等人的工作[26],上述研究都没有试图量化影响轴向声学的阻尼机制。在这项研究中,准一维的结果与线性理论和轴对称模拟的预测有很大不同。

毫无疑问,在过去的研究中使用的Riemman求解器[23-28]可以清晰地捕捉到激波,但保留这一特征以及处理非光滑区域变化的能力,尽管有可能(特别是使用ENO方案),但还没有得到证明。Riemann求解器捕捉激波的能力伴随着数值耗散,这对于正确解决波动力学问题是不可取的。各种求解器/方案中动能的数值耗散已被广泛研究,但声学的数值耗散得到的关注却少得多。包括耗散在内的数值学所造成的误差一般在小尺度上比较大。除非小尺度的动力学影响到大尺度的动力学,否则各种求解器预测的大尺度动力学应该是非常相似的。Shima和Kitamura[29]提出了相反的证据,他们表明,即使是大尺度声学的数值耗散水平也会因不同的黎曼求解器而有所不同。

本文记录了开发一个用于弹道和稳定性预测的非稳态准一维模型的工作,要求如下。首先,该方法应该是确定性的(与RCM不同),并且基于一个单一的方案而不是不同方案的组合。为了简单起见,与ENO方案相比,MUSCL方法更可取。由于非光滑面积变化和激波引起的数值振荡应该是最小的。与物理阻尼相比,声学的数值阻尼应该是可以忽略的。为此,测试了两种不同的方案。第一个是带有人工耗散的二阶MacCormack方案,第二个是SLAU2激波捕捉方案。对于稳态弹道预测,与射击方法的预测进行了比较。通过模拟压缩波崩溃成小激波,以及预测与圆柱形发动机中不同物理过程相关的阻尼率,对声学的零燃烧响应,测试了非稳态能力。

显然,仅用准一维方程无法捕捉到与分段发动机中涡流脱落和顶层涡流脱落[30]相关的旋转动力学。但是,与轴向压力振荡产生的不稳定涡度场相关的旋转项,甚至在理论分析中出现的简单几何形状(如圆柱形装药)[31],使人们对准一维模型预测阻尼率的使用产生了疑问,因此,其用于稳定性预测。在稳定性理论中,流动转弯损失似乎是解释实验结果所必需的[32],但据作者所知,使用准一维模型对其进行预测,从未被验证过。

2. 控制方程

用于固体火箭发动机建模的准一维方程中包括的项因研究重点不同而不同[18,23-26,33]。这里使用的方程只包括说明空间和时间上横截面积变化的项,以及在每个轴向位置上对当地装药边界的法线质量加入:

与其他项相比,由于推进剂的火焰速度而增加的能量(当装药边界不垂直于轴向流动方向时,还有动量)被忽略了。

前面的集合对应于方程的最守恒的形式,在文献中报道的大多数基于ENO的方案[24,25,27]中使用。公式(2)右侧的第一个项使得完全的守恒形式成为不可能,并给高阶空间离散化带来了问题。然而,这一组被称为守恒形式。其他方案(RCM[18]和DG[26])是基于另一种形式,其中面积变化的影响包括使用一维欧拉方程的源项。这套方案被称为非守恒形式。

公式(2)中的下划线项意味着边界处的传入质量以当地的轴向速度进入控制体。这相当于一个滑移的边界条件。没有它的计算对应于无滑移边界条件。

3. 数值计算方法

二阶显式Runge-Kutta方案被用于时间积分。空间离散采用了有限体积方案。单元面上的通量使用MacCormack方案[34]或SLAU2激波捕捉方案[35]进行计算。

MacCormack格式是用前面列出的方程的守恒形式实现的。一个面上的通量是这样定义的,在奇数步中,在预测步中使用以单元为中心向右向右方向的值(向前差分),在校正步中使用以单元为中心向左方向的值(向后差分)。在偶数的时间步骤中遵循相反的偏置顺序,以防止整体方向性的偏置。这个方案的细节和性能已经被很好地记录下来了[34]。人工耗散项被加入以克服声障并抑制激波附近的数值振荡。人工扩散性D2与压力开关量(使其仅在激波周围不为零)和系统中最快波的速度成正比。

虽然该方案可以捕捉到激波,但其解决方案表现出激波后的数值振荡。虽然在准一维模型的背景下从未报道过,但使用四阶超扩散项可以抑制它们[13]。这里考虑了带有和不带有这些超扩散项的人工粘性方法。

由Shima和Kitamura[29,36]开发的简单的低耗散AUSM方案(SLAU)在本工作中被考虑作为MacCormack方案的替代。它的设计是为了在通量分裂框架内保持动能的数值阻尼最小,但其对声学的数值阻尼也被证明比Roe方案[29]小得多,而许多正在使用的Riemann求解器都是基于该方案的。这里考虑的是SLAU2,它是SLAU的一个改进版本,具有更强的激波预测能力[35]。

公式(A1-A3)左边的项使用SLAU2方案的有限体积版本进行离散,而右边的项则作为源项处理。在最初的SLAU2中,ρA和pA代替了ρ和p。

面的左边和右边的属性分别用下标L和R来表示。用于计算通量的面的压力是用以下公式得到的。

其中

M代表马赫数,用面的法线速度和界面声速c1/2计算,K代表比动能。Kitamura和Shima[35]指出,SLAU2对临界声速的规格不是很敏感,因此为了简单起见,采用了两边数值的几何平均值。

跨越面的质量通量是用以下公式计算的。

△(pA)表示压力和横截面积在单元表面的跳跃。来自上游一侧的速度矢量和总比焓以及前述方程中的质量通量(由\(\dot m\)决定)被用来计算对流通量。

在一阶版本中,一个给定面的左边和右边的值对应于它两侧的单元中心的值。为了实现空间的二阶精度,使用MUSCL方法从单元中心的高阶插值来计算左和右的值。

由于基于方程守恒形式的MacCormack和SLAU2方案都不能满足必要的要求(如后面几节所示),所以还必须考虑基于控制方程非守恒形式的SLAU2方案[18]。在这些方程中,左侧与一维欧拉方程完全相同。左边的项用文献[29,35,36]中提出的有限体积SLAU2方案进行离散化;而右边的项,包括面积变化的影响,被当作源项处理。

控制方程和解决它们的SLAU2方案列在附录中。

4. 结果

本研究报告中的模拟是针对不变的装药几何形状,但所有流动量的非稳态动力学是完全时间精确的。在固体火箭发动机中,面积随时间的变化通常是非常缓慢的,这就是为什么许多弹道学程序中的准稳态假设被证明是合理的。

表1列出了用于模拟的参数的默认值。这些值是使用AP/HTPB复合推进剂的小型火箭发动机的典型值。除非另有说明,所有的模拟都使用这些值。由于燃烧反应没有被模拟,燃烧率保持不变。默认的空间分辨率是5毫米,时间步长是通过假设库伦特-弗里德里希-卢伊(CFL)数为0.4(二阶方案的稳定性极限为0.5)来确定。

为了测试数值方法的稳态性能,考虑了四种燃烧器的几何形状,如图1所示。在头端和喷管内没有推进剂燃烧。后锥体和前锥体装药几何形状是轴对称的finocyl和反向finocyl装药配置的代表,通常用于实现时间上的恒定推力。除圆柱形外,所有装药的几何形状都具有相同的燃烧面积,因此应导致相同的燃烧室压力(当使用零维模型计算时)。在研究中使用的锥角与实际火箭发动机中遇到的数值相比,故意选择得相当高,以便测试程序在处理快速和突扩的横截面积变化方面的能力。它们通常被保持在较低的水平,以避免涡流脱落和相关问题。

Fig. 1 模拟中使用的内部推进剂装药的几何形状。从上到下,这些几何形状被称为a)圆柱形,b)后锥体,c)前锥体,和d)中锥体。

4.1 稳态仿真

中锥体装药的几何形状具有所有的复杂性,如横截面积的急剧变化和源项的不连续变化(即,当流量进入喷管时,质量增加突然停止),这可能给数值建模带来问题。它被用来比较各种技术,而不是模拟所有的几何形状。使用打靶方法预测的解决方案被用作参考。仅在燃烧室(图1中0.0 m < x < 2.0 m)进行比较,而不是在喷管(图1中x > 0.5 m)进行比较,有两个原因。首先,如果包括喷管部分,由于那里的压力变化陡峭,各种方法的预测之间的差异变得不那么明显。第二,射击方法没有延伸到喷管(其中使用了短喷管的近似值)。

图2显示了使用各种方案预测的轴向压力曲线。为了说明各种方案预测的差异,图3显示了在其中一个区域突然变化的位置周围的预测的放大图。MacCormack方案似乎工作得非常好,即使在亚声速区域没有任何人工粘性。人工粘性使其能够处理喷管中的超音速流动。激波可以用人工粘性来捕捉,但激波后的数值振荡是不可取的。四阶耗散项被加入以抑制激波后振荡,当只使用二阶项时就会出现这种振荡。然而,当这些额外的四阶耗散项被添加时,在面积变化不平滑的地方会产生数值振荡。

Fig. 2 使用各种方案预测的中锥体装药几何形状的压力曲线比较。

Fig. 3使用各种方案预测的中锥体装药几何形状的压力曲线的放大图。图例与图2相同。

Fig. 4 对前锥体装药几何形状预测的压力曲线的比较。

令人惊讶的是,基于方程守恒形式的SLAU2方案非常相似。在这些位置附近产生的不是震荡,而是大的尖峰。而基于非守恒方程形式的SLAU2方案则不存在这个问题。总的来说,带有二阶人工耗散的MacCormack方案和基于非守恒方程形式的SLAU2方案预测的结果相似,与射击法的预测结果相比,似乎也可以接受。尽管轮廓非常相似,但与射击法的预测相比,由于大约0.6%的质量误差,预测结果向上偏移了。这个误差可能是MacCormack中的人工粘性造成的,而在SLAU2的预测中,它是由于使用了控制方程的非守恒形式造成的。在本文的其余部分,SLAU2指的是该方案的这个版本。

图4和图5分别显示了SLAU2对前锥体和后锥体几何形状的预测和射击法的相应预测。在后一种情况下,SLAU2的预测在圆锥体开始的那一点附近有明显的小幅振荡。两点振荡的振幅是平均压力的0.02%。这些数值振荡在给定的网格上具有最高的可能波数,不太可能影响具有最低波数的大规模声学模式。

Fig. 5 对后锥体颗粒几何形状预测的压力曲线的比较.

Fig. 6 只在前圆柱形部分燃烧的中锥体装药几何形状预测的压力曲线比较

(0.0 m < x < 0.5 m).

为了进一步测试,我们考虑了中锥体的几何形状,其中只有靠近头端的圆柱形部分(图1中0.0米<x<0.5米)有燃烧的推进剂。在这种情况下,由于面积也以无差别的方式变化,质量增加量突然下降到零。图6中比较了使用SLAU2和带有二阶人工耗散的MacCormack方案对这种情况的预测结果。在这两种情况下,在燃烧停止的角落里,数值振荡的振幅小于平均压力的0.02%。

对于圆柱形装药发动机,准一维预测与使用FLUENT(一种商业计算流体力学软件[37])进行的轴对称模拟预测的参考方案进行了比较,而不是使用弹道学程序。轴对称预测是基于一个结构化的网格,该网格在轴向方向上是均匀的,在壁面附近有足够的集群来解决层流边界层。如图7所示,整体压力水平高出0.3%,可能是由于喉部边界层的额外拥塞效应。这个轴对称案例的程序设置也被用于非稳态验证。

4.2 非稳态验证

选择两类问题进行非稳态验证。鉴于使用许多基于准稳态假设的预测在点火瞬时后被证明是相当合理的,在没有燃烧不稳定性的情况下,平均流动的不稳定性可以被忽略。然而,对于燃烧不稳定性的研究,声学预测的时间精度是非常重要的。建立程序模型不稳定性能力的第一步是测试它对声学阻尼项的预测。在此基础上,验证了在不发生数值振荡的情况下,在触发不稳定性开始时捕获尖锐前沿波的能力。

Fig. 7 对圆柱形装药的几何形状,准1-D和轴对称模拟预测的压力曲线的比较。

4.2.1 声学阻尼

除了与凝相粒子相关的声学阻尼项,所有较低轴向模态的声学阻尼项都可以在一维框架内以合理的精度计算出来,前提是:1)数值阻尼可以忽略不计;2)轴向和横向模态之间没有能量交换(即没有混合模态)。

Shimada等人[26]研究了准一维方法对轴向声学建模的适用性,观察到预测的阻尼率与使用轴对称模拟得到的阻尼率不同。数值阻尼被认为是一种可能的解释。也有人推测能量转移到径向模式,但推断其可能性不大。

尽管有数值阻尼,一维方法足以捕捉到与喷管有关的损失。这些损失主要取决于喷管入口马赫数。在准一维方法中,这个数量并不含糊,而在轴对称模拟的情况下,它不清楚应该如何计算。 非均匀性也是必要的,但不清楚。这种不均匀性的影响可能会导致一些错误,这一点还没有得到系统的分析。

与喷管损失相比,流量转弯损失一直是争论不休的话题。它在Culick的一维分析中得到了体现[38],但在多维公式中却没有体现[8]。Flandro澄清说[39],这是由于没有使用装药边界的无滑移条件。鉴于Shimada等人[26]在多维模拟中使用了允许滑移条件的欧拉求解器,而不是NS求解器,流动转弯的贡献应该是没有的。然而,多维模拟的总体损失比准一维模拟要高。可能的解释是数值阻尼或能量转移到径向模式,正如他们所推测的那样[26]。

正如Flandro[39]所指出的,流动可以在任何指定的角度被注入。但是,在准一维计算中,由于滑移和非滑移条件而导致的损失有多准确,仍然没有得到证实。有理由预计,一些多维效应,如由涡度场引起的效应,不能在一维求解器中准确捕捉[39]。在滑移边界下,非稳态波动速度的相位不随径向坐标变化。然而,在粘性计算的情况下,情况并非如此。轴向波动压力梯度被中心线上的非稳态项所平衡,而它被壁面上的非稳态粘性项所平衡。在装药边界本身,没有对轴向速度使用滑移条件,但是靠近壁面的轴向波动速度与中心线上的相同数量不在一个相位上。

Fig. 8 在封闭的管子里进行声学模拟,以估计数值阻尼

这个相位差已经在理论上[31,40]、实验上[41,42]以及本研究的轴对称计算中得到证实。Majdalani[43]最近的分析解决方案表明,壁面的剪切力与压力梯度同相,就像前面的方程式一样,但在更高的频率下,滞后接近π/4。由于这里的重点是较低的轴向模式,这种相位滞后无法得到证实。这种沿径向坐标的声速相位变化的影响还不清楚。这部分研究试图解决数值阻尼、能量转移到径向模式和边界条件的影响,所有这些都会影响一维预测的准确性。

一个简单的测试有助于估计与数值有关的声学阻尼。第一谐波被叠加在长度为2米、直径为123毫米的圆柱形室内的停滞状态(类似于图1中圆柱形情况下的燃烧区),并使用Runge-Kutta时间积分方案和SLAU2通量分裂空间离散方案预测几个周期的时间演变。该领域使用400个点进行离散,这意味着空间分辨率与用于火箭发动机几何形状的空间分辨率相同。图8显示了其中一个端点的压力。图中还显示了与预测值相适应的衰减正弦波。数值阻尼率常数原来是0.008/s,与大约87Hz的特征频率相比,可以忽略不计。与典型的火箭发动机中的喷管和流转损失造成的物理速率常数相比,它也是可以忽略的。

对于图1所示的圆柱形发动机几何形状,在达到稳定状态后,在头端加入随机速度波动,在一段时间内收集几个位置的压力数据。这模拟了扬声器产生的白噪声。在压力数据的傅里叶变换中,对应于特征模式的频率很突出。如果选择一个点来计算频谱,那么在该位置有节点的特征模式将被遗漏。因此,为了避免这种可能性,要计算多个点的频谱。所有峰值位置的超集提供了特征频率的集合。尽管由于喷管的存在而不是封闭端,可以预期会有小的偏差,但圆柱形室内的声学大致对应于管道谐波,节点上的误差很容易被避免。在x=0.065米处(距离头端第13个网格点)使用压力计算出的圆柱形几何的频谱如图9所示。频谱中存在大量由随机强迫产生的噪声,因此特性频率的精确值很难确定。在该图确定的最低特征频率下,在头端施加一个正弦波强制,并收集压力数据。在这个频率上的强迫只激活了第一个模式。非稳态模拟几乎可以在任何时候停止。压力的空间变化对应于一个稳态轮廓,在一定的振幅下增加了第一模式。以此为初始条件,停止头端强迫,并注意压力与时间的关系。压力信号正好是一个衰减的(单色)正弦波,正如预期的那样。如果衰减过程中的频率与用于强制的频率相差很大,则要重复单色强制和衰减步骤。如果频率规格与特征频率相差甚远,或者初始模式形状不正确,就会出现压力干扰的多模态衰减,正如Javed和Chakraborty[44]所观察到的。

Fig. 9 圆柱形发动机随机脉动时的压力频谱。

这里使用的方法在两个方面有帮助。对于圆柱形装药的几何形状,基波的频率大致计算为c/2L,其中L是腔体的长度。目前还不清楚它是否应该包括喷管的收敛段部分的长度。对于模拟的几何形状,喷管大约是12厘米,这大约是燃烧室长度的6%。所以,在估计频率方面有6%的不确定性。对于模拟的圆柱形装药,频率在258赫兹(包括喷管的收敛段部分的计算)和273赫兹(通过排除喷管部分的计算)之间。通过前面详述的测试,这个几何形状的实际频率可以被准确地确定。

同样的程序也适用于轴对称程序,以产生用于验证准一维预测的参考值。与稳态模拟的情况一样,轴向速度采用无滑移边界条件。装药边界的径向速度在空间和时间上是固定的,以模拟对声学无反应的燃烧。使用一维程序计算了四种不同内径(100、123、160和200毫米)的圆柱形发动机的第一纵向模式的衰减率。对于其中两个配置,也使用轴对称模拟计算了衰减率。

图10显示了一维(1-D)和轴对称程序在直径为123毫米的发动机中第一模式衰减期间预测的压力曲线。衰减的正弦波已被拟合以提取衰减率。在这种情况下,频率预测相差约0.12%,而衰减率预测相差约0.67%。图11显示了直径为200毫米的情况下的类似比较图。频率和阻尼率的差异分别为0.8和0.4%。这验证了没有能量损失给径向模式。只要一个模态的波长比径向尺寸大得多,这就可能是真的。对于波长较短的高模,可能确实存在对径向模的损失,这取决于喷管收敛段部分的形状。

Fig. 10 直径为0.1234米的圆柱形装药几何形状的第一模式压力脉冲的衰减。

图11直径为0.2 m的圆柱形装药的第一模态压力脉冲衰减。

基于Culick和Yang的线性分析[45]或Vuillot和Casalis[46]提供的表达式,对于圆柱形的几何形状,可以分析计算由喷管引起的衰变率贡献。对流和辐射分量的计算方法如下:

流量转折的贡献估计为Vinj/R。因此,总衰减率变成了(γ+2)Vinj/R。在准一维程序中保持所有其他参数不变的情况下,通过改变燃烧室圆柱形部分的内径进行的预测与图12中的分析表达式进行了比较。

如果假定一维程序的预测是正确的(因为它们与多维程序的预测相比很好),分析表达式导致对衰减率的预测不足。误差百分比随着直径的增加而增加,在这里考虑的最高腔体直径(对应于长度与直径之比为10)时,误差约为6%。随着发动机直径的增加(喉部直径固定),喷管入口马赫数减少,衰减率的线性近似值应该更加有效。如果能用一维程序分别计算各个贡献,就能更好地理解这种差异。

Fig. 12 圆柱形发动机的总阻尼随端口半径的变化。

Fig. 13 具有滑移边界条件的第一模式压力脉冲在直径为0.2米的圆柱形装药几何形状的发动机中的衰减。

为了从喷管阻尼的贡献中分离出流动转向损失的贡献,在稍作修改后重复了上述练习。修改了动量和能量方程,假设由于推进剂燃烧而增加的质量有一个轴向速度分量。在原来的方程中,增加的质量是在流动的法线方向,因此对动量方程没有贡献。质量增加的轴向分量可以指定为局部瞬时流速(这与多维模拟中使用滑移边界条件相对应)或局部波动(局部瞬时流速减去局部时间平均流速)。使用这两种不同的规格,得到的平均压力曲线略有不同。当添加的质量与系统的波动同相时,就不会有流转损失。对于一个具有200毫米直径圆柱形部分的发动机,图13显示了对第一种模式的预测衰减的曲线拟合;而衰减率原来是19.97/s。通过从28.18/s减去,在有流量转弯损失的情况下得到的流量转弯损失贡献值为8.21/s,估计这个几何形状。这与理论值8.26/s(Vinj/R)相差约0.6%。通过简单地考虑无滑移条件,流动转向损失被准一维程序所捕获,这意味着其大小与室中流动转向的速度(即壁附近的平均涡度场)或与声速相关的非稳态涡度动力学无关。

由于流转损失被预测为非常接近其理论值(Vinj/R),所以总衰减率的预测值和理论值之间的差异是由于与喷管有关的损失。减少喷管收敛段部分的长度对总衰减率的影响可以忽略不计。这再次证实了在轴对称或准一维模拟中没有径向模式的损失。方程(18)是基于使用短喷管近似法确定的导纳值。研究表明,这种近似的结果是预测不足。这可能是理论表达导致整体阻尼值较低的原因。

鉴于各种阻尼贡献的轴对称预测与圆柱形装药几何形状的理论表达式一致[45],没有证据表明与声学边界层中的非稳态涡度场有关的额外旋转机制可以抵消Flandro和Majdalani[31]提出的流动转向损失。

4.2.2 非线性效应:波的陡峭化

随着压力扰动振幅的增加,非线性效应变得明显。使用离散模式叠加的线性分析不再有效。压缩波变陡,而膨胀波则相反。这导致了在实验中观察到的锯齿状的压力空间分布的产生[9]。这种行为被认为会引起触发的不稳定性和直流偏移[7,9,12]。

在达到稳定状态后,通过在发动机头端引入高频(3600赫兹)和高振幅的正弦波强迫,来测试捕获尖锋波的能力。高振幅的压缩波会自我深化为几乎类似激波的结构。为了对这些结构提供更好的分辨率,网格间距从5毫米减少到1毫米。

空间压力分布图见图14。带有二阶人工耗散的守恒的MacCormack方案在尖锐锋面后面留下了非物理的振荡。Baum等人[16]也用这种方案观察到这种振荡。这些不是两点振荡,所以谐波扩散算子无法摆脱它们。这些振荡通常是通过使用一个额外的四阶人工耗散来控制的,这个人工耗散模拟了一个双谐波算子[13]。这种做法对有静止激波的稳态问题效果很好,但在激波移动的情况下可能效果不佳。另一点需要注意的是,它在一维和多维模拟中的效果是已知的;对于准一维流动(具有不可区分的面积变化),它的效果如何,至今还没有报道。除了紧跟在激波后面的小的非物理性结点,振荡大多被抑制。Baum等人[12,16]并没有试图控制MacCormack方案产生的振荡。相反,他们转向了一个完全不同的方案。这里也是这样做的,因为四阶耗散在处理非平滑面积变化时导致了问题(如图2所示)。SLAU2,作为一个激波捕捉方案,在不产生数值振荡的情况下清晰地捕捉到了小激波,并且比等人构建的方案组合要简单得多[12,16]。

5. 结论

开发了一个准一维程序来模拟固体火箭发动机内部的流动。在横截面积的轴向变化为非光滑的位置附近,其解决方案中的数值振荡的振幅足够小,在弹道或线性稳定性预测中没有任何影响。声学的数值阻尼足够低,以至于可以准确计算出特征模的线性阻尼率。在触发的不稳定性开始时形成的类似激波的结构被捕获,没有任何数值振荡。

装药边界的轴向速度的边界条件的类型决定了由于流动转向而产生的声能损失。滑移和非滑移的边界条件都可以在准一维框架中进行模拟。通过比较使用两种类型的边界条件所做的衰减率预测,可以准确地估计出流动转向损失的贡献。尽管在多维模拟和实验中,已知声速波动的相位会沿径向变化。

预计该程序对弹道以及稳定性预测都很有用。目前正在研究装药的几何形状对线性稳定性的影响,并将在未来报告。目前正在努力将其与装药演变程序结合起来,以跟踪几何形状随时间的变化。

附录:非守恒型方程的SLAU2格式

准一维流动的控制方程可以改写为一维欧拉方程,其中源项包括横截面面积变化的影响,以及横流质量和壁面上的加热:

参考文献

[1] Hermsen, R. W., Lamberty, J. T., and McCormick, R. E., “A Computer Program for the Prediction of Solid Propellant Rocket Motor Performance (SPP),“ Air Force Rocket Propulsion Lab. Rept. TR-84- 036, Edwards AFB, CA, 1984.

[2] Blomshield, F., “Lessons Learned in Solid Rocket Combustion Instability,” AIAA Paper 2007-5803, July 2007.

[3] Hodge, B. K., and Koenig, K., “SOLROC, A Solid Rocket Motor Internal Ballistics Software Element,” AIAA Paper 1994-3115, July 1994.

[4] Culick, F. E. C., “Unsteady Motions in Combustion Chambers for Propulsion Systems,” NATO-AGARD Rept. AG-AVT-039, 2006.

[5] Hughes, P., and Saber, A., “Nonlinear Combustion Instability in a Solid Propellant Two-Dimensional Window Motor,” AIAA Paper 1978­1008, July 1978.

[6] Levine, J. N., and Baum, J. D., “A Numerical Study of Nonlinear Instability Phenomena in Solid Rocketmotors,” AIAA Journal, Vol. 21, No. 4, 1983, pp. 557-564.

[7] Malhotra, S., Flandro, G., Malhotra, S., and Flandro, G., “On the Origin of the DC Shift,” Joint Propulsion Conferences, AIAA Paper 1998­3249, 1997.

[8] Nickerson, G. R., Culick, F. E. C., and Dang, L. D., “Standard Stability Prediction Program for Solid Rocket Motors,” Air Force Rocket Propulsion Lab. Rept. TR-83-017, Edwards AFB, CA, 1983.

[9] Flandro, G. A., Fischbach, S. R., and Majdalani, J., “Nonlinear Rocket Motor Stability Prediction: Limit Amplitude, Triggering, and Mean Pressure Shift,” Physics ofFluids, Vol. 19, No. 9, 2007, Paper 094101.

[10] French, J. C., “Nonlinear Combustion Stability Prediction of SRMs Using SPP/SPP,” AIAA Paper 2003-4668, July 2003.

[11] Baum, J. D., Levine, J. N., and Lovine, R. L., “Pulse Triggered Instability in Solid Rocket Motors,” AIAA Journal, Vol. 22, No. 10, 1984, pp. 1413-1419.

[12] Baum, J. D., Levine, J. N., and Lovine, R. L., “Pulse Triggered Instability in Rocket Motors: A Comparison Between Predictions and Experiments,” Journal of Propulsion and Power, Vol. 4, No. 4, 1988, pp. 308-316.

[13] Jameson, A., Schmidt, W., and Turkel, E., “Numerical Simulation of Euler Equations by Finite Volume Methods Using Runge-Kutta Time Stepping Schemes,” AIAA Paper 1981-1250, July 1981.

[14] Coats, D. E., Dang, L., and Nickerson, G. R., “Internal Ballistics Calculations for Nozzleless Solid Propellant Rocket Motors,” AIAA Paper 1982-1199,July 1982.

[15] Wornom, S. F., and Hafez, M. M., “Calculation of Quasi-One- Dimensional Flows with Shocks,” Computers and Fluids, Vol. 14, No. 2, 1986, pp. 131-140.

[16] Baum, J. D., Levine, J. N., and Lovine, R. L., “Pulsed Instability in Rocket Motors-A Comparison Between Predictions and Experiments,” Journal of Propulsion and Power, Vol. 4, No. 4, 1988, pp. 308-316. doi:10.2514/3.23068

[17] Levine, J. N., and Baum, J. D., “A Numerical Study of Nonlinear Instability Phenomena in Solid Rocket Motors,” AIAA Journal, Vol. 21, No. 4, 1983, pp. 557-564.

[18] Loncaric, S., Greatrix, D. R., and Fawaz, Z., “Star-Grain Rocket Motor Nonsteady Internal Ballistics,” Aerospace Science and Technology, Vol. 8, No. 1, 2004, pp. 47-55.

[19] Montesano, J., Behdinan, K., Greatrix, D. R., and Fawaz, Z., “Internal Chamber Modeling of a Solid Rocket Motor: Effects of Coupled Structural and Acoustic Oscillations on Combustion,” Journal of Sound and Vibration, Vol. 311, Nos. 1-2, 2008, pp. 20-38. doi:10.1016/j.jsv.2007.08.030

[20] Baczynski, C., and Greatrix, D. R., “Steepness of Grain Geometry Transitions on Instability Symptom Suppression in Solid Rocket Motor,” AIAA Paper 2009-5177, July 2009.

[21] Greatrix, D. R., “Parametric Evaluation of Solid Rocket Combustion Instability Behavior,” AIAA Paper 2009-5251, July 2009.

[22] Montesano, J., Greatrix, D., Behdinan, K., and Fawaz, Z., “Prediction of Unsteady Non-Linear Combustion Instability in Solid Rocket Motors,” Journal of Aerospace Engineering, Vol. 223, No. 7, 2009, pp. 901-

[23] Pekkan, K., and Ucer, A., “One-Dimensional Combustion Instability Studies with Moving Boundaries in an End-Burning Test Motor,” AIAA Paper 2002-3608, July 2002.

[24] d,Agostino, L., and Andernucci, M., “Unsteady Ballistic Code for Performance Prediction of Solid Propellant Rocket Engines,” AIAA Paper 2007-3131,July2007.

[25] Willcox, M. A., Brewster, M. Q., Tang, K. C., Stewart, D. S., and Kuznetsov, I., “Solid Rocket Motor Internal Ballistics Simulation Using Three-Dimensional Grain Burnback,” Journal of Propulsion and Power, Vol. 23, No. 3, 2007, pp. 575-584.

[26] Shimada, T., Hanzawa, M., Morita, T., Yoashikawa, T., and Wada, Y., “Stability Analysis of Solid Rocket Motor Combustion by Computational Fluid Dynamics,” AIAA Journal, Vol. 46, No. 4, 2008, pp. 947-957.

[27] Ferretti, V., “Numerical Simulations of Acoustic Resonance of Solid Rocket Motors,” D. Thesis, Univ. of Rome, “La Sapienza,” Sapienza, Italy, 2009.

[28] Cavallini, E., Ferretti, V., Favini, B., Di Giacinto, M., and Serraglia, F., “Pressure Oscillations Numerical Simulation in Solid Rocket Motors,” AIAA Paper 2012-3445,July2012.

[29] Shima, E., and Kitamura, K., “On New Simple Low-Dissipation Scheme of AUSM-Family for All Speeds,” AIAA Paper 2009-0136, July 2009.

[30] Boyer, G., Casalis, G., and Estivalalzes, J. L., “Stability Analysis and Numerical Simulation of Simplified Solid Rocket Motors,” Physics of Fluids, Vol. 25, No. 8, Aug. 2013, Paper 084109.

[31] Flandro, G. A., andMajdalani, J., “Aeroacoustic Instability in Rockets,” AIAA Journal, Vol. 41, No. 3, 2003, pp. 485-497.

[32] Blomshield, F., Crump, J., Mathes, H., Stalnaker, R. A., and Beckstead, M., “Stability Testing of Full-Scale Tactical Motors,” Journal of Propulsion and Power, Vol. 13, No. 3, 1997, pp. 349-355.

[33] Smith, R.,Ellis, M., Xia, G., Sankaran, V., Anderson, W., and Merkle, C. L., “Computational Investigation of Acoustics and Instabilities in a Longitudinal-Mode Rocket Combustor,” AIAA Journal, Vol. 46, No. 11, 2008, pp. 2659-2673.

[34] Anderson, J. D., Computational Fluid Dynamics: The Basics with Applications, 1st ed., McGraw-Hill, New York, 1995.

[35] Kitamura, K., and Shima, E., “Improvements of Simple Low- Dissipation AUSM Against Shock Instabilities in Consideration of Interfacial Speed of Sound,” European Conference on Computational FluidDynamics, ECCOMAS CFD Paper 1283, Portugal, 2010.

[36] Shima, E., and Kitamura, K., “Parameter-Free Simple Low-Dissipation AUSM-Family Scheme for All Speeds,” AIAA Journal, Vol. 49, No. 8, 2011, pp. 1693-1709.

[37]ANSYS Fluent Users Guide, Fluent, Inc., Lebanon, NH, Nov. 2011.

[38] Culick, F. E. C., “The Stability of One-Dimensional Motions in a Rocket Motor,” Combustion Science and Technology, Vol. 7, No. 4, 1973, pp. 165-175.

[39] Flandro, G. A., “OnFlowTurning,”AIAAPaper 1995-2730,July 1995.

[40] Majdalani, J., “Boundary-Layer Structure in Cylindrical Rocket Mo­Tors,” AIAA Journal, Vol. 37, No. 4, 1999, pp. 505-508.

[41] Matta, L., and Zinn, B., “Investigation of Flow Turning Loss in a Simulated Unstable Solid Propellant Rocket Motor,” AIAA Paper 1993-0115,July 1993.

[42] Brown, R. S., Blackner, A. M., Willoughby, P. G., and Dunlap, R., “Coupling Between Acoustic Velocity Oscillations and Solid Propellant Combustion,” Journal of Propulsion and Power, Vol. 2, No. 5, 1986, pp. 428-437.

[43] Majdalani, J., “Exact Navier-Stokes Solution for Pulsatory Viscous Channel Flow with Arbitrary Pressure Gradient,” Journal of Propulsion andPower, Vol. 24, No. 6, 2008,pp. 1412-1423. doi:10.2514/1.37815

[44] Javed, A., and Chakraborty, D., “Damping Coefficient Prediction of Solid Rocket Motor Nozzle Using Computational Fluid Dynamics,” Journal of Propulsion and Power, Vol. 30, No. 1, 2013, pp. 19-23.

[45] Culick, F. E. C., and Yang, V., “Stability Predictions in Rockets,” Nonsteady Burning and Combustion Stabiliiy of Solid Propellants, Progress in Astronautics and Aeronautics, AIAA, Washington, D.C., 1992.

[46] Vuillot, F., and Casalis, G., “Motor Flow Instabilities Part 1, Internal Aerodynamics in Solid Rocket Propulsion,” NATO Rept. AVT-096, Brussels, Belgium, 2002.

C. Oefelein Associate Editor

译自:

Quasi-One-Dimensional Modeling of Internal Ballistics and Axial Acoustics in Solid Rocket Motors

Kalyana Chakravarthy,* Arvind S. Iyer,and Debasis Chakraborty*
Defence Research and Development Laboratory, Hyderabad 500 058, India

Published Online:https://doi.org/10.2514/1.B35754

发表回复

您的电子邮箱地址不会被公开。 必填项已用*标注

You cannot copy content of this page