Skip to content

Commit

Permalink
Rename lud to lu. Update docs.
Browse files Browse the repository at this point in the history
  • Loading branch information
Viral B. Shah committed Feb 7, 2013
1 parent b96e9c2 commit c5b0fe0
Show file tree
Hide file tree
Showing 4 changed files with 23 additions and 26 deletions.
3 changes: 1 addition & 2 deletions base/exports.jl
Expand Up @@ -553,8 +553,7 @@ export
ldltd,
linreg,
lu,
lud,
lud!,
lu!,
norm,
normfro,
null,
Expand Down
22 changes: 9 additions & 13 deletions base/linalg_dense.jl
Expand Up @@ -536,17 +536,13 @@ function factors{T}(lu::LUDense{T})
L, U, P
end

function lud!{T<:BlasFloat}(A::Matrix{T})
function lu!{T<:BlasFloat}(A::Matrix{T})
lu, ipiv, info = LAPACK.getrf!(A)
LUDense{T}(lu, ipiv, info)
end

lud{T<:BlasFloat}(A::Matrix{T}) = lud!(copy(A))
lud{T<:Number}(A::Matrix{T}) = lud(float64(A))

## Matlab-compatible
lu{T<:Number}(A::Matrix{T}) = factors(lud(A))
lu(x::Number) = (one(x), x, [1])
lu{T<:BlasFloat}(A::Matrix{T}) = lu!(copy(A))
lu{T<:Number}(A::Matrix{T}) = lu(float64(A))

function det{T}(lu::LUDense{T})
m, n = size(lu)
Expand Down Expand Up @@ -680,7 +676,7 @@ function det(A::Matrix)
m, n = size(A)
if m != n; error("det only defined for square matrices"); end
if istriu(A) | istril(A); return prod(diag(A)); end
return det(lud(A))
return det(lu(A))
end
det(x::Number) = x

Expand Down Expand Up @@ -894,7 +890,7 @@ function cond(A::StridedMatrix, p)
elseif p == 1 || p == Inf
m, n = size(A)
if m != n; error("Use 2-norm for non-square matrices"); end
cnd = 1 / LAPACK.gecon!(p == 1 ? '1' : 'I', lud(A).lu, norm(A, p))
cnd = 1 / LAPACK.gecon!(p == 1 ? '1' : 'I', lu(A).lu, norm(A, p))
else
error("Norm type must be 1, 2 or Inf")
end
Expand Down Expand Up @@ -1190,16 +1186,16 @@ end

#show(io, lu::LUTridiagonal) = print(io, "LU decomposition of ", summary(lu.lu))

lud!{T}(A::Tridiagonal{T}) = LUTridiagonal{T}(LAPACK.gttrf!(A.dl,A.d,A.du)...)
lud{T}(A::Tridiagonal{T}) = LUTridiagonal{T}(LAPACK.gttrf!(copy(A.dl),copy(A.d),copy(A.du))...)
lu(A::Tridiagonal) = factors(lud(A))
lu!{T}(A::Tridiagonal{T}) = LUTridiagonal{T}(LAPACK.gttrf!(A.dl,A.d,A.du)...)
lu{T}(A::Tridiagonal{T}) = LUTridiagonal{T}(LAPACK.gttrf!(copy(A.dl),copy(A.d),copy(A.du))...)
lu(A::Tridiagonal) = factors(lu(A))

function det{T}(lu::LUTridiagonal{T})
n = length(lu.d)
prod(lu.d) * (bool(sum(lu.ipiv .!= 1:n) % 2) ? -one(T) : one(T))
end

det(A::Tridiagonal) = det(lud(A))
det(A::Tridiagonal) = det(lu(A))

(\){T<:BlasFloat}(lu::LUTridiagonal{T}, B::StridedVecOrMat{T}) =
LAPACK.gttrs!('N', lu.dl, lu.d, lu.du, lu.du2, lu.ipiv, copy(B))
Expand Down
16 changes: 10 additions & 6 deletions doc/stdlib/base.rst
Expand Up @@ -1716,21 +1716,25 @@ Linear algebra functions in Julia are largely implemented by calling functions f

Compute the norm of a ``Vector`` or a ``Matrix``

.. function:: lu(A) -> LU
.. function:: lu(A)

Compute LU factorization. LU is an "LU factorization" type that can be used as an ordinary matrix.
Compute the LU factorization of ``A`` and return a ``LUDense`` object. ``factors(lu(A))`` returns the triangular matrices containing the factorization. The following functions are available for ``LUDense`` objects: ``size``, ``factors``, ``\``, ``inv``, ``det``.

.. function:: lu!(A)

``lu!`` is the same as ``lu``, but overwrites the input matrix A with the factorization.

.. function:: chol(A, [LU])

Compute the Cholesky factorization of ``A`` and return a ``CholeskyDense`` object. ``LU`` may be 'L' for using the lower part or 'U' for the upper part. The default is to use 'U'. The ``factors(chol(A))`` returns the triangular matrix containing the factorization. The following functions are available for ``CholeskyDense`` objects: ``size``, ``factors``, ``\``, ``inv``, ``det``.
Compute the Cholesky factorization of ``A`` and return a ``CholeskyDense`` object. ``LU`` may be 'L' for using the lower part or 'U' for the upper part. The default is to use 'U'. ``factors(chol(A))`` returns the triangular matrix containing the factorization. The following functions are available for ``CholeskyDense`` objects: ``size``, ``factors``, ``\``, ``inv``, ``det``. A ``LAPACK.PosDefException`` error is thrown in case the matrix is not positive definite.

.. function:: chol!(A, [LU])

``chol!`` is the same as ``chol``, but overwrites the input matrix A with the factorization.

.. function:: cholpivot(A, [LU])

Compute the pivoted Cholesky factorization of ``A`` and return a ``CholeskyDensePivoted`` object. ``LU`` may be 'L' for using the lower part or 'U' for the upper part. The default is to use 'U'. The ``factors(cholpivot(A))`` returns the triangular matrix containing the factorization. The following functions are available for ``CholeskyDensePivoted`` objects: ``size``, ``factors``, ``\``, ``inv``, ``det``.
Compute the pivoted Cholesky factorization of ``A`` and return a ``CholeskyDensePivoted`` object. ``LU`` may be 'L' for using the lower part or 'U' for the upper part. The default is to use 'U'. ``factors(cholpivot(A))`` returns the triangular matrix containing the factorization. The following functions are available for ``CholeskyDensePivoted`` objects: ``size``, ``factors``, ``\``, ``inv``, ``det``. A ``LAPACK.RankDeficientException`` error is thrown in case the matrix is rank deficient.

.. function:: cholpivot!(A, [LU])

Expand All @@ -1744,9 +1748,9 @@ Linear algebra functions in Julia are largely implemented by calling functions f

Compute QR factorization with pivoting

.. function:: factors(D)
.. function:: factors(F)

Return the factors of a decomposition D. For an LU decomposition, factors(LU) -> L, U, p
Return the factors from a ``Factorization`` object.

.. function:: eig(A) -> D, V

Expand Down
8 changes: 3 additions & 5 deletions test/linalg.jl
Expand Up @@ -32,10 +32,8 @@ for elty in (Float32, Float64, Complex64, Complex128)
@assert_approx_eq inv(bc2) * apd eye(elty, n)
@assert_approx_eq apd * (bc2\b) b

lua = lud(a) # LU decomposition
l,u,p = lu(a)
L,U,P = factors(lua)
@test l == L && u == U && p == P
lua = lu(a) # LU decomposition
l,u,p = factors(lua)
@assert_approx_eq l*u a[p,:]
@assert_approx_eq l[invperm(p),:]*u a
@assert_approx_eq a * inv(lua) eye(elty, n)
Expand Down Expand Up @@ -290,7 +288,7 @@ for elty in (Float32, Float64, Complex64, Complex128)
@assert_approx_eq solve(T,v) invFv
B = convert(Matrix{elty}, B)
@assert_approx_eq solve(T, B) F\B
Tlu = lud(T)
Tlu = lu(T)
x = Tlu\v
@assert_approx_eq x invFv
@assert_approx_eq det(T) det(F)
Expand Down

0 comments on commit c5b0fe0

Please sign in to comment.