Class LUDecomp

java.lang.Object
org.bzdev.math.LUDecomp
All Implemented Interfaces:
TriangularDecomp

public class LUDecomp extends Object implements TriangularDecomp
LU Decomposition class. This class implements LU decomposition, which for a matrix A and permutation matrix P, factors PA so that PA = LU, where L is a lower triangular matrix and U is an upper triangular matrix.

The implementation is based on the one used in the Jama package (a public-domain package provided via a collaboration between MathWorks and the National Institute of Standards and Technology). Jama defines a matrix class whereas LUDecomp treats matrices as two-dimensional arrays, and was written so that the org.bzdev library is self-contained: a few classes in that library have to solve systems of linear equations. It is intended for incidental use where performance is important but where a comprehensive linear algebra package is overkill.

LUDecomp uses the class Permutation to represent permutations, so the use of the permutation matrix P is implicit. One can also provide the matrix used to store the triangular matrices, which avoids additional memory allocation when multiple sets of equations have to be solved sequentially, and store the triangular matrices in the same arrays used to store the matrix being decomposed, which saves space when if the matrix being decomposed will no longer be needed.

Note: this is a correct implementation of the algorithm. Sample implementations shown at various web sites sometimes miss important details. For example, as of late 2013, the implementation shown at rosettacode.org had (and maybe still has) an implementation that could fail in some circumstances. It computes the permutation matrix so as to maximize the value of elements of the original matrix on the diagonal, which does not necessarily prevent the u_{jj} terms from vanishing. The Jama implementation, by contrast, maximizes the absolute value of u_{jj} when computing the permutation, which works in all cases.

  • Constructor Summary

    Constructors
    Constructor
    Description
    LUDecomp(double[][] A)
    Constructor.
    LUDecomp(double[][] A, double[][] LU)
    Constructor providing an LU matrix.
    LUDecomp(double[][] A, int m, int n, boolean strict)
    Constructor specifying matrix dimensions.
    LUDecomp(double[] flatMatrix, boolean columnMajorOrder, int m, int n)
    Constructor given a flat representation of a matrix.
  • Method Summary

    Modifier and Type
    Method
    Description
    double
    det()
    Get the determinate of the matrix to which this decomposition applies.
    double[][]
    Get the inverse of a matrix.
    void
    getInverse(double[][] result)
    Get the inverse of a matrix and store the inverse.
    void
    getInverse(double[] result, boolean columnMajorOrder)
    Get the inverse of a matrix and store the inverse in a flat matrix.
    double[][]
    Get the lower triangular matrix associated with the decomposition.
    int
    Get the number of columns in the triangular matrices associated with an instance of this class.
    int
    Get the number of rows in the triangular matrices associated with an instance of this class.
    Get the permutation to which this decomposition applies.
    double[][]
    Get the upper triangular matrix associated with the decomposition.
    boolean
    Determine if the matrix to which this decomposition applies is nonsingular.
    double[]
    solve(double[] b)
    Solve the system of linear equations Ax = b where x and b are vectors and A is the matrix to which this decomposition applies.
    void
    solve(double[][] x, double[][] b)
    Solve the system of linear equations AX = B where X and B are matrices and A is the matrix to which this decomposition applies.
    void
    solve(double[] x, double[] b)
    Solve the system of linear equations Ax = b where x and b are vectors and A is the matrix to which this decomposition applies.

    Methods inherited from class java.lang.Object

    clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
  • Constructor Details

    • LUDecomp

      public LUDecomp(double[][] A) throws IllegalArgumentException
      Constructor.
      Parameters:
      A - the matrix whose LU decomposition is to be computed
      Throws:
      IllegalArgumentException - the input matrix is ill formed
    • LUDecomp

      public LUDecomp(double[] flatMatrix, boolean columnMajorOrder, int m, int n) throws IllegalArgumentException
      Constructor given a flat representation of a matrix. The elements of the matrix can be stored in one of two orders. Given a matrix Aij, if the row index varies the fastest, the flat representation will consist of the sequence A00, A10, .... Otherwise the sequence is A00, A 01, ....
      Parameters:
      flatMatrix - a one-dimensional array representing a two dimensional matrix
      columnMajorOrder - true if the inverse is stored in column major order (e.g., the order used by Fortran); false if the inverse is stored in row-major-order (e.g., the order used by C)
      m - the number of rows
      n - the number of columns
      Throws:
      IllegalArgumentException - the flat matrix is too small to represent an m by n matrix
    • LUDecomp

      public LUDecomp(double[][] A, int m, int n, boolean strict) throws IllegalArgumentException
      Constructor specifying matrix dimensions. This is provided to allow a single array to be used multiple times in order to avoid unnecessary array allocations.
      Parameters:
      A - the matrix whose LU decomposition is to be computed
      m - the number of rows
      n - the number of columns
      strict - true if matrixA must have exactly the number of rows and columns required; false if A may have rows or columns longer than needed, in which case the additional entries will be ignored
      Throws:
      IllegalArgumentException - the matrix and other parameters are not consistent with each other
    • LUDecomp

      public LUDecomp(double[][] A, double[][] LU) throws IllegalArgumentException
      Constructor providing an LU matrix. The dimensions of matrixLU are used to compute the number of rows and columns. matrixLU must not be modified after the constructor is called until the instance created is no longer used.
      Parameters:
      A - the matrix whose LU decomposition is to be computed.
      LU - the matrix used to store the triangular matrices L and U in a compact form
      Throws:
      IllegalArgumentException - the matrices are not consistent with each other
  • Method Details

    • getNumberOfRows

      public int getNumberOfRows()
      Description copied from interface: TriangularDecomp
      Get the number of rows in the triangular matrices associated with an instance of this class.
      Specified by:
      getNumberOfRows in interface TriangularDecomp
      Returns:
      the number of rows.
    • getNumberOfColumns

      public int getNumberOfColumns()
      Description copied from interface: TriangularDecomp
      Get the number of columns in the triangular matrices associated with an instance of this class.
      Specified by:
      getNumberOfColumns in interface TriangularDecomp
      Returns:
      the number of columns.
    • getP

      public Permutation getP()
      Description copied from interface: TriangularDecomp
      Get the permutation to which this decomposition applies. The permutation returned has the property that its getMatrix() method returns the permutation matrix.
      Specified by:
      getP in interface TriangularDecomp
      Returns:
      the permutation
    • getL

      public double[][] getL()
      Description copied from interface: TriangularDecomp
      Get the lower triangular matrix associated with the decomposition.
      Specified by:
      getL in interface TriangularDecomp
      Returns:
      the lower triangular matrix
    • getU

      public double[][] getU()
      Description copied from interface: TriangularDecomp
      Get the upper triangular matrix associated with the decomposition. For Cholesky decomposition, U is the transpose of the lower triangular matrix returned by getL().
      Specified by:
      getU in interface TriangularDecomp
      Returns:
      the upper triangular matrix
    • det

      public double det() throws IllegalStateException
      Description copied from interface: TriangularDecomp
      Get the determinate of the matrix to which this decomposition applies.
      Specified by:
      det in interface TriangularDecomp
      Returns:
      the determinate
      Throws:
      IllegalStateException - - if the matrix to which this decomposition applies is not a square matrix
    • isNonsingular

      public boolean isNonsingular()
      Description copied from interface: TriangularDecomp
      Determine if the matrix to which this decomposition applies is nonsingular.
      Specified by:
      isNonsingular in interface TriangularDecomp
      Returns:
      true if the matrix is nonsingular; false otherwise
    • solve

      public double[] solve(double[] b) throws IllegalArgumentException, IllegalStateException
      Description copied from interface: TriangularDecomp
      Solve the system of linear equations Ax = b where x and b are vectors and A is the matrix to which this decomposition applies.
      Specified by:
      solve in interface TriangularDecomp
      Parameters:
      b - the vector b in the equation Ax = b
      Returns:
      the vector x that satisfies the equation Ax = b
      Throws:
      IllegalArgumentException - the argument has the wrong size
      IllegalStateException - the matrix this object decomposed is a singular matrix
    • solve

      public void solve(double[] x, double[] b) throws IllegalArgumentException, IllegalStateException
      Description copied from interface: TriangularDecomp
      Solve the system of linear equations Ax = b where x and b are vectors and A is the matrix to which this decomposition applies. The vector x will be modified.
      Specified by:
      solve in interface TriangularDecomp
      Parameters:
      x - the vector x that satisfies the equation Ax = b
      b - the vector b in the equation Ax = b
      Throws:
      IllegalArgumentException - the argument has the wrong size
      IllegalStateException - the matrix this object decomposed is a singular matrix
    • solve

      public void solve(double[][] x, double[][] b) throws IllegalArgumentException, IllegalStateException
      Description copied from interface: TriangularDecomp
      Solve the system of linear equations AX = B where X and B are matrices and A is the matrix to which this decomposition applies. The matrix X will be modified.
      Specified by:
      solve in interface TriangularDecomp
      Parameters:
      x - the matrix X that satisfies the equation AX = B
      b - the matrix B in the equation AX = B
      Throws:
      IllegalArgumentException - the arguments have the wrong size
      IllegalStateException - the matrix this object decomposed is a singular matrix
    • getInverse

      public double[][] getInverse()
      Description copied from interface: TriangularDecomp
      Get the inverse of a matrix. The matrix A whose inverse is computed satisfies PA = LU, where P is the permutation matrix associated with this instance, L is the lower triangular matrix associated with this instance, and U is the upper triangular matrix associated with this instance.

      For some cases (e.g., Cholesky decomposition), the permutation matrix is the identity matrix.

      Specified by:
      getInverse in interface TriangularDecomp
      Returns:
      the inverse of the matrix that this instance decomposed
    • getInverse

      public void getInverse(double[][] result) throws IllegalStateException, IllegalArgumentException
      Description copied from interface: TriangularDecomp
      Get the inverse of a matrix and store the inverse. The matrix A whose inverse is computed satisfies PA = LU, where P is the permutation matrix associated with this instance, L is the lower triangular matrix associated with this instance, and U is the upper triangular matrix associated with this instance.
      Specified by:
      getInverse in interface TriangularDecomp
      Parameters:
      result - the matrix holding the inverse
      Throws:
      IllegalStateException - the matrix is singular
      IllegalArgumentException - - if result has the wrong dimensions
    • getInverse

      public void getInverse(double[] result, boolean columnMajorOrder) throws IllegalStateException, IllegalArgumentException
      Description copied from interface: TriangularDecomp
      Get the inverse of a matrix and store the inverse in a flat matrix. The matrix A whose inverse is computed satisfies A = LLT, where L is the lower triangular matrix associated with this instance.
      Specified by:
      getInverse in interface TriangularDecomp
      Parameters:
      result - the matrix holding the inverse
      columnMajorOrder - true if the inverse is stored in column major order (e.g., the order used by Fortran); false if the inverse is stored in row-major-order (e.g., the order used by C)
      Throws:
      IllegalStateException - the matrix is singular
      IllegalArgumentException - - if result has the wrong dimensions