cov of diagonal normal dist should be diagonal

[oschulz]

Currently we have

julia> d = product_distribution(Normal.(randn(3), rand(3)))
DiagNormal(...)

julia> cov(d) isa Matrix
true

But we'd like to have cov(d) isa Diagonal, of course, so that cov doesn't blow up our memory for (very) large products of Normal (resp. a large diagonal MvNormal).

This happens because we currently define

cov(d::MvNormal) = Matrix(d.Σ)

(so cov will annoyingly also silently convert a GPU array to a CPU array). Do we have a better generic mechanism to turn an AbstractPDMat into a "normal matrix" that's not automatically a Matrix? We could introduce something like

without_chol(A::AbstractMatrix) = A
without_chol(A::PDMat) = A.mat
without_chol(A::PDiagMat) = Diagonal(diag(A))
without_chol(A::PDSparseMat) = A.mat
without_chol(A::ScalMat) = Diagonal(Fill(A.value, A.dim))

[devmotion]

Fixed by #1373 (I assume).

[devmotion]

Although that won't give you a Diagonal but whatever the covariance matrix is. I.e., an AbstractPDMat for MvNormal. IMO this is fine though since they are also optimized AbstractArrays.

[oschulz]

IMO this is fine though since they are also optimized AbstractArrays.

Oh yes, that's much better, since it preserves the Cholesky decomposition (for non-diagonals). I'll close this, then.