同济大学数值分析报告matlab编程题总汇编

MATLAB 编程题库 1.下面的数据表近似地满足函数2

1cx

b

ax y ++=，请适当变换成为线性最小二乘问题，编程求最好的系数c b a ,,，并在同一个图上画出所有数据和函数图像.

625

.0718.0801.0823.0802

.0687

.0606

.0356

.0995.0628.0544.0008.0213.0362.0586.0931.0i

i y x ----

解：

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=A\y;

a=z(1); b=z(2); c=z(3); xh=-1:0.1:1;

yh=(a.*xh+b)./(1+c.*xh.^2); plot(x,y,'r+',xh,yh,'b*')

2.若在Matlab工作目录下已经有如下两个函数文件，写一个割线法程序，求出这两个函数

10 的近似根，并写出调用方式：

精度为10

解：

>> edit gexianfa.m

function [x iter]=gexianfa(f,x0,x1,tol)

iter=0;

while(norm(x1-x0)>tol)

iter=iter+1;

x=x1-feval(f,x1).*(x1-x0)./(feval(f,x1)-feval(f,x0));

x0=x1;x1=x;

end

>> edit f.m

function v=f(x)

v=x.*log(x)-1;

>> edit g.m

function z=g(y)

z=y.^5+y-1;

>> [x1 iter1]=gexianfa('f',1,3,1e-10)

x1 =

1.7632

iter1 =

6

>> [x2 iter2]=gexianfa('g',0,1,1e-10)

x2 =

0.7549

iter2 =

8

3.使用GS 迭代求解下述线性代数方程组：

1231231

2

3

521242103103

x x x x x x x x x

解：

>> edit gsdiedai.m

function [x iter]=gsdiedai(A,x0,b,tol) 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); >>x x =

-3.0910 1.2372 0.9802 >>iter iter = 6

4.用四阶Range-kutta 方法求解下述常微分方程初值问题（取步长h=0.01）

,

(1)2

x

dy

y e xy dx

y

解：

>> edit ksf2.m

function 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; end >>y

调用函数方法

>> edit Rangekutta.m

function [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); >>y

5.取0.2h =，请编写Matlab 程序，分别用欧拉方法、改进欧拉方法在12x ≤≤上求解初值问题。

3,(1)0.4

dy y

x dx x y

解：

>> edit Euler.m

function [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.m

function[x 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=y(i-1)+h*feval(f,x(i),y1); y(i)=(y1+y2)./2; end

>> edit ksf3.m

function 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.0000 y =

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.0000 y =

0.4000 0.5895 0.9278 1.4615 2.2464 3.3466

6.请编写复合梯形积分公式的Matlab程序，计算下面积分的近似值，区间等分20

n=。编写辛普森积分公式的Matlab程序，计算下面积分的近似值，区间等分10

n=。

1

2 01

1

dx

x 、

1

1

sin x

dx

x

解：

>> edit tixingjifen.m

function 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.m

function 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.m

function v=ksf4(x)

v=1./(x.^2+1);

>>tixingjifen('ksf4',0,1,20)

ans =

0.7853

>>simpson('ksf4',0,1,10)

ans =

0.7854

>> edit ksf5.m

function v=ksf5(x)

if(x==0)

v=1;

else

v=sin(x)./x;

end

(第二个函数‘ksf5’调用求积函数时，总显示有错误：“NaN”，还没调试好。见谅！)