王能超 计算方法——算法设计及MATLAB实现课后代码 联系客服

发布时间 : 星期一 文章王能超 计算方法——算法设计及MATLAB实现课后代码更新完毕开始阅读2b5c7cbbad02de80d5d84004

为衡量数值解的精度,我们求出该方程的解析解为y?e?x?x.在此也以文件的形式表示如下:

function y=solvef3(x) y=exp(-x)+x;

令f=@f3; a=0; b=1; N=10; ya=1; 运行:A2=Adams2PC(f,a,b,N,ya); A4=Adams4PC(f,a,b,N,ya); CA4=CAdams4PC(f,a,b,N,ya); y=solvef3(a:(b-a)/N:b); m=[A2,A4(:,2),CA4(:.2),y’]

其中m共有5列数据:左→右依次为离散节点值、二阶Adams预报校正系统所求解、四阶Adams预报校正系统所求解、改进四阶Adams预报校正系统所求解和准确解,精度依次递增.

第四章 方程求根 4.1 二分法

用二分法求解非线性方程f(x)?0在区间[a,b]内的根. MATLAB文件:(文件名:demimethod.m) function [x,k]=demimethod(a,b,f,emg) %a,b表示求解区间[a,b]的端点 %f表示所求解方程德函数名 %emg是精度指标 %x表示所求近似解 %k表示循环次数 fa=feval(f,a);

fab=feval(f,(a+b)/2); k=0;

while abs(b-a)>emg if fab==0 x=(a+b)/2; return;

elseif fa*fab<0 b=(a+b)/2; else

a=(a+b)/2; end

fa=feval(f,a);

fab=feval(f,(a+b)/2); k=k+1; end

x=(a+b)/2;

13

算例9:求解方程f(x)?解:

function f=func2(x)

f=sqrt(x^2+1)-tan(x);

x2?1?tanx在区间[0,?/2]内的实根,使精度达到10?5.

输入:f=@@func2; a=0; b=pi/2; emg=10^-5; [x0,k]=demimethod(a,b,f,emg)

4.2 开方法

求实数的开方运算. MATLAB文件:(文件名:Kaifang.m) function y=Kaifang(a,eps,x0) %a是被开方数 %eps是精度指标 %x0表示初值 %y是a的开方 x(1)=x0;

x(2)=(x(1)+a/x(1))/2; k=2;

while abs(x(k)-x(k-1))>eps x(k+1)=(x(k)+a/x(k))/2; k=k+1; end y=x’;

算例10:用开方法求2,设取x0=1. 解:令a=2; eps=10^-6; x0=1; 运行 y=Kaifang(a,eps,x0)

4.3 Newton下山法

用牛顿下山法求解非线性方程f(x)?0的根.

MATLAB文件:(文件名:Mendnewton.m) function [x,k]=Mendnewton(f,x0,emg) %用牛顿下山法求解非线性方程 %f表示非线性方程

%x0是迭代初值,此方法是局部收敛,初值选择要恰当 %emg是精度指标

%k,u分别表示迭代次数和下山因子

[f1,d1]=feval(f,x0); ?表示非线性方程f在x0处的导数值,以下类同 k=1; x(1)=x0;

14

x(2)=x(1)-f1/d1;

while abs(f1)>emg %控制精度 u=1; k=k+1;

[f1,d1]=feval(f,x(k));

x(k+1)=x(k)-u*f1/d1; %牛顿下山迭代 [f2,d2]=feval(f,x(k+1));

while abs(f2)>abs(f1) %保证迭代后的函数值比迭代前的函数值小 u=u/2;

x(k+1)=x(k)-u*f1/d1; %牛顿下山迭代 [f2,d2]=feval(f,x(k+1)); end end

算例11:用牛顿下山法求方程f(x)?(1)x0=-1.2; (2)x0=2.0;

解:编写func3.m. function [f,d]=func3(x) f=sqrt(x^2+1)-tan(x); d1=’sqrt(x^2+1)-tan(x)’;

d=subs(diff(d1)); %对函数f求一次导数 在命令窗口中输入:

f=@func3; [x,k]=Mendnewton(f,x0,10^-6) 其中,x0需要自己输入.

4.4 快速弦截法

用快速弦截法求解非线性方程f(x)?0的根. MATLAB文件:(文件名:Fast_chord.m) function [x,k]=Fast_chord(f,x1,x2,emg)

%用快速弦截法求解非线性方程 %f表示非线性方程函数 %x1,x2是迭代初值 %emg是精度指标 %k表示循环次数 k=1;

y1=feval(f,x1); y2=feval(f,x2);

x(k)=x2-(x2-x1)*y2/(y2-y1); %用快速弦截法进行迭代求解 y(k)=feval(f,x(k)); k=k+1;

x(k)=x(k-1)-(x(k-1)-x2)*y(k-1)/(y(k-1)-y2);

15

x2?1?tanx的根,使精度达到10?6.初值分别为:

while abs(x(k)-x(k-1))>emg %控制精度 y(k)=feval(f,x(k));

x(k+1)=x(k)-(x(k)-x(k-1))*y(k)/(y(k)-y(k-1)); k=k+1; end

算例12:用快速弦截法求解非线性方程4cosx?e的根,要求精度为??10,初值为:

x?6x1??/4,x2??/2.

解:编写函数文件func4.m function f=func4(x)

f=exp(x)-4*cos(x)

在命令窗口输入:

f=@func4; [x,k]=Fast_chord(f,pi/4,pi/2,10^-6)

第五章 线性方程组的迭代法 5.1 Jacobi迭代

用Jacobi迭代法求解线性方程组. MATLAB文件:(文件名:Jacobimethod.m) function [x,k]=Jacobimethod(A,b,x0,N,emg)

%A是线性方程组的左端矩阵 %b是右端向量 %x0是迭代初始值

%N表示迭代次数上限,若迭代次数大于N,则迭代失败 %emg表示控制精度

%用Jacobi迭代法求线性方程组A*x=b的解 %k表示迭代次数

%x表示用迭代法求得的线性方程组的近似解 n=length(A);

x1=zeros(n,1); x2=zeros(n,1); x1=x0;

r=max(abs(b-A*x1)); k=0;

while r>emg for i=1:n sum=0; for j=1:n

if i~=j

sum=sum+A(i,j)*x1(j); end end

x2(i)=(b(i)-sum)/A(i,i);

16