Skip to content

Commit 17026e0

Browse files
committed
Deprecate keyword argument "thin" in favor of "square" in lq methods.
1 parent bee25f0 commit 17026e0

File tree

4 files changed

+33
-23
lines changed

4 files changed

+33
-23
lines changed

NEWS.md

+3
Original file line numberDiff line numberDiff line change
@@ -363,6 +363,9 @@ Deprecated or removed
363363
* The `cholfact`/`cholfact!` methods that accepted an `uplo` symbol have been deprecated
364364
in favor of using `Hermitian` (or `Symmetric`) views ([#22187], [#22188]).
365365

366+
* The `thin` keyword argument for orthogonal decomposition methods has
367+
been deprecated in favor of `square` ([#WHATAREBIRDS]).
368+
366369
* `isposdef(A::AbstractMatrix, UL::Symbol)` and `isposdef!(A::AbstractMatrix, UL::Symbol)`
367370
have been deprecated in favor of `isposdef(Hermitian(A, UL))` and `isposdef!(Hermitian(A, UL))`
368371
respectively ([#22245]).

base/deprecated.jl

+1
Original file line numberDiff line numberDiff line change
@@ -1889,6 +1889,7 @@ end
18891889
# TODO: after 0.7, remove thin keyword argument and associated logic from...
18901890
# (1) base/linalg/svd.jl
18911891
# (2) base/linalg/qr.jl
1892+
# (3) base/linalg/lq.jl
18921893

18931894
@deprecate find(x::Number) find(!iszero, x)
18941895
@deprecate findnext(A, v, i::Integer) findnext(equalto(v), A, i)

base/linalg/lq.jl

+17-11
Original file line numberDiff line numberDiff line change
@@ -33,24 +33,30 @@ lqfact(A::StridedMatrix{<:BlasFloat}) = lqfact!(copy(A))
3333
lqfact(x::Number) = lqfact(fill(x,1,1))
3434

3535
"""
36-
lq(A; [thin=true]) -> L, Q
36+
lq(A; square = false) -> L, Q
3737
38-
Perform an LQ factorization of `A` such that `A = L*Q`. The default is to compute
39-
a "thin" factorization. The LQ factorization is the QR factorization of `A.'`.
40-
`L` is not extendedwith zeros if the explicit, square form of `Q`
41-
is requested via `thin = false`.
38+
Perform an LQ factorization of `A` such that `A = L*Q`. The default (`square = false`)
39+
computes a factorization with possibly-rectangular `L` and `Q`, commonly the "thin"
40+
factorization. The LQ factorization is the QR factorization of `A.'`. If the explicit,
41+
square form of `Q` is requested via `square = true`, `L` is not extended with zeros.
4242
4343
!!! note
4444
While in QR factorization the "thin" factorization is so named due to yielding
45-
either a square or "tall"/"thin" factor `Q`, in LQ factorization the "thin"
46-
factorization somewhat confusingly produces either a square or "short"/"wide"
47-
factor `Q`. "Thin" factorizations more broadly are also (more descriptively)
48-
referred to as "truncated" or "reduced" factorizatons.
45+
either a square or "tall"/"thin" rectangular factor `Q`, in LQ factorization the
46+
"thin" factorization somewhat confusingly produces either a square or "short"/"wide"
47+
rectangular factor `Q`. "Thin" factorizations more broadly are also
48+
referred to as "reduced" factorizatons.
4949
"""
50-
function lq(A::Union{Number,AbstractMatrix}; thin::Bool = true)
50+
function lq(A::Union{Number,AbstractMatrix}; square::Bool = false, thin::Union{Bool,Void} = nothing)
51+
# DEPRECATION TODO: remove deprecated thin argument and associated logic after 0.7
52+
if thin != nothing
53+
Base.depwarn(string("the `thin` keyword argument in `lq(A; thin = $(thin))` has ",
54+
"been deprecated in favor of `square`, i.e. `lq(A; square = $(!thin))`."), :lq)
55+
square::Bool = !thin
56+
end
5157
F = lqfact(A)
5258
L, Q = F[:L], F[:Q]
53-
return L, thin ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 2)))
59+
return L, !square ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 2)))
5460
end
5561

5662
copy(A::LQ) = LQ(copy(A.factors), copy(A.τ))

test/linalg/lq.jl

+12-12
Original file line numberDiff line numberDiff line change
@@ -111,35 +111,35 @@ end
111111
@testset "correct form of Q from lq(...) (#23729)" begin
112112
# where the original matrix (say A) is square or has more rows than columns,
113113
# then A's factorization's triangular factor (say L) should have the same shape
114-
# as A independent of factorization form (truncated, square), and A's factorization's
114+
# as A independent of factorization form (square, rectangular/"thin"), and A's factorization's
115115
# orthogonal factor (say Q) should be a square matrix of order of A's number of
116-
# columns independent of factorization form (truncated, square), and L and Q
116+
# columns independent of factorization form (square, rectangular/"thin"), and L and Q
117117
# should have multiplication-compatible shapes.
118118
m, n = 4, 2
119119
A = randn(m, n)
120-
for thin in (true, false)
121-
L, Q = lq(A, thin = thin)
120+
for square in (false, true)
121+
L, Q = lq(A, square = square)
122122
@test size(L) == (m, n)
123123
@test size(Q) == (n, n)
124124
@test isapprox(A, L*Q)
125125
end
126126
# where the original matrix has strictly fewer rows than columns ...
127127
m, n = 2, 4
128128
A = randn(m, n)
129-
# ... then, for a truncated factorization of A, L should be a square matrix
129+
# ... then, for a rectangular/"thin" factorization of A, L should be a square matrix
130130
# of order of A's number of rows, Q should have the same shape as A,
131131
# and L and Q should have multiplication-compatible shapes
132-
Lthin, Qthin = lq(A, thin = true)
133-
@test size(Lthin) == (m, m)
134-
@test size(Qthin) == (m, n)
135-
@test isapprox(A, Lthin * Qthin)
136-
# ... and, for a square/non-truncated factorization of A, L should have the
132+
Lrect, Qrect = lq(A, square = false)
133+
@test size(Lrect) == (m, m)
134+
@test size(Qrect) == (m, n)
135+
@test isapprox(A, Lrect * Qrect)
136+
# ... and, for a square factorization of A, L should have the
137137
# same shape as A, Q should be a square matrix of order of A's number of columns,
138138
# and L and Q should have multiplication-compatible shape. but instead the L returned
139-
# has no zero-padding on the right / is L for the truncated factorization,
139+
# has no zero-padding on the right / is L for the rectangular/"thin" factorization,
140140
# so for L and Q to have multiplication-compatible shapes, L must be zero-padded
141141
# to have the shape of A.
142-
Lsquare, Qsquare = lq(A, thin = false)
142+
Lsquare, Qsquare = lq(A, square = true)
143143
@test size(Lsquare) == (m, m)
144144
@test size(Qsquare) == (n, n)
145145
@test isapprox(A, [Lsquare zeros(m, n - m)] * Qsquare)

0 commit comments

Comments
 (0)