Eigen使用

Eigen的安装参见博客.

1.开始小例子

#include <iostream>
#include <Dense>
using Eigen::MatrixXd;
int main()
{
  MatrixXd m(2,2);
  m(0,0) = 3;
  m(1,0) = 2.5;
  m(0,1) = -1;
  m(1,1) = m(1,0) + m(0,1);
  std::cout << m << std::endl;
}

在testeigen文件夹下面建立main.cpp文件,将以上内容考入其中,并将已经解压缩的eigen-3.3.7文件夹下面的Eigen文件夹整个的考到testeigen文件夹下面,最后进行编译:

g++ -I ./Eigen  main.cpp -o main.o

执行

./main.o

可以输出

  3  -1
2.5 1.5

2.Eigen与Cmake

文件结构为

|--testEigen
|  -- src
|     -- main.cpp
|  -- Eigen
|  -- CMakeLists.txt
|  -- build

其中CMakeLists.txt中内容为

cmake_minimum_required (VERSION 2.6)
project (TEST)
set (TEST_VERSION 0.1)
set(CMAKE_BUILD_TYPE "Debug")
set(CMAKE_CXX_FLAGS_DEBUG "$ENV{CXXFLAGS} -O0 -Wall -g -ggdb -DDEBUG")
set(CMAKE_CXX_FLAGS_RELEASE "$ENV{CXXFLAGS} -O3 -Wall")
set(EXECUTABLE_OUTPUT_PATH ${PROJECT_BINARY_DIR}/../bin)
aux_source_directory(${PROJECT_SOURCE_DIR}/src DIR_SRC)
include_directories(${PROJECT_SOURCE_DIR}/src)
include_directories(${PROJECT_SOURCE_DIR}/Eigen)
#link_directories(${PROJECT_SOURCE_DIR}/lib)
add_executable(main ${DIR_SRC})

3.使用压缩列存储直接初始化矩阵

#include <Sparse>
#include <vector>
#include <iostream>
typedef Eigen::SparseMatrix<double> SpMat; // declares a column-major sparse matrix type of double
int main(int argc, char** argv)
{
    int rows = 5,cols = 5,nnz = 8;
    double values[8] = {22,  7, 3, 5, 14,   1, 17,   8};
    int InnerIndices[8] = {1,   2, 0, 2, 4, 2, 1, 4};
    int OuterStarts[6] = {0, 2, 4, 5, 6, 8};
    Eigen::Map<SpMat> sm1(rows,cols,nnz,OuterStarts,InnerIndices,values);
    std::cout<< sm1;
    return 0;
}

输出结果:

0 3 0 0 0 
22 0 0 0 17 
7 5 0 1 0 
0 0 0 0 0 
0 0 14 0 8 

4.Eigen求解系数为对称矩阵的线性方程组

#include <Sparse>
#include <vector>
#include <iostream>
using namespace Eigen;
using namespace std;
typedef Eigen::SparseMatrix<double> SpMat; // declares a column-major sparse matrix type of double
int main(int argc, char** argv)
{
  const int N = 10,ANZ = 19;
  int Ap[N+1] = {0, 1, 2, 3, 4, 6, 7, 9, 11, 15, ANZ};
  int Ai[ANZ] = {0, 1, 2, 3, 1,4, 5, 4,6, 4,7, 0,4,7,8, 1,4,6,9};
  double Ax[ANZ] = {1.7, 1., 1.5, 1.1, .02,2.6, 1.2, 
              .16,1.3, .09,1.6,.13,.52,.11,1.4, .01,.53,.56,3.1};
  VectorXd b(N), x(N);
  b(0) = .287;b(1) =.22;b(2)=.45;b(3)=.44;b(4)=2.486;
  b(5)=.72;b(6)=1.55;b(7)=1.424;b(8)=1.621;b(9)=3.759; 
  Eigen::Map<SpMat,RowMajor> sm2(N,N,ANZ,Ap,Ai,Ax);//注意这里是压缩行存储,因为是对称矩阵,只存了上半部分
  std::cout<< sm2;
  Eigen::SimplicialLLT<SpMat,Upper > solver;//声明这是矩阵的上半部分
  solver.analyzePattern(sm2);   
  solver.factorize(sm2);
  if(solver.info()!=Success) {
     cout<<" decomposition failed!"<<endl;
      return -1;
  }
  x =  solver.solve(b);
  if(solver.info()!=Success) {
     cout<<" solve failed!"<<endl;
      return -1;
  }
  cout<< x;
  return 0;
}

这里我们使用了求解器SimplicialLLT,我们也可以使用其他的求解器,例如PastixLLT,PastixLDLT,PastixLU等,但是需要包含相应的包与头文件。输出结果如下:

1.7 0 0 0 0 0 0 0 0.13 0 
0 1 0 0 0.02 0 0 0 0 0.01 
0 0 1.5 0 0 0 0 0 0 0 
0 0 0 1.1 0 0 0 0 0 0 
0 0 0 0 2.6 0 0.16 0.09 0.52 0.53 
0 0 0 0 0 1.2 0 0 0 0 
0 0 0 0 0 0 1.3 0 0 0.56 
0 0 0 0 0 0 0 1.6 0.11 0 
0 0 0 0 0 0 0 0 1.4 0 
0 0 0 0 0 0 0 0 0 3.1 
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1

利用压缩存储构造稀疏矩阵也可以采用如下方式:

//typedef Eigen::SparseMatrix<double> SpMat;
SpMat A;
At = Map<SpMat,RowMajor>(n,m,nnz,Ap,Ai,Ax);

5.Eigen中的对角矩阵

参考链接
对角矩阵的使用有两种方式。第一种是声明一个对角矩阵,然后用向量来赋值,第二种是把向量看做对角矩阵
首先在上面文件开头加上

typedef Eigen::DiagonalMatrix<double,Dynamic> DiagMat;

结尾处加上

DiagMat S(N);
VectorXd c(N);
for(int i=0;i<N;i++)
    c(i) = i+1;
S.diagonal() = c;
cout<< S*sm2<<endl;//方式一
cout<< c.asDiagonal()*sm2<<endl;//方式二

输出结果

1.7 0 0 0 0 0 0 0 0.13 0 
0 2 0 0 0.04 0 0 0 0 0.02 
0 0 4.5 0 0 0 0 0 0 0 
0 0 0 4.4 0 0 0 0 0 0 
0 0 0 0 13 0 0.8 0.45 2.6 2.65 
0 0 0 0 0 7.2 0 0 0 0 
0 0 0 0 0 0 9.1 0 0 3.92 
0 0 0 0 0 0 0 12.8 0.88 0 
0 0 0 0 0 0 0 0 12.6 0 
0 0 0 0 0 0 0 0 0 31 

1.7 0 0 0 0 0 0 0 0.13 0 
0 2 0 0 0.04 0 0 0 0 0.02 
0 0 4.5 0 0 0 0 0 0 0 
0 0 0 4.4 0 0 0 0 0 0 
0 0 0 0 13 0 0.8 0.45 2.6 2.65 
0 0 0 0 0 7.2 0 0 0 0 
0 0 0 0 0 0 9.1 0 0 3.92 
0 0 0 0 0 0 0 12.8 0.88 0 
0 0 0 0 0 0 0 0 12.6 0 
0 0 0 0 0 0 0 0 0 31 

Author: zrjer
Reprint policy: All articles in this blog are used except for special statements CC BY 4.0 reprint polocy. If reproduced, please indicate source zrjer !
 Previous
开源求解器汇总 开源求解器汇总
线性规划问题 name 官网 github地址 语言 方法 HiGHS 链接 链接 c++ Simplex/Interior Point method PCx 链接 链接 c Interior Point method C
2020-02-26 zrjer
Current 
Eigen使用 Eigen使用
Eigen的安装参见博客. 1.开始小例子#include <iostream> #include <Dense> using Eigen::MatrixXd; int main() { MatrixXd m(2
2020-02-26 zrjer
  TOC