- #1
iwonde
- 31
- 0
Homework Statement
I need to write a Matlab program that solves the following:
Use the Runge-Kutta Fehlberg Algorithm with tolerance TOL= 10^(-4) to approximate the solution to the initial-value problem:
y'=(y/t)^2 + y/t, 1<=t<=1.2, y(1)=1, with hmax=0.05 and hmin=0.02.
Homework Equations
The Attempt at a Solution
%This is the source code for the function rkf.
% w denotes the approximation to y
% h=step size
function [w,t,h]=rkf(a,b,w0,TOL,hmax,hmin)
%Step 1
h=hmax; t=a; w=w0; FLAG=1;
%Step 2
while(FLAG==1)
%Step 3
F1 = h*f(t,w);
F2 = h*f(t+h/4 ,w+F1/4);
F3 = h*f(t+3*h/8,w+ 3*F1/32+9*F2/32);
F4 = h*f(t+12*h/13,w+932*F1/2197-7200*F2/2197+7296*F3/2197);
F5 = h*f(t+h,w+(439*F1)/216-(8*F2)+3680*F3/513-(845*F4)/4104);
F6 = h*f(t+h/2,w-8*F1/27+2*F2-(3544*F3)/2565+1859*F4/4104-11*F5/40);
%Step 4
R = abs(F1/360-128*F3/4275-2197*F4/75240 + F5/50 + 2*F6/55)/h;
%Step 5
if (R<=TOL)
%Step 6
t= t+h;
w= w + 25*F1/216 + 1408*F3/2565 + 2197*F4/4104 - F5/5;
end
%Step 8
delta=0.84*(TOL/R)^(1/4);
%Step 9
if delta <= 0.1
h=0.1*h;
elseif delta>=4
h=4*h;
else
h=delta*h;
end
%Step 10
if h>hmax
h=hmax;
end
%Step 11
if t>=b
FLAG=0;
elseif (t+h)>b
h=b-t;
elseif h<hmin
FLAG=0;
fprintf('minimum h exceeded\n');
return;
end
end
end
=====================================
%I've tried compiling the code with:
a=1; b=1.2; TOL=10^(-4); w0=1; hmax=0.05; hmin=0.02;
f= @(t,w) (w/t)^2 + w/t;
[w,t,h]=rkf(a,b,w0,f,TOL,hmax,hmin);
w
t
h
======================================
The above code didn't work. I'm attempting at a code that outputs something like:
i t_i w_i h_i
---------------------------------
1 1.05
2 1.10
3 1.15
4 1.20
I'm not really sure of how to add in the index for t,w, h to obtain a value at every step. Any help would be greatly appreciated. Thanks!