[理學(xué)]CP030-計算物理數(shù)據(jù)插值0.ppt

1、實驗數(shù)據(jù)的插值,根據(jù)這些數(shù)據(jù)，希望合理地估計出其它溫度（如25攝氏度，40攝氏度）時的電阻,舉例,這就是本章要討論的“插值問題”,當(dāng)精確函數(shù) y = f(x) 非常復(fù)雜或未知時，在區(qū)間a,b上一系列節(jié)點 x0 xm處測得函數(shù)值 y0 = f(x0), , ym = f(xm)，由此構(gòu)造一個簡單易算的 近似函數(shù) g(x) f(x)，滿足條件 g(xj) = f(xj) (j = 0, m) (*) 這個問題稱為“插值問題”,插值問題的定義,這里的 g(x) 稱為f(x) 的插值函數(shù)。,節(jié)點 x0 xm稱為插值節(jié)點,條件（*)稱為插值條件，區(qū)間a,b稱為插值區(qū)間,f(x),g(x),最常用的插值函

2、數(shù)是 ?,代數(shù)多項式,用代數(shù)多項式作插值函數(shù)的插值稱為代數(shù)插值,本章主要討論的內(nèi)容,插值函數(shù)的類型,插值問題,插值法,插值函數(shù),代數(shù)插值中的三個問題,一、插值問題解的存在唯一性？ 二、插值多項式的常用構(gòu)造方法？ 三、插值函數(shù)的誤差如何估計？,代數(shù)插值,代數(shù)插值問題解的存在惟一性,給定區(qū)間a,b上互異的n+1個點xjnj=0的一 組函數(shù)值f(xj)，j =0，, n，求一個n次多項式pn(x)Pn，使得 pn(xj)=f(xj)，j=0,1,，n. . (1) 令 pn(x)=a0+a1x+anxn, . (2) 只要證明Pn(x)的系數(shù)a0 ,a1, an存在唯一即可,(3),為此由插值條件(

3、1)知Pn(x)的系數(shù)滿足下列n+1個代數(shù)方程構(gòu)成的線性方程組 a0+a1x0+anx0n=f(x0) a0+a1x1+anx1n= f(x1) . a0+a1xn+anxnn= f(xn),代數(shù)插值問題解的存在惟一性,而ai(i=0,1,2,n)的系數(shù)行列式是Vandermonde行列式,由于xi互異，所以(4)右端不為零，從而方程組(3)的解 a0 ,a1 ,an 存在且唯一。,代數(shù)插值問題解的存在惟一性,通過解上述方程組(3)求得插值多項式pn(x)的方法并不可取.這是因為當(dāng)n較大時解方程組的計算量較大，而且方程組系數(shù)矩陣的條件數(shù)一般較大（可能是病態(tài)方程組）,當(dāng)階數(shù)n越高時，病態(tài)越重。,

4、為此我們必須從其它途徑來求Pn(x)：不通過求解方程組而獲得插值多項式,代數(shù)插值問題解的存在惟一性,基本思想：在n次多項式空間Pn中找一組合適的基函數(shù) 0(x),1(x), n（x),使,pn(x)=a0 0(x) +a1 1(x) +an n（x),不同的基函數(shù)的選取導(dǎo)致不同的插值方法,Lagrange插值,Newton插值,代數(shù)插值問題解的存在惟一性,兩節(jié)點-一次（線性）插值,兩節(jié)點-一次（線性）插值,兩節(jié)點-一次（線性）插值,兩節(jié)點-一次（線性）插值,兩節(jié)點-一次（線性）插值,Function y=chazhi1(x0,y0,x1,y1,x)%拉格朗日插值法 y=(x-x1)/(x0-x

5、1)*y0+( x-x0)/(x1-x0)*y1 Function y=chazhi1(x0,y0,x1,y1,x)%牛頓插值法 y=y0+(y1-y0)/(x1-x0)*( x-x0),Matlab 編程實現(xiàn),三節(jié)點-二次插值,三節(jié)點-二次插值,三節(jié)點-二次插值,Function y=chazhi1(x0,y0,x1,y1,x2,y2,x)%拉格朗日插值法 y=(x-x1)*(x-x2)/(x0-x1)*(x0-x2)*y0+ (x-x0)*(x-x2)/(x1-x0)*(x1-x2)*y1+ (x-x0)*(x-x1)/(x2-x0)*(x2-x1)*y2 Function y=chazh

6、i1(x0,y0,x1,y1, x2,y2,x)%牛頓插值法 y=y0+(y1-y0)/(x1-x0)*( x-x0)+ (x-x0)*(x-x1)*(y2-y1)/(x2-x1)-(y1-y0)/(x1-x0)/(x2-x0),Matlab 編程實現(xiàn),三節(jié)點-二次插值,三節(jié)點-逐次線性插值,function chazhi1 %一次插值 x0=pi/6;y0=1/2; x1=pi/4;y1=1/sqrt(2); x=pi/180*50; y=chazhi12(x0,y0,x1,y1,x) error=y-sin(x) function y=chazhi11(x0,y0,x1,y1,x)%拉格朗

7、日插值法 y=(x-x1)/(x0-x1)*y0+( x-x0)/(x1-x0)*y1; function y=chazhi12(x0,y0,x1,y1,x)%牛頓法 y=y0+(y1-y0)/(x1-x0)*( x-x0);,解:(一次插值),function chazhi2 %二次插值 x0=pi/6;y0=1/2; x1=pi/4;y1=1/sqrt(2); x2=pi/3;y2=sqrt(3)/2; x=pi/180*50; y=chazhi22(x0,y0,x1,y1,x2,y2,x)%選擇插值函數(shù) error=y-sin(x) function y=chazhi21(x0,y0,x

8、1,y1,x2,y2,x)%拉格朗日插值法 y=(x-x1)*(x-x2)/(x0-x1)*(x0-x2)*y0+(x-x0)*(x-x2)/. (x1-x0)*(x1-x2)*y1+(x-x0)*(x-x1)/(x2-x0)*(x2-x1)*y2; function y=chazhi22(x0,y0,x1,y1,x2,y2,x)%牛頓插值法 y=y0+(y1-y0)/(x1-x0)*( x-x0)+(x-x0)*(x-x1)*. (y2-y1)/(x2-x1)-(y1-y0)/(x1-x0)/(x2-x0);,解:(二次插值),N+1節(jié)點-N次插值,Lagrange插值基函數(shù),（2）,與 節(jié)

9、點有關(guān)，而與f 無關(guān),這里每個lj(x)都是n次多項式，且由(1)式容易驗證lj(x）滿足,j=0,1,，n （1）,對任意的pn(x)Pn，都有pn(x)=c0 l0(x)+c1 l1(x)+cnln(x)其中c0,c1,cn為組合系數(shù)。,可以證明函數(shù)組l0(x)，l1(x)，, ln(x) 在插值區(qū)間a,b上線性無關(guān)，所以這n+1個函數(shù)可作為Pn的一組基函數(shù)，稱為Lagrange插值基函數(shù)。,Lagrange插值基函數(shù),由Lagrange插值基函數(shù)滿足(2)式可知，方程組變成,因此得到插值多項式 pn(x)= f(x0)l0(x)+f(x1) l1(x)+ f(xn) ln(x),記為Ln

10、(x) f(xj)lj(x),稱Ln(x)為n次Lagrange插值多項式,插值余項 /* Remainder */,定理 若,在a , b內(nèi)存在, 則在a , b上,的n+1個互異的點，對 f（x)所作的n次Lagrange插值多項式Ln (x) 有誤差估計,Rolles Theorem的推論: 若 充分光滑，且,證明：由于Rn(xi) 0 ，i=0,1,，n,任意固定 x xi (i = 0, , n), 考察,(t)有 n+2 個不同的根 x0 xn x,解：,n = 1,分別利用x0, x1 以及 x1, x2 計算,利用,sin 50 = 0.7660444,利用x0, x1 作為插

11、值節(jié)點的實際誤差 0.01001,利用,計算得：sin 50 0.76008,利用x1, x2作為插值節(jié)點的實際誤差 0.00596,n = 2,sin 50 = 0.7660444,2次插值的實際誤差 0.00061,function main %N次插值 x0=pi/6;y0=1/2; x1=pi/4;y1=1/sqrt(2); x2=pi/3;y2=sqrt(3)/2; xy=x0 x1 x2;y0 y1 y2; x=pi/180*50; y=chazhiN1(xy,x) error=y-sin(x),Matlab 編程實現(xiàn),N+1節(jié)點-N次插值,function y=chazhiN1(

12、xy,x) %拉格朗日插值法 y=0; N=size(xy,2); for i=1:N L=1; for j=1:N if i=j L=L*(x-xy(1,j)/(xy(1,i)-xy(1,j); end end y=L*xy(2,i)+y; end,牛頓插值(Newtons Interpolation ),Lagrange 插值雖然易算，但若要增加一個節(jié)點時，全部基函數(shù) li(x) 都需要重新計算。,能否重新在Pn中尋找新的基函數(shù) ？,希望每加一個節(jié)點時，只附加一項上去即可。,1,x - x0, (x - x0)(x - x1) ，(x-x0)(x-x1) (x-xn-1)是否構(gòu)成Pn的一組

13、基函數(shù)？,利用插值條件Nn(xj)=f(xj), j=0,1，n代入上式，得關(guān)于Ak (k=0,1，n)的線性代數(shù)方程組,牛頓插值法的基函數(shù),當(dāng)xj 互異時，系數(shù)矩陣非奇異，且容易求解,基函數(shù),How complex the expression are!, 差商(亦稱均差) /* divided difference */,稱為在xi,xj處的1階差商,稱為在xi,xj,xk處的2階差商,k階差商：,利用插值條件和差商,可求出Nn(x)的系數(shù) Ai :,Newton插值多項式,因此，每增加一個結(jié)點，Newton插值多項式只增加一項，克服了Lagrange插值的缺點。,Newton插值多項式,

14、. xk f(xk) 一階差商 二階差商 三階差商 n 階差商,差商表,例1:給定f(x)=lnx的數(shù)據(jù)表 xi 2.20 2.40 2.60 2.80 3.00 f(xi) 0.78846 0.87547 0.95551 1.02962 1.09861 1.構(gòu)造差商表 2.分別寫出二次、四次Newton插值多項式 解:差商表,N2(x)=0.78846 +0.43505(x-2.20) - 0.087375(x-2.20) (x-2.40） N4(x)= 0.78846 +0.43505(x-2.20） - 0.087375(x-2.20)(x-2.40) +0.0225(x-2.20)(x-2.40)(x-2.60) -0.00755(x-2.20)(x-2.40)(x-2.60)(x-2.80),function main %N次插值 x0=pi/6;y0=1/2; x1=pi/4;y1=1/sqrt(2); x2=pi/3;y2=sqrt(3)/2; xy=x0 x1 x2;y0 y1 y2; x=pi/180*50; y=chazhiN1(xy,x) error=y-sin(x),Matlab 編程實現(xiàn),N+1節(jié)點-N次插值,function y=chazhiN2(xy,x)%牛頓法 N