Projectiles - Physics for Game Programmers - page 126

// The getRightHandSide() method returns the right-hand
// sides of the six first-order projectile ODEs.
// q[0] = vx = dxdt
// q[1] = x
// q[2] = vy = dydt
// q[3] = y
// q[4] = vz = dzdt
// q[5] = z
public double[] getRightHandSide(double s, double q[],
double deltaQ[], double ds,
double qScale) {
double dQ[] = new double[6];
double newQ[] = new double[6];
// Compute the intermediate values of the
// dependent variables.
for(int i=0; i<6; ++i) {
newQ[i] = q[i] + qScale*deltaQ[i];
}
// Declare some convenience variables representing
// the intermediate values of velocity.
double vx = newQ[0];
double vy = newQ[2];
double vz = newQ[4];
// Compute the apparent velocities by subtracting
// the wind velocity components from the projectile
// velocity components.
double vax = vx - getWindVx();
double vay = vy - getWindVy();
double vaz = vz;
// Compute the apparent velocity magnitude. The 1.0e-8 term
// ensures there won't be a divide by zero later on
// if all of the velocity components are zero.
double va = Math.sqrt(vax*vax + vay*vay + vaz*vaz) + 1.0e-8;
The directional components of drag force are computed separately, which will clean up
the right-hand-side code a little later in the method.
// Compute the total drag force and the dirctional
// drag components.
double Fd = 0.5*getDensity()*getArea()*getCd()*va*va;
double Fdx = -Fd*vax/va;
double Fdy = -Fd*vay/va;
double Fdz = -Fd*vaz/va;