Integration Basics[Integration#Physics#game#flash#EN]
Integration Basics
Introduction
Hi, I’m Glenn Fiedler and welcome to the first article in my series on Game Physics.
If you have ever wondered how the physics simulation in a computer game works then this series of articles will explain it for you. I assume you are proficient with C++ and have a basic grasp of physics and mathematics. Nothing else will be required if you pay attention and study the example source code.
A game physics simulation works by making many small predictions based on the laws of physics. These predictions are actually quite simple, and basically boil down to something like “given that the object is here, and is traveling this fast, in a short amount of time from now should be over there”. We perform these predictions using a mathematical technique called integration.
Exactly how to implement this integration is the subject of this article.
Integrating the Equations of Motion
You should remember from high school or university physics that force equals mass times acceleration. We can switch this equation around to see that acceleration equals force divided by mass. This makes sense because the more an object weighs, the less acceleration it will receive from the same amount of force.
f = ma
a = f/m
Since acceleration is the rate of change in velocity over time, and acceleration is force divided by mass, we can say that force divided by mass is the rate of change in velocity:
dv/dt = a = F/m
Similarly, velocity is the rate of change in position over time:
dp/dt = v
This means that if we know the current position and velocity of an object, and the forces that will be applied to it, we can integrate to find its position and velocity at some point in the future.
Numerical Integration
For those who have not formally studied differential equations at university, take heart for you are in just as good a position as those who have! This is because we are not going to solve the differential equations as you would normally do in first year mathematics, instead we are going to numerically integrate to find the solution.
Here is how numerical integration works. First, start at a certain initial position and velocity, then take a small step forward in time to find the velocity and position at the next time value. Do this repeatedly to go forward in time in small increments, each time taking the results of the previous integration as the starting point for the next.
For simplicity, lets assume we are stepping forward one second at a time. First we start at our initial velocity and position, then work out how much the velocity and position change over the course of one second. Then, we add this change to the velocity and position and repeat the process. Starting with the velocity values we now have at 1 second, we work out how much the values will now change over the next second to find the values at 2 seconds and so on for 3, 4, 5, 6… etc.
But how do we find the amount of change in velocity and position each step? The answer lies in the equations of motion presented in the paragraph above. We work out how much acceleration is being applied by dividing force by mass. Since acceleration is literally change in velocity per second, we realize that the change in velocity is acceleration * (time in seconds). Similarly, the change in position is velocity * (time in seconds).
We will call our current time ‘t’ and the change in time between each step ‘dt’ (delta time). Now we can present the equations of motion in a form that should make sense to anyone:
- acceleration = force divided by mass
- change in velocity = acceleration times delta time
- change in position = velocity times delta time
It makes sense because we know if you are traveling 60 kilometers per hour in a car, then in one hour we’ll be exactly 60 kilometers further down the road. Similarly, if a car is accelerating at 10 kilometers per hour per second from a standing start, then in 10 seconds it must be traveling at exactly 100 kilometers per hour.
Lets put this into code. Assuming that our initial position and velocity are zero, we can integrate the equations of motion over time simply by adding the change in position and velocity to the previous values each timestep:
float t = 0;
float dt = 1;
float velocity = 0;
float position = 0;
float force = 10;
float mass = 1;
while ( t <= 10 )
{
position = position + velocity * dt;
velocity = velocity + ( force / mass ) * dt;
t = t + dt;
}
Notice how the integrated velocity feeds into the integration of the position next frame. This uses the most recently integrated velocity value to update the position of the car over time so that at each second we know how far the car has travelled down the road and how fast it is traveling.
Why Euler is never good enough
What we just did above is called Euler Integration. To save you any future embarrassment, I must point out now that Euler is pronounced “Oiler” not “yew-ler” as it is the last name of the Swiss mathematician Leonhard Euler who first discovered the technique.
Euler integration is the most basic of numerical integration techniques. It is only 100% accurate if the rate of change is constant over the timestep. For example it is accurate in integrating velocity over time in the example above because acceleration is a constant 10 meters per second per second. But when we integrate velocity to get the position in the example, it is not completely accurate. This is because velocity is not constant, it is changing over time due to the acceleration.
This means that the values for x(t) and v(t) that the integrator produces over time are different from the true values for the functions of position and velocity over time. In the example above the distance travelled down the road is somewhat inaccurate because of the integrator.
Just how inaccurate is the Euler integrator? From basic physics we know that:
s = ut + 0.5at^2
That is, given a constant acceleration a, and an initial velocity u, we can work out how much distance (s) is travelled over the time t.
In our case, u=0 because the car is at a standing start, so:
s = 0.5at^2
s = 0.5(10)(10)^2
s = 0.5(10)(100)
s = 500 meters
If we use the Euler integrator above instead of solving exactly we get:
t=0: position = 0 velocity = 0 t=1: position = 0 velocity = 10 t=2: position = 10 velocity = 20 t=3: position = 30 velocity = 30 t=4: position = 60 velocity = 40 t=5: position = 100 velocity = 50 t=6: position = 150 velocity = 60 t=7: position = 210 velocity = 70 t=8: position = 280 velocity = 80 t=9: position = 360 velocity = 90 t=10: position = 450 velocity = 100
So when using the Euler integrator we have an inaccuracy of 50 meters for this simple calculation after only 10 seconds! The bad news is that this error keeps getting larger as time goes on.
One solution to this problem is to decrease the time step. For example we can reduce the time step to dt = 1/100th of second each Euler integration instead of one second. Doing this will get us a result much closer to the exact value.
However, no matter how much you reduce your timestep, the simple truth is that the error will always be there and that it will keep increasing over time. Given that this is an extremely simple simulation, imagine something as simple as adding a torque curve to the car, or adding gears. Suddenly the car’s acceleration is not a constant, it changes over time. Now there is error when integrating the velocity, and these errors are magnified even further when integrating the position.
If you use Euler then you are a bloody idiot
The bottom line is that the Euler Integrator is basically crap and you should never use it. We need a better alternative. One that is able to detect curvature in the derivatives when integrating a value, so it is reasonably accurate even when the derivatives change over the timestep.
Rather than introduce you to the vast array of different integrators that exist, I’m going to cut to the chase and go straight to the best. This integrator is called the Runge Kutta order 4 integrator aka RK4. This is the standard integrator used for numerical integration these days and is sufficiently accurate for just about anything required in game physics, given an appropriate timestep. As far as I know its safe to just pronounce “Runge Kutta” as it reads, but its quite ugly sounding so it is usually just referred to as RK4 when speaking.
Technically RK4 has error O(5) (read “order 5″) in the Taylor’s Series expansion of the solution to the differential equation. I’m not going to go into the details of deriving or explaining how an RK4 integrator is derived, but I will say intuitively that what it is doing is detecting curvature (change over time) when integrating down to the fourth derivative. RK4 is not completely accurate, but its order of accuracy is extremely good and this is what counts.
An RK4 integrator works by evaluating the derivatives at several points in the timestep to detect this curvature, then it combines these sampled derivatives with a weighted average to get the single best approximation of the derivative that it can provide over the timestep.
Implementing RK4 is not as hard as it looks
I could present to you RK4 in form of general mathematical equations which you could use to integrate some function using its derivative f(x,t) as is usually done, but seeing as my target audience for this article are programmers, not mathematicians, I will explain using code instead. From my point of view, processes like this are best described in code anyway.
So before we go any further I will define the state of an object as a struct in C++ so that we have both position and velocity values conveniently stored in one place.
struct State
{
float x; // position
float v; // velocity
};
We’ll also need a struct to store the derivatives of the state values so we can easily pass them around to functions. The derivatives we store are velocity and acceleration:
struct Derivative
{
float dx; // derivative of position: velocity
float dv; // derivative of velocity: acceleration
};
The final piece of the puzzle that we need to implement the RK4 integrator is a function that can advance the physics state ahead from t to t+dt using one set of derivatives, then once there, recalculate the derivatives using this new state. This routine is the heart of the RK4 integrator and when implemented in C++ it looks like this:
Derivative evaluate(const State &initial, float t, float dt, const Derivative &d)
{
State state;
state.x = initial.x + d.dx*dt;
state.v = initial.v + d.dv*dt;
Derivative output;
output.dx = state.v;
output.dv = acceleration(state, t+dt);
return output;
}
It is absolutely critical that you understand what this method is doing. First it takes the current state of the object (its position and velocity) and advances it ahead dt seconds using an Euler Integration step with the derivatives that were passed in (velocity and acceleration). Once this new position and velocity are calculated, it calculates new derivatives at this point in time using the integrated state. These derivatives will be different to the derivatives that were initially passed in to the method if the derivatives are not constant over the timestep.
In order to calculate the derivatives it copies the current state velocity into the derivative struct (this is the trick to doing position and velocity integration simultaneously) then it calls the acceleration function to calculate the acceleration for the current state at the time t+dt. The acceleration function is what drives the entire simulation and in the example source code for this article I define it as follows:
float acceleration(const State &state, float t)
{
const float k = 10;
const float b = 1;
return -k * state.x - b*state.v;
}
This method calculates a spring and damper force and returns it as the acceleration assuming unit mass. What you write here of course is completely simulation dependent, but it is crucial that you structure your simulation in such a way that it can calculate the acceleration or force derivatives completely from inside this method given the current state and time, otherwise your simulation cannot work with the RK4 integrator.
Finally we get to the integration routine itself which integrates the state ahead from t to t+dt using RK4:
void integrate(State &state, float t, float dt)
{
Derivative a = evaluate(state, t, 0.0f, Derivative());
Derivative b = evaluate(state, t+dt*0.5f, dt*0.5f, a);
Derivative c = evaluate(state, t+dt*0.5f, dt*0.5f, b);
Derivative d = evaluate(state, t+dt, dt, c);
const float dxdt = 1.0f/6.0f * (a.dx + 2.0f*(b.dx + c.dx) + d.dx);
const float dvdt = 1.0f/6.0f * (a.dv + 2.0f*(b.dv + c.dv) + d.dv)
state.x = state.x + dxdt * dt;
state.v = state.v + dvdt * dt;
}
Notice the multiple calls to evaluate that make up this routine. RK4 samples derivatives four times to detect curvature instead of just once in Euler integration. The important thing to understand is how it actually does this sampling.
If you look at the code above it should be clear what is going on. Notice how it uses the previous derivative when calculating the next one. When calculating derivative b it steps ahead from t to t+dt*0.5 using derivative a, then to calculate c it steps ahead using derivative b, and finally d is calculated by stepping ahead with c. This feedback of the current derivative into the calculation of the next one is what gives the RK4 integrator its accuracy.
Once the four derivative samples have been evaluated, the best overall derivative is calculated using a weighted sum of the derivatives which comes from the Taylor Series expansion. This single value can then be used to advance the position and velocity over dt as we did before in the Euler integrator.
Notice that even when using a complicated integrator such as RK4, it all boils down into something = something + change in something * change in time. This is because differentiation and integration are fundamentally linear operations. For now we are just integrating single values, but rest assured that it all ends up like this when integrating vectors, or even quaternions and matrices for rotational dynamics.
Conclusion
I have demonstrated how to implement an RK4 integrator for a basic physics simulation. At this point if you are serious about learning game physics you should study the example source code that comes with this article and play around with it.
Here are some ideas for study:
- Switch from integrating velocity directly from acceleration to integrating momentum from force instead (the derivative of momentum is force). You will need to add “mass” and “inverseMass” to the State struct and I recommend adding a method called “recalculate” which updates velocity = momentum * inverseMass whenever it is called. Every time you modify the momentum value you need to recalculate the velocity. You should also rename the acceleration method to “force”. Why do this? Later on when working with rotational dynamics you will need to work with angular momentum directly instead of angular velocity, so you might as well start now.
- Try modifying the integrate method to implement an Euler integrator. Compare the results of the simulation against the RK4 integrator. How much can you increase the spring constant k before the simulation explodes with Euler? How large can you make it with RK4?
- Extend position, velocity and force to 3D quantities using vectors. If you use your intuition you should easily be able to extend the RK4 integrator to do this.
Once you have finished playing around with the example, lets continue to the next article where we discuss how to display your physics simulation on the screen.
浙公网安备 33010602011771号