ubuntu 14.04 下安装Armadillo库 进行稀疏矩阵特征分解

图像处理中往往会遇到大型稀疏矩阵的特征值分解问题，如果采用一般的密集矩阵求解方法，对内存要求比较大。所以在稀疏矩阵的基础上进行特征值分解，更加合理。opencv2.49版本中还没有对稀疏矩阵的特征值分解函数，所以自己对此类问题进行了网上搜索。调研了ietl、lapack、Armadillo、Eigen这类库。文章先对各种库进行总结，后面在详细的介绍Armadillo库在ubuntu14.04系统上的安装方法，最后给出一个利用Armadillo进行稀疏矩阵特征值分解的例子。

1. 各种库的总结

ietl

网站地址：

http://www.comp-phys.org/software/ietl/

该库可以解决稀疏矩阵的特征值分解问题，但是其目前其不经已经没有更新，而且在求解多个特征值时，只有第一个值是正确的，其余的均为0. 这个是其自身算法稳定性问题，所以自己舍弃了这一个工具包

lapack

网站地址：

http://www.netlib.org/lapack/#_lapack_version_3_5_0

该库解决了密集矩阵的分解问题，但是稀疏矩阵的分解并不是它的优势。

Armadillo

网站地址：

http://arma.sourceforge.net/docs.html

该库目前还在更新，并且函数仿照matlab使用习惯，使用方便，容易上手，并且其对稀疏矩阵的运算支持的很全面。它既支持密集矩阵的求解也支持稀疏矩阵的求解。密集矩阵求解主要基于lapack和blas等库，稀疏矩阵求解主要基于superlu和arpack等库。

Eigen

网站地址：

http://eigen.tuxfamily.org/dox/index.html

该库目前也在更新中，其对密集矩阵的支持比较全面，但是对稀疏矩阵的支持比较欠缺，其没有针对稀疏矩阵特征值分解的函数，有关于LU、和QR分解的功能。

2.Armadillo 库的安装

Armadillo库依赖LAPACK、 BLAS、ARPACK 、SuperLU、ATLAS、OpenBLAS等库。如果只需要对密集矩阵进行处理，安装LAPACK、 BLAS库即可，如果需要对稀疏矩阵进行求解则需要包括ARPACK 、SuperLU。OpenBLAS是BLAS的进化版本，它支持多线程，所以为了提高速度可以采用OpenBLAS代替BLAS。

1) LAPACK包和BLAS包的安装

因为LAPACK包安装时可以同时生成BLAS包，所以这两个一起安装。这里使用的是3.5版本。

将make.inc.example文件改名为make.inc，对该文件进行修改，修改如下，加上lapackelib这一项，这个库后面可能会用到。修改后会生成，执行make 命令生成libblas.a liblapack.a liblapacke.a等库文件，然后将生成的所有库都复制到/usr/local/lib中。

2)ATLAS包的安装

阅读INSTALL.txt文件，进行软件安装

一)关闭CPU THROTTLING功能

切换到root模式下，执行下面的命令，如果提示缺少命令包，可以根据提示安装。然后在root模式下执行下面的命令，关闭CPU THROTTLING功能。

cpufreq-selector -g performance

二).编译库文件

 在解压后的安装包根目录下建立building文件夹，cd到此目录下，然后执行../configure

安装的时候 要最后执行make install 把产生的头文件，库文件分别复制到/usr/local/lib 和/usr/local/include中。 安装的时候 要最后执行make install将头文件和库文件进行打包，并通过命令行 把产生的头文件，库文件分别复制到/usr/local/lib 和/usr/local/include中

3)ARPACK包的安装

按照README内容对ARmake.inc进行修改
这是对用到的源文件路径进行设置。

编译器部分设置如下，主要是“FFLAGS”和“MAKE”选项的设置。

修改完成后可以直接make lib进行编译
编译完成后，将libarpack.a复制到 /usr/local/lib目录中

4)openblas包安装

按照 readme.md文件进行安装。
直接make编译，make install安装。安装文件默认在./OpenBLAS-0.2.14/lib目录中.
./OpenBLAS-0.2.14/lib/lib 中所有库文件复制到/usr/local/lib目录中，./OpenBLAS-0.2.14/lib/include中的所有文件，复制到/usr/local/include/openblas目录下

5)superlu包安装

根据readme文档进行安装，这里采用4.3版本，因为Armadillo说明文档中提示只支持这一版本。

复制MAKE_INC文件夹下的make.linux文件到SuperLU_4.3目录下，修改名称make.inc，并对其进行修改，因为我的blas库在/usr/loval/lib目录中所以BLASLIB设置为该值。blas库可以在安装lapack时得到。

修改完成后，执行make文件。

superlu库函数，安装后将 superlu.a文件复制到/usr/local/lib下，并且把SRC中的.h文件复制到/usr/local/include/superlu目录下，这样后面 Armadillo 整体安装时才能成功检测到该库函数。
因为暂时用不到matlab相关的接口，所以暂时到此结束。

6) libf2c.a库的安装

如果系统中没有libf2c.a这一个库函数，需要自己编译，并将其复制到/usr/local/lib目录中。在下面的网址中下载libf2c.zip,然后按照readme文件进行安装。

http://www.netlib.org/f2c/

3. 使用Armadillo例子

1) 基本例子

这个是Armadillo自带的一个例子。利用下面的命令进行编译，我的头文件放在目录/home/xuehen/armadillo-6.100.0/include中，这里只是对密集矩阵进行求解工作，所以只需要liblapack.a libblas.a 库文件即可。-lgfortran是加上gortran的库。

g++ example1.cpp -o example1 -O2 -I/home/xuehen/armadillo-6.100.0/include -DARMA_DONT_USE_WRAPPER -llapack -lblas -lgfortran

#include <iostream>
#include <armadillo>

using namespace std;
using namespace arma;

// Armadillo documentation is available at:
// http://arma.sourceforge.net/docs.html

int
main(int argc, char** argv)
  {
  cout << "Armadillo version: " << arma_version::as_string() << endl;
  
  mat A(2,3);  // directly specify the matrix size (elements are uninitialised)
  
  cout << "A.n_rows: " << A.n_rows << endl;  // .n_rows and .n_cols are read only
  cout << "A.n_cols: " << A.n_cols << endl;
  
  A(1,2) = 456.0;  // directly access an element (indexing starts at 0)
  A.print("A:");
  
  A = 5.0;         // scalars are treated as a 1x1 matrix
  A.print("A:");
  
  A.set_size(4,5); // change the size (data is not preserved)
  
  A.fill(5.0);     // set all elements to a particular value
  A.print("A:");
  
  // endr indicates "end of row"
  A << 0.165300 << 0.454037 << 0.995795 << 0.124098 << 0.047084 << endr
    << 0.688782 << 0.036549 << 0.552848 << 0.937664 << 0.866401 << endr
    << 0.348740 << 0.479388 << 0.506228 << 0.145673 << 0.491547 << endr
    << 0.148678 << 0.682258 << 0.571154 << 0.874724 << 0.444632 << endr
    << 0.245726 << 0.595218 << 0.409327 << 0.367827 << 0.385736 << endr;
  
  A.print("A:");
  
  // determinant
  cout << "det(A): " << det(A) << endl;
  
  // inverse
  cout << "inv(A): " << endl << inv(A) << endl;
  
  // save matrix as a text file
  A.save("A.txt", raw_ascii);
  
  // load from file
  mat B;
  B.load("A.txt");
  
  // submatrices
  cout << "B( span(0,2), span(3,4) ):" << endl << B( span(0,2), span(3,4) ) << endl;
  
  cout << "B.row(0): " << endl << B.row(0) << endl;
  
  cout << "B.col(1): " << endl << B.col(1) << endl;
  
  // transpose
  cout << "B.t(): " << endl << B.t() << endl;
  
  // maximum from each column (traverse along rows)
  cout << "max(B): " << endl << max(B) << endl;
  
  // maximum from each row (traverse along columns)
  cout << "max(B,1): " << endl << max(B,1) << endl;
  
  // maximum value in B
  cout << "max(max(B)) = " << max(max(B)) << endl;
  
  // sum of each column (traverse along rows)
  cout << "sum(B): " << endl << sum(B) << endl;
  
  // sum of each row (traverse along columns)
  cout << "sum(B,1) =" << endl << sum(B,1) << endl;
  
  // sum of all elements
  cout << "accu(B): " << accu(B) << endl;
  
  // trace = sum along diagonal
  cout << "trace(B): " << trace(B) << endl;
  
  // generate the identity matrix
  mat C = eye<mat>(4,4);
  
  // random matrix with values uniformly distributed in the [0,1] interval
  mat D = randu<mat>(4,4);
  D.print("D:");
  
  // row vectors are treated like a matrix with one row
  rowvec r;
  r << 0.59119 << 0.77321 << 0.60275 << 0.35887 << 0.51683;
  r.print("r:");
  
  // column vectors are treated like a matrix with one column
  vec q;
  q << 0.14333 << 0.59478 << 0.14481 << 0.58558 << 0.60809;
  q.print("q:");
  
  // convert matrix to vector; data in matrices is stored column-by-column
  vec v = vectorise(A);
  v.print("v:");
  
  // dot or inner product
  cout << "as_scalar(r*q): " << as_scalar(r*q) << endl;
  
  // outer product
  cout << "q*r: " << endl << q*r << endl;
  
  // multiply-and-accumulate operation (no temporary matrices are created)
  cout << "accu(A % B) = " << accu(A % B) << endl;
  
  // example of a compound operation
  B += 2.0 * A.t();
  B.print("B:");
  
  // imat specifies an integer matrix
  imat AA;
  imat BB;
  
  AA << 1 << 2 << 3 << endr << 4 << 5 << 6 << endr << 7 << 8 << 9;
  BB << 3 << 2 << 1 << endr << 6 << 5 << 4 << endr << 9 << 8 << 7;
  
  // comparison of matrices (element-wise); output of a relational operator is a umat
  umat ZZ = (AA >= BB);
  ZZ.print("ZZ:");
  
  // cubes ("3D matrices")
  cube Q( B.n_rows, B.n_cols, 2 );
  
  Q.slice(0) = B;
  Q.slice(1) = 2.0 * B;
  
  Q.print("Q:");
  
  // 2D field of matrices; 3D fields are also supported
  field<mat> F(4,3); 
  
  for(uword col=0; col < F.n_cols; ++col)
  for(uword row=0; row < F.n_rows; ++row)
    {
    F(row,col) = randu<mat>(2,3);  // each element in field<mat> is a matrix
    }
  
  F.print("F:");
  
  return 0;
  }

2). 求解对称稀疏实矩阵的例子

g++ example.cpp -o example -O2 -I/home/xuehen/armadillo-6.100.0/include /usr/local/lib/libarpack.a /usr/local/lib/liblapack.a /usr/local/lib/libopenblas.a /usr/local/lib/libf2c.a  -lgfortran

或者：
g++ example.cpp -o example -O2 -I/home/xuehen/armadillo-6.100.0/include -L/usr/local/lib -larpack -llapack -lopenblas -lf2c  -lgfortran

对其进行编译，由于是稀疏矩阵求解所以要加上 libarpack.a库。 这里才用了libopenblas.a替代了libblas.a, libf2c.a这一个库是必须的否则会出现链接错误。

#include <iostream>
#include <armadillo>

using namespace std;
using namespace arma;

// Armadillo documentation is available at:
// http://arma.sourceforge.net/docs.html

int
main(int argc, char** argv)
{
// generate sparse matrix
sp_fmat A = sprandu<sp_fmat>(1000, 1000, 0.1);
sp_fmat B = A.t()*A;

fvec eigval;
fmat eigvec;

eigs_sym(eigval, eigvec, B, 1);  // find 5 eigenvalues/eigenvectors
}

4.总结

由于尝试对Armadillo打包时出现链接错误，所以这里才用了不安装的模式调用Armadillo。