///////////////////////////////////////////////////////////////////////////
//                                                                       //
// Program file name: Jump.java                                          //
//                                                                       //
// © Tao Pang 2006                                                       //
//                                                                       //
// Last modified: January 18, 2006                                       //
//                                                                       //
// (1) This Java program is part of the book, "An Introduction to        //
//     Computational Physics, 2nd Edition," written by Tao Pang and      //
//     published by Cambridge University Press on January 19, 2006.      //
//                                                                       //
// (2) No warranties, express or implied, are made for this program.     //
//                                                                       //
///////////////////////////////////////////////////////////////////////////

// An example of modeling a motorcycle jump through the
// two-point predictor-corrector scheme.

import java.lang.*;
public class Jump {
  static final int n = 100, j = 2;
  public static void main(String argv[]) {
    double x[] = new double[n+1];
    double y[] = new double[n+1];
    double vx[] = new double[n+1];
    double vy[] = new double[n+1];
    double ax[] = new double[n+1];
    double ay[] = new double[n+1];

 // Assign all the constants involved
    double g = 9.80;
    double angle = 42.5*Math.PI/180; 
    double speed = 67;
    double mass = 250;
    double area = 0.93;
    double density = 1.2;
    double k = area*density/(2*mass);
    double dt = 2*speed*Math.sin(angle)/(g*n);
    double d = dt*dt/2;

 // Assign the quantities for the first two points
    x[0] = y[0] = 0;
    vx[0] = speed*Math.cos(angle);
    vy[0] = speed*Math.sin(angle);
    double v = Math.sqrt(vx[0]*vx[0]+vy[0]*vy[0]);
    ax[0] = -k*v*vx[0];
    ay[0] = -g-k*v*vy[0];
    double p = vx[0]*ax[0]+vy[0]*ay[0];
    x[1] = x[0]+dt*vx[0]+d*ax[0];
    y[1] = y[0]+dt*vy[0]+d*ay[0];
    vx[1] = vx[0]+dt*ax[0]-d*k*(v*ax[0]+p*vx[0]/v);
    vy[1] = vy[0]+dt*ay[0]-d*k*(v*ay[0]+p*vy[0]/v);
    v = Math.sqrt(vx[1]*vx[1]+vy[1]*vy[1]);
    ax[1] = -k*v*vx[1];
    ay[1] = -g-k*v*vy[1];

 // Calculate other position and velocity recursively
    double d2 = 2*dt;
    double d3 = dt/3;
    for (int i=0; i<n-1; ++i) {

   // Predict the next position and velocity
      x[i+2] = x[i]+d2*vx[i+1];
      y[i+2] = y[i]+d2*vy[i+1];
      vx[i+2] = vx[i]+d2*ax[i+1];
      vy[i+2] = vy[i]+d2*ay[i+1];
      v = Math.sqrt(vx[i+2]*vx[i+2]+vy[i+2]*vy[i+2]);
      ax[i+2] = -k*v*vx[i+2];
      ay[i+2] = -g-k*v*vy[i+2];
      
   // Correct the new position and velocity
      x[i+2] = x[i]+d3*(vx[i+2]+4*vx[i+1]+vx[i]);
      y[i+2] = y[i]+d3*(vy[i+2]+4*vy[i+1]+vy[i]);
      vx[i+2] = vx[i]+d3*(ax[i+2]+4*ax[i+1]+ax[i]);
      vy[i+2] = vy[i]+d3*(ay[i+2]+4*ay[i+1]+ay[i]);
    }

 // Output the result in every j time steps
    for (int i=0; i<=n; i+=j)
      System.out.println(x[i] +" " + y[i]);
  }
}