# Mathematica code to find R (upper triangular matrix) in QR decomposition using reflection method

I am trying to write Mathematica code to find $R$ (upper triangular matrix) in $A=QR$ $A \in \Bbb{R}^{n\times n}$decomposition using reflection method. I wrote for $3 \times 3$ matrix. Any tips for write it for $n$ dimension.

• How did you implement your Householder reflections? Dec 20, 2016 at 10:18
• Finding the rotation angle and Q[b_] := Module[{[Theta]}, [Theta] = Theta[b]; R[[Theta]].( { {1, 0}, {0, -1} } ).R[-[Theta]]] Dec 20, 2016 at 10:21
• Why not try to write a method for the Householder reflection for $n$-dimensions before anything else? Dec 20, 2016 at 10:23
• Have you see the Method option of Orthogonalize and the Matrix Decompositions guide? Dec 20, 2016 at 13:55

Here's a pretty basic (and inefficient) routine for getting the QR decomposition of a matrix via Householder reflection:

qrd[mat_?MatrixQ] := Module[{r = mat, h, m, n, q, v, v2},
{m, n} = Dimensions[r]; q = IdentityMatrix[m];
Do[v = PadLeft[r[[k ;;, k]], m];
v2 = v - SparseArray[{k -> Norm[v]}, m];
h = If[! TrueQ[Norm[v2, ∞] == 0], ReflectionMatrix[v2], IdentityMatrix[m]];
q = q.h; r = h.r,
{k, n}];
{q, r}]


which returns the $$\mathbf Q$$ and $$\mathbf R$$ factors in a list.

It could be made more efficient by exploiting the structure of the Householder update (see e.g. Golub and Van Loan or Stewart for how to do that), but I'll leave that for the OP or someone else to do.

• Thank you and I am trying to do so. Dec 20, 2016 at 12:20
• Of course, since QRDecomposition[] is actually built-in, the code here is more of a toy than anything else. Dec 20, 2016 at 14:42