Skip to content

Commit

Permalink
Allow vectors in reflectorApply! (JuliaLang#42874)
Browse files Browse the repository at this point in the history
Co-authored-by: Sheehan Olver <[email protected]>
  • Loading branch information
2 people authored and LilithHafner committed Mar 8, 2022
1 parent 380b891 commit fbcbbf3
Show file tree
Hide file tree
Showing 3 changed files with 17 additions and 5 deletions.
4 changes: 2 additions & 2 deletions stdlib/LinearAlgebra/src/generic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1516,9 +1516,9 @@ end
end

# apply reflector from left
@inline function reflectorApply!(x::AbstractVector, τ::Number, A::AbstractMatrix)
@inline function reflectorApply!(x::AbstractVector, τ::Number, A::AbstractVecOrMat)
require_one_based_indexing(x)
m, n = size(A)
m, n = size(A, 1), size(A, 2)
if length(x) != m
throw(DimensionMismatch("reflector has length $(length(x)), which must match the first dimension of matrix A, $m"))
end
Expand Down
6 changes: 3 additions & 3 deletions stdlib/LinearAlgebra/src/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -270,13 +270,13 @@ function qrfactPivotedUnblocked!(A::AbstractMatrix)

# Compute reflector of columns j
x = view(A, j:m, j)
τj = LinearAlgebra.reflector!(x)
τj = reflector!(x)
τ[j] = τj

# Update trailing submatrix with reflector
LinearAlgebra.reflectorApply!(x, τj, view(A, j:m, j+1:n))
reflectorApply!(x, τj, view(A, j:m, j+1:n))
end
return LinearAlgebra.QRPivoted{eltype(A), typeof(A)}(A, τ, piv)
return QRPivoted{eltype(A), typeof(A)}(A, τ, piv)
end

# LAPACK version
Expand Down
12 changes: 12 additions & 0 deletions stdlib/LinearAlgebra/test/generic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -256,6 +256,18 @@ end
@test_throws DimensionMismatch reflect!([x; x], y, c, s)
end

@testset "LinearAlgebra.reflectorApply!" begin
for T in (Float64, ComplexF64)
x = rand(T, 6)
τ = rand(T)
A = rand(T, 6)
B = LinearAlgebra.reflectorApply!(x, τ, copy(A))
C = LinearAlgebra.reflectorApply!(x, τ, reshape(copy(A), (length(A), 1)))
@test B[1] C[1] A[1] - conj(τ)*(A[1] + dot(x[2:end], A[2:end]))
@test B[2:end] C[2:end] A[2:end] - conj(τ)*(A[1] + dot(x[2:end], A[2:end]))*x[2:end]
end
end

@testset "LinearAlgebra.axp(b)y! for element type without commutative multiplication" begin
α = [1 2; 3 4]
β = [5 6; 7 8]
Expand Down

0 comments on commit fbcbbf3

Please sign in to comment.