# 4th Order Runge Kutta

This topic is 4740 days old which is more than the 365 day threshold we allow for new replies. Please post a new topic.

## Recommended Posts

I found some simple C code for Runge Kutta from here: http://www.physics.orst.edu/~rubin/nacphy/ComPhys/DIFFEQ/PRO/srk_c/srk_c.c But I feel it's wrong. Mainly cause f(x, y) isn't doing anything with the time. All it returns is a negative y. I converted most of it by change x into t for time, and got rid of saving it to file, but need some proper code for the f(t,y), or any other corrections that may be needed:
/* Runge Kutta algorithm for first-order differential equations*/

#include <stdio.h>

#define dist 0.5                        /* stepsize in t */
#define MAX 4                           /* max for t */

double runge4(double t, double y, double step); /* Runge-Kutta function */

double f(double t, double y);           /* function for derivative */

double  f(double t, double y)
{
return(-y);
}

double runge4(double t, double y, double step)
{
double h=step/2.0, k1, k2, k3, k4;

k1=step*f(t, y);
k2=step*f(t+h, y+k1/2.0);
k3=step*f(t+h, y+k2/2.0);
k4=step*f(t+step, y+k3);

printf("%f\n",k1);
printf("%f\n",k2);
printf("%f\n",k3);
printf("%f\n",k4);

return(y+(k1+2*k2+2*k3+k4)/6.0);
}

main()
{
double t, y;
int j;

y=1;                                    /* initial position */

//printf("0\t%f\n", y);

y = runge4(1, y, 10);

/*
for (j=1; j*dist<=MAX; j++)
{
t=j*dist;
y=runge4(t, y, dist);

printf("%f\t%f\n", t, y);
}
*/
}



##### Share on other sites
It isn't wrong. It's using the runge-kutta formula to integrate dy/dx = -y.

The output should be an exponential curve (it should exponentially decrease towards 0). If it isn't, then something is wrong.

You need to replace f with whatever you're integrating.

##### Share on other sites
y' = f(t, y) is the general form for an ODE, but many common applications--especially dynamics--won't use t at all. Consider, for instance, an equation describing the acceleration of a pendulum. The acceleration will vary given the bob's position and velocity, but not based on whether it's 3:00 or 5:00.

##### Share on other sites
That's just it though. I want to use t, and it's not using it. What can I do with t? Should I multiply it with y?

##### Share on other sites
Well, I don't know. Should your f(t,y) be t+y? Should it be 4t^2+7y? This has nothing to do with RK4 and everything to do with what it is you are actually simulating. f(t,y) is a function YOU will write, which represents your system. f(t,y) should be the derivative dy/dt at the point (t,y).

##### Share on other sites
Ok, well let's say I wanted y = y0 + a * t, which is a basic physics equation. Would it be ok for the function f(t,y) to return y0 + a * t, assuming I had y0 as an arguement?

##### Share on other sites
Quote:
 Original post by Jacob RomanOk, well let's say I wanted y = y0 + a * t, which is a basic physics equation. Would it be ok for the function f(t,y) to return y0 + a * t, assuming I had y0 as an arguement?

Quote:
 Original post by Sneftelf(t,y) should be the derivative dy/dt at the point (t,y).

ergo, your function f should simply return a.

and yes, since 'a' is a constant, 4-th order runge kutta is overkill for this simple problem.

##### Share on other sites
Actually I planned on using this for Rigid Body physics, so it shouldn't be overkill. Never worked with rigid bodies before.

##### Share on other sites
Just to make sure that the issue is clear...

A set of first-order ordinary differential equations (ODE) can be written in explicit form as such:

(1) dy/dt = f(t, y),
where y is a vector and f is a vector function.

This means that dy/dt (or y') can be given as a function of t and y.

For example, let's look at Newton's equation of motion.
If the force acting on a point mass can be expressed as a function of the points position, we can write the equation of motion like this:

(2) dx/dt = v
(3) dv/dt = m-1F(x).

Now, if we define:
y = (x, v) and f(t, y) = (v, m-1F(x)),

then (1) is equivalent to (2) and (3).

For instance, if we have a satellite orbiting a planet which is located at the origin, we can write F(x) like this:
F(x) = -GMmx / |x|3

Here's a more empirical example:
(4) y(t) = t*e2t

The derivative of y w.r.t t is:
(5) dy/dt = e2t + 2*t*e2t = 1/t*y(t) + 2*y(t)

So if we define:
f(t,y) = 1/t*y(t) + 2*y(t),

we can write (5) like this:
(6) dy/dt = f(t,y).

One of the solutions to the ODE (6) is the equation (4).

PS. Most commerical engines use semi-implicit integration for rigid body simulation.

##### Share on other sites
Quote:
 Original post by Jacob RomanActually I planned on using this for Rigid Body physics, so it shouldn't be overkill. Never worked with rigid bodies before.

rigid body physics can mean a lot of things, but regardless, i doubt runge kutta is the best solution for any rigid body problem i can imagine.

• ### What is your GameDev Story?

In 2019 we are celebrating 20 years of GameDev.net! Share your GameDev Story with us.

(You must login to your GameDev.net account.)

• 14
• 14
• 45
• 22
• 27
• ### Forum Statistics

• Total Topics
634044
• Total Posts
3015211
×