for (j=0;j<i;j++) {
a[i][j]=myran.doub();
a[j][i]=a[i][j];
}
}
Unsymmeig usym = new Unsymmeig(a,true,false);
// for (i=0;i<N;i++) System.out.printf(usym.wri[i]);
// Test that all eigenvalues are real for symmetric matrix
for (i=0;i<N;i++)
localflag = localflag || (usym.wri[i].im() != 0);
globalflag = globalflag || localflag;
if (localflag) {
fail("*** Unsymmeig, symmetric, interface1: Symmetric matrix gave an eigenvalue that was not real");
}
// Test eigenvector/eigenvalue pairs
for (i=0;i<N;i++) { // for each eigenvector
for (j=0;j<N;j++) vec[j]=usym.zz[j][i];
res=matmul(a,vec);
for (j=0;j<N;j++) vec[j] *= usym.wri[i].re();
// System.out.printf(maxel(vecsub(res,vec)));
localflag = localflag || (maxel(vecsub(res,vec)) > sbeps);
}
globalflag = globalflag || localflag;
if (localflag) {
fail("*** Unsymmeig, symmetric, interface1: Matrix times eigenvector was not the same as lambda*eigenvector");
}
// Test the sorting of the eigenvalues
for (i=1;i<N;i++)
localflag = localflag || (usym.wri[i].re() > usym.wri[i-1].re());
globalflag = globalflag || localflag;
if (localflag) {
fail("*** Unsymmeig, symmetric, interface1: Eigenvalues not sorted in high-to-low order");
}
//------------------------------------------------------------------------
// Test Unsymmeig, non-symmetric, interface1
System.out.println("Testing Unsymmeig, non-symmetric, interface1");
ranmat(a);
Unsymmeig usym2 = new Unsymmeig(a,true,false);
// for (i=0;i<N;i++) System.out.printf(usym2.wri[i]);
// Test eigenvector/eigenvalue pairs
for (i=0;i<N;i++) { // for each eigenvector
if (usym2.wri[i].im() == 0.0) {
for (j=0;j<N;j++) vec[j]=usym2.zz[j][i];
res=matmul(a,vec);
for (j=0;j<N;j++) vec[j] *= usym2.wri[i].re();
// System.out.println("real eigenvalue " << maxel(vecsub(res,vec)));
localflag = localflag || (maxel(vecsub(res,vec)) > sbeps);
} else {
if (usym2.wri[i].im() > 0.0)
for (j=0;j<N;j++) zvec[j]=new Complex(usym2.zz[j][i],usym2.zz[j][i+1]);
else
for (j=0;j<N;j++) zvec[j]=new Complex(usym2.zz[j][i-1],-usym2.zz[j][i]);
for (j=0;j<N;j++) {
zres[j]=new Complex(0.,0.);
for (k=0;k<N;k++) zres[j] = zres[j].add(zvec[k].mul(a[j][k]));
}
for (j=0;j<N;j++) zvec[j] = zvec[j].mul(usym2.wri[i]);
max=0.;
for (j=0;j<N;j++) max = (max > zres[j].sub(zvec[j]).abs() ? max : zres[j].sub(zvec[j]).abs());
// System.out.println("imag eigenvalue " << max);
localflag = localflag || (max > sbeps);
}
max=0.;
for (j=0;j<N;j++) max = (max > zres[j].sub(zvec[j]).abs() ? max : zres[j].sub(zvec[j]).abs());
localflag = localflag || (max > sbeps);
}
globalflag = globalflag || localflag;
if (localflag) {
fail("*** Unsymmeig, non-symmetric, interface1: Matrix times eigenvector was not the same as lambda*eigenvector");
}
// Test the sorting of the eigenvalues
for (i=1;i<N;i++)
localflag = localflag || (usym2.wri[i].re() > usym2.wri[i-1].re());
globalflag = globalflag || localflag;
if (localflag) {
fail("*** Unsymmeig, non-symmetric, interface1: Eigenvalues not sorted in high-to-low order of real part");
}
for (i=1;i<N;i++)
localflag = localflag || (usym2.wri[i].im() > 0.0) && (!usym2.wri[i].equals(usym2.wri[i+1].conj()));
globalflag = globalflag || localflag;
if (localflag) {
fail("*** Unsymmeig, non-symmetric, interface1: A complex eigenvalue with positive imag part is not followed by its conjugate");
}
//------------------------------------------------------------------------
// Test Unsymmeig, non-symmetric, interface2
System.out.println("Testing Unsymmeig, non-symmetric, interface2");
for (i=0;i<N;i++)
for (j=0;j<N;j++)
a[i][j]= (i > j+1 ? 0.0 : myran.doub());
Unsymmeig usym3 = new Unsymmeig(a,true,true);
// for (i=0;i<N;i++) System.out.printf(usym3.wri[i]);
// Test eigenvector/eigenvalue pairs
for (i=0;i<N;i++) { // for each eigenvector