物质热力学函数温度系数的通用计算方法

刘 洌 梁国柱

(北京航空航天大学 宇航学院，北京 100191)

在火箭发动机的热力特性分析中，需要使用推进剂及其燃烧产物在一定温度下的热力函数(摩尔定压热容、焓、熵)值.由于热力学性质表大多只给出的整百K下的标准热力函数值，当计算某一非整百K温度的热力函数值时，就必须利用2个或若干个相邻温度的标准数据进行插值，导致较大误差产生.因此，1961年，美国的研究人员提出，考虑到物质的摩尔定压热容与焓、熵的相互耦合关系，可用最小二乘法拟合出计算热力函数的多项式温度系数，达到方便求解某一温度范围内，任意温度下热力函数的目的［1］.20世纪七八十年代，利用相似的方法，美国刘易斯研究中心通过计算机程序计算得到常用化学物质的温度系数，并提高了计算热力函数值的精度［2－3］.20 世纪90 年代，我国的科技工作者通过精确计算，获得了Cu及其氧化物、氟化物的温度系数［4］.最近20年来，美国相关研究中心进一步扩充了具有广泛用途(包括火星大气性质分析［5－6］)的求解物质热力函数的温度系数数据库［7－10］.近年来，火箭发动机推进剂配方不断改进，以及类似于光谱测量、原子振动测试等实验技术的发展，不仅使火箭发动机推进剂与燃烧产物中可能存在的物质数量大幅增加，而且通过实验测定得到的热力函数值也更加精确，因此应该对现有物质温度系数表进行进一步的更新和扩充.

本文根据标准热力学数据［11－12］计算现代火箭发动机中常见的135种单质、化合物以及带电离子在300～5 000 K温度范围内的温度系数，并对计算得到的物质摩尔定压热容误差较大的温度系数进行单独修正，为工程计算提供更准确的热力函数温度系数值.

1 计算方法

1.1 热力函数的温度系数表达

根据文献［3］，并考虑工程中实际应用要求，可采用含7个相同温度系数的经验公式计算上述3种热力函数公式，具体公式表达如下.

摩尔定压热容:

焓:

熵:

由热力学微分关系形式可知，上述3种热力函数具有关联性，应作为整体考虑.热力函数都具有相同的温度系数a1～a5，焓、熵具有各自独立的系数a6，a7;JANAF标准热力函数表(NIST-JANAF thermochemical tables)的实验测得值与经验公式计算值之间具有一定的误差.为了减少经验公式计算值与标准值间的相对误差，可以采用最小二乘法进行拟合计算，求出温度系数a1～a7.

实际上，随着温度不断升高，分子振动加剧，摩尔定压热容、焓、熵等热力函数值发生剧烈变化，尤其要指出的是摩尔定压热容，对某些物质，甚至会改变摩尔定压热容以温度为自变量的单调性.因此可人为将温度分为若干区间，对于相同物质的不同的温度区间，采用不同的温度系数来计算热力函数.对大多数物质而言，可分为高、低温2个区间，对于部分固体和液体物质，需要将温度区间更细划分，以满足精度要求.

1.2 温度系数的求解方法

选取温度区间［T0，T1］，［T1，T2］，…，［Tn－1，Tn］，共有对应的 a0，i，a1，i，…，an－1，i(i=1，2，…，7)，7n个待求温度系数.对某种确定的物质而言，需依据不同温度区间计算不同的温度系数.

此外，为保证函数的精确性，要求通过相邻温度区间的不同温度系数在温度交界点计算得到的热力函数值相等，即使分段函数两端取得的热力函数值相等，则必须满足如下3(n－1)个方程:

为方便表达，用x表示式(1)～式(3)中与温度系 数 相 乘 的 各 组 合 量，如 xCp，m，6=0，xH，2=0.5RT2，xS，1=R ln T，本文下标表示热力函数和与之对应的温度系数下标.则在温度分段交接点上，利用上述热力函数值相等的条件，可将方程变为

对式(4)～ 式(6)而言，xCp，m，i，xH，i，xS，i(i=1，2，…，7)均是交接温度下的已知参数.

对于一个确定的j，式(4)～式(6)是关于温度系数 aj－1，1，aj－1，2和 aj－1，3的三元一次方程组.当j=n－1时，求出其左端关于待求温度系数的已知系数矩阵的逆矩阵B(3，3)，将其左乘于方程右端各项，获得 an－2，k(k=1，2，3)的表达式:

由于在各温度区间交接点上，矩阵A内均为与交接温度相关的可计算参数.至此，如果已知式(7)右端中的4n+3个温度系数，则由该方程即可求出左端的3(n－1)个余下的温度系数.因此，下面将讨论由求解7n个待求温度系数变成求解4n+3个温度系数的求解过程.

当 T ∈［Tn－1，Tn］时:

热力函数计算值对温度系数的偏导数为

式中xm，i仅是温度的函数.

由式(7)，当 T ∈ ［T0，T1］∪ … ∪ ［Tn－2，Tn－1］，对应区间的温度系数为 aj－1，i(i=1，2，…，7)时:

热力函数计算值对温度系数的偏导数为

由式(8)、式(9)，可得

根据偏导数关系式(11)、式(12)、式(14)和式(15)，结合式(16)、式(17)，可得到一个4n+3阶的正规线性方程组，用选主元的Doolittle方法，精确地解出各个温度区间内的4n+3个温度系数aj，i，然后根据确定的相互关系式(7)，解出剩余的3(n－1)个未知数.

本方法能使3种热力函数同时得到很好的拟合，又保证了在若干个不同温度区间的交界处计算的函数值相等，误差小，精度高.

1.3 摩尔定压热容的单独修正

用1.2节中方法得到的温度系数可以计算出对应温度下的摩尔定压热容、焓和熵.但根据热力学分析可知，作为以温度为单自变量的函数，物质的焓、熵具有单调性，但对于某些特殊物质，摩尔定压热容的单调性不成立，以致使用该温度系数计算得到的摩尔定压热容值误差相对较大.其次，某些物质的摩尔定压热容与焓、熵数量级相差较大，也会对摩尔定压热容的计算精度造成影响.

用计算得到的温度系数计算热力函数值，与标准值进行比较，可知摩尔定压热容的相对误差均比焓、熵的相对误差大一个数量级以上，特别是在函数的单调性发生变化的时候，误差更为明显.为了保证数据的精确性，本文进一步对摩尔定压热容进行单独修正.

为了保证拟合摩尔定压热容函数的精确性，可重新确定一组新的温度系数来单独计算物质的摩尔定压热容.

有最小，就可以达到修正的目的.

由式(1)可知，存在5n个待求温度系数，根据1.2节中类似方法，通过函数变换，化简为(4n+1)×(4n+1)的标准矩阵方程形式，从而求得待求温度系数值.用新求解的温度系数计算摩尔定压热容值，可进一步提高摩尔定压热容值的数据精度.

2 计算结果与分析

根据135种火箭发动机常见物质在300～5000K上的最新标准热力函数数据，统一将温度划分为低温区间(300～1 000 K)和高温区间(1000～5000K)，采用上述计算方法，用 FORTRAN90语言编写计算热力函数温度系数的程序，计算得到火箭发动机135种常见物质的温度系数值，以及相应的摩尔定压热容、焓、熵的计算值，并将其与已知的热力函数标准值进行比较，获得上述每一种物质在各计算温度下的热力函数相对误差及标准差.

现代的固体火箭发动机常采用硼和铜元素作燃烧添加剂，以B2O2与CuO气体为例，其温度系数计算结果如表 1、表2所示.利用式(1)～式(3)计算在上述各温度区间上B2O2和CuO的热力函数值，热力函数标准值、计算值结果如图1所示.

表1 B2 O2(气)的温度系数表

表2 CuO(气)的温度系数表

将得到的函数的计算值与标准值比较，可得如图2所示的计算区间上各计算温度点的相对误差，并得到表3、表4所示的热力函数计算值与标准值的误差分析表.

表3 B2O2(气)热力函数计算值与标准值误差分析表

表4 CuO(气)热力函数计算值与标准值误差分析表

由图2可知，对B2O2与CuO而言，摩尔定压热容函数值在300，1000以及5000K(2个区间端点以及交接点)上相对误差较大，B2O2在300K时取得最大相对误差4.13%，CuO在5000K时取得最大相对误差3.07%，而在其他温度下误差均在1%以内，2种物质的平均相对误差均在0.5%左右;而就焓、熵而言，计算值和标准值相比，在误差较大的区间端点处相对误差也小于1‰，计算值与标准值几乎完全一致.

对135种物质分析计算得到的误差结果如表5所示，可以看出利用温度系数计算得到的计算值相对误差均较小，完全适用于工程应用.

图1 B2 O2(气)、CuO(气)热力函数计算值与标准值比较

表5 135种物质热力函数计算值与标准值误差分析表

135种物质中的26种物质(氮原子、液体铅、固体硅等)的摩尔定压热容计算值，从平均相对误差而言，误差小于5%，但最大相对误差大于10%，有必要进行进一步修正.

图2 热力函数计算值与标准值相对误差

以氮原子为例，在温度小于1 900 K时，摩尔定压热容数值基本不发生变化，但随着温度的继续升高，摩尔定压热容迅速增大.在这种条件下，若不对计算摩尔定压热容的温度系数进行修正，计算中会产生很大误差，但焓、熵的相对误差仍然保持较低值.修正后与修正前得到的摩尔定压热容值相比，可以得到图3所示结果.利用重新修正的计算摩尔定压热容的温度系数计算摩尔定压热容，与标准摩尔定压热容值相比较，可以得到表6所示结果.

图3 氮原子摩尔定压热容修正前后比较

表6 氮原子摩尔定压热容修正前后误差比较

由表6可知，与修正前相比，整个温度区间上摩尔定压热容的平均相对误差减小了250倍，达到10－5数量级，与焓、熵精度相同，最大相对误差也大大减小至3‱，进一步提高了摩尔定压热容温度系数的精确度与计算上的适应性，完全满足工程计算需要.

当然，经验公式计算的准确性还依赖于热力函数标准数据的准确性.根据不同物质的特有化学性质，对单质、化合物等还可以采用3系数、9系数等各种形式的经验公式进行热力函数值的计算，以便进一步提高精确性.此外，根据物质不同的物理性质，划分不同的温度区间，也可以进一步减小计算误差.

3 结论

为了达到方便准确计算物质热力函数的目的，采用含7个温度系数的多项式作为计算摩尔定压热容、焓、熵3种热力函数的经验公式，得到如下成果:

1)利用最小二乘法，计算了火箭发动机研究中常见的包括单质、化合物以及带电离子在内的135种物质在300～5 000K温度范围内的温度系数.

2)根据得到的温度系数计算了135种物质的摩尔定压热容、焓、熵.与标准值相比，在任意温度区间下，焓与熵的相对误差均小于10－3，除极个别温度外，摩尔定压热容的计算值相对误差均小于1%.

3)对氮原子等26种摩尔定压热容值随温度变化较大的物质进行摩尔定压热容的温度系数的单独修正，修正后摩尔定压热容相对误差接近10－5.

4)建立了更准确的热力函数数据库及温度系数数据库.所获得的温度系数不仅可运用于火箭发动机的热力分析，还可以应用于热能动力、化工冶金等其他学科，具有广泛的应用价值.

References)

［1］ Zelenznik F J，Gordon S.Simultaneous least-squares approximation of a function and its first integrals with application to thermodynamic data［R］.NASA TN D-767，1961

［2］ McBride B J，Heimel S，Ehlers JG，et al.Thermodynamic properties to6000K for210 substances involving the first18 elements［R］.NASA SP-3001，1963

［3］ McBride B J，Gordon S.Computer program for calculation of complex chemical equilibrium compositions rocket performance incident and reflected shocks and Chapman-Jouguet detonations［R］.NASA SP-273，1971

［4］田德余，翁武军，刘振源.热力学函数温度系数的计算［J］.航空动力学报，1991，6(2):101 －104 Tian Deyu，Weng Wujun，Liu Zhenyuan.Calculation of the thermodynamic functions’temperature coefficient［J］.Journal of Aerospace Power，1991，6(2):101 －104(in Chinese)

［5］ Captelli M，Colonna G，Giordano D.High-temperature thermodynamic properties of Mars-atmosphere components［J］.Journal of Spacecraft and Rockets，2005，42(6):980 －989

［6］ Nelson H F.Thermodynamic properties of the low Martian atmosphere［J］.Journal of Thermophysics and Heat Transfer，2005，19(4):592－594

［7］ McBride B J，Gordon S.Computer program for calculating and fitting thermodynamic functions［R］.NASA RP-1271，1992

［8］ Gordon S，McBride B J.Thermodynamic data to 20 000 K for monoatomic gases［R］.NASA TP-1999-0208523，1999

［9］ Zehe M J，Gordon S，McBride B J.CAP:a computer code for generating tabular thermodynamic functions from NASA Lewis coefficients［R］.NASA/TP-2001-210959/REV1，2001

［10］ McBride B J，Zehe M J，Gordon S.NASA glenn coefficients for calculating thermodynamic properties of individual species［R］.NASA TP-2002-211556，2002

［11］ Chase MW，Davies C A，Downey JR.JANAF thermochemical tables［M］.3 rd ed.Washington:American Chemical Society，1985:1－1856

［12］ Chase MW.NIST-JANAF thermochemical tables［M］.4 th ed.Berlin:Springer Verlag，2000:1 －1952