skprocrustes.py 文件源码

def debug_bidiag(i, s, inds, A, B, U, V, T):
    print("\n       ********* DEBUGGING BLOCKBIDIAG: ************\n")
    # We will check the recurrence relations listed in Karimi, Toutounian 
    print("\n        Iteration i = {}, inds = {}\n".format(i, inds))
    E1 = np.zeros((inds+s, s))
    E1[0:s, :] = np.eye(s,s)
    errorRecurrence1 = sp.norm(B-np.dot(U[:,0:inds+s], np.dot(E1, B1)))
    print("\n        B - UU(i+1)*E1*B1 = {}\n".format(errorRecurrence1))
    #
    # AVk = Ukp1Tk
    errorRecurrence2 = sp.norm(np.dot(A, V[:, 0:inds]) - np.dot(U[:, 0:inds+s], T[0:inds+s, 0:inds]))
    print("\n        A*VV(i) - UU(i+1)T(i) = {}\n".format(errorRecurrence2))
    #
    # ATUkp1 = VkTkT + Vkp1Akp1Ekp1T
    Eip1 = np.zeros((inds+s, s))
    Eip1[inds:inds+s, :] = np.eye(s,s)
    errorRecurrence3 = sp.norm(np.dot(A.T, U[:, 0:inds+s]) - np.dot(V[:, 0:inds], T[0:inds+s, 0:inds].T) - np.dot(V[:, inds:inds+s], np.dot(Aip1, Eip1.T)))
    print("\n        A.T*UU(i+1) - VV(i)*T(i).T - V(i+1)*A(i+1)*E(i+1).T = {}\n".format(errorRecurrence3))