数值计算方法matlab程序

  • 格式:doc
  • 大小:80.50 KB
  • 文档页数:32

下载文档原格式

  / 32
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

function [x0,k]=bisect1(fun1,a,b,ep)

if nargin<4

ep=1e-5;

end

fa=feval(fun1,a);

fb=feval(fun1,b);

if fa*fb>0

x0=[fa,fb];

k=0;

return;

end

k=1;

while abs(b-a)/2>ep

x=(a+b)/2;

fx=feval(fun1,x);

if fx*fa<0

b=x;

fb=fx;

else

a=x;

fa=fx;

end

end

x0=(a+b)/2;

>> fun1=inline('x^3-x-1');

>> [x0,k]=bisect1(fun1,1.3,1.4,1e-4)

x0 =

1.3247

k =

7

>>

简单迭代法

function [x0,k]=iterate1(fun1,x0,ep,N)

if nargin<4

N=500;

end

if nargin<3

ep=1e-5;

end

x=x0;

x0=x+2*ep;

while abs(x-x0)>ep & k

x0=x;

x=feval(fun1,x0);

k=k+1;

end

x0=x;

if k==N

warning('已达最大迭代次数')

end

>> fun1=inline('(x+1)^(1/3)'); >> [x0,k]=iterate1(fun1,1.5)

x0 =

1.3247

k =

7

>> fun1=inline('x^3-1'); >> [x0,k]=iterate1(fun1,1.5)

x0 =

Inf

k =

>>

Steffesen加速迭代(简单迭代法的加速)function [x0,k]=steffesen1(fun1,x0,ep,N)

if nargin<4

N=500;

end

if nargin<3

ep=1e-5;

end

x=x0;

x0=x+2*ep;

k=0;

while abs(x-x0)>ep & k

x0=x;

y=feval(fun1,x0);

z=feval(fun1,y);

x=x0-(y-x0)^2/(z-2*y+x0);

k=k+1;

end

x0=x;

if k==N

warning('已达最大迭代次数')

end

>> fun1=inline('(x+1)^(1/3)');

>> [x0,k]=steffesen1(fun1,1.5)

x0 =

1.3247

k =

3

>> fun1=inline('x^3-1');

>> [x0,k]=steffesen1(fun1,1.5)

x0 =

1.3247

k =

6

Newton迭代

function [x0,k]=Newton7(fname,dfname,x0,ep,N)

if nargin<5

N=500;

end

if nargin<4

ep=1e-5;

end

x=x0;

x0=x+2*ep;

k=0;

while abs(x-x0)>ep & k

x0=x;

x=x0-feval(fname,x0)/feval(dfname,x0);

k=k+1;

end

x0=x;

if k==N

warning('已达最大迭代次数')

end

>> fname=inline('x-cos(x)');

>> dfname=inline('1+sin(x)');

>> [x0,k]=Newton7(fname,dfname,pi/4,1e-8)

x0 =

0.7391

k =4

非线性方程求根的Matlab函数调用举例:

1.求多项式的根:求f(x)=x^3-x-1=0的根:

>> roots([1 0 -1 -1])

ans =

1.3247

-0.6624 + 0.5623i

-0.6624 - 0.5623i

2.求一般函数的根

>> fun=inline('x*sin(x^2-x-1)','x')

fun =

Inline function:

fun(x) = x*sin(x^2-x-1)

>> fplot(fun,[-2 0.1]);grid on

>> x=fzero(fun,[-2,-1]) x =

-1.5956

>> x=fzero(fun,[-1 -0.1]) x =

-0.6180

[x,f,h]=fsolve(fun,-1.6) x =

-1.5956

f =

1.4909e-009

h =