# RQ decomposition

**URL:** <https://discourse.julialang.org/t/rq-decomposition/112795>\
**Category:** Numerics\
**Created:** [April 10, 2024, 6:20pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795 "2024-04-10T18:20:26Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![RJDennis](https://avatars.discourse-cdn.com/v4/letter/r/90db22/32.png) [@RJDennis](https://discourse.julialang.org/u/RJDennis)\
**Post date:** [April 10, 2024, 6:20pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/1 "2024-04-10T18:20:26Z")

</div>

Hi,

can anyone point me to a Julia implementation of the RQ decomposition.

Thanks.

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [April 10, 2024, 6:29pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/2 "2024-04-10T18:29:01Z")

</div>

MatrixFactorizations.jl

I’m curious why you want it?

---

<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:** [April 10, 2024, 6:47pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/3 "2024-04-10T18:47:47Z")

</div>

> [@dlfivefifty](#):
>
> I’m curious why you want it?

I’ve never needed the RQ factorization by itself, but Givens RQ is a useful starting point to perform in-place `ldiv!` solves or `logdet` calculations with Hessenberg matrices ([(H+­μI) \ x solvers for Hessenberg factorizations by stevengj · Pull Request #31853 · JuliaLang/julia · GitHub](https://github.com/JuliaLang/julia/pull/31853)), following [Henry (1994)](https://dl.acm.org/doi/book/10.5555/866641).

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [April 10, 2024, 6:53pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/4 "2024-04-10T18:53:58Z")

</div>

Cool. I don’t understand your first point: I thought QR can be used in place?

I also don’t understand your second point but I’ll have to read the reference first 😅

I’m really interested in the ∞ dimensional case where QR and QL are fundamentally different. Eg QR is unique but you can have multiple QL factorisations depending on the freehold index. I haven’t thought much about RQ and LQ as they are equivalent to QL and QR applied to the adjoint but perhaps I’ve overlooked a nuance.

---

<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:** [April 10, 2024, 7:05pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/5 "2024-04-10T19:05:46Z")

</div>

> [@dlfivefifty](#):
>
> Cool. I don’t understand your first point: I thought QR can be used in place? I also don’t understand your second point but I’ll have to read the reference first 😅

TLDR: basically, you carry out the Givens QR algorithm _at the same time_ as you apply backsubstitution to solve Hx = b \Longleftrightarrow Rx = Q^\* b, or more generally a shifted system (H - \mu I) x = b, not separately. (Both Givens RQ and backsubstitution iterate over the rows in reverse order.) This allows you to leave H _unmodified_ (which is especially useful if you are going to repeat this process for multiple different shifts \mu) rather than overwriting it, and instead you only need auxiliary storage for _one row_ of R (the current row for backsubstitution). (Similarly, Q is not stored, but applied one Givens rotation at a time to b and to construct the preceding row of R.)

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [April 10, 2024, 7:12pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/6 "2024-04-10T19:12:02Z")

</div>

And the reason you need Givens instead of Householder is to avoid creating a householder vector? Interesting!

---

<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:** [April 10, 2024, 7:14pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/7 "2024-04-10T19:14:40Z")

</div>

> [@dlfivefifty](#):
>
> And the reason you need Givens instead of Householder is to avoid creating a householder vector?

You use Givens instead of Householder in order to exploit the sparsity pattern of an upper Hessenberg matrix H, which gives you O(m^2) complexity instead of O(m^3).

It’s for the same reason that [Givens QR is used for QR iteration on Hessenberg matrices](https://people.inf.ethz.ch/arbenz/ewp/Lnotes/chapter4.pdf), in order that the QR iteration can be O(m^2) by exploiting the Hessenberg sparsity.

(Similarly, we could implement fast `qr` methods for `Tridiagonal`, `UpperHessenberg`, etc. using Givens rotations, but currently it seems to just fall back to the O(m^3) dense Householder algorithm.)

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [April 10, 2024, 7:23pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/8 "2024-04-10T19:23:54Z")

</div>

I don’t understand this. Givens is just a 2x2 Householder with a sign change… You can easily do banded Householder and the factors matrix is banded… I never understood the logic for using Givens for banded matrices (and I wrote the only extant banded QR code publicly available that I know of!)

---

<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:** [April 10, 2024, 7:27pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/9 "2024-04-10T19:27:05Z")

</div>

> [@dlfivefifty](#):
>
> I don’t understand this. Givens is just a 2x2 Householder with a sign change… You can easily do banded Householder and the factors matrix is banded…

It sounds like you understand it perfectly — you need to use some form of 2x2 banded Householder to exploit the sparsity of an upper-Hessenberg matrix, and the conventional name for that algorithm is “Givens”.

(Whereas the unqualified term “Householder QR” usually refers to the non-banded version.)

---

<div class="post-metadata">

**Author:** ![dlfivefifty](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dlfivefifty/32/1959_2.png) [@dlfivefifty](https://discourse.julialang.org/u/dlfivefifty)\
**Post date:** [April 10, 2024, 7:30pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/10 "2024-04-10T19:30:49Z")

</div>

Ah ok! Alex Townsend insisted Givens is better for banded matrices and I never understood why 😅 I guess you then agree with me that if the matrix has more bands than tridiagonal/Hessenberg then “banded Householder” is more sensible then Givens

---

<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:** [April 10, 2024, 7:34pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/11 "2024-04-10T19:34:59Z")

</div>

> [@dlfivefifty](#):
>
> Alex Townsend insisted Givens is better for banded matrices and I never understood why 😅 I guess you then agree with me that if the matrix has more bands than tridiagonal/Hessenberg then “banded Householder” is more sensible then Givens

Givens (skipping zero blocks) and banded Householder (matched to bandwidth & skipping zero blocks) seem like they should have the same arithmetic complexity for banded matrices, no? (Givens is certainly better than “default” non-banded Householder, maybe that’s what Alex meant?)

But banded Householder might be better than Givens for large bandwidths for the same reasons that BLAS3 is better than BLAS2 (bigger blocks = better cache locality, wider SIMD, better parallelization, etcetera).

(I’m afraid we’ve hijacked the discussion from the OP!)

---

<div class="post-metadata">

**Author:** ![RJDennis](https://avatars.discourse-cdn.com/v4/letter/r/90db22/32.png) [@RJDennis](https://discourse.julialang.org/u/RJDennis)\
**Post date:** [April 10, 2024, 10:01pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/12 "2024-04-10T22:01:05Z")

</div>

Thanks.

I was trying to perform a particular orthogonalisation of a tensor decomposition.

---

<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:** [April 10, 2024, 11:13pm UTC](https://discourse.julialang.org/t/rq-decomposition/112795/13 "2024-04-10T23:13:04Z")

</div>

> [@RJDennis](#):
>
> can anyone point me to a Julia implementation of the RQ decomposition.

Note that RQ is [equivalent to a QR decomposition](https://math.stackexchange.com/questions/1640695/rq-decomposition) of the transpose after reversing the order of the rows (i.e. `qr(reverse(A,dims=1)')`), or more simply to the LQ decomposition after reversing the rows (i.e. to `lq(reverse(A,dims=1))`).

```julia
julia> using LinearAlgebra

julia> A = rand(3,5);

julia> F = lq(reverse(A, dims=1));

julia> R, Q = reverse(F.L), reverse(Matrix(F.Q), dims=1);

julia> A ≈ R * Q # check A = RQ
true

julia> istriu(R) # check upper triangular
true

julia> Q*Q' ≈ I # check orthonormal rows
true

```

> [@RJDennis](#):
>
> I was trying to perform a particular orthogonalisation of a tensor decomposition.

This is a circular answer to the question — you are saying you want the RQ factorization because you want that particular orthogonalization. But for what _reason_ (in what context) do you need this _particular_ orthogonalization, as opposed to e.g. LQ?