Skip to content

Commit

Permalink
Fix Matrix(QRSparseQ) constructor
Browse files Browse the repository at this point in the history
Fixes #26367
  • Loading branch information
andreasnoack committed Mar 11, 2018
1 parent 6f2c369 commit 05f7f25
Show file tree
Hide file tree
Showing 4 changed files with 21 additions and 14 deletions.
21 changes: 10 additions & 11 deletions stdlib/LinearAlgebra/src/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -527,21 +527,20 @@ AbstractMatrix{T}(Q::QRPackedQ) where {T} = QRPackedQ{T}(Q)
QRCompactWYQ{S}(Q::QRCompactWYQ) where {S} = QRCompactWYQ(convert(AbstractMatrix{S}, Q.factors), convert(AbstractMatrix{S}, Q.T))
AbstractMatrix{S}(Q::QRCompactWYQ{S}) where {S} = Q
AbstractMatrix{S}(Q::QRCompactWYQ) where {S} = QRCompactWYQ{S}(Q)
Matrix(A::AbstractQ{T}) where {T} = lmul!(A, Matrix{T}(I, size(A.factors, 1), min(size(A.factors)...)))
Array(A::AbstractQ) = Matrix(A)
Matrix(Q::Union{QRCompactWYQ{T},QRPackedQ{T}}) where {T} = lmul!(Q, Matrix{T}(I, size(Q, 1), min(size(Q.factors)...)))
Array(Q::AbstractQ) = Matrix(Q)

size(A::Union{QR,QRCompactWY,QRPivoted}, dim::Integer) = size(getfield(A, :factors), dim)
size(A::Union{QR,QRCompactWY,QRPivoted}) = size(getfield(A, :factors))
size(A::AbstractQ, dim::Integer) = 0 < dim ? (dim <= 2 ? size(getfield(A, :factors), 1) : 1) : throw(BoundsError())
size(A::AbstractQ) = size(A, 1), size(A, 2)
size(F::Union{QR,QRCompactWY,QRPivoted}, dim::Integer) = size(getfield(F, :factors), dim)
size(F::Union{QR,QRCompactWY,QRPivoted}) = size(getfield(F, :factors))
size(Q::Union{QRCompactWYQ,QRPackedQ}, dim::Integer) = size(getfield(Q, :factors), dim == 2 ? 1 : dim)
size(Q::AbstractQ) = size(Q, 1), size(Q, 2)


function getindex(A::AbstractQ, i::Integer, j::Integer)
x = zeros(eltype(A), size(A, 1))
function getindex(Q::AbstractQ, i::Integer, j::Integer)
x = zeros(eltype(Q), size(Q, 1))
x[i] = 1
y = zeros(eltype(A), size(A, 2))
y = zeros(eltype(Q), size(Q, 2))
y[j] = 1
return dot(x, lmul!(A, y))
return dot(x, lmul!(Q, y))
end

## Multiplication by Q
Expand Down
4 changes: 2 additions & 2 deletions stdlib/LinearAlgebra/test/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -146,7 +146,7 @@ rectangularQ(Q::LinearAlgebra.AbstractQ) = convert(Array, Q)
@test_throws DimensionMismatch rmul!(Matrix{eltya}(I, n+1, n+1),q)
@test rmul!(squareQ(q), adjoint(q)) Matrix(I, n, n)
@test_throws DimensionMismatch rmul!(Matrix{eltya}(I, n+1, n+1), adjoint(q))
@test_throws BoundsError size(q,-1)
@test_throws ErrorException size(q,-1)
@test_throws DimensionMismatch LinearAlgebra.lmul!(q,zeros(eltya,n1+1))
@test_throws DimensionMismatch LinearAlgebra.lmul!(adjoint(q), zeros(eltya,n1+1))

Expand All @@ -156,7 +156,7 @@ rectangularQ(Q::LinearAlgebra.AbstractQ) = convert(Array, Q)
@test_throws DimensionMismatch rmul!(Matrix{eltya}(I, n+1, n+1),q)
@test rmul!(squareQ(q), adjoint(q)) Matrix(I, n, n)
@test_throws DimensionMismatch rmul!(Matrix{eltya}(I, n+1, n+1),adjoint(q))
@test_throws BoundsError size(q,-1)
@test_throws ErrorException size(q,-1)
@test_throws DimensionMismatch q * Matrix{Int8}(I, n+4, n+4)
end
end
Expand Down
5 changes: 4 additions & 1 deletion stdlib/SuiteSparse/src/spqr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -128,10 +128,13 @@ end
struct QRSparseQ{Tv<:CHOLMOD.VTypes,Ti<:Integer} <: LinearAlgebra.AbstractQ{Tv}
factors::SparseMatrixCSC{Tv,Ti}
τ::Vector{Tv}
n::Int # Number of columns in original matrix
end

Base.size(Q::QRSparseQ) = (size(Q.factors, 1), size(Q.factors, 1))

Matrix(Q::QRSparseQ) where {T} = lmul!(A, Matrix{T}(I, size(Q, 1), min(size(Q, 1), Q.n)))

# From SPQR manual p. 6
_default_tol(A::SparseMatrixCSC) =
20*sum(size(A))*eps(real(eltype(A)))*maximum(norm(view(A, :, i))^2 for i in 1:size(A, 2))
Expand Down Expand Up @@ -303,7 +306,7 @@ julia> F.pcol
"""
@inline function Base.getproperty(F::QRSparse, d::Symbol)
if d == :Q
return QRSparseQ(F.factors, F.τ)
return QRSparseQ(F.factors, F.τ, size(F, 2))
elseif d == :prow
return invperm(F.rpivinv)
elseif d == :pcol
Expand Down
5 changes: 5 additions & 0 deletions stdlib/SuiteSparse/test/spqr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -81,4 +81,9 @@ end
@test A*xs A*xd
end

@testset "Issue 26367" begin
A = sparse([0.0 1 0 0; 0 0 0 0])
@test Matrix(qrfact(A).Q) == Matrix(qrfact(Matrix(A)).Q) == Matrix(I, 2, 2)
end

end

0 comments on commit 05f7f25

Please sign in to comment.