Print Python Recursive Array-CodePudding

assume there are two A & R matrices which are 4x4 and 4x1:

import numpy as np
from numpy.linalg import *
A = np.array([[0,1,1,1],
              [1,0,0,1],
              [0,0,0,1],
              [0,1,1,0]])

R = np.ones((4,1))/4

and there is a function to transpose A matrix after normalizing, which is

def normalized_Adj_matrix_transpose(x):
    D = np.sum(x, axis = 1) # degree matrix
    D_reciprocal = np.array([1/i for i in D]).reshape(len(D),1)
    return np.transpose(x*D_reciprocal)

My goal is to keep updating R matrix until it converged, where R = ((A*D^-1)^T)^n.dot(R) for n is a real number.

R1 = normalized_Adj_matrix_transpose(A).dot(R)
R2 = normalized_Adj_matrix_transpose(A).dot(R1)
R3 = normalized_Adj_matrix_transpose(A).dot(R2)

print(R1.reshape((1,4)))
print(R2.reshape((1,4)))
print(R3.reshape((1,4)))

output:

[[0.125      0.20833333 0.20833333 0.45833333]]
[[0.10416667 0.27083333 0.27083333 0.35416667]]
[[0.13541667 0.21180556 0.21180556 0.44097222]]

My question is: how to use for loop to recursively update the R matrix and print each iteration result? many thanks

such as

iteration 1: [[0.125      0.20833333 0.20833333 0.45833333]]
iteration 2: [[0.10416667 0.27083333 0.27083333 0.35416667]]
iteration 3: [[0.13541667 0.21180556 0.21180556 0.44097222]]

CodePudding user response：

Why don't you put everything in a while loop with an if-condition in it that checks if R1 is converging. Select a small enough e to check for convergence:

x = normalized_Adj_matrix_transpose(A)
R1 = R.copy()
e = np.array([10**(-9)]*(len(R))).reshape(-1,1)
R1_list = []
while True:
    R2 = x.dot(R1)
#    print(R1.reshape(1,4))
#    print(R2.reshape(1,4))
#    R1_list.append(R1.reshape(1,4))
    if ((R2 - R1) < e).all():
        break
    else:
        R1 = R2

Output:

print(R1.reshape(1,4))  #  [[0.12 0.24 0.24 0.4 ]]
print(R2.reshape(1,4))  #  [[0.12 0.24 0.24 0.4 ]]

You can uncomment print inside the loop to see how R1 is converging but printing is not interesting in my opinion. If you want to maybe plot it later, it might be better to keep them in a list e.g. R1_list (you can also uncomment that and use it).

You can also optimize your normalized_Adj_matrix_transpose as:

def normalized_Adj_matrix_transpose(x):
    return A.T /A.sum(axis=1)