1、1MATLAB 编程题库1.下面的数据表近似地满足函数 ,请适当变换成为线性最小二乘问题,编程求21cxbay最好的系数 ,并在同一个图上画出所有数据和函数图像.cba, 625.07180823.0.687.0.356. 954391iyx解:x=-0.931 -0.586 -0.362 -0.213 0.008 0.544 0.628 0.995;y=0.356 0.606 0.687 0.802 0.823 0.801 0.718 0.625;A=x ones(8,1) -x.2.*y;z=Ay;a=z(1); b=z(2); c=z(3);xh=-1:0.1:1;yh=(a.*xh+b
2、)./(1+c.*xh.2);plot(x,y,r+,xh,yh,b*)22.若在 Matlab 工作目录下已经有如下两个函数文件,写一个割线法程序,求出这两个函数精度为 的近似根,并写出调用方式:10文件一 文件二function v = f(x)v = x .* log(x) - 1;function z = g(y)z = y.5 + y - 1;解: edit gexianfa.mfunction x iter=gexianfa(f,x0,x1,tol)iter=0;while(norm(x1-x0)tol) iter=iter+1;x=x1-feval(f,x1).*(x1-x0).
3、/(feval(f,x1)-feval(f,x0);x0=x1;x1=x;end edit f.mfunction v=f(x)v=x.*log(x)-1; edit g.mfunction z=g(y)z=y.5+y-1; x1 iter1=gexianfa(f,1,3,1e-10)x1 =1.7632iter1 =6 x2 iter2=gexianfa(g,0,1,1e-10)x2 =0.7549iter2 =833.使用 GS 迭代求解下述线性代数方程组: 1235140xx+=-解: edit gsdiedai.mfunction x iter=gsdiedai(A,x0,b,tol)
4、D=diag(diag(A);L=D-tril(A);U=D-triu(A);iter=0;x=x0;while(norm(b-A*x)./norm(b)tol)iter=iter+1;x0=x;x=(D-L)(U*x0+b);end A=5 2 1;-1 4 2;1 -3 10; b=-12 10 3;tol=1e-4;x0=0 0 0; x iter=gsdiedai(A,x0,b,tol);xx =-3.09101.23720.9802iteriter =644.用四阶 Range-kutta 方法求解下述常微分方程初值问题(取步长 h=0.01),(1)2xdyey=+解: edit
5、ksf2.mfunction v=ksf2(x,y)v=y+exp(x)+x.*y; a=1;b=2;h=0.01; n=(b-a)./h; x=1:0.01:2;y(1)=2;fori=2:(n+1)k1=h*ksf2(x(i-1),y(i-1);k2=h*ksf2(x(i-1)+0.5*h,y(i-1)+0.5*k1);k3=h*ksf2(x(i-1)+0.5*h,y(i-1)+0.5*k2);k4=h*ksf2(x(i-1)+h,y(i-1)+k3);y(i)=y(i-1)+(k1+2*k2+2*k3+k4)./6;endy调用函数方法 edit Rangekutta.mfunction
6、 x y=Rangekutta(f,a,b,h,y0)x=a:h:b;n=(b-a)/h;y(1)=y0;fori=2:(n+1)k1=h*(feval(f,x(i-1),y(i-1);k2=h*(feval(f,x(i-1)+0.5*h,y(i-1)+0.5*k1);k3=h*(feval(f,x(i-1)+0.5*h,y(i-1)+0.5*k2);k4=h*(feval(f,x(i-1)+h,y(i-1)+k3);y(i)=y(i-1)+(k1+2*k2+2*k3+k4)./6;end x y=Rangekutta(ksf2,1,2,0.01,2);y55.取 ,请编写 Matlab 程序
7、,分别用欧拉方法、改进欧拉方法在 上求解初0.2h 12x值问题。 3,(1)0.4dyx=-解: edit Euler.mfunction x y=Euler(f,a,b,h,y0)x=a:h:b;n=(b-a)./h;y(1)=y0;fori=2:(n+1)y(i)=y(i-1)+h*feval(f,x(i-1),y(i-1);end edit gaijinEuler.mfunctionx y=gaijinEuler(f,a,b,h,y0)x=a:h:b;n=(b-a)./h;y(1)=y0;fori=2:(n+1)y1=y(i-1)+h*feval(f,x(i-1),y(i-1);y2=
8、y(i-1)+h*feval(f,x(i),y1);y(i)=(y1+y2)./2;end edit ksf3.mfunction v=ksf3(x,y)v=x.3-y./x;x y=Euler(ksf3,1,2,0.2,0.4)x =1.0000 1.2000 1.4000 1.6000 1.8000 2.0000y =0.4000 0.5200 0.7789 1.2165 1.8836 2.8407 x y=gaijinEuler(ksf3,1,2,0.2,0.4)x =1.0000 1.2000 1.4000 1.6000 1.8000 2.0000y =0.4000 0.5895 0.
9、9278 1.4615 2.2464 3.346666.请编写复合梯形积分公式的 Matlab 程序,计算下面积分的近似值,区间等分 。20n编写辛普森积分公式的 Matlab 程序,计算下面积分的近似值,区间等分 。1、120dx+1sindx-解: edit tixingjifen.mfunction s=tixingjifen(f,a,b,n)x=linspace(a,b,(n+1);y=zeros(1,length(x);y=feval(f,x)h=(b-a)./n;s=0.5*h*(y(1)+2*sum(y(2:n)+y(n+1);end edit simpson.mfunction
10、 I=simpson(f,a,b,n)h=(b-a)/n;x=linspace(a,b,2*n+1);y=feval(f,x);I=(h/6)*(y(1)+2*sum(y(3:2:2*n-1)+4*sum(y(2:2:2*n)+y(2*n+1); edit ksf4.mfunction v=ksf4(x)v=1./(x.2+1);tixingjifen(ksf4,0,1,20)ans =0.7853simpson(ksf4,0,1,10)ans =0.7854 edit ksf5.mfunction v=ksf5(x)if(x=0)v=1;elsev=sin(x)./x;end(第二个函数ks
11、f5调用求积函数时,总显示有错误:“NaN”,还没调试好。见谅!)77.用 迭代方法对下面方程组求解,取初始向量 。Jacobi (0)3,21)Tx123243x解:edit Jacobi.mfunctionx iter=Jacobi(A,x0,b,tol)D=diag(diag(A);L=D-tril(A);U=D-triu(A);x=x0;iter=0;while(norm(A*x-b)/norm(b)tol)iter=iter+1;x0=x;x=D(L+U)*x0+b);end A=2 4 -4;3 3 3;4 4 2; b=2 -3 -2;x0=3 2 -1; x,iter=Jaco
12、bi(A,x0,b,1e-4)x =1-1-1iter =388.用牛顿法求解方程 在 附近的根。cos20xx解: edit Newton.mfunction x iter=Newton(f,g,x0,tol)iter=0;done=0while donex=x0-feval(f,x0)/feval(g,x0);done=norm(x-x0) edit ksf6.mfunction v=ksf6(x)v=x*cos(x)+2; edit ksg6.mfunction z=ksg(y)z=y.5+y-1; x iter=Newton(ksf6,ksg6,2,1e-4)x =2.4988iter
13、 =399.分别用改进乘幂法、反幂法计算矩阵 A 的按模最大特征值及其对应的特征向量、按模最小特征值及其对应的特征向量。 126解: edit ep.mfunction t,x=ep(A,x0,tol)tv0 ti0=max(abs(x0);lam0=x0(ti0);x0=x0./lam0;x1=A*x0;tv1 ti1=max(abs(x1);lam1=x1(ti1);x1=x1./lam1;while(abs(lam0-lam1)tol)x0=x1;lam0=lam1;x1=A*x0;tv1 ti1=max(abs(x1);lam1=x1(ti1);x1=x1./lam1;endt=lam
14、1;x=x1; edit fanep.mfunctiont,x=fanep(A,x0,tol)tv0 ti0=max(abs(x0);lam0=x0(ti0);x0=x0./lam0;x1=Ax0;tv1 ti1=max(abs(x1);lam1=x1(ti1);x1=x1./lam1;while(abs(1/lam0-1/lam1)tol)x0=x1;lam0=lam1;x1=Ax0;tv1 ti1=max(abs(x1);lam1=x1(ti1);x1=x1./lam1;endt=1/lam1;x=x1; A=12 6 -6;6 16 2;-6 2 16;x0=1 0.5 -0.5;tol
15、=1e-4;t,x=ep(A,x0,tol)t =21.5440x =1.00000.7953-0.7953 A=12 6 -6;6 16 2;-6 2 16;x0=1 0.5 -0.5;tol=1e-4;t,x=fanep(A,x0,tol)t =4.4560x =1.0000-0.62870.62871010.将积分区间 n 等分,用复合梯形求积公式计算定积分 ,比较不同 值时的误321dxn差(画出平面上 log(n)-log(Error)图) 解: edit ksf7.mfunction v=ksf7(x)v=sqrt(1+x.2); I=quad(ksf7,1,3); n=1:100
16、;fori=1:100x(i)=tixingjifen(ksf7,1,3,i);error(i)=abs(I-x(i);endplot(log10(n),log10(error(n)1111.用 迭代方法对下面方程组求解,取初始向量 。SOR(0)1,)Tx12321805x解: edit sor.mfunction x,iter=sor(A,x0,b,omega,tol)D=diag(diag(A);L=D-tril(A);U=D-triu(A);x=x0;iter=0;while(norm(A*x-b)/norm(b)tol)iter=iter+1;x0=x;x=(D-omega*L)(o
17、mega*b+(1-omega)*D*x+omega*U*x);end A=2 -1 0;-1 3 -1;0 -1 2; b=1 8 -5;x0=1 0 -1; x,iter=sor(A,x0,b,1.1,1e-4)x =1.99993.0000-1.0000iter =51212.编程实现求解满足下列条件的区间-1,2上的三次样条函数 S(x),并画出此样条函数的图形:xi -1 0 1 2f(xi) -1 0 1 0f(xi) 0 -1function splx=-1 0 1 2y=0 -1 0 1 0 -1pp=csape(x,y,complete)breaks,coefs,npolys
18、,ncoefs,dim=unmkpp(pp)xh=-1:0.1:2if -1tolx=(a+b)/2fx=feval(f,x)if sign(fx)=sign(fa)a=xfa=fxelseif sign(fx)=sign(fb)b=xfb=fxelsereturnendend14其他 A=1 2 3 4;5 6 7 8; AA =1 2 3 45 6 7 8 m n=size(A)m =2n =4 x=1 2 3 4 5;xx =1 2 3 4 5length(x)ans =5 A=2 3 4;1 1 9; 1 2 -6; AA =2 3 41 1 91 2 -6 L U=lu(A); LL
19、 =1.0000 0 00.5000 1.0000 00.5000 -1.0000 1.0000 UU =2.0000 3.0000 4.00000 -0.5000 7.00000 0 -1.000015 x=linspace(0,1,11);xans =00.10000.20000.30000.40000.50000.60000.70000.80000.90001.0000 x=1 2 3 4; y=6 11 18 27; p=polyfit(x,y,2)p =1.0000 2.0000 3.0000diag(ones(4,1),1)ans =0 1 0 0 00 0 1 0 00 0 0
20、1 00 0 0 0 10 0 0 0 01612.编程实现求解满足下列条件的区间-1,2上的三次样条函数 S(x),并画出此样条函数的图形:xi -1 0 1 2f(xi) -1 0 1 0f(xi) 0 -1function splx=-1 0 1 2y=0 -1 0 1 0 -1pp=csape(x,y,complete)breaks,coefs,npolys,ncoefs,dim=unmkpp(pp)xh=-1:0.1:2if -1tolx=(a+b)/2fx=feval(f,x)if sign(fx)=sign(fa)a=xfa=fxelseif sign(fx)=sign(fb)b=xfb=fxelsereturnendend