犰狳:用 sp_mat 解决

Posted

技术标签:

【中文标题】犰狳:用 sp_mat 解决【英文标题】:Armadillo: solve with sp_mat 【发布时间】:2014-03-07 13:18:29 【问题描述】:

我正在使用 Armadillo 编写一个基本的 FEM 程序。我使用sp_matvec 作为矩阵和向量类型。问题是,当我做solve(X, b) 时,我得到一个错误。难道solve不支持sp_mat。除了使用密集矩阵之外还有其他选择吗?下面是代码,sp_mat 无法编译。如果我使用注释行 mat A 代替它工作正常。

    int N = 3;
double h = 1./N;

//mat A = zeros<mat>(N+1,N+1);
sp_mat A(N+1,N+1); 
for(int i=0;i<=N;i++) 
    if(i>0) A(i,i-1)=-1.;
    A(i,i)=2.;
    if(i<N) A(i,i+1)=-1.;

A(N,N)=1;

vec b = zeros(N+1);
for(int i=0;i<=N;i++) 
    b(i)=h; 


vec zeta = solve(A,b);
cout << zeta;

错误:

make all 
Building file: ../src/FEM.cpp
Invoking: GCC C++ Compiler
g++ -O0 -g3 -Wall -c -fmessage-length=0 -MMD -MP -MF"src/FEM.d" -MT"src/FEM.d" -o "src/FEM.o" "../src/FEM.cpp"
../src/FEM.cpp: In function ‘int main()’:
../src/FEM.cpp:37:22: error: no matching function for call to ‘solve(arma::sp_mat&, arma::vec&)’
  vec zeta = solve(A,b);
                      ^
../src/FEM.cpp:37:22: note: candidates are:
In file included from /usr/include/armadillo:397:0,
                 from ../src/FEM.cpp:10:
/usr/include/armadillo_bits/fn_solve.hpp:25:1: note: template<class T1, class T2> const arma::Glue<T1, T2, arma::glue_solve> arma::solve(const arma::Base<typename T1::elem_type, T1>&, const arma::Base<typename T1::elem_type, T2>&, bool, const typename arma::arma_blas_type_only<typename T1::elem_type>::result*)
 solve
 ^
/usr/include/armadillo_bits/fn_solve.hpp:25:1: note:   template argument deduction/substitution failed:
../src/FEM.cpp:37:22: note:   ‘arma::sp_mat aka arma::SpMat<double>’ is not derived from ‘const arma::Base<typename T1::elem_type, T1>’
  vec zeta = solve(A,b);
                      ^
In file included from /usr/include/armadillo:397:0,
                 from ../src/FEM.cpp:10:
/usr/include/armadillo_bits/fn_solve.hpp:44:1: note: template<class T1, class T2> const arma::Glue<T1, T2, arma::glue_solve> arma::solve(const arma::Base<typename T1::elem_type, T1>&, const arma::Base<typename T1::elem_type, T2>&, const char*, const typename arma::arma_blas_type_only<typename T1::elem_type>::result*)
 solve
 ^
/usr/include/armadillo_bits/fn_solve.hpp:44:1: note:   template argument deduction/substitution failed:
../src/FEM.cpp:37:22: note:   ‘arma::sp_mat aka arma::SpMat<double>’ is not derived from ‘const arma::Base<typename T1::elem_type, T1>’
  vec zeta = solve(A,b);
                      ^
In file included from /usr/include/armadillo:397:0,
                 from ../src/FEM.cpp:10:
/usr/include/armadillo_bits/fn_solve.hpp:67:1: note: template<class T1, class T2> const arma::Glue<T1, T2, arma::glue_solve_tr> arma::solve(const arma::Op<T1, arma::op_trimat>&, const arma::Base<typename T1::elem_type, T2>&, bool, const typename arma::arma_blas_type_only<typename T1::elem_type>::result*)
 solve
 ^
/usr/include/armadillo_bits/fn_solve.hpp:67:1: note:   template argument deduction/substitution failed:
../src/FEM.cpp:37:22: note:   ‘arma::sp_mat aka arma::SpMat<double>’ is not derived from ‘const arma::Op<T1, arma::op_trimat>’
  vec zeta = solve(A,b);
                      ^
In file included from /usr/include/armadillo:397:0,
                 from ../src/FEM.cpp:10:
/usr/include/armadillo_bits/fn_solve.hpp:87:1: note: template<class T1, class T2> const arma::Glue<T1, T2, arma::glue_solve_tr> arma::solve(const arma::Op<T1, arma::op_trimat>&, const arma::Base<typename T1::elem_type, T2>&, const char*, const typename arma::arma_blas_type_only<typename T1::elem_type>::result*)
 solve
 ^
/usr/include/armadillo_bits/fn_solve.hpp:87:1: note:   template argument deduction/substitution failed:
../src/FEM.cpp:37:22: note:   ‘arma::sp_mat aka arma::SpMat<double>’ is not derived from ‘const arma::Op<T1, arma::op_trimat>’
  vec zeta = solve(A,b);
                      ^
In file included from /usr/include/armadillo:397:0,
                 from ../src/FEM.cpp:10:
/usr/include/armadillo_bits/fn_solve.hpp:110:1: note: template<class T1, class T2> bool arma::solve(arma::Mat<typename T1::elem_type>&, const arma::Base<typename T1::elem_type, T1>&, const arma::Base<typename T1::elem_type, T2>&, bool, const typename arma::arma_blas_type_only<typename T1::elem_type>::result*)
 solve
 ^
/usr/include/armadillo_bits/fn_solve.hpp:110:1: note:   template argument deduction/substitution failed:
../src/FEM.cpp:37:22: note:   ‘arma::sp_mat aka arma::SpMat<double>’ is not derived from ‘arma::Mat<typename T1::elem_type>’
  vec zeta = solve(A,b);
                      ^
In file included from /usr/include/armadillo:397:0,
                 from ../src/FEM.cpp:10:
/usr/include/armadillo_bits/fn_solve.hpp:139:1: note: template<class T1, class T2> bool arma::solve(arma::Mat<typename T1::elem_type>&, const arma::Base<typename T1::elem_type, T1>&, const arma::Base<typename T1::elem_type, T2>&, const char*, const typename arma::arma_blas_type_only<typename T1::elem_type>::result*)
 solve
 ^
/usr/include/armadillo_bits/fn_solve.hpp:139:1: note:   template argument deduction/substitution failed:
../src/FEM.cpp:37:22: note:   ‘arma::sp_mat aka arma::SpMat<double>’ is not derived from ‘arma::Mat<typename T1::elem_type>’
  vec zeta = solve(A,b);
                      ^
make: *** [src/FEM.o] Error 1

【问题讨论】:

【参考方案1】:

Armadillo 中的稀疏矩阵支持尚未完成。

您可以使用 ARPACK 对稀疏矩阵进行(特征)分解。稀疏矩阵求解器可能会在下一个版本中出现,它可能会使用 SuiteSparse 项目中的 CXSparse 库。

【讨论】:

【参考方案2】:

从 5.0 版开始,Armadillo 具有用于求解稀疏系统的 spsolve() 函数。

【讨论】:

以上是关于犰狳:用 sp_mat 解决的主要内容,如果未能解决你的问题,请参考以下文章

将犰狳中的矩阵从稀疏转换为密集(spmat 到 mat)

使用 Rcpp 通过引用传递犰狳稀疏矩阵

使用 Rcpp 代码访问和修改 arma::sp_mat 类稀疏矩阵的非零元素

如何在犰狳中按元素划分2个稀疏矩阵?

犰狳错误:没有匹配函数调用‘inv(arma::SpMat<double>&)’

犰狳:解决 Ax=b 分配堆?