Decompositions

Methods

================================ ============ Method Description ================================ ============ lu_ (A) LU decomposition luFactor_ (A) LU factors for solving equations luSolve_ (lu, piv, b) Solve equations using LU decomposition svd_ (a[, full_matrices]) Singular value decomposition svdvals_ (a) Compute singular values of a matrix diagsvd_ (s, M, N) Construct sigma matrix from singular values orth_ (a) Construct orthonormal basis for range of A using SVD nullSpace_ (a) Construct orthonormal basis for null space of A using SVD cholesky_ (a[, lower]) Cholesky decomposition choleskyBanded_ (a) Cholesky decomposition for banded matrices choFactor_ (a[, lower]) Cholesky factor for solving equations choSolve_ (C, b[, lower]) Solve equations using Cholesky decomposition choSolveBanded_ (C, b) Solve equations using Cholesky decomposition for banded matrices polar_ (a) Compute the polar decomposition qr_ (a) QR decomposition qrMultiply_ (a, c) QR decomposition, with Q multiplied by a matrix qrDelete_ (Q, R, k, p, alter) QR downdate on row or column deletions rq_ (a) RQ decomposition schur_ (a) Schur decomposition ================================ ============

Description

.. _lu:

.. function:: lu(A)

:param A: Matrix<Float64>
:return: (Matrix<Float64>, Matrix<Float64>, Matrix<Float64>)

Compute the LU decomposition :math:`A=PLU`. The three returned matrices
are :math:`P, L, U`, respectively.

.. _luFactor:

.. function:: luFactor(A)

:param A: Matrix<Float64>
:return: (Matrix<Float64>, Vector<Int32>)

Compute the LU decomposition. The first return value stores both :math:`L`
and :math:`U`. The second matrix stores the pivot information. The return
values can be sent to :ref:`luSolve <luSolve>` for equation solving.

.. _luSolve:

.. function:: luSolve(lu, piv, b)

:param lu: Matrix<Float64>
:param piv: Vector<Int32>
:param b: Matrix<Float64>
:return: Matrix<Float64>

Solve the equation :math:`Ax=b` using LU decomposition. The input arguments
should be obtained from :ref:`luFactor <luFactor>`.

.. _svd:

.. function:: svd(A[, full_matrices])

:param A: Matrix<Float64>
:param full_matrices: Bool
:return: (Matrix<Float64>, Vector<Float64>, Matrix<Float64>), the triple :math:`(U,s,V^T)`
    for the singular value decomposition.

Compute the singular value decomposition :math:`A=USV^T` of :math:`A`.

Let :math:`A` has dimension :math:`M\times N` and let :math:`K=\min(M,N)`.
If ``full_matrices`` is true, then both :math:`U` and :math:`V` are unitary,
with dimensions :math:`M\times M` and :math:`N\times N`, respectively.
If ``full_matrices`` is false, then :math:`U` has dimension :math:`M\times K`,
with columns orthonormal, and :math:`V` has dimension :math:`K\times N`, with
rows orthonormal. The singular values real and non-positive.

.. _svdvals:

.. function:: svdvals(A)

:param A: Matrix<Float64>
:return: Vector<Float64>, singular values of the matrix.

Return the singular values of the matrix :math:`A`.

.. _diagsvd:

.. function:: diagsvd(s, m, n)

:param s: Vector<Float64>, the singular values
:param m: Int64, number of rows
:param n: Int64, number of columns
:return: Matrix<Float64>, sigma matrix in SVD

Construct the sigma matrix in SVD from the vector of singular values, as well
as the dimensions :math:`m` and :math:`n`.

This is simply an :math:`m\times n` matrix with singular values on the diagonal.

.. _orth:

.. function:: orth(A[, rcond])

:param A: Matrix<Float64>
:param rcond: Float64, the relative condition number. Optional, default
    value 0.0 to be explained below.
:return: Matrix<Float64>, the columns form an orthonormal basis for the
    range of :math:`A`.

Compute an orthonormal basis for the range of :math:`A`.

This is done by first taking singular value decomposition :math:`A=USV^T`,
then take the columns of :math:`U` corresponding to "zero" singular
values. Singular values ``s`` smaller than ``rcond * max(s)`` are considered
zero. If ``rcond`` is set to zero, it is considered to be floating point eps
multiplied by the maximum dimension of :math:`A`.

.. _nullSpace:

.. function:: nullSpace(A[, rcond])

:param A: Matrix<Float64>
:param rcond: Float64, the relative condition number. Optional, default
    value 0.0 to be explained below.
:return: Matrix<Float64>, the columns form an orthonormal basis for the
    null space of :math:`A`.

Compute an orthonormal basis for the null space of :math:`A`.

This is done by first taking singular value decomposition :math:`A=USV^T`,
then take the rows of :math:`V` corresponding to "zero" singular
values. Singular values ``s`` smaller than ``rcond * max(s)`` are considered
zero. If ``rcond`` is set to zero, it is considered to be floating point eps
multiplied by the maximum dimension of :math:`A`.

.. _cholesky:

.. function:: cholesky(A[, lower])

:param A: Matrix<Float64>
:param lower: Bool, whether to compute upper- or lower-triangular Cholesky
    decomposition. Default to false for upper-triangular.
:return: Matrix<Float64>

Compute the Cholesky factorization of a symmetric/Hermitian positive-definite
matrix. If ``lower`` is true, the decomposition is :math:`A=LL^*`, Otherwise
the decomposition is :math:`A=U^*U`.

.. _choleskyBanded:

.. function:: choleskyBanded(A)

:param A: Matrix<Complex64>
:return: Matrix<Complex64>

Compute the Cholesky factorization of a banded matrix.

.. _choFactor:

.. function:: choFactor(A[, lower])

:param A: Matrix<Float64>
:param lower: Bool, whether to compute upper- or lower-triangular Cholesky
    decomposition. Default to false for upper-triangular.
:return: Matrix<Float64>

Compute the Cholesky factorization of a symmetric/Hermitian positive-definite
matrix. If ``lower`` is true, the decomposition is :math:`A=LL^*`, Otherwise
the decomposition is :math:`A=U^*U`.

This function *does not* zero-out the part of the result that is unused. That part
may contain arbitrary data. Hence the result should only be used as input to
:ref:`choSolve <choSolve>`.

.. _choSolve:

.. function:: choSolve(C, b, [, lower])

:param C: Matrix<Float64>, the Cholesky factor obtained from
    :ref:`choFactor <choFactor>`.
:param b: Matrix<Float64>, right side of the equation.
:param lower: Bool, whether the Cholesky factor is lower-triangular or
    upper triangular.
:return: Matrix<Float64>

Solve equation :math:`Ax = b` using the Cholesky factor :math:`C` of :math:`A`.

.. _choSolveBanded:

.. function:: choSolveBanded(C, b)

:param C: Matrix<Complex64>
:param b: Matrix<Complex64>
:return: Matrix<Complex64>

Solve equation :math:`Ax = b` using the Cholesky factor :math:`C` of banded matrix :math:`A`.

.. _polar:

.. function:: polar(A)

:param A: Matrix<Float64>
:return: (Matrix<Float64>, Matrix<Float64>), the pair of matrices U and P.

Compute the polar decomposition of a matrix :math:`A=UP`. Suppose :math:`A`
has dimension :math:`m\times n`, then :math:`U` has dimension :math:`m\times n`
and :math:`P` has dimension :math:`n\times n`. If :math:`A` is a square matrix,
then :math:`U` is unitary. If :math:`m > n`, then the columns of :math:`U` are
orthonormal. If :math:`m < n`, then the rows of :math:`U` are orthonormal. The
matrix :math:`P` is positive semi-definite. If :math:`A` is nonsingular, then
:math:`P` is positive definite.

Let :math:`A=WSV^T` be the singular value decomposition of :math:`A`, then the
polar decomposition is given by :math:`U=WV^T` and :math:`P=VSV^T`.

Example:

1. Compute polar decomposition of the following square matrix.

.. math::
    
    A = \begin{bmatrix} 1 & -1 \\ 2 & 4 \end{bmatrix}

::

    var A = matrix<Float64>([[1.0, -1.0], [2.0, 4.0]])
    var (U, P) = polar(A)

The result is:

.. math::

    U = \begin{bmatrix} 0.857493 & -0.514496 \\ 0.514496 & 0.857493 \end{bmatrix}

    P = \begin{bmatrix} 1.886484 & 1.200490 \\ 1.200490 & 3.944467 \end{bmatrix}

2. Case when the input matrix has :math:`m < n`.

3. Case when the input matrix has :math:`m > n`.

.. _qr:

.. function:: qr(a)

:param a: Matrix<Float64>
:return: (Matrix<Float64>, Matrix<Float64>)

Compute the QR decomposition :math:`A=QR`, where :math:`Q` is unitary/orthogonal,
and :math:`R` is upper triangular. The two matrices :math:`Q` and :math:`R` are
returned.

.. _qrMultiply:

.. function:: qrMultiply(a, c)

:param a: Matrix<Float64>
:param c: Matrix<Float64>
:return: (Matrix<Float64>, Matrix<Float64>)

Compute the QR decomposition :math:`A=QR`, where :math:`Q` multiplied by
another matrix :math:`C` is returned. The two returned matrices are :math:`QC` and
:math:`R`.

.. _qrDelete:

.. function:: qrDelete(Q, R, k, p, alter)

:param Q: Matrix<Float64>
:param R: Matrix<Float64>
:param k: Int64
:param p: Int64
:param alter: String, one of "row" and "col"
:return: (Matrix<Float64>, Matrix<Float64>)

Given QR decomposition :math:`A=QR` of a matrix :math:`A`, downdate the
decomposition when :math:`p` rows or columns have been removed from :math:`A`
starting at row or column :math:`k`.

The method implemented is an :math:`O(n^2)` algorithm according to
Section 6.5 of *Matrix Computations, 4th edition*.

.. _rq:

.. function:: rq(a)

:param a: Matrix<Float64>
:return: (Matrix<Float64>, Matrix<Float64>)

Compute the RQ decomposition of a matrix, :math:`A = RQ` where :math:`R`
is upper triangular and :math:`Q` is unitary/orthogonal.

.. _schur:

.. function:: schur(a)

:param a: Matrix<Float64>
:return: (Matrix<Float64>, Matrix<Float64>)

Compute the Schur decomposition of a matrix, :math:`A = ZTZ^H` where
:math:`Z` is unitary, and :math:`T` is upper triangular in the complex case,
and quasi-triangular in the real case, where there may be 2-by-2 blocks
on the diagonal for complex eigenvalue pairs.

**Note.** Currently only the real case is implemented.