r - R 中线性系统迭代求解器的典型用例是什么？

Question

我想知道迭代求解器是否是求解线性系统(非稀疏、对称、正定)的更快方法。

我尝试了 R 包中的共轭梯度方法 Rlinsolve和 cPCG ，但与非迭代 base::solve() 相比，两者似乎都不太准确且速度较慢.

library(Rlinsolve)
library(cPCG)
library(microbenchmark)

n <- 2000
A <- tcrossprod(matrix(rnorm(n^2),nrow=n) + diag(rep(10,n)))
x <- rnorm(n)
b <- A%*%x

mean(abs(x - solve(A,b)))
## [1] 1.158205e-08
mean(abs(x - lsolve.cg(A,b)$x))
## [1] 0.03836865
mean(abs(x - cgsolve(A,b)))
## [1] 0.02642611
mean(abs(x - pcgsolve(A,b)))
## [1] 0.02638013

microbenchmark(solve(A, b), lsolve.cg(A, b),
               cgsolve(A, b), pcgsolve(A, b), times=5)
## Unit: milliseconds
##            expr       min        lq      mean    median         uq       max neval  cld
##     solve(A, b)  183.3039  188.6678  189.7362  188.8665   189.8514  197.9914     5 a   
##  lsolve.cg(A, b) 7178.7477 7784.7646 7934.8406 8114.5838 8218.7356 8377.3714     5    d
##   cgsolve(A, b) 1907.0940 2020.8368 2226.0513 2121.2917  2483.1947 2597.8393     5  b  
##  pcgsolve(A, b) 4059.5856 4109.0319 4203.4093 4242.7750  4275.9537 4329.7005     5   c 

(带有 OpenBLAS 和 4 个内核的 R 版本 3.6.1。)

我错过了什么吗？这种迭代方法的典型用例是什么？

什么是非稀疏线性系统的好的 R 示例，可以显示迭代求解器的优点？

Answer 1

作为Rlinsolve的创建者包，我有点不同意这个问题隐含的争论。如果你有密集矩阵 A ，使用定制计算存储稀疏矩阵的所有好处立即消失。当高斯模型下的协方差结构带状时，我已经看到了稀疏求解器的一些不错的用法，但此类文献非常少。

请记住 没有一种工具可以解决所有问题 .如果您有对称的正定矩阵，cholesky 或基于 EVD 的解决方案是一个很好的工具。

仅供引用，我见过 Rlinsolve用于比较基于 EM 的迭代求解器的性能的统计计算论文，这是他们的创造，与我的包提供的方法。我相信这在某种程度上起到了很好的作用。