实验2非线性方程f(x)0的解法.doc

实验2非线性方程f(x)0的解法.doc

《实验2非线性方程f(x)0的解法.doc》由会员分享，可在线阅读，更多相关《实验2非线性方程f(x)0的解法.doc（21页珍藏版）》请在冰点文库上搜索。

百度文库-让每个人平等地提升自我

SolutionofNonlinearEquationsf（x）=0

1.实验描述

1．:

参照程序求解出单调收敛的不动点。

2．:

已知初值，时间和末值，求解汇率。

3.：

已知函数求函数在区间内的极值和根。

4.：

已知运动方程求解运动时间和距离。

2.实验内容

1．使用程序求解下面每个函数的不动点（尽肯能多）近似值，答案精确到小数点后12位。

同时，构造每个函数的图和直线y=x来显示所有不动点。

（a）

（b）

（c）

（d）

2．如果在240个月内每月付款300美元，求解满足全部年金A为500000美元的汇率I的近似值（精确到小数点后10位）。

3.一个计算机程序使用点（x0,y0），（x1,y1），（xn,yn）可画出函数y=f（x）的图形，通常还标记出图形的纵向高度，而且必须写出一个子程序来确定函数f在区间[a,b]内的最大值和最小值。

（a）构造一个算法寻找值Ymax=maxk{yk}和Ymin=mink{yk}.

（b）写一个MATLAB程序寻找函数f（x）在区间[a,b]内根的近似值位置和极值。

（c）使用（b）中的程序寻找第1题和第2题中根的位置和极值，并用真实值进行比较。

4.设投射方的运动程为

（a）求当撞击地面时经过的时间，精确到小数点后10位。

（b）求水平飞行行程，精确到小数点后10位。

3.实验结果及分析

1.:

算法：

（1）输入函数g，p0,tol,max1,令k=2。

（2）判断k>max1是否成立，如成立输出结果，如不成立，执行（3）。

（3）令p（k）=g（p（k-1））,err=|p（k）-p（k-1）|。

（4）判断err

（2）。

start

Inputg,p0,tol;k=2

k>max1

Y

N

Letp（k）=g（p（k-1））;err=|p（k）-p（k-1）|

k=k+1

err

output

Y

N

output

end

对（a）令，简单画出与的图像，可令x=0,,.

图1.与

运行程序后，当输入[k,p,err,P]=fixpt（@g,-1,1e-12,100）

输出k=3

p=2

err=0

P=-122

通过图形我们知道不动点有3个，X=2是吸引不动点，其余2个为排斥不动点。

对（b）令，简单画出与的图像，可令x=1,,.

图2.与

运行程序后，当输入[k,p,err]=fixpt（@g,1,1e-12,100）

输出：

k=37p=err=

通过图形我们知道不动点有1个，X=是吸引不动点。

因而函数只有一个吸引不动点

对（c）令，简单画出与的图像，可令x=1,,

图3.与

运行程序后，当输入[k,p,err]=fixpt（@g,0,1e-12,100）

输出k=38p=err=

运行程序后，当输入[k,p,err]=fixpt（@g,,1e-12,100）

输出k=38p=err=

通过图形我们知道不动点有2个，X=,X=是吸引不动点。

因而函数有两个吸引迭代点。

令（4）令，简单画出与的图像，可令x=1,,。

图4.与

当输入[k,p,err,P]=fixpt（@g,1,1e-12,100）

输出k=2p=1err=0

P=

1

1

通过图形我们知道不动点有2个，X=1是吸引不动点.因而函数只有一个吸引迭代点。

算法：

（1）输入f,a,b,delta,max1,令k=1。

（2）f（a）=feval（f,a）;f（b）=feval（f,b）;k=1.

（3）判断f（a）*f（b）>0是否成立，若不成立输出结果，若成立执行步骤（4）。

（4）判断k>max1是否成立,若成立输出结果，若不成立执行步骤（5）。

（5）令dx=yb*（b-a）/（yb-ya）,c=b-x,yc=feval（f,c）.

（6）判断yc=0是否成立，若成立输出结果，若不成立执行步骤（7）。

（7）判断yc*yb>0是否成立，若成立，令b=c,yb=yc.若不成立,令a=c,ya=yc.

（8）令err=|b-a|.

（9）判断err

若不成立，令k=k+1,执行（4）。

流程图：

start

Inputf,a,b,delta,max1，，

k=k+1

Y

f（a）*f（b）>0

N

output

Y

k>max1

N

dx=yb*（b-a）/（b-a）

c=b-x

yc=feval（f,c）

Y

yc=0

N

Y

N

yc*yb>0

b=c

yb=yc

a=c

ya=yc

err=|b-a|

N

err

Y

output

end

当输入[c,err,yc]=regula（@f,,,1e-10,100）

输出c=

err=

yc=

3.（a）

算法:

（1）输入f,a,b,N,令k=1.

（2）判断k>N+1是否成立，若成立，输出并执行步骤（3）；若不成立，令p（k）=feval（f,b*（k-1）/（a*N））,k=k+1,再执行步骤

（2）。

（3）令m=p

（1）,k=2.

（4）判断k>N+1是否成立，若成立，输出并执行步骤（6）；若不成立，执行步骤（5）。

（5）判断m>p（k）是否成立，若成立，m=p（k）,执行步骤（4）；若不成立，执行步骤（4）。

（6）令n=p

（1）,k=2.

（7）判断k>N+1是否成立，若成立，输出；若不成立，执行步骤（8）。

（8）判断n>p（k）是否成立，若成立，n=p（k）,输出；若不成立，执行步骤（7）。

流程图：

start

Inputf,a,b,Nk=1

k=k+1

Y

k>N+1

N

output

p（k）=feval（f,b*（k-1）/（a*N））

m=p

（1）,k=2

k=k+1

Y

k>N+1

N

N

output

m>p（k）

Y

m=p（k）

n=p

（1）,k=2

k=k+1

Y

k>N+1

N

N

n>p（k）

output

Y

end

n=p（k）

（b）

算法：

（1）输入f,X,delta.

（2）令Y=f（X）,n=length（X）,m=0,l=0,X（n+1）=X（n）,Y（n+1）=Y（n）,k=2.

（3）判断k>n+1是否成立，若成立，输出；若不成立，执行步骤（4）。

（4）判断Y（k-1）*Y（k）<=0是否成立，若成立，m=m+1，r（m）=（X（k-1）+X（k））/2;若不成立，执行步骤（5）。

（5）令s=（Y（k）-Y（k-1））*（Y（k+1）-Y（k））。

（6）判断（abs（Y（k））

（7）判断s<=1e-6是否成立，若成立，l=l+1，v（l）=Y（k）;若不成立，执行步骤（3）.

start

流程图：

Inputf,X,delta

Y=f（X）,n=length（X）,m=0,l=0,X（n+1）=X（n）,Y（n+1）=Y（n）,k=2

k=k+1

Y

k>n+1

N

N

Y（k-1）*Y（k）<=0

Y

m=m+1，r（m）=（X（k-1）+X（k））/2

s=（Y（k）-Y（k-1））*（Y（k+1）-Y（k））

N

（abs（Y（k））

Y

m=m+1，r（m）=X（k）

N

s<=1e-6

Y

l=l+1，v（l）=Y（k）

output

end

（c）用（b）中程序求第1题

输入>>X=-2:

:

2;

>>[r,v]=chao（@k,X,1e-12）

输出

r=

v=+006*

用程序求得c=

用（b）中程序求的根的位置的近似值中有3个，真实值却只在c=左右，说明程序（b）计算结果不过精确。

用（b）中程序求第2题

输入>>X=-2:

:

2;

>>[r,v]=chao（@o,X,1e-12）

输出r=

v=+004*

用程序求得c=

用（b）中程序求的根的位置的近似值中有1个，与真实值c=

接近，说明程序（b）计算结果较为精确。

算法：

（1）输入h,r,p0,p1,delta,max1.令k=1.

（2）判断k>max1是否成立，若成立，输出；若不成立，执行（3）。

（3）令p2=p1-feval（f,p1）*（p1-p0）/（feval（f,p1）-feval（f,p0））,err1=|p1-p0|,

err2=|feval（r,p1）-feval（r,p0）|,p0=p1,p1=p2

（4）判断err1

（2）。

start

流程图：

Inputf,r,p0,p1,deltamax1k=1

k>max1

P2=p1-feval（f,p1）*（p1-p0）/（feval（h,p1）-feval（h,p0））

Err1=|p1-p0|err2=|feval（r,p1）-feval（r,p0）|

C=feval（r,p1）p0=p1p1=p2

Err2

output

end

输入[p1,err1,err2,k,c]=secant（@h,@r,9,,1e-10,100）

输出p1=

err1=

err2=

k=6

c=+003

4.结论

1.采用的不动点迭代法，只能计算出吸引不动点，不能算出排斥迭代点。

2.计算得所求的利率为。

3.用（b）中程序求的根的位置的近似值中有合理的，也有不合理的。

4.计算得投射体撞击地面经过的时间p1=

飞行行程c=+003

21

附件（代码）：

1.（a）%plotthefigureofg=x.^5-3.*x.^3-2.*x.^2+2withy=x

x=-2:

:

;

y=x;

g=x.^5-3.*x.^3-2.*x.^2+2;

plot（x,g,x,y）

xlabel（'x'）;

ylabel（'g=x.^5-3.*x.^3-2.*x.^2+2y=x'）;

functionb=g（x）

b=x^5-3*x^3-2*x^2+2;

%definitefunctiong（x）

（b）%plotthefigureofcos（sin（x））withy=x

x=-2:

:

;

y=x;

g=cos（sin（x））;

plot（x,g,x,y）

xlabel（'x'）;

ylabel（'cos（sin（x））y=x'）;

%definitefunctiong（x）

functionb=g（x）

b=cos（sin（x））;

（c）%plotthefigureofg=x.^2-sin（x+withy=x

x=-2:

:

;

y=x;

g=x.^2-sin（x+;

plot（x,g,x,y）

xlabel（'x'）;

ylabel（'g=x.^2-sin（x+y=x'）;

%definitefunctiong（x）

functionb=g（x）

b=x.^2-sin（x+;

（d）plotthefigureofx.^（x-cos（x））withy=x

x=-2:

:

;

y=x;

g=x.^（x-cos（x））;

plot（x,g,x,y）

xlabel（'x'）;

ylabel（'x.^（x-cos（x））y=x'）;

%definitefunctiong（x）

functionb=g（x）

b=x.^（x-cos（x））;

function[k,p,err,P]=fixpt（g,p0,tol,max1）

%input--gisiteratinfunctioninputasstring'g'

%-p0isinitialguessforthefixedpoint

%-tolisthetolerance

%-max1isthemaximumnumberofiterations

%output-kisthenumberofiteratinsthatwerecarriedout

%-pisapproximationtothefixedpoint

%-erristheerrorintheapproximation

%-Pcontainsthesequence{pn}

P

（1）=p0;

fork=2:

max1

P（k）=feval（g,P（k-1））;

err=abs（P（k）-P（k-1））;

relerr=err/（abs（P（k））+eps）;

p=P（k）;

if（err

end

ifk==max1

disp（'maximumnumberofiterationsexceeded'）

end

P=P';

2.%definitethefunctionf（x）

functionA=f（x）

p=300;

n=240;

A=12*p*（（1+x/12）^n-1）/x-500000;

function[c,err,yc]=regula（f,a,b,delta,max1）

%input--fisthefunctioninputasastring'f'

%-aandbaretheleftandrightpoints

%-deltaisthetoleranceforthezero

%-max1isthemaximumnumberofiterations

%output-cisthezero

%-yc=f（c）

%-erristheerrorestimateforc

ya=feval（f,a）;

yb=feval（f,b）;

ifya*yb>0

disp（'ya,ybarenotsuitable'）

return

end

fork=1:

max1

dx=yb*（b-a）/（yb-ya）;

c=b-dx;

yc=feval（f,c）;

ifyc==0,break;

elseifyc*yb<0

a=c;ya=yc;

else

b=c;yb=yc;

end

err=abs（b-a）;

iferr

end

3

%definitek（x）

functiony=k（x）

y=1000000.*x.^3-111000.*x.^2+1110.*x-1;

function[r,v]=chao（f,X,delta）

%input--fistheobjectfunction

%-Xisthevectorofabscissas

%-deltaisthetolerance

%output-risthevectorofapproximateroots

%-visthevectorofdfapproximateroots

Y=f（X）;

n=length（X）;

m=0;l=0;

X（n+1）=X（n）;

Y（n+1）=Y（n）;

fork=2:

n

ifY（k-1）*Y（k）<=0

m=m+1;

r（m）=（X（k-1）+X（k））/2;

end

s=（Y（k）-Y（k-1））*（Y（k+1）-Y（k））;

if（abs（Y（k））

m=m+1;

r（m）=X（k）;

end

ifs<=1e-6

l=l+1;

v（l）=Y（k）;

end

end

%definitefunctiono（x）

functiony=o（x）

y=5.*x.^10-38.*x.^9+21.*x.^8-5.*pi.*x.^6-3.*pi.*x.^5-5.*x.^2+8.*x-3;

function[r,v]=chao（f,X,delta）

%input--fistheobjectfunction

%-Xisthevectorofabscissas

%-deltaisthetolerance

%output-risthevectorofapproximateroots

%-visthevectorofdfapproximateroots

Y=f（X）;

n=length（X）;

m=0;l=0;

X（n+1）=X（n）;

Y（n+1）=Y（n）;

fork=2:

n

ifY（k-1）*Y（k）<=0

m=m+1;

r（m）=（X（k-1）+X（k））/2;

end

s=（Y（k）-Y（k-1））*（Y（k+1）-Y（k））;

if（abs（Y（k））

m=m+1;

r（m）=X（k）;

end

ifs<=1e-12

l=l+1;

v（l）=Y（k）;

end

end

4.

%definitefunctionr（t）

functionx=r（t）

x=2400*（1-exp（-t/15））;

%definitefunctionh（t）

functiony=h（t）

y=9600*（1-exp（-t/15））-480*t;

function[p1,err1,err2,k,c]=secant（h,r,p0,p1,delta,max1）

%input--histhefunctionofaltitudeinpputasstring'h'

%--risthefunctionofpathinputasastring'r'

%--p0andp1istheinitialapproximationstoazero

%--deltaisthetoleranceforp1andy

%--max1isthemaximumnumberofiterations

%output-p1isthesecantmethodapproximationtothezero

%-err1isestimateforp1

%-err2isestimatefory

%-kisthenumberofiterations

%-yisthefunctionvaluer（p1）

fork=1:

max1

p2=p1-feval（h,p1）*（p1-p0）/（feval（h,p1）-feval（h,p0））;

err1=abs（p1-p0）;

err2=abs（feval（r,p1）-feval（r,p0））;

c=feval（r,p1）;

p0=p1;

p1=p2;

iferr1

end