Package org.geoforge.alg.jama

Source Code of org.geoforge.alg.jama.JamMatrix

package org.geoforge.alg.jama;

import java.text.NumberFormat;
import java.text.DecimalFormat;
import java.text.DecimalFormatSymbols;
import java.util.Locale;
import java.io.PrintWriter;
import java.io.BufferedReader;
import java.io.StreamTokenizer;
import org.geoforge.alg.jama.utils.JamMaths;

/**
Jama = Java JamMatrix class.
<P>
The Java JamMatrix Class provides the fundamental operations of numerical
linear algebra.  Various constructors create Matrices from two dimensional
arrays of double precision floating point numbers.  Various "gets" and
"sets" provide access to submatrices and matrix elements.  Several methods
implement basic matrix arithmetic, including matrix addition and
multiplication, matrix norms, and element-by-element array operations.
Methods for reading and printing matrices are also included.  All the
operations in this version of the JamMatrix Class involve real matrices.
Complex matrices may be handled in a future version.
<P>
Five fundamental matrix decompositions, which consist of pairs or triples
of matrices, permutation vectors, and the like, produce results in five
decomposition classes.  These decompositions are accessed by the JamMatrix
class to compute solutions of simultaneous linear equations, determinants,
inverses and other matrix functions.  The five decompositions are:
<P><UL>
<LI>Cholesky Decomposition of symmetric, positive definite matrices.
<LI>LU Decomposition of rectangular matrices.
<LI>QR Decomposition of rectangular matrices.
<LI>Singular Value Decomposition of rectangular matrices.
<LI>Eigenvalue Decomposition of both symmetric and nonsymmetric square matrices.
</UL>
<DL>
<DT><B>Example of use:</B></DT>
<P>
<DD>Solve a linear system A x = b and compute the residual norm, ||b - A x||.
<P><PRE>
double[][] vals = {{1.,2.,3},{4.,5.,6.},{7.,8.,10.}};
JamMatrix A = new JamMatrix(vals);
JamMatrix b = JamMatrix.random(3,1);
JamMatrix x = A.solve(b);
JamMatrix r = A.times(x).minus(b);
double rnorm = r.normInf();
</PRE></DD>
</DL>

@author The MathWorks, Inc. and the National Institute of Standards and Technology.
@version 5 August 1998
*/
public class JamMatrix implements Cloneable, java.io.Serializable
{

   /* ------------------------
   Class variables
    * ------------------------ */
   /** Array for internal storage of elements.
   @serial internal array storage.
    */
   private double[][] A;

   /** Row dimensions.
   @serial row dimension.
    */
   private int m;

   /** Column dimensions.
   @serial column dimension.
    */
   private int n;

   /* ------------------------
   Constructors
    * ------------------------ */
   /** Construct an m-by-n matrix of zeros.
   @param m    Number of rows.
   @param n    Number of colums.
    */
   public JamMatrix(int m, int n)
   {
      this.m = m;
      this.n = n;
      A = new double[m][n];
   }

   /** Construct an m-by-n constant matrix.
   @param m    Number of rows.
   @param n    Number of colums.
   @param s    Fill the matrix with this scalar value.
    */
   public JamMatrix(int m, int n, double s)
   {
      this.m = m;
      this.n = n;
      A = new double[m][n];
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = s;
         }
      }
   }

   /** Construct a matrix from a 2-D array.
   @param A    Two-dimensional array of doubles.
   @exception  IllegalArgumentException All rows must have the same length
   @see        #constructWithCopy
    */
   public JamMatrix(double[][] A)
   {
      m = A.length;
      n = A[0].length;
      for (int i = 0; i < m; i++)
      {
         if (A[i].length != n)
         {
            throw new IllegalArgumentException("All rows must have the same length.");
         }
      }
      this.A = A;
   }

   /** Construct a matrix quickly without checking arguments.
   @param A    Two-dimensional array of doubles.
   @param m    Number of rows.
   @param n    Number of colums.
    */
   public JamMatrix(double[][] A, int m, int n)
   {
      this.A = A;
      this.m = m;
      this.n = n;
   }

   /** Construct a matrix from a one-dimensional packed array
   @param vals One-dimensional array of doubles, packed by columns (ala Fortran).
   @param m    Number of rows.
   @exception  IllegalArgumentException Array length must be a multiple of m.
    */
   public JamMatrix(double vals[], int m)
   {
      this.m = m;
      n = (m != 0 ? vals.length / m : 0);
      if (m * n != vals.length)
      {
         throw new IllegalArgumentException("Array length must be a multiple of m.");
      }
      A = new double[m][n];
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = vals[i + j * m];
         }
      }
   }

   /* ------------------------
   Public Methods
    * ------------------------ */
   /** Construct a matrix from a copy of a 2-D array.
   @param A    Two-dimensional array of doubles.
   @exception  IllegalArgumentException All rows must have the same length
    */
   public static JamMatrix constructWithCopy(double[][] A)
   {
      int m = A.length;
      int n = A[0].length;
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         if (A[i].length != n)
         {
            throw new IllegalArgumentException("All rows must have the same length.");
         }
         for (int j = 0; j < n; j++)
         {
            C[i][j] = A[i][j];
         }
      }
      return X;
   }

   /** Make a deep copy of a matrix
    */
   public JamMatrix copy()
   {
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = A[i][j];
         }
      }
      return X;
   }

   /** Clone the JamMatrix object.
    */
   public Object clone()
   {
      return this.copy();
   }

   /** Access the internal two-dimensional array.
   @return     Pointer to the two-dimensional array of matrix elements.
    */
   public double[][] getArray()
   {
      return A;
   }

   /** Copy the internal two-dimensional array.
   @return     Two-dimensional array copy of matrix elements.
    */
   public double[][] getArrayCopy()
   {
      double[][] C = new double[m][n];
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = A[i][j];
         }
      }
      return C;
   }

   /** Make a one-dimensional column packed copy of the internal array.
   @return     JamMatrix elements packed in a one-dimensional array by columns.
    */
   public double[] getColumnPackedCopy()
   {
      double[] vals = new double[m * n];
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            vals[i + j * m] = A[i][j];
         }
      }
      return vals;
   }

   /** Make a one-dimensional row packed copy of the internal array.
   @return     JamMatrix elements packed in a one-dimensional array by rows.
    */
   public double[] getRowPackedCopy()
   {
      double[] vals = new double[m * n];
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            vals[i * n + j] = A[i][j];
         }
      }
      return vals;
   }

   /** Get row dimension.
   @return     m, the number of rows.
    */
   public int getRowDimension()
   {
      return m;
   }

   /** Get column dimension.
   @return     n, the number of columns.
    */
   public int getColumnDimension()
   {
      return n;
   }

   /** Get a single element.
   @param i    Row index.
   @param j    Column index.
   @return     A(i,j)
   @exception  ArrayIndexOutOfBoundsException
    */
   public double get(int i, int j)
   {
      return A[i][j];
   }

   /** Get a submatrix.
   @param i0   Initial row index
   @param i1   Final row index
   @param j0   Initial column index
   @param j1   Final column index
   @return     A(i0:i1,j0:j1)
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public JamMatrix getMatrix(int i0, int i1, int j0, int j1)
   {
      JamMatrix X = new JamMatrix(i1 - i0 + 1, j1 - j0 + 1);
      double[][] B = X.getArray();
      try
      {
         for (int i = i0; i <= i1; i++)
         {
            for (int j = j0; j <= j1; j++)
            {
               B[i - i0][j - j0] = A[i][j];
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
      return X;
   }

   /** Get a submatrix.
   @param r    Array of row indices.
   @param c    Array of column indices.
   @return     A(r(:),c(:))
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public JamMatrix getMatrix(int[] r, int[] c)
   {
      JamMatrix X = new JamMatrix(r.length, c.length);
      double[][] B = X.getArray();
      try
      {
         for (int i = 0; i < r.length; i++)
         {
            for (int j = 0; j < c.length; j++)
            {
               B[i][j] = A[r[i]][c[j]];
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
      return X;
   }

   /** Get a submatrix.
   @param i0   Initial row index
   @param i1   Final row index
   @param c    Array of column indices.
   @return     A(i0:i1,c(:))
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public JamMatrix getMatrix(int i0, int i1, int[] c)
   {
      JamMatrix X = new JamMatrix(i1 - i0 + 1, c.length);
      double[][] B = X.getArray();
      try
      {
         for (int i = i0; i <= i1; i++)
         {
            for (int j = 0; j < c.length; j++)
            {
               B[i - i0][j] = A[i][c[j]];
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
      return X;
   }

   /** Get a submatrix.
   @param r    Array of row indices.
   @param j0   Initial column index
   @param j1   Final column index
   @return     A(r(:),j0:j1)
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public JamMatrix getMatrix(int[] r, int j0, int j1)
   {
      JamMatrix X = new JamMatrix(r.length, j1 - j0 + 1);
      double[][] B = X.getArray();
      try
      {
         for (int i = 0; i < r.length; i++)
         {
            for (int j = j0; j <= j1; j++)
            {
               B[i][j - j0] = A[r[i]][j];
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
      return X;
   }

   /** Set a single element.
   @param i    Row index.
   @param j    Column index.
   @param s    A(i,j).
   @exception  ArrayIndexOutOfBoundsException
    */
   public void set(int i, int j, double s)
   {
      A[i][j] = s;
   }

   /** Set a submatrix.
   @param i0   Initial row index
   @param i1   Final row index
   @param j0   Initial column index
   @param j1   Final column index
   @param X    A(i0:i1,j0:j1)
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public void setMatrix(int i0, int i1, int j0, int j1, JamMatrix X)
   {
      try
      {
         for (int i = i0; i <= i1; i++)
         {
            for (int j = j0; j <= j1; j++)
            {
               A[i][j] = X.get(i - i0, j - j0);
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
   }

   /** Set a submatrix.
   @param r    Array of row indices.
   @param c    Array of column indices.
   @param X    A(r(:),c(:))
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public void setMatrix(int[] r, int[] c, JamMatrix X)
   {
      try
      {
         for (int i = 0; i < r.length; i++)
         {
            for (int j = 0; j < c.length; j++)
            {
               A[r[i]][c[j]] = X.get(i, j);
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
   }

   /** Set a submatrix.
   @param r    Array of row indices.
   @param j0   Initial column index
   @param j1   Final column index
   @param X    A(r(:),j0:j1)
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public void setMatrix(int[] r, int j0, int j1, JamMatrix X)
   {
      try
      {
         for (int i = 0; i < r.length; i++)
         {
            for (int j = j0; j <= j1; j++)
            {
               A[r[i]][j] = X.get(i, j - j0);
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
   }

   /** Set a submatrix.
   @param i0   Initial row index
   @param i1   Final row index
   @param c    Array of column indices.
   @param X    A(i0:i1,c(:))
   @exception  ArrayIndexOutOfBoundsException Submatrix indices
    */
   public void setMatrix(int i0, int i1, int[] c, JamMatrix X)
   {
      try
      {
         for (int i = i0; i <= i1; i++)
         {
            for (int j = 0; j < c.length; j++)
            {
               A[i][c[j]] = X.get(i - i0, j);
            }
         }
      }
      catch (ArrayIndexOutOfBoundsException e)
      {
         throw new ArrayIndexOutOfBoundsException("Submatrix indices");
      }
   }

   /** JamMatrix transpose.
   @return    A'
    */
   public JamMatrix transpose()
   {
      JamMatrix X = new JamMatrix(n, m);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[j][i] = A[i][j];
         }
      }
      return X;
   }

   /** One norm
   @return    maximum column sum.
    */
   public double norm1()
   {
      double f = 0;
      for (int j = 0; j < n; j++)
      {
         double s = 0;
         for (int i = 0; i < m; i++)
         {
            s += Math.abs(A[i][j]);
         }
         f = Math.max(f, s);
      }
      return f;
   }

   /** Two norm
   @return    maximum singular value.
    */
   public double norm2()
   {
      return (new JamSingularValueDecomposition(this).norm2());
   }

   /** Infinity norm
   @return    maximum row sum.
    */
   public double normInf()
   {
      double f = 0;
      for (int i = 0; i < m; i++)
      {
         double s = 0;
         for (int j = 0; j < n; j++)
         {
            s += Math.abs(A[i][j]);
         }
         f = Math.max(f, s);
      }
      return f;
   }

   /** Frobenius norm
   @return    sqrt of sum of squares of all elements.
    */
   public double normF()
   {
      double f = 0;
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            f = JamMaths.s_getHypot(f, A[i][j]);
         }
      }
      return f;
   }

   /**  Unary minus
   @return    -A
    */
   public JamMatrix uminus()
   {
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = -A[i][j];
         }
      }
      return X;
   }

   /** C = A + B
   @param B    another matrix
   @return     A + B
    */
   public JamMatrix plus(JamMatrix B)
   {
      checkMatrixDimensions(B);
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = A[i][j] + B.A[i][j];
         }
      }
      return X;
   }

   /** A = A + B
   @param B    another matrix
   @return     A + B
    */
   public JamMatrix plusEquals(JamMatrix B)
   {
      checkMatrixDimensions(B);
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = A[i][j] + B.A[i][j];
         }
      }
      return this;
   }

   /** C = A - B
   @param B    another matrix
   @return     A - B
    */
   public JamMatrix minus(JamMatrix B)
   {
      checkMatrixDimensions(B);
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = A[i][j] - B.A[i][j];
         }
      }
      return X;
   }

   /** A = A - B
   @param B    another matrix
   @return     A - B
    */
   public JamMatrix minusEquals(JamMatrix B)
   {
      checkMatrixDimensions(B);
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = A[i][j] - B.A[i][j];
         }
      }
      return this;
   }

   /** Element-by-element multiplication, C = A.*B
   @param B    another matrix
   @return     A.*B
    */
   public JamMatrix arrayTimes(JamMatrix B)
   {
      checkMatrixDimensions(B);
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = A[i][j] * B.A[i][j];
         }
      }
      return X;
   }

   /** Element-by-element multiplication in place, A = A.*B
   @param B    another matrix
   @return     A.*B
    */
   public JamMatrix arrayTimesEquals(JamMatrix B)
   {
      checkMatrixDimensions(B);
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = A[i][j] * B.A[i][j];
         }
      }
      return this;
   }

   /** Element-by-element right division, C = A./B
   @param B    another matrix
   @return     A./B
    */
   public JamMatrix arrayRightDivide(JamMatrix B)
   {
      checkMatrixDimensions(B);
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = A[i][j] / B.A[i][j];
         }
      }
      return X;
   }

   /** Element-by-element right division in place, A = A./B
   @param B    another matrix
   @return     A./B
    */
   public JamMatrix arrayRightDivideEquals(JamMatrix B)
   {
      checkMatrixDimensions(B);
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = A[i][j] / B.A[i][j];
         }
      }
      return this;
   }

   /** Element-by-element left division, C = A.\B
   @param B    another matrix
   @return     A.\B
    */
   public JamMatrix arrayLeftDivide(JamMatrix B)
   {
      checkMatrixDimensions(B);
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = B.A[i][j] / A[i][j];
         }
      }
      return X;
   }

   /** Element-by-element left division in place, A = A.\B
   @param B    another matrix
   @return     A.\B
    */
   public JamMatrix arrayLeftDivideEquals(JamMatrix B)
   {
      checkMatrixDimensions(B);
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = B.A[i][j] / A[i][j];
         }
      }
      return this;
   }

   /** Multiply a matrix by a scalar, C = s*A
   @param s    scalar
   @return     s*A
    */
   public JamMatrix times(double s)
   {
      JamMatrix X = new JamMatrix(m, n);
      double[][] C = X.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            C[i][j] = s * A[i][j];
         }
      }
      return X;
   }

   /** Multiply a matrix by a scalar in place, A = s*A
   @param s    scalar
   @return     replace A by s*A
    */
   public JamMatrix timesEquals(double s)
   {
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            A[i][j] = s * A[i][j];
         }
      }
      return this;
   }

   /** Linear algebraic matrix multiplication, A * B
   @param B    another matrix
   @return     JamMatrix product, A * B
   @exception  IllegalArgumentException JamMatrix inner dimensions must agree.
    */
   public JamMatrix times(JamMatrix B)
   {
      if (B.m != n)
      {
         throw new IllegalArgumentException("Matrix inner dimensions must agree.");
      }
      JamMatrix X = new JamMatrix(m, B.n);
      double[][] C = X.getArray();
      double[] Bcolj = new double[n];
      for (int j = 0; j < B.n; j++)
      {
         for (int k = 0; k < n; k++)
         {
            Bcolj[k] = B.A[k][j];
         }
         for (int i = 0; i < m; i++)
         {
            double[] Arowi = A[i];
            double s = 0;
            for (int k = 0; k < n; k++)
            {
               s += Arowi[k] * Bcolj[k];
            }
            C[i][j] = s;
         }
      }
      return X;
   }

   /** LU Decomposition
   @return     JamLUDecomposition
   @see JamLUDecomposition
    */
   public JamLUDecomposition lu()
   {
      return new JamLUDecomposition(this);
   }

   /** QR Decomposition
   @return     JamQRDecomposition
   @see JamQRDecomposition
    */
   public JamQRDecomposition qr()
   {
      return new JamQRDecomposition(this);
   }

   /** Cholesky Decomposition
   @return     JamCholeskyDecomposition
   @see JamCholeskyDecomposition
    */
   public JamCholeskyDecomposition chol()
   {
      return new JamCholeskyDecomposition(this);
   }

   /** Singular Value Decomposition
   @return     JamSingularValueDecomposition
   @see JamSingularValueDecomposition
    */
   public JamSingularValueDecomposition svd()
   {
      return new JamSingularValueDecomposition(this);
   }

   /** Eigenvalue Decomposition
   @return     JamEigenvalueDecomposition
   @see JamEigenvalueDecomposition
    */
   public JamEigenvalueDecomposition eig()
   {
      return new JamEigenvalueDecomposition(this);
   }

   /** Solve A*X = B
   @param B    right hand side
   @return     solution if A is square, least squares solution otherwise
    */
   public JamMatrix solve(JamMatrix B)
   {
      return (m == n ? (new JamLUDecomposition(this)).solve(B)
            : (new JamQRDecomposition(this)).solve(B));
   }

   /** Solve X*A = B, which is also A'*X' = B'
   @param B    right hand side
   @return     solution if A is square, least squares solution otherwise.
    */
   public JamMatrix solveTranspose(JamMatrix B)
   {
      return transpose().solve(B.transpose());
   }

   /** JamMatrix inverse or pseudoinverse
   @return     inverse(A) if A is square, pseudoinverse otherwise.
    */
   public JamMatrix inverse()
   {
      return solve(identity(m, m));
   }

   /** JamMatrix determinant
   @return     determinant
    */
   public double det()
   {
      return new JamLUDecomposition(this).det();
   }

   /** JamMatrix rank
   @return     effective numerical rank, obtained from SVD.
    */
   public int rank()
   {
      return new JamSingularValueDecomposition(this).rank();
   }

   /** JamMatrix condition (2 norm)
   @return     ratio of largest to smallest singular value.
    */
   public double cond()
   {
      return new JamSingularValueDecomposition(this).cond();
   }

   /** JamMatrix trace.
   @return     sum of the diagonal elements.
    */
   public double trace()
   {
      double t = 0;
      for (int i = 0; i < Math.min(m, n); i++)
      {
         t += A[i][i];
      }
      return t;
   }

   /** Generate matrix with random elements
   @param m    Number of rows.
   @param n    Number of colums.
   @return     An m-by-n matrix with uniformly distributed random elements.
    */
   public static JamMatrix random(int m, int n)
   {
      JamMatrix A = new JamMatrix(m, n);
      double[][] X = A.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            X[i][j] = Math.random();
         }
      }
      return A;
   }

   /** Generate identity matrix
   @param m    Number of rows.
   @param n    Number of colums.
   @return     An m-by-n matrix with ones on the diagonal and zeros elsewhere.
    */
   public static JamMatrix identity(int m, int n)
   {
      JamMatrix A = new JamMatrix(m, n);
      double[][] X = A.getArray();
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            X[i][j] = (i == j ? 1.0 : 0.0);
         }
      }
      return A;
   }

   /** Print the matrix to stdout.   Line the elements up in columns
    * with a Fortran-like 'Fw.d' style format.
   @param w    Column width.
   @param d    Number of digits after the decimal.
    */
   public void print(int w, int d)
   {
      print(new PrintWriter(System.out, true), w, d);
   }

   /** Print the matrix to the output stream.   Line the elements up in
    * columns with a Fortran-like 'Fw.d' style format.
   @param output Output stream.
   @param w      Column width.
   @param d      Number of digits after the decimal.
    */
   public void print(PrintWriter output, int w, int d)
   {
      DecimalFormat format = new DecimalFormat();
      format.setDecimalFormatSymbols(new DecimalFormatSymbols(Locale.US));
      format.setMinimumIntegerDigits(1);
      format.setMaximumFractionDigits(d);
      format.setMinimumFractionDigits(d);
      format.setGroupingUsed(false);
      print(output, format, w + 2);
   }

   /** Print the matrix to stdout.  Line the elements up in columns.
    * Use the format object, and right justify within columns of width
    * characters.
    * Note that is the matrix is to be read back in, you probably will want
    * to use a NumberFormat that is set to US Locale.
   @param format A  Formatting object for individual elements.
   @param width     Field width for each column.
   @see java.text.DecimalFormat#setDecimalFormatSymbols
    */
   public void print(NumberFormat format, int width)
   {
      print(new PrintWriter(System.out, true), format, width);
   }

   // DecimalFormat is a little disappointing coming from Fortran or C's printf.
   // Since it doesn't pad on the left, the elements will come out different
   // widths.  Consequently, we'll pass the desired column width in as an
   // argument and do the extra padding ourselves.
   /** Print the matrix to the output stream.  Line the elements up in columns.
    * Use the format object, and right justify within columns of width
    * characters.
    * Note that is the matrix is to be read back in, you probably will want
    * to use a NumberFormat that is set to US Locale.
   @param output the output stream.
   @param format A formatting object to format the matrix elements
   @param width  Column width.
   @see java.text.DecimalFormat#setDecimalFormatSymbols
    */
   public void print(PrintWriter output, NumberFormat format, int width)
   {
      output.println()// start on new line.
      for (int i = 0; i < m; i++)
      {
         for (int j = 0; j < n; j++)
         {
            String s = format.format(A[i][j]); // format the number
            int padding = Math.max(1, width - s.length()); // At _least_ 1 space
            for (int k = 0; k < padding; k++)
               output.print(' ');
            output.print(s);
         }
         output.println();
      }
      output.println();   // end with blank line.
   }

   /** Read a matrix from a stream.  The format is the same the print method,
    * so printed matrices can be read back in (provided they were printed using
    * US Locale).  Elements are separated by
    * whitespace, all the elements for each row appear on a single line,
    * the last row is followed by a blank line.
   @param input the input stream.
    */
   public static JamMatrix read(BufferedReader input) throws java.io.IOException
   {
      StreamTokenizer tokenizer = new StreamTokenizer(input);

      // Although StreamTokenizer will parse numbers, it doesn't recognize
      // scientific notation (E or D); however, Double.valueOf does.
      // The strategy here is to disable StreamTokenizer's number parsing.
      // We'll only get whitespace delimited words, EOL's and EOF's.
      // These words should all be numbers, for Double.valueOf to parse.

      tokenizer.resetSyntax();
      tokenizer.wordChars(0, 255);
      tokenizer.whitespaceChars(0, ' ');
      tokenizer.eolIsSignificant(true);
      // TODO: replace by eg:
     
      java.util.Vector<Double> vD = new java.util.Vector<Double>();

      // Ignore initial empty lines
      while (tokenizer.nextToken() == StreamTokenizer.TT_EOL);
      if (tokenizer.ttype == StreamTokenizer.TT_EOF)
         throw new java.io.IOException("Unexpected EOF on matrix read.");
      do
      {
         vD.addElement(Double.valueOf(tokenizer.sval)); // Read & store 1st row.
      }
      while (tokenizer.nextToken() == StreamTokenizer.TT_WORD);

      int n = vD.size()// Now we've got the number of columns!
      double row[] = new double[n];
      for (int j = 0; j < n; j++// extract the elements of the 1st row.
         row[j] = vD.elementAt(j).doubleValue();
     
      // TODO: replace obsolete Vector by ArrayList
      java.util.Vector<double[]> v = new java.util.Vector<double[]>();
     
     
      v.addElement(row)// Start storing rows instead of columns.
      while (tokenizer.nextToken() == StreamTokenizer.TT_WORD)
      {
         // While non-empty lines
         v.addElement(row = new double[n]);
         int j = 0;
         do
         {
            if (j >= n)
               throw new java.io.IOException("Row " + v.size() + " is too long.");
            row[j++] = Double.valueOf(tokenizer.sval).doubleValue();
         }
         while (tokenizer.nextToken() == StreamTokenizer.TT_WORD);
         if (j < n)
            throw new java.io.IOException("Row " + v.size() + " is too short.");
      }
      int m = v.size()// Now we've got the number of rows.
      double[][] A = new double[m][];
      v.copyInto(A)// copy the rows out of the vector
      return new JamMatrix(A);
   }


   /* ------------------------
   Private Methods
    * ------------------------ */
   /** Check if size(A) == size(B) **/
   private void checkMatrixDimensions(JamMatrix B)
   {
      if (B.m != m || B.n != n)
      {
         throw new IllegalArgumentException("Matrix dimensions must agree.");
      }
   }

   @Override
   public String toString()
   {
      String str = "";
     
      for (int i = 0; i < this.m; i++)
      {
         str+= "[";
         for (int j = 0; j < this.n; j++)
         {
            str+= this.get(i, j);
            str+= " ";
         }
         str+= "]\n";
      }
     
      return str;
   }
   private static final long serialVersionUID = 1;

}
TOP

Related Classes of org.geoforge.alg.jama.JamMatrix

TOP
Copyright © 2018 www.massapi.com. All rights reserved.
All source code are property of their respective owners. Java is a trademark of Sun Microsystems, Inc and owned by ORACLE Inc. Contact coftware#gmail.com.