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:
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.
What is actually written in the matrix A
when one calls qrfact!(A)
?
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 collect
ed Array
.