實驗數據的處理及模型參數的確定第一頁，共五十一頁，2022年，8月28日●數學模型中各參數的確定例：鎳硅藻土上苯加氫合成環(huán)己烷是表面反應控制的固體催化劑上的氣相反應。在160oC，微分反應器中的初始反應速率方程為模型參數

ka

──表觀速率常數

bH

──H2的吸附系數

bB──

C6H6的吸附系數利用實驗得到的全部信息，確定數學模型中的待定參數第二頁，共五十一頁，2022年，8月28日

●線性插值

●

Lagrange插值

●埃米爾特插值

●一元線性回歸

●線性模型的推廣

●

多元回歸

可化為多元線性回歸的問題

●

多項式擬合簡介

●逐次回歸分析函數關系插值法回歸分析相關關系數值微分★★★★★引言：2.常用的數學方法第三頁，共五十一頁，2022年，8月28日例：72型分光光度計測得某試樣的吸收值如下：2-1-1–1線性插值——問題的提出λ/nm430440450460470480A0.4100.3750.3250.2800.2400.205希望：根據給定的函數表作一個既能反應f(x)的特性，又便于計算的簡單函數p(x)，用p(x)近似f(x)，計算出任意

x對應的y值求在435，445，455，465，475nm處的吸收值。第四頁，共五十一頁，2022年，8月28日

定義：設y=f(x)在區(qū)間[a,b]上有意義，且已知在點a<x0<x1<…<xn<b上的值y0,y1,…,yn,若存在一簡單函數pn(x),使

pn(xi)=yi(i=0,1,…,n)成立，則稱pn(x)為f(x)的插值函數，x0，x1，…，xn為插值節(jié)點區(qū)間[a,b]為插值區(qū)間，求pn(x)的方法稱為插值法

xyy=f(x)y=p(x)x1y1xn

yn幾何意義：2-1-1–2線性插值——方法原理ab第五頁，共五十一頁，2022年，8月28日線性插值原理：兩點間直線方程：xyy=f(x)y=p(x)xi-1yi-1xi

yi2-1-1–2線性插值——方法原理第六頁，共五十一頁，2022年，8月28日分段線性插值：實驗點個數為n時，求插值結點x的函數值。首先確定x在哪兩點間

2-1-1–2線性插值——方法原理第七頁，共五十一頁，2022年，8月28日LINEPLOT(N,X,Y,X0,Y0)DOJ=1,N-1J1=J+1X0<=X(J1)CONTINUEJ=J-1T=(X0-X(J))/(x(J1)-x(J))Y0=Y(J)+T*(Y(J1)-Y(J))RETURNnoyes2-1-1–3線性插值——程序框圖第八頁，共五十一頁，2022年，8月28日開始輸入：數據點X(I),Y(I)，未知點X0調用線性插值子程序求未知點X0對應的函數值Y0輸出：X0,Y0值結束2-1-1–4線性插值——應用示例顯示程序顯示輸入顯示輸出第九頁，共五十一頁，2022年，8月28日2-1-2-1一元三點Lagrange插值——問題的提出例：計算乙醇的平均摩爾體積實驗測得25℃時乙醇溶液的平均摩爾體積（cm2mol-1）與乙醇的物質的量分數的關系如下x/cm2mol-10.089121.220.115322.160.143523.180.173924.320.206825.570.242426.950.281128.470.323430.150.369732.010.420734.070.477136.37計算x=0.1,0.2,0.3,0.4時的。第十頁，共五十一頁，2022年，8月28日線性插值公式：二點(xi-1,yi-1),(xi,yi)(兩點式)Lagrange插值（三點插值,拋物線插值）：xi-1xixi+1即2-1-2-2一元三點Lagrange插值——方法原理第十一頁，共五十一頁，2022年，8月28日xyy=f(x)y=p(x)xi-1yi-1xi+1

yi+1xiyi編程難點：如何確定使用哪三個結點進行插值xj-1xj+1xjxj-2xj+22-1-2-2一元三點Lagrange插值——方法原理第十二頁，共五十一頁，2022年，8月28日LGRG2(X,Y,N,T,Z)DoJ=3,N-1I=JT>X(I)CONTINUEP=(T-X(I))*(T-X(I+1))/(X(I-1)-X(I))/(X(I-1)-X(I+1))Q=(T-X(I-1))*(T-X(I+1))/(X(I)-X(I-1))/(X(I)-X(I+1))R=(T-X(I-1))*(T-X(I))/(X(I+1)-X(I-1))/(X(I+1)-X(I))Z=P*Y(I-1)+Q*Y(I)+R*Y(I+1)RETURNnoyesI=I-1|T-X(I-1)|<=|T-X(I)|yesno2-1-2-3一元三點Lagrange插值——程序框圖第十三頁，共五十一頁，2022年，8月28日開始輸入：數據點X(I),Y(I)，未知點X0調用lagrange插值子程序求未知點X0對應的函數值Y0輸出：X0,Y0值結束2-1-2-4一元三點Lagrange插值——應用示例顯示程序顯示輸入顯示輸出第十四頁，共五十一頁，2022年，8月28日2-2-1-1一元線性回歸——問題的提出例:銅鉬礦中鉬對銅含量的線性依賴關系No123456789x(Cu)285290300303310318325335338y(Mo)4.64.74.74.94.95.15.05.35.4一元線性回歸的數學模型：y=ax+b+εyi=axi+b+εin個實驗點回歸直線：y=ax+b殘差：εi=yi-（axi+b）第十五頁，共五十一頁，2022年，8月28日最小二乘法：第i點殘差：εi=yi-（axi+b）當殘差的平方和為最小時，對應的a、b值是最佳值。（正規(guī)方程組）2-2-1-2一元線性回歸——方法原理第十六頁，共五十一頁，2022年，8月28日令平均值離差平方和2-2-1-2一元線性回歸——方法原理第十七頁，共五十一頁，2022年，8月28日1.線性相關系數R—衡量回歸方程式與數據相符合的程度。若R1，則數據點落在直線上。注意：2.加權最小二乘法3.剔除可疑數據2-2-1-2一元線性回歸——方法原理第十八頁，共五十一頁，2022年，8月28日PK(N,X,Y,A,B,R)SX,SY,SXX,SXY,SYY=0DoI=1,NX1=X(I),Y1=Y(I)SX=SX+X1,SY=SY+Y1SXX=SXX+X1*X1SYY=SYY+Y1*Y1,SXY=SXY+X1

LXX=SXX-SX*SX/n,LYY=SYY-SY*SY/n,LXY=SXY-SY*SX/NB=LXY/LXX,A=(SY-B*SX)/N,R=LXY/SQRT(LXX*LYY)

RETURN2-2-1-3一元線性回歸——程序框圖第十九頁，共五十一頁，2022年，8月28日開始輸入：數據點數N

銅與鉬的實驗數據X(I),Y(I)（I=1，N）調用一元線性回歸子程序計算A，B，R輸出：A，B，R結束2-2-1-4一元線性回歸——應用示例顯示程序顯示輸入顯示輸出第二十頁，共五十一頁，2022年，8月28日2-2-2-1線性模型的推廣——方法原理變量x與y之間存在某種非線性關系確定曲線類型（非線性關系）實際經驗散點圖形狀線性關系最小二乘法確定系數非線性關系第二十一頁，共五十一頁，2022年，8月28日曲線類型及變換公式雙曲線型冪指數型指數型S型對數型平方根曲線2-2-2-1線性模型的推廣——方法原理第二十二頁，共五十一頁，2022年，8月28日例1：Arrhenius公式的應用根據k和T數據，可確定指前因子A和活化能Ea。2-2-2-2線性模型的推廣——應用示例第二十三頁，共五十一頁，2022年，8月28日例2：Clausius-Clapryron方程式的應用純組分氣-液（氣-固）兩相平衡的方程式：上式中：p:T/K時液（固）飽和蒸氣壓；ΔH：相變熱不定積分：測定不同溫度下的飽和蒸氣壓，將lnp~1/T進行線性回歸，可算出ΔH，并計算其它溫度下的蒸氣壓2-2-2-2線性模型的推廣——應用示例第二十四頁，共五十一頁，2022年，8月28日開始輸入：數據點數N

溫度T與蒸氣壓p的實驗數據T(I),P(I)（I=1，N）輸出：A,B,ΔH,Ti,pi結束調用線性回歸子程序計算A，B（相變熱ΔH=-8.314E-3*B）T(I)=T(I)+273.15X(I)=1/T(I)Y(I)=ln[P(I)](I=1,N)計算其它溫度Ti下的蒸氣壓pi2-2-2-2線性模型的推廣——應用示例第二十五頁，共五十一頁，2022年，8月28日例3：幾種常用的吸附等溫式的計算（氣固吸附）2-2-2-2線性模型的推廣——應用示例Freundlich經驗式：Langmuir方程：B.E.T.方程：第二十六頁，共五十一頁，2022年，8月28日2-2-3-1多元線性回歸——方法原理數學模型：函數y與多個自變量間x1，x2，…，xm的線性相關關系設共進行了n次測定（i=1，2，…，n）的自變量取值分別為xi1，xi2，…，xim，函數值的測定值為yi，所得的值如表：第二十七頁，共五十一頁，2022年，8月28日選擇回歸系數，以使殘差的平方和最小，即殘差殘差的平方和Q殘差的平方和最小，相當于求解以下方程組（i=1，2，…，n）（k=，0，1，…，m）最小二乘法2-2-3-1多元線性回歸——方法原理第二十八頁，共五十一頁，2022年，8月28日正規(guī)方程組：全相關系數R式中（m+1）元線性方程組，Guass消去法求解b0，b1，…，bm2-2-3-1多元線性回歸——方法原理第二十九頁，共五十一頁，2022年，8月28日MLR(M,N,XV,S,B0,R)輸入：X（I，J）I=1,2,…,NJ=1,2,…,M+1CALLGS(A,M,M1,EPS)RETURN程序框圖DOI=1,MA(I,I)=1.0DOJ=I+1,M1A(I,J)=0DOK=1,NA(I,J)=A(I,J)+X(K,I)*X(K,J)DOI=1,M-1J=I+1,MA(J,I)=A(I,J)2-2-3-2多元線性回歸第三十頁，共五十一頁，2022年，8月28日變量xi(i=1,2,…,m)與y之間存在某種非線性關系確定曲線類型（非線性關系）實際經驗散點圖形狀線性關系最小二乘法確定系數非線性關系2-2-3-3多元線性回歸——可化為多元線性回歸的問題變量代換變量代換第三十一頁，共五十一頁，2022年，8月28日函數式通式：y’=b0+b1x1+b2x2變換變量變換常數y’=x1=x2=b0=b1=b2=y=a+bx+cx2yxx2abcy=alnx+bcosx+cylnxcosxcaby2

=a+bx+cx2y2xx2abcy＝x/

(a+bx+cx2)x/yxx2abcy=a+b/x+c/x2y1/x1/x2abclnyxx2lnabc2-2-3-3多元線性回歸——可化為多元線性回歸的問題第三十二頁，共五十一頁，2022年，8月28日例1：Antoine方程式的應用（p：蒸氣壓，T：溫度）令Tlgp=y,T=x1,lgp=x2,b0=ac+b,b1=a,b2=c2-2-3-4多元線性回歸——應用示例第三十三頁，共五十一頁，2022年，8月28日例2：用鎳硅藻土作催化劑，苯加氫合成環(huán)己烷。用微分反應器測定和分析得到160oC的初始反應速率以及相應的氫和苯的分壓值pH/kPa，pB/Kpa和r0/(mol/Kg·h)的數據初始反應速率方程為：請利用上述實驗數據擬合出參數ka、bH及bB。74.9467.2157.7650.7592.5626.7034.2443.4250.4310.20218.2220.8223.5189.2117.692.6687.6675.6456.1752.419.9714.7126.0745.0148.77115.1147.2217.8212.2202.4pHpH

pB

pBr0r02-2-3-4多元線性回歸——應用示例第三十四頁，共五十一頁，2022年，8月28日實驗測得pH,pB和r0數據取倒數移項，開1/4方（線性化）2-2-3-4多元線性回歸——應用示例第三十五頁，共五十一頁，2022年，8月28日開始輸入：自變量個數M，數據點個數N,反應溫度T，氣體分壓pH，pB，R0的實驗數據X(I,j)輸出：bH、bB和k結束調用多元線性回歸子程序計算B0，B1和B2M1=M+1，X1(I,j)＝X(I,j)X1(I,M1)＝DSQRT（DSQRT（X(I,1)**3*X(I,2)/X(I,M1)））（I=1,N）

X(I,J)=X1(I,J)（I=1,N，J=1,M1）計算：bH=A(1,M1)/B0；bB=A(2,M1)/B0；k=1/((B0**4)*(KH**3)*KB)顯示程序

顯示輸入

顯示輸出2-2-3-4多元線性回歸——應用示例第三十六頁，共五十一頁，2022年，8月28日例3：化學反應動力學方程的總級數n設：某氣體反應A+BC化學反應速率方程式中：ka—反應表觀速率常數；

pA，pB

和pC—參加反應各氣體A，B和C的氣相分壓；

，和—參加反應各物質在化學反應速率方程中的分級數。此反應的總級數n=++取對數n=++（線性化）2-2-3-4多元線性回歸——應用示例第三十七頁，共五十一頁，2022年，8月28日開始輸入：數據點數N,反應溫度T，氣體分壓p的實驗數據,pa(I),pb(I),pc(I)（I=1，N）輸出：n,ka，，和,結束調用多元線性回歸子程序計算ka，，和X1(I)=ln(pa(I)),x2(I)=ln(pb(I)),x3(I)=ln(pb(I))Y(I)=lnv(I)(I=1,N)計算反應的總級數n=++2-2-3-4多元線性回歸——應用示例第三十八頁，共五十一頁，2022年，8月28日例4分子結構-性能的多元線性回歸Hammet方程：苯環(huán)上間位或對位取代基對反應速率的影響式中：k——反應速率常數

k0——未取代時母體的反應速率常數ρ——反映取代基電子效應的結構常數σ——與反應類型有關的結構常數2-2-3-4多元線性回歸——應用示例第三十九頁，共五十一頁，2022年，8月28日推廣：考慮其他多種因素對分子性能的影響分子具有的性質：“應答”具有特定功能的結構參數π——反映取代基親油/親水能力E——反映取代基的空間效應腎上腺阻斷劑：N-N-二甲基-2-溴苯乙胺衍生物結構—性能之間的關系數學模型：2-2-3-4多元線性回歸——應用示例第四十頁，共五十一頁，2022年，8月28日

2-2-4多項式擬合簡介——方法原理數學模型：殘差的平方和Q殘差的平方和最小，相當于求解以下方程組最小二乘法第四十一頁，共五十一頁，2022年，8月28日2-3-1數值微分——問題的提出實驗測定的一批離散點t/[t]t1t2……tncA/[cA]cA1cA2……cAn有化學反應

aA產物化學反應速率方程：第四十二頁，共五十一頁，2022年，8月28日在恒容反應中，化學反應速率r可用反應物A濃度隨時間變化率來表示式中，dcA/dt即為反應物A的濃度cA~時間t曲線上某點的斜率，數學上若cA~t是連續(xù)單值函數，則dcA/dt是在t點的一階導數。t/[t]cA/[cA]t2cA22-3-1數值微分——問題的提出cA1t1第四十三頁，共五十一頁，2022年，8月28日2-3-2數值微分——方法原理在微積分中，函數的導數是通過極限來定義的，有完整的計算方法。而在實際工作中，常常需要求列表函數在節(jié)點和非節(jié)點處的導數值，這就是數值微分要解決的問題?；瘜W化工中有不少實際問題都需用數值微分求導來解決。

t/[t]t1t2……tncA/[cA]cA1cA2……cAn……第四十四頁，共五十一頁，2022年，8月28日2-3-2數值微分——方法原理設有若干等距離的節(jié)點，欲求節(jié)點處的一階導數，可用差商代替微商進行計算考慮函數f(x)在點（x0+h）處的泰勒展開式（向前差商式）CB（向后差商式）AC（中心差商式）ABxyx0-hx0x0+hABCy＝f(x)第四十五頁，共五十一頁，2022年，8月28日中心差商定義：（1）若有函數y=f(t)，按等間隔在s個t點(t1,t2,…,ts-1,ts)測出相應的y值（y1

，y2

，…，ys

），除兩端t1與ts以外，任意點的導數yi’=f’(ti)可按下式計算：（2）（i=2，3，…，s-1）若y＝f(x)為連續(xù)有界單值函數，客觀存在在x點的一階導數,可用較為精確的“中心差商”來計算：用插值法求2-3-2數值微分——方法原理第四十六頁，共五十一頁，2022年，8月28日兩端點的導數值用以下兩式計算：（3）（4）式中t是獨立變量t的間隔。：用插值法求2-3-2數值微分——方法原理第四十七頁，共五十一頁，2022年，8月28日CF(N,X,Y,H,R)T=X(1)+H調用插值子程序計算FYA=FT=X(1)+2H調用插值子程序計算FYB=FW=2T=X(W)-HR(1)=[4YA-3Y(1)-YB]/2HAB