.. _qrFactor:
.. _qrcode:

QR Decomposition
================

The QR decomposition is used to solve systems of linear equations
(:math:`\m{A}\bm{x} = \m{Q\, R}\bm{x} = \bm{b}`) and in the
eigendecomposition algorithm described in :ref:`QReig`, The MATLAB
command, ``[Q, R] = qr(A)`` finds the factors :math:`\bf{Q}` and
:math:`\bf{R}`. As shown in equation :eq:`eq-QR`, :math:`\bf{Q}` is
unitary, and :math:`\bf{R}` is upper triangular. The ``qrFactor``
function is a simple algorithm implementation. It has an *economy*
option to remove rows of zeros from :math:`\bf{R}` and unused columns
from :math:`\bf{Q}` for over-determined matrices. The ``qrFactor``
function makes use of the ``householder`` function to compute the
Householder reflection matrices that, with multiplication, convert
:math:`\bf{R}` to upper triangular form. Then, the :math:`\bf{Q}` matrix
is built as an orthogonal product of transposed reflection matrices. The
code is based on pseudo-code in Demmel’s *Applied Numerical Linear
Algebra* text [DEMMEL97]_ and Golub and Van Loan’s
*Matrix Computations* text [GOLUB13]_.

.. math::
   :label: eq-QR

   \m{A} = \m{Q\, R} = \begin{bmatrix}
                          \vertbar{}  & \vertbar{} & {} & \vertbar{} \\ 
                          \bm{q}_1  & \bm{q}_2 & \cdots{} & \bm{q}_n \\ 
                          \vertbar{}  & \vertbar{} & {} & \vertbar{}
                       \end{bmatrix}
   \begin{bmatrix}
                    r_{11} & r_{12} & \cdots{} & r_{1n} \\
                         0 & r_{22} & \cdots{} & r_{2n} \\
   \vdots{} & \vdots{} & \ddots{} &  \vdots{} \\
                         0 &  0     & \cdots{} & r_{nn}
   \end{bmatrix}

.. index:: qrFactor

::

   function [Q, R] = qrFactor(A, econ)
   % QRFACTOR - Simple QR factorization using Householder reflections.
   % [Q, R] = QRFACTOR(A), yields Q and R factors such that A = Q*R.
   %
   % A - input real or complex m-by-n matrix, (m >= n).
   % econ - (optional) enter 0 (zero) to get the economy QR,
   %        where Q is m-by-n and R is n-by-n.
   % Q - Unitary basis vectors of A (m-by-m)
   % R - Upper triangular matrix (m-by-n)
   %
   % See also: QR, EIGQR, HOUSEHOLDER
   %
       [m, n] = size(A); R = A; Q = eye(m);
       for k = 1:min((m-1), n)
           H = householder(R(k:m, k), m, k-1);
           R = H*R;
           R(k+1:m,k) = zeros(m-k,1); % no round-off error
           Q = Q*H';
           % For real numbers, H = H', but not in the complex case
       end
       % The econ option is useful for m > n matrices.
       if exist('econ','var') && econ == 0 && m > n
           Q = Q(:,1:n); R = R(1:n, :);
       end
   end

.. index:: QR decomposition