How to write general function of 4th Order Runge-Kutta Method?
182 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
I am trying to develop a Matlab function for the 4th Order Runge-Kutta Method. It needs to be able to work with any function for given initial conditions, step size, etc. and then plot the results afterwards. Here is the code that I have so far. There are no error in the function itself. However, when I try to use it, I get a couple of error messages that I have also shown below the mfile. I am not sure what the error messages mean or how I would go about correcting them. Any help would be much appreciated!
FUNCTION:
function [tp,yp] = rk4(dydt,tspan,y0,h)
if nargin < 4, error('at least 4 input arguments required'), end
if any(diff(tspan)<=0), error('tspan not ascending order'), end
n = length(tspan);
ti = tspan(1); tf = tspan(n);
tp = ti:h:tf;
yp = zeros(ti,length(tspan));
for n = ti:tf
tp = ti;
k1 = dydt(tp(n),y0(n));
k2 = dydt(tp(n) + (1/2)*h, y0(n) + (1/2)*k1*h);
k3 = dydt(tp(n) + (1/2)*h, y0(n) + (1/2)*k2*h);
k4 = dydt(tp(n) + h, y0(n) + k3*h);
yp = y0(n) + ((1/6)*(k1 + (2*k2) + (2*k3) + k4)*h);
ti = ti + h;
y0 = yp;
end
plot(tp,yp)
COMMAND WINDOW:
>> [tp,yp] = rk4(@(CA) ((1/5)*(20-CA))-(0.12*CA),15,0,1)
Attempted to access tp(15); index out of bounds because
numel(tp)=1.
Error in rk4 (line 13)
k1 = dydt(tp(n),y0(n));
1 commentaire
Réponse acceptée
Star Strider
le 6 Déc 2014
I didn’t run your code, but see if changing the initial ‘tp’ assignment to:
tv = ti:h:tf;
and then in your loop, changing:
tp = ti;
to:
tp = tv(n);
and changing all the subsequent references to ‘tp(n)’ to simply ‘tp’.
Also, you will need change this line (adding a subscript to ‘yp’) to save ‘yp’ at each iteration:
yp(n) = y0(n) + ((1/6)*(k1 + (2*k2) + (2*k3) + k4)*h);
See if that improves things. You may have to make some other changes as well to accommodate these, and there may be other problems in your code, but this should keep you progressing towards a successful conclusion.
7 commentaires
Steve Chuks
le 6 Mar 2015
Hi all... @StarStrider & Alicia.
I have a similar work as to the Runge-Kutta method to solve for ODE.
Can you give me a hint on how the plot will look like and the error estimation of that? Just thinking of how to implement that on matlab.
I appreciate your assistance. Thanks.
Plus de réponses (1)
SHIVANI TIWARI
le 26 Avr 2019
hi everyone.. can u guys help me to generate the code for solving 1st order differential equations using improved runge kutta method.
2 commentaires
ahti Maître
le 6 Sep 2020
Modifié(e) : ahti Maître
le 6 Sep 2020
this worked for me at least, with y1 the derivative of your function,(yo,ti) the starting point, tf the final abscissa, and h the step size:
function [tp,yp] = rk4_1(y1,ti,tf,h,y0)
n=fix((tf-ti)/h)+1;
tv = ti:h:tf;
yp = zeros(1,n);
for j = 1:n
tp = tv(j);
k1 = y1(tp,y0);
k2 = y1(tp + (1/2)*h, y0 + (1/2)*k1*h);
k3 = y1(tp + (1/2)*h, y0 + (1/2)*k2*h);
k4 = y1(tp + h, y0 + k3*h);
yp(j) = y0+((1/6)*(k1 + (2*k2) + (2*k3) + k4)*h);
ti = ti + h;
y0 = yp(j);
end
plot(tv,yp)
Voir également
Catégories
En savoir plus sur Ordinary Differential Equations dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!