# Symmetric matrices

**URL:** https://discourse.julialang.org/t/symmetric-matrices/4086
**Category:** General Usage
**Tags:** proposal, performance, linearalgebra, matrices
**Created:** [June 4, 2017, 11:33pm UTC](https://discourse.julialang.org/t/symmetric-matrices/4086 "2017-06-04T23:33:18Z")
**Posts on this page:** 14
**Page:** 1

<div class="post-metadata">

### Author: ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [June 4, 2017, 11:33pm UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/1 "2017-06-04T23:33:18Z")

</div>

According to the docs, we are able to create a Symmetric matrix from another matrix that is already stored in memory:

```julia
A = eye(3)
B = Symmetric(A)
sizeof(A) # 72 bytes
sizeof(B) # 16 bytes

```

How to populate a Symmetric matrix directly without having to store a full matrix? In particular, if we could index a Symmetric matrix and have the expected behavior that whenever `B[i,j]` is set to `b`, `B[j,i]` is automatically set to `b`, that would be very useful.

I would like to save memory while working with big covariance matrices and also be able to enforce symmetry in the type to dispatch the appropriate LinAlg solvers. Is this a feature that you consider relevant?

```julia
B[1,1] # ERROR: not supported

```

---

<div class="post-metadata">

### Author: ![jebej](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jebej/32/1784_2.png) [@jebej](https://discourse.julialang.org/u/jebej)
#### Post date: [June 5, 2017, 12:30am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/2 "2017-06-05T00:30:02Z")

</div>

> [@juliohm](#):
>
> How to populate a Symmetric matrix directly without having to store a full matrix? In particular, if we could index a Symmetric matrix and have the expected behavior that whenever B[i,j] is set to b, B[j,i] is automatically set to b, that would be very useful.

If you are using dense matrices (or sparse for that matter, although see [#22200](https://github.com/JuliaLang/julia/pull/22200)), you can’t. Most algorithms depend on matrix multiplication, diagonalisation etc which are BLAS or LAPACK routines that take in the full dense array, despite not needing it all. Even if your application does not depend on those operations, the Symmetric type in julia is just a wrapper for the full array. You would need to make your own type that stores only the entries you want. Examples of similar types are `Diagonal` or `SymmetricTridiagonal`, which only store the required entries.

---

<div class="post-metadata">

### Author: ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [June 5, 2017, 12:43am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/3 "2017-06-05T00:43:11Z")

</div>

@jebej, thank you for sharing the issue on GitHub. I tried to digest the information there, but couldn’t understand why `Symmetric` cannot be made to store only half of the elements like `Diagonal` and `SymmetricTridiagonal`?

Couldn’t a Symmetric matrix be converted to a dense representation on the last minute when BLAS/LAPACK routines are called? We would still save memory, correct?

---

<div class="post-metadata">

### Author: ![jebej](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jebej/32/1784_2.png) [@jebej](https://discourse.julialang.org/u/jebej)
#### Post date: [June 5, 2017, 12:56am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/4 "2017-06-05T00:56:29Z")

</div>

> [@juliohm](#):
>
> Couldn’t a Symmetric matrix be converted to a dense representation on the last minute when BLAS/LAPACK routines are called? We would still save memory, correct?

It could, but I assume the reason it’s not done is because the benefits do not outweigh the disadvantages (i.e. it is pretty expensive to do this conversion). Plus, if you do end up making the matrix for a BLAS call, then you would need both the hypothetical `SymmetricWhichSavesMemory` matrix stored, plus the BLAS-required matrix, so it would end up using more memory.

---

<div class="post-metadata">

### Author: ![jebej](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jebej/32/1784_2.png) [@jebej](https://discourse.julialang.org/u/jebej)
#### Post date: [June 5, 2017, 1:01am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/5 "2017-06-05T01:01:28Z")

</div>

Note that if you are only looking to save memory for long-term storage, you could always store only the upper triangle in a sparse matrix. To make it dense again you could then call `Symmetric(full(A))`.

---

<div class="post-metadata">

### Author: ![nalimilan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nalimilan/32/147_2.png) [@nalimilan](https://discourse.julialang.org/u/nalimilan)
#### Post date: [June 5, 2017, 6:38pm UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/6 "2017-06-05T18:38:53Z")

</div>

LAPACK supports a [packed storage layout](http://www.netlib.org/lapack/lug/node123.html) for symmetric matrices. That’s not what `Symmetric` uses, I think for performance reasons. But it would make sense to have a type supporting it (or maybe even to do that by default) to save memory.

---

<div class="post-metadata">

### Author: ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [June 5, 2017, 11:38pm UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/7 "2017-06-05T23:38:20Z")

</div>

I would love to have these memory savings with big symmetric matrices, it sounds much smarter and efficient. I don’t have a clear picture of why this is not feasible, and thus I will leave that question to the LAPACK experts.

---

<div class="post-metadata">

### Author: ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [June 6, 2017, 1:28am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/8 "2017-06-06T01:28:48Z")

</div>

> [@juliohm](#):
>
> I would love to have these memory savings with big symmetric matrices, it sounds much smarter and efficient. I don’t have a clear picture of why this is not feasible

It’s only a factor of two in memory, which is usually not the difference between feasible and infeasible, and the price you pay is that operations on the packed format are much less efficient (you probably lose much more than a factor of two in performance for most matrix operations).

---

<div class="post-metadata">

### Author: ![Ralph\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ralph_smith/32/10344_2.png) [@Ralph\_Smith](https://discourse.julialang.org/u/Ralph_Smith)
#### Post date: [June 6, 2017, 2:52am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/9 "2017-06-06T02:52:10Z")

</div>

For some operations (esp. Cholesky and triangular solves) LAPACK has routines for **rectangular full packed** format (newer than the format @nalimilan and @stevengj address above), which

> enables the use of half the full storage while maintaining efficiency by using Level 3 BLAS/LAPACK kernels.

(detaiils [here](http://www.netlib.org/lapack/lapack-3.2.html), about haflway down the page).

If these capabilities fit the bill, the price to be paid is in developer time.

---

<div class="post-metadata">

### Author: ![nalimilan](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nalimilan/32/147_2.png) [@nalimilan](https://discourse.julialang.org/u/nalimilan)
#### Post date: [June 6, 2017, 9:09am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/10 "2017-06-06T09:09:09Z")

</div>

@Ralph_Smith Sounds very interesting. Could you file an issue describing the new format a bit and (if possible) how it could be implemented in Julia?

---

<div class="post-metadata">

### Author: ![chriscoey](https://avatars.discourse-cdn.com/v4/letter/c/3ab097/32.png) [@chriscoey](https://discourse.julialang.org/u/chriscoey)
#### Post date: [October 8, 2018, 10:21pm UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/11 "2018-10-08T22:21:16Z")

</div>

Have there been any developments or further discussions on this? Symmetric matrix operations are ubiquitous in optimization, for which Julia has great potential.

---

<div class="post-metadata">

### Author: ![Ralph\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ralph_smith/32/10344_2.png) [@Ralph\_Smith](https://discourse.julialang.org/u/Ralph_Smith)
#### Post date: [October 9, 2018, 2:13am UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/12 "2018-10-09T02:13:22Z")

</div>

cf. discussion here:  
[https://github.com/JuliaLang/julia/issues/22259](https://github.com/JuliaLang/julia/issues/22259)

Note that Andreas renamed his package with wrappers for RFP routines to [GenericLinearAlgebra](https://github.com/andreasnoack/GenericLinearAlgebra.jl)

---

<div class="post-metadata">

### Author: ![bjarthur](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bjarthur/32/9638_2.png) [@bjarthur](https://discourse.julialang.org/u/bjarthur)
#### Post date: [January 20, 2022, 2:55pm UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/13 "2022-01-20T14:55:34Z")

</div>

i’ve just made public [SymmetricFormats.jl](https://github.com/JaneliaSciComp/SymmetricFormats.jl), which provides the [packed storage format](http://www.netlib.org/lapack/lug/node123.html) for symmetric matrices that @nalimilan mentioned above. should be in the [GeneralRegistry](https://github.com/JuliaRegistries/General/pull/52804) in a couple days.

my use case is a simulation where a symmetric corellation matrix consumes the vast majority of memory. saving 2x on memory permits me to train models 2x bigger. moreover, i find that using packed arrays for this particular codebase is also almost 2x _faster_, contrary to @stevengj comment above. perhaps because i only rely on two BLAS routines ([spmv](http://www.netlib.org/lapack/explore-html/d6/d30/group __single__ blas__level2_gad1af0d0777da05d1c27ea99a69c8017c.html) and [syr](http://www.netlib.org/lapack/explore-html/d6/d30/group__single __blas__ level2_ga7b8a99048765ed2bf7c1e770bff0b622.html)). syr computes A += alpha \* x \* x’, and so the iteration over A is halved.

i’m particularly interested in feedback on an [open PR](https://github.com/JaneliaSciComp/SymmetricFormats.jl/pull/4) to this package which provides an optional way to set off-diagonal elements, which is a [frequently requested](https://github.com/JuliaLang/julia/pull/43818) feature.

---

<div class="post-metadata">

### Author: ![bjarthur](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bjarthur/32/9638_2.png) [@bjarthur](https://discourse.julialang.org/u/bjarthur)
#### Post date: [February 22, 2022, 11:04pm UTC](https://discourse.julialang.org/t/symmetric-matrices/4086/14 "2022-02-22T23:04:12Z")

</div>

i’ve also just made public [BatchedBLAS.jl](https://github.com/JaneliaSciComp/BatchedBLAS.jl), which provides additional tooling for symmetric packed arrays on GPUs. specifically, `batched_{syr,spr,symv,spmv}` and for comparison `batched_{ger,gemv,dot}`. it was surprisingly easy to get within 5% of CUDA code with pure julia kernels!
