财富算法特征求解

问：

是否存在财富特征求解算法的实现？

目前，我正在尝试将其移植到 python，但到目前为止它没有提供正确的结果。 我唯一找到的是他在这里描述算法的论文，我不明白我的算法哪里出了问题。

import numpy as np
from sympy import Poly, symbols
lamda = symbols('lamda')

#Creates generalized companion matrix
def gen_companion(p,S):
    #create lagrange coeffs
    L_i = []
    for i in range(len(S)):
        den_s = 1
        for j in range(len(S)):
            if (not (i == j)):
                den_s *= p(S[i]-S[j])
        p_s = p(S[i])/den_s
        L_i.append(p_s)
    #create lagrange matrix L 
    dim = S.shape[0]
    L = np.zeros((dim,dim),np.complex64)
    for i in range(dim):
        for j in range(dim):
            L[j,i] = np.abs(L_i[i])
    B = np.zeros((dim,dim),np.complex64)
    #Create Matrix S
    for i in range(dim):
        B[i,i] = S[i]
    return B-L
    


def eigensolve(p,n_iterations):
    A = np.polynomial.polynomial.polycompanion(coeffs).astype(np.complex64)
    S = qr_algorithm(A,n_iterations)
    for _ in range(n_iterations):
        L = gen_companion(p,S)
        S = qr_algorithm(L,n_iterations)
    return S[0]


def qr_algorithm(A,n_iterations):
    for _ in range(n_iterations):
        Q,R = np.linalg.qr(A)
        A = R @ Q
    return np.array([A[i,i] for i in range(A.shape[0])])

p = Poly((lamda - 10)*(lamda-20)*(lamda-4)*(lamda-5),lamda)
coeffs = np.array(list(reversed(p.all_coeffs())))


n_iterations = 10
v = eigensolve(p,n_iterations)
print("")

我只得到了我在算法的第二部分输入的相同特征值：

(20.153511+0j)
(20.153511+0j)
(20.153511+0j)
(20.153511+0j)
(20.153511+0j)