How to most efficiently use QR-decomposition in Julia?

Question

Avoiding array allocations is good for performance. However, I have yet to understand what is the most possible efficient way one can perform a QR decomposition of a matrix A. (note: both Q and R matrices are needed)

Simply using Q, R = qr(A) is probably not the best idea, since it allocates both Q and R, where both could be re-allocated.

The function qrfact allows one to store factorization in a packed format. However, I would still write afterwards: F = qrfact(A); Q = F[:Q]; R = F[:R] once again allocating new arrays for Q and R. Finally, the documentation also suggests the qrfact! function, which saves space by overwriting the input A, instead of creating a copy. However, if one uses F = qrfact!(A) the over-written A is not useful in the sense that it is not either Q or R, which one (specifically, I) would need.

So my two questions are:

1. What is the best/most efficient way to perform a QR decomposition if you only care about the matrices Q and R and you have no problem re-allocating them.

2. What is actually written in the matrix A when one calls qrfact!(A) ?

Answer 1

In

F = qrfact!(A)

or

F = qrfact(A)

F[:Q] and F[:R] do not allocate new dense arrays; they are simply views over the packed format from which Q and R are easily computed. This means that qrfact!(A) doesn't need to allocate arrays for Q and R, it simply computes the packed format in place for A.

However, that also means that F[:Q] and F[:R] cannot be mutated. If you need to modify one of them for whatever reason, you will need to collect it into a mutable Array, and this will certainly allocate. It will still be more efficient to use qrfact!(A) instead of qrfact(A), because the latter will allocate space for the packed QR factorization as well as for the collected Array.