样条插值及应用深入研究

学院： 研究生学院 专业： 机械工程 组号： 39 成绩：

值设计车门的曲线，其中一小段的数据如下：

x f(x) f'(x) 0 0.0 0.8 1 0.79 2 1.53 3 2.19 4 2.71 5 3.03 6 3.27 7 2.89 8 3.06 9 3.19 10 3.29 0.2

Mathematics程序： Clear[x,y,a,b,c,M] x[i_]:=i-1; y[1]=0; y[2]=0.79; y[3]=1.53; y[4]=2.19; y[5]=2.71; y[6]=3.03; y[7]=3.27; y[8]=2.89; y[9]=3.06; y[10]=3.19; y[11]=3.29;

B=Table[{x[i],y[i]},{i,1,11}] y'[1]=0.8; y'[11]=0.2; h[j_]:=1;

a[j_]:=h[j-1]/(h[j-1]+h[j]); a[11]=1; b[1]=1; b[j_]:=1-a[j];

c[1]=6/h[1]((y[2]-y[1])/h[1]-y'[1]);

-45-

学院： 研究生学院 专业： 机械工程 组号： 39 成绩：

c[j_]:=6((y[j+1]-y[j])/h[j]-(y[j]-y[j-1])/h[j-1])/(h[j-1]+h[j]); c[11]=6/h[11-1](y'[11]-(y[11]-y[11-1])/h[11-1]);

A=Table[Switch[i-j,-1,b[j-1],0,2,1,a[j+1],_,0],{i,1,11},{j,1,11}]; MatrixForm[%] CC=Table[c[j],{j,1,11}]; MatrixForm[%] LinearSolve[A,CC]; MatrixForm[%];

M[j_]:=LinearSolve[A,CC][[j]] Table[M[j],{j,1,11}]

S[j_]:=M[j+1](x-x[j])^3/(6h[j])-M[j](x-x[j+1])^3/(6h[j])+ (y[j+1]-M[j+1]h[j]^2/6)(x-x[j])/h[j]- (y[j]-M[j]h[j]^2/6)(x-x[j+1])/h[j] Table[S[j],{j,1,10}]; Expand[%]; MatrixForm[%] g1=Plot[%[[1]],{x,0,1}] g2=Plot[%%[[2]],{x,1,2}] g3=Plot[%%%[[3]],{x,2,3}] g4=Plot[%%%%[[4]],{x,3,4}] g5=Plot[%%%%%[[5]],{x,4,5}] g6=Plot[%%%%%%[[6]],{x,5,6}] g7=Plot[%%%%%%%[[7]],{x,6,7}] g8=Plot[%%%%%%%%[[8]],{x,7,8}] g9=Plot[%%%%%%%%%[[9]],{x,8,9}] g10=Plot[%%%%%%%%%%[[10]],{x,9,10}] g11=ListPlot[B,Prolog->AbsolutePointSize[15]]

Show[g1,g2,g3,g4,g5,g6,g7,g8,g9,g10,Prolog->AbsolutePointSize[15]]

运行结果如下：

-46-

学院： 研究生学院 专业： 机械工程 组号： 39 成绩：

4、已知函数在下列各点的值为

xi f(xi) 0.2 0．98 0.4 0.92 0.6 0.81 0.8 0.64 1.0 0.38 使用三次样条函数(自然边界条件)对数据进行插值。

Mathematics程序： Clear[x,y,a,b,c,M] x[i_]:=0.2(i+1); y[0]=0.98; y[1]=0.92; y[2]=0.81; y[3]=0.64; y[4]=0.38;

B=Table[{x[i],y[i]},{i,0,4}]; y''[0]=0; y''[4]=0; h[j_]:=0.2;

a[j_]:=h[j-1]/(h[j-1]+h[j]); b[j_]:=1-a[j];

-47-

学院： 研究生学院 专业： 机械工程 组号： 39 成绩：

c[1]=6((y[2]-y[1])/h[1]-(y[1]-y[0])/h[0])/(h[0]+h[1])-h[0]y''[0]/(h[0]+h[1]); c[j_]:=6((y[j+1]-y[j])/h[j]-(y[j]-y[j-1])/h[j-1])/(h[j-1]+h[j]);

c[3]=6((y[4]-y[3])/h[3]-(y[3]-y[2])/h[2])/(h[2]+h[3])-h[3]y''[4]/(h[2]+h[3]); A=Table[Switch[i-j,-1,b[j-1],0,2,1,a[j+1],_,0],{i,1,3},{j,1,3}]; MatrixForm[%] CC=Table[c[j],{j,1,3}]; MatrixForm[%] LinearSolve[A,CC]; MatrixForm[%];

M[j_]:=LinearSolve[A,CC][[j]] Table[M[j],{j,1,3}] M[0]=0; M[4]=0;

S[j_]:=M[j+1](x-x[j])^3/(6h[j])-M[j](x-x[j+1])^3/(6h[j])+ (y[j+1]-M[j+1]h[j]^2/6)(x-x[j])/h[j]- (y[j]-M[j]h[j]^2/6)(x-x[j+1])/h[j] Table[S[j],{j,0,3}]; Expand[%]; MatrixForm[%]

g1=Plot[%[[1]],{x,0.2,0.4}] g2=Plot[%%[[2]],{x,0.4,0.6}] g3=Plot[%%%[[3]],{x,0.6,0.8}] g4=Plot[%%%%[[4]],{x,0.8,1.0}]

g5=ListPlot[B,Prolog->AbsolutePointSize[15]] Show[g1,g2,g3,g4,g5,Prolog->AbsolutePointSize[15]]

运行结果如下：

-48-