Lanczos过程 Python实现

import numpy as np
from sklearn.datasets import make_spd_matrix
def Lanczos(A,r,m):
    '''
    当A时对称正定矩阵时，上Hessenberg矩阵化为三对角阵
    :param r: 列向量 
    :param A: 矩阵
    :param m: 迭代次数
    :return: 返回一个标准正交基
    '''
    v = []
    n1 = A.shape[0]
    v0 = np.zeros(n1)
    v.append(v0)
    beta0 = 0
    beta = []
    beta.append(beta0)
    v1 = r / np.linalg.norm(r,2)
    v.append(v1)
    for j in range(m-1):
        z = np.dot(A,v[j+1])
        alpha = np.inner(v[j+1],z)
        z = z - alpha*v[j+1] - beta[j]*v[j]
        '''修正'''
        alpha = np.inner(v[j+1],z)
        z = z - alpha * v[j+1]
        betaj = np.linalg.norm(z,2)
        beta.append(betaj)
        if betaj == 0:
            break
        v2 = z / betaj
        v.append(v2)
    return v[1::]


A =make_spd_matrix(10, random_state=42)#随机生成对阵正定矩阵
m = 10
r = [3,9,1,9,7,6,5,9,10,7]
print(A)
print(Lanczos(A,r,m))

当A为对称正定矩阵时，Arnoldi过程可化简为Lanczos过程