Initial Value Problems - Numerical Methods in Engineering with MATLAB - page 280

the errorbehaves erraticallyat the pointofdiscontinuity, the program can get stuck
in an infinite loop trying to find the appropriate valueof h .We wouldalso use a
nonadaptive methodif the output istohaveevenly spacedvalues of x .
runKut5
The adaptive Runge-Kuttamethodis implementedinthe function runKut5 listed
below. The input argument h is the trial value of the incrementfor the first integration
step.
function [xSol,ySol] = runKut5(dEqs,x,y,xStop,h,eTol)
% 5th-order Runge-Kutta integration.
% USAGE: [xSol,ySol] = runKut5(dEqs,x,y,xStop,h,eTol)
% INPUT:
% dEqs = handle of function that specifyies the
% 1st-order differential equations
% F(x,y) = [dy1/dx dy2/dx dy3/dx ...].
% x,y = initial values; y must be row vector.
% xStop = terminal value of x.
% h = trial value of increment of x.
% eTol = per-step error tolerance (default = 1.0e-6).
% OUTPUT:
% xSol = x-values at which solution is computed.
% ySol = values of y corresponding to the x-values.
ifsize(y,1)>1;y=y';end %ymustberowvector
if nargin < 6; eTol = 1.0e-6; end
n = length(y);
A=[01/53/103/517/8];
B=[ 0 0 0 0 0
1/5 0 0 0 0
3/40 9/40 0 0 0
3/10 -9/10 6/5 0 0
-11/54 5/2 -70/27 35/27 0
1631/55296 175/512 575/13824 44275/110592 253/4096];
C = [37/378 0 250/621 125/594 0 512/1771];
D = [2825/27648 0 18575/48384 13525/55296 277/14336 1/4];
% Initialize solution
xSol = zeros(2,1); ySol = zeros(2,n);
xSol(1) = x; ySol(1,:) = y;
stopper=0;k=1;