#QR-Zerlegung mit Householder

import numpy as np

def mysign(x): #Konvention sign(0)=1
    y=np.sign(x)
    if abs(x)<1e-15:
        y=1
    return y

def qr(A):
    n_row, n_col = A.shape
    P = np.eye(n_row,n_row)
    for k in range(n_col):
        T1= A[k:,k]
        alpha = -mysign(T1[0])*np.linalg.norm(T1)
        v=np.array(T1.copy())
        v[0] = v[0]-alpha
        c = 2/np.linalg.norm(v)**2
        P[k:,:] = P[k:,:]-np.outer(c*v,v@P[k:,:])
        A[k:,k:] = A[k:,k:]-np.outer(c*v,v@A[k:,k:])
    return P.transpose(), A
 
A = 1.0 * np.array( [ [0,0,-6],[1,2,3],[0,4,5] ] )

Q,R = qr(A)

print("Q=",Q)
print("R=",R)
