計(jì)算機(jī)化工應(yīng)用w08插值與擬合

通下周課程安排上機(jī)實(shí)習(xí)，請(qǐng)準(zhǔn)時(shí)

插值與擬合LagrangeHermite分段三次樣條

一般采用函數(shù)表示自然現(xiàn)象、過(guò)程中各參數(shù)的變化規(guī)yF(x)解析 列表 圖示要的函數(shù)近方法4.1

(x)F(x)

(φ(x)具有簡(jiǎn)單的形式函數(shù)溫密粘函數(shù)插值的定插值就是構(gòu)造一個(gè)函數(shù)在指定點(diǎn)取給定值的x0，x1，…，xn處的函數(shù)值y0，y1，…，yn插值原①f(x)(①f(x)(x),x[a,b];②(xi)y,i0,1,i,當(dāng)選用不同的插值函數(shù)時(shí)，相應(yīng)的得到不同的插值方 代數(shù)多項(xiàng)式插值代數(shù)多

P(x)P(x)aaxax2

ann代數(shù)多項(xiàng)式插值基本原 構(gòu)造多次插值多項(xiàng))n+1 P(x)P(x)aaxax2

ai如果以上的多項(xiàng)式滿足插值原則②，就稱(chēng)P(x)F(x)的插值多項(xiàng)式代數(shù)多項(xiàng)式插值基本原 確定ｎ次多項(xiàng)式系aaxax2 確定ｎ次多項(xiàng)式系 1

aaxax2

axn 1 2 axa

n

1

axn 000,xn)1x111 0 000,xn)1x111 0i nnn由于節(jié)點(diǎn)互異，行列式不等于0，方程多項(xiàng)式一些簡(jiǎn)單的方法進(jìn)行構(gòu)造， 日(Lagrange)線性插值插值多項(xiàng) 1(x)a0y0

y0y1

y0a0 y a , ya y x

x 日線性插1ay0y1, ya 1

0y0y10x0

x01(x)a0P(x)

xx1

x

直線的兩點(diǎn)式表達(dá)L(x)xx1 LL(x)xx1 L(x)x0x101x10 分別稱(chēng)為節(jié)點(diǎn)x0和x1的一次插值日二次插二次插值多項(xiàng)P2(xy0L0xy1L1(xy2L2x)基函數(shù)如何確定

L(x)1,(xxii0,x在節(jié)點(diǎn)上，但xxii對(duì)于i=0，可假定L0(x)＝A(x-x1)(x-

(xx)(xxL0(x) (x0x1)(x0x2 (xx)(xx1L(x) 1 (x1x2)(x1x0 (xx)(xx)2L(x) 2 (x2x0)(x2日n次插nnLn1xlkxknln

x

xx1xx2 xxk1xxk1 xxnk xk

xxxx x x

xx

k

k 1 lx1

iki例例 求二次插值多項(xiàng)式

nLn1xlkxykn

L0(x)

(xx1)(xx2)(x0x1)(x0(xx)(xxk

L1(x) (x1x2)(x1x0) (xx)(xxL2(x) 解 日方法，有

(x2x0)(x2L(x)(x2)(x4)8(x1)(x4)1(x1)(x2)(12)(1 (21)(2 (41)(43x216x例題171512]；試編寫(xiě)一個(gè) 例題functionCha4demo2x=1:6;y=[161821171512];fork=1:nforj=[1:k-1k+1:n]w=(xi-x(j))./(x(k)-

plot(x,y,'o',xi,v,'-日插值法的特LagrangeLagrange7現(xiàn) （Runge）就給出了一個(gè)例子f(x)1(125x2 1xi10取等矩節(jié)點(diǎn) 12i/n(i0,1,,n)， 日插值多項(xiàng)i10當(dāng)n10

yf(x)及插值多項(xiàng)式y(tǒng) (x

現(xiàn)y

f

L(x)0-0-

(x)比較接近f 在區(qū)間[-1，1]兩端則誤差很大。當(dāng)n繼續(xù)增大時(shí)，部分區(qū)上插值多項(xiàng)式誤差偏大的現(xiàn)象更嚴(yán)重。這種現(xiàn)象稱(chēng) 現(xiàn)。分段插 分段線性插值在每個(gè)分段插值區(qū)間[xi,xi+1]

L(i)(x)yl(x) (x)

x

x i i1 x

L(0)

x0x1 xx

L1(xi)

x分段三次Hermite許多最有效的插值技術(shù)都是基于分段三次多項(xiàng)

H(x) axax2ax3,i0,1,,① i ,①, H(xi)yi,(i,, H'(xi)mi,(i,三次樣條插值三次樣條插值函三次樣條插值原在n個(gè)小區(qū)間構(gòu)造S(x)，共有n個(gè)三次多項(xiàng)式，需確定4n個(gè)參在所有節(jié)

S(xi)yi,i0,1,,

n+1個(gè)方limS(x)limS(x)S(xi)yi,i1,,nlimS(x)S(xi),i1,,nlimS(x)S(xi)Mi i1,,n13(n-1)3(n-1)方邊界條

S”(x0)=S”(xn)=0，自然邊界條已知

2個(gè)方一維插值中的函數(shù)interp1可實(shí)現(xiàn)一維插值運(yùn)調(diào)用格式y(tǒng)i=x、y為樣本點(diǎn)，yi為插值點(diǎn)自變量值xi對(duì)應(yīng)的函數(shù)值xi‘nearest’——最近插例已知插值節(jié)點(diǎn)x=1:6，節(jié)點(diǎn)處函數(shù)值y=[1618211512]；試采用interp1函數(shù)計(jì)算最近插值、線性插 例functionCha4demo3x=1:6;y=[161821171512];plot(x,y,'o',xi,y4,'-例不同插值方法的效果例刻度校正刻度校正例3已知某轉(zhuǎn)子流量刻度校正刻度校正例plot(X,Y,'b-o'),holdon執(zhí)行結(jié)果：Yk逆插已知函數(shù)值f(x)，求自變量x的值現(xiàn)需調(diào)整實(shí)際流量為850ml/min，則流量計(jì)逆插plot(Y,X,'b-o'),holdonf=@(x)Yk-interp1(X,Y,x);pchip和 插 樣條插yi=pchip(x,y,xi)

yi=spline(x,y,xi)法與interp1函數(shù)功能相同第二種方法將三次多項(xiàng)式以結(jié)構(gòu)體pp于利用其它函數(shù)進(jìn)行微分、積分等運(yùn)例題利用點(diǎn)(x=sin(kπ/6),y=cos(kπ/6)，其中k=[0123]來(lái)近單位圓的前四分之一圓周。比較pchip與gridon,axisequallegend('Orignal例題10

pchip與spline的區(qū)樣條插值函數(shù)和pchip插值函數(shù)，以及它們的一、二、三保形分段三次Hermitepchip(piecewisecHermitegpolynomial):(shap-preserving端點(diǎn)處另spline。例利用(1,0(0,11,0),(0,-1),(1,0)xx=y=[10-10010-1pp=yy=ppval(pp,axis例xx=y=10pp=yy=ppval(pp,axis例表乙醇－水二元溶液的氣液平衡曲線實(shí)驗(yàn)x0y0xy例解 計(jì)算程序如下這這里：x和y分別表示樣本點(diǎn)的液相和氣相中乙醇摩爾分?jǐn)?shù)xk和yk分別表示插值點(diǎn)的液相和氣相中乙醇摩爾分?jǐn)?shù)例樣本點(diǎn)三樣本點(diǎn)三次樣條0.80.70.60.50.40.30.20.10 0.8 0.9三次樣條插值乙醇－水二元溶液的氣液平衡曲插值函數(shù)的選初始數(shù)據(jù)對(duì)于插值的影已知x=[0:2:10]，y=[00.9093- -0.27940.9894-(y＝sinx)，比較一維線性、線性最近、立方和三次樣條插值所ezplot(@sin,[010])plot(xi,yi,'r+'),text(0.5876,0.2537,'\leftarrow線性插值

10

三次樣條插三次插線性插

最近插

插值函Zi=VI=例題：精餾塔理論板例題：精餾塔理論板xy溫xy溫 xy溫xy溫例題：精餾塔理論板McCabe-Thiele法求二元精餾的理論塔板數(shù)實(shí)際就是逐板計(jì)yn1

R

xn

R1yq

q1

q1例題：精餾塔理論板function計(jì)算所需q=1.103;R=3;xD=0.78;xW=0.02;xF=%相平衡數(shù)x0=[0.05.10.20.30.40.50.60.70.80.894.95y0=[0.31.43.52.575.614.657.698.755.82.894Yr@(xR/(R+1).*xxD/(R+1);%精餾段操作線fun=@(x)(q-1)*(R/(R+1).*x+xD/(R+1))-(q*x-xQfzero(fun,0.5%求操作線交yQ=xOP=[xW,xQ,yOP=[xW,yQ,yfit=xfit=例題：精餾塔理論板%續(xù)上繪制圖holdonboxonplot([01],[0k=1;yn(1)=xn(1)=interp1(y0,x0,yn(1),'pchip');text(xn(1),yn(1),num2str(1),...例題：精餾塔理論板%續(xù)上whilexn(k)>yn(k+1)=interp1(xOP,yOP,xn(k));k=k+1;xn(k)=interp1(y0,x0,yn(k),'pchip');text(xn(k),yn(k),num2str(k),...N=text(xn(N),yn(N),num2str(N),...text(0.5,0.5,strcat('所需理論板：',num2str(N)),…例題：精餾塔理論板1 3 yy6

所需理論板：進(jìn)料板位置：70 x

插值與擬合線性最小二非線性最小二乘

最小二乘法曲線擬0實(shí)曲線f(x)系數(shù)f(x) 地近（i,i 據(jù)數(shù)稱(chēng)。最小二乘法曲線擬合的常用常用最小二乘法進(jìn)行曲線擬

y k若Qe2 最k稱(chēng)此曲線擬合法為稱(chēng)此曲線擬合法為最小二乘法按計(jì)算方法特點(diǎn)又分為線性和非線性最小二乘法一元線性最小二乘 y

yaxk(k1,2, ,myk(k1,2, ,m)ky*(k1,k (ab 一元線性最小二乘由最小二乘法：設(shè)Q

e2

( (abx QQ 求解該法方程即可得求解該法方程即可得

Q2m

a

mabxk

yk

ax x2 x

abxkxk

k k

多元線性最小二乘yb1x1b2x2b3x3

bnxnbixinin共n個(gè)影響因子，有mmxkyk代表第k

y*ky*b b b bk f(x1k,x2k,x3k ,xnk;1,b2,

多元線性最小二乘

e2

kk(yy*kk

k k

b2(b1x1kxikb2x2kxik bnxnkxikx1kyk)

,Sijxikxjk(ij, ,, x y(is11b1

s1nbn

s22b2

2nbn

s b b 非線性最小二乘在應(yīng)用最小二乘法曲線擬合時(shí)，通常遇 的是非線性函數(shù)對(duì)比線性模型擬合，非線性模型擬合 的多 通過(guò)代換轉(zhuǎn)化為線性關(guān)

直接采用非線性擬非線性函數(shù)的線性化處化工中常見(jiàn)的函指數(shù)函負(fù)指數(shù)對(duì)數(shù)S型曲n次多項(xiàng)

均可以線性化處非線性函數(shù)的線性化 雙曲線 a YYabXY1 X 冪函數(shù)y

lnylncblnx1xYA令：Yln XlnYA非線性函數(shù)的線性化處YAYA

ycebxlnylnc

a

x令：Yln A