invalid value encountered in true_divide for Gram_Schimdt algorithm

Viewed 41

I was trying to implement Gram-Schimdt using the following function. However, it seems that there is a zero divider in the code. But I am not sure how to solve it.

def gs(A):
    x = A.shape[0]
    y = A.shape[1]
    U = np.zeros((x,y))
    U[:,1] = A[:,1]/np.linalg.norm(A[:,1])
    for i in range(1,y):
        U[:,i] = A[:,i]
        for j in range(0,i):
            U[:,i] = U[:,i] - U[:,j].dot(U[:,i])*U[:,j]/(U[:,j].dot(U[:,j]))
        U[:,i] = U[:,i]/np.linalg.norm(U[:,i])
    return U

error message:

/var/folders/tr/ssszhd_x03jdz162mg0j01vm0000gn/T/ipykernel_39952/685295545.py:9: RuntimeWarning: invalid value encountered in true_divide
  U[:,i] = U[:,i] - U[:,j].dot(U[:,i])*U[:,j]/(U[:,j].dot(U[:,j]))
1 Answers

The problem was in the initialization

U[:,1] = A[:,1]/np.linalg.norm(A[:,1])

The first column is U[:,0], since you are asigning to U[:,1], the U[:,0] will remain all zeros, and the loop will fill the matrix with zero-divisions.

Fixing that line it works

def gs(A):
    x = A.shape[0]
    y = A.shape[1]
    U = np.zeros((x,y))
    U[:,0] = A[:,0]/np.linalg.norm(A[:,0])
    for i in range(1,y):
        U[:,i] = A[:,i]
        for j in range(0,i):
            U[:,i] = U[:,i] - U[:,j].dot(U[:,i])*U[:,j]/(U[:,j].dot(U[:,j]))
        U[:,i] = U[:,i]/np.linalg.norm(U[:,i])
    return U

# test
A = gs(np.random.rand(100, 50))
assert np.allclose(A.T @ A, np.eye(50)) # check the results
Related