数值分析作业

1、数值分析实验作业 指导老师：院 系： 姓 名： 学 号： 第二章 用matlab编写的程序如下：function t_charpt2 %数值实验二：含“实验2.1：多项式插值的震荡现象”和“实验2.2：样条插值的收敛”%输入：实验选择，函数式选择，插值结点数%输出：拟合函数及原函数的图形result=inputdlg('请选择实验，若选2.1，请输入1，否则输入2：','charpt2',1,'1');Nb=str2num(char(result);if(Nb=1)&(Nb=2)errordlg('实验选择错误！');re

2、turn;end promps='请选择实验函数，若选f(x),请输入f，若选h(x),请输入h,若选g(x),请输入g:' titles='charpt2' result=inputdlg(promps,'charpt2',1,'f'); Nb_f=char(result); if(Nb_f='f'&Nb_f='h'&Nb_f='g')errordlg('实验函数选择错误！');return;endresult=inputdlg('请输入插值

3、结点数N：','charpt2',1,'10');Nd=str2num(char(result);if(Nd<1)errordlg('结点输入错误！');return;end switch Nb_f case 'f' f=inline('1./(1+25*x.2)');a=-1;b=1; case 'h' f=inline('x./(1+x.4)');a=-5;b=5; case 'g' f=inline('atan(x)');a=-5;

4、b=5; end if(Nb=1) x0=linspace(a,b,Nd+1);y0=feval(f,x0); x=a:0.1:b;y=Lagrange(x0,y0,x); fplot(f,a b,'co'); hold on; plot(x,y,'b-'); xlabel('x');ylabel('y=f(x) o and y=Ln(x)-'); elseif(Nb=2) x0=linspace(a,b,Nd+1);y0=feval(f,x0); x=a:0.1:b; cs=spline(x0,y0);y=ppval(cs,x)

5、; plot(x0,y0,'o');hold on ;plot(x,y,'k-'); xlabel('x');ylabel('y=f(x) o and y=Spline(x)-'); end % function y=Lagrange(x0,y0,x) n=length(x0);m=length(x); for i=1:m; z=x(i); s=0.0; for k=1:n p=1.0; for j=1:n if(j=k) p=p*(z-x0(j)/(x0(k)-x0(j); end end s=s+p*y0(k); end y(

6、i)=s;end若选择实验2.1，实验函数为f，插值结点数为6，则结果为：若选择实验2.1，实验函数为f，插值结点数为20，则结果为：若选择实验2.1，实验函数为f，插值结点数为27，则结果为：选择其它的函数重复上述的实验，在这里我选择的是h函数，具体结果如下：若选择实验2.1，实验函数为h，插值结点数为6，则结果为：若选择实验2.1，实验函数为h，插值结点数为20，则结果为：若选择实验2.1，实验函数为h，插值结点数为27，则结果为：实验2.2输入2，实验函数为f，插值节点数为20，则结果为：输入2，实验函数为h，插值节点数为20，则结果为：分析：1.对于n次Lagrange插值多项式Pn（

7、x）近似逼近函数f（x）时，Pn（x）的次数并非越高，越逼近目标函数f（x）。在实验中选取一个较小数6，可以看到Pn（x）图像与f（x）吻合情况不是很好，当选取n=20时，Pn（x）与f（x）在区间-1,1上的吻合情况较好；但是当选取较大数27，Pn（x）开始发生振荡，吻合情况较差，也即“龙格”现象。在选取的函数h（x）不会发生类似的现象，精度会随着次数的增加而增加，不会出现所谓的“龙格”现象。2.从样条插值的实验结果可以看到：样条插值的逼近情况比Lagrange插值要好，并且Lagrange插值在区间-5，5上端点处会有一定的波动，但是样条插值结果与精确值吻合较好。第三章 实验3用matla

8、b编写的程序如下：function charpt3 %数值实验三：含“实验3.1”和“实验3.2”%子函数调用：dlsa%输入：实验选择%输出：原函数及求得的相应插值多项式的函数的图像以及参数alph和误差rresult=inputdlg('请选择实验，若选3.1，请输入1，否则输入2：','charpt3',1,'1');Nb=str2num(char(result);if(Nb=1)&(Nb=2)errordlg('实验选择错误！');return;end x0=-1:0.5:2; y0=-4.447 -0.452 0

9、.551 0.048 -0.447 0.549 4.552; if(Nb=1) n=3; %n为拟合阶次 alph=polyfit(x0,y0,n); y=polyval(alph,x0); r=(y0-y)*(y0-y)' x=-1:0.01:2; y=polyval(alph,x); plot(x,y,'k-'); xlabel('x');ylabel('y0 * and ployfit.y-'); hold on plot(x0,y0,'*') grid on; else result=inputdlg('请

10、输入权向量w:','charpt3',1,'1 1 1 1 1 1 1'); w=str2num(char(result); n=3; a,b,c,alph,r=dlsa(x0,y0,w,n); end disp('平方误差:',num2str(r) disp('参数alph:',num2str(alph) function a,b,c,alph,r=dlsa(x,y,w,n) %功能：用正交化方法对离散数据作多项式最小二乘拟合。 %输入：m+1个离散点（x,y,w),x,y,w分别用行向量给出。 % 拟合多项式的次数n,

11、0<n<m. % 平方误差r=(y-s,y-s),并作离散点列和拟合曲线的图形 m=length(x)-1; if(n<1|n>=m)errordlg('错误：n<1或者n>=m!');return;end%求三项递推公式的参数a,b,拟合多项式s(x)的系数c,其中d(k)=(y,sk); s1=0;s2=ones(1,m+1);v2=sum(w); d(1)=y*w'c(1)=d(1)/v2; for k=1:n xs=x.*s2.2*w'a(k)=xs/v2; if k=1 b(k)=0; else b(k)=v2/v1

12、; end s3=(x-a(k).*s2-b(k)*s1; v3=s3.2*w' d(k+1)=y.*s3*w'c(k+1)=d(k+1)/v3; s1=s2;s2=s3;v1=v2;v2=v3; end %求平方误差r r=y.*y*w'-c*d'%求拟合多项式s(x)的降幂系数alph alph=zeros(1,n+1);T=zeros(n+1,n+2); T(:,2)=ones(n+1,1);T(2,3)=-a(1); if n>=2 for k=3:n+1 for i=3:k+1 T(k,i)=T(k-1,i)-a(k-1)*T(k-1,i-1)-

13、b(k-1)*T(k-2,i-2); end end end for i=1:n+1 for k=i:n+1 alph(n+2-i)=alph(n+2-i)+c(k)*T(k,k+2-i); end end %用秦九韶方法计算s(t)的输出序列（t,s) xmin=min(x);xmax=max(x);dx=(xmax-xmin)/(25*m); t=(xmin-dx):dx:(xmax+dx); s=alph(1); for k=2:n+1 s=s.*t+alph(k); end%输出点列x-y和拟合曲线t-s的图形 plot(x,y,'*',t,s,'-')

14、; title('离散数据的多项式拟合'); xlabel('x');ylabel('y');grid on;实验3.1，输入1，权函数为1 1 1 1 1 1 1，结果为：>> 平方误差:2.1762e-005参数alph:1.9991 -2.9977-3.9683e-005 0.54912图形为：实验3.2：输入2，权函数为1 1 1 1 1 1 1，结果为：>> 平方误差:2.1762e-005参数alph:1.9991 -2.9977-3.9683e-005 0.54912图形为：分析：利用最小二乘法作曲线的拟合，

15、对实验3.1给出的数据作的三次多项式的图形和数据节点拟合得较好。正交化多项式，取w=1 1 1 1 1 1 1时，曲线的拟合也是十分吻合的。由于没有出现病态法方程组的情况，二者图形结果及平方误差也相差不大。第四章 实验4用matlab编写的程序如下：function t_charpt4%数值实验四：含“实验4.1：复化求积公式计算定积分”和“实验4.2：高斯数值积分法用于积分方程求解”%子函数调用：CG_L_I function(复化Gauss_Legendre I)公式、CTrapezia(复化梯形公式）、CSimpson(复化Simpson公式）%输入：实验选择、积分法选择、积分式题号选择

16、%输出：积分（实验4.1）或方程解（实验4.2）的精确值和数值解及误差 result=inputdlg('请选择实验,若选4.1,请输入1,否则输入2:','charpt4',1,'1'); Nt=str2num(char(result); if(Nt=1)&(Nt=2)errordlg('实验选择错误！');return;end promps='请选择积分法,若用复化梯形,输入T,用复化Simpson,输入S,用复化Gauss_Legendre,输入GL:' result=inputdlg(promps,

17、'charpt4',1,'T'); Nb=char(result); if(Nb='T' & Nb='S' & Nb='GL')errordlg('积分公式选择错误！');return;end if(Nt=1) result=inputdlg('请输入积分式题号14:','实验4.1',1,'1'); Nb_f=str2num(char(result); if(Nb_f<1)|(Nb_f>4)errordlg('没有

18、该积分式!');return;end switch Nb_f case 1 fun=inline('-2./(x.2-1)');a=2;b=3; case 2 fun=inline('4./(x.2+1)');a=0;b=1; case 3 fun=inline('3.x');a=0;b=1; case 4 fun=inline('x.*exp(x)');a=1;b=2; end tol=0.5e-7;h=0.01; if(Nb='T')%用复化梯形公式 t=(fun(a)+fun(b)*(b-a)/2;

19、k=1;t0=0; while(abs(t-t0)>=tol*3) t0=t;h=(b-a)/2k; t=t0/2+h*sum(fun(a+h:2*h:b-h); k=k+1; end elseif(Nb='S')%用复化Simpson公式 t=quad(fun,a,b,tol); elseif(Nb='GL')%用复化Gauss_Legendre I N=floor(b-a)/h);t=0;xk=0; for k=0:N xk=a+k*h+h/2; t=t+fun(xk-h/(2*sqrt(3)+fun(xk+h/(2*sqrt(3); end t=t*

20、h/2; end elseif(Nt=2) result=inputdlg('请输入方程式题号1或2:','实验4.2',1,'1'); Nb_f=str2num(char(result); if(Nb_f=1 & Nb_f=2)errordlg('没有该方程式!');return;end result=inputdlg('请输入步长:','实验4.2',1,'0.01'); h=str2num(char(result); if(h<=0)errordlg('请

21、输入正确的步长!');return;end if(Nb='T')%用复化梯形公式 x,t=CTrapezia(0,1,h,Nb_f); elseif(Nb='S')%用复化Simpson公式 x,t=CSimpson(0,1,h,Nb_f); elseif(Nb='GL')%用复化Gauss_Legendre I公式 x,t=CG_L_I(0,1,h,Nb_f); end plot(x,t,'g-'); xlabel('x');ylabel('y'); title('积分方程求解&#

22、39;) hold on disp('实验4.2(',num2str(Nb_f),')的计算结果:',num2str(t'); if(Nb_f=1) fplot('exp(x)',0 1,'*'); hold off disp('实验4.2题（1）的精确解:',num2str(exp(x'); disp('绝对误差和:',num2str(sum(abs(t'-exp(x'); else fplot('1/(1+t)2',0 1,'*')

23、; hold off y=1./(x+1).2; disp('实验4.2题（2）的精确解:',num2str(exp(y'); disp('绝对误差和:',num2str(sum(abs(t'-y'); end end if(Nt=1) disp('实验4.1题(',num2str(Nb_f),')的计算结果:',num2str(t); switch Nb_f case 1 disp('精确解:ln2-ln3=-0.4054651081') disp('绝对误差:',num2

24、str(abs(t+0.4054651081); case 2 disp('精确解:pi=3.14159265358979') disp('绝对误差:',num2str(abs(t-pi); case 3 disp('精确解:2/ln3=1.82047845325368') disp('绝对误差:',num2str(abs(t-1.82047845325368); case 4 disp('精确解:e2=7.38905609893065') disp('绝对误差:',num2str(abs(t-7.

25、38905609893065); endend %function x,y=CG_L_I(a,b,h,N)%复化Gauss_Legendre I,用于积分方程求解%输入：a、b分别为求积下、上限，h为步长，N为方程式题号实验选择、积分法选择、积分式题号选择%输出：x,y为方程离散解y(i)=f(x(i) n=floor(b-a)/h); A=zeros(2*n+2,2*n+2);b=zeros(2*n+2,1);x=b; for i=0:n t1=a+h*(i+0.5-sqrt(3)/6);t2=a+h*(i+0.5+sqrt(3)/6); x(2*i+1)=t1;x(2*i+2)=t2; i

26、f(N=1) b(2*i+1)=exp(t1);b(2*i+2)=exp(t2); else b(2*i+1)=-(4*t1.3+5*t1.2-2*t1+5)/(8*(t1+1).2); b(2*i+2)=-(4*t2.3+5*t2.2-2*t2+5)/(8*(t2+1).2); end for j=0:n x1=a+h*(j+0.5-sqrt(3)/6);x2=a+h*(j+0.5+sqrt(3)/6); if(N=1) A(2*i+1,2*j+1)=exp(t1)*2/(exp(1)-1); A(2*i+1,2*j+2)=exp(t1)*2/(exp(1)-1); A(2*i+2,2*j+

27、1)=exp(t2)*2/(exp(1)-1); A(2*i+2,2*j+2)=exp(t2)*2/(exp(1)-1); else A(2*i+1,2*j+1)=1/(1+x1)-t1; A(2*i+1,2*j+2)=1/(1+x2)-t1; A(2*i+2,2*j+1)=1/(1+x1)-t2; A(2*i+2,2*j+2)=1/(1+x2)-t2; end end endA=h*A/2-eye(2*n+2);y=inv(A)*b; %function x,y=CTrapezia(a,b,h,N)%复化Gauss_%复化梯形公式，用于积分方程求解%输入：a、b分别为求积下、上限，h为步长，

28、N为方程式题号实验选择、积分法选择、积分式题号选择%输出：x,y为方程离散解y(i)=f(x(i) n=floor(b-a)/h); A=zeros(n+1,n+1);b=zeros(n+1,1);x=b; for i=0:n t=a+i*h;x(i+1)=t; if(N=1) b(i+1)=-exp(t); else b(i+1)=(4*t.3+5*t.2-2*t+5)/(8*(t+1).2); end for j=0:n s=a+j*h; if(j=0)|(j=n) if(N=1) A(i+1,j+1)=exp(t)*2/(exp(1)-1); else A(i+1,j+1)=1/(1+s

29、)-t; end else if(N=1) A(i+1,j+1)=exp(t)*2/(exp(1)-1)*2; else A(i+1,j+1)=(1/(1+s)-t)*2; end end end endA=A-eye(n+1)*2/h;b=-2*b/h;y=inv(A)*b; %function x,y=CSimpson(a,b,h,N)%复化Simpson公式，用于积分方程求解%输入：a、b分别为求积下、上限，h为步长，N为方程式题号实验选择、积分法选择、积分式题号选择%输出：x,y为方程离散解y(i)=f(x(i) n=floor(b-a)/h);h1=h/2; A=zeros(2*n+

30、1,2*n+1);b=zeros(2*n+1,1);x=b; for i=0:n*2 t=a+i*h1;x(i+1)=t; if(N=1) b(i+1)=-exp(t); else b(i+1)=(4*t.3+5*t.2-2*t+5)/(8*(t+1).2); end for j=0:n*2 s=a+j*h1; if(j=0)|(j=n*2) if(N=1) A(i+1,j+1)=exp(t)*2/(exp(1)-1); else A(i+1,j+1)=1/(1+s)-t; end elseif(mod(j,2)=0) if(N=1) A(i+1,j+1)=exp(t)*2/(exp(1)-1

31、)*2; else A(i+1,j+1)=(1/(1+s)-t)*2; end elseif(mod(j,2)=0) if(N=1) A(i+1,j+1)=exp(t)*2/(exp(1)-1)*4; else A(i+1,j+1)=(1/(1+s)-t)*4; end end end end A=A-eye(2*n+1)*6/h;b=-6*b/h; y=inv(A)*b;1.若选择实验4.1，1）积分法选择T(即复化梯形公式)：>> 实验4.1题(1)的计算结果:-0.40547精确解:ln2-ln3=-0.4054651081绝对误差:1.3944e-008>> 实

32、验4.1题(2)的计算结果:3.1416精确解:pi=3.14159265358979绝对误差:3.9736e-008>> 实验4.1题(3)的计算结果:1.8205精确解:2/ln3=1.82047845325368绝对误差:4.3655e-008>> 实验4.1题(4)的计算结果:7.3891精确解:e2=7.38905609893065绝对误差:2.0775e-0082）积分法选择S(即复化Simpson公式)：>> 实验4.1题(1)的计算结果:-0.40547精确解:ln2-ln3=-0.4054651081绝对误差:1.2625e-009>

33、> 实验4.1题(2)的计算结果:3.1416精确解:pi=3.14159265358979绝对误差:2.7517e-010>> 实验4.1题(3)的计算结果:1.8205精确解:2/ln3=1.82047845325368绝对误差:1.0877e-010>> 实验4.1题(4)的计算结果:7.3891精确解:e2=7.38905609893065绝对误差:8.2168e-0113）积分法选择GL(即复化Gauss_Legendre公式):>> 实验4.1题(1)的计算结果:-0.40796精确解:ln2-ln3=-0.4054651081绝对误差:0

34、.0024907>> 实验4.1题(2)的计算结果:3.1615精确解:pi=3.14159265358979绝对误差:0.0199>> 实验4.1题(3)的计算结果:1.8506精确解:2/ln3=1.82047845325368绝对误差:0.030165>> 实验4.1题(4)的计算结果:7.538精确解:e2=7.38905609893065绝对误差:0.148892.若选择实验4.2，1）积分法选择T(即复化梯形公式),题号选择1，积分步长为0.01，则结果如下：>> 实验4.2(1)的计算结果:0.99998 1.01 1.0202 1

35、.0304 1.0408 1.0513 1.0618 1.0725 1.0833 1.0942 1.1052 1.1163 1.1275 1.1388 1.1503 1.1618 1.1735 1.1853 1.1972 1.2092 1.2214 1.2337 1.2461 1.2586 1.2712 1.284 1.2969 1.3099 1.3231 1.3364 1.3498 1.3634 1.3771 1.3909 1.4049 1.419 1.4333 1.4477 1.4623 1.477 1.4918 1.5068 1.5219 1.5372 1.5527 1.5683 1.5

36、84 1.6 1.616 1.6323 1.6487 1.6653 1.682 1.6989 1.716 1.7332 1.7506 1.7682 1.786 1.804 1.8221 1.8404 1.8589 1.8776 1.8964 1.9155 1.9348 1.9542 1.9738 1.9937 2.0137 2.034 2.0544 2.075 2.0959 2.117 2.1382 2.1597 2.1814 2.2034 2.2255 2.2479 2.2705 2.2933 2.3163 2.3396 2.3631 2.3869 2.4109 2.4351 2.4596

37、2.4843 2.5092 2.5345 2.5599 2.5857 2.6117 2.6379 2.6644 2.6912 2.7182实验4.2题（1）的精确解:1 1.0101 1.0202 1.0305 1.0408 1.0513 1.0618 1.0725 1.0833 1.0942 1.1052 1.1163 1.1275 1.1388 1.1503 1.1618 1.1735 1.1853 1.1972 1.2092 1.2214 1.2337 1.2461 1.2586 1.2712 1.284 1.2969 1.31 1.3231 1.3364 1.3499 1.3634 1

38、.3771 1.391 1.4049 1.4191 1.4333 1.4477 1.4623 1.477 1.4918 1.5068 1.522 1.5373 1.5527 1.5683 1.5841 1.6 1.6161 1.6323 1.6487 1.6653 1.682 1.6989 1.716 1.7333 1.7507 1.7683 1.786 1.804 1.8221 1.8404 1.8589 1.8776 1.8965 1.9155 1.9348 1.9542 1.9739 1.9937 2.0138 2.034 2.0544 2.0751 2.0959 2.117 2.138

39、3 2.1598 2.1815 2.2034 2.2255 2.2479 2.2705 2.2933 2.3164 2.3396 2.3632 2.3869 2.4109 2.4351 2.4596 2.4843 2.5093 2.5345 2.56 2.5857 2.6117 2.6379 2.6645 2.6912 2.7183绝对误差和:0.00289482）积分法选择S(即复化Simpson公式),题号选择1，积分步长为0.01，则结果如下>> 实验4.2(1)的计算结果:1 1.005 1.0101 1.0151 1.0202 1.0253 1.0305 1.0356 1

40、.0408 1.046 1.0513 1.0565 1.0618 1.0672 1.0725 1.0779 1.0833 1.0887 1.0942 1.0997 1.1052 1.1107 1.1163 1.1219 1.1275 1.1331 1.1388 1.1445 1.1503 1.156 1.1618 1.1677 1.1735 1.1794 1.1853 1.1912 1.1972 1.2032 1.2092 1.2153 1.2214 1.2275 1.2337 1.2399 1.2461 1.2523 1.2586 1.2649 1.2712 1.2776 1.284 1.2

41、905 1.2969 1.3034 1.31 1.3165 1.3231 1.3298 1.3364 1.3431 1.3499 1.3566 1.3634 1.3703 1.3771 1.384 1.391 1.3979 1.4049 1.412 1.4191 1.4262 1.4333 1.4405 1.4477 1.455 1.4623 1.4696 1.477 1.4844 1.4918 1.4993 1.5068 1.5144 1.522 1.5296 1.5373 1.545 1.5527 1.5605 1.5683 1.5762 1.5841 1.592 1.6 1.608 1.

42、6161 1.6242 1.6323 1.6405 1.6487 1.657 1.6653 1.6736 1.682 1.6905 1.6989 1.7074 1.716 1.7246 1.7333 1.7419 1.7507 1.7594 1.7683 1.7771 1.786 1.795 1.804 1.813 1.8221 1.8313 1.8404 1.8497 1.8589 1.8682 1.8776 1.887 1.8965 1.906 1.9155 1.9251 1.9348 1.9445 1.9542 1.964 1.9739 1.9838 1.9937 2.0037 2.0138 2.0238 2.034 2.0442 2.0544 2.0647 2.0751 2.0855 2.0959 2.1064 2.117 2.1276 2.1383 2.149 2.1598 2.1706 2.1815 2.1924 2.2034 2.2144 2.2255 2.2367 2.2479 2.2592 2.2705 2.2819 2.2933 2.3048 2.3164 2.328 2.3396 2.3514 2.3632 2.375 2.3869 2.3989 2.4109 2.423