Skip to content

Commit 5a55cc8

Browse files
lxvmstevengj
andauthored
Fix overflow and type issues with generic det of a triangular matrix (#1418)
Fixes #1415 --------- Co-authored-by: Steven G. Johnson <[email protected]>
1 parent 6c92d28 commit 5a55cc8

File tree

2 files changed

+14
-1
lines changed

2 files changed

+14
-1
lines changed

src/generic.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1846,7 +1846,7 @@ julia> det(BigInt[1 0; 2 2]) # exact integer determinant
18461846
function det(A::AbstractMatrix{T}) where {T}
18471847
if istriu(A) || istril(A)
18481848
S = promote_type(T, typeof((one(T)*zero(T) + zero(T))/one(T)))
1849-
return convert(S, det(UpperTriangular(A)))
1849+
return prod(Base.Fix1(convert, S), @view A[diagind(A)]; init=one(S))
18501850
end
18511851
return det(lu(A; check = false))
18521852
end

test/generic.jl

Lines changed: 13 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,7 @@ using Main.LinearAlgebraTestHelpers.OffsetArrays
1717
using Main.LinearAlgebraTestHelpers.DualNumbers
1818
using Main.LinearAlgebraTestHelpers.FillArrays
1919
using Main.LinearAlgebraTestHelpers.SizedArrays
20+
using Main.LinearAlgebraTestHelpers.Furlongs
2021

2122
Random.seed!(123)
2223

@@ -96,6 +97,18 @@ n = 5 # should be odd
9697
@testset "det with nonstandard Number type" begin
9798
elty <: Real && @test det(Dual.(triu(A), zero(A))) isa Dual
9899
end
100+
if elty <: Int
101+
@testset "det no overflow - triangular" begin
102+
A = diagm([typemax(elty), typemax(elty)])
103+
@test det(A) == det(float(A))
104+
end
105+
end
106+
@testset "det with units - triangular" begin
107+
for dim in 0:4
108+
A = diagm(Furlong.(ones(elty, dim)))
109+
@test det(A) == Furlong{dim}(one(elty))
110+
end
111+
end
99112
end
100113

101114
@testset "diag" begin

0 commit comments

Comments
 (0)