Parallel Colt 0.7.2

cern.colt.matrix.tfloat.algo
Class FloatAlgebra

java.lang.Object
  extended by cern.colt.PersistentObject
      extended by cern.colt.matrix.tfloat.algo.FloatAlgebra
All Implemented Interfaces:
Serializable, Cloneable

public class FloatAlgebra
extends PersistentObject

Linear algebraic matrix operations operating on FloatMatrix2D; concentrates most functionality of this package.

Version:
1.0, 09/24/99, 1.1 10/20/2007
Author:
wolfgang.hoschek@cern.ch, Piotr Wendykier (piotr.wendykier@gmail.com)
See Also:
Serialized Form

Field Summary
static FloatAlgebra DEFAULT
          A default Algebra object; has FloatProperty.DEFAULT attached for tolerance.
static FloatAlgebra ZERO
          A default Algebra object; has FloatProperty.ZERO attached for tolerance.
 
Constructor Summary
FloatAlgebra()
          Constructs a new instance with an equality tolerance given by Property.DEFAULT.tolerance().
FloatAlgebra(float tolerance)
          Constructs a new instance with the given equality tolerance.
 
Method Summary
 FloatMatrix1D backwardSolve(FloatMatrix2D U, FloatMatrix1D b)
           
 FloatCholeskyDecomposition chol(FloatMatrix2D matrix)
          Constructs and returns the cholesky-decomposition of the given matrix.
 Object clone()
          Returns a copy of the receiver.
 float cond(FloatMatrix2D A)
          Returns the condition of matrix A, which is the ratio of largest to smallest singular value.
 float det(FloatMatrix2D A)
          Returns the determinant of matrix A.
 FloatEigenvalueDecomposition eig(FloatMatrix2D matrix)
          Constructs and returns the Eigenvalue-decomposition of the given matrix.
 FloatMatrix1D forwardSolve(FloatMatrix2D L, FloatMatrix1D b)
           
static float hypot(float a, float b)
          Returns sqrt(a^2 + b^2) without under/overflow.
static FloatFloatFunction hypotFunction()
          Returns sqrt(a^2 + b^2) without under/overflow.
 FloatMatrix2D inverse(FloatMatrix2D A)
          Returns the inverse or pseudo-inverse of matrix A.
 FComplexMatrix1D kron(FComplexMatrix1D x, FComplexMatrix1D y)
          Computes the Kronecker product of two complex matrices.
 FloatMatrix1D kron(float[] A, float[] B)
          Computes the Kronecker product of two arrays.
 FloatMatrix1D kron(FloatMatrix1D x, FloatMatrix1D y)
          Computes the Kronecker product of two real matrices.
 FloatLUDecomposition lu(FloatMatrix2D matrix)
          Constructs and returns the LU-decomposition of the given matrix.
 float mult(FloatMatrix1D x, FloatMatrix1D y)
          Inner product of two vectors; Sum(x[i] * y[i]).
 FloatMatrix1D mult(FloatMatrix2D A, FloatMatrix1D y)
          Linear algebraic matrix-vector multiplication; z = A * y.
 FloatMatrix2D mult(FloatMatrix2D A, FloatMatrix2D B)
          Linear algebraic matrix-matrix multiplication; C = A x B.
 FloatMatrix2D multOuter(FloatMatrix1D x, FloatMatrix1D y, FloatMatrix2D A)
          Outer product of two vectors; Sets A[i,j] = x[i] * y[j].
 float norm(FloatMatrix1D x, Norm type)
           
 float norm(FloatMatrix2D A, Norm type)
           
 float norm1(FloatMatrix1D x)
          Returns the one-norm of vector x, which is Sum(abs(x[i])).
 float norm1(FloatMatrix2D A)
          Returns the one-norm of matrix A, which is the maximum absolute column sum.
 float norm2(FloatMatrix1D x)
          Returns the two-norm (aka euclidean norm) of vector x; equivalent to Sqrt(mult(x,x)).
 float norm2(FloatMatrix2D A)
          Returns the two-norm of matrix A, which is the maximum singular value; obtained from SVD.
 float normF(FComplexMatrix2D A)
          Returns the Frobenius norm of matrix A, which is Sqrt(Sum(A[i,j]2)).
 float normF(FloatMatrix1D A)
          Returns the Frobenius norm of matrix A, which is Sqrt(Sum(A[i]2)).
 float normF(FloatMatrix2D A)
          Returns the Frobenius norm of matrix A, which is Sqrt(Sum(A[i,j]2)).
 float normInfinity(FloatMatrix1D x)
          Returns the infinity norm of vector x, which is Max(abs(x[i])).
 float normInfinity(FloatMatrix2D A)
          Returns the infinity norm of matrix A, which is the maximum absolute row sum.
 FloatMatrix1D permute(FloatMatrix1D A, int[] indexes, float[] work)
          Modifies the given vector A such that it is permuted as specified; Useful for pivoting.
 FloatMatrix2D permute(FloatMatrix2D A, int[] rowIndexes, int[] columnIndexes)
          Constructs and returns a new row and column permuted selection view of matrix A; equivalent to FloatMatrix2D.viewSelection(int[],int[]).
 FloatMatrix2D permuteColumns(FloatMatrix2D A, int[] indexes, int[] work)
          Modifies the given matrix A such that it's columns are permuted as specified; Useful for pivoting.
 FloatMatrix2D permuteRows(FloatMatrix2D A, int[] indexes, int[] work)
          Modifies the given matrix A such that it's rows are permuted as specified; Useful for pivoting.
 FloatMatrix2D pow(FloatMatrix2D A, int p)
          Linear algebraic matrix power; B = Ak <==> B = A*A*...*A.
 FloatProperty property()
          Returns the property object attached to this Algebra, defining tolerance.
 FloatQRDecomposition qr(FloatMatrix2D matrix)
          Constructs and returns the QR-decomposition of the given matrix.
 int rank(FloatMatrix2D A)
          Returns the effective numerical rank of matrix A, obtained from Singular Value Decomposition.
 void setProperty(FloatProperty property)
          Attaches the given property object to this Algebra, defining tolerance.
 FloatMatrix1D solve(FloatMatrix2D A, FloatMatrix1D b)
          Solves A*x = b.
 FloatMatrix2D solve(FloatMatrix2D A, FloatMatrix2D B)
          Solves A*X = B.
 FloatMatrix2D solveTranspose(FloatMatrix2D A, FloatMatrix2D B)
          Solves X*A = B, which is also A'*X' = B'.
 FloatMatrix2D subMatrix(FloatMatrix2D A, int[] rowIndexes, int columnFrom, int columnTo)
          Copies the columns of the indicated rows into a new sub matrix.
 FloatMatrix2D subMatrix(FloatMatrix2D A, int rowFrom, int rowTo, int[] columnIndexes)
          Copies the rows of the indicated columns into a new sub matrix.
 FloatMatrix2D subMatrix(FloatMatrix2D A, int fromRow, int toRow, int fromColumn, int toColumn)
          Constructs and returns a new sub-range view which is the sub matrix A[fromRow..toRow,fromColumn..toColumn].
 FloatSingularValueDecomposition svd(FloatMatrix2D matrix)
          Constructs and returns the SingularValue-decomposition of the given matrix.
 FloatSingularValueDecompositionDC svdDC(FloatMatrix2D matrix)
          Constructs and returns the SingularValue-decomposition of the given matrix.
 String toString(FloatMatrix2D matrix)
          Returns a String with (propertyName, propertyValue) pairs.
 String toVerboseString(FloatMatrix2D matrix)
          Returns the results of toString(A) and additionally the results of all sorts of decompositions applied to the given matrix.
 float trace(FloatMatrix2D A)
          Returns the sum of the diagonal elements of matrix A; Sum(A[i,i]).
 FloatMatrix2D transpose(FloatMatrix2D A)
          Constructs and returns a new view which is the transposition of the given matrix A.
 FloatMatrix2D trapezoidalLower(FloatMatrix2D A)
          Modifies the matrix to be a lower trapezoidal matrix.
 float vectorNorm2(FloatMatrix2D X)
          Returns the two-norm (aka euclidean norm) of vector X.vectorize();
 float vectorNorm2(FloatMatrix3D X)
          Returns the two-norm (aka euclidean norm) of vector X.vectorize();
 FloatMatrix2D xmultOuter(FloatMatrix1D x, FloatMatrix1D y)
          Outer product of two vectors; Returns a matrix with A[i,j] = x[i] * y[j].
 FloatMatrix2D xpowSlow(FloatMatrix2D A, int k)
          Linear algebraic matrix power; B = Ak <==> B = A*A*...*A.
 
Methods inherited from class java.lang.Object
equals, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
 

Field Detail

DEFAULT

public static final FloatAlgebra DEFAULT
A default Algebra object; has FloatProperty.DEFAULT attached for tolerance. Allows ommiting to construct an Algebra object time and again. Note that this Algebra object is immutable. Any attempt to assign a new Property object to it (via method setProperty), or to alter the tolerance of its property object (via property().setTolerance(...)) will throw an exception.


ZERO

public static final FloatAlgebra ZERO
A default Algebra object; has FloatProperty.ZERO attached for tolerance. Allows ommiting to construct an Algebra object time and again. Note that this Algebra object is immutable. Any attempt to assign a new Property object to it (via method setProperty), or to alter the tolerance of its property object (via property().setTolerance(...)) will throw an exception.

Constructor Detail

FloatAlgebra

public FloatAlgebra()
Constructs a new instance with an equality tolerance given by Property.DEFAULT.tolerance().


FloatAlgebra

public FloatAlgebra(float tolerance)
Constructs a new instance with the given equality tolerance.

Parameters:
tolerance - the tolerance to be used for equality operations.
Method Detail

chol

public FloatCholeskyDecomposition chol(FloatMatrix2D matrix)
Constructs and returns the cholesky-decomposition of the given matrix.


clone

public Object clone()
Returns a copy of the receiver. The attached property object is also copied. Hence, the property object of the copy is mutable.

Overrides:
clone in class PersistentObject
Returns:
a copy of the receiver.

cond

public float cond(FloatMatrix2D A)
Returns the condition of matrix A, which is the ratio of largest to smallest singular value.


det

public float det(FloatMatrix2D A)
Returns the determinant of matrix A.

Returns:
the determinant.

eig

public FloatEigenvalueDecomposition eig(FloatMatrix2D matrix)
Constructs and returns the Eigenvalue-decomposition of the given matrix.


hypot

public static float hypot(float a,
                          float b)
Returns sqrt(a^2 + b^2) without under/overflow.


hypotFunction

public static FloatFloatFunction hypotFunction()
Returns sqrt(a^2 + b^2) without under/overflow.


inverse

public FloatMatrix2D inverse(FloatMatrix2D A)
Returns the inverse or pseudo-inverse of matrix A.

Returns:
a new independent matrix; inverse(matrix) if the matrix is square, pseudoinverse otherwise.

lu

public FloatLUDecomposition lu(FloatMatrix2D matrix)
Constructs and returns the LU-decomposition of the given matrix.


kron

public FComplexMatrix1D kron(FComplexMatrix1D x,
                             FComplexMatrix1D y)
Computes the Kronecker product of two complex matrices.

Parameters:
x -
y -
Returns:
the Kronecker product of two complex matrices

kron

public FloatMatrix1D kron(float[] A,
                          float[] B)
Computes the Kronecker product of two arrays.

Parameters:
A -
B -
Returns:
the Kronecker product of two arrays

kron

public FloatMatrix1D kron(FloatMatrix1D x,
                          FloatMatrix1D y)
Computes the Kronecker product of two real matrices.

Parameters:
x -
y -
Returns:
the Kronecker product of two real matrices

mult

public float mult(FloatMatrix1D x,
                  FloatMatrix1D y)
Inner product of two vectors; Sum(x[i] * y[i]). Also known as dot product.
Equivalent to x.zDotProduct(y).

Parameters:
x - the first source vector.
y - the second source matrix.
Returns:
the inner product.
Throws:
IllegalArgumentException - if x.size() != y.size().

mult

public FloatMatrix1D mult(FloatMatrix2D A,
                          FloatMatrix1D y)
Linear algebraic matrix-vector multiplication; z = A * y. z[i] = Sum(A[i,j] * y[j]), i=0..A.rows()-1, j=0..y.size()-1.

Parameters:
A - the source matrix.
y - the source vector.
Returns:
z; a new vector with z.size()==A.rows().
Throws:
IllegalArgumentException - if A.columns() != y.size().

mult

public FloatMatrix2D mult(FloatMatrix2D A,
                          FloatMatrix2D B)
Linear algebraic matrix-matrix multiplication; C = A x B. C[i,j] = Sum(A[i,k] * B[k,j]), k=0..n-1.
Matrix shapes: A(m x n), B(n x p), C(m x p).

Parameters:
A - the first source matrix.
B - the second source matrix.
Returns:
C; a new matrix holding the results, with C.rows()=A.rows(), C.columns()==B.columns().
Throws:
IllegalArgumentException - if B.rows() != A.columns().

multOuter

public FloatMatrix2D multOuter(FloatMatrix1D x,
                               FloatMatrix1D y,
                               FloatMatrix2D A)
Outer product of two vectors; Sets A[i,j] = x[i] * y[j].

Parameters:
x - the first source vector.
y - the second source vector.
A - the matrix to hold the results. Set this parameter to null to indicate that a new result matrix shall be constructed.
Returns:
A (for convenience only).
Throws:
IllegalArgumentException - if A.rows() != x.size() || A.columns() != y.size().

norm1

public float norm1(FloatMatrix1D x)
Returns the one-norm of vector x, which is Sum(abs(x[i])).


norm1

public float norm1(FloatMatrix2D A)
Returns the one-norm of matrix A, which is the maximum absolute column sum.


norm2

public float norm2(FloatMatrix1D x)
Returns the two-norm (aka euclidean norm) of vector x; equivalent to Sqrt(mult(x,x)).


vectorNorm2

public float vectorNorm2(FloatMatrix2D X)
Returns the two-norm (aka euclidean norm) of vector X.vectorize();


vectorNorm2

public float vectorNorm2(FloatMatrix3D X)
Returns the two-norm (aka euclidean norm) of vector X.vectorize();


norm

public float norm(FloatMatrix2D A,
                  Norm type)

norm

public float norm(FloatMatrix1D x,
                  Norm type)

norm2

public float norm2(FloatMatrix2D A)
Returns the two-norm of matrix A, which is the maximum singular value; obtained from SVD.


normF

public float normF(FloatMatrix2D A)
Returns the Frobenius norm of matrix A, which is Sqrt(Sum(A[i,j]2)).


normF

public float normF(FloatMatrix1D A)
Returns the Frobenius norm of matrix A, which is Sqrt(Sum(A[i]2)).


normF

public float normF(FComplexMatrix2D A)
Returns the Frobenius norm of matrix A, which is Sqrt(Sum(A[i,j]2)).


normInfinity

public float normInfinity(FloatMatrix1D x)
Returns the infinity norm of vector x, which is Max(abs(x[i])).


normInfinity

public float normInfinity(FloatMatrix2D A)
Returns the infinity norm of matrix A, which is the maximum absolute row sum.


permute

public FloatMatrix1D permute(FloatMatrix1D A,
                             int[] indexes,
                             float[] work)
Modifies the given vector A such that it is permuted as specified; Useful for pivoting. Cell A[i] will go into cell A[indexes[i]].

Example:

   Reordering
   [A,B,C,D,E] with indexes [0,4,2,3,1] yields 
   [A,E,C,D,B]
   In other words A[0]<--A[0], A[1]<--A[4], A[2]<--A[2], A[3]<--A[3], A[4]<--A[1].
 
   Reordering
   [A,B,C,D,E] with indexes [0,4,1,2,3] yields 
   [A,E,B,C,D]
   In other words A[0]<--A[0], A[1]<--A[4], A[2]<--A[1], A[3]<--A[2], A[4]<--A[3].
 
 

Parameters:
A - the vector to permute.
indexes - the permutation indexes, must satisfy indexes.length==A.size() && indexes[i] >= 0 && indexes[i] < A.size() ;
work - the working storage, must satisfy work.length >= A.size(); set work==null if you don't care about performance.
Returns:
the modified A (for convenience only).
Throws:
IndexOutOfBoundsException - if indexes.length != A.size().

permute

public FloatMatrix2D permute(FloatMatrix2D A,
                             int[] rowIndexes,
                             int[] columnIndexes)
Constructs and returns a new row and column permuted selection view of matrix A; equivalent to FloatMatrix2D.viewSelection(int[],int[]). The returned matrix is backed by this matrix, so changes in the returned matrix are reflected in this matrix, and vice-versa. Use idioms like result = permute(...).copy() to generate an independent sub matrix.

Returns:
the new permuted selection view.

permuteColumns

public FloatMatrix2D permuteColumns(FloatMatrix2D A,
                                    int[] indexes,
                                    int[] work)
Modifies the given matrix A such that it's columns are permuted as specified; Useful for pivoting. Column A[i] will go into column A[indexes[i]]. Equivalent to permuteRows(transpose(A), indexes, work).

Parameters:
A - the matrix to permute.
indexes - the permutation indexes, must satisfy indexes.length==A.columns() && indexes[i] >= 0 && indexes[i] < A.columns() ;
work - the working storage, must satisfy work.length >= A.columns(); set work==null if you don't care about performance.
Returns:
the modified A (for convenience only).
Throws:
IndexOutOfBoundsException - if indexes.length != A.columns().

permuteRows

public FloatMatrix2D permuteRows(FloatMatrix2D A,
                                 int[] indexes,
                                 int[] work)
Modifies the given matrix A such that it's rows are permuted as specified; Useful for pivoting. Row A[i] will go into row A[indexes[i]].

Example:

   Reordering
   [A,B,C,D,E] with indexes [0,4,2,3,1] yields 
   [A,E,C,D,B]
   In other words A[0]<--A[0], A[1]<--A[4], A[2]<--A[2], A[3]<--A[3], A[4]<--A[1].
 
   Reordering
   [A,B,C,D,E] with indexes [0,4,1,2,3] yields 
   [A,E,B,C,D]
   In other words A[0]<--A[0], A[1]<--A[4], A[2]<--A[1], A[3]<--A[2], A[4]<--A[3].
 
 

Parameters:
A - the matrix to permute.
indexes - the permutation indexes, must satisfy indexes.length==A.rows() && indexes[i] >= 0 && indexes[i] < A.rows() ;
work - the working storage, must satisfy work.length >= A.rows(); set work==null if you don't care about performance.
Returns:
the modified A (for convenience only).
Throws:
IndexOutOfBoundsException - if indexes.length != A.rows().

pow

public FloatMatrix2D pow(FloatMatrix2D A,
                         int p)
Linear algebraic matrix power; B = Ak <==> B = A*A*...*A. Implementation: Based on logarithms of 2, memory usage minimized.

Parameters:
A - the source matrix; must be square; stays unaffected by this operation.
p - the exponent, can be any number.
Returns:
B, a newly constructed result matrix; storage-independent of A.
Throws:
IllegalArgumentException - if !property().isSquare(A).

property

public FloatProperty property()
Returns the property object attached to this Algebra, defining tolerance.

Returns:
the Property object.
See Also:
setProperty(FloatProperty)

qr

public FloatQRDecomposition qr(FloatMatrix2D matrix)
Constructs and returns the QR-decomposition of the given matrix.


rank

public int rank(FloatMatrix2D A)
Returns the effective numerical rank of matrix A, obtained from Singular Value Decomposition.


setProperty

public void setProperty(FloatProperty property)
Attaches the given property object to this Algebra, defining tolerance.

Parameters:
property - the Property object to be attached.
Throws:
UnsupportedOperationException - if this==DEFAULT && property!=this.property() - The DEFAULT Algebra object is immutable.
UnsupportedOperationException - if this==ZERO && property!=this.property() - The ZERO Algebra object is immutable.
See Also:
property

backwardSolve

public FloatMatrix1D backwardSolve(FloatMatrix2D U,
                                   FloatMatrix1D b)

forwardSolve

public FloatMatrix1D forwardSolve(FloatMatrix2D L,
                                  FloatMatrix1D b)

solve

public FloatMatrix1D solve(FloatMatrix2D A,
                           FloatMatrix1D b)
Solves A*x = b.

Returns:
x; a new independent matrix; solution if A is square, least squares solution otherwise.

solve

public FloatMatrix2D solve(FloatMatrix2D A,
                           FloatMatrix2D B)
Solves A*X = B.

Returns:
X; a new independent matrix; solution if A is square, least squares solution otherwise.

solveTranspose

public FloatMatrix2D solveTranspose(FloatMatrix2D A,
                                    FloatMatrix2D B)
Solves X*A = B, which is also A'*X' = B'.

Returns:
X; a new independent matrix; solution if A is square, least squares solution otherwise.

subMatrix

public FloatMatrix2D subMatrix(FloatMatrix2D A,
                               int[] rowIndexes,
                               int columnFrom,
                               int columnTo)
Copies the columns of the indicated rows into a new sub matrix. sub[0..rowIndexes.length-1,0..columnTo-columnFrom] = A[rowIndexes(:),columnFrom..columnTo] ; The returned matrix is not backed by this matrix, so changes in the returned matrix are not reflected in this matrix, and vice-versa.

Parameters:
A - the source matrix to copy from.
rowIndexes - the indexes of the rows to copy. May be unsorted.
columnFrom - the index of the first column to copy (inclusive).
columnTo - the index of the last column to copy (inclusive).
Returns:
a new sub matrix; with sub.rows()==rowIndexes.length; sub.columns()==columnTo-columnFrom+1 .
Throws:
IndexOutOfBoundsException - if columnFrom<0 || columnTo-columnFrom+1<0 || columnTo+1>matrix.columns() || for any row=rowIndexes[i]: row < 0 || row >= matrix.rows() .

subMatrix

public FloatMatrix2D subMatrix(FloatMatrix2D A,
                               int rowFrom,
                               int rowTo,
                               int[] columnIndexes)
Copies the rows of the indicated columns into a new sub matrix. sub[0..rowTo-rowFrom,0..columnIndexes.length-1] = A[rowFrom..rowTo,columnIndexes(:)] ; The returned matrix is not backed by this matrix, so changes in the returned matrix are not reflected in this matrix, and vice-versa.

Parameters:
A - the source matrix to copy from.
rowFrom - the index of the first row to copy (inclusive).
rowTo - the index of the last row to copy (inclusive).
columnIndexes - the indexes of the columns to copy. May be unsorted.
Returns:
a new sub matrix; with sub.rows()==rowTo-rowFrom+1; sub.columns()==columnIndexes.length .
Throws:
IndexOutOfBoundsException - if rowFrom<0 || rowTo-rowFrom+1<0 || rowTo+1>matrix.rows() || for any col=columnIndexes[i]: col < 0 || col >= matrix.columns() .

subMatrix

public FloatMatrix2D subMatrix(FloatMatrix2D A,
                               int fromRow,
                               int toRow,
                               int fromColumn,
                               int toColumn)
Constructs and returns a new sub-range view which is the sub matrix A[fromRow..toRow,fromColumn..toColumn]. The returned matrix is backed by this matrix, so changes in the returned matrix are reflected in this matrix, and vice-versa. Use idioms like result = subMatrix(...).copy() to generate an independent sub matrix.

Parameters:
A - the source matrix.
fromRow - The index of the first row (inclusive).
toRow - The index of the last row (inclusive).
fromColumn - The index of the first column (inclusive).
toColumn - The index of the last column (inclusive).
Returns:
a new sub-range view.
Throws:
IndexOutOfBoundsException - if fromColumn<0 || toColumn-fromColumn+1<0 || toColumn>=A.columns() || fromRow<0 || toRow-fromRow+1<0 || toRow>=A.rows()

svd

public FloatSingularValueDecomposition svd(FloatMatrix2D matrix)
Constructs and returns the SingularValue-decomposition of the given matrix.


svdDC

public FloatSingularValueDecompositionDC svdDC(FloatMatrix2D matrix)
Constructs and returns the SingularValue-decomposition of the given matrix. This is a divide-and-conquer version.


toString

public String toString(FloatMatrix2D matrix)
Returns a String with (propertyName, propertyValue) pairs. Useful for debugging or to quickly get the rough picture. For example,
   cond          : 14.073264490042144
   det           : Illegal operation or error: Matrix must be square.
   norm1         : 0.9620244354009628
   norm2         : 3.0
   normF         : 1.304841791648992
   normInfinity  : 1.5406551198102534
   rank          : 3
   trace         : 0
 
 


toVerboseString

public String toVerboseString(FloatMatrix2D matrix)
Returns the results of toString(A) and additionally the results of all sorts of decompositions applied to the given matrix. Useful for debugging or to quickly get the rough picture. For example,
   A = 3 x 3 matrix
   249  66  68
   104 214 108
   144 146 293
 
   cond         : 3.931600417472078
   det          : 9638870.0
   norm1        : 497.0
   norm2        : 473.34508217011404
   normF        : 516.873292016525
   normInfinity : 583.0
   rank         : 3
   trace        : 756.0
 
   density                      : 1.0
   isDiagonal                   : false
   isDiagonallyDominantByColumn : true
   isDiagonallyDominantByRow    : true
   isIdentity                   : false
   isLowerBidiagonal            : false
   isLowerTriangular            : false
   isNonNegative                : true
   isOrthogonal                 : false
   isPositive                   : true
   isSingular                   : false
   isSkewSymmetric              : false
   isSquare                     : true
   isStrictlyLowerTriangular    : false
   isStrictlyTriangular         : false
   isStrictlyUpperTriangular    : false
   isSymmetric                  : false
   isTriangular                 : false
   isTridiagonal                : false
   isUnitTriangular             : false
   isUpperBidiagonal            : false
   isUpperTriangular            : false
   isZero                       : false
   lowerBandwidth               : 2
   semiBandwidth                : 3
   upperBandwidth               : 2
 
   -----------------------------------------------------------------------------
   LUDecompositionQuick(A) --> isNonSingular(A), det(A), pivot, L, U, inverse(A)
   -----------------------------------------------------------------------------
   isNonSingular = true
   det = 9638870.0
   pivot = [0, 1, 2]
 
   L = 3 x 3 matrix
   1        0       0
   0.417671 1       0
   0.578313 0.57839 1
 
   U = 3 x 3 matrix
   249  66         68       
   0 186.433735  79.598394
   0   0        207.635819
 
   inverse(A) = 3 x 3 matrix
   0.004869 -0.000976 -0.00077 
   -0.001548  0.006553 -0.002056
   -0.001622 -0.002786  0.004816
 
   -----------------------------------------------------------------
   QRDecomposition(A) --> hasFullRank(A), H, Q, R, pseudo inverse(A)
   -----------------------------------------------------------------
   hasFullRank = true
 
   H = 3 x 3 matrix
   1.814086 0        0
   0.34002  1.903675 0
   0.470797 0.428218 2
 
   Q = 3 x 3 matrix
   -0.814086  0.508871  0.279845
   -0.34002  -0.808296  0.48067 
   -0.470797 -0.296154 -0.831049
 
   R = 3 x 3 matrix
   -305.864349 -195.230337 -230.023539
   0        -182.628353  467.703164
   0           0        -309.13388 
 
   pseudo inverse(A) = 3 x 3 matrix
   0.006601  0.001998 -0.005912
   -0.005105  0.000444  0.008506
   -0.000905 -0.001555  0.002688
 
   --------------------------------------------------------------------------
   CholeskyDecomposition(A) --> isSymmetricPositiveDefinite(A), L, inverse(A)
   --------------------------------------------------------------------------
   isSymmetricPositiveDefinite = false
 
   L = 3 x 3 matrix
   15.779734  0         0       
   6.590732 13.059948  0       
   9.125629  6.573948 12.903724
 
   inverse(A) = Illegal operation or error: Matrix is not symmetric positive definite.
 
   ---------------------------------------------------------------------
   EigenvalueDecomposition(A) --> D, V, realEigenvalues, imagEigenvalues
   ---------------------------------------------------------------------
   realEigenvalues = 1 x 3 matrix
   462.796507 172.382058 120.821435
   imagEigenvalues = 1 x 3 matrix
   0 0 0
 
   D = 3 x 3 matrix
   462.796507   0          0       
   0        172.382058   0       
   0          0        120.821435
 
   V = 3 x 3 matrix
   -0.398877 -0.778282  0.094294
   -0.500327  0.217793 -0.806319
   -0.768485  0.66553   0.604862
 
   ---------------------------------------------------------------------
   SingularValueDecomposition(A) --> cond(A), rank(A), norm2(A), U, S, V
   ---------------------------------------------------------------------
   cond = 3.931600417472078
   rank = 3
   norm2 = 473.34508217011404
 
   U = 3 x 3 matrix
   0.46657  -0.877519  0.110777
   0.50486   0.161382 -0.847982
   0.726243  0.45157   0.51832 
 
   S = 3 x 3 matrix
   473.345082   0          0       
   0        169.137441   0       
   0          0        120.395013
 
   V = 3 x 3 matrix
   0.577296 -0.808174  0.116546
   0.517308  0.251562 -0.817991
   0.631761  0.532513  0.563301
 
 


trace

public float trace(FloatMatrix2D A)
Returns the sum of the diagonal elements of matrix A; Sum(A[i,i]).


transpose

public FloatMatrix2D transpose(FloatMatrix2D A)
Constructs and returns a new view which is the transposition of the given matrix A. Equivalent to A.viewDice(). This is a zero-copy transposition, taking O(1), i.e. constant time. The returned view is backed by this matrix, so changes in the returned view are reflected in this matrix, and vice-versa. Use idioms like result = transpose(A).copy() to generate an independent matrix.

Example:

2 x 3 matrix:
1, 2, 3
4, 5, 6
transpose ==> 3 x 2 matrix:
1, 4
2, 5
3, 6
transpose ==> 2 x 3 matrix:
1, 2, 3
4, 5, 6

Returns:
a new transposed view.

trapezoidalLower

public FloatMatrix2D trapezoidalLower(FloatMatrix2D A)
Modifies the matrix to be a lower trapezoidal matrix.

Returns:
A (for convenience only).

xmultOuter

public FloatMatrix2D xmultOuter(FloatMatrix1D x,
                                FloatMatrix1D y)
Outer product of two vectors; Returns a matrix with A[i,j] = x[i] * y[j].

Parameters:
x - the first source vector.
y - the second source vector.
Returns:
the outer product A.

xpowSlow

public FloatMatrix2D xpowSlow(FloatMatrix2D A,
                              int k)
Linear algebraic matrix power; B = Ak <==> B = A*A*...*A.

Parameters:
A - the source matrix; must be square.
k - the exponent, can be any number.
Returns:
a new result matrix.
Throws:
IllegalArgumentException - if !Testing.isSquare(A).

Parallel Colt 0.7.2

Jump to the Parallel Colt Homepage