Fix logdet(PDiagMat(zeros(0))) (#139)

* fix #138

* Update src/pdiagmat.jl

Co-authored-by: David Widmann <[email protected]>

* make _randPDMat more reliable

Co-authored-by: David Widmann <[email protected]>

Author: st--

src/pdiagmat.jl

@@ -58,7 +58,10 @@ Base.kron(A::PDiagMat, B::PDiagMat) = PDiagMat( vcat([A.diag[i] * B.diag for i i
 ### Algebra
 
 Base.inv(a::PDiagMat) = PDiagMat(a.inv_diag, a.diag)
-LinearAlgebra.logdet(a::PDiagMat) = sum(log, a.diag)
+function LinearAlgebra.logdet(a::PDiagMat)
+    diag = a.diag
+    return isempty(diag) ? zero(log(zero(eltype(diag)))) : sum(log, diag)
+end
 LinearAlgebra.eigmax(a::PDiagMat) = maximum(a.diag)
 LinearAlgebra.eigmin(a::PDiagMat) = minimum(a.diag)
 

test/kron.jl

@@ -1,7 +1,8 @@
 using PDMats
 using Test
+using LinearAlgebra: LinearAlgebra
 
-_randPDMat(T, n) = (X = randn(T, n, n); PDMat(X * X'))
+_randPDMat(T, n) = (X = randn(T, n, n); PDMat(X * X' + LinearAlgebra.I))
 _randPDiagMat(T, n) = PDiagMat(rand(T, n))
 _randScalMat(T, n) = ScalMat(n, rand(T))
 function _pd_compare(A::AbstractPDMat, B::AbstractPDMat)

test/pdmtypes.jl

@@ -38,6 +38,10 @@ x = one(Float32); d = 4
 s = SparseMatrixCSC{Float32}(I, 2, 2)
 @test convert(PDSparseMat{Float64}, PDSparseMat(s)).mat == PDSparseMat(convert(SparseMatrixCSC{Float64}, s)).mat
 
+Z = zeros(0, 0)
+test_pdmat(PDMat(Z), Z; t_eig=false)
+test_pdmat(PDiagMat(diag(Z)), Z; t_eig=false)
+
 # no-op conversion with correct eltype (#101)
 X = PDMat((Y->Y'Y)(randn(Float32, 4, 4)))
 @test convert(AbstractArray{Float32}, X) === X
@@ -51,4 +55,4 @@ a = PDMat([1 0.5; 0.5 1])
 @inferred unwhiten(a, ones(2))
 
 # convert Matrix type to the same Cholesky type (#117)
-@test PDMat([1 0; 0 1]) == [1.0 0.0; 0.0 1.0]
+@test PDMat([1 0; 0 1]) == [1.0 0.0; 0.0 1.0]