return a qrdecomposition of the matrix a using householder reflection the qrdecomposition decomposes the matrix a of shape m n into an orthogonal matrix q of shape m m and an upper triangular matrix r of shape m n note that the matrix a does not have to be square this method of decomposing a uses the householder reflection which is numerically stable and of complexity on3 https en wikipedia orgwikiqrdecompositionusinghouseholderreflections arguments a a numpy ndarray of shape m n note several optimizations can be made for numeric efficiency but this is intended to demonstrate how it would be represented in a mathematics textbook in cases where efficiency is particularly important an optimized version from blas should be used a np array12 51 4 6 167 68 4 24 41 dtypefloat q r qrhouseholdera check that the decomposition is correct np allcloseqr a true check that q is orthogonal np allcloseqq t np eyea shape0 true np allcloseq tq np eyea shape0 true check that r is upper triangular np allclosenp triur r true select a column of modified matrix a construct first basis vector determine scaling factor construct vector v for householder reflection construct the householder matrix pad with ones and zeros as necessary return a qr decomposition of the matrix a using householder reflection the qr decomposition decomposes the matrix a of shape m n into an orthogonal matrix q of shape m m and an upper triangular matrix r of shape m n note that the matrix a does not have to be square this method of decomposing a uses the householder reflection which is numerically stable and of complexity o n 3 https en wikipedia org wiki qr_decomposition using_householder_reflections arguments a a numpy ndarray of shape m n note several optimizations can be made for numeric efficiency but this is intended to demonstrate how it would be represented in a mathematics textbook in cases where efficiency is particularly important an optimized version from blas should be used a np array 12 51 4 6 167 68 4 24 41 dtype float q r qr_householder a check that the decomposition is correct np allclose q r a true check that q is orthogonal np allclose q q t np eye a shape 0 true np allclose q t q np eye a shape 0 true check that r is upper triangular np allclose np triu r r true select a column of modified matrix a construct first basis vector determine scaling factor construct vector v for householder reflection construct the householder matrix pad with ones and zeros as necessary
import numpy as np


def qr_householder(a: np.ndarray):
    
    m, n = a.shape
    t = min(m, n)
    q = np.eye(m)
    r = a.copy()

    for k in range(t - 1):
        
        x = r[k:, [k]]
        
        e1 = np.zeros_like(x)
        e1[0] = 1.0
        
        alpha = np.linalg.norm(x)
        
        v = x + np.sign(x[0]) * alpha * e1
        v /= np.linalg.norm(v)

        
        q_k = np.eye(m - k) - 2.0 * v @ v.T
        
        q_k = np.block([[np.eye(k), np.zeros((k, m - k))], [np.zeros((m - k, k)), q_k]])

        q = q @ q_k.T
        r = q_k @ r

    return q, r


if __name__ == "__main__":
    import doctest

    doctest.testmod()