I implemented the Householder transformation in Python, so that I can later use it in a QR decomposition. Unfortunately I haven't found a good concise source for reading up on the algorithm.

I am not really satisfied with my code - I mostly dislike its readability due to the number of various computations. I also find my usage of the identity matrices \$I\$ and \$J\$ inefficient. So I would be happy about any suggestions on improving this implementation.

def houseHolder2(m):    
    l = len(m)
    for i in range(0, l - 1):
        I = np.eye(l - i)
        ai = m[:, i] 
        abs = np.linalg.norm(ai[i:l]) 
        sign = np.sign(m[i, i])
        vi = np.array([ai[i:l]]).T + sign * abs * np.array([I[:, 0]]).T
        Qi = I - 2 * (vi @ vi.T) / (vi.T @ vi)
        J = np.eye(l)
        J[i:l, i:l] = Qi
        m = J@m
    return m  
  • 2
    \$\begingroup\$ Are you aware of numpy.linalg.qr? \$\endgroup\$ – Grajdeanu Alex. Jul 29 '17 at 12:29
  • \$\begingroup\$ @MrGrj, ye, I think I have seen something of sorts. However it is not my goal to solve a specific problem or computation, but rather just improve my understanding how algorithm works. Being able to write them down in clean code should be also part of the process imo. \$\endgroup\$ – Imago Jul 29 '17 at 14:47
  • \$\begingroup\$ What's a typical m? From an efficiency view point, the real question is how many times are you repeating the loop, a few times, or hundreds. \$\endgroup\$ – hpaulj Jul 29 '17 at 23:58
  • 1
    \$\begingroup\$ I think you can simplify the indexing of I with I[:,[0]]; Or the whole line with vi=(ai[i:l] + sign * abs * I[:,0])[:,None] \$\endgroup\$ – hpaulj Jul 30 '17 at 0:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.