import java.io.*;
import javax.swing.*;
//import Numeric;

class SCO13A
{ 

public static double[][] tridiagonal(double b[][], double d[], double e[])
{
//reduces matrix to tridiaonal form by using HOUSEHOLDER TRANSFORMATION
//this method is used by QL method to calculate eigen values
//and eigen  vectors of a symmetric matrix
	int l,k,j,i;
    int n=b.length;
	double scale,hh,h,g,f;
    double a[][]=new double[n+1][n+1];
    double c[][]=new double[n][n];
    for(i=0;i<n;i++)
      for(j=0;j<n;j++)
         a[i][j]=b[i][j];
	for (i=n;i>=2;i--) {
		l=i-1;
		h=scale=0.0;
		if (l > 1) {
			for (k=1;k<=l;k++)
				scale += Math.abs(a[i-1][k-1]);
			if (scale == 0.0)
				e[i-1]=a[i-1][l-1];
			else {
				for (k=1;k<=l;k++) {
					a[i-1][k-1] /= scale;
					h += a[i-1][k-1]*a[i-1][k-1];
				}
				f=a[i-1][l-1];
				g=(f >= 0.0 ? -Math.sqrt(h) : Math.sqrt(h));
				e[i-1]=scale*g;
				h -= f*g;
				a[i-1][l-1]=f-g;
				f=0.0;
				for (j=1;j<=l;j++) {
					a[j-1][i-1]=a[i-1][j-1]/h;
					g=0.0;
					for (k=1;k<=j;k++)
						g += a[j-1][k-1]*a[i-1][k-1];
					for (k=j+1;k<=l;k++)
						g += a[k-1][j-1]*a[i-1][k-1];
					e[j-1]=g/h;
					f += e[j-1]*a[i-1][j-1];
				}
				hh=f/(h+h);
				for (j=1;j<=l;j++) {
					f=a[i-1][j-1];
					e[j-1]=g=e[j-1]-hh*f;
					for (k=1;k<=j;k++)
						a[j-1][k-1] -= (f*e[k-1]+g*a[i-1][k-1]);
				}
			}
		} else
			e[i-1]=a[i-1][l-1];
		d[i-1]=h;
	}
	d[1-1]=0.0;
	e[1-1]=0.0;
	/* Contents of this loop can be omitted if eigenvectors not
			wanted except for statement d[i-1]=a[i-1][i-1]; */
	for (i=1;i<=n;i++) {
		l=i-1;
		if (d[i-1] != 0) {
			for (j=1;j<=l;j++) {
				g=0.0;
				for (k=1;k<=l;k++)
					g += a[i-1][k-1]*a[k-1][j-1];
				for (k=1;k<=l;k++)
					a[k-1][j-1] -= g*a[k-1][i-1];
			}
		}
		d[i-1]=a[i-1][i-1];
		a[i-1][i-1]=1.0;
		for (j=1;j<=l;j++) a[j-1][i-1]=a[i-1][j-1]=0.0;
	}
return  a;
}
	
public static double pythag(double a, double b)
{
//this method is used by QL method
	double absa,absb;
	absa=Math.abs(a);
	absb=Math.abs(b);
	if (absa > absb) return absa*Math.sqrt(1.0+(absb/absa)*(absb/absa));
	else return (absb==0.0 ? 0.0 : absb*Math.sqrt(1.0+(absa/absb)*(absa/absb)));
}

public static double[][] QL(double d[], double e[], double a[][])
{
// QL algorithm : eigenvalues of a symmetric matrix reduced to tridiagonal
// form by using method tridiagonal
    int n=d.length;
	int m,l,iter,i,j,k;
	double s,r,p,g,f,dd,c,b;
	for (i=2;i<=n;i++) e[i-2]=e[i-1];
	e[n-1]=0.0;
    double z[][]=new double[n][n];
    for(i=0;i<n;i++)
      for(j=0;j<n;j++)
        z[i][j]=a[i][j];
	for (l=1;l<=n;l++) {
		iter=0;
		do {
			for (m=l;m<=n-1;m++) {
				dd=Math.abs(d[m-1])+Math.abs(d[m]);
				if ((double)(Math.abs(e[m-1])+dd) == dd) break;
			}
			if (m != l) {
              if (iter++ == 30) System.out.println("Too many iterations in QL");
              g=(d[l]-d[l-1])/(2.0*e[l-1]);
              r=Matrix.pythag(g,1.0);
              g=d[m-1]-d[l-1]+e[l-1]/(g+Matrix.SIGN(r,g));
              s=c=1.0;
              p=0.0;
              for (i=m-1;i>=l;i--) {
					f=s*e[i-1];
					b=c*e[i-1];
					e[i]=(r=Matrix.pythag(f,g));
					if (r == 0.0) {
						d[i] -= p;
						e[m-1]=0.0;
						break;
					}
					s=f/r;
					c=g/r;
					g=d[i ]-p;
					r=(d[i-1]-g)*s+2.0*c*b;
					d[i ]=g+(p=s*r);
					g=c*r-b;
					for (k=1;k<=n;k++) {
						f=z[k-1][i ];
						z[k-1][i ]=s*z[k-1][i-1]+c*f;
						z[k-1][i-1]=c*z[k-1][i-1]-s*f;
					}
				}
				if (r == 0.0 && i >= l) continue;
				d[l-1] -= p;
				e[l-1]=g;
				e[m-1]=0.0;
			}
		} while (m != l);
	}
    return z;
}

public static double[][] eigenQL(double a[][])
{
 // QL algoritm to solve eigen value problems
 // symmetric matrices only (real eigen values)
 // first column of the matrix returns eigen values
 // second..n+1 column returns eigen vectors.
 // Note : If matrix is not symmetric DO NOT use
 // this method use eigenValue method (a QR algorithm)
 int i,j;
 int n=a.length;
 double sum[]=new double[n];;
 double d[]=new double[n];
 double b[][]=new double[n][n];
 double e[]=new double[n];
 double z[][]=new double[n+1][n];
 b=tridiagonal(a,d,e);
 b=QL(d,e,b);
 for(j=0;j<n;j++)
    {
    z[0][j]=d[j];
    for(i=0;i<n;i++)
       {
       z[i+1][j]=b[i][j]/b[0][j];
       if(z[i+1][j]<1e-13) z[i+1][i]=0;
       }
    }
 // result: first row is eigenvalues, the rest is eigenvector
 return z;
}
	
  public static void main(String args[]) throws IOException
  {
  //Eigen Value calculations
  double c[][];
  c=new double[3][3];
  double a[][]={{3.556,-1.778,0},{-1.778,3.556,-1.778},{0,-1.778,3.556}};
  //double a[][]={{-19,-4.7,6.45},{-4.7,4.6,11.8},  {6.45,11.8,-8.3}};
  //double a[][]={{3,-1,0},{-1,2,-1},{0,-1,3}};
  String s="Orijinal matris : \n";
  s+=Matrix.toString(a);
  c=new double[3][3];
  c=eigenQL(a);
  s+="sınır değer matrisi (birinci kolon sınır değerler) \n";
  s+=Matrix.toString(c);
  JOptionPane.showMessageDialog(null,s,
     "QL eigenvalue metodu : ",JOptionPane.PLAIN_MESSAGE);
  }
}
