Sign in to follow this  
orbano

ODE Solver / Adaptive stepsize ... how?

Recommended Posts

Hi there! I have finally implemented a 4th order Runge-Kutta approximation in my physics simulation. My step function looks something like this: Step(ObjectArray,ObjectCount,Stepsize) { For each object { Collect forces acting on the actual object (from global variables); Calculate displacement with RK4 (or something similar); Set new object state; } } ---------- The probles is, i really have no idea, how adaptive stepsize comes in. Now in my simulation i use no framerate limitig, physical rendering is called each time graphics rendering is called. I want my simulation to run deltaT = 1.0 in each second. If i use framerate limiting, its quite easy to achieve this, but how should i do it with adaptive stepsizing? Could anyone help? some of my questions: -should i use different stepsize for each object to reach de desired final deltaT (and set the desired deltaT according to the best approximation error)? -should i compute the worst approximation error, and adapt the whole simulation to that? please help :)

Share this post


Link to post
Share on other sites
1. If your force is fixed during a step you don't need RK4. What do you mean by collect forces? Values or functions?
2. When to do collision detection?

You need some kind of error metric which you calculate every ode-step and from which you decide to adapt your stepsize. You store the stepsize between gameloop-steps.

If you adapt your stepsize per object:
instead of Calculate displacement with RK4 (or something similar);
you do something like:
while(rk4.stepsize()<dt) {
dt-=rk4.step(); // returns actual stepsize, updates predicted maximum stepsize for next step
}
rk4.step(dt);

Share this post


Link to post
Share on other sites
You need an adaptive stepsize to get a solution with an error under some threshold while reducing calculations.

Perform this algorithm to find the best stepsize:


1)
x0 = ?
y0 = ?
H = ?

2)
calculate:
e(x0 + H, H)
e(x0 + H, H/2)
h

3)
if H/h is much greater than 2, goto 4)
else goto 5)

4)
H = 2*h
goto 2)

5)
x0 = x0 + H
y0 = e(x0 + H, H/2)
H = 2*h
finished this step


The function e(x0 + H, H) is an estimation for y, where the second parameter is the current stepsize. So e(x0 + H, H) is one big step, while e(x0 + H, H/2) are two half steps. It is used to approximate the error. You can use RK4 for those functions.

While H is the currently selected stepsize, h is used to decrease it until the error goes below some threashold:


H / h = ( 2^p/(2^p - 1) * abs( e(x0 + H, H) - e(x0 + H, H/2) ) / epsilon )^(1/(p+1))




Here epsilon should not be smaller than K*eps, where eps is the machine accuracy and K is:
K = max { abs(y(x)) | x0 <= x <= x0+h }


Good luck, its not as complicated as it looks.



Edit:
You may also want to look for Runge-Kutta-Fehlberg.

Share this post


Link to post
Share on other sites
nmi! thank you for your description, but i know, how to approximate the error, the only thing i dont know is how to deal with it :)
Im still trying to understan what Trap is saying :) My english is not so good...

Share this post


Link to post
Share on other sites
1. yes, i need RK4 even with fixed forces/step. its much more accurate. ive tested it.
2. collision detection is after the Step. it doesnt interfere with ODE solver (yet).

I can calculate the error metric you said, but dont know what to do with it :)

Bit i dont really understand this:

while(rk4.stepsize()<dt) {
dt-=rk4.step();
}
rk4.step(dt);

step(); - steps with its own estimated required stepsize to maintain accuracy?
step(dt); - what does this do? it steps back if stepsize exceeds dt?

BTW: how does collision detection come into play? now i use a very simple collision detection:

for each object(a)
{
for each object(b)
{
if a intersects b and a and b are closing, then do the collision.
}
}

//of course im not doing O(n^2) check, im using AABB trees...//

Share this post


Link to post
Share on other sites
If forces are fixed during a step you can explicitly solve the ode: x_new=x_old+v_old*dt+0.5*a*dt*dt, v_new=v_old+a*dt
Way better than rk4 :)

If your adaptive stepsize is bigger than dt, you just do a step of dt length. If it is smaller than dt you step forward with the adaptive stepsize till the remaining step is smaller than your adaptive stepsize and do a final step of dt length.

Share this post


Link to post
Share on other sites
oh yes got it :)) just didnt look at the loop...

i think your equation is worse (its something like RK2). but tonight i will sit down and calculate the equations on paper... but the tests showed, that RK4 need ten times larger stepsize than Midpoint method, to achieve same accuracy... we only have to sit down and calculate the errors...

Share this post


Link to post
Share on other sites

Create an account or sign in to comment

You need to be a member in order to leave a comment

Create an account

Sign up for a new account in our community. It's easy!

Register a new account

Sign in

Already have an account? Sign in here.

Sign In Now

Sign in to follow this