基于遗传算法与能量特性图形的固体推进剂配方优化方法.pdf

本发明公开了一种基于遗传算法与能量特性图形的固体推进剂配方优化方法，包括以下步骤：步骤一、建模：根据最小自由能原理，建立固体推进剂的四个能量特性计算模型；步骤二、初始参数设定与存储；步骤三、采用遗传算法进行配方优化设计；步骤四、图形绘制用取值范围设定及各原料类别的质量含量计算；步骤五、能量特性图形绘制，其绘制过程如下：能量特性参数输入、能量特性曲线方程拟合、能量特性图形绘制和能量特性图形同步显示。本发明方法步骤简单、设计合理且实现方便、使用效果好，将遗传算法与能量特性图形相结合进行固体推进剂配方优化设计，能有效解决现有固体推进剂配方设计过程中存在的能量特性试验成本高、周期长、试验量大等缺陷。

权利要求书一种基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于该方法包括以下步骤：
步骤一：建模：根据最小自由能原理，建立固体推进剂的四个能量特性计算模型；其中，四个所述能量特性计算模型分别为根据固体推进剂的质量配比对该固体推进剂的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量，分别进行计算的比冲计算模型、特征速度计算模型、燃烧室温度计算模型和燃烧产物平均相对分子质量计算模型；
步骤二、初始参数设定与存储：首先，通过与数据处理器相接的参数输入单元，输入配制所设计固体推进剂所用的组分种数N和各组分的组分信息，并将所输入的信息同步存储至与所述数据处理器相接的数据存储单元内；其中，各组分的组分信息均包括化学式、所属原料类别和质量配比范围（mi0～miz），i=1、2、…、N，N为配制所设计固体推进剂所用的组分数量，mi0为配制所设计固体推进剂所用第i种组分的质量配比下限值且miz为第i种组分的质量配比上限值；其中，mi0和miz均为实数，且0＜mi0＜100，0＜miz＜100，N≥2；其中，所属原料类别为粘合剂、氧化剂或添加剂；待N种组分的组分信息均输入完毕后，再采用所述参数输入单元对所设计固体推进剂的原料类别总数Q进行设定，且Q=2或3；
之后，通过所述参数输入单元自预先建立的燃烧产物数据库中选出所设计固体推进剂燃烧后所产生的所有燃烧产物；所述燃烧产物数据库中存储有多种燃烧产物的属性信息；其中，各燃烧产物的属性信息均包括化学式、相对分子质量和相态，其中相态为气相或凝聚相；
步骤三、采用遗传算法进行配方优化设计：所述数据处理器调用遗传算法模块，对所设计固体推进剂的质量配比进行优化设计，其优化设计过程如下：
步骤301、种群初始化：将配制所设计固体推进剂所用的N种组分中各组分的一个质量配比取值mix作为一个个体，并将多个个体集合为一个种群，同时所述种群中的所有个体均进行二进制编码后形成初始化种群，且将所形成的初始化种群同步存储至所述数据存储单元内；其中，每一个个体均包括N个质量配比取值mix，i=1、2…、N，且m1x+m2x+…+mNx=100，m1x、m2x、…、mNx分别为N种组分的质量配比取值；mix为第i种组分的质量配比取值且其为自该组分的质量配比范围（mi0～miz）中随机抽取一个数值；
步骤302、初始化种群中各个体适应度值计算：初始化种群中所有个体的适应度值计算方法均相同；对于所述初始化种群中的任一个个体来说，均先利用步骤一中所建立的四个所述能量特性计算模型对该个体的四个能量特性参数进行计算并对计算结果进行同步存储，并将计算得出的该个体的比冲值作为该个体的适应度值；其中，四个所述能量特性参数分别为该个体的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量；
待所述初始化种群中所有个体的四个所述能量特性参数均计算完毕后，所述初始化种群中所有个体的适应度值均计算出来；之后，再相应计算出所述初始化种群的种群平均适应度值；
步骤303、选择操作：根据步骤302中计算得出的所述初始化种群中所有个体的适应度值，选出所述初始化种群中适应度值高的多个个体作为子代群体；
步骤304、交叉操作与变异操作：对步骤303中选出的子代群体进行交叉操作与变异操作，获得新一代的子代群体，并将所获得的子代种群同步存储至所述数据存储单元内；
步骤305、子代群体中各个体适应度值计算：对步骤304中所获得的子代群体中各个体的适应度值进行计算，且所述子代群体中所有个体的适应度值计算方法均相同；对于所述子代群体中的任一个个体来说，均先利用步骤一中所建立的四个所述能量特性计算模型对该个体的四个能量特性参数进行计算并对计算结果进行同步存储，之后将计算得出的该个体的比冲值作为该个体的适应度值；其中，四个所述能量特性参数分别为该个体的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量；
待所述子代群体中所有个体的四个所述能量特性参数均计算完毕后，所述子代群体中所有个体的适应度值均计算出来；之后，再相应计算出所述子代群体的种群平均适应度值；
步骤306、选择操作：根据步骤305中计算得出的所述子代群体中所有个体的适应度值，选出所述子代群体中适应度值高的多个个体作为子代群体；
步骤307、判断是否满足终止条件：当进化代数超过预先设定的最大进化代数itex或者子代群体中个体的最大适应度值大于或等于预先设定的适应度设定值时，遗传算法终止并输出当前所获得所述子代群体中适应度值最高的个体；否则，返回步骤304；
步骤307中所输出的适应度值最高的个体，便为通过遗传算法所设计的固体推进剂质量配比m1、m2、…、mN；
步骤四、图形绘制用取值范围设定及各原料类别的质量含量计算，其计算过程如下：
步骤401、原料类别质量含量取值范围及能量特性参数取值范围设定：当Q=2时，采用所述参数输入单元输入配制所设计固体推进剂所用粘合剂的质量含量取值范围（m粘0～m粘z）和所用氧化剂的质量含量取值范围（m氧0～m氧z），并同步存储至所述数据存储单元内；当Q=3时，采用所述参数输入单元输入配制所设计固体推进剂所用粘合剂的质量含量取值范围（m'粘0～m'粘z）、所用氧化剂的质量含量取值范围（m'氧0～m'氧z）和所用添加剂的质量含量取值范围（m'添0～m'添z），并同步存储至所述数据存储单元内；同时，采用所述参数输入单元分别输入四个所述能量特性参数的取值范围，并同步存储至所述数据存储单元内；
步骤402、原料类别质量含量计算：所述数据处理器根据步骤二中所输入的各组分的原料类别，调用原料类别质量含量计算模块，对步骤301中所述初始化种群中和步骤304中所获得的子代群体中各个体的原料类别质量含量分别进行计算，并对计算结果进行同步存储，且所有个体的原料类别质量配比计算方法均相同；
对于步骤301中所述初始化种群中和步骤304中所获得的子代群体中任一个个体来说，当Q=2时，根据该个体中所包含的N种组分的质量配比取值，先计算得出该个体中所属原料类别为粘合剂的所有组分的质量配比取值之和m粘，之后根据公式m粘=100‑m氧，计算得出该个体中所属原料类别为氧化剂的所有组分的质量配比取值之和m氧；其中，m粘为该个体的粘结剂的质量含量，且m氧为该个体的氧化剂的质量含量；
对于步骤301中所述初始化种群中和步骤304中所获得的子代群体中任一个个体来说，当Q=3时，根据该个体中所包含的N种组分的质量配比取值，先计算得出该个体中所属原料类别为粘合剂的所有组分的质量配比取值之和m’粘，再计算得出该个体中所属原料类别为氧化剂的所有组分的质量配比取值之和m’氧，之后再根据公式m’添=100‑m’粘‑m’氧，计算得出该个体中所属原料类别为添加剂的所有组分的质量配比取值之和m’添；其中，m’粘为该个体的粘结剂的质量含量，m’氧为该个体的氧化剂的质量含量，且m’添为该个体的添加剂的质量含量；
步骤五、能量特性图形绘制，其绘制过程如下：
步骤501、能量特性参数输入：采用所述参数输入单元输入需通过能量特性图形进行分析的能量特性参数y，其中y为比冲、特征速度、燃烧室温度或燃烧产物平均相对分子质量；
步骤502、能量特性曲线方程拟合：所述数据处理器调用曲线方程拟合模块，根据步骤302中计算得出的所述初始化种群中各个体的四个能量特性参数和步骤305中所述子代群体中各个体的四个能量特性参数，并结合步骤四中计算得出的所述初始化种群中和步骤304中所获得的子代群体中各个体的原料类别质量含量，拟合出能量特性参数y的能量特性曲线方程；
其中，当Q=2时，拟合出的能量特性曲线方程为y=a·x1+b·x2+c·x1·x2+d（4），式（4）中a、b、c和d均为常数，x1为配制所设计固体推进剂所用粘合剂的质量含量，x2为配制所设计固体推进剂所用氧化剂的质量含量，且x1+x2=100；
当Q=3时，拟合出的能量特性曲线方程为y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5），式（5）中a'、b'、c'、d'、e'、f'、g'和h'均为常数，且x'1为配制所设计固体推进剂所用粘合剂的质量含量，x'2为配制所设计固体推进剂所用氧化剂的质量含量，x'3为配制所设计固体推进剂所用添加剂的质量含量，且x'1+x'2+x'3=100；
步骤503、能量特性图形绘制：所述数据处理器调用能量特性图形绘制模块，且根据步骤502中拟合出的能量特性参数y的能量特性曲线方程，对能量特性参数y的能量特性图形进行绘制；
当Q=2时，所述数据处理器调用能量特性图形绘制模块，绘制出能量特性参数y的等高线图，所绘制的等高线图上包括能量特性参数y取不同值时的多条等高线，每一条所述等高线均为由能量特性参数y的取值相等的各点所连成的曲线，且多条所述等高线上均标注有该条等高线所对应的能量特性参数y的取值；实际对能量特性参数y的等高线图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x1，且纵坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x2的平面直角坐标系，再根据步骤401中所输入的（m粘0～m粘z）和（m氧0～m氧z），对所建立平面直角坐标系中横坐标轴和纵坐标轴的刻度分别进行设置；之后，根据公式y=a·x1+b·x2+c·x1·x2+d（4），在所建立的平面直角坐标系上绘制出能量特性参数y取不同值时的多条等高线；
当Q=3时，所述数据处理器调用能量特性图形绘制模块，绘制出能量特性参数y的等性能三角图，所绘制的等性能三角图上包括能量特性参数y取不同值时的多条等高线，每一条所述等高线均为由能量特性参数y的取值相等的各点所连成的曲线，且多条所述等高线上均标注有该条等高线所对应的能量特性参数y的取值；实际对能量特性参数y的等性能三角图进行绘制时，先建立一个等边三角形，并将所建立等边三角形的三条边分别定义为配制所设计固体推进剂所用粘合剂的质量含量x'1、配制所设计固体推进剂所用氧化剂的质量含量x'2和配制所设计固体推进剂所用添加剂的质量含量x'3；之后，根据步骤401中所输入的（m'粘0～m'粘z）、（m'氧0～m'氧z）和（m'添0～m'添z），在所建立三角形的三条边上分别设置刻度；然后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x1'·x'2·x'3+h'（5），在所建立的等边三角形上绘制出能量特性参数y取不同值时的多条等高线；
步骤504、能量特性图形同步显示：通过与所述数据处理器相接的显示单元，对步骤503中绘制出的能量特性图形进行同步显示。
按照权利要求1所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤504中对绘制出的能量特性图形进行同步显示后，还需人为对所显示的能量特性图形进行分析，并根据分析结果对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
按照权利要求1或2所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的等高线图；
实际对能量特性参数y的等高线图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x'1，且纵坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x'2的平面直角坐标系；之后，根据步骤401中所输入的（m'粘0～m'粘z）和（m'氧0～m'氧z），对所建立平面直角坐标系中横坐标轴和纵坐标轴的刻度分别进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，在所建立的平面直角坐标系上绘制出能量特性参数y取不同值时的多条等高线。
按照权利要求1或2所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤503中当Q=2时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的二维能量特性图；实际对二维能量特性图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x1，且纵坐标轴为所设计固体推进剂的能量特性参数y的平面直角坐标系；之后，再根据步骤401中所输入的（m粘0～m粘z），对所建立平面直角坐标系中横坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立平面直角坐标系中纵坐标轴的刻度进行设置；之后，根据公式y=a·x1+b·x2+c·x1·x2+d（4）和公式x1+x2=100，在所建立的平面直角坐标系上绘制出能量特性参数y的取值随x1的取值变化而变化的曲线；
步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的二维能量特性图；实际对二维能量特性图进行绘制时，先通过参数输入单元将x'1、x'2和x'3三个参数中的一个参数设定为固定值，之后再自x'1、x'2和x'3中另外两个参数中选取一个参数作为横坐标轴，并以所设计固体推进剂的能量特性参数y作为纵坐标轴，建立平面直角坐标系；之后，再根据步骤401中所输入的原料类别质量含量取值范围，对所建立平面直角坐标系中横坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立平面直角坐标系中纵坐标轴的刻度进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，在所建立的平面直角坐标系上绘制出一个参数固定条件下能量特性参数y的取值随另一个参数的取值变化而变化的曲线。
按照权利要求1或2所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的三维能量特性图；所绘制的三维能量特性图上包括多条能量特性曲线一和多条能量特性曲线二；
实际对三维能量特性图进行绘制时，先建立一个X坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x'1、Y坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x'2且Z坐标轴为所设计固体推进剂的能量特性参数y的三维空间直角坐标系；其次，根据步骤401中所输入的原料类别质量含量取值范围，对所建立三维空间直角坐标系中X坐标轴和Y坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立三维空间直角坐标系中Z坐标轴的刻度进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，沿X坐标轴由左至右绘制出多条x'1取不同值时能量特性参数y的取值随x'2的取值变化而变化的能量特性曲线一，每一条所述能量特性曲线一上x'1的取值均相等；然后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x1'·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，沿Y坐标轴由前至后绘制出多条x'2取不同值时能量特性参数y的取值随x'1的取值变化而变化的能量特性曲线二，每一条所述能量特性曲线二上x'2的取值均相等。
按照权利要求1或2所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤一中建立所述比冲计算模型、特征速度计算模型、燃烧室温度计算模型和燃烧产物平均相对分子质量计算模型之前，还需根据最小自由能原理建立最小自由能数学模型；
所建立的最小自由能数学模型为式（2）中：j=1、2、…、A，A为固体推进剂所含化学元素的种类；s=1、2、…、m，m为固体推进剂燃烧后处于化学平衡状态时所含燃烧产物的种数；其中μs为预先通过所述参数输入单元输入的第s种燃烧产物的化学势（KJ/mol），ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg）且ns≥0，asj为1mol第s种燃烧产物中所含第j种化学元素的原子数；bj为1000g固体推进剂中所含第j种化学元素的原子数，πj为拉格朗日乘子；
步骤一中所建立的燃烧产物平均相对分子质量计算模型为（6），式（6）中n为1000g气相燃烧产物的摩尔数，且n=1000g固体推进剂燃烧后处于化学平衡状态时所含的所有气相燃烧产物的摩尔数（mol/Kg）之和；
所建立的燃烧室温度计算模型为绝热温度计算模型，且所建立的绝热温度计算模型为式（3）中Hc1为绝热温度T=T1时的1000g燃烧产物的焓，Hc2为绝热温度T=T2时的1000g燃烧产物的焓，Hc1＜Hp＜Hc2，并且Tc1和Tc2均根据公式其中ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg），Hcs为1mol第s种燃烧产物在绝热温度为T时的焓，Hcs=RT(αs1+αs2T/2+αs3T2/3)+αs4T3/4+αs5T4/5+αs6T5/6），其中，R为通用气体常数（Kg·m/mol·K），T为绝热温度，αs1、αs2、αs3、αs4、αs5和αs6为预先通过所述参数输入单元输入的第s种燃烧产物的热力学函数温度系数；其中Mi为配制固体推进剂所用第i种组分的相对分子质量，Hi为预先通过所述参数输入单元输入的1mol第i种组分的焓，Wi为配制固体推进剂所用第i种组分的质量百分数；
所建立的比冲计算模型为理论比冲计算模型，且所建立的理论比冲计算模型为式（1）中：Isp为理论比冲（N·s/Kg），Tc为燃烧室温度（K），Pe为发动机喷管出口处压强（Pa），Pc为燃烧室内压强（Pa），为燃烧产物平均相对分子质量，R为通用气体常数（Kg·m/mol·K），k为绝热指数且其定压比热与定容比热之比，其中，Pe和Pc为预先通过所述参数输入单元输入的发动机设计参数；相应地，步骤302中对于所述初始化种群中的任一个个体来说，均以利用步骤一中所建立的理论比冲计算模型计算得出的该个体的理论比冲值作为该个体的适应度值；步骤305中对于所述子代群体中的任一个个体来说，均以利用步骤一中所建立的理论比冲计算模型计算得出的该个体的理论比冲值作为该个体的适应度值；
所建立的特征速度计算模型为式（5）中n为1000g气相燃烧产物的摩尔数，R为通用气体常数（Kg·m/mol·K），Tc为燃烧室温度（K），k为绝热指数且其定压比热与定容比热之比，
按照权利要求6所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤302和步骤305中利用步骤一中所建立的四个所述能量特性计算模型，对任一个个体的四个能量特性参数进行计算时，其计算过程如下：
步骤Ⅰ、平衡组成计算：所述数据处理器调用参数计算模块，根据当前所计算的个体中各组分的化学式和质量配比取值m1x、m2x…mNx，并结合各组分的相对分子质量，对asj和bj进行计算；之后，所述数据处理器结合预先输入的μs和πj，且根据公式计算得出ns，此时便获得当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成；
步骤Ⅱ、燃烧室温度计算：所述数据处理器结合步骤Ⅰ中计算得出的当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成，并根据公式计算得出燃烧室温度Tc；其中，其中Mi为当前所计算个体所对应固体推进剂配方中第i种组分的相对分子质量，Hi为1mol第i种组分的焓，Wi为当前所计算个体所对应固体推进剂配方中第i种组分的质量百分数；
步骤Ⅲ、燃烧产物平均相对分子质量计算：先根据步骤一中推算得出的当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成，计算得出当前所计算个体所对应固体推进剂配方燃烧后所产生的1000g气相燃烧产物的摩尔数n，所述数据处理器公式并结合计算得出当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平均相对分子质量；
步骤Ⅳ、理论比冲值计算：所述数据处理器根据公式并结合步骤Ⅱ中计算得出的燃烧室温度Tc、步骤Ⅲ中计算得出的燃烧产物平均相对分子质量以及预先通过所述参数输入单元输入的发动机设计参数Pe和Pc，计算得出当前所计算个体所对应固体推进剂配方的理论比冲值Isp；
步骤Ⅴ、特征速度计算：所述数据处理器根据公式并结合步骤Ⅱ中计算得出的燃烧室温度Tc和步骤Ⅲ中计算得出的1000g气相燃烧产物的摩尔数n，计算得出当前所计算个体所对应固体推进剂配方的特征速度C*。
按照权利要求2所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤504中对绘制出的能量特性图形进行同步显示后，人为对所显示的等高线图或等性能三角图进行分析时，先自所述等高线图或等性能三角图中找出能量特性参数y取值最佳的区域，并根据该区域所对应x1与x2的取值范围或者x'1、x'2和x'3的取值范围，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
按照权利要求1或2所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤307中最大进化代数itex=20～100，步骤301中所述初始化种群中所包含个体的数量为20～30。
按照权利要求1或2所述的基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于：步骤304中进行交叉操作与变异操作时，按照交叉概率pc进行交叉操作，且按照变异概率pm进行变异操作；其中，交叉概率pc=0.6～0.9，变异概率pm=0.1～0.2；且步骤303选择操作中，采用最优个体保留策略来获得最优个体。

说明书基于遗传算法与能量特性图形的固体推进剂配方优化方法
技术领域
本发明属于固体推进剂配方优化设计技术领域，尤其是涉及一种基于遗传算法与能量特性图形的固体推进剂配方优化方法。
背景技术
随着国防（航天、兵器等）科技的发展对推进剂的要求越来越严格，要求其能量高，性能优良。众所周知，在保持推进剂基本性能的条件下，推进剂的能量越高，射程越远。如ICBM洲际导弹（9260千米）所用的推进剂的比冲增加1%，射程增加7.3%（676千米）；当推进剂的比冲增加5%时，射程增加45%（4167千米）。原来为9260千米射程的洲际导弹，只要改进的推进剂比冲提高5%，即由2601N·S/Kg提高至2731N·S/Kg，射程就将提高到13427千米，这在战略和战术上都有极为重要的意义。为提高推进剂的比冲等性能，需要做很多实验，花费成千上亿元的资金，而用优化设计方法可以少做许多试验，节省大量的人力物力和研制经费，能缩短推进剂配方研制周期，而高性能推进剂的研制成功对提高我国的航天技术水平，对提高国威军威、振兴中华有重要意义。在相同的装药条件下，推进剂的能量(即比冲)提高直接影响到火箭和导弹的射程。综上，现如今急需一种方法步骤简单、实现方便、投入成本低且使用操作简便、所需时间短、能快速最高比冲的固体推进剂最佳配比，从而大幅度缩短推进剂配方的研制周期。
发明内容
本发明所要解决的技术问题在于针对上述现有技术中的不足，提供一种基于遗传算法与能量特性图形的固体推进剂配方优化方法，其方法步骤简单、设计合理且实现方便、使用效果好，能有效解决现有固体推进剂配方设计过程中存在的能量特性试验成本高、周期长、试验量大等缺陷。
为解决上述技术问题，本发明采用的技术方案是：一种基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征在于该方法包括以下步骤：
步骤一：建模：根据最小自由能原理，建立固体推进剂的四个能量特性计算模型；其中，四个所述能量特性计算模型分别为根据固体推进剂的质量配比对该固体推进剂的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量，分别进行计算的比冲计算模型、特征速度计算模型、燃烧室温度计算模型和燃烧产物平均相对分子质量计算模型；
步骤二、初始参数设定与存储：首先，通过与数据处理器相接的参数输入单元，输入配制所设计固体推进剂所用的组分种数N和各组分的组分信息，并将所输入的信息同步存储至与所述数据处理器相接的数据存储单元内；其中，各组分的组分信息均包括化学式、所属原料类别和质量配比范围（mi0～miz），i=1、2、…、N，N为配制所设计固体推进剂所用的组分数量，mi0为配制所设计固体推进剂所用第i种组分的质量配比下限值且miz为第i种组分的质量配比上限值；其中，mi0和miz均为实数，且0＜mi0＜100，0＜miz＜100，N≥2；其中，所属原料类别为粘合剂、氧化剂或添加剂；待N种组分的组分信息均输入完毕后，再采用所述参数输入单元对所设计固体推进剂的原料类别总数Q进行设定，且Q=2或3；
之后，通过所述参数输入单元自预先建立的燃烧产物数据库中选出所设计固体推进剂燃烧后所产生的所有燃烧产物；所述燃烧产物数据库中存储有多种燃烧产物的属性信息；其中，各燃烧产物的属性信息均包括化学式、相对分子质量和相态，其中相态为气相或凝聚相；
步骤三、采用遗传算法进行配方优化设计：所述数据处理器调用遗传算法模块，对所设计固体推进剂的质量配比进行优化设计，其优化设计过程如下：
步骤301、种群初始化：将配制所设计固体推进剂所用的N种组分中各组分的一个质量配比取值mix作为一个个体，并将多个个体集合为一个种群，同时所述种群中的所有个体均进行二进制编码后形成初始化种群，且将所形成的初始化种群同步存储至所述数据存储单元内；其中，每一个个体均包括N个质量配比取值mix，i=1、2…、N，且m1x+m2x+…+mNx=100，m1x、m2x、…、mNx分别为N种组分的质量配比取值；mix为第i种组分的质量配比取值且其为自该组分的质量配比范围（mi0～miz）中随机抽取一个数值；
步骤302、初始化种群中各个体适应度值计算：初始化种群中所有个体的适应度值计算方法均相同；对于所述初始化种群中的任一个个体来说，均先利用步骤一中所建立的四个所述能量特性计算模型对该个体的四个能量特性参数进行计算并对计算结果进行同步存储，并将计算得出的该个体的比冲值作为该个体的适应度值；其中，四个所述能量特性参数分别为该个体的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量；
待所述初始化种群中所有个体的四个所述能量特性参数均计算完毕后，所述初始化种群中所有个体的适应度值均计算出来；之后，再相应计算出所述初始化种群的种群平均适应度值；
步骤303、选择操作：根据步骤302中计算得出的所述初始化种群中所有个体的适应度值，选出所述初始化种群中适应度值高的多个个体作为子代群体；
步骤304、交叉操作与变异操作：对步骤303中选出的子代群体进行交叉操作与变异操作，获得新一代的子代群体，并将所获得的子代种群同步存储至所述数据存储单元内；
步骤305、子代群体中各个体适应度值计算：对步骤304中所获得的子代群体中各个体的适应度值进行计算，且所述子代群体中所有个体的适应度值计算方法均相同；对于所述子代群体中的任一个个体来说，均先利用步骤一中所建立的四个所述能量特性计算模型对该个体的四个能量特性参数进行计算并对计算结果进行同步存储，之后将计算得出的该个体的比冲值作为该个体的适应度值；其中，四个所述能量特性参数分别为该个体的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量；
待所述子代群体中所有个体的四个所述能量特性参数均计算完毕后，所述子代群体中所有个体的适应度值均计算出来；之后，再相应计算出所述子代群体的种群平均适应度值；
步骤306、选择操作：根据步骤305中计算得出的所述子代群体中所有个体的适应度值，选出所述子代群体中适应度值高的多个个体作为子代群体；
步骤307、判断是否满足终止条件：当进化代数超过预先设定的最大进化代数itex或者子代群体中个体的最大适应度值大于或等于预先设定的适应度设定值时，遗传算法终止并输出当前所获得所述子代群体中适应度值最高的个体；否则，返回步骤304；
步骤307中所输出的适应度值最高的个体，便为通过遗传算法所设计的固体推进剂质量配比m1、m2、…、mN；
步骤四、图形绘制用取值范围设定及各原料类别的质量含量计算，其计算过程如下：
步骤401、原料类别质量含量取值范围及能量特性参数取值范围设定：当Q=2时，采用所述参数输入单元输入配制所设计固体推进剂所用粘合剂的质量含量取值范围（m粘0～m粘z）和所用氧化剂的质量含量取值范围（m氧0～m氧z），并同步存储至所述数据存储单元内；当Q=3时，采用所述参数输入单元输入配制所设计固体推进剂所用粘合剂的质量含量取值范围（m'粘0～m'粘z）、所用氧化剂的质量含量取值范围（m'氧0～m'氧z）和所用添加剂的质量含量取值范围（m'添0～m'添z），并同步存储至所述数据存储单元内；同时，采用所述参数输入单元分别输入四个所述能量特性参数的取值范围，并同步存储至所述数据存储单元内；
步骤402、原料类别质量含量计算：所述数据处理器根据步骤二中所输入的各组分的原料类别，调用原料类别质量含量计算模块，对步骤301中所述初始化种群中和步骤304中所获得的子代群体中各个体的原料类别质量含量分别进行计算，并对计算结果进行同步存储，且所有个体的原料类别质量配比计算方法均相同；
对于步骤301中所述初始化种群中和步骤304中所获得的子代群体中任一个个体来说，当Q=2时，根据该个体中所包含的N种组分的质量配比取值，先计算得出该个体中所属原料类别为粘合剂的所有组分的质量配比取值之和m粘，之后根据公式m粘=100‑m氧，计算得出该个体中所属原料类别为氧化剂的所有组分的质量配比取值之和m氧；其中，m粘为该个体的粘结剂的质量含量，且m氧为该个体的氧化剂的质量含量；
对于步骤301中所述初始化种群中和步骤304中所获得的子代群体中任一个个体来说，当Q=3时，根据该个体中所包含的N种组分的质量配比取值，先计算得出该个体中所属原料类别为粘合剂的所有组分的质量配比取值之和m’粘，再计算得出该个体中所属原料类别为氧化剂的所有组分的质量配比取值之和m’氧，之后再根据公式m’添=100‑m’粘‑m’氧，计算得出该个体中所属原料类别为添加剂的所有组分的质量配比取值之和m’添；其中，m’粘为该个体的粘结剂的质量含量，m’氧为该个体的氧化剂的质量含量，且m’添为该个体的添加剂的质量含量；
步骤五、能量特性图形绘制，其绘制过程如下：
步骤501、能量特性参数输入：采用所述参数输入单元输入需通过能量特性图形进行分析的能量特性参数y，其中y为比冲、特征速度、燃烧室温度或燃烧产物平均相对分子质量；
步骤502、能量特性曲线方程拟合：所述数据处理器调用曲线方程拟合模块，根据步骤302中计算得出的所述初始化种群中各个体的四个能量特性参数和步骤305中所述子代群体中各个体的四个能量特性参数，并结合步骤四中计算得出的所述初始化种群中和步骤304中所获得的子代群体中各个体的原料类别质量含量，拟合出能量特性参数y的能量特性曲线方程；
其中，当Q=2时，拟合出的能量特性曲线方程为y=a·x1+b·x2+c·x1·x2+d（4），式（4）中a、b、c和d均为常数，x1为配制所设计固体推进剂所用粘合剂的质量含量，x2为配制所设计固体推进剂所用氧化剂的质量含量，且x1+x2=100；
当Q=3时，拟合出的能量特性曲线方程为y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5），式（5）中a'、b'、c'、d'、e'、f'、g'和h'均为常数，且x'1为配制所设计固体推进剂所用粘合剂的质量含量，x'2为配制所设计固体推进剂所用氧化剂的质量含量，x'3为配制所设计固体推进剂所用添加剂的质量含量，且x'1+x'2+x'3=100；
步骤503、能量特性图形绘制：所述数据处理器调用能量特性图形绘制模块，且根据步骤502中拟合出的能量特性参数y的能量特性曲线方程，对能量特性参数y的能量特性图形进行绘制；
当Q=2时，所述数据处理器调用能量特性图形绘制模块，绘制出能量特性参数y的等高线图，所绘制的等高线图上包括能量特性参数y取不同值时的多条等高线，每一条所述等高线均为由能量特性参数y的取值相等的各点所连成的曲线，且多条所述等高线上均标注有该条等高线所对应的能量特性参数y的取值；实际对能量特性参数y的等高线图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x1，且纵坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x2的平面直角坐标系，再根据步骤401中所输入的（m粘0～m粘z）和（m氧0～m氧z），对所建立平面直角坐标系中横坐标轴和纵坐标轴的刻度分别进行设置；之后，根据公式y=a·x1+b·x2+c·x1·x2+d（4），在所建立的平面直角坐标系上绘制出能量特性参数y取不同值时的多条等高线；
当Q=3时，所述数据处理器调用能量特性图形绘制模块，绘制出能量特性参数y的等性能三角图，所绘制的等性能三角图上包括能量特性参数y取不同值时的多条等高线，每一条所述等高线均为由能量特性参数y的取值相等的各点所连成的曲线，且多条所述等高线上均标注有该条等高线所对应的能量特性参数y的取值；实际对能量特性参数y的等性能三角图进行绘制时，先建立一个等边三角形，并将所建立等边三角形的三条边分别定义为配制所设计固体推进剂所用粘合剂的质量含量x'1、配制所设计固体推进剂所用氧化剂的质量含量x'2和配制所设计固体推进剂所用添加剂的质量含量x'3；之后，根据步骤401中所输入的（m'粘0～m'粘z）、（m'氧0～m'氧z）和（m'添0～m'添z），在所建立三角形的三条边上分别设置刻度；然后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5），在所建立的等边三角形上绘制出能量特性参数y取不同值时的多条等高线；
步骤504、能量特性图形同步显示：通过与所述数据处理器相接的显示单元，对步骤503中绘制出的能量特性图形进行同步显示。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤504中对绘制出的能量特性图形进行同步显示后，还需人为对所显示的能量特性图形进行分析，并根据分析结果对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的等高线图；
实际对能量特性参数y的等高线图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x'1，且纵坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x'2的平面直角坐标系；之后，根据步骤401中所输入的（m'粘0～m'粘z）和（m'氧0～m'氧z），对所建立平面直角坐标系中横坐标轴和纵坐标轴的刻度分别进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，在所建立的平面直角坐标系上绘制出能量特性参数y取不同值时的多条等高线。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤503中当Q=2时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的二维能量特性图；实际对二维能量特性图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x1，且纵坐标轴为所设计固体推进剂的能量特性参数y的平面直角坐标系；之后，再根据步骤401中所输入的（m粘0～m粘z），对所建立平面直角坐标系中横坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立平面直角坐标系中纵坐标轴的刻度进行设置；之后，根据公式y=a·x1+b·x2+c·x1·x2+d（4）和公式x1+x2=100，在所建立的平面直角坐标系上绘制出能量特性参数y的取值随x1的取值变化而变化的曲线；
步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的二维能量特性图；实际对二维能量特性图进行绘制时，先通过参数输入单元将x'1、x'2和x'3三个参数中的一个参数设定为固定值，之后再自x'1、x'2和x'3中另外两个参数中选取一个参数作为横坐标轴，并以所设计固体推进剂的能量特性参数y作为纵坐标轴，建立平面直角坐标系；之后，再根据步骤401中所输入的原料类别质量含量取值范围，对所建立平面直角坐标系中横坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立平面直角坐标系中纵坐标轴的刻度进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，在所建立的平面直角坐标系上绘制出一个参数固定条件下能量特性参数y的取值随另一个参数的取值变化而变化的曲线。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的三维能量特性图；所绘制的三维能量特性图上包括多条能量特性曲线一和多条能量特性曲线二；
实际对三维能量特性图进行绘制时，先建立一个X坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x'1、Y坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x'2且Z坐标轴为所设计固体推进剂的能量特性参数y的三维空间直角坐标系；其次，根据步骤401中所输入的原料类别质量含量取值范围，对所建立三维空间直角坐标系中X坐标轴和Y坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立三维空间直角坐标系中Z坐标轴的刻度进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，沿X坐标轴由左至右绘制出多条x'1取不同值时能量特性参数y的取值随x'2的取值变化而变化的能量特性曲线一，每一条所述能量特性曲线一上x'1的取值均相等；然后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，沿Y坐标轴由前至后绘制出多条x'2取不同值时能量特性参数y的取值随x'1的取值变化而变化的能量特性曲线二，每一条所述能量特性曲线二上x'2的取值均相等。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤一中建立所述比冲计算模型、特征速度计算模型、燃烧室温度计算模型和燃烧产物平均相对分子质量计算模型之前，还需根据最小自由能原理建立最小自由能数学模型；
所建立的最小自由能数学模型为式（2）中：j=1、2、…、A，A为固体推进剂所含化学元素的种类；s=1、2、…、m，m为固体推进剂燃烧后处于化学平衡状态时所含燃烧产物的种数；其中μs为预先通过所述参数输入单元输入的第s种燃烧产物的化学势（KJ/mol），ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg）且ns≥0，asj为1mol第s种燃烧产物中所含第j种化学元素的原子数；bj为1000g固体推进剂中所含第j种化学元素的原子数，πj为拉格朗日乘子；
步骤一中所建立的燃烧产物平均相对分子质量计算模型为（6），式（6）中n为1000g气相燃烧产物的摩尔数，且n=1000g固体推进剂燃烧后处于化学平衡状态时所含的所有气相燃烧产物的摩尔数（mol/Kg）之和；
所建立的燃烧室温度计算模型为绝热温度计算模型，且所建立的绝热温度计算模型为式（3）中Hc1为绝热温度T=T1时的1000g燃烧产物的焓，Hc2为绝热温度T=T2时的1000g燃烧产物的焓，Hc1＜Hp＜Hc2，并且Tc1和Tc2均根据公式其中ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg），Hcs为1mol第s种燃烧产物在绝热温度为T时的焓，Hcs=RT(αs1+αs2T/2+αs3T2/3)+αs4T3/4+αs5T4/5+αs6T5/6），其中，R为通用气体常数（Kg·m/mol·K），T为绝热温度，αs1、αs2、αs3、αs4、αs5和αs6为预先通过所述参数输入单元输入的第s种燃烧产物的热力学函数温度系数；其中Mi为配制固体推进剂所用第i种组分的相对分子质量，Hi为预先通过所述参数输入单元输入的1mol第i种组分的焓，Wi为配制固体推进剂所用第i种组分的质量百分数；
所建立的比冲计算模型为理论比冲计算模型，且所建立的理论比冲计算模型为式（1）中：Isp为理论比冲（N·s/Kg），Tc为燃烧室温度（K），Pe为发动机喷管出口处压强（Pa），Pc为燃烧室内压强（Pa），为燃烧产物平均相对分子质量，R为通用气体常数（Kg·m/mol·K），k为绝热指数且其定压比热与定容比热之比，其中，Pe和Pc为预先通过所述参数输入单元输入的发动机设计参数；相应地，步骤302中对于所述初始化种群中的任一个个体来说，均以利用步骤一中所建立的理论比冲计算模型计算得出的该个体的理论比冲值作为该个体的适应度值；步骤305中对于所述子代群体中的任一个个体来说，均以利用步骤一中所建立的理论比冲计算模型计算得出的该个体的理论比冲值作为该个体的适应度值；
所建立的特征速度计算模型为式（5）中n为1000g气相燃烧产物的摩尔数，R为通用气体常数（Kg·m/mol·K），Tc为燃烧室温度（K），k为绝热指数且其定压比热与定容比热之比，
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤302和步骤305中利用步骤一中所建立的四个所述能量特性计算模型，对任一个个体的四个能量特性参数进行计算时，其计算过程如下：
步骤Ⅰ、平衡组成计算：所述数据处理器调用参数计算模块，根据当前所计算的个体中各组分的化学式和质量配比取值m1x、m2x…mNx，并结合各组分的相对分子质量，对asj和bj进行计算；之后，所述数据处理器结合预先输入的μs和πj，且根据公式计算得出ns，此时便获得当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成；
步骤Ⅱ、燃烧室温度计算：所述数据处理器结合步骤Ⅰ中计算得出的当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成，并根据公式计算得出燃烧室温度Tc；其中，其中Mi为当前所计算个体所对应固体推进剂配方中第i种组分的相对分子质量，Hi为1mol第i种组分的焓，Wi为当前所计算个体所对应固体推进剂配方中第i种组分的质量百分数；
步骤Ⅲ、燃烧产物平均相对分子质量计算：先根据步骤一中推算得出的当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成，计算得出当前所计算个体所对应固体推进剂配方燃烧后所产生的1000g气相燃烧产物的摩尔数n，所述数据处理器公式并结合计算得出当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平均相对分子质量
步骤Ⅳ、理论比冲值计算：所述数据处理器根据公式并结合步骤Ⅱ中计算得出的燃烧室温度Tc、步骤Ⅲ中计算得出的燃烧产物平均相对分子质量以及预先通过所述参数输入单元输入的发动机设计参数Pe和Pc，计算得出当前所计算个体所对应固体推进剂配方的理论比冲值Isp；
步骤Ⅴ、特征速度计算：所述数据处理器根据公式并结合步骤Ⅱ中计算得出的燃烧室温度Tc和步骤Ⅲ中计算得出的1000g气相燃烧产物的摩尔数n，计算得出当前所计算个体所对应固体推进剂配方的特征速度C*。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤504中对绘制出的能量特性图形进行同步显示后，人为对所显示的等高线图或等性能三角图进行分析时，先自所述等高线图或等性能三角图中找出能量特性参数y取值最佳的区域，并根据该区域所对应x1与x2的取值范围或者x'1、x'2和x'3的取值范围，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤307中最大进化代数itex=20～100，步骤301中所述初始化种群中所包含个体的数量为20～30。
上述基于遗传算法与能量特性图形的固体推进剂配方优化方法，其特征是：步骤304中进行交叉操作与变异操作时，按照交叉概率pc进行交叉操作，且按照变异概率pm进行变异操作；其中，交叉概率pc=0.6～0.9，变异概率pm=0.1～0.2；且步骤303选择操作中，采用最优个体保留策略来获得最优个体。
本发明与现有技术相比具有以下优点：
1、方法步骤简单、设计合理且实现方便。
2、投入成本低且使用操作简便，明显简化了固体推进剂的配方研制过程，大大缩短了配方研制周期，大幅减小了配方研制成本。
3、采用遗传算法（Genetic Algorithm）进行固体推进剂配方优化设计，不仅计算过程简单，并且计算量小，计算结果的准确度易于保证，所获得推进剂配方的优化结果准确度较高。遗传算法是一种通过模拟自然进化过程搜索最优解的方法，其优点是擅长全局搜索，它本身并不要求对优化问题的性质作一些深入的数学分析。由于进行配方设计工作时，总是力求从各种可能的方案中选择最优方案，推进剂配方优化设计就是基于这种思想，建立在最优化数学理论和现代计算技术的基础上，运用计算机自动计算和设计得出最优方案。遗传算法通过对解空间进行进化从而求得最优方案，其好处在于即使算法中的某些参数不起作用了，整个算法还是可以正常地工作，也就是说，整体种群的走向是越来越好的。
4、遗传算法优化设计方法设计合理、模拟效果好且操作简便，在怀特(White)的最小自由能原理基础上，将固体推进剂配方的能量特性计算与遗传优化算法、图形绘制功能相结合，能根据推进剂配方组分配比迅速计算出能量特性，并能快速优化设计出最高比冲下的最佳配比，同时能将推进剂配方中各组分对能量特性的贡献大小直观表示出来。因而，本发明能有效解决现有固体推进剂配方设计过程中存在的能量特性试验成本高、周期长、试验量大等缺陷，能大幅度减少固体推进剂配方设计中的试验量，节省大量的人力、物力，缩短推进剂配方的研制周期，进一步提高推进剂的能量性能，可用于指导新型高能固体推进剂的配方研制。
5、能够绘制出多种能量特性图形，包括等性能三角图、等高线图、二维能量特性图和三维能量特性图，通过上述能量特性图形能对配制所设计固体推进剂所用的各组分对各能量特性贡献的大小进行直观、全面表示出来，从而能将从而为固体推进剂的配方设计提供了极大的方便。
6、遗传算法优化设计方法的使用效果好、实用价值高且推广应用前景广泛、适用面广，能简便投入推进剂配方的优化设计过程中，可经济、高效地完成推进剂配方优化设计过程，并且具有优化结果准确度高、优化速度快、能结合“所设计固体推进剂的比冲值随第i种组分的质量配比取值变化的变化曲线”对优化结果进行相应调整等诸多优点，能快速获得固体推进剂具有最高比冲的最佳配比，并且所设计固体推进剂的组分可达3～8种，能有效满足固体推进剂配方的设计要求，实用性非常强。
7、将遗传算法优化方法与能量特性图形绘制进行有效结合，在通过遗传算法对所设计固体推进剂进行快速、准确地最高比冲优化设计的同时，利用遗传算法进行配方优化设计过程中所用的所有个体，绘制出所设计固体推进剂的多种能量特性图形。
8、所绘制的能量特性图形的用途主要体现在以下两个方面：第一方面，在采用遗传算法对某一种固体推进剂进行配方优化设计的同时，可以利用遗传算法优化设计过程中所产生的个体（即多个配方），自动绘制出该固体推进剂的多种能量特性图形，能量特性图形的绘制方法非常简单且易于操作，以便技术人员对将推进剂配方中各组分对各能量特性的贡献大小进行直观了解，同时也可以作为后续研究工作的参考资料，实用价值非常高，做到一举两得，也就是说，所绘制的能量特性图形本身的实用价值非常大；第二方面，对所绘制的能量特性图形进行分析，具体是通过能量特性图形对推进剂配方中各组分对各能量特性的贡献大小进行分析，并根据分析结果且结合具体的设计需求，对上述遗传算法中所设计的推进剂配方进行人为微调，相应对遗传算法中所设计的推进剂配方进行进一步优化，以获得适宜当前设计需求的最佳配方。
综上所述，本发明方法步骤简单、设计合理且实现方便、使用效果好，将遗传算法与能量特性图形相结合进行固体推进剂配方优化设计，能有效解决现有固体推进剂配方设计过程中存在的能量特性试验成本高、周期长、试验量大等缺陷。
下面通过附图和实施例，对本发明的技术方案做进一步的详细描述。
附图说明
图1为本发明的优化方法流程图。
图2为本发明所采用遗传算法的方法流程框图。
图3为本发明所绘制的等比冲三角图。
图4为本发明所绘制的等特征速度三角图。
图5为本发明所绘制的比冲等高线图。
图6为本发明所绘制的特征速度等高线图。
图7为本发明所绘制的二维能量特性总图。
图8为本发明所绘制的比冲三维能量特性图.
具体实施方式
如图1所示的一种基于遗传算法与能量特性图形的固体推进剂配方优化方法，包括以下步骤：
步骤一：建模：根据最小自由能原理（即最小吉布斯自由能原理），建立固体推进剂的四个能量特性计算模型；其中，四个所述能量特性计算模型分别为根据固体推进剂的质量配比对该固体推进剂的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量，分别进行计算的比冲计算模型、特征速度计算模型、燃烧室温度计算模型和燃烧产物平均相对分子质量计算模型。
步骤二、初始参数设定与存储：首先，通过与数据处理器相接的参数输入单元，输入配制所设计固体推进剂所用的组分种数N和各组分的组分信息，并将所输入的信息同步存储至与所述数据处理器相接的数据存储单元内；其中，各组分的组分信息均包括化学式、所属原料类别和质量配比范围（mi0～miz），i=1、2、…、N，N为配制所设计固体推进剂所用的组分数量，mi0为配制所设计固体推进剂所用第i种组分的质量配比下限值且miz为第i种组分的质量配比上限值；其中，mi0和miz均为实数，且0＜mi0＜100，0＜miz＜100，N≥2；其中，所属原料类别为粘合剂、氧化剂或添加剂；待N种组分的组分信息均输入完毕后，再采用所述参数输入单元对所设计固体推进剂的原料类别总数Q进行设定，且Q=2或3；
之后，通过所述参数输入单元自预先建立的燃烧产物数据库中选出所设计固体推进剂燃烧后所产生的所有燃烧产物；所述燃烧产物数据库中存储有多种燃烧产物的属性信息；其中，各燃烧产物的属性信息均包括化学式、相对分子质量和相态，其中相态为气相或凝聚相。
本实施例中，Q=3。
步骤三、采用遗传算法进行配方优化设计：结合图2，所述数据处理器调用遗传算法模块，对所设计固体推进剂的质量配比进行优化设计，其优化设计过程如下：
步骤301、种群初始化：将配制所设计固体推进剂所用的N种组分中各组分的一个质量配比取值mix作为一个个体，并将多个个体集合为一个种群，同时所述种群中的所有个体均进行二进制编码后形成初始化种群，且将所形成的初始化种群同步存储至所述数据存储单元内；其中，每一个个体均包括N个质量配比取值mix，i=1、2…、N，且m1x+m2x+…+mNx=100，m1x、m2x、…、mNx分别为N种组分的质量配比取值；mix为第i种组分的质量配比取值且其为自该组分的质量配比范围（mi0～miz）中随机抽取一个数值；
步骤302、初始化种群中各个体适应度值计算：初始化种群中所有个体的适应度值计算方法均相同；对于所述初始化种群中的任一个个体来说，均先利用步骤一中所建立的四个所述能量特性计算模型对该个体的四个能量特性参数进行计算并对计算结果进行同步存储，并将计算得出的该个体的比冲值作为该个体的适应度值；其中，四个所述能量特性参数分别为该个体的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量；
待所述初始化种群中所有个体的四个所述能量特性参数均计算完毕后，所述初始化种群中所有个体的适应度值均计算出来；之后，再相应计算出所述初始化种群的种群平均适应度值；
步骤303、选择操作：根据步骤302中计算得出的所述初始化种群中所有个体的适应度值，选出所述初始化种群中适应度值高的多个个体作为子代群体；
步骤304、交叉操作与变异操作：对步骤303中选出的子代群体进行交叉操作与变异操作，获得新一代的子代群体，并将所获得的子代种群同步存储至所述数据存储单元内；
步骤305、子代群体中各个体适应度值计算：对步骤304中所获得的子代群体中各个体的适应度值进行计算，且所述子代群体中所有个体的适应度值计算方法均相同；对于所述子代群体中的任一个个体来说，均先利用步骤一中所建立的四个所述能量特性计算模型对该个体的四个能量特性参数进行计算并对计算结果进行同步存储，之后将计算得出的该个体的比冲值作为该个体的适应度值；其中，四个所述能量特性参数分别为该个体的比冲、特征速度、燃烧室温度和燃烧产物平均相对分子质量；
待所述子代群体中所有个体的四个所述能量特性参数均计算完毕后，所述子代群体中所有个体的适应度值均计算出来；之后，再相应计算出所述子代群体的种群平均适应度值；
步骤306、选择操作：根据步骤305中计算得出的所述子代群体中所有个体的适应度值，选出所述子代群体中适应度值高的多个个体作为子代群体；
步骤307、判断是否满足终止条件：当进化代数超过预先设定的最大进化代数itex或者子代群体中个体的最大适应度值大于或等于预先设定的适应度设定值时，遗传算法终止并输出当前所获得所述子代群体中适应度值最高的个体；否则，返回步骤304；
步骤307中所输出的适应度值最高的个体，便为通过遗传算法所设计的固体推进剂质量配比m1、m2、…、mN。
本实施例中，步骤一中建立所述比冲计算模型、特征速度计算模型、燃烧室温度计算模型和燃烧产物平均相对分子质量计算模型之前，还需根据最小自由能原理建立最小自由能数学模型。
根据热力学原理，在高温条件下固体推进剂的燃烧产物可视为理想气体，则整个系统的自由能就等于该系统各组分自由能的总和，己知物质的自由能是压力、温度和浓度的函数，当该体系达到化学平衡时，体系的自由能最小。因此，在一定的压力和温度条件下，求出既能使物系自由能最小又符合质量守恒定律的一组分值，则这组组分值即为该条件下物系的燃烧产物平衡组成。根据系统达到化学平衡时其自由能函数之总和为最小的原理，采用迅速收敛的数学方法，通过重复迭代可解任何复杂系统的化学平衡组成。
本实施例中，所建立的最小自由能数学模型为式（2）中：j=1、2、…、A，A为固体推进剂所含化学元素的种类；s=1、2、…、m，m为固体推进剂燃烧后处于化学平衡状态时所含燃烧产物的种数；其中μs为预先通过所述参数输入单元输入的第s种燃烧产物的化学势（KJ/mol），ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg）且ns≥0，asj为1mol第s种燃烧产物中所含第j种化学元素的原子数；bj为1000g固体推进剂中所含第j种化学元素的原子数，πj为拉格朗日乘子。
实际使用时，最小自由能法就是求出在满足式（2）条件下的一组ns值使得体系自由能最小，其中s=1、2、…、m，且ns≥0，这是一个多元函数的条件极值问题，用拉格朗日乘数法则可求解。
实际对化学平衡组成进行求解时，也可以采用河南科学技术出版社于1999年出版的由田德余、刘剑洪著的《化学推进剂计算能量学》一书中，第六章“能量特性计算的基本方程”中第6.3节“含凝聚相燃烧产物的化学平衡组成”中所记载的化学平衡组成计算方法进行计算。
本实施例中，步骤一中所建立的燃烧产物平均相对分子质量计算模型为式（6）中n为1000g气相燃烧产物的摩尔数，且n=1000g固体推进剂燃烧后处于化学平衡状态时所含的所有气相燃烧产物的摩尔数（mol/Kg）之和。
所建立的燃烧室温度计算模型为绝热温度计算模型，且所建立的绝热温度计算模型为式（3）中Hc1为绝热温度T=T1时的1000g燃烧产物的焓，Hc2为绝热温度T=T2时的1000g燃烧产物的焓，Hc1＜Hp＜Hc2，并且Tc1和Tc2均根据公式其中ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg），Hcs为1mol第s种燃烧产物在绝热温度为T时的焓，Hcs=RT(αs1+αs2T/2+αs3T2/3)+αs4T3/4+αs5T4/5+αs6T5/6），其中，R为通用气体常数（Kg·m/mol·K），T为绝热温度，αs1、αs2、αs3、αs4、αs5和αs6为预先通过所述参数输入单元输入的第s种燃烧产物的热力学函数温度系数；其中Mi为配制固体推进剂所用第i种组分的相对分子质量，Hi为预先通过所述参数输入单元输入的1mol第i种组分的焓，Wi为配制固体推进剂所用第i种组分的质量百分数（即第i种组分的质量配比）。其中，Mi和Hi均预先通过所述参数输入单元输入并存储至所述数据存储单元内。
实际根据公式对Hc1进行计算时，ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg），Hcs为1mol第s种燃烧产物在绝热温度为T1时的焓，并且Hcs=RT(αs1+αs2T/2+αs3T2/3)+αs4T3/4+αs5T4/5+αs6T5/6）中的T=T1。而实际根据公式对Hc2进行计算时，ns为1000g固体推进剂燃烧后处于化学平衡状态时所含第s种燃烧产物的摩尔数（mol/Kg），Hc2为1mol第s种燃烧产物在绝热温度为T2时的焓，并且Hcs=RT(αs1+αs2T/2+αs3T2/3)+αs4T3/4+αs5T4/5+αs6T5/6）中的T=T2。
所建立的比冲计算模型为理论比冲计算模型，且所建立的理论比冲计算模型为式（1）中：Isp为理论比冲（N·s/Kg），Tc为燃烧室温度（K），Pe为发动机喷管出口处压强（Pa），Pc为燃烧室内压强（Pa），为燃烧产物平均相对分子质量，R为通用气体常数（Kg·m/mol·K），k为绝热指数且其定压比热与定容比热之比，其中，Pe和Pc为预先通过所述参数输入单元输入的发动机设计参数；相应地，步骤302中对于所述初始化种群中的任一个个体来说，均以利用步骤一中所建立的理论比冲计算模型计算得出的该个体的理论比冲值作为该个体的适应度值；步骤305中对于所述子代群体中的任一个个体来说，均以利用步骤一中所建立的理论比冲计算模型计算得出的该个体的理论比冲值作为该个体的适应度值。
本实施例中，公式中：R=8.3144Kg·m/mol·K，k=1.1～1.3。
实际使用时，可以根据实际需要对R和k的取值进行相应调整。
步骤一中所建立的特征速度计算模型为式（5）中n为1000g气相燃烧产物的摩尔数，R为通用气体常数（Kg·m/mol·K），Tc为燃烧室温度（K），k为绝热指数且其定压比热与定容比热之比，。
本实施例中，步骤302和步骤305中利用步骤一中所建立的四个所述能量特性计算模型，对任一个个体的四个能量特性参数进行计算时，其计算过程如下：
步骤Ⅰ、平衡组成计算：所述数据处理器调用参数计算模块，根据当前所计算的个体中各组分的化学式和质量配比取值m1x、m2x…mNx，并结合各组分的相对分子质量，对asj和bj进行计算；之后，所述数据处理器结合预先输入的μs和πj，且根据公式计算得出ns，此时便获得当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成；
步骤Ⅱ、燃烧室温度计算：所述数据处理器结合步骤Ⅰ中计算得出的当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成，并根据公式计算得出燃烧室温度Tc；其中，其中Mi为当前所计算个体所对应固体推进剂配方中第i种组分的相对分子质量，Hi为1mol第i种组分的焓，Wi为当前所计算个体所对应固体推进剂配方中第i种组分的质量百分数；
步骤Ⅲ、燃烧产物平均相对分子质量计算：先根据步骤一中推算得出的当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平衡组成，计算得出当前所计算个体所对应固体推进剂配方燃烧后所产生的1000g气相燃烧产物的摩尔数n，所述数据处理器公式并结合计算得出当前所计算个体所对应固体推进剂配方燃烧后的燃烧产物平均相对分子质量；
步骤Ⅳ、理论比冲值计算：所述数据处理器根据公式并结合步骤Ⅱ中计算得出的燃烧室温度Tc、步骤Ⅲ中计算得出的燃烧产物平均相对分子质量以及预先通过所述参数输入单元输入的发动机设计参数Pe和Pc，计算得出当前所计算个体所对应固体推进剂配方的理论比冲值Isp；
步骤Ⅴ、特征速度计算：所述数据处理器根据公式并结合步骤Ⅱ中计算得出的燃烧室温度Tc和步骤Ⅲ中计算得出的1000g气相燃烧产物的摩尔数n，计算得出当前所计算个体所对应固体推进剂配方的特征速度C*。
实际操作过程中，步骤307中最大进化代数itex=20～100。步骤301中所述初始化种群中所包含个体的数量为20～30。步骤304中进行交叉操作与变异操作时，按照交叉概率pc进行交叉操作，且按照变异概率pm进行变异操作；其中，交叉概率pc=0.6～0.9，变异概率pm=0.1～0.2。
本实施例中，步骤307中最大进化代数itex=100，步骤301中所述初始化种群中所包含个体的数量为25；步骤304中交叉概率pc=0.8，变异概率pm=0.2。具体操作过程中，可根据实际需要，对步骤307中最大进化代数itex、步骤301中所述初始化种群中所包含个体的数量、步骤304中的交叉概率pc和变异概率pm的取值进行相应调整。
本实施例中，步骤303选择操作中，采用最优个体保留策略来获得最优个体。最优保留策略是在进行交叉和变异后比较产生的新个体适应度值是增加还是减小，如果产生的新个体适应度值增加就保留该新个体，否则保留原个体。该改进策略可以有效地保证遗传操作所产生的都是优良新个体，确定了进化的方向，避免了个别个体在进化过程中退化现象，增强了算法的收敛性能。
本实施例中，采用本发明所公开的基于遗传算法的固体推进剂配方优化设计方法，对NEPE推进剂配方进行优化设计，优化设计结果见表1。也就是说，所设计的固体推进剂为NEPE推进剂。配制NEPE推进剂用以下6种组分，其中，粘合剂为PEG、NG和BTTN，氧化剂为高氯酸铵（AP）和HMX，添加剂为铝粉，HTPB推进剂的质量配比设置范围、优化结果及真实配方对比表，见表1：
表1HTPB推进剂的质量配比设置范围、优化结果及真实配方对比表
组分PEGAPALHMXNGBTTN质量配比设置范围5‑1010‑2010‑2035‑457‑117‑11优化结果7.115.7319.3340.188.758.90真实配方715174299
该真实配方为国內研制多年的高能固体推进剂配方，按常规的实践经验看，要设计、调试，定型一个实用的配方，需要几年到十几年的时间，并且设计一个实用的配方，也需要很长的时间，而利用本发明所采用的配方优化设计方法只需要几分到几十分钟即可搞定，大大缩短了推进剂配方的研制时间，节省了大量的人力、物力和财力。
步骤四、图形绘制用取值范围设定及各原料类别的质量含量计算，其计算过程如下：
步骤401、原料类别质量含量取值范围及能量特性参数取值范围设定：当Q=2时，采用所述参数输入单元输入配制所设计固体推进剂所用粘合剂的质量含量取值范围（m粘0～m粘z）和所用氧化剂的质量含量取值范围（m氧0～m氧z），并同步存储至所述数据存储单元内；当Q=3时，采用所述参数输入单元输入配制所设计固体推进剂所用粘合剂的质量含量取值范围（m'粘0～m'粘z）、所用氧化剂的质量含量取值范围（m'氧0～m'氧z）和所用添加剂的质量含量取值范围（m'添0～m'添z），并同步存储至所述数据存储单元内。同时，采用所述参数输入单元分别输入四个所述能量特性参数的取值范围，并同步存储至所述数据存储单元内。
步骤402、原料类别质量含量计算：所述数据处理器根据步骤二中所输入的各组分的原料类别，调用原料类别质量含量计算模块，对步骤301中所述初始化种群中和步骤304中所获得的子代群体中各个体的原料类别质量含量分别进行计算，并对计算结果进行同步存储，且所有个体的原料类别质量配比计算方法均相同；
对于步骤301中所述初始化种群中和步骤304中所获得的子代群体中任一个个体来说，当Q=2时，根据该个体中所包含的N种组分的质量配比取值，先计算得出该个体中所属原料类别为粘合剂的所有组分的质量配比取值之和m粘，之后根据公式m粘=100‑m氧，计算得出该个体中所属原料类别为氧化剂的所有组分的质量配比取值之和m氧；其中，m粘为该个体的粘结剂的质量含量，且m氧为该个体的添加剂的质量含量；
对于步骤301中所述初始化种群中和步骤304中所获得的子代群体中任一个个体来说，当Q=3时，根据该个体中所包含的N种组分的质量配比取值，先计算得出该个体中所属原料类别为粘合剂的所有组分的质量配比取值之和m’粘，再计算得出该个体中所属原料类别为氧化剂的所有组分的质量配比取值之和m’氧，之后再根据公式m’添=100‑m’粘‑m’氧，计算得出该个体中所属原料类别为添加剂的所有组分的质量配比取值之和m’添；其中，m’粘为该个体的粘结剂的质量含量，m’氧为该个体的氧化剂的质量含量，且m’添为该个体的氧化剂的质量含量。
步骤五、能量特性图形绘制，其绘制过程如下：
步骤501、能量特性参数输入：采用所述参数输入单元输入需通过能量特性图形进行分析的能量特性参数y，其中y为比冲、特征速度、燃烧室温度或燃烧产物平均相对分子质量。
步骤502、能量特性曲线方程拟合：所述数据处理器调用曲线方程拟合模块，根据步骤302中计算得出的所述初始化种群中各个体的四个能量特性参数和步骤305中所述子代群体中各个体的四个能量特性参数，并结合步骤四中计算得出的所述初始化种群中和步骤304中所获得的子代群体中各个体的原料类别质量含量，拟合出能量特性参数y的能量特性曲线方程；
其中，当Q=2时，拟合出的能量特性曲线方程为y=a·x1+b·x2+c·x1·x2+d（4），式（4）中a、b、c和d均为常数，x1为配制所设计固体推进剂所用粘合剂的质量含量，x2为配制所设计固体推进剂所用氧化剂的质量含量，且x1+x2=100；
当Q=3时，拟合出的能量特性曲线方程为y=a'·x'1+b'·x'2+c'·x'3+d'·x'1?x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5），式（5）中a'、b'、c'、d'、e'、f'、g'和h'均为常数，且x'1为配制所设计固体推进剂所用粘合剂的质量含量，x'2为配制所设计固体推进剂所用氧化剂的质量含量，x'3为配制所设计固体推进剂所用添加剂的质量含量，且x'1+x'2+x'3=100；
步骤503、能量特性图形绘制：所述数据处理器调用能量特性图形绘制模块，且根据步骤502中拟合出的能量特性参数y的能量特性曲线方程，对能量特性参数y的能量特性图形进行绘制。
当Q=2时，所述数据处理器调用能量特性图形绘制模块，绘制出能量特性参数y的等高线图，所绘制的等高线图上包括能量特性参数y取不同值时的多条等高线，每一条所述等高线均为由能量特性参数y的取值相等的各点所连成的曲线，且多条所述等高线上均标注有该条等高线所对应的能量特性参数y的取值；实际对能量特性参数y的等高线图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x1，且纵坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x2的平面直角坐标系，再根据步骤401中所输入的（m粘0～m粘z）和（m氧0～m氧z），对所建立平面直角坐标系中横坐标轴和纵坐标轴的刻度分别进行设置；之后，根据公式y=a·x1+b·x2+c·x1·x2+d（4），在所建立的平面直角坐标系上绘制出能量特性参数y取不同值时的多条等高线。
当Q=3时，所述数据处理器调用能量特性图形绘制模块，绘制出能量特性参数y的等性能三角图，所绘制的等性能三角图上包括能量特性参数y取不同值时的多条等高线，每一条所述等高线均为由能量特性参数y的取值相等的各点所连成的曲线，且多条所述等高线上均标注有该条等高线所对应的能量特性参数y的取值；实际对能量特性参数y的等性能三角图进行绘制时，先建立一个等边三角形，并将所建立等边三角形的三条边分别定义为配制所设计固体推进剂所用粘合剂的质量含量x'1、配制所设计固体推进剂所用氧化剂的质量含量x'2和配制所设计固体推进剂所用添加剂的质量含量x'3；之后，根据步骤401中所输入的（m'粘0～m'粘z）、（m'氧0～m'氧z）和（m'添0～m'添z），在所建立三角形的三条边上分别设置刻度；然后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5），在所建立的等边三角形上绘制出能量特性参数y取不同值时的多条等高线。
本实施例中，Q=3。其中，当能量特性参数y为比冲值Isp时，所绘制出的等性能三角图为等比冲三角图，详见图3；当能量特性参数y为特征速度C*时，所绘制出的等性能三角图为等特征速度三角图，详见图4。同理。也可以绘制出等平均相对分子质量三角图和等燃烧室温度三角图。
步骤504、能量特性图形同步显示：通过与所述数据处理器相接的显示单元，对步骤503中绘制出的能量特性图形进行同步显示。
本实施例中，步骤504中对绘制出的能量特性图形进行同步显示后，还需人为对所显示的能量特性图形进行分析，并根据分析结果对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
本实施例中，步骤504中对绘制出的能量特性图形进行同步显示后，还需人为对所显示的能量特性图形进行分析时，先自所述能量特性图形中找出能量特性参数y取值最佳的区域，并根据该区域所对应x1与x2的取值范围或者x'1、x'2和x'3的取值范围，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
本实施例中，人为对等性能三角图进行分析时，先自等性能三角图中找出能量特性参数y取值最佳的区域，并根据该区域所对应x1与x2的取值范围或者x'1、x'2和x'3的取值范围，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
如图3所示的等比冲三角图中，由于能量特性参数y越大越好，因而分析得出：能量特性参数y取值最佳的区域为比冲值Isp在2635Ns/kg‑2646Ns/kg之间的区域，此区域对应的AP‑HMX氧化剂的质量含量在64%～76%，添加剂的质量含量在14%～24%，其余为粘合剂PEG‑NG‑BTTN。之后，根据上述分析结果，并结合所设计固体推进剂的设计需求，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调，以获得最佳配比。
如图4所示的等特征速度三角图中，由于能量特性参数y越大越好，因而分析得出：能量特性参数y取值最佳的区域为特征速度C*在1600m/s～1613m/s之间的区域，此区域对应的AP‑HMX氧化剂的质量含量在64%～78%，添加剂的质量含量在8%～24%之间，其余为粘合剂PEG‑NG‑BTTN。之后，根据上述分析结果，并结合所设计固体推进剂的设计需求，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调，以获得最佳配比。
同理，等燃烧室温度三角图的分析方法，与等比冲三角图和等特征速度三角图的方法相同。而等平均相对分子质量三角图的分析方法，与等比冲三角图和等特征速度三角图的分析方法区别仅在于，由于燃烧产物平均相对分子质量越小越好，因而能量特性参数y取值最佳的区域为能量特性参数y取值较小的区域。
实际对等三角形图进行绘制时，先建立一个等边三角形，并将所建立等边三角形的三条边分别定义为x'1、x'2和x'3三个参数，且相应在所建立三角形的三条边上分别设置刻度；之后，根据鼠标的点击位置来绘制图形，当用户点击了鼠标的时候，首先取出鼠标所处位置的平面坐标，并判断该位置是否落在所建立等边三角形所在的区域之内：如果点击区域落在等边三角形的区域内，就把鼠标的位置转换为所建立等边三角图上三个边上对应坐标的值（即获得与该位置对应的x'1、x'2和x'3三个参数），然后代入所得的拟合方程（5）中,计算出该位置所对应的能量特性参数y的取值；最后，搜索整个计算空间把相同性能值（即能量特性参数y取值相等）的所有配方绘制出来，绘制出来的图形会形成一个闭合的曲线；同时，还可以根据鼠标所在的位置显示或绘出当前鼠标位置处的能量特性参数y的取值。
本实施例中，步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的等高线图；
实际对能量特性参数y的等高线图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x'1，且纵坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x'2的平面直角坐标系；之后，根据步骤401中所输入的（m'粘0～m'粘z）和（m'氧0～m'氧z），对所建立平面直角坐标系中横坐标轴和纵坐标轴的刻度分别进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，在所建立的平面直角坐标系上绘制出能量特性参数y取不同值时的多条等高线。
本实施例中，当能量特性参数y为比冲值Isp时，所绘制出的等高线图为比冲等高线图，详见图5；当能量特性参数y为特征速度C*时，所绘制出的等性能三角图为特征速度等高线图，详见图6。同理。也可以绘制出平均相对分子质量等高线图和等燃烧室温度等高线图。
本实施例中，步骤504中还需通过与所述数据处理器相接的显示单元，对所绘制出的能量特性参数y的等高线图进行同步显示。并且，还需人为对所显示的等高线图进行分析，并根据分析结果对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。并且，人为对等高线图进行分析时，先自等高线图中找出能量特性参数y取值最佳的区域，并根据该区域所对应x1与x2的取值范围或者x'1、x'2和x'3的取值范围，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
如图5所示的比冲等高线图中，由于能量特性参数y越大越好，因而分析得出：能量特性参数y取值最佳的区域为比冲值Isp在2635Ns/kg～2646Ns/kg之间的区域，此区域对应的AP‑HMX氧化剂的质量含量在65%～76%，粘合剂PEG‑NG‑BTTN的质量含量在11%～22%，其余为添加剂AL粉。之后，根据上述分析结果，并结合所设计固体推进剂的设计需求，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调，以获得最佳配比。
如图6所示的特征速度等高线图中，由于能量特性参数y越大越好，因而分析得出：能量特性参数y取值最佳的区域为特征速度C*在1600m/s～1613m/s之间，此区间对应的AP‑HMX氧化剂的质量含量在64%～76%，添加剂的质量含量在13%～28%之间，其余为粘合剂PEG‑NG‑BTTN。之后，根据上述分析结果，并结合所设计固体推进剂的设计需求，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调，以获得最佳配比。
同理，燃烧室温度等高线图的分析方法，与比冲等高线图和特征速度等高线图的方法相同。而平均相对分子质量等高线图的分析方法，与比冲等高线图和特征速度等高线图的分析方法区别仅在于，由于燃烧产物平均相对分子质量越小越好，因而能量特性参数y取值最佳的区域为能量特性参数y取值较小的区域。
实际对等高线图进行绘制时，先建立一个平面直角坐标系，并将所建立平面直角坐标系的横坐标轴和纵坐标轴分别定义为x'1和x'2两个参数，且相应在横坐标轴和纵坐标轴分别设置刻度；之后，根据鼠标的点击位置来绘制图形，当用户点击了鼠标的时候，先将当前鼠标的位置坐标转换为所建立平面直角坐标系的两个坐标轴上的坐标，也就是拟合方程（5）中的两个输入参数的值，而x'3根据x'1+x'2+x'3=100即可求得，因此获得当前鼠标的位置对应的x'1、x'2和x'3三个参数，然后代入所得的拟合方程（5）中,计算出该位置所对应的能量特性参数y的取值；最后，搜索整个计算空间把相同性能值（即能量特性参数y取值相等）的所有配方绘制出来，绘制出来的图形会形成一个闭合的曲线；同时，还可以根据鼠标所在的位置显示或绘出当前鼠标位置处的能量特性参数y的取值。
当Q=3时，与等性能三角图所不同的是：等性能三角图中有三个输入参数，而等高线图中只有两个参数。其实，二者也并没有本质上的不同，只不过，等高线图的第三个参数（而是通过公式x'1+x'2+x'3=100求得）没有显示出来。
另外，当Q=2时，绘制等高线图的结构与当Q=3时的图形结构相同；且所采用的绘制方法与当Q=3时的绘制方法区别仅在于，当Q=2时只有x1和x2两个参数，并且根据拟合方程（4）进行绘制。
另外，步骤503中当Q=2时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的二维能量特性图；实际对二维能量特性图进行绘制时，先建立一个横坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x1，且纵坐标轴为所设计固体推进剂的能量特性参数y的平面直角坐标系；之后，再根据步骤401中所输入的（m粘0～m粘z），对所建立平面直角坐标系中横坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立平面直角坐标系中纵坐标轴的刻度进行设置；之后，根据公式y=a·x1+b·x2+c·x1·x2+d（4）和公式x1+x2=100，在所建立的平面直角坐标系上绘制出能量特性参数y的取值随x1的取值变化而变化的曲线；
步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的二维能量特性图；实际对二维能量特性图进行绘制时，先通过参数输入单元将x'1、x'2和x'3三个参数中的一个参数设定为固定值，之后再自x'1、x'2和x'3中另外两个参数中选取一个参数作为横坐标轴，并以所设计固体推进剂的能量特性参数y作为纵坐标轴，建立平面直角坐标系；之后，再根据步骤401中所输入的原料类别质量含量取值范围，对所建立平面直角坐标系中横坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立平面直角坐标系中纵坐标轴的刻度进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，在所建立的平面直角坐标系上绘制出一个参数固定条件下能量特性参数y的取值随另一个参数的取值变化而变化的曲线。
本实施例中，如图7所示，绘制出当能量特性参数y为比冲值Isp、特征速度C*、燃烧室温度Tc和燃烧产物平均相对分子质量时的二维能量特性总图。
本实施例中，步骤504中还需通过与所述数据处理器相接的显示单元，对所绘制出的二维能量特性总图进行同步显示。并且，还需人为对所显示的二维能量特性总图进行分析，并根据分析结果对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。并且，人为对二维能量特性总图进行分析时，先自二维能量特性总图中分别找出比冲值Isp、特征速度C*、燃烧室温度Tc和燃烧产物平均相对分子质量的取值最佳的位置点，再根据所找出各位置点所对应x1与x2的取值范围或者x'1、x'2和x'3的取值范围，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
如图7所示的二维能量特性总图中，由于比冲值Isp、特征速度C*和燃烧室温度Tc的取值越大越好，且燃烧产物平均相对分子质量的取值越小越好，因而分析得出：当x'2的取值固定（即氧化剂AP‑HMX的质量含量固定）时，随粘合剂PEG‑NG‑BTTN与添加剂AL粉的质量含量变化，固体推进剂的比冲值Isp、特征速度C*和燃烧室温度Tc的取值越大越好，且燃烧产物平均相对分子质量随之变化，且各自变化规律的曲线一目了然，最高比冲值Isp约为2660Ns/kg，最高特征速度C*约为1620m/s，最高燃烧温度Tc约为4200K。之后，根据上述分析结果，并结合所设计固体推进剂的设计需求，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调，以获得最佳配比。
实际对二维能量特性图进行绘制时，先将x'1、x'2和x'3三个参数中的一个参数设定为固定值，再建立一个平面直角坐标系并将x'1、x'2和x'3中另外两个参数中选取一个参数作为横坐标轴，且以所设计固体推进剂的能量特性参数y作为纵坐标轴，同时对所建立平面直角坐标系中横坐标轴的刻度进行设置；之后，根据鼠标的点击位置来绘制图形，当用户点击了鼠标的时候，先将当前鼠标的位置坐标转换为所建立平面直角坐标系的两个坐标轴上的坐标，也就是拟合方程（5）中的函数值和一个输入参数的值，由于x'1、x'2和x'3三个参数中的一个参数为固定值，因而根据x'1+x'2+x'3=100可以求得第三个参数值，因此获得当前鼠标的位置对应的x'1、x'2和x'3三个参数；最后，利用拟合方程（5）求出不同能量特性参数y所对应的x'1、x'2和x'3三个参数。在二维能量特性图，拟合方程（5）中的参数只有一个。
另外，当Q=2时，绘制二维能量特性图的结构与当Q=3时的图形结构相同；且所采用的绘制方法与当Q=3时的绘制方法区别仅在于，当Q=2时只有x1和x2两个参数，因而无需对参数进行固定，并且根据拟合方程（4）进行绘制。
本实施例中，步骤503中当Q=3时，所述数据处理器还需调用能量特性图形绘制模块，绘制出能量特性参数y的三维能量特性图；所绘制的三维能量特性图上包括多条能量特性曲线一和多条能量特性曲线二；
实际对三维能量特性图进行绘制时，先建立一个X坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x'1、Y坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x'2且Z坐标轴为所设计固体推进剂的能量特性参数y的三维空间直角坐标系；其次，根据步骤401中所输入的原料类别质量含量取值范围，对所建立三维空间直角坐标系中X坐标轴和Y坐标轴的刻度进行设置；并根据步骤401中所输入的四个能量特性参数的取值范围，对所建立三维空间直角坐标系中Z坐标轴的刻度分别进行设置；之后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，沿X坐标轴由左至右绘制出多条x'1取不同值时能量特性参数y的取值随x'2的取值变化而变化的能量特性曲线一，每一条所述能量特性曲线一上x'1的取值均相等；然后，根据公式y=a'·x'1+b'·x'2+c'·x'3+d'·x'1·x'2+e'·x'1·x'3+f'·x'2·x'3+g'·x'1·x'2·x'3+h'（5）和公式x'1+x'2+x'3=100，沿Y坐标轴由前至后绘制出多条x'2取不同值时能量特性参数y的取值随x'1的取值变化而变化的能量特性曲线二，每一条所述能量特性曲线二上x'2的取值均相等。
本实施例中，步骤504中还需通过与所述数据处理器相接的显示单元，对所绘制出的三维能量特性图进行同步显示。并且，还需人为对所显示的三维能量特性图进行分析，并根据分析结果对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。并且，人为对三维能量特性图进行分析时，先自三维能量特性图找出能量特性参数y取值最佳的区域，并根据该区域所对应x'1、x'2和x'3的取值范围，对步骤307中所输出的固体推进剂质量配比m1、m2、…、mN进行微调。
实际对三维能量特性图进行绘制时，先建立一个X坐标轴为配制所设计固体推进剂所用粘合剂的质量含量x'1、Y坐标轴为配制所设计固体推进剂所用氧化剂的质量含量x'2且Z坐标轴为所设计固体推进剂的能量特性参数y的三维空间直角坐标系，并X、Y和Z三个坐标轴的刻度分别进行设置；之后，根据鼠标的点击位置来绘制图形，当用户点击了鼠标的时候，先将当前鼠标的位置坐标转换为所建立三维空间直角坐标系中的三个坐标轴上的坐标，也就是拟合方程（5）中的函数值和两个输入参数的值，柄根据x'1+x'2+x'3=100求得第三个参数值，因此获得当前鼠标的位置对应的x'1、x'2和x'3三个参数；最后，利用拟合方程（5）求出所设定区域内不同能量特性参数y所对应的x'1、x'2和x'3三个参数。不同能量特性参数y的选取，按照一定间隔进行选取。
而当采用本发明所述的配方设计方法对HTPB推进剂配方进行优化设计时，先按照遗传算法进行优化设计且遗传算法的优化设计结果见表2。也就是说，所设计的固体推进剂为HTPB推进剂。配制HTPB推进剂用以下三种组分，其中，粘合剂为丁羟粘合系统(含固化剂)即HTPB，氧化剂为高氯酸铵（AP），添加剂为铝粉，HTPB推进剂的质量配比设置范围、优化结果及真实配方对比表，见表2：
表2HTPB推进剂的质量配比设置范围、优化结果及真实配方对比表
组分HTPBAPAL质量配比设置范围8‑1560‑8012‑18优化结果12.8869.117.74真实配方13.167115
该真实配方为国外巳使用，按常规的实践经验看，要设计、调试，定型一个实用的配方，需要几年到十几年的时间，并且设计一个实用的配方，也需要很长的时间，而利用本发明所采用的配方优化设计方法只需要几分到几十分钟即可搞定，大大缩短了推进剂配方的研制时间，节省了大量的人力、物力和财力。
同时，按照本发明步骤四和步骤五中所述的方法进行能量特性图形绘制。本实施例中，Q=3。其中，当能量特性参数y为比冲值Isp时，所绘制出的三维能量特性图为比冲三维能量特性图，详见图8；当能量特性参数y为特征速度C*时，所绘制出的三维能量特性图为特征速度三维能量特性图。同理。也可以绘制出平均相对分子质量三维能量特性图和燃烧室温度三维能量特性图。如图8所示，比冲三维能量特性图能直观反映出固体推进剂的比冲值Isp与所用氧化剂和粘合剂的定量立体关系。
以上所述，仅是本发明的较佳实施例，并非对本发明作任何限制，凡是根据本发明技术实质对以上实施例所作的任何简单修改、变更以及等效结构变化，均仍属于本发明技术方案的保护范围内。