Maple紧致ludecomposition A.x=B，x为矩阵答案

【问题标题】：Maple compact ludecomposition A.x=B, x is matrixMaple紧致ludecomposition A.x=B，x为矩阵
【发布时间】：2016-10-10 12:37:21
【问题描述】：

我最近实现了来自 Numerical Recipe 2nd edition 的 LUDecomposition 例程，称为 bandec()。它使用紧凑形式的带状矩阵创建 lu 分解（它还以紧凑形式返回 L 和 U）。我的问题是，我如何求解方程组

A.x=B，如果 x 是矩阵？

有没有我可以使用的例程？

【问题讨论】：

标签： matrix maple

【解决方案1】：

你写了“稀疏带状”。虽然没有针对带状情况的特殊稀疏求解器，但可以使用 LAPACK 的（密集）带状求解器。我写这篇文章的机会很小，可能实际上是你的情况。

例如，

restart;
with(LinearAlgebra):
N:=10000:
M:=RandomMatrix(N,datatype=float[8]):
V:=RandomVector(N,datatype=float[8]):
infolevel[LinearAlgebra]:=1:
X:=CodeTools:-Usage(LinearSolve(M,V)):
  LinearSolve: using method LU
  LinearSolve: calling external function
  LinearSolve: NAG hw_f07adf
  LinearSolve: NAG hw_f07aef
  memory used=0.75GiB, alloc change=0.78GiB, cpu time=16.24s, real time=4.30s, gc time=8.00ms
Norm(M.X-V);
  unknown: NAG hw_f06paf
  unknown: NAG hw_f06paf
  Norm: calling external function
  Norm: NAG:  hw_f06raf
                                      -8
                1.06381179421077832 10  

restart;
with(LinearAlgebra):
N:=10000:
B:=max(1,floor(0.005*N)):
2*B+1;
                          101
M:=RandomMatrix(N,shape=band[B,B],datatype=float[8]):
V:=RandomVector(N,datatype=float[8]):
infolevel[LinearAlgebra]:=1:
X:=CodeTools:-Usage(LinearSolve(M,V)):
  LinearSolve: using method LU
  LinearSolve: calling external function
  LinearSolve: CLAPACK hw_dgbtrf_
  LinearSolve: CLAPACK hw_dgbtrs_
  memory used=13.09MiB, alloc change=11.52MiB, cpu time=20.00ms, real time=23.00ms, gc time=0ns
Norm(M.X-V);
  Multiply: calling external function
  Multiply: NAG hw_f06pbf
  Multiply: calling external function
  Multiply: NAG hw_f06pbf
Norm: calling external function
  Norm: NAG:  hw_f06raf
                                     -11
               8.03126454229641240 10

【讨论】：

嗨！是的，这是我的情况，谢谢！我正在尝试 1000x1000 的矩阵，它没有及时给出任何速度。它实际上变得更慢了。只有当大小增加到 10000 时，带状的速度才会更快。你知道为什么较小的矩阵会变慢吗？矩阵真的很稀疏，所以它不应该总是更快地求解吗？
上面的代码是否显示了相似的时间，比较，完全按照写的？而当你只改变大小到1000以上？在矩阵/向量上使用 datatype=float[8] 是关键。该 infolevel 设置显示使用了哪些 clapack 函数，并且在 band 情况下显示它使用来自 LAPACK 的专用 band-LU（打包存储）。如果您发布为 www.mapleprimes.com 并将实际工作表示例上传到您的问题中，则调试会更容易。

【解决方案2】：

常规命令 LinearAlgebra:-LinearSolve 处理该问题。如果 B 是一个矩阵，那么 x 也是。

【讨论】：

是的，效果很好，但我一直在尝试实现一个巨大的稀疏带状矩阵。 Maple 提供的所有“稀疏”方法都变慢了（这很奇怪！）。这就是为什么我试图实现一个利用带状矩阵非常稀疏这一事实的例程。
我看不到您的评论和您的问题之间的联系。 A的稀疏性与x是否为矩阵有什么关系？
是的，对不起，我不清楚。我有一个稀疏带状矩阵 J，我想求解系统 J.X=B，其中 X 是一个矩阵。这可以使用 LinearSolve 完成，但速度非常慢。所以我一直在尝试实现一个利用稀疏性的 LUdecomposition 代码，但现在我需要一个执行 LinearSolve 所做的例程（对于上面的 J.X=B），但这与 bandec() 的输出兼容，即紧凑形式的 PLU 矩阵。