# Generate a positive definite matrix

**URL:** <https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582>\
**Category:** General Usage\
**Tags:** question, matrices\
**Created:** [October 18, 2020, 12:43pm UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582 "2020-10-18T12:43:16Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![F-YF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/f-yf/32/17363_2.png) [@F-YF](https://discourse.julialang.org/u/F-YF)\
**Post date:** [October 18, 2020, 12:43pm UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582/1 "2020-10-18T12:43:16Z")

</div>

How do I randomly generate a positive definite matrix？ What functions can directly implement this step？

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [October 18, 2020, 1:06pm UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582/2 "2020-10-18T13:06:40Z")

</div>

Depends on exactly what you want. A trivial way is

```julia
A = randn(n,n); A = A'*A; A = (A + A')/2

```

If you want to control the eigenvalues, you can use something like

```julia
Q, _ = qr(randn(n, n)); D = Diagonal(eigvals); A = Q*D*Q'

```

where `eigvals` is a vector of eigenvalues.

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [October 18, 2020, 1:08pm UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582/3 "2020-10-18T13:08:12Z")

</div>

Aside: having gotten used to JuliaMono (thanks, @cormullion!), the typesetting of that `*` looks really weird and confusing…

---

<div class="post-metadata">

**Author:** ![cormullion](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cormullion/32/49131_2.png) [@cormullion](https://discourse.julialang.org/u/cormullion)\
**Post date:** [October 18, 2020, 1:39pm UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582/4 "2020-10-18T13:39:24Z")

</div>

Marginal note: I’m very happy if I’ve contributed even a tiny amount to your legendary productivity, Tim. (Even though I’ve seen a complaint about the asterisk. But it was Hacker News, so to be expected 😀.)

---

<div class="post-metadata">

**Author:** ![yha](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/yha/32/3502_2.png) [@yha](https://discourse.julialang.org/u/yha)\
**Post date:** [October 18, 2020, 4:58pm UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582/5 "2020-10-18T16:58:50Z")

</div>

> [@tim.holy](#):
>
> Depends on exactly what you want. A trivial way is
> 
> ```julia
> A = randn(n,n); A = A'*A; A = (A + A')/2
> 
> ```

Is that `A = (A + A')/2` for numerical reasons? It doesn’t seem to be necessary:

```julia
As = [(A = randn(n,n); A'*A) for _=1:10_000]
symm_err(A) = maximum(abs.(A-A'))
maximum(symm_err.(As)) # == 0.0

```

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [October 19, 2020, 9:25am UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582/6 "2020-10-19T09:25:33Z")

</div>

> It doesn’t seem to be necessary

Hmm, interesting. I _know_ I’ve found cases where it is (where the result has ulp-level errors), so I just add it by default now. But I get the same thing you do from your test case.

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [October 19, 2020, 10:32am UTC](https://discourse.julialang.org/t/generate-a-positive-definite-matrix/48582/7 "2020-10-19T10:32:56Z")

</div>

You can map the Cholesky factor of a PD matrix into a vector of reals. The code below uses [`TransformVariables.CorrCholeskyFactor`](https://tamaspapp.eu/TransformVariables.jl/dev/#TransformVariables.CorrCholeskyFactor), which maps to a Cholesky factor of the _correlation_ matrix (the diagonal of `U'*U` is ones):

```julia
using TransformVariables
n = 10
C = CorrCholeskyFactor(n)
U = C(rand(dimension(C)))
σ = abs.(randn(n))
L = U * Diagonal(σ)
A = L' * L

```

then `A` is guaranteed to be PD, also, _all_ PD matrices of this size will be generated for some random values.

In practice, you would want to work with the factors (eg `L`) for almost all numerical applications.
