multistep method 求解数值常微分方程——matlab程序

博大广阔

%the local error of the adams bashforth computer by a method is o(h^5)
%compare with fourth_order runge-kutaa method,a multistep method requires
%only one new function evalution for each step .this can lead to great
%saving in time and expense.
function multistep_method=multistep_method()
clear all
close all
%%%%%%%%%%%%%%%%%%%%%%%%%algorithm parameter and the initial value
x10=[1];                                   %state vector of initial value
                                              %should be modify for different system   
u=0;                                       %control law   
t=0; h=0.2; tmax=10;                       %simulation time and step length
%%%%%%%%%%%%%%%%%%%%%四阶龙格计算初始四个值%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  XX(:,1)=x10; time(1)=0;
for n=2:1:4   
    [u,dX]=derives(t,x10,u);  if n==2; yy(:,1)=dX;end
       m1=h*dX;              t=t+0.5*h;X=x10+0.5*m1;               %更新导数     
     [u,dX]=derives(t,X,u);
      m2=h*dX;            X=x10+0.5*m2;                             %更新导数
    [u,dX]=derives(t,X,u);
     m3=h*dX;          t=t+0.5*h; X=x10+m3;                         %更新导数
    [u,dX]=derives(t,X,u);
    m4=h*dX;   
     x10=x10+(1/6)*(m1+2*m2+2*m3+m4);
      XX(:,n)=x10;  time(n)=t;                     %存输出值
     [a,dXa]=derives(t,x10,u);     yy(:,n)=dXa;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%multistep initial value %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
YY(1:4)=XX(1,:);
Xold=XX(:,4);             %state vector
y0=yy(:,1);y1=yy(:,2);y2=yy(:,3);y3=yy(:,4);                       %存导数   

n=5;          %loop needed paramter use as save the output result
while t<tmax   
%%%%%%%%%%%%%%%%%%%%%%%%%multistep method%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    aa=Xold+h*(55*y3-59*y2+37*y1-9*y0)/24; t=t+h;
     [u,y4]=derives(t,aa,u);
    Xold=Xold+h*(9*y4+19*y3-5*y2+y1)/24;
   
     y0=y1;y1=y2;y2=y3;  y3=y4;    %rewrite the parameter ！！！！！认真点，，，，
   
    YY(n)=Xold(1)                        %output which should be modify for different system
    time(n)=t; n=n+1;         
end

  plot(time,YY,'r') ;grid on  
   
end

function [u,dX]=derives(t,X,u)
%fi=cos(3*t);u=0;
%A=[0 1;-10 -3]; B=[0;1];
dX=t+X-1;
end

fudy10

这个程序的用处是？

zhong124

这个是常微分方程里面的多重打靶法吗？