BUG: sparse.linalg.lsqr: 返回的方差显示出系统性错误

import numpy as np import matplotlib.pyplot as plt from scipy.sparse import linalg from scipy.linalg import cho_factor, cho_solve np.random.seed(0) fig2, ax_metas = plt.subplots(2, 4) fig, ax = plt.subplots() for n in (250, 501, 1001): ms = [] ratios = [] metas = [] pylops_ratios = [] for m in np.unique(np.logspace(0, np.log10(n), 100).astype(int)): if m > n: continue A = np.random.randn(n, m) x_true = np.random.randn(m) #noise, var = (0, np.ones(n)) var = np.random.uniform(0.5, 1.5, n) noise = np.random.randn(n) * 0.1 y = A @ x_true + noise C = np.diag(var) Cinv = np.diag(1/var) # Method 1: Direct solve direct_cov = np.linalg.solve(A.T @ Cinv @ A, np.eye(m)) direct_var = np.diag(direct_cov) # Method 2: LSQR with Cholesky L = np.linalg.cholesky(C) Linv = np.linalg.inv(L) A_white = Linv @ A y_white = Linv @ y x_chol, *meta, var_chol = linalg.lsqr(A_white, y_white, calc_var=True, show=False, atol=1e-12, btol=1e-12, damp=1e-6, iter_lim=10000000) ms.append(m) ratios.append(np.mean(var_chol/direct_var)) …