// test_lsfit.java

import java.io.*;
import java.text.*;

public class test_lsfit
{
  int nx = 5;
  int ny = 5;
  int nz = 5;
  int dof = (nx-2)*(ny-2)*(nz-2);  // Degrees of Freedom                        
  double hx, hy, hz;
  double xg[] = new double[nx]; // X grid, does not need to be uniform          
  double yg[] = new double[ny];
  double zg[] = new double[nz];
  double us[] = new double[dof];
  double ug[][][] = new double[nx][ny][nz]; 
  double uA[][][] = new double[nx][ny][nz]; 
  double xmin = -1.0;
  double xmax =  1.0;
  double ymin = -1.0;
  double ymax =  1.0;
  double zmin = -1.0;
  double zmax =  1.0;
  double x, y, z;
  DecimalFormat f4 = new DecimalFormat("0.000");

  public test_lsfit()
  {
    System.out.println("test_lsfit.java");
    System.out.println("fitting");
    System.out.println("x*x*x*x+2.0*y*y*y*y+3.0*z*z*z*z+4.0*x*y*z+5.0*x*y+z+2.0");
    System.out.println(" ");
    // initialize hx, hy, hy, xg, yg, zg                
    hx = (xmax-xmin)/(double)(nx-1);
    hy = (ymax-ymin)/(double)(ny-1);
    hz = (zmax-zmin)/(double)(nz-1);
    for(int i=0; i<nx; i++) xg[i] = xmin+hx*(double)i;
    for(int ii=0; ii<ny; ii++) yg[ii] = ymin+hy*(double)ii;
    for(int iii=0; iii<nz; iii++) zg[iii] = zmin+hz*(double)iii;

    for(int i=0; i<nx; i++)
    {
      x = xg[i];
      for(int ii=0; ii<ny; ii++)
      {
	y = yg[ii];
	for(int iii=0; iii<nz; iii++)
	{
	  z = zg[iii];
	  uA[i][ii][iii] = ub(x,y,z);
	  System.out.println("xg="+f4.format(x)+", yg="+f4.format(y)+
			     ", zg="+f4.format(z)+
                             ", uA["+i+"]["+ii+"]["+iii+"]="+
			     f4.format(uA[i][ii][iii]));
	  if(i==0 || ii==0 || i==nx-1 || ii==ny-1 || iii==0 || iii==nz-1)
	  {
	    ug[i][ii][iii] = ub(x,y,z);
	  }
	  else // not boundary                                                  
	  {
	    ug[i][ii][iii] = 0.0;
	  }
	}
      }
    }

    lsfit LS = new lsfit(3, 4, dof); // dimension, power, points                
    LS.fit3_bpoly(nx, ny, nz, xg, yg, zg, ug);
    for(int i=1; i<nx-1; i++)
    {
      x = xg[i];
      for(int ii=1; ii<ny-1; ii++)
      {
	y = yg[ii];
	for(int iii=1; iii<nz-1; iii++)
	{
	  z = zg[iii];
	  us[S(i,ii,iii)] = LS.eval3_poly(x, y, z);
          System.out.println("want="+f4.format(uA[i][ii][iii])+
                             "  fit="+f4.format(us[S(i,ii,iii)])+
			     "  err="+(us[S(i,ii,iii)]-uA[i][ii][iii]));
	} // iii
      } // ii
    } // i
    System.out.println("test_lsfit.java ");

    System.out.println("finished test_lsfit.java");
  } // end test_lsfit

  double ub(double x, double y, double z)
  {
    return x*x*x*x+2.0*y*y*y*y+3.0*z*z*z*z+4.0*x*y*z+5.0*x*y+z+2.0;
  }

  int S(int i, int ii, int iii) // linear index of non boundary dof, zero base  
  {
    return (i-1)*(ny-2)*(nz-2)+(ii-1)*(nz-2)+(iii-1);
  }

  public static void main (String[] args)
  {
    new test_lsfit();
  } // end main of test_lsfit.java
} // end class test_lsfit