// Newton metodu ile optimizasyon
// 
import visad.*;
import visad.java2d.DisplayImplJ2D;
import visad.java3d.DisplayImplJ3D;
import visad.util.*;
import java.rmi.RemoteException;
import java.awt.*;
import javax.swing.*;
import java.awt.event.*;


class f1 extends f_xj
{
public double func(double x[]) 
{
//çözümü istenen fonksiyon	
double ff;
ff=2.0*x[0]*x[1]+2.0*x[0]-x[0]*x[0]-2.0*x[1]*x[1];
return ff; //minimum değil maxsimum istendiğinden - ile çarptık
}
}


public class OPO6
{
// 1. Newton metodu
//________________________________________________________________
public static double[] newtond( f_xj f,double x[])
{
// lineer olmıyan denklem siteminin kökleri
// fj(xi)=0 lineer olmıyan denklem sistemi için çözüm matrisi 
// türev matrisi dfj/dxi turev alt programıyla hesaplanır.
// Newton-Raphson metodu ile lineer olmıyan denklem sistemi çözümü 
//ti :önem katsayısı
//x bağımsız değişken vektörü
//y bağımlı değişken vektörü
//dy bağımlı değişkenin türev fonksiyonları vektörü
double ti=1.0;
int i,ii,jj;
int nmax=500;
double tolerance=1.0e-10;
int n=x.length;
double b[];
b=new double[n];
double dy[][]=H(f,x);
//dy=new double[n][n];
i=0;
for(i=0;i<n;i++)
{
b[i]=1.0;
}
while( i++ < nmax && Matrix.abs(b) > tolerance )
{
    b=Matrix.multiply(Matrix.divide(df(f,x),dy),-ti);
    x=Matrix.add(x,b);
    System.out.println("i= "+i+"\n" );
    for(int kk=0;kk<x.length;kk++)
        System.out.println(x[kk]+" ");
}
if(i >= nmax) JOptionPane.showMessageDialog(null,"Uyarı maximum iterasyon sayısı aşıldı \n"+
  " çözüm geçerli olmıyabilir","MAKSİMUM ITERASYON SAYISI UYARISI",JOptionPane.WARNING_MESSAGE);
return x;
}
public static double turev(f_xj f,double x[],int denklem_ref,int x_ref)
{
// 
// verilen fonksiyonun ikinci türevlerinin matrisi
// fonksiyon f in denklem_ref sayılı fonksiyonunun(sayılar 0 dan başlar)
// x[x_ref] değişkenine göre türevi (sayılar 0 dan başlar) 
// df_denklem_ref(x)/d_x_ref
// bu metod newtond metodu içinde kullanılmak içindir.
double h0=0.256808;
int i,m;
int n=7;
double f1[];
f1=new double[x.length];
double f2[];
f2=new double[x.length];
double x1[];
x1=new double[x.length];
double x2[];
x2=new double[x.length];
for(i=0;i<x.length;i++)
{
x1[i]=x[i];
x2[i]=x[i];
}
//turev of a simple function
double T[][];
T=new double[n][n];
double h[];
h=new double[n];
//vector<double> h(n,0);
for(i=0;i<n;i++)
{
  h[i]=0;
  for(int j=0;j<n;j++)
    T[i][j]=0;
}
h[0]=h0;
double r=0.5;
for( i=1;i<n;i++)
{
h[i]=h0*Math.pow(r,i);
}

for(i=0;i<n;i++)
{
x1[x_ref]+=h[i];
x2[x_ref]-=h[i]*1.0094847;
f1=df(f,x1);
f2=df(f,x2);
T[i][0]=( f1[denklem_ref] - f2[denklem_ref])/(2.0*h[i]);
x1[x_ref]=x[x_ref];
x2[x_ref]=x[x_ref];
}
for(m=1;m<n;m++)
{
  for(i=0;i<n-m;i++)
  {
  T[i][m]=(h[i]*h[i]*T[i+1][m-1] - h[i+m]*h[i+m]*T[i][m-1])/(h[i]*h[i]
  - h[i+m]*h[i+m]);
  }
}
double xx=T[0][n-1];
return xx;
}

public static double[][]H(f_xj f,double x[])
{
  int n=x.length;	
  double dy[][]=new double[n][n];;
  for(int ii=0;ii<n;ii++)
  {
    for(int jj=0;jj<n;jj++)
    {   //Hessian matrisi
        dy[ii][jj]=turev(f,x,ii,jj);
    }
  }
  return dy;
}

public static double turev(f_xj f,double xxi[],int x_ref)
{
// df/dxj
double h0=0.256808;
int i,m;
int n=7;
double x1[];
x1=new double[xxi.length];
double x2[];
x2=new double[xxi.length];
for(i=0;i<xxi.length;i++)
{
x1[i]=xxi[i];
x2[i]=xxi[i];
}
//turev of a simple function
double T[][];
T=new double[n][n];
double h[];
h=new double[n];
for(i=0;i<n;i++)
{
  h[i]=0;
  for(int j=0;j<n;j++)
    T[i][j]=0;
}
h[0]=h0;
double r=0.5;
for( i=1;i<n;i++)
{
//h[i]=h0*Math.pow(r,i);
h0/=2.0;
h[i]=h0;
}
//first turev (difference formula)
for(i=0;i<n;i++)
{
x1[x_ref]+=h[i];
x2[x_ref]-=h[i];
T[i][0]=(f.func(x1)-f.func(x2))/(2.0*h[i]);
x1[x_ref]=xxi[x_ref];
x2[x_ref]=xxi[x_ref];
}

for(m=1;m<n;m++)
{
  for(i=0;i<n-m;i++)
  {
  T[i][m]=(h[i]*h[i]*T[i+1][m-1] -
  h[i+m]*h[i+m]*T[i][m-1])/(h[i]*h[i] - h[i+m]*h[i+m]);
  }
}
double xx=T[0][n-1];
return xx;
}
public static double df(f_x f,double x)
{
// türev alt programı	
// df fonksiyonunun türevini hesaplar
double h0=0.0256;
int i,m;
int n=7;
//turev of a simple function
double T[][];
T=new double[n][n];
double h[];
h=new double[n];
//vector<double> h(n,0);
for(i=0;i<n;i++)
{
  h[i]=0;
  for(int j=0;j<n;j++)
    T[i][j]=0;
}
h[0]=h0;
double r=0.5;
for( i=1;i<n;i++)
{
h[i]=h0*Math.pow(r,i);
}

for(i=0;i<n;i++)
{
T[i][0]=( f.func(x + h[i]) - f.func( x - h[i]))/(2.0*h[i]);
}
for(m=1;m<n;m++)
{
  for(i=0;i<n-m;i++)
  {
  T[i][m]=(h[i]*h[i]*T[i+1][m-1] - h[i+m]*h[i+m]*T[i][m-1])/
          (h[i]*h[i]- h[i+m]*h[i+m]);
  }
}
double xx=T[0][n-1];
return xx;
}
public static double[] df(f_xj f,double xxi[])
{
//Gredyen vektörü
// f(x0,x1,x2,...xn) founksiyonunun türev vektörü
// [df/dx0, df/dx1,df/dx2,....,df/dxn]
// df/dxj  j=0...x.length
//
// bu fonksiyon bilgisayar tarafından hesaplanmaktadır
// kullanıcı tarafından tanımlanan f(x) fonksiyonunun
// sayısal olarak alınan türevidir.	

// 
double a[]=new double[xxi.length];
for(int i=0;i<xxi.length;i++)
   {
   a[i]=turev(f,xxi,i);
   }
return a;   	   
}
     public static void main (String args[]) throws RemoteException, VisADException 
     {  
	 String s="kök tahmin değerlerini arada bir boşluk bırakarak giriniz: ";
     double [] x0=opt.verigir(s);
     f1 f=new f1();
     double [] r1= newtond(f,x0);
     s=" çözüm seti  : \n"+Matrix.toStringT(r1);
     String s2="Newton metodu ile lineer olmıyan çok değişkenli fonksiyon optimizasyonu:";
     JOptionPane.showMessageDialog(null,s,s2,JOptionPane.PLAIN_MESSAGE);
     double xmin[]=new double[2];
	double xmax[]=new double[2];
	int NROWS=200;
	int NCOLS=200;
	xmin[0]=-2.0;
	xmin[1]=-2.0;
	xmax[0]= 2.0;
	xmax[1]= 2.0;
	f2 fa=new f2();
    Plot3D.func(f,0.0,10.0,0.0,10.0);
	//Plot_3D.contour_Plot(f,xmin,xmax,NROWS,NCOLS);
     //System.exit(0);
     } 
}