//
//	Collapse2.java:
//
//	Simulation of Twin Tower collapse times using a model that embodies
//	the "pile driver" concept embraced by Bazant and Zhou
//	(despite the clear evidence the tops disintegrated)
//
//	This simulation generalizes the earlier simulation, Collapse1.java,
//	by adding parameters to model the following:
//
//	* A specified fraction of the mass falling within the Tower's
//	  profile moves to outside the profile during each story collapse,
//	  and thereafter does not participate in accelerating mass downward.
//	* The increasing mass of stories lower in the Tower can be specified.
//
//	With those generatlizations, Collapse2.java retains the
//	following assumptions of Collapse1.java
//
//	* Each floor is an infinately thin slab, and all the mass
//	  of a story is concentrated in the slab.
//	* Mass is uniformally distributed among the stories.
//	* The overhanging portion (eg: 14 floors in the North Tower)
//	  falls as a block, with its bottom floor accumulating pancaked slabs
//	  of the once-intact floors as it encounters them.
//	* Once the bottom of the block reaches the ground, the floors in it
//	  start to pancake from bottom to top, the roof of the tower
//	  falling at freefall at that point.
//	* The falling block remains perfectly centered over he intact portion.
//	* The accumulation of floors is inelastic.
//	* Each floor's support vanishes when touched by the falling block.
//	* Momentum is conserved.
//	* None of the kinetic energy of the falling mass is diverted to other
//	  sinks (concrete pulverization, steel bending, etc.)
//
//
//	Equation for elapsed time, given initial velocity v, height h,
//	and gravitational acceleration g:
//
//	t = - 2*v/g + sqrt(v^2/g^2 + 2*h/g)
//

class Collapse2 {
				// floors
  int numberFloors = 110;	// number of floors
  int startFloor;		// floor at which collapse begins
  int stopFloor;		// floor to collapse to
  int currentFloor;		// current topmost floor of intact portion
  double thicken = 0d;		// mass of top floor as fraction of bottom
  double dispersion = 0d;	// fraction of mass dispersed per floor
				// constants
  static double g = 32;		// gravitational acceleration in feet/sec^2
  double h = 12.4;		// distance between floors
  //
  //	utility functions
  //
  public static double sqr(double x) {
    return x*x;
  }
  public static double sqrt(double x) {
    return Math.pow(x,0.5);
  }
  //
  //	compute time of one-story fall
  //	given initial velocity v, and height of story h
  //
  static double timeFall(double v, double h) {
    return - v/g + sqrt(sqr(v)/sqr(g) + 2*h/g);
  }
  //
  //	compute velocity at end of fall
  //	given initial velocity v, and elapsed time of fall t
  //
  static double velocFall(double v, double t) {
    return v + t*g;
  }
  //
  //	compute distance of fall
  //	given an elapsed time t
  //
  static double distanceFall(double t) {
    return 0.5*sqr(t)*g;
  }
  //
  //	constructor creates an instance
  //	given the floor at which the collapse starts
  //
  Collapse2(int start) {
    startFloor = start;
  }
  //
  //
  //
  double floorPosition(int floor) {
    return (double)floor/(double)numberFloors;
  }
  double floorMass(int floor) {
    return 1d + thicken * (1 - floorPosition(floor));
  }
  double floorsMass(int top, int bottom) {
    double m = 0d;
    for (int i = top; i >= bottom; i--)
      m += floorMass(i);
    return m;
  }
  //
  //	timeCollapse returns the time that it takes the roof
  //	to reach the ground.
  //
  double timeCollapse() {	// time from start until roof reaches ground
    double v = 0d;		// set initial velocity to zero
    double m;			// mass of falling portion over footprint
    double e = 0d;		// mass of rubble outside footprint
    double ts = 0d;		// set elapsed time to zero
    double t;			// time
    double mc;

    m = floorsMass(numberFloors,startFloor); // mass of block
    for (currentFloor = startFloor;
	 currentFloor >= stopFloor;
	 currentFloor--) {
      t = timeFall(v,h);	// time to fall story
      v = velocFall(v,t);	// velocity after falling story
      mc = floorMass(currentFloor); // mass of current story
      v = v*m/(m + mc);		// velocity after inelastic impact
      ts += t;			// update total elapsed time.
      m *= 1d - dispersion;	// remove mass dispersed from pile-driver
      e += m*dispersion;	// add to dispersed mass
      m += mc;			// increment falling mass to include story
    }
    if (stopFloor != 0) {
      System.out.println("crash zone to floor "+ stopFloor +": "+ ts);
      return ts;
    }
				// time for bottom of falling portion of
				// building to reach the ground
    System.out.println("crash zone to ground: "+ ts);

    t = timeFall(v,h*(numberFloors - startFloor));
    ts += t;
				// time for roof to reach ground
    System.out.println("roof to ground: "+ ts);

    System.out.println("mass dispersion ratio: "+ e/(e + m));

    return ts;
  }
  //
  //
  void readArgs(String[] argv) {
    int narg = argv.length;
    for (int iarg = 2; iarg < narg; iarg += 2) {
	String flag = argv[iarg];
	String val = argv[iarg + 1];
	if        (flag.equals("--stop")) {
	  stopFloor = Integer.valueOf(val).intValue();
	} else if (flag.equals("--thicken")) {
	  thicken = Double.valueOf(val).doubleValue();
	} else if (flag.equals("--dispersion")) {
	  dispersion = Double.valueOf(val).doubleValue();
	}
    }
  }
  //
  //
  public static void main(String[] argv) {
    if (argv.length < 1) {
      System.out.println("syntax: \n"+
	"\t java Collapse2 -t <height> <velocity>\n"+
	"\t\t prints time an object in freefall takes to reach the ground\n"+
	"\t\t <height>: height in feet at start\n"+
	"\t\t <velocity>: velocity in feet/second at start\n"+
	"\t java Collapse2 -d <time>\n"+
	"\t\t prints distance an object in freefall drops in a time interval\n"+
	"\t\t <time>: time interval in seconds\n"+
	"\t java Collapse2 -c <start-floor> <options>\n"+
	"\t\t prints time to total collapse of a WTC tower\n"+
	"\t\t first prints time bottom of falling top hits ground\n"+
	"\t\t next prints time roof hits ground\n"+
	"");
      System.out.println("\noptions: \n"+
	"\t --start <start-floor>\n"+
	"\t\t floor at which collapse starts, eg: 82, 96\n"+
	"\t --stop <stop-floor>\n"+
	"\t\t floor at which collapse stops, eg: 0\n"+
	"\t --thicken <mass-thicken>\n"+
	"\t\t fractional mass of bottom floor beyond top, eg: 0.5\n"+
	"\t --dispersion <mass-dispersion>\n"+
	"\t\t fractional dispersion of mass outside footprint, eg: 0.1\n"+
	"");
    } else {
      Collapse2 that;
      char c;
      c = argv[0].charAt(1);
      switch (c) {
       case 't':
	double h = Double.valueOf(argv[1]).doubleValue();
	double v = Double.valueOf(argv[2]).doubleValue();
	System.out.println(
	  "time "+ timeFall(v,h) +" for height "+ h +" starting "+ v);
       break;
       case 'd':
        double t = Double.valueOf(argv[1]).doubleValue();
	System.out.println(
	  "distance "+ distanceFall(t) +" for time "+ t);
       break;
       case 'c':
	int start = Integer.valueOf(argv[1]).intValue();
	that = new Collapse2(start);
	that.readArgs(argv);
	that.timeCollapse();
       break;
       default:
	System.out.println(
	  "invalid option: "+ c);
      }
    }
  }
}