The following code trows a `DomainError`

after `cholesky!()`

was called:

```
mat = [2.4 2.05
2.05 2.4]
logdet(mat)
logdet(cholesky(mat))
cholesky!(mat)
logdet(mat)
```

Why is this the case?

The following code trows a `DomainError`

after `cholesky!()`

was called:

```
mat = [2.4 2.05
2.05 2.4]
logdet(mat)
logdet(cholesky(mat))
cholesky!(mat)
logdet(mat)
```

Why is this the case?

`cholesky!`

modifies the matrix in place and returns a very thin wrapper excluding the ignored data. The original binding still has that data though:

```
julia> using LinearAlgebra
julia> mat = [2.4 2.05
2.05 2.4]
2×2 Array{Float64,2}:
2.4 2.05
2.05 2.4
julia> a = cholesky!(mat)
Cholesky{Float64,Array{Float64,2}}
U factor:
2×2 UpperTriangular{Float64,Array{Float64,2}}:
1.54919 1.32327
⋅ 0.80558
julia> logdet(a)
0.4430819716794719
julia> mat
2×2 Array{Float64,2}:
1.54919 1.32327
2.05 0.80558
```

Since `cholesky!`

might be called with very large matrices, preserving the existing memory is vital for good performance. Dispatch takes it from there and efficiently works with the wrapper instead of the underlying data.

Thanks, so the most performant solution would be (like you did) calling `a = cholesky!(mat)`

and then working with `a`

?

Also, if I am calling `isposdef!(mat)`

to first check whether `mat`

is positive definite, is there a way to further work with the resulting cholesky factorization?

Because in my code, I would like to check whether `mat`

is positive definite and if that is the case, I would like to perform a cholesky factorization.

Okay, looking at LinearAlgebra I found that

```
isposdef!(A::AbstractMatrix) =
ishermitian(A) && isposdef(cholesky!(Hermitian(A); check = false))
```

since I know that my matrix is hermitian, I guess the best solution in my case is

```
mat = [2.4 2.05
2.05 2.4]
a = cholesky!(Hermitian(mat); check = false)
if isposdef(a)
... here I can work with a
else
...
```

1 Like