发射药氧平衡对燃烧产物浓度影响计算.doc

中北大学2012届毕业设计说明书目 录1 引言11.1 选题的背景和意义11.2 国内外研究概况11.2.1 国外研究概况11.2.2 国内研究现状21.3 主要研究内容32 发射药燃烧产物的计算原理及方法42.1 发射药的氧系数42.2 发射药的氧平衡42.3 迭代法计算火药燃烧产物浓度52.4最小吉布斯能法72.5 列主元法解线性方程组113 最小自由能法的建模及结果与讨论133.1最小自由能法的建模133.2 氧平衡对火药燃烧产物影响153.3 小结184 结论19附录20参考文献27致谢29 第 30 页 共 32 页1 引言1.1 选题的背景和意义火药是指在适当的外界能量作用下，自身能进行迅速而有规律的燃烧，同时生成大量高温气体的物质。火药的燃烧是一种氧化还原反应，火药自身包含氧化元素与还原元素，从元素组成来看，火药通常由碳（C）、氢（H）、氧（O）、氮（N）四种元素组成，其中碳、氢是可燃元素，氧是助燃氧化剂，氮是载气体。它可以在隔绝外界环境的条件下在适当的外界作用时发生自身的氧化还原反应。人们利用火药燃烧化学反应将化学能转变成热能。为了充分利用此热能，在火药设计中会尽可能增加燃烧放热量及比容。火药的燃烧放热量及比容都与火药组成有很大关系，在不同氧平衡下火药被氧化的程度不同，燃烧产物浓度自然也就不同。由此可见，要研究火药的能量性能和燃烧性能，就需要研究火药燃烧产物的生成规律，通过最小吉布斯自由能法可以求出火药燃烧产物浓度，进而研究不同组分的加入对火药燃烧产物的影响。这不仅可以帮助人们判断火药的实际燃烧程度及对人体的危害和对环境的污染等，对于指导火药配方的理论设计和实验研究也有着非常重要的意义1-3。1.2 国内外研究概况1.2.1 国外研究概况为了了解火药燃烧过程的物理化学本质，从理论上控制和调节火药的燃烧性能参数，许多学者从40年代以后就对火药的燃烧机理进行了研究4。另外国外许多著名的烟火学家也建立了若干烟火药燃烧理论模型，如希洛夫的烟火药燃烧火焰结构物理模型5、埃利弗里曼稳态燃烧模型、希特洛夫斯基模型6以及针对Mg与聚四氟乙烯(PTEE)体系的Ladouceur和Kubota燃烧数值模型7。至于火药燃烧产物浓度的计算，国内外最常采用的都是最小自由能法。另外国外也经常会采用布莱克林法。布莱克林法是从基本方程出发推导出的，能够较好地避免平衡组分法在计算中初始值难以设定及在要求较高精度时代次数高的弊病。最小吉布斯自由能法则是将火药的燃烧产物组成问题转化为以元素守恒原理为约束条件的平衡状态下自由能最小的极值问题。根据热力学原理，在高温条件下，燃烧气体产物可视为理想气体，这时物系的自由能就等于组成该体系各组分的自由能之和，G(n)=Gi(Xi)。一种物质的自由能是温度，压力和浓度的函数，当体系达到化学平衡时，体系的自由能就达到最小。因此在一定温度和压力下，求出既能使体系的自由能为最小，又能满足体系质量守恒的一组组分值，即为该条件下的平衡组分，整个过程的计算采用电子计算机编程，迭代至要求精度为止，求得火药燃烧最优平衡组成。这种方法具有广泛的应用意义8。1.2.2 国内研究现状对于火药燃烧产物浓度的计算目前国内常用的方法有三种：第一种是利用化学平衡原理计算；第二种是最小吉布斯自由能法计算；第三种是采用图解法得到的经验公式。最小吉布斯自由能法是目前用得较多的方法之一。在化学平衡原理方法中，首先假设主要的燃烧产物和主要的解离产物，根据原子守恒与解离平衡常数列等式，将各产物浓度以各元素数表示，通过计算机编程计算，适合含化学元素较少的系统的简化计算8。而随着武器的不断发展，对火药的要求越来越高，火药的组分和组成火药的化学元素越来越多，火药的能量越来越高，各种解离产物随温度升高而迅速增大，接下来的计算会变得非常繁琐。由此可见，化学平衡原理的方法具有明显的局限性。在经验公式法中，工程上通常都是采用直线拟合方法通过变量代换把曲线回归转化为直线回归问题，利用回归分析方法来求得经验公式（回归分析是指应用数学的方法，对大量的观察数据进行处理，从而得出比较符合事物内部规律的数学表达式）。火箭火药燃烧产物成分随温度的变化近似为幂函数，属于一元非线性回归问题。大量的计算表明对于双基火箭火药的燃烧产物组分浓度的计算，本方法的精确度很高。可见当无法求得理论公式或理论公式太繁琐时，经验公式能较客观地反映事物的内在规律性，也能使计算更为简洁，便于应用。但缺点便是适用范围比较狭窄9-10。希洛夫的烟火药燃烧火焰结构物理模型、埃利弗里曼稳态燃烧模型、希特洛夫斯基模型等的研究仅限于局部参数在有限条件下的模拟运算，理论模型大都建立在各种假设和约束条件之上，涉及到的因素较少，不足以指导火药的合成和燃烧性能预测。故而北京理工大学的焦庆忠、焦清介等人研究了神经网络技术，通过模拟人类认知的微结构，避开了燃烧理论的建模与运算，选用多层前馈型神经网络结构,模拟黑火药反应并预测燃烧体系产物种类、数量及特性参数，并对仿真结果与试验结果进行对比研究，误差低于7，实现对期望函数(试验数据)的无限逼近11。对于火药燃烧产物的实际测定，西安近代化学研究所在国内率先开展了比较系统的发射药烟焰特性以及推进剂烟雾测试研究工作，通过气象色谱、红外光谱完成气体成分定量，它利用气体采集装置与密闭爆发器直接相连，电化学传感器实时采样读数的方法，对发射药燃烧产物进行了分析。准确性和重复性都很好12。另外西安近代化学研究所还公开了一种真空条件下火药燃烧特性测试装置，主要包括燃烧室、真空泵、电加热器、热电偶、高速摄像仪及数据采集处理系统。本发明的突出优点是，集成度高，一次测量可获得多种参数，试验费用低13。最小自由能法的应用十分广泛。北京理工大学的陈明华、温玉全等人在对MTV烟火剂的燃烧性能计算中采用最小自由能法计算了不同配比的MT烟火剂燃烧后的温度、反应热、析碳量和氟化镁等参数，确定了燃烧后产生最大温度和最大反应热的配比，并由此对在红外辐射中起重要作用的析碳和氟化镁的量进行了计算和分析14。而在单基发射药氧平衡调解的课题中，中北大学的吕智星、贺增弟等人采用最小自由能法计算了单基发射药的爆温、爆热、火药力和燃气组成，探讨了硝酸铵含量对于发射药炮口焰的影响15。另外，南京理工大学的赵志建、刘明伟等人在对于多基硝铵火药稳态燃烧模型的研究中也采用了最小自由能法，提出了描述多基硝铵火药燃烧过程中个含能组分热分解动力学、主要分解和燃烧产物特性，燃烧比表面变化，火药中的空隙率等因素的数学模型16。1.3 主要研究内容 应用最小自由能法原理，编辑C语言程序，求解火药燃烧产物浓度。最小自由能法计算过程中的初始值的求解应用平衡迭代原理；最小自由能法计算过程中的线性方程组采用列主元法计算。改变火药组分比例，计算改变后的燃烧产物浓度，由此探讨氧平衡对火药燃烧产物浓度的影响。2 发射药燃烧产物的计算原理及方法 火药燃烧产物的组成主要取决于火药自身的元素组成和燃烧时的条件等。火药燃烧时放出大量的热，因而燃烧产物所能达到的温度是很高的。高温时燃烧产物间要发生化学平衡反应和燃烧产物的解离反应。应用化学平衡法求得初始值，结合最小自由能法即可求解火药燃烧产物浓度。计算过程中的线性方程组采用列主元法求解。2.1 发射药的氧系数氧系数即火药中所含的氧化元素物质的量与火药中所含可燃元素完全氧化所需氧的物质的量之比。即 （2.1）式中 氧系数 no,i火药中所含氧化元素的物质的量，mol no,i火药中所含可燃元素完全氧化所需氧化元素的物质的量，mol当=1时，物质中氧化元素物质的量与可燃元素完全氧化所需氧化元素物质的量相等，符合等物质的量规则，此时进行燃烧化学反应所放出的热量最大；当1时，属富氧化合物；当1时，属缺氧物质。在火药设计中，为提高火药的做功能力，除考虑燃烧时放热量大以外，还要考虑比容的大小。2.2 发射药的氧平衡所谓氧平衡是指：火炸药或组分中所含氧的质量与所含可燃元素完全氧化所需氧的质量之差值。计算时常以百分数表示。 （2.2）式中 Ob氧平衡 no,i火药或组分中所含氧化元素原子的物质的量之和，mol/mol nf,i火药或组分中所含可燃元素完全氧化所需氧元素原子的物质的量之和，mol/mol M火药或组分的摩尔质量 当Ob0时，该物质为正氧平衡。发射药成分通常包括硝化纤维素、溶剂、安定剂等，而这些组分大多说为负氧平衡，故而发射药也大多为负氧平衡。2.3 迭代法计算火药燃烧产物浓度在高温下燃烧产物间要发生化学平衡反应和燃烧产物的解离反应。例如，含C、H、O、N四种元素的一般发射药燃烧时的主要产物有CO、CO2、H2、H2O、N2。这些燃烧产物中存在水煤气平衡反应，即CO2+H2CO+H2O这个平衡反应随温度变化而向不同方向移动，温度升高，平衡向右移动，产物之间含量比例发生变化。在高温下，燃烧产物间要发生解离，主要解离产物有OH,H,O2,O,N,NO等。有如下的解离反应：H22HH2OOH+HN22NCO2CO+OCO2CO+1/2O2。 。解离反应是一种吸热反应，所以随着燃烧温度的升高，解离程度愈大8。对于一般能量不高的发射药，燃烧温度不高，解离很少，可以忽略。但对于高能的推进剂，由于燃烧温度在3500K以上，甚至4000K以上，这时，燃烧产物的解离反应很严重，燃烧产物种类明显增多。解离反应一方面使火药的化学潜能释放不完全，另一方面又增加火药的比容，对于高能火药可以增加比容而弥补解离消耗的热能。压力对于燃烧产物组成有重要的影响，压力增加，解离反应就困难，对于发射药在枪炮膛压力条件下燃烧，压力高达数百兆帕时，解离程度不高；推进剂在火箭燃烧室压力条件下燃烧，压力为数十兆帕至数百兆帕，则解离反应很明显，燃烧产物的浓度要受到相互间一系列化学反应平衡的制约。对此，只考虑燃烧产物均为气态时的平衡组成。火药假定化学式为CaHbOcNd，对于一般火药，氧系数小于1，但假设其均能使C原子完全气化。那么根据化学反应原理，可将火药燃烧的化学反应写为： 根据原子守恒定律，可列下列方程组， （2.3） 由于含有CO、CO2、H2、H2O,存在水煤气平衡反应：H2+CO2CO+H2O其平衡常数为 （2.4） 其平衡常数Kw是随温度的改变而改变，在给定温度下可热力学数据手册查得。在第一次近似求解时，先假设次要组分OH,H,O2,O,N,NO为0,。这样可将式（2.3）和式（2.4）简化、合并得 （2.5）由上式可以求出各燃烧产物的质量摩尔浓度。将上式转化为只有的一元二次方程，求解，即 （Kw-1)+Kw(1/2nH-nO+nC)+nO+nC(nc-no)=0 （2.6）由式（2.6）求出，代入式（2.5）得主要燃烧产物质量摩尔浓度的一次近似解 （2.7）式中 （2.8） 求值的过程中，根号前的正负号一般取正，有时也取负，以得正值为合理。由此即可求出5种主要组分的平衡组成。2.4最小吉布斯能法随着火药中的高能添加成分越来越复杂，燃烧产物种类越来越多。人们越来越倾向于应用最小自由能法计算系统化学平衡。一种物质的自由能是温度、压力和浓度的函数。当体系达到化学平衡时，体系的自由能为最小。因此，在一定的温度和压力条件下，求出既能使体系的自由能最小，又能满足体系质量守恒的一组组分值，即为该条件下的平衡组分。这就是最小自由能法计算复杂系统平衡组成的基本原理18。根据热力学原理，在高温条件下，火药的燃烧气体产物可视为理想气体，这时物系的自由能就等于组成该体系各组分的自由能之和。 （2.9）设一个系统由l种化学元素组成，燃烧后该系统生成n种气态产物,则系统的自由能函数之和可表示为： （2.10）对于每一组分的自由能函数为： （2.11） （2.12） （2.13） （2.14）则系统之总自由能函数 （2.15）式中 G(n)系统之总吉布斯自由能函数 第i种气态组分之吉布斯自由能函数 物质的标准吉布斯自由能 第i种气态组分物质的量 n系统组分物质的量之和 p系统的压力 T系统的温度 R摩尔气体常数对于复杂体系的元素原子守恒方程为 （2.16）式中 nj系统中j元素的物质的量 aij 第i种气体组分中j元素原子物质的量以系统自由能函数方程式（2.15）和元素的原子质量守恒方程式（2.16）为基础，进行一系列的数学推导，即可求出火药燃烧产物浓度。迭代计算公式的推导如下：假设一组正值（）作为平衡组分的近似值，则由式（2.10）可知系统的总自由能为 （2.17）式中 （2.18） （2.19）那么可在x=y处用泰勒级数展开式的前三项近似值Q(x)来表示原来的自由能函数 （2.20）由式（2.11）知 （2.21） （2.22）将式（2.21）、式（2.22）代入式（2.20），经整理后得 （2.23）下一步是要求满足元素原子物质的量守恒条件下，求Q(x)是极小值，即转化成了求式（2.23）极小值。求一个函数的条件极值，可用拉格朗日乘数法，即使Q(x)函数与元素原子物质的量守恒方程合并，则 （2.24）式中 j拉格朗日乘数（j=1,2,.l）求出G(x)的最小值即Q(x)的最小值为系统的自由能函数最小。对式（2.25）求导，并令其等于0，即： （2.25）则 （2.26） 由式（2.26）得 （2.27）式中 （2.28）将式（2.27）各项累加得 （2.29） （2.30）式中 （2.31）所以 （2.32）令，则式（2.27）又可写为 （2.33）将式（2.33）代入式（2.16）则 （2.34）令 （2.35）所以式（2.34）可以改写为 （2.36）将方程式（2.32）和（2.36）联立，构成了（n+1+l）个线性方程组 （2.37）解出未知数1、2、j及（1+u），将j值代入式（2.33）得到一组改善了的，若所得一组组分全部为正值，则可直接作为下一次迭代计算的起始值；如果出现负值，使迭代过程无法进行下去，则通常采用减幅因子的办法来进行校正，以保证迭代计算能顺利进行。 （2.38）则 （2.39）Zi为下一次迭代计算的初值。如此反复迭代下去，直至最后两次计算结果达到精度要求且自由能为最小时为止。2.5 列主元法解线性方程组列主元素消去法具有很好的数值稳定性，在实际应用中它与全主元素消去法精度差不多，但计算过程比较简单，故应用的更广泛。对于线性方程组 (2.40)列主元素消去法首先在增广矩阵A,b的第一列n个元素中选取最大的一个作为主元素，并把此主元素所在的行与第一行交换，然后通过初等行变换把第一列的后n-1个元素消为0，得增广矩阵A(1),b(1)。其次在A(1),b(1)的第二列后n-1个元素中选取绝对值最大的一个作为主元素，并把它所在行与第二行交换，然后通过初等行变换，把第二列的后n-2个元素消为0，依此类推。只要|A|0，消元过程就能进行到底18。列主元素消去法的一般计算过程如下：消元过程：对于k=1,2，n-1，执行（1） 选行号ik,使。（2） 交换第ik行与第k行（3） 对于i=k+1,k+2,n,计算 回代过程：3 最小自由能法的建模及结果与讨论3.1最小自由能法的建模假设火药组成为CaHbOcNd,火药燃烧产物有CO2,CO,H2O,H2,N2,H,OH,O2,O,N,NO，由迭代法求出第一组值，应用最小自由能法进行编程计算，当精度小于10-5时停止迭代。 本程序的流程图如下：输入火药假定化学式CaHbOcNd由平衡迭代法求取燃烧产物浓度初始值xi将最小自由能方程展开成泰勒级数，取前三项列原子守恒方程直接带入进行下一次计算以列主元法求解方程组得出新一组产物浓度值yiY若yi此均大于0N采用减幅因子修正，得到一组正值Y若xi,yi相差精度10-5N结束图3.1 最小自由能法的流程图最小吉布斯能法求火药燃烧产物浓度的C语言程序代码详见附录。3.2 SF-2火药燃烧产物计算已知SF-2火药，其中硝化纤维素含氮量为12%，其假定化学式为C25.5243H33.1277O33.0757N8.3393，1kgSF-2火药中各元素物质的量及各组分氧平衡值具体如下表3.1所示，NC表示硝化纤维素,NG表示硝化甘油,DNT表示二硝基甲苯,DNP表示邻苯二甲酸二丁酯,C-2表示二号中定剂，V表示凡士林：表3.1 1kgSF-2火药中各元素物质的量及各组分氧平衡组分名称Ob100i100 C H O Nn(C)iin(C)in(H)iin(H)in(O)iin(O)in(N)iin(N)iNC-38.55622.7412.734429.3316.424836.0820.20488.574.7992NG+3.52513.213.302522.025.505039.639.907513.213.3025DNT-114.2938.433.458732.942.964621.961.976410.980.9882DNP-224.5657.483.448879.044.742414.370.8622C-2-246.7362.421.872666.591.99774.160.12488.320.2496V-347.0170.730.7073149.321.4932总计10025.524333.127733.07579.3395 此时火药氧平衡Ob=-55.304%,根据最小自由能法，求得SF-2火药在3000K,300MPa下燃烧产物浓度，结果如下表3.2所示表3.2 SF-2火药氧平衡及燃烧产物浓度产物名称产物初始浓度第一组迭代值产物最终浓度CO21.6759701.6352191.635178CO23.84833023.88904623.889179H2O5.8754295.9146135.915847H210.68842110.64836010.647875N24.6697504.6695304.669707H0.0000100.0009640.000111OH0.0000100.0000990.000098O20.0000100.0007320.000039O0.0000100.0002010.000054N0.0000100.0000490.000028NO0.0000100.0003980.000074sum46.7579646.75921146.75819现行的火药为负氧平衡，并且其所含的氧不足以将C元素完全氧化成CO2，故而实际上只有少数的C得到了充分的氧化。在计算火药燃烧产物浓度时，假设火药中所含的氧首先将碳原子氧化成CO，然后剩余的氧再均匀地将CO氧化成CO2、将H2氧化成H2O。由表3.2与表3.1数据对比可知燃烧产物符合原子守恒，可见用最小自由能法计算火药燃烧产物浓度是可行的。另外，CO与H2的生成量明显较CO2和H2O多，而并且的H2O生成量略微多于CO2。因为这两个反应的热效应接近，氢对氧的亲和力稍大于碳，故而生成H2O的能力会略微强于生成CO2。由于火药在高温下燃烧，因而产物之间发生解离反应。但由于压力较高，解离程度并不大。3.2 氧平衡对火药燃烧产物影响为了研究火药的氧平衡对火药燃烧产物浓度的影响，现在不改变火药其它组成的基础上改变其中某一组份的含量，从而达到改变火药整体的氧平衡的目的。（一）作为双基药的主要能量组分，硝化甘油的氧平衡为+3.5，很明显增加火药中硝化甘油的含量可以有效提高发射药的氧平衡，在不改变其它组分比例的条件下增加硝化甘油含量，具体如下表3.3所示：表3.3 提高硝化甘油比例后各元素物质含量组分名称i100 C H O Nn(C)iin(C)in(H)iin(H)in(O)iin(O)in(N)iin(N)iNC52.2722.7411.886229.3315.330836.0818.8598.574.4795NG3013.213.963022.026.606039.6311.88913.213.9630DNT8.438.433.228132.942.767021.961.844610.980.9223DNP5.657.483.218979.044.426214.370.8047C-22.862.421.747766.591.86454.160.11658.320.2330V0.9370.730.6578149.321.3887总计10024.701732.383233.51389.5978 改变组分组成后火药假定化学式为C24。7017H32.3832O33.5138N9.5978，根据氧平衡的加和性可得改变组分后其氧平衡为-51.373%。根据最小自由能法，火药在3000K,300MPa下燃烧产物浓度，结果如下表3.4所示：表3.4 提高硝化甘油比例后氧平衡及燃烧产物浓度产物名称产物初始浓度产物最终浓度CO22.0665112.019482CO22.63518922.682178H2O6.7455906.792343H29.4460119.399177N24.7989004.798847H0.0000100.000110OH0.0000100.000099O20.0000100.000044O0.0000100.000057N0.0000100.000028NO0.0000100.000076sum45.69226145.692441表3.4结果表明，增加硝化甘油含量，氧平衡提高。对比表3.2可见，当火药氧平衡增加时，由于CO2和H2O标准生成焓较小，故而CO2和H2O的生成量增加，根据原子守恒原理，CO和H2的生成量则相应降低。因为温度、压力均没有发生变化，火药组成的变化也不大，故而解离产物浓度没有太大变化。C被完全氧化的比例提高，但比容却有所降低。 （二）作为双基药的助溶剂，二硝基甲苯的氧平衡为-114.2，可见增加火药中二硝基甲苯的含量可以显著降低发射药的氧平衡，在不改变其它组分比例的条件下增加二硝基甲苯含量，具体如下表3.5所示：表3.5 提高二硝基甲苯比例后各元素物质含量组分名称i100 C H O Nn(C)iin(C)in(H)iin(H)in(O)iin(O)in(N)iin(N)iNC52.923122.7412.034729.3315.522336.0819.47578.574.5355NG23.626313.213.121022.025.202539.639.363113.213.1210DNT1438.435.380232.944.611621.963.074410.981.5372DNP5.670357.483.259379.044.481814.370.8148C-22.835262.421.769766.591.88794.160.11798.320.2359V0.945170.730.6685149.321.4112总计10026.233433.117332.84599.4296改变组分组成后火药假定化学式为C26.2334H33.1173O32.8459N9.4296，根据氧平衡的加和性可得改变组分后其氧平衡为-58.54%。根据最小吉布斯自由能法可知，火药在3000K,300MPa下燃烧产物浓度，结果如下表3.6所示：表3.6 提高二硝基甲苯比例后氧平衡及燃烧产物浓度产物名称产物初始浓度产物最终浓度CO21.4399541.403724CO24.79344624.829660H2O5.1725475.208431H211.38610211.350094N24.7148004.714759H0.0000100.000111OH0.0000100.000096O20.0000100.000035O0.0000100.000052N0.0000100.000028NO0.0000100.000072sum47.50690947.490966由表3.6计算结果可知，增加二硝基甲苯含量，氧平衡降低。对比表3.2可见，当火药氧平衡降低时，火药的不完全氧化现象表现的更为明显，生成焓较小的CO2和H2O的生成量降低，根据原子守恒原理，CO和H2的生成量则增加。由于N为载气体，只有在下面的解离反应中才会与O由关联H2O+1/2N2H2+NO解离程度不大，因此N2的生成量不会随氧平衡的改变而产生较大浮动。此时完全氧化程度降低，但比容有所增加。3.3 小结火药的氧平衡直接影响着火药的燃烧产物浓度及比容。氧平衡的提高，可以使可燃元素燃烧的更完全，火药的完全氧化程度越大。另外，氧化元素O的摩尔物质的量较可燃元素C、H均大，因此氧平衡的提高会使火药的火药中C、H元素相对减少，生成气体也随之减少，故而比容通常会降低。4 结论本文根据最小吉布斯自由能法建立数学模型求火药燃烧产物浓度。在求解过程中，以平衡迭代法计算初始值，以列主元法计算求解线性方程组。以SF-2火药为例，计算其燃烧产物浓度，并改变SF-2火药中的各组分比例，使其氧平衡发生变化，前后的计算结果对比表明：当改变火药组分比使火药氧平衡提高时，标准生成焓较大的CO、H2更易充分燃烧生成标准生成焓较小的CO2和H2O，致使CO、H2含量相对降低，CO2、H2O含量相对增加；相应地，改变火药组分比使氧平衡降低则会使产物发生与之相反的变化。因此在火药设计中，要综合考虑武器对火药的性能要求，合理选取氧平衡，调节放热量与比容，减少有毒有害气体的产生，保障工作人员的安全。由平衡迭代法求解的初始值与燃烧产物实际浓度相差不多，因而减少了后续工作的计算量。从上面例子可以看出, 利用最小自由能法可以计算火药燃烧产物浓度，进而为研究火药燃烧性能与能量性能奠定基础，使火药的配方设计更具针对性。 但是, 该计算的不足之处在于只是考虑了理想状态下的气体产物。 在实际计算中还应考虑固态产物和实际气体, 计算量会加大。随着计算机技术的发展, 这个问题会得到很好的解决。 附录C语言程序代码#include #include #define n 5#define P 300000#define T 3500#define R 8.314#define K 7.822float ini(float a,float b,float c,float d,int i)float x5,A;A=(K*(a+b/2-c)+c)/(2*(K-1);x0=-A+sqrt(A*A+a*(c-a)/(K-1);x1=a-x0;x2=c-a-x0;x3=b/2-x2;x4=d/2; return(xi);float cal(float Bnn+1,int t)float m,p,x5;int i,j,k,q;for(i=0;i=n-2;i+)m=fabs(Bii);j=i;for(k=i+1;km)m=fabs(Bki);j=k;for(q=i;q=n;q+)p=Biq;Biq=Bjq;Bjq=p;for(k=i+1;k=n-1;k+)m=-Bki/Bii;for(j=i;j=0;i-)p=0.0;for(j=i+1;j=n-1;j+)p=p+Bij*xj;xi=(Bin-p)/Bii; return xt;void main()float temp,min,di,x11,y11,C11,z11,s5=0.0,B56=0,Gm11=-1412297.8,-975783.55,-1119951.55,-622513.9,-842503.9,-295786.7,-771454.8,-898911.4,-428669.8,-175946.3,-824347.15,G11,sum5=0.0,PI5,a,b,c,d,yave=0.0,u,subi,v,t11=0.0;double A114=1,0,2,0,1,0,1,0,0,2,1,0,0,2,0,0,0,0,0,2,0,1,0,0,0,1,1,0,0,0,2,0,0,0,1,0,0,0,0,1,0,0,1,1,r44=0;int i,j,k;printf(input:a,b,c,dn);scanf(%f,%f,%f,%f,&a,&b,&c,&d);printf(燃烧产物初始浓度为:n);for(k=0;k=4;k+)xk=ini(a,b,c,d,k);for(i=5;i=10;i+)xi=1e-5;do for(i=0;i=10;i+)yi=xi;for(i=0;i=10;i+) yave+=yi;for(i=0;i=10;i+) Ci=(Gmi)/(R*T)+log(P/101.325);for(i=0;i=10;i+)Gi=yi*(Ci+log(yi)-log(yave);for(j=0;j=3;j+)for(i=0;i=10;i+) sumj+=Aij*Gi;for(i=0;i=10;i+)sum4+=Gi;for(j=0;j=3;j+)for(k=j;k=3;k+) for(i=0;i=10;i+) rjk+=Aij*Aik*yi; rkj=rjk;for(j=0;j=3;j+)for(i=0;i=10;i+)sj+=Aij*yi;for(i=0;i=3;i+)for(j=0;j=3;j+)Bij=rij;for(i=0;i=3;i+)Bi4=si;for(j=0;j=3;j+)B4j=sj;for(i=0;i=4;i+)Bi5=sumi;for(j=0;j=4;j+)PIj=cal(B,j);for(i=0;i=10;i+)for(j=0;j=3;j+)ti+=Aij*PIj; for(i=0;i=10;i+)xi=yi*(1+PI4)+ti*yi-Gi;for(i=0;i=10;i+)printf(x%d=%fn,i,xi); float u11=0; j=0;for(i=0;i=10;i+)if(xi0) temp=yi/(yi-xi); uj+=temp; v=u0; for(i=0;i=10;i+) if(uiv) v=ui; for(i=0;i=10;i+) subi=xi-yi; zi=yi+subi*v; xi=zi;float oo11;for(i=0;i=10;i+)ooi=(xi-yi)/yi;min=oo0;for(i=0;iooi)min=ooi;while(min=1e-5);参考文献1吕秉峰,刘幼平,董凤云等.定容条件下火药实际燃烧规律的数值模拟.火炸药学报,2007,12(6)2徐复明.低敏感高能发射药新配方、装药、点火和理论模拟技术.火炸药学报,2003,26(4)3郑东升，梁磊.提高发射药氧平衡的途径.山西化工,2009,29(4)4黄人骏.火药设计基础.北京：北京理工大学出版社.19975田丰,杨慧群,张俊.密闭爆发器实验中热散失修正方法研究.山西化工.2011，31（1）6Hagan M T，Demuth H BNeural network design.PWS Publishing Company.2006,47Conkling J AChemistry of pyrotechnicsMarcel Dekker，Inc2005,68刘继华.火炸药物理化学性能.北京：北京理工大学出版社.19979万俊华,钟宏伟.