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 的迭代求解器（这是他们的创建）与我的包中提供的方法进行了性能比较。我相信这在某种程度上起到了很好的作用。