Is there a clever way to get rid of these loops using numpy?

I'm reaching the maximum recursion depth and I've been trying to use np.tensordot() I couldn't really get an insight into how to use it in this case.

def stopping_condtion(a,V,V_old,eps):
    return  np.max(la.norm(V - V_old)) < ((1 - a) * eps)  / a

def value_iteration(net_profit, a, P,V,k = 0):
    eps = 0.1
    m = len(net_profit)
    n = len(V)

    A = np.zeros((n,m))

    for i in range(0,n):
        for j in range(0,m):
            A[i,j] = net_profit[j,i] + a * np.sum(P[j,:,i]) * V[j] 
            
    V_new = np.max(A,axis = 1)

    if stopping_condtion(a,V_new,V,eps):
        print(f'a* =  {np.argmax(A,axis = 1)} with alpha = {a} after n = {k} iterations ')
        return np.argmax(A,axis = 1)
    
    return value_iteration(net_profit, a, P,V_new,k+1)

These are the inputs

profit = np.array([900, 800 , 600 , 400, 100])
cost   = np.array([0  , 80  , 800])

net_profit = (np.tile(profit,(3,1)).transpose() - cost).transpose()
alpha  = np.array([0.3, 0.6 , 0.9])


P  = np.array([ [[0.6, 0.4  , 0   , 0  , 0  ],
                 [0  , 0.5  , 0.3 , 0.2, 0  ],
                 [0  , 0    , 0.4 , 0.3, 0.3],
                 [0  , 0    , 0   , 0.5, 0.5],
                 [0  , 0    , 0   , 0  , 1  ]],
               
                [[0.8, 0.2  , 0   , 0  , 0  ],
                 [0  , 0.8  , 0.2 , 0  , 0  ],
                 [0  , 0.2  , 0.6 , 0.2, 0  ],
                 [0  , 0    , 0.3 , 0.6, 0.1],
                 [0  , 0    , 0   , 0.5, 0.5]],
               
                [[1  , 0    , 0   , 0  , 0  ],
                 [1  , 0    , 0   , 0  , 0  ],
                 [1  , 0    , 0   , 0  , 0  ],
                 [1  , 0    , 0   , 0  , 0  ],
                 [1  , 0    , 0   , 0  , 0  ]] ])

V = np.zeros(len(P[0,0]))
value_iteration(net_profit,alpha[0],P,V)

I was wondering if there is a way to get rid of the loop and only use the Numpy operations for efficiency sake.

You can make use of a transposition and broadcasting like the following (untested) code.

A = net_profit.T + a * np.sum(P, axis=1).T * V[:,None]