matlab数值分析程序--高等数学,数值代数的matlab实现-文字版, matlab电子书, 和matlab 有关的电子书:

10.4.4 改进的欧拉公式的MATLAB程序

10.4.4  改进的欧拉公式的MATLAB程序

用自适应改进的欧拉公式求解常微分方程初值问题(10.5)的数值解的MATLAB主程序

function [H,X,Y,k,h,P]=gaiEuler(funfcn,x0,b,y0,tol)

%初始化.

pow=1/3;

if nargin < 5 | isempty(tol),

tol = 1.e-6;

end;

if nargin < 6 | isempty(trace),

 trace = 0;

end;

x=x0; h=0.0078125*(b-x);

y=y0(:);p=128; n=fix((b-x0)/h);

H=zeros(p,1); X=zeros(p,1);

Y=zeros(p,length(y)); k=1;

X(k)=x; Y(k,:)=y';

% 绘图.

clc,x,h,y

end

% 主循环.

while (x<b)&(x+h>x)

  if x+h>b

h=b-x;

  end

   %计算斜率.

   fxy=feval(funfcn,x,y); fxy=fxy(:);

%计算误差,设定可接受误差.

delta=norm(h*fxy,'inf'); wucha=tol*max(norm(y,'inf'),1.0);

% 当误差可接受时重写解.

if delta<=wucha

x=x+h; y1=y+h*fxy; fxy1=feval(funfcn,x,y1);

fxy=fxy(:);y2=(fxy+fxy1)/2; y=y+h*y2; k=k+1;

if k>length(X)

X=[X;zeros(p,1)]; Y=[Y;zeros(p,length(y))];

H=[H;zeros(p,1)];

        end

H(k)=h;X(k)=x;Y(k,:)=y'; plot(X,Y,'mh'), grid

xlabel('自变量 X'), ylabel('因变量 Y')

title('用改进的欧拉公式计算dy/dx=f(x,y),y(x0)=y0在[x0,b]上的数值解')

end

%更新步长.

if delta~=0.0

h=min(h*8,0.9*h*(wucha/delta)^pow);

end

end

if (x<b)

disp('Singularity likely.'), x

  end

H=H(1:k);  X=X(1:k); Y=Y(1:k,:); n=1:k; P=[n',H,X,Y]

 

例10.4.5  分别用梯形公式和改进的欧拉公式求解区间上的初值问题,取精度为,与精确解作比较,在同一个坐标系中作出图形.

解  输入程序

>> x0=0;y0=1; b=2;tol=1.e-1;

[Ht,X,Yt,k,h,Pt]=odtixing2(@funfcn,x0,b,y0,tol)

hold on

S1 = 8/3*X-29/9+38/9*exp(-3*X),

plot(X,S1,'b-')

hold off

hold on

[H,X,Y,k,h,P]=gaiEuler(@funfcn,x0,b,y0,tol)

hold off

legend('用梯形公式计算dy/dx=8x-3y-7,y(0)=1在[0,2]上的数值解','dy/dx= 8x-3y-7,y(0)=1在[0,2]上的精确解', '用改进的欧拉公式计算dy/dx=8x-3y-7,y(0)=1在[0,2]上的数值解')

,juwY=S1-Y; xiwY=juwY./Y;

L=[P,S1,juwY,xiwY]

运行后屏幕显示取精度为,分别用梯形公式和改进的欧拉公式求此解初值问题在区间上的自变量X处数值解Yt, Y和精确解S1及其图形,步长H, Y的相对误差xiwY和绝对误差juwY(略).

 

欢迎转载,转载请注明来自一手册:http://yishouce.com/book/3/3100404.html
友情链接It题库(ittiku.com)| 版权归yishouce.com所有| 友链等可联系 admin#yishouce.com|粤ICP备16001685号-1