//======================================================
// Numerical Analysis package in java
// gauss integration of two dimensional function
// and differentiation (derivative) functions
// Dr. Turhan Coban
// =====================================================
import java.io.*;

abstract class fdi
{
public double a,b;
public fdi(double ai,double bi)
{a=ai;b=bi;}
public fdi()
{a=0;b=1;}
public void setab(double ai,double bi)
{a=ai;b=bi;}	
abstract public  double c(double x);
abstract public  double d(double x);
abstract public  double f(double x[]); 
}

class f1 extends fdi 
{
  public f1(double ai,double bi)
  {super(ai,bi);}
  
  public double c(double x)
  {return 1.0;}
  
  public double d(double x)
  {return 1.5;}
  
  public double f(double x[]) 
  {
  //çözümü istenen fonksiyon	
  double ff;
  ff=Math.log(x[0]+2.0*x[1]);
  return ff; 
  }
}  


class integ_gauss3
{
	
public static double[][] gauss_legendre_coefficients(double x1,double x2,int n)
{
	//calculates legendre gauss-coefficients as coefficients of the integral
	//for n terms
	double EPS=3.0e-15;
	int m,j,i;
	double z1,z,xm,xl,pp,p3,p2,p1;
    double a[][]=new double[2][n];//a[0][i]=x[i]   a[1][i]=w[i]
	m=(n+1)/2;
	xm=0.5*(x2+x1);
	xl=0.5*(x2-x1);
	for (i=1;i<=m;i++)  {
		z=Math.cos(Math.PI*((i-0.25)/(n+0.5)));
		do {
			p1=1.0;
			p2=0.0;
			for (j=1;j<=n;j++) {
				p3=p2;
				p2=p1;
				p1=((2.0*j-1.0)*z*p2-(j-1.0)*p3)/j;
			}
			pp=n*(z*p1-p2)/(z*z-1.0);
			z1=z;
			z=z1-p1/pp;
		} while (Math.abs(z-z1) > EPS);
		a[0][i-1]=xm-xl*z;
		a[0][n-i]=xm+xl*z;
		a[1][i-1]=2.0*xl/((1.0-z*z)*pp*pp);
		a[1][n-i]=a[1][i-1];
			}
return a;
}	

public static double I(fdi f,int n,int m) 
{
//integral f(x)dx
//integral of a function by using gauss-legendre quadrature
//coefficients are pre-calculated for n,m terms for [-1,1]
//band then utilises variable transform
 int i,j;
 double rn[],cn[];
 double rm[],cm[];
 rn=new double[n];
 rm=new double[m];
 cn=new double[n];
 double a[][]=new double[2][n];
 a=gauss_legendre_coefficients(-1.0,1.0,n);
 rn=a[0];
 cn=a[1];
 a=gauss_legendre_coefficients(-1.0,1.0,m);
 rm=a[0];
 cm=a[1];
double x[],h1,h2,J,a1,b1,c1,d1,k1,k2,JX,Q;
x=new double[2];
a1=f.a;
b1=f.b;
h1=(b1-a1)/2.0;
h2=(b1+a1)/2.0;
J=0;
for(i=0;i<m;i++)
{
JX=0;
x[0]=h1*rm[i]+h2;
d1=f.d(x[0]);
c1=f.c(x[0]);
k1=(d1-c1)/2.0;
k2=(d1+c1)/2.0;
//4=========================
    for(j=0;j<n;j++)
    {
	x[1]=k1*rn[j]+k2;
	Q=f.f(x);
	JX+=cn[j]*Q;
	}
//System.out.println("JX="+JX);	
//5=========================	
J=J+cm[i]*k1*JX;	
}
J=h1*J;
return J;
}

public static double I(fdi f,int n)
{return I(f,n,n);}

public static double I(fdi f)
{return I(f,60,60);}

  public static void main(String args[]) throws IOException
  {
  f1 ff=new f1(1.4,2.0);	  
  System.out.println("integral of class         f1 : "+I(ff));
}
}
