制陶材料强度设计.
制陶材料优化设计
摘要:本文主要解决了要使硅酸盐(SiN)制陶材料烧制出来的陶瓷制品的强34
度能达到最大值,应如何选择适当的工艺条件的问题。首先,我们建立了一个线性回归模型,利用残差来检验可靠性,剔除了一个残差点。然后利用剩下的合理
用该模型算出了强度的最大的数据重新构建了一个较为精确的线性回归模型,并
值。为了减小误差,增加精度,我们在线性模型的基础上进行了非线性修正——为每一个变量增加了一个平方项,使得精度与效率达到一种合理的、符合实际的平衡状态。解出该模型后,通过残差的标准差和相关系数对两个模型的结果进行了分析和比较。我们还利用减少因素的
提出了更合理的试验设计
。
问题回答:
四种烧
烧结温烧结YO:MgO YO:AlO232323 加热结添 Ca0含量 强度最 度 时间
加剂的(mol%) 大值 0(C) (h)
总量
线14
一步 2.0 1:1 1:4 1900 1 1164.6 性 mol%
非
14
线一步 2.0 1:6 1:4 1900 1 1252.7
mol%
性
关键词:陶瓷强度 多元线性回归 非线性修正
残差的标准差 假设检验
制陶材料优化设计
一(问题重述:
陶瓷制品是一种耐磨、抗氧化和耐高温的材料,能够用于许多有这些特殊需要的领域。但是,陶瓷制品的低强度严重阻碍了它自身的广泛应用。随着人类科技水平的不断发展和社会对高强度、耐磨、抗氧化和耐高温的材料需求量的不断增大,硅酸盐(SiN)制陶材料应运而生。 34
硅酸盐(SiN)制陶材料是一种强度高、耐磨、抗氧化和耐高温的材料,34
它广泛应用于高温结构的材料中,如切割工具、齿轮、内燃机部件及航空、航天飞行器的有关部件等。
影响这种材料强度的因素有:
A: 加热方案。A=两步,A=一步; 12
B: 四种烧结添加剂CaO,YO,MgO和AlO的总量。B=14摩尔%,23231
B=16摩尔%,B=18摩尔%。 23
C: CaO的含量。C=0.0摩尔%,C=1.0摩尔%,C=2.0摩尔%。 123
D: YO的摩尔%与MgO的摩尔的比率。D=1:1, D=1:2, D=1:6。 23123
E: YO的摩尔%与AlO的摩尔%的比率。E=2:1, E=1:1, E=1:4。 2323123
oooF: 烧结温度。F=1800C, F=1850C, F=1900C。 123
G: 烧结时间。G=1h, G=2h, G=3h。 123
为了使强度达到最高工艺水平,经过试验,测量数据见附录一: 要求:
1. 根据附录一表中的测量数据,建立合理的数学模型,并对试验结果
进行分析;
2. 寻找使得强度最大的最优工艺条件;
3. 对所建立的模型进行误差分析并做出
;
4. 提出一种较有效的试验计划及试验结果的分析方法,
就此研究向有关部门写一份申报科技进步奖的报告。
2
二( 基本假设
1. 除了表中所列因素外,不考虑其他因素对陶瓷强度的影响;
2. 表中的各种因素水平只能取题中已经给出的值;
3. 在相同因素下所得到的强度值服从正态分布;
4. 各因素对强度的影响是相互独立的。
三( 符号说明
ˆY 求出的陶瓷强度
Yi(1,2,,7), 给出的陶瓷强度均值 i
ˆ 残差 YY,ii
X(i=1,2,…,m) 与因素相对应的自变量 i
him(0,1,), 自变量的系数 i
, 随机误差
K 残差标准差
回归平方和 SR
残差平方和 Se
H 检验假设的原假设 0
H 检验假设的备择假设 1
α 显著性水平
四( 问题分析与模型建立
4.1 问题分析
“强度达到最高工艺水平”是以硅酸盐(SiN)制陶材料的强度达到最大34值为标准。由题目中所给出的列表可知,该材料的强度受以下七个因素所制约:
3
?加热方案;?四种烧结添加剂CaO,YO,MgO和AlO的总量;?CaO的2323
含量;? YO的摩尔%与MgO%的摩尔的比率;?YO的摩尔%与AlO的摩232323尔%的比率;?烧结温度;?烧结时间。这些因素对强度的影响程度各不相同。
为了简便,我们先用多元线性回归来建立模型,以强度为应变量,七个因素为七个自变量,算出各个因素对强度影响的加权系数,进而求出强度关于这七个因素的线性 关系式。
另外,由于仪器误差、环境影响等不确定因素,可能会导致列表中给出的数据不够准确,使得原本可以被线性化的模型变得不满足线性化条件。残差能反映数据的合理性及偏差度,所以在初步求出关系式后,若是能通过线性假设检验,为了保证结果的精确性,必须对各个因素的残差进行检验,剔除残差极大的样本点;若是不能通过检验,仍有可能是因为坏点的影响使其不能线性化了,所以需检验残差,剔除坏点,并对剔除坏点后的线性函数关系式做线性假设检验。
对剔除坏点后的线性函数关系式进行线性假设检验:如果检验结果在置信度为95%的置信区间,则说明此线性函数关系式成立;否则需要试用非线性回归或是更高级的算法。
线性假设检验后用得出的线性函数关系式求出强度的最大值。 通过
4.2模型建立
非确定性关系是指变量之间存在某种相关关系,但不存在确定的函数关系。而线性回归就是寻找这类数量关系式并进行统计推断的一种方法。题中给出表格中所反映的强度与七个因素之间的关系就是这种非确定性关系,所以可先建立多
ˆY元线性回归模型,设强度为,第i个因素为X得: i
,,,h1,,,
,,,h2ˆ,,,,YXXX,h,,,,n012,, ? ,,,
,,,h,,n,,,
,2N(0,),,,
其中是不可观测的随机误差,可以忽略,而h和是未知参数。所以实,,i
4
际使用的方程可简化为:
h,,1,,h2ˆ,, ? Y,XXXh,,,n120,,
,,hn,,
上式中各个系数及常数项都是利用最小二乘法拟合得到的。在得到这个模
ˆ型后,我们开始对数据进行残差检验。设18个样本点的残差为,残差标准YY,ii差为K:
18112ˆ2KYY,,() ? ,,,iin,1i,1
如果某个样本点的残差超过了残差标准差的两倍,而且在正态概率纸上该样本点也明显偏低,则舍去该样本点。去除该点后,我们用剩余的样本点重新做最小二乘拟合,得到经过修正后的模型。之后,对该模型采用逐步回归的计算方法求系数。
接下来,我们对该关系式进行线性假设检验。为检验应变量Y和自变量Xi之间是否存在线性关系,检验假设如下:
Hhhhh:0,,,,,,0123m ? ,Hhim:0,1,2,3,,至少有一个,,1i,
若拒绝H,则说明应变量与自变量之间存在线性关系,线性回归法可行;若拒0
绝H,则说明应变量与自变量之间不存在线性关系,必须放弃线性回归法。 1
在m元线性回归模型中,当?中的假设H为真时,可得0
SS与相互独立(): Re
SmRFFmnm,,,(,1) ? (1)Snm,,e
Fmnm(,1),,对于给定的显著水平(取0.05)可查表得。对于给出回归观,,
测值可算得数值F,若
FFmnm,,,(,1) ? ,
ˆY则拒绝H,即认为应变量与自变量X之间确有线性关系;否则就接受H,即0i1
5
认为应变量与自变量之间没有线性关系。
五( 模型求解与结果分析
5.1 模型求解
在Matlab中,先用regstats(y,x)函数第一次求出?式:
67.673,,
,,,58.058,,
,,52.476
,,ˆ ? YXXXXXXX,,,1.4812563.74,,1234567,,
,,111.06,,,5.778,,
,,,33.346,,
接下来,将十八个样本点中的自变量的取值一一代入?式,算出一组强度
ˆ值,利用前面提到的残差的计算公式求出每一样本点的残差值: Yi
k=54.025 k= -2.9903 110
k=-63.419 k=26.718 211
k=84.345 k=-65.798 312
k=119.5 k=29.443 413
k=-4.1604 k=-143.89 514
k=-76.22 k=58.878 615
k=25.894 k=211.92 716
k=-30.691 k=-3.4617 817
k=-109.28 k=-110.83 918
18112ˆ2KYY,,()=89.1175 ,,,iin,1i,1
因为k>2K,所以第十六个样本点是应该删除的残差点Matlab模拟图如下: 16
6
将该点删除后,把剩下的十七个样本点的内容输入Matlab,仍然用最小二乘法拟合得到下面的式子:
95.167,,
,,,82.039,,
,,66.97
,,ˆYXXXXXXX,,16.442443.75 ? ,,1234567,,
,,124.36,,0.80486,,
,,,13.182,,
接下来,对?式进行线性假设检验,以判断数据是否适用于线性模型。
由已知条件和?和?式可求得:
1725ˆSYY,,,()2.3291*10 ,Rii,1
172ˆSYY,,,()52680, ,eiii,1
进而算出
7
SmS7RR FF,,,4.2117(7,9) SnmS(1)9,,ee
FF,,2.710.05
可知应变量与自变量有线性关系,多元线性回归模型可用。
最后通过仔细观察?式的线性特征,容易得到当系数为负值的变量取其最小值,系数为正的变量取其最大值时,强度达到最大。将最佳工艺条件符号(2,1,3,3,3,3,1)代入得:
ˆ=1164.59 Ymax
即硅酸盐(SiN)制陶材料的强度在题中所给出的七大条件的影响下所能达到34
的最大值为1164.59。
5.2结果分析
我们使用残差和相关系数平方这两个指标对结果进行分析。
5.2.1残差标准差分析
残差的概念已经在上文中作过较为详细地介绍了,它实际上是反映了通过线性回归所得出的强度和相同的因素下给定强度的偏差量,能够在一定程度上说明建立的模型是否合理。
算出?式所引出的残差值和残差标准差:
k=27.1401 k=20.0219 110
k=-87.0122 k=39.5797 211
k=107.9382 k=-5.75 312
k=88.6648 k=17.5494 413
k=-31.2941 k=-92.4123 514
k=-76.7451 k=43.2783 615
k=47.4956 k=76.0861 716
k=-5.3832 k=-98.3532 817
k=-70.8042 9
18112ˆ2KYY,,(),,,iin,1i,1 又由算出K=66.2245。用Matlab作图如下:
8
通过上图可以清楚地看出某些点,诸如1,8,10,12,13五个点,残差很小,基本与给出数据相吻合;但是多数点的残差都很大,有的甚至达到了无法接受的程度,例如3、4点。显然,线性模型不够合理,无法如实地反映出强度随七个因素变化的真正趋势。但是,我们也必须看到,残差点剔除后,残差和残差的标准差明显的减小了,这说明了我们正确地剔除了残差点。
5.2.2相关系数分析
相关系数平方能够体现模型的拟合程度,即从点的分布趋势上反映了模型与实际情况的吻合程度。数值上越接近1,它的拟合程度就越好。它的计算公式为:
cov(,)XY22,,(,)()XY 6,DXDY
2,(,)XY由上式可以算出?式所对应的相关系数平方=0.5785,而去除残差点后5
2,(,)XY的?式所对应的相关系数平方=0.7661。由此可见,在取出了残差点后,6
整个模型的拟合程度得到了很大的改善。
9
本题的每一个样本点中都包含一组强度值,为了简单,我们将它们的平均值作为计算中的强度值。而由于这样的强度值应该是服从于正态分布的,用平均值计算必然会造成一定的误差。
在用线性回归之前,我们假定了各个因素对强度的影响是相对独立的。而事实上,这些因素在现实的操作中,总是互相影响,互相制约的。因此,当我们将各个因素主观地隔离后,就会产生一些不可预测的误差。
六( 模型优化
虽然残差点被删除后的数据集合符合线性回归的条件,可是线性回归本身就是一种不精确的拟合方法。相对于非线性拟合等方法来说,它的误差是比较大的,很多时候是不可以忽略的。所以,对于一些对精度要求比较高的应用领域,线性回归法是行不通的。这个时候,我们可以在线性模型的基础上进行非线性修正,提高精度,减少误差。
多项式拟合法一般取最高次幂小于或等于,。一旦超过,,精度增加得很小,计算量却会以几何级数增大。在这里,因为制陶产业里对精度的要求不是特别高,所以我们取最高次幂为,。在关系式中增加交叉项时,会极大的增加计算难度,而对精度的帮助却不大。权衡再三,我们决定只在优化模型中增加平方项,不加交叉项。牺牲一点精度,换来较高的效率和较好的可用性。
ˆY 首先建立非线性回归的数学模型。仍然设强度为,第i个因素为X,在此i基础上,建立平方项非线性回归的数学模型,也就是增加所有X的平方项,关i系式如下所示:
h,,01,,h022222222,,ˆ,, YXXXXXXX,,1234567,,,,
,,h,,07
10
h,,11,,h12,,,, ? ,XXXXXXX,1234567,,,,
,,h,,17
2N(0,) ,,
其中是不可观测的随机误差,可以忽略,而hij和是未知参数。所以实,,
际使用的方程可简化为:
h,,01,,h022222222,,ˆ,,+ YXXXXXXX,1234567,,,,
,,h07,,
h,,11,,h12,,,, ? XXXXXXX1234567,,,,
,,h17,,
下面求?式的未知数的系数:由最小二乘法推导可得:
TT()XXHXY, ?
T,1()XX当存在时,可以导出,的计算关系式如下:
T,1HXXXY,() ? 由于插入了平方项,所以X矩阵为一个17*14的矩阵,行向量是由Xi(i=1,2,…,7)
2和Xi(i=1,2,…,7)组成的;Y矩阵是一个17*1的矩阵,由17个强度值组成;H
2是系数矩阵,由Xi(i=1,2,…,7)和Xi(i=1,2,…,7)的所有系数构成。依据?式,通过矩阵运算可以求出H:(具体运算过程见附录二中的非线性模型程序。) H=[2301 -215.12 -375.8 174.13 -96.344 -386.17 -173.92 -733.43
33.098 106.36 -41.854 55.616 98.916 34.986]
然后通过?式求出当因素标号为(2,1,3,1,3,3,1)时达到最大值:
ˆ ,1252.7 Ymax
接下来,我们计算非线性模型的的残差和残差标准差:
k=2.3584 k=-4.8211 110
k=4.0161 k=-1.9842 211
11
k=3.2156 k=-2.7847 312
k=0.5163 k=3.552 413
k=-3.1778 k=0.7718 514
k=0.4596 k=-2.1218 615
k=-1.6055 k=5.6737 716
k=-5.2996 k=1.7143 817
k=-0.4830 9
又根据?式算出K=3.1749
最后,我们作出实验值、线性计算值和非线性计算值的折线比较图(见彩页)。从图中可以明显看出非线性计算值与实验值的折线图每个点都几乎完全吻合,远比线性计算值折线的拟合度好。非线性的残差和残差标准差远比线性的要小,也就是说的变化也可以看出,这些都说明了非线性模型的精度远远超过了线性模型的精度,可以应用于一些对精度要求较高的领域。
七(更合理试验设计计划
7.1优化设计的切入点
我们从减少约束因数的个数入手,削弱弱约束对于达到同样强度要求的限制。 7.2优化设计的具体方案
在我们做完线性模型的基础之上,调用Matlab中的逐步回归函数stepwise(x,y),通过比较均方差和回归系数来评价回归模型的优劣,即较高的回归系数和较小的均方差。而做到这一点就是通过剔除变量来实现的。最终我们得到在剔除烧结温度因数(F)下的模型Y。
222222Y=-551.59X+24.992X+72.826X-20.218X+40.417X+56.622X 123457
+1745.6X-176.01X-238.32X+90.505X-43.488X-250.21X123457
求解得到当
X=(2,1,1,2,3, ,1)时,Ymax=1102.2 在这种情况下,生产过程中烧结温度对陶瓷材料强度的影响变的不那么苛刻,从而降低了工艺的难度和生产的成本。虽然精度有所降低,但综合考虑精度和成本,根据实际的工程需要进行选择。
12
逐步回归分析图
13
参考文献
[1] 姜启源等 《数学模型》 高等教育出版社
[2] 张肇炽 《代数与几何》 高等教育出版社
[3] 田铮等 《随机数学基础》 西工大教学
[4] 张航等 《精通MATLAB 6.X》 清华大学出版社 [5] 苏金明等 《MATLAB 6.1 实用指南》 电子工业出版社
14
附录一:
陶瓷试验方案及强度数据表
试验号 因素 强度
A B C D E F G
1 1 2 2 1 3 1 3 996.8 783.6 796.9 2 1 2 1 2 2 3 1 843.8 816.2 714.3 824.4 3 1 2 3 3 1 2 2 647.1 667.9 534.3 617.7 4 1 3 2 1 2 3 2 616.3 552.3 552.6 596.0 5 1 3 1 2 1 2 3 517.8 526.1 498.1 499.5 6 1 3 3 3 3 1 1 1002.0 1097.0 882.9 940.1 7 1 1 2 2 3 2 1 806.5 933.5 964.9 1046.0 8 1 1 1 3 2 1 2 801.5 803.2 846.2 756.4 9 1 1 3 1 1 3 3 739.2 863.3 797.0 929.6 10 2 2 2 3 1 3 1 615.0 627.5 583.9 597.1 563.9 11 2 2 1 1 3 2 2 795.9 854.0 937.0 999.2 724.8 12 2 2 3 2 2 1 3 850.9 921.8 990.6 943.5 840.9 13 2 3 2 2 1 1 2 513.0 665.9 718.9 646.4 14 2 3 1 3 3 3 3 831.3 981.4 912.5 950.7 987.3 15 2 3 3 1 2 2 1 806.1 908.1 627.6 855.0 16 2 1 2 3 2 2 3 727.3 643.9 584.0 643.4 602.1 17 2 1 1 3 2 2 3 836.8 716.3 862.9 796.2 18 2 1 3 1 1 1 1 1001.0 937.6 955.3 995.8 1009.0
15
附录二:
初等模型函数:
%初等模型的建立
%输入自变量和0因变量 y_chu=[
858.8
799.67
616.75
579.3
510.38
980.5
937.73
801.73
832.27
746.85
862.18
919.54
636.05
932.64
799.2
640.14
803.05
979.74];
x_chu=[
1 2 2 1 3 1 1
1 2 1 2 2 3 1
1 2 3 3 1 2 2
16
1 3 2 1 2 3 2
1 3 1 2 1 2 3
1 3 3 3 3 1 1
1 1 2 2 3 2 1
1 1 1 3 2 1 2
1 1 3 1 1 3 3
2 2 2 3 1 3 1
2 2 1 1 3 2 2
2 2 3 2 2 1 3
2 3 2 2 1 1 2
2 3 1 3 3 3 3
2 3 3 1 2 2 1
2 1 2 3 2 2 3
2 1 1 3 2 2 3
2 1 3 1 1 1 1
];
%进行一般多元线性回归
[b,bint,r,rint,stats]=regress(y_chu,x_chu)%由regress()函数获取自变量的系数向量 %回归诊断
regstats(y_chu,x_chu)
one=[1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1]' x_chu1=[one,x_chu];
beta_chu %各自变量对应的回归系数
y_chu_mo=x_chu1*beta_chu
%相关系数的平方R-square--说明模型的拟合程度 R-square=0.5785 %F--因变量与自变量之间线性关系判断 F= 2.2168
F(7,9)=2.72 满足线性假设条件
%显著性概率p--回归的显著程度判断 p=0.0220<0.05 回归显著
19
非线性模型函数:
%非线性优化模型
x=[
1 2 2 1 3 1 1 1 2 1 2 2 3 1
3 1 2 2 1 2 3
1 3 2 1 2 3 2 1 3 1 2 1 2 3
3 3 1 1 1 3 3
1 1 2 2 3 2 1 1 1 1 3 2 1 2 1 1 3 1 1 3 3 2 2 2 3 1 3 1 2 2 1 1 3 2 2 2 2 3 2 2 1 3 2 3 2 2 1 1 2 2 3 1 3 3 3 3 2 3 3 1 2 2 1 2 1 1 3 2 2 3 2 1 3 1 1 1 1 ];
y=[858.8
799.67
616.75
579.3
510.38
980.5
20
937.73
801.73
832.27
746.85
862.18
919.54
636.05
932.64
799.2
803.05
979.74
];
w=x(:,[1:7]);
v=[w w.^2] %矩阵x的扩展 xishu=inv(v'*v)*v'*y result=v*sishu
n=0:18;
plot(y,'r');hold on; plot(result,'b'); axis([0,18,450,1000]) grid on
21
附录三:
科技进步申请报告
申请人:
申请项目:材料强度最大化工艺的快速定位
申请奖项:
问题背景:在寻找使硅酸盐(SiN)制陶材料具有最大强度值的最佳工艺时,34
我们需要先在不同的工艺条件下测出其强度值,然后对所得数据进行处理,从而得到最佳工艺。然而,在众多的试验数据中,要想找到某种材料能形成最大强度的精确加工条件是一件费时又费力的事情。传统的流程分为两步:首先要做足够多的试验,得到尽可能多的数据;然后对这些数据进行无差别的拟合或是不加判断地综合分析。这种方法是极耗时间的。假如有,个加工条件,每个因素(加工
nn条件)有,个水平。若是进行无差别分析时,则至少需做m次试验,对m个数据进行分析处理。而实际的加工过程当中,影响强度的因素是极多的,这种方法几乎是不可能实现的。因此,处理数据的方法显得极为重要。我们在本试验的基础上提出了一种非常有效快捷的数据处理方法。
基本思想:这种处理方法核心是剔除对强度影响可忽略的数据量。面对浩如烟海的数据,我们的想法是从中剔除一些严重不符实际的的数据,然后再依据一些条件和要求,找出一些对强度影响很大的数据,仅对这些数据进行分析,初步求出结果,再利用一些影响很小,但又不可忽略的数据对结果进行修正。
核心技术:这种方法的关键是找到剔除阀值和正确挑出影响很大的数据。我们在做寻找使硅酸盐(SiN)制陶材料具有最大强度值的最佳工艺的试验时,具34
22
体做法是:
? 由于数据比较多,我们首先对已有的数据初步建立一个线性回归模型,利用这个模型引出残差和残差标准差,由随机数学中的定理可知,残差大于残差标准差两倍的数据是严重不合理的,必须剔除。经我们验证,剔除不合理的数据后,计算量减少得不明显,可是偏差大大地减小了。
? 我们在确立了较为精确的模型后,利用Matlab中的一个函数,将各个因素对强度的影响用可调的图形表现出来。在图中可以暂时删除某个因素对强度的影响,这样对精度的影响若是在工程允许的范围内(比如说相关系数平方大于某个阀值或是均方差小于某个阀值),我们就可以不考虑它们的影响,只考虑剩下的重要因素的影响。这样便大大减少了计算量。
在该试验的基础上我们提出了一种处理大量数据的有效快捷的方法:对需要处理的数据,首先为它找一个合适的模型,利用这个模型和上面所说的残差理论剔除不合理的数据。然后,利用工程问题对参量特定的要求,剔除对精度影响较小的因素。这样,在工程允许的情况下,计算量和偏差值都大大减少了。
应用前景:材料强度最大化的快速定位法能够广泛的应用于工程实践当中。传统的方法有可能会得出比快速定位法好得多的结果,但却是建立在高成本和高劳动量的基础上的。快速定位法为了提高效率和速度牺牲了一些精度。然而,这在工程上是允许的。因为工程所关注的是材料所具有的强度能否满足要求以及收益有多大。在同时满足条件的情况下,工程人员必然会选用快速定位法。所以,快速定位法是一种高效实用的方法,将会为我们赢取很多的时间和极大的效益。
23
24