// integrate3d.java
//                         compile javac -cp . integrate3d.java
//                         run     java  -cp . integrate3d
//
// First by calculus and calculator
// Then using gauss legendre integration, integrate function f :
//
//  gaulegf(ax, bx, xx, wx, nx);   x's adjusted for x and weight
//  gaulegf(ay, by, yy, wy, ny);   y's adjusted for y and weight
//  gaulegf(az, bz, zz, wz, nz);   z's adjusted for z and weight
//  volume_int = 0.0;
//  for(i=1; i<=nx; i++)
//    for(j=1; j<=ny; j++)
//      for(k=1; k<=nz; k++)
//        volume_int = volume_int + wx[i]*wy[j]*wz[k]*f(xx[i],yy[j],zz[k]);
//
//
// given f(x,y,z) = x*y*z
//  
//    /               /3 /3 /3  
// 3D | f(x,y,z)dv =  |  |  | f(x,y,z)dx dy dz
//    /               /1 /1 /1        
//     volume          z  y  x        
// 
//                   do by parts
//             /3
// fx(x,y,z) = | f(x,y,z)dx  y and z constant
//             /1
//
// integrate f(x,y,z) dx = x^2/2 *y*z
// evaluate  fx(y,z) = 9/2 *y*z - 1/2 *y*z = 4*y*z
//
//
//                  /3
//        fy(y,z) = | fx(y,z)dy  x gone z constant
//                  /1
//
// integrate fx(y,z) dy = 4y^2/2 *z
// evaluate  fy(z) = 36/2 *z - 4/2 *z = 16*z
//
//                         /3
//                  F(z) = | fy(z)dz  x and y gone
//                         /1
// integrate fy(z) dz = 16z^2/2
// evaluate  fy(z) = 144/2 - 16/2 = 64
//
//
//            /3
// F(x,y,z) = | f(x,y,z)dv = 64  done
//            /1
//

public class integrate3d
{
  double volume_int = 0.0;
  int nx = 7;
  int ny = 7;
  int nz = 7;
  double xx[] = new double[nx+1];
  double wx[] = new double[nx+1];
  double yy[] = new double[ny+1];
  double wy[] = new double[ny+1];
  double zz[] = new double[nz+1];
  double wz[] = new double[nz+1];
    
  public integrate3d()
  {
    System.out.println("integrate3d.java running");
    System.out.println("given f(x,y,z) = x*y*z");
    System.out.println(" /               /3 /3 /3");  
    System.out.println(" | f(x,y,z)dv =  |  |  | f(x,y,z)dx dy dz");
    System.out.println(" /               /1 /1 /1");
    System.out.println("  volume          z  y  x");        
    System.out.println(" ");

    new gaulegf(1.0, 3.0, xx, wx, nx);  // x's adjusted for x and weight
    new gaulegf(1.0, 3.0, yy, wy, ny);  // y's adjusted for y and weight
    new gaulegf(1.0, 3.0, zz, wz, nz);  // z's adjusted for z and weight
    volume_int = 0.0;
    for(int i=1; i<=nx; i++)
      for(int j=1; j<=ny; j++)
        for(int k=1; k<=nz; k++)
          volume_int = volume_int + wx[i]*wy[j]*wz[k]*f(xx[i],yy[j],zz[k]);

    System.out.println("volume_int ="+volume_int);
    System.out.println("using numerical integration");
    System.out.println(" ");
    System.out.println("close enough to calculus with calculator 64");
    System.out.println("integrate3d.java finished");    
  }

  double f(double x, double y, double z)
  {
    return x*y*z;
  }
    
  public static void main (String[] args)
  {
    new integrate3d();
  }
} // end integrate3d.java