Skip to content

Commit

Permalink
Deprecate keyword argument "thin" in favor of "square" in qr methods.
Browse files Browse the repository at this point in the history
  • Loading branch information
Sacha0 committed Oct 23, 2017
1 parent b46e0ba commit bee25f0
Show file tree
Hide file tree
Showing 3 changed files with 23 additions and 14 deletions.
1 change: 1 addition & 0 deletions base/deprecated.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1888,6 +1888,7 @@ end

# TODO: after 0.7, remove thin keyword argument and associated logic from...
# (1) base/linalg/svd.jl
# (2) base/linalg/qr.jl

@deprecate find(x::Number) find(!iszero, x)
@deprecate findnext(A, v, i::Integer) findnext(equalto(v), A, i)
Expand Down
28 changes: 18 additions & 10 deletions base/linalg/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -265,7 +265,7 @@ The following functions are available for the `QR` objects: [`inv`](@ref), [`siz
and [`\\`](@ref). When `A` is rectangular, `\\` will return a least squares
solution and if the solution is not unique, the one with smallest norm is returned.
Multiplication with respect to either thin or full `Q` is allowed, i.e. both `F[:Q]*F[:R]`
Multiplication with respect to either square or non-square `Q` is allowed, i.e. both `F[:Q]*F[:R]`
and `F[:Q]*A` are supported. A `Q` matrix can be converted into a regular matrix with
[`full`](@ref) which has a named argument `thin`.
Expand Down Expand Up @@ -305,24 +305,32 @@ end
qrfact(x::Number) = qrfact(fill(x,1,1))

"""
qr(A, pivot=Val(false); thin::Bool=true) -> Q, R, [p]
qr(A, pivot=Val(false); square::Bool = false) -> Q, R, [p]
Compute the (pivoted) QR factorization of `A` such that either `A = Q*R` or `A[:,p] = Q*R`.
Also see [`qrfact`](@ref).
The default is to compute a thin factorization. Note that `R` is not
extended with zeros when the full `Q` is requested.
The default is to compute a "thin" factorization. Note that `R` is not
extended with zeros when a square orthogonal factor `Q` is requested.
"""
qr(A::Union{Number, AbstractMatrix}, pivot::Union{Val{false}, Val{true}}=Val(false); thin::Bool=true) =
_qr(A, pivot, thin=thin)
function _qr(A::Union{Number, AbstractMatrix}, ::Val{false}; thin::Bool=true)
function qr(A::Union{Number,AbstractMatrix}, pivot::Union{Val{false},Val{true}} = Val(false);
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 `qr(A, pivot; thin = $(thin))` has ",
"been deprecated in favor of `square`, i.e. `qr(A, pivot; square = $(!thin))`."), :qr)
square::Bool = !thin
end
return _qr(A, pivot, square = square)
end
function _qr(A::Union{Number,AbstractMatrix}, ::Val{false}; square::Bool = false)
F = qrfact(A, Val(false))
Q, R = getq(F), F[:R]::Matrix{eltype(F)}
return (thin ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 1)))), R
return (!square ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 1)))), R
end
function _qr(A::Union{Number, AbstractMatrix}, ::Val{true}; thin::Bool=true)
function _qr(A::Union{Number, AbstractMatrix}, ::Val{true}; square::Bool = false)
F = qrfact(A, Val(true))
Q, R, p = getq(F), F[:R]::Matrix{eltype(F)}, F[:p]::Vector{BlasInt}
return (thin ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 1)))), R, p
return (!square ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 1)))), R, p
end

"""
Expand Down
8 changes: 4 additions & 4 deletions test/linalg/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@ bimg = randn(n,2)/2

# helper functions to unambiguously recover explicit forms of an implicit QR Q
squareQ(Q::LinAlg.AbstractQ) = A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 1)))
truncatedQ(Q::LinAlg.AbstractQ) = convert(Array, Q)
nonsquareQ(Q::LinAlg.AbstractQ) = convert(Array, Q)

@testset for eltya in (Float32, Float64, Complex64, Complex128, BigFloat, Int)
raw_a = eltya == Int ? rand(1:7, n, n) : convert(Matrix{eltya}, eltya <: Complex ? complex.(areal, aimg) : areal)
Expand Down Expand Up @@ -75,9 +75,9 @@ truncatedQ(Q::LinAlg.AbstractQ) = convert(Array, Q)
q,r = qra[:Q], qra[:R]
@test_throws KeyError qra[:Z]
@test q'*squareQ(q) eye(a_1)
@test q'*truncatedQ(q) eye(a_1, n1)
@test q'*nonsquareQ(q) eye(a_1, n1)
@test q*r a[:, 1:n1]
@test q*b[1:n1] truncatedQ(q)*b[1:n1] atol=100ε
@test q*b[1:n1] nonsquareQ(q)*b[1:n1] atol=100ε
@test q*b squareQ(q)*b atol=100ε
@test_throws DimensionMismatch q*b[1:n1 + 1]
@test_throws DimensionMismatch b[1:n1 + 1]*q'
Expand Down Expand Up @@ -163,7 +163,7 @@ end

@testset "Issue 7304" begin
A = [-√.5 -√.5; -√.5 .5]
Q = truncatedQ(qrfact(A)[:Q])
Q = nonsquareQ(qrfact(A)[:Q])
@test vecnorm(A-Q) < eps()
end

Expand Down

0 comments on commit bee25f0

Please sign in to comment.