#include <fstream>
#include "tops.h"

using namespace std;

void testTop( Top & rotor, const char* filename )
{
  Vector omega0 ( 1.0, -1.0, 0.5 );
  Matrix A0 ( 1,0,0,
 	      0,1,0,
	      0,0,1 );

  rotor.Initialization( omega0, A0 );
  
  ofstream f( filename );
  
  for ( double t = 0.0; t < 12.0; t += 0.1 ) {
    Vector omega;
    Matrix A;
    rotor.Evolution( t, omega, A );
    f << t 
      << ' ' << A.row(0)
      << ' ' << A.row(1)
      << ' ' << A.row(2)
      << ' ' << omega 
      << endl;
  }
}

int main()
{
  double Ix = 1.0, Iy = 1.5, Iz = 2.0;

  TopSpherical  sphere (Ix);          // Ix = Iy = Iz
  TopProlate    prolate(Ix, Iz);      // Ix < Iy = Iz
  TopOblate     oblate (Ix, Iz);      // Ix = Iy < Iz
  TopAsymmetric asymtop(Ix, Iy, Iz);  // Ix < Iy < Iz
  TopRecur      recurse(Ix, Iy, Iz);

  testTop( sphere,  "sphere.dat" );
  testTop( prolate, "prolate.dat" );
  testTop( oblate,  "oblate.dat" );
  testTop( asymtop, "asymtop.dat" );
  testTop( recurse, "recurse.dat" );

  return 0;
}