Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

deprecate orthogonal decomposition keyword "thin" to "full" #24279

Merged
merged 3 commits into from
Oct 28, 2017
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 4 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
@@ -381,6 +381,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 `full`, which has the opposite meaning:
`thin == true` if and only if `full == 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]).
5 changes: 5 additions & 0 deletions base/deprecated.jl
Original file line number Diff line number Diff line change
@@ -2031,6 +2031,11 @@ end
@deprecate strwidth textwidth
@deprecate charwidth textwidth

# 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)
@deprecate findfirst(A, v) findfirst(equalto(v), A)
20 changes: 18 additions & 2 deletions base/linalg/bidiag.jl
Original file line number Diff line number Diff line change
@@ -182,11 +182,27 @@ similar(B::Bidiagonal, ::Type{T}, dims::Union{Dims{1},Dims{2}}) where {T} = spze

#Singular values
svdvals!(M::Bidiagonal{<:BlasReal}) = LAPACK.bdsdc!(M.uplo, 'N', M.dv, M.ev)[1]
function svdfact!(M::Bidiagonal{<:BlasReal}; thin::Bool=true)
function svdfact!(M::Bidiagonal{<:BlasReal}; full::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 `svdfact!(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g. `svdfact!(A; full = $(!thin))`."), :svdfact!)
full::Bool = !thin
end
d, e, U, Vt, Q, iQ = LAPACK.bdsdc!(M.uplo, 'I', M.dv, M.ev)
SVD(U, d, Vt)
end
svdfact(M::Bidiagonal; thin::Bool=true) = svdfact!(copy(M),thin=thin)
function svdfact(M::Bidiagonal; full::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 `svdfact(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g. `svdfact(A; full = $(!thin))`."), :svdfact)
full::Bool = !thin
end
return svdfact!(copy(M), full = full)
end

####################
# Generic routines #
4 changes: 2 additions & 2 deletions base/linalg/dense.jl
Original file line number Diff line number Diff line change
@@ -1263,7 +1263,7 @@ function pinv(A::StridedMatrix{T}, tol::Real) where T
return B
end
end
SVD = svdfact(A, thin=true)
SVD = svdfact(A, full = false)
Stype = eltype(SVD.S)
Sinv = zeros(Stype, length(SVD.S))
index = SVD.S .> tol*maximum(SVD.S)
@@ -1305,7 +1305,7 @@ julia> nullspace(M)
function nullspace(A::StridedMatrix{T}) where T
m, n = size(A)
(m == 0 || n == 0) && return eye(T, n)
SVD = svdfact(A, thin = false)
SVD = svdfact(A, full = true)
indstart = sum(SVD.S .> max(m,n)*maximum(SVD.S)*eps(eltype(SVD.S))) + 1
return SVD.Vt[indstart:end,:]'
end
29 changes: 18 additions & 11 deletions base/linalg/lq.jl
Original file line number Diff line number Diff line change
@@ -33,24 +33,31 @@ lqfact(A::StridedMatrix{<:BlasFloat}) = lqfact!(copy(A))
lqfact(x::Number) = lqfact(fill(x,1,1))

"""
lq(A; [thin=true]) -> L, Q
lq(A; full = 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 (`full = 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,
full/square form of `Q` is requested via `full = 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}; full::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 `full`, which has the opposite meaning, ",
"e.g. `lq(A; full = $(!thin))`."), :lq)
full::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, !full ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 2)))
end

copy(A::LQ) = LQ(copy(A.factors), copy(A.τ))
29 changes: 19 additions & 10 deletions base/linalg/qr.jl
Original file line number Diff line number Diff line change
@@ -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 full/square or non-full/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
[`Matrix`](@ref).
@@ -305,24 +305,33 @@ 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); full::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 full/square orthogonal factor `Q` is requested (via `full = true`).
"""
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);
full::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 `full`, which has the opposite meaning, ",
"e.g. `qr(A, pivot; full = $(!thin))`."), :qr)
full::Bool = !thin
end
return _qr(A, pivot, full = full)
end
function _qr(A::Union{Number,AbstractMatrix}, ::Val{false}; full::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 (!full ? 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}; full::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 (!full ? Array(Q) : A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 1)))), R, p
end

"""
86 changes: 67 additions & 19 deletions base/linalg/svd.jl
Original file line number Diff line number Diff line change
@@ -11,23 +11,30 @@ end
SVD(U::AbstractArray{T}, S::Vector{Tr}, Vt::AbstractArray{T}) where {T,Tr} = SVD{T,Tr,typeof(U)}(U, S, Vt)

"""
svdfact!(A, thin::Bool=true) -> SVD
svdfact!(A; full::Bool = false) -> SVD
`svdfact!` is the same as [`svdfact`](@ref), but saves space by
overwriting the input `A`, instead of creating a copy.
"""
function svdfact!(A::StridedMatrix{T}; thin::Bool=true) where T<:BlasFloat
function svdfact!(A::StridedMatrix{T}; full::Bool = false, thin::Union{Bool,Void} = nothing) where T<:BlasFloat
# DEPRECATION TODO: remove deprecated thin argument and associated logic after 0.7
if thin != nothing
Base.depwarn(string("the `thin` keyword argument in `svdfact!(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g. `svdfact!(A; full = $(!thin))`."), :svdfact!)
full::Bool = !thin
end
m,n = size(A)
if m == 0 || n == 0
u,s,vt = (eye(T, m, thin ? n : m), real(zeros(T,0)), eye(T,n,n))
u,s,vt = (eye(T, m, full ? m : n), real(zeros(T,0)), eye(T,n,n))
else
u,s,vt = LAPACK.gesdd!(thin ? 'S' : 'A', A)
u,s,vt = LAPACK.gesdd!(full ? 'A' : 'S', A)
end
SVD(u,s,vt)
end

"""
svdfact(A; thin::Bool=true) -> SVD
svdfact(A; full::Bool = false) -> SVD
Compute the singular value decomposition (SVD) of `A` and return an `SVD` object.
@@ -36,9 +43,9 @@ Compute the singular value decomposition (SVD) of `A` and return an `SVD` object
The algorithm produces `Vt` and hence `Vt` is more efficient to extract than `V`.
The singular values in `S` are sorted in descending order.
If `thin=true` (default), a thin SVD is returned. For a ``M \\times N`` matrix
`A`, `U` is ``M \\times M`` for a full SVD (`thin=false`) and
``M \\times \\min(M, N)`` for a thin SVD.
If `full = false` (default), a "thin" SVD is returned. For a
``M \\times N`` matrix `A`, `U` is ``M \\times M`` for a "full" SVD (`full = true`) and
``M \\times \\min(M, N)`` for a "thin" SVD.
# Examples
```jldoctest
@@ -59,22 +66,47 @@ julia> F[:U] * Diagonal(F[:S]) * F[:Vt]
0.0 2.0 0.0 0.0 0.0
```
"""
function svdfact(A::StridedVecOrMat{T}; thin::Bool = true) where T
function svdfact(A::StridedVecOrMat{T}; full::Bool = false, thin::Union{Bool,Void} = nothing) where T
# DEPRECATION TODO: remove deprecated thin argument and associated logic after 0.7
if thin != nothing
Base.depwarn(string("the `thin` keyword argument in `svdfact(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g. `svdfact(A; full = $(!thin))`."), :svdfact)
full::Bool = !thin
end
S = promote_type(Float32, typeof(one(T)/norm(one(T))))
svdfact!(copy_oftype(A, S), thin = thin)
svdfact!(copy_oftype(A, S), full = full)
end
function svdfact(x::Number; full::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 `svdfact(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g. `svdfact(A; full = $(!thin))`."), :svdfact)
full::Bool = !thin
end
return SVD(x == 0 ? fill(one(x), 1, 1) : fill(x/abs(x), 1, 1), [abs(x)], fill(one(x), 1, 1))
end
function svdfact(x::Integer; full::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 `svdfact(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g. `svdfact(A; full = $(!thin))`."), :svdfact)
full::Bool = !thin
end
return svdfact(float(x), full = full)
end
svdfact(x::Number; thin::Bool=true) = SVD(x == 0 ? fill(one(x), 1, 1) : fill(x/abs(x), 1, 1), [abs(x)], fill(one(x), 1, 1))
svdfact(x::Integer; thin::Bool=true) = svdfact(float(x), thin=thin)

"""
svd(A; thin::Bool=true) -> U, S, V
svd(A; full::Bool = false) -> U, S, V
Computes the SVD of `A`, returning `U`, vector `S`, and `V` such that
`A == U * Diagonal(S) * V'`. The singular values in `S` are sorted in descending order.
If `thin = true` (default), a thin SVD is returned. For a ``M \\times N`` matrix
`A`, `U` is ``M \\times M`` for a full SVD (`thin=false`) and
``M \\times \\min(M, N)`` for a thin SVD.
If `full = false` (default), a "thin" SVD is returned. For a ``M \\times N`` matrix
`A`, `U` is ``M \\times M`` for a "full" SVD (`full = true`) and
``M \\times \\min(M, N)`` for a "thin" SVD.
`svd` is a wrapper around [`svdfact`](@ref), extracting all parts
of the `SVD` factorization to a tuple. Direct use of `svdfact` is therefore more
@@ -99,11 +131,27 @@ julia> U * Diagonal(S) * V'
0.0 2.0 0.0 0.0 0.0
```
"""
function svd(A::AbstractArray; thin::Bool=true)
F = svdfact(A, thin=thin)
function svd(A::AbstractArray; full::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 `svd(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g `svd(A; full = $(!thin))`."), :svd)
full::Bool = !thin
end
F = svdfact(A, full = full)
F.U, F.S, F.Vt'
end
svd(x::Number; thin::Bool=true) = first.(svd(fill(x, 1, 1)))
function svd(x::Number; full::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 `svd(A; thin = $(thin))` has ",
"been deprecated in favor of `full`, which has the opposite meaning, ",
"e.g. `svd(A; full = $(!thin))`."), :svd)
full::Bool = !thin
end
return first.(svd(fill(x, 1, 1)))
end

function getindex(F::SVD, d::Symbol)
if d == :U
54 changes: 27 additions & 27 deletions test/linalg/lq.jl
Original file line number Diff line number Diff line change
@@ -22,7 +22,7 @@ bimg = randn(n,2)/2

# helper functions to unambiguously recover explicit forms of an LQPackedQ
squareQ(Q::LinAlg.LQPackedQ) = A_mul_B!(Q, eye(eltype(Q), size(Q.factors, 2)))
truncatedQ(Q::LinAlg.LQPackedQ) = convert(Array, Q)
rectangularQ(Q::LinAlg.LQPackedQ) = convert(Array, Q)

@testset for eltya in (Float32, Float64, Complex64, Complex128)
a = eltya == Int ? rand(1:7, n, n) : convert(Matrix{eltya}, eltya <: Complex ? complex.(areal, aimg) : areal)
@@ -98,10 +98,10 @@ truncatedQ(Q::LinAlg.LQPackedQ) = convert(Array, Q)
@testset "Matmul with LQ factorizations" begin
lqa = lqfact(a[:,1:n1])
l,q = lqa[:L], lqa[:Q]
@test truncatedQ(q)*truncatedQ(q)' eye(eltya,n1)
@test rectangularQ(q)*rectangularQ(q)' eye(eltya,n1)
@test squareQ(q)'*squareQ(q) eye(eltya, n1)
@test_throws DimensionMismatch A_mul_B!(eye(eltya,n+1),q)
@test Ac_mul_B!(q, truncatedQ(q)) eye(eltya,n1)
@test Ac_mul_B!(q, rectangularQ(q)) eye(eltya,n1)
@test_throws DimensionMismatch A_mul_Bc!(eye(eltya,n+1),q)
@test_throws BoundsError size(q,-1)
end
@@ -111,38 +111,38 @@ 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 ("full", "reduced"/"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 ("full", "reduced"/"thin"), and L and Q
# should have multiplication-compatible shapes.
local m, n = 4, 2
A = randn(m, n)
for thin in (true, false)
L, Q = lq(A, thin = thin)
for full in (false, true)
L, Q = lq(A, full = full)
@test size(L) == (m, n)
@test size(Q) == (n, n)
@test isapprox(A, L*Q)
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, full = false)
@test size(Lrect) == (m, m)
@test size(Qrect) == (m, n)
@test isapprox(A, Lrect * Qrect)
# ... and, for a full 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)
@test size(Lsquare) == (m, m)
@test size(Qsquare) == (n, n)
@test isapprox(A, [Lsquare zeros(m, n - m)] * Qsquare)
Lfull, Qfull = lq(A, full = true)
@test size(Lfull) == (m, m)
@test size(Qfull) == (n, n)
@test isapprox(A, [Lfull zeros(m, n - m)] * Qfull)
end

@testset "getindex on LQPackedQ (#23733)" begin
@@ -153,14 +153,14 @@ end
return implicitQ, explicitQ
end

m, n = 3, 3 # truncated Q 3-by-3, square Q 3-by-3
m, n = 3, 3 # reduced Q 3-by-3, full Q 3-by-3
implicitQ, explicitQ = getqs(lqfact(randn(m, n)))
@test implicitQ[1, 1] == explicitQ[1, 1]
@test implicitQ[m, 1] == explicitQ[m, 1]
@test implicitQ[1, n] == explicitQ[1, n]
@test implicitQ[m, n] == explicitQ[m, n]

m, n = 3, 4 # truncated Q 3-by-4, square Q 4-by-4
m, n = 3, 4 # reduced Q 3-by-4, full Q 4-by-4
implicitQ, explicitQ = getqs(lqfact(randn(m, n)))
@test implicitQ[1, 1] == explicitQ[1, 1]
@test implicitQ[m, 1] == explicitQ[m, 1]
@@ -169,7 +169,7 @@ end
@test implicitQ[m+1, 1] == explicitQ[m+1, 1]
@test implicitQ[m+1, n] == explicitQ[m+1, n]

m, n = 4, 3 # truncated Q 3-by-3, square Q 3-by-3
m, n = 4, 3 # reduced Q 3-by-3, full Q 3-by-3
implicitQ, explicitQ = getqs(lqfact(randn(m, n)))
@test implicitQ[1, 1] == explicitQ[1, 1]
@test implicitQ[n, 1] == explicitQ[n, 1]
@@ -178,11 +178,11 @@ end
end

@testset "size on LQPackedQ (#23780)" begin
# size(Q::LQPackedQ) yields the shape of Q's square form
# size(Q::LQPackedQ) yields the shape of Q's full/square form
for ((mA, nA), nQ) in (
((3, 3), 3), # A 3-by-3 => square Q 3-by-3
((3, 4), 4), # A 3-by-4 => square Q 4-by-4
((4, 3), 3) )# A 4-by-3 => square Q 3-by-3
((3, 3), 3), # A 3-by-3 => full/square Q 3-by-3
((3, 4), 4), # A 3-by-4 => full/square Q 4-by-4
((4, 3), 3) )# A 4-by-3 => full/square Q 3-by-3
@test size(lqfact(randn(mA, nA))[:Q]) == (nQ, nQ)
end
end
@@ -205,12 +205,12 @@ end
@test Ac_mul_Bc(C, implicitQ) Ac_mul_Bc(C, explicitQ)
end
# where the LQ-factored matrix has at least as many rows m as columns n,
# Q's square and truncated forms have the same shape (n-by-n). hence we expect
# Q's full/square and reduced/rectangular forms have the same shape (n-by-n). hence we expect
# _only_ *-by-n (n-by-*) C to work in A_mul_B*(C, Q) (Ac_mul_B*(C, Q)) ops.
# and hence the n-by-n C tests above suffice.
#
# where the LQ-factored matrix has more columns n than rows m,
# Q's square form is n-by-n whereas its truncated form is m-by-n.
# Q's full/square form is n-by-n whereas its reduced/rectangular form is m-by-n.
# hence we need also test *-by-m C with
# A*_mul_B(C, Q) ops, as below via m-by-m C.
mA, nA = 3, 4
8 changes: 4 additions & 4 deletions test/linalg/qr.jl
Original file line number Diff line number Diff line change
@@ -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)
rectangularQ(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)
@@ -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'*rectangularQ(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] rectangularQ(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'
@@ -163,7 +163,7 @@ end

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