# Initialization of vector of matrix decompositions

**URL:** <https://discourse.julialang.org/t/initialization-of-vector-of-matrix-decompositions/58909>\
**Category:** Performance\
**Tags:** question\
**Created:** [April 9, 2021, 8:51am UTC](https://discourse.julialang.org/t/initialization-of-vector-of-matrix-decompositions/58909 "2021-04-09T08:51:54Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![BambOoxX](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bambooxx/32/22179_2.png) [@BambOoxX](https://discourse.julialang.org/u/BambOoxX)\
**Post date:** [April 9, 2021, 8:51am UTC](https://discourse.julialang.org/t/initialization-of-vector-of-matrix-decompositions/58909/1 "2021-04-09T08:51:54Z")

</div>

For one of my projects, I need to perform multiple QR (up to multiple hundreds) decompositions of dense tall complex matrices. These decompositions are later used in some linear algebra operations, where the operations are performed on views of subblocks.

After trying `qr`, `qrfact` and their **bang** versions, I am facing the issue of the storage.  
The packed QR format returning objects of type e.g. `LinearAlgebra.QRCompactWY` seems to be efficient both in terms of allocation and for later manipulation when called directly, but in my case, i will have to use subarrays of such objects.  
Is it therefore still better to store the decompositions this way of to store separately the`Q and R matrices ?  
Also,my initial guess to store the decompositions was to initialize a vector of those.  
However, while on can initialize a vector of matrices using a comprehension e.g.

```julia
[Matrix{ComplexF64}(undef, i, j) for _ in 1:k]

```

it does not seem to be possible to preallocate for

```julia
[LinearAlgebra.QRCompactWY{Float64,Array{Float64,2}}(undef, i, j) for _ in 1:k]

```

How should this storage be pre-allocated ?

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [April 9, 2021, 9:03am UTC](https://discourse.julialang.org/t/initialization-of-vector-of-matrix-decompositions/58909/2 "2021-04-09T09:03:59Z")

</div>

You can’t fully preallocate a QR factorization using the julia wrappers anyway: eg try

```julia
a = randn(2,2)
b = qr!(a)
a .= 0
b.Q

```

---

<div class="post-metadata">

**Author:** ![BambOoxX](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bambooxx/32/22179_2.png) [@BambOoxX](https://discourse.julialang.org/u/BambOoxX)\
**Post date:** [April 9, 2021, 9:09am UTC](https://discourse.julialang.org/t/initialization-of-vector-of-matrix-decompositions/58909/3 "2021-04-09T09:09:47Z")

</div>

From there I see only two alternatives :

- iteratively `push!` the decompositions into the vector
- Pre-allocate for dense Q and R matrices

Do you concur ?

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [April 9, 2021, 9:14am UTC](https://discourse.julialang.org/t/initialization-of-vector-of-matrix-decompositions/58909/4 "2021-04-09T09:14:31Z")

</div>

My point is that in terms of pre-allocations, you can’t do better than `ret = [qr!(arr[i]) for i=1:N]` (which will overwrite `a[i]` and allocate memory for the Q factors), unless you bypass julia’s qr wrapper and go straight for the LAPACK calls.

---

<div class="post-metadata">

**Author:** ![BambOoxX](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/bambooxx/32/22179_2.png) [@BambOoxX](https://discourse.julialang.org/u/BambOoxX)\
**Post date:** [April 9, 2021, 2:08pm UTC](https://discourse.julialang.org/t/initialization-of-vector-of-matrix-decompositions/58909/5 "2021-04-09T14:08:19Z")

</div>

I tested the method you proposed, but I is overall slower than by a direct allocation of the matrices, seemingly due to the access on views afterwards.