
2.2 常用的数学建模方法
2.2.1 机理分析法
应用自然科学中的定理和定律,对被研究系统的有关因素进行分析、演绎、归纳,从而建立系统的数学模型。机理分析法是人们在一切科学研究中广泛使用的方法。
【例2-1】 在渗碳工艺过程中通过平衡理论找出控制参量与炉气碳势之间的机理关系式。甲醇加煤油渗碳气氛中,描述炉气碳势与CO2含量关系的实际数据,见表2-1。
表2-1 甲醇加煤油渗碳气氛(930℃)

【解】 渗碳过程中的炉气化学反应式如下:
CFe+CO22CO (2-1)
由式(2-1)可得:
(2-2)
其中,p为总压,设p=1atm(101.325kPa),pCO、分别为CO、CO2气体的分压,φCO、
分别为CO、CO2所占的体积百分数。K2为平衡常数,αC为碳的活度。
(2-3)
(2-4)
其中,CC表示平衡碳浓度,即炉气碳势。CCA表示加热到温度T时奥氏体中的饱和碳浓度。
同样,可得:
(2-5)
对式(2-5)取对数,可得:
(2-6)
由于在温度一定时,CCA和K2均为常数,式(2-6)右边前两项也应为常数。因此,可设lgCCA-lgK2=a。而对于这项,由于φCO、
与CC有关,且要建立CC和
之间的数学模型,于是有
(2-7)
设lgCC=Y,lg=x,可得:
Y=a+bx (2-8)
以上就是利用试验数据进行最小二乘法拟合,拟合过程中的数据见表2-2。
表2-2 碳势控制单参数数学模型最小二乘法拟合过程中的各计算值

求出a=-0.2336,b=-0.41077,于是方程为Y=-0.2336-0.41077x。即
CC=0.5839()-0.41077 (2-9)
式(2-9)即为碳势控制的单参数数学模型。
2.2.2 模拟方法
模型的结构及性质已知,但其数量描述及求解都相当麻烦。如果有另一种系统,结构和性质与其相同,而且构造出的模型也类似,就可以把后一种模型看成是原来模型的模拟,对后一个模型去分析或实验并求得其结果。
例如,研究钢铁材料中裂纹尖端在外载荷作用下的应力、应变分布,可以通过弹塑性力学及断裂力学知识进行分析计算,但求解非常麻烦。此时可以借助试验光测力学的手段来完成分析。首先,根据一定比例,采用模具将环氧树脂制成具有同样结构的模型,并根据钢铁材料中裂纹形式在环氧树脂模型中加工出裂纹。随后,将环氧树脂模型放入恒温箱内,对环氧树脂模型在冻结应力的温度下加载,并在载荷不变的条件下缓缓冷却到室温卸载。已冻结应力的环氧树脂模型放在平面偏振光场或圆偏振光场下观察,环氧树脂模型中将出现一定分布的条纹,这些条纹反映了模型在受载时的应力、应变情况,用照相法将条纹记录下来并确定条纹级数,再根据条纹级数计算应力。最后,根据相似原理、材料等因素确定一定的比例系数,将计算出的应力换算成钢铁材料中的应力,从而获得了裂纹尖端的应力、应变分布。
以上是用试验模型来模拟理论模型,分析时也可用相对简单理论模型来模拟、分析较复杂理论模型,或用可求解的理论模型来分析尚不可求解的理论模型。
【例2-2】 经试验获得低碳钢的屈服点σs与晶粒直径d的对应关系见表2-3,用最小二乘法建立起d与σs之间关系的数学模型(霍尔-配奇公式)。
表2-3 低碳钢屈服点与晶粒直径的对应关系

【解】 以作为x,σs作为y,取y=a+bx,为一直线。设试验数据点为
,一般来说,直线并不通过其中任一试验数据点。因为每点均有偶然误差ei,有
e1=a+bXi-Yi
所有试验数据点误差的平方和为
(2-10)
按照最小二乘法原理,误差平方和最小的直线为最佳直线,求最小值的条件是
(2-11)
得
(2-12)
将计算结果代入式(2-12)联立解得:
(2-13)
取a=σ0,b=K,得到公式
σ=σ0=K=64.09+393.69
(2-14)
这是典型的霍尔-配奇公式。
以上是用试验模型来模拟理论模型,分析时也可用相对简单的理论模型来模拟、分析较复杂的理论模型,或用可求解的理论模型来分析尚不可求解的理论模型。
2.2.3 类比分析法
若两个不同的系统可以用同一形式的数学模型来描述,则此两个系统就可以互相类比。类比分析法是根据两个(或两类)系统某些属性或关系的相似,去猜想两者的其他属性或关系也可能相似的一种方法。
【例2-3】 在聚合物的结晶过程中,结晶度随时间的延续不断增加,最后趋于该结晶条件下的极限结晶度。现期望在理论上描述这一动力学过程(即推导Avrami方程)。
【解】 采用类比分析法。聚合物的结晶过程包括成核和晶体生长两个阶段,这与下雨时雨滴落在水面上生成一个个圆形水波向外扩展的情形相类似,因此可通过水波扩散模型来推导聚合物结晶时的结晶度与时间的关系。
在水面上任选一参考点,根据概率分析,在时间0-t时刻范围内通过该点的水波数为m的概率P(m)为Poisson分布(假设落下的雨滴数大于m,t时刻通过任意点p的水波数的平均值为量E)。
P(m)=e-E(m=0,1,2,3,…) (2-15)
显然有
(2-16)
(2-17)
把水波扩散模型作为结晶前期的模拟来讨论薄层熔体形成“二维球晶”的情况。雨滴接触水面相当于形成晶核,水波相当于二维球晶的生长表面,当m=0时,意味着所有的球晶面都不经过p点,即p点仍处于非晶态。根据式(2-15)可知其概率为
P(0)=e-E (2-18)
设此时球晶部分占有的体积分数为φc,则有
1-φc=P(0)=e-E (2-19)
下面求平均值E,它应为时间的函数。先考虑一次性同时成核的情况,对应所有雨滴同时落入水面,到t时刻,雨滴所产生的水波都将通过p点(见图2-1),把这个面积称为有效面积,通过p点的水波数等于这个有效面积内落入的雨滴数。设单位面积内的平均雨滴数为N,当时间由t增加到t+dt时,有效面积的增量即图中阴影部分的面积为2πrdr,平均值E的增量为
dE=N2πrdr (2-20)

图2-1 有效面积示意
若水波前进速度即球晶径向生长速度为v,则r=vt,对式(2-20)作积分得平均值相同t的关系为
(2-21)
代入式(2-19),得
1-φc= (2-22)
式(2-22)表示晶核密度为N,一次性成核时体系中的非晶部分与时间的关系。
如果晶核是不断形成的,相当于不断下雨的情况,设单位时间内单位面积上平均产生的晶核数即晶核生成速度为I,到t时刻产生的晶核数(相当于生成的水波)则为h。时间增加dt,有效面积的增量仍为2πrdr,其中,只有满足t>r/v的条件下产生的水波才是有效的,因此有
(2-23)
积分得
(2-24)
代入可得
(2-25)
同样的方法可以用来处理三维晶球,这时把圆环确定的有效面积增量用球壳确定的有效体积增量来代替,对于同时成核体系(N为单位体积的晶核数),则
(2-26)
对于不断成核体系,定义I为单位时间、单位体积中产生的晶核数,则
(2-27)
将上述情况归纳起来,可用一个通式表示:
1-φc= (2-28)
式中,k是同核密度及晶体一维生长速度有关的常数,称为结晶速度倍数;n是与成核方式及核结晶生长方式有关的常数。该式称为Avrami方程。下面对所建模型进行检验。
图2-2为尼龙1010等温结晶体数据的Avrami处理结果,可见在结晶前期试验同理论相符,在结晶的最后部分同理论发生了偏离。

图2-2 尼龙1010等温结晶的Avrami作图
分析Avrami方程的推导过程,这种后期的偏离是可以理解的,因为生长着的球晶面相互接触后,接触区的增长即刻停止。在结晶前期球晶尺寸较小,非晶部分很多,球晶之间不发生接触,可以由式(2-23)来描述,随着时间的延长,球晶增长到满足相互接触的体积时,总体的结晶速度就要降低,Avrami方程将出现偏差。
2.2.4 数据分析法
当系统的结构性质不大清楚,无法从理论分析中得到系统的规律,也不便于类比分析,但有若干能表征系统规律、描述系统状态的数据可利用时,就可以通过描述系统功能的数据分析来连接系统的结构模型。回归分析是处理这类问题的有力工具。
求一条通过或接近两组数据点的曲线,这一过程叫曲线拟合,而表示曲线的数学式称为回归方程。求系统回归方程的一般方法如下。
【例2-4】 设有一未知系统,已测得该系统有n个输入、输出数据点,为(xi,yi),i=1,2,3,…,n。现寻求其函数关系Y=f(x)或F(x,y)=0。
【解】 无论x,y为什么函数关系,假设用一多项式
(2-29)
作为对输出(观测值)y的估计(用表示)。若能确定其阶数及系数b0,b1,…,bm,则所得到的就是回归方程-数学模型。各项系数即回归系数。
当输入为xi,输出为yi时,多项式拟合曲线相应于xi的估计值为
(2-30)
现在要使多项式估计值与观测值yi的差的平方和
(2-31)
为最小,这就是最小二乘法,令
(2-32)
得到下列正规方程组:
(2-33)
一般数据点个数n大于多项式阶数m,m取决于残差的大小,这样从式(2-33)可求出回归系数b0,b1,…,bm,从而建立回归方程数学模型。