Skip to content

Commit

Permalink
Allow rank computation for QRPivoted matrices (#54283)
Browse files Browse the repository at this point in the history
  • Loading branch information
danilonumeroso authored Apr 30, 2024
1 parent cc26f4a commit 831ebe0
Show file tree
Hide file tree
Showing 3 changed files with 21 additions and 0 deletions.
2 changes: 2 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -87,6 +87,8 @@ Standard library changes

#### LinearAlgebra

* `rank` can now take a `QRPivoted` matrix to allow rank estimation via QR factorization ([#54283]).

#### Logging

#### Printf
Expand Down
7 changes: 7 additions & 0 deletions stdlib/LinearAlgebra/src/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -535,6 +535,13 @@ function ldiv!(A::QRCompactWY{T}, B::AbstractMatrix{T}) where {T}
return B
end

function rank(A::QRPivoted; atol::Real=0, rtol::Real=min(size(A)...) * eps(real(float(one(eltype(A.Q))))) * iszero(atol))
m = min(size(A)...)
m == 0 && return 0
tol = max(atol, rtol*abs(A.R[1,1]))
return something(findfirst(i -> abs(A.R[i,i]) <= tol, 1:m), m+1) - 1
end

# Julia implementation similar to xgelsy
function ldiv!(A::QRPivoted{T,<:StridedMatrix}, B::AbstractMatrix{T}, rcond::Real) where {T<:BlasFloat}
require_one_based_indexing(B)
Expand Down
12 changes: 12 additions & 0 deletions stdlib/LinearAlgebra/test/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -528,4 +528,16 @@ end
end
end

@testset "issue #53214" begin
# Test that the rank of a QRPivoted matrix is computed correctly
@test rank(qr([1.0 0.0; 0.0 1.0], ColumnNorm())) == 2
@test rank(qr([1.0 0.0; 0.0 0.9], ColumnNorm()), rtol=0.95) == 1
@test rank(qr([1.0 0.0; 0.0 0.9], ColumnNorm()), atol=0.95) == 1
@test rank(qr([1.0 0.0; 0.0 1.0], ColumnNorm()), rtol=1.01) == 0
@test rank(qr([1.0 0.0; 0.0 1.0], ColumnNorm()), atol=1.01) == 0

@test rank(qr([1.0 2.0; 2.0 4.0], ColumnNorm())) == 1
@test rank(qr([1.0 2.0 3.0; 4.0 5.0 6.0 ; 7.0 8.0 9.0], ColumnNorm())) == 2
end

end # module TestQR

0 comments on commit 831ebe0

Please sign in to comment.