From b43fd703d8f41ca5258859ffd837a40a231ccd95 Mon Sep 17 00:00:00 2001 From: Sacha Verweij Date: Sun, 22 Oct 2017 15:41:43 -0700 Subject: [PATCH] Deprecate keyword argument "thin" in favor of "square" in lq methods. --- NEWS.md | 4 ++++ base/deprecated.jl | 1 + base/linalg/lq.jl | 28 +++++++++++++++++----------- test/linalg/lq.jl | 24 ++++++++++++------------ 4 files changed, 34 insertions(+), 23 deletions(-) diff --git a/NEWS.md b/NEWS.md index 899ce39c96cdd..b998273a52cfe 100644 --- a/NEWS.md +++ b/NEWS.md @@ -363,6 +363,10 @@ Deprecated or removed * The `cholfact`/`cholfact!` methods that accepted an `uplo` symbol have been deprecated in favor of using `Hermitian` (or `Symmetric`) views ([#22187], [#22188]). + * The `thin` keyword argument for orthogonal decomposition methods has + been deprecated in favor of `square`, which has the opposite meaning: + `thin == true` if and only if `square == false` ([#24279]). + * `isposdef(A::AbstractMatrix, UL::Symbol)` and `isposdef!(A::AbstractMatrix, UL::Symbol)` have been deprecated in favor of `isposdef(Hermitian(A, UL))` and `isposdef!(Hermitian(A, UL))` respectively ([#22245]). diff --git a/base/deprecated.jl b/base/deprecated.jl index a218ea172224f..2153f0e948284 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -1889,6 +1889,7 @@ end # TODO: after 0.7, remove thin keyword argument and associated logic from... # (1) base/linalg/svd.jl # (2) base/linalg/qr.jl +# (3) base/linalg/lq.jl @deprecate find(x::Number) find(!iszero, x) @deprecate findnext(A, v, i::Integer) findnext(equalto(v), A, i) diff --git a/base/linalg/lq.jl b/base/linalg/lq.jl index 356b3e99b9228..4c36ad41bbca7 100644 --- a/base/linalg/lq.jl +++ b/base/linalg/lq.jl @@ -33,24 +33,30 @@ lqfact(A::StridedMatrix{<:BlasFloat}) = lqfact!(copy(A)) lqfact(x::Number) = lqfact(fill(x,1,1)) """ - lq(A; [thin=true]) -> L, Q + lq(A; square = false) -> L, Q -Perform an LQ factorization of `A` such that `A = L*Q`. The default is to compute -a "thin" factorization. The LQ factorization is the QR factorization of `A.'`. -`L` is not extendedwith zeros if the explicit, square form of `Q` -is requested via `thin = false`. +Perform an LQ factorization of `A` such that `A = L*Q`. The default (`square = false`) +computes a factorization with possibly-rectangular `L` and `Q`, commonly the "thin" +factorization. The LQ factorization is the QR factorization of `A.'`. If the explicit, +square form of `Q` is requested via `square = true`, `L` is not extended with zeros. !!! note While in QR factorization the "thin" factorization is so named due to yielding - either a square or "tall"/"thin" factor `Q`, in LQ factorization the "thin" - factorization somewhat confusingly produces either a square or "short"/"wide" - factor `Q`. "Thin" factorizations more broadly are also (more descriptively) - referred to as "truncated" or "reduced" factorizatons. + either a square or "tall"/"thin" rectangular factor `Q`, in LQ factorization the + "thin" factorization somewhat confusingly produces either a square or "short"/"wide" + rectangular factor `Q`. "Thin" factorizations more broadly are also + referred to as "reduced" factorizatons. """ -function lq(A::Union{Number,AbstractMatrix}; thin::Bool = true) +function lq(A::Union{Number,AbstractMatrix}; square::Bool = false, thin::Union{Bool,Void} = nothing) + # DEPRECATION TODO: remove deprecated thin argument and associated logic after 0.7 + if thin != nothing + Base.depwarn(string("the `thin` keyword argument in `lq(A; thin = $(thin))` has ", + "been deprecated in favor of `square`, i.e. `lq(A; square = $(!thin))`."), :lq) + square::Bool = !thin + end F = lqfact(A) L, Q = F[:L], F[:Q] - return L, thin ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 2))) + return L, !square ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 2))) end copy(A::LQ) = LQ(copy(A.factors), copy(A.τ)) diff --git a/test/linalg/lq.jl b/test/linalg/lq.jl index 2ae6c634a7f2f..d3c7f19939525 100644 --- a/test/linalg/lq.jl +++ b/test/linalg/lq.jl @@ -111,14 +111,14 @@ end @testset "correct form of Q from lq(...) (#23729)" begin # where the original matrix (say A) is square or has more rows than columns, # then A's factorization's triangular factor (say L) should have the same shape - # as A independent of factorization form (truncated, square), and A's factorization's + # as A independent of factorization form (square, rectangular/"thin"), and A's factorization's # orthogonal factor (say Q) should be a square matrix of order of A's number of - # columns independent of factorization form (truncated, square), and L and Q + # columns independent of factorization form (square, rectangular/"thin"), and L and Q # should have multiplication-compatible shapes. m, n = 4, 2 A = randn(m, n) - for thin in (true, false) - L, Q = lq(A, thin = thin) + for square in (false, true) + L, Q = lq(A, square = square) @test size(L) == (m, n) @test size(Q) == (n, n) @test isapprox(A, L*Q) @@ -126,20 +126,20 @@ end # where the original matrix has strictly fewer rows than columns ... m, n = 2, 4 A = randn(m, n) - # ... then, for a truncated factorization of A, L should be a square matrix + # ... then, for a rectangular/"thin" factorization of A, L should be a square matrix # of order of A's number of rows, Q should have the same shape as A, # and L and Q should have multiplication-compatible shapes - Lthin, Qthin = lq(A, thin = true) - @test size(Lthin) == (m, m) - @test size(Qthin) == (m, n) - @test isapprox(A, Lthin * Qthin) - # ... and, for a square/non-truncated factorization of A, L should have the + Lrect, Qrect = lq(A, square = false) + @test size(Lrect) == (m, m) + @test size(Qrect) == (m, n) + @test isapprox(A, Lrect * Qrect) + # ... and, for a square factorization of A, L should have the # same shape as A, Q should be a square matrix of order of A's number of columns, # and L and Q should have multiplication-compatible shape. but instead the L returned - # has no zero-padding on the right / is L for the truncated factorization, + # has no zero-padding on the right / is L for the rectangular/"thin" factorization, # so for L and Q to have multiplication-compatible shapes, L must be zero-padded # to have the shape of A. - Lsquare, Qsquare = lq(A, thin = false) + Lsquare, Qsquare = lq(A, square = true) @test size(Lsquare) == (m, m) @test size(Qsquare) == (n, n) @test isapprox(A, [Lsquare zeros(m, n - m)] * Qsquare)