// tetra_cube.java  develope 5 tetra make a cube

import java.io.*;

public class tetra_cube
{
  boolean debug = true;
  int pts, pt, np, j;
  int k = 0; // index into newpts k++
  datread DR; // = new datread(filename); // defines data structures
  int new_points = 0; // for each x,y,z point kept
  int new_polys = 0; // k/4 for datwrite
  double xyz[][] = new double[10000][3]; // x, y, z for output
  int newpts[] = new int[40000]; // 3, xyz1, xyz2, xyz3 for output
  String fileout = "tetra_cube.dat";
  double x, y, z, dz;
  int box = 1;
  double vol = 0.0;
    
  public  tetra_cube()
  {
    System.out.println("tetra_cube.java running");
    new_polys = 0;
    new_points = 0;

    System.out.println(" ");
    System.out.println("generate vertex");
    xyz[new_points][0] =-1.0; // 1
    xyz[new_points][1] =-1.0;
    xyz[new_points][2] =-1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;
    xyz[new_points][0] =-1.0; // 2
    xyz[new_points][1] =-1.0;
    xyz[new_points][2] = 1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;
    xyz[new_points][0] =-1.0; // 3
    xyz[new_points][1] = 1.0;
    xyz[new_points][2] =-1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;
    xyz[new_points][0] =-1.0; // 4
    xyz[new_points][1] = 1.0;
    xyz[new_points][2] = 1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;
    xyz[new_points][0] = 1.0; // 5
    xyz[new_points][1] =-1.0;
    xyz[new_points][2] =-1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;
    xyz[new_points][0] = 1.0; // 6
    xyz[new_points][1] = 1.0;
    xyz[new_points][2] =-1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;
    xyz[new_points][0] = 1.0; // 7
    xyz[new_points][1] =-1.0;
    xyz[new_points][2] = 1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;
    xyz[new_points][0] = 1.0; // 8
    xyz[new_points][1] = 1.0;
    xyz[new_points][2] = 1.0;
    System.out.println("pt "+(new_points+1)+" x="+xyz[new_points][0]+" y="+xyz[new_points][1]+" z="+xyz[new_points][2]);
    new_points = new_points +1;

    System.out.println("tetra 1"); // 1, 2, 3, 5
    newpts[k++] = 3;
    newpts[k++] = 1;
    newpts[k++] = 5;
    newpts[k++] = 3;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;  // 1b
    newpts[k++] = 1;
    newpts[k++] = 2;
    newpts[k++] = 3;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 1c
    newpts[k++] = 1;
    newpts[k++] = 5;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 1d
    newpts[k++] = 3;
    newpts[k++] = 5;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);

    System.out.println("tetra 2"); // 2, 5, 7, 8
    newpts[k++] = 3;
    newpts[k++] = 7;
    newpts[k++] = 5;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;  // 2b
    newpts[k++] = 7;
    newpts[k++] = 2;
    newpts[k++] = 8;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 2c
    newpts[k++] = 7;
    newpts[k++] = 5;
    newpts[k++] = 8;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 2d
    newpts[k++] = 2;
    newpts[k++] = 5;
    newpts[k++] = 8;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);

    System.out.println("tetra 3"); // 3, 5, 6, 8
    newpts[k++] = 3;
    newpts[k++] = 6;
    newpts[k++] = 3;
    newpts[k++] = 8;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;  // 3b
    newpts[k++] = 6;
    newpts[k++] = 3;
    newpts[k++] = 5;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 3c
    newpts[k++] = 6;
    newpts[k++] = 8;
    newpts[k++] = 5;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 3d
    newpts[k++] = 3;
    newpts[k++] = 8;
    newpts[k++] = 5;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);

    System.out.println("tetra 4"); // 2, 3, 4, 8
    newpts[k++] = 3;
    newpts[k++] = 4;
    newpts[k++] = 3;
    newpts[k++] = 8;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;  // 4b
    newpts[k++] = 4;
    newpts[k++] = 8;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 4c
    newpts[k++] = 4;
    newpts[k++] = 3;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 4d
    newpts[k++] = 3;
    newpts[k++] = 8;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);

    System.out.println("tetra 5 inside"); // 2, 3, 5, 8
    newpts[k++] = 3;
    newpts[k++] = 2;
    newpts[k++] = 3;
    newpts[k++] = 8;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;  // 5b
    newpts[k++] = 2;
    newpts[k++] = 3;
    newpts[k++] = 5;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 5c
    newpts[k++] = 5;
    newpts[k++] = 8;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3; // 5d
    newpts[k++] = 5;
    newpts[k++] = 8;
    newpts[k++] = 3;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);

    System.out.println("new_points="+new_points);
    System.out.println("new_polys="+new_polys);

    new datwrite(fileout, xyz, new_points, newpts, new_polys);
    System.out.println("file "+fileout+" written");

    System.out.println("check that volume of 5 tetrahedron = cube volume");
    vol = volume(1, 2, 3, 5);
    System.out.println("tetra 1 volume(1, 2, 3, 5)="+volume(1, 2, 3, 5));
    vol = vol + volume(2, 5, 7, 8);
    System.out.println("tetra 2 volume(2, 5, 7, 8)="+volume(2, 5, 7, 8));
    vol = vol + volume(3, 5, 6, 8);
    System.out.println("tetra 3 volume(3, 5, 6, 8)="+volume(3, 5, 6, 8));
    vol = vol + volume(3, 4, 5, 8);
    System.out.println("tetra 4 volume(2, 3, 4, 8)="+volume(2, 3, 4, 8));
    vol = vol + volume(2, 3, 5, 8);
    System.out.println("tetra 5 volume(2, 3, 5, 8)="+volume(2, 3, 5, 8));
    System.out.println("sum volume of 5 tetrahedron ="+vol);
    System.out.println("tetra 5 could be two tetra with half volume");
    				    
    k = 0;
    new_polys = 0;
    System.out.println("inside, tetra 5");
    newpts[k++] = 3;
    newpts[k++] = 3;
    newpts[k++] = 2;
    newpts[k++] = 8;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;
    newpts[k++] = 3;
    newpts[k++] = 2;
    newpts[k++] = 5;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;
    newpts[k++] = 5;
    newpts[k++] = 8;
    newpts[k++] = 2;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    newpts[k++] = 3;
    newpts[k++] = 5;
    newpts[k++] = 8;
    newpts[k++] = 3;
    new_polys = new_polys + 1;
    System.out.println("tetra "+(new_polys)+"   "+newpts[k-4]+", "+newpts[k-3]+" ,"+newpts[k-2]+" ,"+newpts[k-1]);
    fileout = "tetra_cube2.dat";
    new datwrite(fileout, xyz, new_points, newpts, new_polys);
    System.out.println("file "+fileout+" written");
    
    System.out.println("tetra_cube.java finished");
    
  } // end tetra_cube

  double volume(int p1, int p2, int p3, int p4)
  {
    double A[][] = new double[4][4];
    double vol1;
    A[0][0] = xyz[p1-1][0]; 
    A[0][1] = xyz[p1-1][1]; 
    A[0][2] = xyz[p1-1][2]; 
    A[0][3] = 1.0; 
    A[1][0] = xyz[p2-1][0]; 
    A[1][1] = xyz[p2-1][1]; 
    A[1][2] = xyz[p2-1][2]; 
    A[1][3] = 1.0; 
    A[2][0] = xyz[p3-1][0]; 
    A[2][1] = xyz[p3-1][1]; 
    A[2][2] = xyz[p3-1][2];
    A[2][3] = 1.0; 
    A[3][0] = xyz[p4-1][0]; 
    A[3][1] = xyz[p4-1][1]; 
    A[3][2] = xyz[p4-1][2];
    A[3][3] = 1.0; 
    vol1 = Math.abs(Matrix.determinant(A)/6.0);
    return vol1;
  } // end volume of tetrahedron
  
  public static void main(String[] args)
  {
    new tetra_cube(); // writes tetra_cube.dat
  }
} // end class  tetra_cube.java