第二章 固体火箭发动机总体优化设计
2.1 固体火箭发动机系统参数优化
发动机性能好、质量轻、成本低和可靠性高等项指标是设计师们历来所追求的目标。过去设计发动机总体时，一般均采用陈旧的参数分析方法。设计师多半凭经验确定发动机结构形式、直径、壳体长度、燃烧室压强、喷管扩张比、药形参数及推进剂配方等。尽管对燃烧室压强、喷管扩张比等参数也采用了单变量函数求极值的方法进行优选，但由于发动机的设计质量指标与许多设计项目和设计变量有关，孤立地仅对某一设计项目和设计变量进行优选，不可能实现最优设计。再加上分析采用手工计算，因而总体设计速度慢、精度低，设计结果只能是可行方案。
利用计算机进行发动机总体优化设计，以寻求价格低廉、结构可靠、性能最优的发动机。为此，引入现代数学规划方法，建立发动机优化数学模型，对发动机性能、内弹道、药形结构和多级固体发动机进行优化设计。
图2-1 整个发动机初步设计过程
固体火箭发动机优化设计问题，属于有约束的非线性规划范畴。其数学描述是：求设计变量X＝(x1,x2,…,xn)T.，使目标函数f(X)为最小或最大，且满足约束条件：
hj(X)≥0 ，j=1,2,…,m
及 ai≤X≤bi，i=1,2,…,n
若满足minf(X) =f(X*) \(X\in R\) （1-34）
则X*为X的最优设计方案。
对于实际问题，优化不是立即进行的，它需要进行数据准备，要进行数据结果的分析，从设计开始到结束都可以认为是优化的组成部分，整个过程如图1-39所示。
由于参数寻优是上述整个过程中的重点，也是计算最关键的一部分，先前的许多计算多为数据准备，或者进行类型选择等等，所以有必要详细地给出该过程的描述，如图1-40所示。
图2-2 优化过程
正确的数学模型应能准确表达设计问题。为此，要正确选择设计变量、目标函数和约束条件，并把它们组合在一起，成为一组能准确地反映最优化设计问题实质的数学表达式；同时，所建立的数学模型要容易计算和处理。
2.1.1 优化准则和目标函数
固体火箭发动机常用的优化准则归纳如下：
（1）对于给定的有效载荷，在保证射程一定的条件下，导弹的总成本尽量最低；
（2）对于给定的有效载荷，在保证射程一定的条件下，导弹的起飞质量尽量最小；
（3）对于给定的有效载荷，在保证发动机容积（直径和长度）一定的条件下，导弹射程或主动段终点速度尽可能大；
（4）在满足发动机战术技术性能要求的条件下，发动机的可靠性最高；
（5）在满足给定的发动机容积（直径和长度）条件下，发动机的质量与总冲之比尽可能小。
发动机是导弹的重要组成部分，从系统工程观点考虑，导弹是系统，发动机是分系统。在优化设计中，一般分系统的优化准则应取系统的优化准则。不同类型的发动机，可以选用不同的优化准则，设计师可以根据自己的设计意图合理地选取，优化准则的数学表达式即为目标函数。一般应用比较多的优化准则为：对于给定的有效载荷，在保证导弹射程一定的条件下，导弹的起飞质量最小，其数学表达式为
\(\left\{ \begin{array}{l}\min {f_s}(\overline \chi )\\s.t.{}_{}^{}\mathop {}\limits_{}^{} L(\overline \chi ) = const\\f(\overline \chi ) = \frac{{{m_0}(\overline \chi )}}{{m{}_{eff}}} = (1 + \frac{{{m_m}(\overline \chi )}}{{{m_{eff}}}})\exp [\frac{1}{{{I_s}}}({V_k} + \overline g {t_a}\sin \overline \theta )]\end{array} \right.\)
（1-35）
式中 s.t——Subject to之略，意即约束条件为；
\({m_0}\left( {{\rm{\bar X}}} \right)\)——导弹起飞质量；
meff——导弹有效载荷（除发动机外的一切质量）；
\({m_m}\left( {{\rm{\bar X}}} \right)\)——发动机结构质量；
Vk——导弹飞行的最大末速度；
ta——发动机工作时间；
\(\overline{g}\)——飞行过程中重力加速度的平均值；
\(\bar{\theta }\)——导弹速度向量与当地水平面夹角的平均值。
上式可以从最大末速度公式推导得出，若忽略空气阻力的作用，导弹的最大末速度为
\({{V}_{k}}={{I}_{s}}\ln (\frac{{{m}_{0}}(\overline{x})}{{{m}_{m}}(\overline{x})+{{m}_{eff}}})-\overline{g}{{t}_{a}}\sin \overline{\theta }\)
（1-36）
由于 \({{m}_{0}}(\overline{x})={{m}_{p}}+{{m}_{m}}(\overline{x})+{{m}_{eff}}\)
所以 \({{m}_{0}}(\overline{x})=[{{m}_{m}}(\overline{x})+{{m}_{eff}}]\exp [\frac{1}{{{I}_{s}}}({{V}_{k}}+\overline{g}{{t}_{a}}\sin \overline{\theta })]\)
\(\frac{{{m}_{0}}(\overline{x})}{{{m}_{eff}}}=[1+\frac{{{m}_{m}}(\overline{x})}{{{m}_{eff}}}]\exp [\frac{1}{{{I}_{s}}}({{V}_{k}}+\overline{g}{{t}_{a}}\sin \overline{\theta })]\)
从上式可以看出，在导弹最大末速度一定的条件下，目标函数值决定于发动机的能量特性（Is）和质量特性（\({{m}_{p}},{{m}_{m}}(\overline{x})\)）。
2.1.2 能量特性
（一）多元回归分析法
迄今为止，已经制造并试验了许多不同推进剂，具有不同设计参数的固体火箭发动机，已经积累了许多有关实际比冲的资料。因此有可能对现有发动机的具体设计参数和实际比冲之间的关系，采用多元回归分析方法进行最佳拟合，得到实际比冲的计算模型。
在进行回归分析时，确定设计变量与曲线造型同时进行。所谓曲线造型就是确定回归方程的形式。最初采用了下述形式的回归方程
ξ0=\(\Sigma {{\alpha }_{i}}\ln {{x}_{i}}\) （1-37）
式中 ξ0——冲量系数，ξ0=Is,ex/Is,th；
xi——设计变量；
αi——待定的回归系数。
因为两相流损失大致与推进剂中铝粉含量成正比，因而将有铝粉含量的一项不用对数形式，而采用原参数。所以回归方程的形式改变为
\({{\xi }_{0}}=\Sigma {{\alpha }_{i}}\ln {{x}_{i}}+{{\alpha }_{j}}{{x}_{j}}\) （1-38）
式中xj——不采用对数形式的设计变量。
在曲线造型和设计变量（自变量）都确定以后，应用回归分析的方法确定经验公式（1-37）的特定系数。最后得到经验公式如下
\({{\xi }_{0}}=1.949+0.0122\ln {{D}_{t}}-0.0465\ln \beta -0.230Al\)
—0.00448lnε+0.0153ln （1-δ）—0.00842lnS
—0.005lnpc+0.00333lnρ （1-39）
式中，\({{D}_{t}}\)——喷管喉径（cm）；
\(\beta \)——反映喷管扩张损失的参数（deg）；
\(\varepsilon \) ——喷管扩张比；
Al——推进剂中铝粉含量的百分比（%）；
\(\delta \)——喷管的潜入分数，=喷管潜入长度/药柱长度；
S——喷管入口面积与喉径面积之比；
\({{p}_{c}}\)——燃烧室压强（kPa）；
\(\rho \)——喷管上游的曲率半径（mm）。
式（1-38）的标准偏差σ=0.00277，相关系数R=0.968，这两个参数都作为回归方程好坏的标准，σ越小则回归方程的预报精度越好，R值越接近1，则说明回归方程的拟合度越好。
如果减少自变量，只取Dt、β、Al、ε四个参数作为设计变量，则对应的回归方程为
\({{\xi }_{0}}=1.0502+0.0111\ln {{D}_{t}}-0.0328\ln \beta \)
\(-0.254Al-0.00617\ln \varepsilon \) （1-40）
\(\sigma =0.00342\)
\(R=0.940\)
该公式对应的稍大一些，作为减少自变量的代价。
（二）美国SPP中采用的能量计算模型
在SPP程序中有两套预估能量损失的方法。一套是在解析方法的基础上进行的，另一套采用经验公式。用解析方法预估比冲损失，计算工作量大，耗费计算机时多，因而不宜于优化设计采用。在经验公式法中，各项比冲损失与发动机设计参数之间的关系式比较直观而且简单，易于计算，因此适宜优化设计采用。
经验方法预估性能损失的计算公式包括两大部分：一是关于燃烧室冲量系数xc与发动机设计参数间的经验关系式；另一部分是喷管冲量系数xN与发动机设计参数间的关系。
\({\xi _0} = {\xi _c} \bullet {\xi _N}\) （1-41）
xN比较详细地反映了发动机设计参数对喷管比冲量损失之间的关系。
\({\xi _N} = 1 – \left( {{\varsigma _{div}} + {\varsigma _{kin}} + {\varsigma _{bl}} + {\varsigma _{tp}} + {\varsigma _{sub}}} \right)\) （1-42）
其中：zdiv喷管扩散损失；
zkin非化学平衡流损失；
zbl 附面层损失；
ztp 两相流损失；
zsub喷管潜入损失。
下面将简要介绍各项损失。
1．燃烧室冲量与发动机设计参数间的关系xc
\({{\xi }_{N}}=\left[ K+\frac{10-a}{10}\bullet \left( 1-K \right) \right]\bullet b\bullet c\) （1-43）
式中：K——与燃速有关的系数（如下表）；
a——推进剂中铝的百分数；
b——与发动机散热状况有关的系数；
c——与粘合剂种类有关的系数。
表2-1 K值与燃速的数据表
| 燃速r（mm/s） | K值 | 燃速r（mm/s） | K值 |
| 2.79 | 0.914 | 12.70 | 0.986 |
| 3.05 | 0.931 | 15.24 | 0.989 |
| 3.30 | 0.940 | 17.78 | 0.991 |
| 3.56 | 0.946 | 20.32 | 0.992 |
| 3.81 | 0.951 | 22.87 | 0.993 |
| 4.06 | 0.956 | 25.40 | 0.994 |
| 4.32 | 0.960 | 30.48 | 0.996 |
| 4.57 | 0.964 | 35.56 | 0.997 |
| 4.83 | 0.967 | 40.64 | 0.998 |
| 5.08 | 0.970 | 45.72 | 0.999 |
| 7.62 | 0.977 | 大于50.8 | 1.000 |
| 10.16 | 0.982 |
2．与喷管设计参数有关的损失
(1) 喷管扩散损失zdiv
\({{\varsigma }_{diiv}}=\frac{1}{2}\cdot \left[ 1-\cos \left( \frac{{{\alpha }_{M}}+{{\alpha }_{e}}}{2} \right) \right]\) （1-44）
在这里aM——喷管初始膨胀半角
ae——喷管出口膨胀半角
(2) 非化学平衡流损失zkin
\({{\varsigma }_{kin}}=\frac{1}{3}\cdot \left( 1-\frac{{{I}_{s,FRZ}}}{{{I}_{s,EQ}}} \right)\cdot \left( \frac{1.38}{P} \right)\) （1-45）
在这里Is,FRZ——冻结流比冲；
Is,EQ ——平衡流比冲;
P ——喷管进口压强(MPa)。
(3) 附面层损失zbl
\({{\varsigma }_{bl}}={{C}_{1}}\cdot \frac{{{P}^{0.8}}}{D_{t}^{0.2}}\cdot \left[ 1+2\cdot \exp \left( -{{C}_{2}}\frac{{{P}^{0.8}}\cdot t}{D_{t}^{0.2}} \right) \right]\cdot \left[ 1+0.016\cdot \left( \varepsilon -9 \right) \right]\) （1-46）
在这里P——喷管进口处压强（MPa）
Dt——喷管喉部直径（cm）
C1,C2——与喷管材料有关的常数，对于一般喷管
C1=0.00236,C2=0.0006052，
对于钢制喷管
C1=0.00327, C2=0；
t——工作时间（s）
e——面积膨胀比。当e<9时，1+0.0016(e–9)=1。
（4）两相流损失ztp
\({{\varsigma }_{tp}}={{C}_{3}}\cdot \frac{{{Z}^{{{C}_{4}}}}\cdot D_{p}^{{{C}_{5}}}}{{{P}^{0.15}}\cdot {{\varepsilon }^{0.08}}\cdot D_{t}^{{{C}_{6}}}}\) （1-47）
在这里C3,C4,C5,C6——计算常数，详见下表;
P ——喷管进口处压强（MPa）；
Dt——喷管喉部直径（cm）；
e ——面积膨胀比；
Z ——凝相摩尔数(mol/100g)；
Dp——凝相颗粒平均直径()，其中Dp的计算有如下三种：
1、Smith和Delaney提出的经验公式：
\(Dp = \left\{ {} \right.\begin{array}{*{20}{c}}{1.983D{t^{0.5}},2.54 \le Dt \le 25.4}\\{10,Dt > 25.4}\end{array}\)
在这里Dp的单位为微米，Dt的单位为厘米。
2、早期的SPP程序中使用的经验公式如下：
\(Dp=2.385\bullet {{P}^{\frac{1}{3}}}\bullet {{Z}^{\frac{1}{3}}}\bullet \left[ 1-\exp (-0.00157\bullet {{L}^{*}}) \right]\bullet (1+0.01772\bullet Dt)\) （1-48）
在这里P——喷管进口处压强（MPa）；
Z——凝相摩尔数(mol/100g)；
L*――发动机的特征长度(cm)
Dt——喷管喉部直径（cm）；
Dp——凝相颗粒平均直径()。
3、新的SPP程序中采用了改进的经验公式，如下所示：
\(Dp=2.7622\bullet D{{t}^{0.2932}}\bullet \left[ 1-\exp (-0.1179\bullet Z\bullet P\bullet \tau ) \right]\) （1-49）
在这里Dp——凝相颗粒的质量加权平均直径()；
Dt——喷管喉部直径（cm）；
Z——凝相摩尔数(mol/100g)；
P ——喷管进口处压强（MPa）；
表2-2 计算常数与凝相摩尔数和喉径的关系
| Z(mol/100g) | Dt(cm) | Dp () | C3 | C4 | C5 | C6 |
| <2.54 | 0.1084 | 0.5 | 1.0 | 1.0 | ||
| 2.545.08 | 0.9000 | 1.0 | 0.8 | |||
| >5.08 | <4 | 0.1339 | 0.8 | 0.8 | ||
| 48 | 0.0702 | 0.8 | 0.4 | |||
| >8 | 0.0489 | 0.8 | 0.33 | |||
| <2.54 | 0.3612 | 1.0 | 1.0 | 1.0 | ||
| 2.545.08 | 0.2998 | 1.0 | 0.8 | |||
| >5.08 | <4 | 0.4456 | 0.8 | 0.8 | ||
| 48 | 0.2340 | 0.8 | 0.4 | |||
| >8 | 0.1612 | 0.8 | 0.33 |
(5) 喷管潜入损失zsub
\({{\varsigma }_{sub}}=0.0442\cdot {{\left( \frac{P\cdot Z}{S} \right)}^{0.8}}\cdot \frac{{{L}_{N}}^{0.4}}{D_{t}^{0.2}}\) （1-50）
在这里P ——喷管进口处压强（Mpa）；
Dt——喷管喉部直径（cm）
Z凝相摩尔数(mol/100g)；
S ——喷管收缩面积比=喷管入口面积/喷管喉部面积，S>1.0；
LN——喷管的潜入深度/发动机壳体长度。
对于有些小型发动机，由于绝热层的烧蚀也将产生一部分附加的燃气，因此，引入一个由于惰性气体的流出而引起的比冲效率因子：
\({{\xi }_{b}}=\sqrt{1+{{\varphi }_{b}}}\) （1-51）
式中：{{\varphi }_{b}}\)为惰性排出物与装药质量之比。对于一些大型的贴壁浇铸的发动机尽管绝热层占据一定的质量，但并非参与燃烧，这种情况下可取\({{\varphi }_{b}}=0.0\)。而对于那些绝热层完全暴露在燃气中的发动机，则该项不能忽略。综上所述，比冲效率可表示为：
\(\xi ={{\xi }_{c}}{{\xi }_{N}}{{\xi }_{b}}\) （1-52）
2.1.3 质量特性
计算发动机质量的数学方程表达了发动机质量与设计参数之间的关系，在设计变量变化过程中，发动机的质量将随着设计变量而变化。发动机质量的大小最终将影响目标函数值，因此通过设计变量和发动机质量之间的关系，可以间接地建立设计变量与目标函数之间的关系。如果描述发动机质量特性的数学方程反映实际的结构情况越好，则计算发动机质量精度越高，从而优化计算所得的结果也越符合实际。
在对发动机质量作粗略估算时，常采用统计数据和经验公式，并且对于同一零件可能有不同的估算壁厚及质量的方法，计算方法本身也在发展和不断完善。
由于固体火箭发动机结构复杂，建立统一的发动机质量计算方程组相当困难。为了简化质量计算，将发动机分为若干个组件，再根据各个组件中零件的形状，受力情况及使用材料等因素，这样可以足够准确地反映参数变化对组件质量的影响，从而反映对整个发动机质量的影响。由此得出发动机质量计算模型主要由五部分组成，即燃烧室、绝热层、喷管、药柱和点火器及附件等。
（一）燃烧室质量mc（由圆柱段、前封头和后封头质量组成）
\({{m}_{c}}=\pi {{\rho }_{c}}{{D}_{c}}{{L}_{c}}{{\delta }_{c}}+\frac{\pi m_{1}^{2}{{f}_{\alpha }}{{f}_{p}}{{p}_{c}}{R}_{c}^{3}{{\rho }_{c}}}{(m_{1}^{2}-1){{\sigma }_{bc}}\cos {{\theta }_{2}}}[1-\frac{1}{1+(m_{1}^{2}-1){{\sin }^{2}}{{\theta }_{2}}}]+\)
\(\frac{\pi m_{1}^{2}{{f}_{\alpha }}{{f}_{p}}{{p}_{c}}R_{c}^{3}{{\rho }_{c}}}{(m_{1}^{2}-1){{\sigma }_{bc}}\cos {{\theta }_{2}}}[\frac{m_{1}^{2}R_{c}^{2}-(m_{1}^{2}-1){{y}_{2}^{2}}}{m_{1}^{2}R_{c}^{2}}-\frac{1}{1+(m_{1}^{2}-1){{\sin }^{2}}{{\theta }_{2}}}]\) （1-53）
式中 \({{\rho }_{c}}\)——壳体材料密度（kg/cm3）；
\({{f}_{\alpha }}\)——壳体安全系数；
\({{f}_{p}}\)——燃烧室最大工作压强与平均压强之比；
\({{p}_{c}}\)——燃烧室平均压强（kPa）；
\({{\sigma }_{\text{bc}}}\)——壳体材料强度极限（kPa）；
\({{\theta }_{2}}\)——从封的转角处至燃烧室轴线间的夹角取值为60°~65°；
\({{m}_{1}}\)——椭球封头长半轴与短半轴之比；
\({{y}_{2}}\)——后开口半径（cm）。
（二）绝热、衬层质量
\({{m}_{in}}=2\pi R_{c}^{3}{{\rho }_{in}}\{\frac{{{{\dot{R}}}_{a}}{{t}_{a}}}{2{{R}_{c}}}+\frac{{{\delta }_{in}}{{L}_{c}}}{R_{c}^{2}}+\frac{{{{\dot{R}}}_{c}}{{t}_{a}}}{2{{R}_{c}}(m_{1}^{2}-1)}[\frac{m_{1}^{2}R_{c}^{2}-(m_{1}^{2}-1)y_{2}^{2}}{R_{c}^{2}}-1]\}\) （1-54）
式中，\({{\rho }_{in}}\)——绝热，衬层材料密度（kg/cm3）
\({{\dot{R}}_{a}}\)——前封头绝热衬层烧蚀率（cm/s）；
\({{t}_{a}}\)——发动机工作时间（s）；
\({{\delta }_{in}}\)——燃烧室绝热衬层厚度；
\({{\dot{R}}_{c}}\)——后开口处绝热衬层炭化率。
（三）喷管质量\({{m}_{n}}\)
\({{m}_{n}}={{f}_{1}}\frac{{{m}_{p}}{{I}_{s}}}{{{\lambda }_{mn}}}\) （1-55）
式中，\({{f}_{1}}\)——比例系数；
\({{m}_{p}}\)——药柱质量；
\({{\lambda }_{mn}}\)——喷管的冲质比。
（四）药柱质量\({{m}_{p}}\)
\({{m}_{p}}=4\pi {{b}^{3}}{{\eta }_{v}}{{\rho }_{p}}/(3{{m}_{1}}+D{{L}_{c}})\) （1-56）
式中，b——药柱外半径，\(b={{R}_{c}}-{{\delta }_{c}}-{{\delta }_{in}}\)；
\({{\eta }_{v}}\)——体积装填分数，为装药容积与总容积之比。
（五）点火器、附件质量\({{m}_{j}}\)
\({{m}_{j}}={{k}_{f}}{{f}_{\alpha }}{{f}_{p}}{{p}_{c}}\) （1-57）
式中，\({{k}_{f}}\)——质量系数，\({{k}_{f}}=(4.4\tilde{\ }4.5)\)×10-3。
2.1.4 设计变量
在设计过程中进行选择，且最后必须确定的各项独立参数，称为设计变量。在选择过程中它们是变量，但这些变量确定后，则设计对象就完全确定。优化设计正是研究怎样合理地优选这些设计变量值的一种现代设计方法。在这些参数中，凡是可以根据设计要求事先给定的，是设计常量；只有那些需要在设计过程中优选的参数，才可看成是优化设计过程中的设计变量。
设计变量应选择对目标函数影响较大的变量，而且它们对目标函数有着矛盾的影响，这样，目标函数将有明显的极值存在。例如，发动机直径Dc、工作压强pc、扩张比ε、发动机工作时间ta、喷管潜入分数δ、推进剂配方和燃速r等。
设计变量应相互独立。例如，比冲Is是推进剂配方、工作压强和扩张比的函数，如果把推进剂配方、工作压强和扩张比作为设计变量，就不能将Is作为设计变量，因为它们不再是独立的。又如，当发动机体积给定时，若将直径列为设计变量，则长度不能再作为设计变量，但可以作为约束条件。
设计变量应尽量采用具有物理意义的无因次量，这样不仅便于计算，更主要的是对同类的工程问题具有通用性。
在满足设计要求条件下，应充分分析各设计变量的主次，减少变量的数目，使优化设计问题简化。此外，经常把设计变量按其对目标函数影响重要程度依次排列，在考虑节省计算机时的情况下，可以选择最重要的变量进行寻优，或作为全空间寻优的起始过程，或根本就不在全空间寻优。
一般情况，可以将设计变量缩减为10个以内，例如
\(\overline{x}={{(D,{{p}_{c}},\varepsilon ,r,\delta ,t)}^{T}}={{({{x}_{1}},{{x}_{2}},{{x}_{3}},{{x}_{4}},{{x}_{5}},{{x}_{6}})}^{T}}\)
2.1.5 约束条件
目标函数取决于设计变量，但在很多实际问题中，设计变量的取值范围是有限制的或必须满足一定的条件。在优化设计中，这种对设计变量取值的限制条件，称为约束条件或设计约束。可分为界限约束和不等式约束两类。
（一）界限约束条件
设计变量通常有一定的允许变化范围。
1．燃烧室直径Dc的界限
燃烧室直径即导弹直径Dc，通常由导弹总体综合考虑导弹的强度、刚度、气动外形、发射条件等要求而限定的。发动机设计只能在允许的变化范围内优选。
Dcmin≤Dc≤Dcmax
2．工作压强Pc的界限
工作压强Pc受到推进剂临界压强Pcr的限制，低于Pcr将出现不正常燃烧。在有些发动机中，还规定了最大压强Pmax的限制。
Pcr≤Pc≤Pmax （1-58）
3．扩张比ε的界限
对于高空工作的发动机，扩张比ε受到喷管出口直径不得大于导弹直径Dc的限制，当然除非采用可延伸扩张段的新技术；对于低空工作的发动机，扩张比还受到\({{p}_{{{c}_{\min }}}}/{{p}_{a}}\)≥0.3~0.4的限制，以免由于过膨胀，喷管内出现冲击波。因此
1≤ε≤\({{({{D}_{c}}/{{d}_{t}})}^{2}}\)
4．燃速r的界限
推进剂燃速r，可以用加入燃速调节剂，改变氧化剂颗粒大小和匹配等方法进行调节，但对于一定的基本配方，其变化范围是受限制的
rmin≤≤rmax
5．工作时间ta的界限
发动机最大工作时间tmax，受到材料耐烧蚀特性和最小推力值Fmin的限制；工作时间最小值又受到导弹允许过载系数的限制，因此
tmin≤ta≤tmax
（二）不等式约束条件
在导弹设计的战术技术要求中，除已作目标函数的设计要求外，其余的设计要求都可以作为不等式约束提出。
例如，可以把体积（长度）作为不等式的约束条件，则
L≤Lmax
可以把导弹起飞质量作为不等式约束条件，则
\({{L}_{0}}\le {{L}_{0\max }}\)
同时，研制费用和研制周期的要求，也可以作为不等式约束条件。
等式约束和不等式约束，在数学模型中可用三种形式表示不同的约束条件
\({{g}_{j}}(\overline{x})\)≥0或\({{g}_{j}}(\overline{x})\)≤0，j=1，2，…，m
\({{h}_{j}}\left( X \right)=0,j=\left( m+1 \right),\left( m+2 \right),…,p\).
\({{a}_{i}}\le {{x}_{i}}\le {{b}_{i}},i=1,2,…,n\)
2.1.6 优化方法
固体火箭发动机总体优化设计是带约束的非线性规划问题，因此只介绍非线性规划的优经方法。
无约束多维函数的优化方法，可分为两类。
一类是目标函数能用数字解析式来表达，并有连续的一阶偏导数或二阶偏导数。在这种情况下，可以利用其一阶段导数或二阶偏导数来确定寻优方向和路线，此法称为解析法，或间接法。多数工程问题用数学解析式来表示比较困难，求其一阶偏导数或二阶偏导数更困难，甚至不可能，因此间接法在工程运用中受到限制。
另一类是目标函数不能用明确的数学解析式表达，或其一阶、二阶偏导数难以求得。在这种情况下，可以借助于直接比较在局部区域内的目标函数值来确定寻优的方向和路线，此法称为直接法。前面指出，发动机总体设计的目标函数是非线性的，因无法用一个明确的数学解析式来表达，而目标函数值需要经过若干次迭代计算才能求得。因此，该问题只能用直接法来寻优。
在带约束优化问题的间接解析中，用一种特定乘子将约束优化问题转换为无约束优化问题，然后求最优解的近似方法，即罚函数法。
间接解法的基本思想是按照一定的原则构造一个包含原目标函数和约束条件的新目标函数，即将约束最优化问题的求解转换成无约束最优化问题求解。
本文用外点罚函数法嵌套Powell法求解固体火箭发动机优化设计问题，其中的一维搜索最优化方法选用黄金分割法。此外本文还提供复合型法和变量轮换法供用户选择。下面简要介绍本系统使用的三种优化方法。
（1）外点罚函数法
罚函数法是用待定乘子将约束最优化问题转换为无约束最优化问题，然后求最优解的近似但很有效的方法。外点罚函数法的特点是将罚函数定义于约束可行域之外，且求解无约束问题的探索点是从可行域外逼近原目标函数的约束最优解，它可用于既含不等式约束，又含等式约束的问题。其基本原理如下：
对于有约束的优化问题：
min f ( X ) \(X \in {R^n}\)
\({g_j}(X) = 0\) j = 1,2, …,m （1-59）
\({h_j}(X) = 0\) j = m+1,m+2,…,p
构造一个罚函数
\(\varphi \left( {X,{R^{\left( k \right)}}} \right) = f\left( X \right) + {R^{\left( k \right)}}\left( {\sum\limits_{j = 1}^m {g_j^2\left( X \right){\delta _j}\left( {{g_j}} \right) + \sum\limits_{j = m + 1}^p {h_j^2\left( X \right)} } } \right)\) （1-60）
式中： \({\delta _j}({g_j}) = 0\) \({\delta _j}(X) \ge 0\)时
\({\delta _j}({g_j}) = 1\) , \({\delta _j}(X) < 0\)时
罚因子 \({R^{(k)}}\) 是一个递增的数列，且递增率c >1，即：
\({R^{(k + 1)}} = {R^{(k)}} \cdot c\) （1-61）
这样就把原来的约束优化问题（1-59）变成了对（1-60）式求无约束极小值的问题。
即： \(\min \varphi (X,{R^{(k)}})\) \(X \in {R^n}\)
当X使全部约束条件满足时，惩罚项为零，此时(1-60 )式的无约束极小值，等价于原函数在已满足全部约束条件下的极小值。当X不满足约束条件时，则惩罚项是一个较大的正值，以起到惩罚的作用。随着罚因子的递增，逐次对j ( X,R(k) )求极小，迫使惩罚项趋于零，最后有：
\(\min f(X,{R^{(k)}}) \approx f({X^*})\) （1-62）
其算法框图见图1-41。
图2-3 外点罚函数法算法框图
用罚函数法实际上是求罚函数的极值，计算的结果并不能严格满足约束条件，会存在一定程度的偏差。
（2）Powell法
外点罚函数法求解优化问题最终还是归结为求罚函数的无约束极小值，因此还必须选定无约束优化方法。无约束优化方法也分直接法和间接法。对固体火箭发动机一体化优化设计问题，用间接法显然是不可取的，这是由于间接法需要对目标函数求一阶或二阶导数的缘故。因此本文选用直接法中的Powell法进行无约束优化。
Powell法是一种改进的共轭方向法。此方法不需对目标函数作求导计算，用于变量n<10~20或目标函数的一阶导数不连续的优化问题，均能得到很好的计算结果。因此Powell法是求无约束最优化问题较为有效的一种方法，它既在每次搜索中避免了搜索方向的线性相关性，所产生的搜索方向又是趋向相互共轭的，所以这个方法不仅对二次函数有效，而且对其它许多函数也颇为有效。
在Powell法中嵌套的一维搜索法采用了黄金分割法。该方法较常用，且对目标函数无连续性要求，适应性强。
Powell法对目标的单峰性有一定的要求，如果是多峰的，求解的结果不能保证是全局最优解，它取决于迭代的初值。初值选择的不同所收敛的区域也不同。
（3）复合型法
本文还采用了复合型法，它是求解约束非线性优化问题的一种应用比较广泛的直接解法，它源于求解无约束非线性优化问题的单纯形法，实际上是单纯形法在约束问题中的发展。
所谓复合型是指在n维设计空间内，选择n+1(或P，而n+1≤P≤2n)个初始点，构造一个初始复合型，只是这些初始点所构成的复合型要位于受约束条件限制的可行域内。对于二维问题，复合型为由三（即n+1）个顶点构成的三角形，或由四（即2n）个顶点构成四边形。对于三维问题，复合型为由四（即n+1）个顶点构成的四面体，或由五个顶点或六（即2n）个顶点构成五边形。对于n维问题，复合型则为n+1~2n个顶点构成一个不规则的多面体。
在复合型中，比较这些顶点的目标函数值，取最大者为坏点，以其余各点的中心为映射的轴心，寻找坏点的映射点，一般说来，此映射点的目标函数值总是小于坏点。这时以映射点替换坏点构成新的复合型。按照这个步骤重复多次使复合型的位置逐步移向最优点附近，同时复合型也不断收缩。当其收缩到满足精度要求时，输出复合型顶点中目标函数值为最小的点最为最优解。这种方法不必保持规则图形，较为灵活，所求的结果也较为可靠，收敛也比较快。
初始复合型的顶点应该都选择在可行域内，而且为了不使问题降维，在P个顶点中至少有n+1个点的n个矢量是线性无关的。
初始复合型顶点的产生有三种：
（1）全部顶点利用随机方法产生；
（2）全部顶点人为地预先按实际情况确定；
（3）人为地预先确定一个顶点，其余顶点按随机方法产生。
图2-4 Powell算法框图
这三种方法中，第一种方法最简便，但收敛性差；第二种方法收敛性好但要在可行域内人为地确定全部顶点比较麻烦；第三种方法介于两者之间，只要在可行域内选择一个顶点，其余各点随机产生，此法比较简单，但收敛性略差一些。对于较复杂问题，人为确定可行点比较困难，所以本文采用第一种方法。
用复合型法计算的结果将严格满足边界条件，只要能收敛，结果是比较可靠的。
（4）变量轮换法
本文还采用了变量轮换法，首选在全范围内利用随机方法产生一个点，该点必须满足所有的约束条件。然后以该点作为起始点，采用变量轮换法轮流在各个方向上搜索，直到搜索出最小目标函数，过程即告结束（如图1-43）。
图2-5 变量轮换法求函数极值
固体火箭发动机优化设计，不论采用什么样的方法，都是以满足某一设计准则为条件，要求某一设计参数为最优为目的的，所以它属于最优准则的优化设计类型。作为优化方法的数学基础已经很成熟，而且有现成的程序，作为应用部门不必要花费太大的精力去研制优化方法程序，而应把注意力集中在所研究问题的数学模型建立和优化准则的确定。
2.1.7 参数分析
参数分析就是暂时固定其余的设计变量（参数），有目的地改变所研究的参数，计算优化设计的目标函数或有关性能，研究这些参数变化对目标函数或有关性能的影响。
参数分析有下列两种主要形式：
1．设计变量偏离最优点的影响
当求出最优点（最优设计方案）之后，以这个最优点为基础，有规律地改变其中某一个设计变量，而暂时固定其余变量，研究该设计变量偏离最优点对目标函数及有关性能的影响。这种研究所观察的方案多数不是最优的。目的在于观察某些设计变量偏离最优点所造成的后果。因此，大多用于方案确定之后，研究深化过程中对控制各参数的影响，并为修改设计时提供参考。
现以某发动机的参数分析结果来说明。如图1-44所示，为发动机直径与目标函数的关系曲线。图1-45为工作时间t与目标函数的关系曲线。
图2-6 Dc与目标函数（obj）的关系
图2-7 t与目标函数（obj）的关系
这种参数分析曲线，又叫做参数敏感性分析曲线。图1-44表示参数偏离最优点时对目标函数值的影响大，很敏感。对这种参数，在整个设计过程中，以及在试制的工艺过程、使用过程中，都应十分注意，使它尽量保持在最佳值上。相反，有的参数敏感曲线很平缓，如图2-7所示。这类参数在设计和制造中都可放宽要求。
2．改变设计要求（约束条件或性能指标）对最优设计的影响
研究方法与上面所讲的方法相反。在这种情况下，改变某一项设计要求，按一定间隔改变它的数值，而每一次改变都可以求得其最优设计方案。比较这些方案的判别，可以看出各个设计要求对最优设计的影响，从而可以进一步论证设计要求的合理性和现实性。
2.2 固体动力系统推力方案优化
随着固体动力系统向着智能化、能量管理、多脉冲、可调可控和吸气式、组合式的方向发展，飞行器总体设计与推进系统设计的耦合程度越来越高。而且，由于固体动力系统推力方案与飞行器外弹道的匹配，更多的时候需要在研制的过程中预先设定。因此，飞行器和发动机研究单位越来越强烈的认识到孤立的动力系统优化设计结果，对于飞行器总体而言，往往会使得飞行器整体性能更为低劣。飞行器总体基于传统总冲、推力方案、质量分数、冲质比等战术性能指标的设计任务要求无法实现精细化的描述其对动力系统优化设计的需求。
为了解决飞行器设计迭代过程中的冲突，优化飞行器总体性能，固体动力系统推力方案优化最大的特征是以飞行器总体性能最优，作为固体动力系统推力方案优化的目标，以此来优化固体动力系统的性能参数和内部细节结构。
图 2-2-1固体动力系统推力方案优化流程
2.2.1 推力方案优化实例
以一个临近空间飞行器的火箭助推为例，下面的流程图展示了该算例的内外弹道一体化计算过程。
图 2-2-2 内外弹道一体计算过程
算例中采用的飞行器（含助推级）外形如图2-2-3示，这里只针对助推级进行发动机推力方案优化计算。
图 2-2-3、算例中导弹的外形结构
飞行器总长度为8.17m，其中助推级的翼展为2.5m，飞行器（含助推级发动机）总质量为1.27t，助推级发动机的总冲1520N·s，其工作时间为8s，发动机总重量为900kg，装药质量为641kg，平均推力为190kN，其工作的推力时间曲线如图2-2-4所示。
图 2-2-4 助推级发动机推力-时间曲线
分别以67deg发射角度进行内外弹道联合计算，进行内外弹道一体化计算，计算出的导弹飞行速度时间及高度时间曲线如图2-2-5所示。
图 5 优化前外弹道计算结果
针对助推级发动机，在保证总冲及装药质量不变的条件下，对助推级进行推力方案优化计算，优化后获得的推力曲线为单室双推力，总工作时间为28s，发动机总重量为900kg，装药质量也为641kg，第一级平均推力为170kN，第二级平均推力为38kN，优化后的推力时间曲线如图2-2-6所示。
图 2-2-6、优化后发动机推力-时间曲线
同样对其在67deg，不考虑风作用的条件下进行外弹道求解。求解后的外弹道计算速度时间及高度时间曲线如图2-2-7所示。
图 2-2-7、优化后发动机外弹道计算结果
对比优化前及优化后的外弹道计算结果，选取工作时间最长的发动机关机时刻(28s)为对比点,可以发现：对于单推力外弹道其速度为839m/s，而双推力外弹道其关机速度可达到1070m/s。同时还选取当导弹飞行分别达到30公里时刻：单推力外弹道其速度为723m/s，而双推力外弹道速度为822m/s。
对比上述计算结果可以看出：对于单室单推力发动机，当达到其最大速度时，导弹仍飞行处于低空飞行，此时由于空气相对密稠，速度大会引起导弹的气动加热大，而对于单室双推却能够较好解决上述问题。同时对达到相同飞行高度时，单室双推能够具有相对较高的速度，因此也可以看出采用双推力发动机相对而言能够提高导弹飞行的续航能力。
2.2.1 结构参数驱动的发动机性能曲线详细优化
影响发动机推力方案最大的因素是发动机内部的装药结构。对于目前发动机研制领域普遍采用的三维装药结构。在设计的过程中有大量的三维结构参数需要设计人员进行协调，才能取得与外弹道匹配性最好的内弹道和推力曲线。例如图2-2- 8所示的翼柱形装药结构就具有17个独立的药柱结构参数。
图 2-2-8 典型翼柱形药柱结构示意
这些参数对设计结果产生的影响难以简单的进行评价，设计人员也很难手动的对如此大量的设计参数进行协调。为了解决这一难题，需要利用计算机自动寻优方法，进行自动的发动机燃面调优。优化设计流程如图 2-2-9所示。
图 2-2-9 燃面优化设计流程
由于固体火箭发动机药柱燃面的优化设计，是一个高度非线性、不可微、不连续、多峰且包括众多约束的问题，一般使用探索性解法进行求解，例如多岛遗传算法、蚁群算法、模拟退火算法等等。
而对发动机推力曲线的评价，则可以采用下面的最小二乘曲线拟合方法来计算优化目标函数。
其中，Pc[i]表示在第i个采样点Pc的数据。显然，较为理想的平直内弹道，应具有尽量小的E和S，以及尽量大的T和Ttail，四者的权重系数设置如表1。
表1 目标函数权重系数
| T | E | S | Ttail |
| 3 | 10 | 15 | 1 |
下面展示一个翼柱形装药的性能曲线优化为例，展示这一优化过程。
该翼柱形装药的初始尺寸为：外轮廓：L =2100mm，r1=225mm，r2=450mm；内孔：L1=70mm，L2=105mm，L3=210mm，L4=80mm，L5 =1780mm，r3=20mm，r4=15mm；翼槽：α=30°，β=17°，γ=60°，R1 =20mm，R2=110mm，R3=10mm，H1=200mm，H2=210mm，T=30mm，L6 =1200mm，翼槽个数为10。发动机壳体壁厚4mm，绝热层厚5mm，头部脱粘深度为药柱外径的50%，尾部脱粘深度为药柱外径的50%。
（a）药型参数
（b）药型三维结构
图 2-2-10 翼柱形装药结构
优化中只将翼槽的10个尺寸作为设计变量，来进行参数化建模，变量取值为整数，范围见下表。
表2 设计变量的优化值
| 设计变量 | 上限 | 下限 | 优化值 |
| α(度) | 17 | 40 | 37 |
| β(度) | 6 | 23 | 21 |
| γ(度) | 46 | 69 | 68 |
| R1 (mm) | 10 | 120 | 65 |
| R2 (mm) | 20 | 130 | 121 |
| R3 (mm) | 5 | 10 | 10 |
| H1 (mm) | 100 | 200 | 110 |
| H2 (mm) | 200 | 400 | 240 |
| T (mm) | 20 | 40 | 38 |
| L6 (mm) | 1000 | 1400 | 1324 |
初始药形的燃面-肉厚曲线如下图，放映曲线平直程度的方差值为0.1634。
图 2-2-11初始燃面-肉厚曲线
采用遗传算法，目标函数的迭代收敛过程如下图所示。
图 2-2-11评价函数收敛过程
优化得到的药形，其燃面-肉厚曲线如图9，方差值为0.0423，优化后药柱的燃面-肉厚曲线值的方差减小了0.1211，发动机燃面-肉厚曲线更加平直，从而获得了性能更好地发动机推力方案曲线。
图 2-2-12优化后的燃面-肉厚曲线
