JSci.maths.matrices.IntegerTridiagonalMatrix Maven / Gradle / Ivy
Show all versions of jsci Show documentation
/* AUTO-GENERATED */
package JSci.maths.matrices;
import JSci.maths.ArrayMath;
import JSci.maths.ExtraMath;
import JSci.maths.Mapping;
import JSci.maths.LinearMath;
import JSci.maths.MaximumIterationsExceededException;
import JSci.maths.DimensionException;
import JSci.maths.vectors.AbstractIntegerVector;
import JSci.maths.vectors.IntegerVector;
import JSci.maths.groups.AbelianGroup;
import JSci.maths.algebras.*;
import JSci.maths.fields.*;
/**
* The IntegerTridiagonalMatrix class provides an object for encapsulating integer tridiagonal matrices.
* @version 2.2
* @author Mark Hale
*/
public class IntegerTridiagonalMatrix extends AbstractIntegerSquareMatrix implements TridiagonalMatrix {
/**
* Tridiagonal data.
*/
protected final int ldiag[];
protected final int diag[];
protected final int udiag[];
/**
* Constructs an empty matrix.
* @param size the number of rows/columns
*/
public IntegerTridiagonalMatrix(final int size) {
super(size);
ldiag=new int[size];
diag=new int[size];
udiag=new int[size];
}
/**
* Constructs a matrix from an array.
* Any non-tridiagonal elements in the array are ignored.
* @param array an assigned value
* @exception MatrixDimensionException If the array is not square.
*/
public IntegerTridiagonalMatrix(final int array[][]) {
this(array.length);
if(!ArrayMath.isSquare(array))
throw new MatrixDimensionException("Array is not square.");
diag[0]=array[0][0];
udiag[0]=array[0][1];
int i=1;
for(;i=0 && i=0 && j=0 && i=0 && j-norm.
* @author Taber Smith
*/
public int infNorm() {
int result=Math.abs(diag[0])+Math.abs(udiag[0]);
int tmpResult;
int i=1;
for(;iresult)
result=tmpResult;
}
tmpResult=Math.abs(ldiag[i])+Math.abs(diag[i]);
if(tmpResult>result)
result=tmpResult;
return result;
}
/**
* Returns the Frobenius (l2) norm.
* @author Taber Smith
*/
public double frobeniusNorm() {
double result=diag[0]*diag[0]+udiag[0]*udiag[0];
int i=1;
for(;i3) {
array[1][0]=ldiag[1]*m.diag[0]+diag[1]*m.ldiag[1];
array[1][1]=ldiag[1]*m.udiag[0]+diag[1]*m.diag[1]+udiag[1]*m.ldiag[2];
array[1][2]=diag[1]*m.udiag[1]+udiag[1]*m.diag[2];
array[1][3]=udiag[1]*m.udiag[2];
}
if(mRow==3) {
array[1][0]=ldiag[1]*m.diag[0]+diag[1]*m.ldiag[1];
array[1][1]=ldiag[1]*m.udiag[0]+diag[1]*m.diag[1]+udiag[1]*m.ldiag[2];
array[1][2]=diag[1]*m.udiag[1]+udiag[1]*m.diag[2];
} else if(mRow>4) {
for(int i=2;i3) {
array[mRow-2][mRow-4]=ldiag[mRow-2]*m.ldiag[mRow-3];
array[mRow-2][mRow-3]=ldiag[mRow-2]*m.diag[mRow-3]+diag[mRow-2]*m.ldiag[mRow-2];
array[mRow-2][mRow-2]=ldiag[mRow-2]*m.udiag[mRow-3]+diag[mRow-2]*m.diag[mRow-2]+udiag[mRow-2]*m.ldiag[mRow-1];
array[mRow-2][mRow-1]=diag[mRow-2]*m.udiag[mRow-2]+udiag[mRow-2]*m.diag[mRow-1];
}
mRow--;
array[mRow][mRow-2]=ldiag[mRow]*m.ldiag[mRow-1];
array[mRow][mRow-1]=ldiag[mRow]*m.diag[mRow-1]+diag[mRow]*m.ldiag[mRow];
array[mRow][mRow]=ldiag[mRow]*m.udiag[mRow-1]+diag[mRow]*m.diag[mRow];
return new IntegerSquareMatrix(array);
} else
throw new MatrixDimensionException("Incompatible matrices.");
}
private IntegerSquareMatrix multiplyTridiagonal(final AbstractIntegerSquareMatrix m) {
int mRow=numRows;
if(numCols==m.rows()) {
final int array[][]=new int[mRow][mRow];
array[0][0]=diag[0]*m.getElement(0,0)+udiag[0]*m.getElement(1,0);
array[0][1]=diag[0]*m.getElement(0,1)+udiag[0]*m.getElement(1,1);
array[0][2]=udiag[0]*m.getElement(1,2);
if(mRow>3) {
array[1][0]=ldiag[1]*m.getElement(0,0)+diag[1]*m.getElement(1,0);
array[1][1]=ldiag[1]*m.getElement(0,1)+diag[1]*m.getElement(1,1)+udiag[1]*m.getElement(2,1);
array[1][2]=diag[1]*m.getElement(1,2)+udiag[1]*m.getElement(2,2);
array[1][3]=udiag[1]*m.getElement(2,3);
}
if(mRow==3) {
array[1][0]=ldiag[1]*m.getElement(0,0)+diag[1]*m.getElement(1,0);
array[1][1]=ldiag[1]*m.getElement(0,1)+diag[1]*m.getElement(1,1)+udiag[1]*m.getElement(2,1);
array[1][2]=diag[1]*m.getElement(1,2)+udiag[1]*m.getElement(2,2);
} else if(mRow>4) {
for(int i=2;i3) {
array[mRow-2][mRow-4]=ldiag[mRow-2]*m.getElement(mRow-3,mRow-4);
array[mRow-2][mRow-3]=ldiag[mRow-2]*m.getElement(mRow-3,mRow-3)+diag[mRow-2]*m.getElement(mRow-2,mRow-3);
array[mRow-2][mRow-2]=ldiag[mRow-2]*m.getElement(mRow-3,mRow-2)+diag[mRow-2]*m.getElement(mRow-2,mRow-2)+udiag[mRow-2]*m.getElement(mRow-1,mRow-2);
array[mRow-2][mRow-1]=diag[mRow-2]*m.getElement(mRow-2,mRow-1)+udiag[mRow-2]*m.getElement(mRow-1,mRow-1);
}
mRow--;
array[mRow][mRow-2]=ldiag[mRow]*m.getElement(mRow-1,mRow-2);
array[mRow][mRow-1]=ldiag[mRow]*m.getElement(mRow-1,mRow-1)+diag[mRow]*m.getElement(mRow,mRow-1);
array[mRow][mRow]=ldiag[mRow]*m.getElement(mRow-1,mRow)+diag[mRow]*m.getElement(mRow,mRow);
return new IntegerSquareMatrix(array);
} else {
throw new MatrixDimensionException("Incompatible matrices.");
}
}
// TRANSPOSE
/**
* Returns the transpose of this matrix.
* @return a int matrix
*/
public Matrix transpose() {
final IntegerTridiagonalMatrix ans=new IntegerTridiagonalMatrix(numRows);
System.arraycopy(ldiag,1,ans.udiag,0,ldiag.length-1);
System.arraycopy(diag,0,ans.diag,0,diag.length);
System.arraycopy(udiag,0,ans.ldiag,1,udiag.length-1);
return ans;
}
// CHOLESKY DECOMPOSITION
/**
* Returns the Cholesky decomposition of this matrix.
* Matrix must be symmetric and positive definite.
* @return an array with [0] containing the L-matrix and [1] containing the U-matrix.
*/
public AbstractDoubleSquareMatrix[] choleskyDecompose() {
final int N=numRows;
final double arrayL[][]=new double[N][N];
final double arrayU[][]=new double[N][N];
double tmp=Math.sqrt(diag[0]);
arrayL[0][0]=arrayU[0][0]=tmp;
arrayL[1][0]=arrayU[0][1]=ldiag[1]/tmp;
for(int j=1;jJAMA (public domain).
* @return an array with [0] containing the Q-matrix and [1] containing the R-matrix.
* @jsci.planetmath QRDecomposition
*/
public AbstractDoubleSquareMatrix[] qrDecompose() {
final int N=numRows;
final double array[][]=new double[N][N];
final double arrayQ[][]=new double[N][N];
final double arrayR[][]=new double[N][N];
// copy matrix
array[0][0]=diag[0];
array[0][1]=udiag[0];
for(int i=1;i=0; k--) {
arrayQ[k][k] = 1.0;
for(int j=k; jJAMA (public domain).
* @return an array with [0] containing the U-matrix, [1] containing the S-matrix and [2] containing the V-matrix.
*/
public AbstractDoubleSquareMatrix[] singularValueDecompose() {
final int N=numRows;
final int Nm1=N-1;
final double array[][]=new double[N][N];
final double arrayU[][]=new double[N][N];
final double arrayS[]=new double[N];
final double arrayV[][]=new double[N][N];
final double e[]=new double[N];
final double work[]=new double[N];
// copy matrix
array[0][0]=diag[0];
array[0][1]=udiag[0];
for(int i=1;i=0;k--) {
if(arrayS[k]!=0.0) {
for(int j=k+1;j=0;k--) {
if(k0) {
int k, action;
// action = 1 if arrayS[p] and e[k-1] are negligible and k=-1;k--) {
if(k==-1)
break;
if(Math.abs(e[k])<=eps*(Math.abs(arrayS[k])+Math.abs(arrayS[k+1]))) {
e[k]=0.0;
break;
}
}
if(k==p-2) {
action=4;
} else {
int ks;
for(ks=p-1;ks>=k;ks--) {
if(ks==k)
break;
double t=(ks!=p ? Math.abs(e[ks]) : 0.0)+(ks!=k+1 ? Math.abs(e[ks-1]) : 0.0);
if(Math.abs(arrayS[ks])<=eps*t) {
arrayS[ks]=0.0;
break;
}
}
if(ks==k) {
action=3;
} else if(ks==p-1) {
action=1;
} else {
action=2;
k=ks;
}
}
k++;
switch(action) {
// deflate negligible arrayS[p]
case 1: {
double f=e[p-2];
e[p-2]=0.0;
for(int j=p-2;j>=k;j--) {
double t=ExtraMath.hypot(arrayS[j],f);
final double cs=arrayS[j]/t;
final double sn=f/t;
arrayS[j]=t;
if(j!=k) {
f=-sn*e[j-1];
e[j-1]*=cs;
}
for(int i=0;i=arrayS[k+1])
break;
double tmp=arrayS[k];
arrayS[k]=arrayS[k+1];
arrayS[k+1]=tmp;
if(k