/* AUTO-GENERATED */
package jsci.maths.matrices;

import jsci.GlobalSettings;
import jsci.maths.ArrayMath;
import jsci.maths.ExtraMath;
import jsci.maths.LinearMath;
import jsci.maths.Mapping;
import jsci.maths.MaximumIterationsExceededException;
import jsci.maths.algebras.*;
import jsci.maths.fields.*;
import jsci.maths.groups.AbelianGroup;
import jsci.maths.vectors.AbstractIntegerVector;

/**
* The AbstractIntegerSquareMatrix class provides an object for encapsulating integer square matrices.
* @version 2.2
* @author Mark Hale
*/
public abstract class AbstractIntegerSquareMatrix extends AbstractIntegerMatrix implements SquareMatrix {
        protected transient AbstractDoubleSquareMatrix[] LU;
        protected transient int[] LUpivot;
        /**
        * Constructs a matrix.
        */
        protected AbstractIntegerSquareMatrix(final int size) {
                super(size, size);
        }
        /**
        * Converts this matrix to a double matrix.
        * @return a double square matrix
        */
        public AbstractDoubleMatrix toDoubleMatrix() {
                final double ans[][]=new double[numRows][numCols];
                for(int i=0;irows()+1
        * to hold the pivot information (null if not interested).
        * The last array element will contain the parity.
        * @return an array with [0] containing the L-matrix
        * and [1] containing the U-matrix.
        * @jsci.planetmath LUDecomposition
        */
        public AbstractDoubleSquareMatrix[] luDecompose(int pivot[]) {
                if(LU!=null) {
                        if(pivot!=null)
                                System.arraycopy(LUpivot,0,pivot,0,pivot.length);
                        return LU;
                }
                int pivotrow;
                final int N=numRows;
                final double arrayL[][]=new double[N][N];
                final double arrayU[][]=new double[N][N];
                if(pivot==null)
                        pivot=new int[N+1];
                for(int i=0;imax) {
                                        max=tmp;
                                        pivotrow=i;
                                }
                        }
                // swap row j with pivotrow
                        if(pivotrow!=j) {
                                double[] tmprow = arrayU[j];
                                arrayU[j] = arrayU[pivotrow];
                                arrayU[pivotrow] = tmprow;
                                int k=pivot[j];
                                pivot[j]=pivot[pivotrow];
                                pivot[pivotrow]=k;
                                // update parity
                                pivot[N]=-pivot[N];
                        }
                // divide by pivot
                        double tmp=arrayU[j][j];
                        for(int i=j+1;iJAMA (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
                for(int i=0;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.
        * @jsci.planetmath SingularValueDecomposition
        */
        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
                for(int i=0;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