目錄
- 介紹
- 安裝
- Demo
- 矩陣、向量初始化
- C++數組和矩陣轉換
- 矩陣基礎操作
- 點積和叉積
- 轉置、伴隨、行列式、逆矩陣
- 計算特征值和特征向量
- 解線性方程
- 最小二乘求解
- 稀疏矩陣
-
-
介紹
Eigen是一個輕量級的矩陣庫,除了稀疏矩陣不成熟(3.1有較大改進)以外,其他的矩陣和向量操作都比較完善,而且速度不錯.
安裝
在eigen 3.1.3下載最新的版本,然后解壓文件,將解壓出來的文件夾下的\\Eigen\文件夾拷貝到程序文件夾下,包括頭文件,即可使用
Demo
eigendemo.zip
示例是vs2010環境下的程序,主要的文件就只有main.cpp和Eigen文件夾。
矩陣、向量初始化
#include <iostream> #include "Eigen/Dense" using namespace Eigen; int main() {MatrixXf m1(3,4); //動態矩陣,建立3行4列。MatrixXf m2(4,3); //4行3列,依此類推。MatrixXf m3(3,3);Vector3f v1; //若是靜態數組,則不用指定行或者列/* 初始化 */Matrix3d m = Matrix3d::Random();m1 = MatrixXf::Zero(3,4); //用0矩陣初始化,要指定行列數m2 = MatrixXf::Zero(4,3);m3 = MatrixXf::Identity(3,3); //用單位矩陣初始化v1 = Vector3f::Zero(); //同理,若是靜態的,不用指定行列數m1 << 1,0,0,1, //也可以以這種方式初始化1,5,0,1,0,0,9,1;m2 << 1,0,0,0,4,0,0,0,7,1,1,1;//向量初始化,與矩陣類似Vector3d v3(1,2,3);VectorXf vx(30); }
C++數組和矩陣轉換
使用Map函數,可以實現Eigen的矩陣和c++中的數組直接轉換,語法如下:
//@param MatrixType 矩陣類型 //@param MapOptions 可選參數,指的是指針是否對齊,Aligned, or Unaligned. The default is Unaligned. //@param StrideType 可選參數,步長 /*Map<typename MatrixType,int MapOptions,typename StrideType> */int i;//數組轉矩陣double *aMat = new double[20];for(i =0;i<20;i++){aMat[i] = rand()%11;}//靜態矩陣,編譯時確定維數 Matrix<double,4,5> Eigen:Map<Matrix<double,4,5> > staMat(aMat);//輸出for (int i = 0; i < staMat.size(); i++)std::cout << *(staMat.data() + i) << " ";std::cout << std::endl << std::endl;//動態矩陣,運行時確定 MatrixXdMap<MatrixXd> dymMat(aMat,4,5);//輸出,應該和上面一致for (int i = 0; i < dymMat.size(); i++)std::cout << *(dymMat.data() + i) << " ";std::cout << std::endl << std::endl;//Matrix為列優先,如下返回指針dymMat.data();
矩陣基礎操作
eigen重載了基礎的+ - * / += -= *= /= *可以表示標量和矩陣或者矩陣和矩陣
#include <iostream> #include <Eigen/Dense> using namespace Eigen; int main() {//單個取值,單個賦值double value00 = staMat(0,0);double value10 = staMat(1,0);staMat(0,0) = 100;std::cout << value00 <<value10<<std::endl;std::cout <<staMat<<std::endl<<std::endl;//加減乘除示例 Matrix2d 等同于 Matrix<double,2,2>Matrix2d a;a << 1, 2,3, 4;MatrixXd b(2,2);b << 2, 3,1, 4;Matrix2d c = a + b;std::cout<< c<<std::endl<<std::endl;c = a - b;std::cout<<c<<std::endl<<std::endl;c = a * 2;std::cout<<c<<std::endl<<std::endl;c = 2.5 * a;std::cout<<c<<std::endl<<std::endl;c = a / 2;std::cout<<c<<std::endl<<std::endl;c = a * b;std::cout<<c<<std::endl<<std::endl;
點積和叉積
#include <iostream> #include <Eigen/Dense> using namespace Eigen; using namespace std; int main() {//點積、叉積(針對向量的)Vector3d v(1,2,3);Vector3d w(0,1,2);std::cout<<v.dot(w)<<std::endl<<std::endl;std::cout<<w.cross(v)<<std::endl<<std::endl; } */
轉置、伴隨、行列式、逆矩陣
小矩陣(4 * 4及以下)eigen會自動優化,默認采用LU分解,效率不高
#include <iostream> #include <Eigen/Dense> using namespace std; using namespace Eigen; int main() {Matrix2d c;c << 1, 2,3, 4;//轉置、伴隨std::cout<<c<<std::endl<<std::endl;std::cout<<"轉置\n"<<c.transpose()<<std::endl<<std::endl;std::cout<<"伴隨\n"<<c.adjoint()<<std::endl<<std::endl;//逆矩陣、行列式std::cout << "行列式: " << c.determinant() << std::endl;std::cout << "逆矩陣\n" << c.inverse() << std::endl; }
計算特征值和特征向量
#include <iostream> #include <Eigen/Dense> using namespace std; using namespace Eigen; int main() {//特征向量、特征值std::cout << "Here is the matrix A:\n" << a << std::endl;SelfAdjointEigenSolver<Matrix2d> eigensolver(a);if (eigensolver.info() != Success) abort();std::cout << "特征值:\n" << eigensolver.eigenvalues() << std::endl;std::cout << "Here's a matrix whose columns are eigenvectors of A \n"<< "corresponding to these eigenvalues:\n"<< eigensolver.eigenvectors() << std::endl; }
解線性方程
#include <iostream> #include <Eigen/Dense> using namespace std; using namespace Eigen; int main() {//線性方程求解 Ax =B;Matrix4d A;A << 2,-1,-1,1,1,1,-2,1,4,-6,2,-2,3,6,-9,7;Vector4d B(2,4,4,9);Vector4d x = A.colPivHouseholderQr().solve(B);Vector4d x2 = A.llt().solve(B);Vector4d x3 = A.ldlt().solve(B); std::cout << "The solution is:\n" << x <<"\n\n"<<x2<<"\n\n"<<x3 <<std::endl; }
除了colPivHouseholderQr、LLT、LDLT,還有以下的函數可以求解線性方程組,請注意精度和速度:?解小矩陣(4*4)基本沒有速度差別
Decomposition | Method | 矩陣特殊要求 | 速度 | 精度 |
PartialPivLU | partialPivLu() | 可逆 | ++ | + |
FullPivLU | fullPivLu() | None | - | +++ |
HouseholderQR | householderQr() | None | ++ | + |
ColPivHouseholderQR | colPivHouseholderQr() | None | + | ++ |
FullPivHouseholderQR | fullPivHouseholderQr() | None | - | +++ |
LLT | llt() | 正定 | +++ | + |
LDLT | ldlt() | 正或負半定 Positive or negative semidefinite | +++ | ++ |
最小二乘求解
最小二乘求解有兩種方式,jacobiSvd或者colPivHouseholderQr,4*4以下的小矩陣速度沒有區別,jacobiSvd可能更快,大矩陣最好用colPivHouseholderQr
#include <iostream> #include <Eigen/Dense> using namespace std; using namespace Eigen; int main() {MatrixXf A1 = MatrixXf::Random(3, 2);std::cout << "Here is the matrix A:\n" << A1 << std::endl;VectorXf b1 = VectorXf::Random(3);std::cout << "Here is the right hand side b:\n" << b1 << std::endl;//jacobiSvd 方式:Slow (but fast for small matrices)std::cout << "The least-squares solution is:\n"<< A1.jacobiSvd(ComputeThinU | ComputeThinV).solve(b1) << std::endl;//colPivHouseholderQr方法:faststd::cout << "The least-squares solution is:\n"<< A1.colPivHouseholderQr().solve(b1) << std::endl; }
稀疏矩陣
稀疏矩陣的頭文件包括:
#include <Eigen/SparseCore>#include <Eigen/SparseCholesky>#include <Eigen/IterativeLinearSolvers>#include <Eigen/Sparse>
初始化有兩種方式: 1.使用三元組插入
typedef Eigen::Triplet<double> T; std::vector<T> tripletList; triplets.reserve(estimation_of_entries); //estimation_of_entries是預估的條目 for(...) {tripletList.push_back(T(i,j,v_ij));//第 i,j個有值的位置的值 } SparseMatrixType mat(rows,cols); mat.setFromTriplets(tripletList.begin(), tripletList.end()); // mat is ready to go!
2.直接將已知的非0值插入
SparseMatrix<double> mat(rows,cols); mat.reserve(VectorXi::Constant(cols,6)); for(...) {// i,j 個非零值 v_ij != 0mat.insert(i,j) = v_ij; } mat.makeCompressed(); // optional
稀疏矩陣支持大部分一元和二元運算:
sm1.real() sm1.imag() -sm1 0.5*sm1 sm1+sm2 sm1-sm2 sm1.cwiseProduct(sm2)
二元運算中,稀疏矩陣和普通矩陣可以混合使用
//dm表示普通矩陣 dm2 = sm1 + dm1;
也支持計算轉置矩陣和伴隨矩陣