rk_wjl.asv

来自「(有源代码)数值分析作业,本文主要包括两个部分,第一部分是常微分方程(ODE)的」· ASV 代码 · 共 19 行

function [xx,yy]=rk_wjl(odefile,xi,xf,yi,h,varargin)

x = [xi:h:xf];		% Vector of x values
if x(end) ~= xf
   x(end+1) = xf;
end

d = diff(x);
% if no this sentence ,you  change d(i) to  h ,is OK
% this is a fine method 
for i = 1:length(x)-1
      k1 = d(i) * feval(odefile,x(i),y(i),varargin{:});
      k2 = d(i) * feval(odefile,x(i)+d(i)/2,y(i)+k1/2, ...
         varargin{:});
      k3 = d(i) * feval(odefile,x(i)+d(i)/2,y(i)+k2/2, ...
         varargin{:});
      k4 = d(i) * feval(odefile,x(i+1),y(i)+k3,varargin{:});
      y(i+1) = y(i) + (k1+2*k2+2*k3+k4)/6;
   end