Eigen 矩阵计算工具

  • 1. 源码
  • 2. 矩阵的定义
    • 2.1. 模板函数
    • 2.2. 动态矩阵
    • 2.3. 静态矩阵
    • 2.4. 构造函数
  • 3. 元素访问和设置
  • 4. 重置大小
  • 5. 矩阵运算
    • 5.1. 加减运算符
    • 5.2. 标量乘除运算符
    • 5.3. 转置、共轭和伴随矩阵
    • 5.4. 矩阵乘法
    • 5.5. 点乘和叉乘运算
      • 5.5.1. 点乘
      • 5.5.2. 叉乘
      • 5.5.3. 点到线的距离
    • 5.6. reduction运算
    • 5.7. 分解函数
      • 5.7.1. QR分解
      • 5.7.2. 矩阵求逆
      • 5.7.3. 计算数值法求解和真实值的残差
      • 5.7.4. 最小二乘解

1. 源码

Eigen是C++开源矩阵计算工具
可以从官网下载的源码包的eigen-x.x.x\Eigen解压到目标目录下就可以了


Eigen用源码的方式提供给用户使用,在使用时只需要包含Eigen的头文件即可进行使用

g++ -I${PATH_TO_EIGEN} example.cpp -o example
g++ -I$./ example.cpp -o example /* 同一目录下时 */

之所以采用这种方式,是因为Eigen采用模板方式实现
由于模板函数不支持分离编译,所以只能提供源码而不是动态库的方式供用户使用


2. 矩阵的定义


2.1. 模板函数

emplate<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
  • _Scalar 矩阵元素的类型
  • _Rows 矩阵行数
  • _Cols 矩阵列数

矩阵定义时可以使用Dynamic 来表示矩阵的行列数是未知,例如:

typedef matrix<double, Dynamic, Dynamic> MatrixXd;

在Eigen中也提供了很多常见的简化定义形式,例如

typedef Matrix<Scalar,3,3> Matrix3;
typedef Matrix<Scalar,3,1> Vector3;

无论是矩阵还是数组、向量,无论是静态矩阵还是动态矩阵都提供默认构造函数
也就是定义这些数据结构时都可以不用提供任何参数,其大小均由运行时来确定


2.2. 动态矩阵

动态矩阵是指其大小在运行时确定 MatrixXd VectorXd

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;int main()
{MatrixXd m = MatrixXd::Random(3, 3);m = (m + MatrixXd::Constant(3, 3, 1.0)) * 50;cout << "m =" << endl<< m << endl;VectorXd v(3);v << 1, 2, 3;cout << "m * v =" << endl<< m * v << endl;
}

2.3. 静态矩阵

静态矩阵是指其大小在编译时确定 Matrix3d Vector3d

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;int main()
{Matrix3d m = Matrix3d::Random();m = (m + Matrix3d::Constant(1.0)) * 50;cout << "m =" << endl<< m << endl;Vector3d v(1, 2, 3);cout << "m * v =" << endl<< m * v << endl;
}

对于小矩阵(一般大小小于16)的使用固定大小的静态矩阵,它可以带来比较高的效率
对于大矩阵(一般大小大于32)建议使用动态矩阵 特别大的矩阵使用了固定大小的静态矩阵则可能造成栈溢出的问题


2.4. 构造函数

Matrix3f a;       /* 不分配内存,不初始化 */
MatrixXf a(1, 2); /* 配内存,不初始化 */Vector2d a(2.0, 3.0); /* 向量,eigen默认列优先 */
RowVector2d a;        /* 定义行向量 */Matrix3f b;
b << 1, 2, 3, 4, 5, 6, 7, 8, 9;/* 常用矩阵初始化 */
MatrixXf::Zero(3, 4); /* 置零 */
MatrixXf::Ones(3, 3); /* 置一 */
Vector3f::Ones();
MatrixXi::Identity(5, 5); /* 单位阵 */
Matrix3d::Random();       /* 随机阵 */

3. 元素访问和设置

在矩阵的访问中,行索引总是作为第一个参数
需注意Eigen中遵循大家的习惯让矩阵、数组、向量的下标都是从0开始
矩阵元素的访问可以通过 () 操作符完成

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;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 << "Here is the matrix m:\n"<< m << std::endl;VectorXd v(2);v(0) = 4;v(1) = v(0) - 1;std::cout << "Here is the vector v:\n"<< v << std::endl;
}// Here is the matrix m:
//   3  -1
// 2.5 1.5
// Here is the vector v:
// 4
// 3

Vector 还可以用x()y() 访问

在Eigen中重载了"<<"操作符,通过该操作符即可以一个一个元素的进行赋值,也可以一块一块的赋值

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{Matrix3f m;m << 1, 2, 3,4, 5, 6,7, 8, 9;std::cout << m;
}
// 1 2 3
// 4 5 6
// 7 8 9

4. 重置大小

当前矩阵的行数、列数、大小可以通过rows()cols()size()来获取
对于动态矩阵可以通过 resize() 函数来动态修改矩阵的大小

需注意:

  1. 固定大小的矩阵是不能使用resize()来修改矩阵的大小
  2. resize()函数会析构掉原来的数据,因此调用resize()函数之后将不能保证元素的值不改变
  3. 使用“=”操作符操作动态矩阵时,如果左右边的矩阵大小不等,则左边的动态矩阵的大小会被修改为右边的大小
#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{MatrixXf a(2, 2);std::cout << "a is of size " << a.rows() << "x" << a.cols() << std::endl;MatrixXf b(3, 3);a = b;std::cout << "a is now of size " << a.rows() << "x" << a.cols() << std::endl;
}// a is of size 2x2
// a is now of size 3x3

5. 矩阵运算


5.1. 加减运算符

根据矩阵加减法的运算规则,知道参与运算的两个矩阵应当具有相同的行数和列数
此外,由于C++的语法要求,两个Matrix对象的元素应当具有相同的数据类型

Eigen提供了如下的运算符:

  • ‘+’, 二元运算符,a + b
  • ‘-’, 二元运算符,a - b
  • ‘-’, 单元运算符,-a
  • ‘+=’ 复合运算符,a += b
  • ‘-=’ 复合运算符,a -= b
#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{Matrix2d a;MatrixXd b(2, 2);a << 1, 2, 3, 4;b << 2, 3, 1, 4;std::cout << "a + b =\n"<< a + b << std::endl;// a + b =
// 3 5
// 4 8std::cout << "a - b =\n"<< a - b << std::endl;// a - b =
// -1 -1
//  2  0std::cout << "Doing a += b;" << std::endl;a += b;std::cout << "Now a =\n"<< a << std::endl;// Now a =
// 3 5
// 4 8Vector3d v(1, 2, 3);Vector3d w(1, 0, 0);std::cout << "-v + w - v =\n"<< -v + w - v << std::endl;// -v + w - v =
// -1
// -4
// -6
}

5.2. 标量乘除运算符

  • ’ * ', 二元运算符,matrix * scalar
  • ’ * ', 二元运算符,scalar * matrix
  • ’ / ', 二元运算符,matrix / scalar
  • ‘*=’, 复合运算符,matrix *= scalar
  • ‘/=’, 复合运算符,matrix /= scalar
#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{Matrix2d a;a << 1, 2, 3, 4;std::cout << "a * 2.5 =\n"<< a * 2.5 << std::endl;Vector3d v(1, 2, 3);std::cout << "0.1 * v =\n"<< 0.1 * v << std::endl;std::cout << "Doing v *= 2;" << std::endl;v *= 2;std::cout << "Now v =\n"<< v << std::endl;
}// a * 2.5 =
// 2.5   5
// 7.5  10// 0.1 * v =
// 0.1
// 0.2
// 0.3// Doing v *= 2;
// Now v =
// 2
// 4
// 6

5.3. 转置、共轭和伴随矩阵

矩阵a的转置aT,通过成员函数transpose()获得
矩阵a的共轭a¯,通过成员函数conjugate()获得
矩阵a的伴随a,通过成员函数adjoint()获得

例如复数矩阵

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{MatrixXcf a = MatrixXcf::Random(2,2);std::cout << "Here is the matrix a\n" << a << std::endl;std::cout << "Here is the matrix a^T\n" << a.transpose() << std::endl;std::cout << "Here is the conjugate of a^¯\n" << a.conjugate() << std::endl;std::cout << "Here is the matrix a^*\n" << a.adjoint() << std::endl;
}// Here is the matrix a
// (-0.211234,0.680375) (-0.604897,0.823295)
//   (0.59688,0.566198) (0.536459,-0.329554)// Here is the matrix a^T
// (-0.211234,0.680375)   (0.59688,0.566198)
// (-0.604897,0.823295) (0.536459,-0.329554)// Here is the conjugate of a^¯
// (-0.211234,-0.680375) (-0.604897,-0.823295)
//   (0.59688,-0.566198)   (0.536459,0.329554)// Here is the matrix a^*
// (-0.211234,-0.680375)   (0.59688,-0.566198)
// (-0.604897,-0.823295)   (0.536459,0.329554)

对于实数矩阵,共轭就是它本身所以 conjugate() 不做任何操作,而 adjoint()transpose() 等价

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{MatrixXd a = MatrixXd::Random(2,2);std::cout << "Here is the matrix a\n" << a << std::endl;std::cout << "Here is the matrix a^T\n" << a.transpose() << std::endl;std::cout << "Here is the conjugate of a^¯\n" << a.conjugate() << std::endl;std::cout << "Here is the matrix a^*\n" << a.adjoint() << std::endl;
}// Here is the matrix a
//  0.680375  0.566198
// -0.211234   0.59688// Here is the matrix a^T
//  0.680375 -0.211234
//  0.566198   0.59688// Here is the conjugate of a^¯
//  0.680375  0.566198
// -0.211234   0.59688// Here is the matrix a^*
//  0.680375 -0.211234
//  0.566198   0.59688

5.4. 矩阵乘法

通过运算符*实现的
因为向量也是一种特殊的矩阵,所以矩阵的乘法对向量同样适用

Eigen提供了两个运算符:

  • ’ * ', 二元运算符,a * b
  • ’ *= ', 复合运算符,a = b(等价于a = ab)
#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{MatrixXd mat(2, 3);mat << 1, 2, 3,4, 5, 6;Vector3d u(-2, 1, 0);std::cout << "Here is mat:\n"<< mat << std::endl;std::cout << "Here is u:\n"<< u << std::endl;std::cout << "Here is mat*u:\n"<< mat * u << std::endl;
}// Here is mat:
// 1 2 3
// 4 5 6
// Here is u:
// -2
//  1
//  0
// Here is mat*u:
//  0
// -3

矩阵乘法,只有在第一个矩阵的列数和第二个矩阵的行数相同时才有意义
第一个矩阵行数为结果的行数,第二个矩阵列数为结果的列数

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;int main()
{Matrix2d mat;mat << 1, 2,3, 4;Vector2d u(-1, 1), v(2, 0);std::cout << "Here is mat*mat:\n"<< mat * mat << std::endl;std::cout << "Here is mat*u:\n"<< mat * u << std::endl;std::cout << "Here is u^T*mat:\n"<< u.transpose() * mat << std::endl;std::cout << "Here is u^T*v:\n"<< u.transpose() * v << std::endl;std::cout << "Here is u*v^T:\n"<< u * v.transpose() << std::endl;std::cout << "Let's multiply mat by itself" << std::endl;mat = mat * mat;std::cout << "Now mat is mat:\n"<< mat << std::endl;
}// Here is mat*mat:
//  7 10
// 15 22
// Here is mat*u:
// 1
// 1
// Here is u^T*mat:
// 2 2
// Here is u^T*v:
// -2
// Here is u*v^T:
// -2 -0
//  2  0
// Let's multiply mat by itself
// Now mat is mat:
//  7 10
// 15 22

5.5. 点乘和叉乘运算

可以通过dot()cross()方法实现点乘和叉乘运算


5.5.1. 点乘

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;int main()
{Vector2d v(1, 2);Vector2d w(3, 4);cout << "Dot product: " << v.dot(w) << endl;
}
// 11

点乘几何意义:投影

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;int main()
{Vector2d v(1, 1);Vector2d w(2, 2);cout << v.dot(w) / w.norm() << endl;
}
// 1.41421

5.5.2. 叉乘

点乘则适用于任何维度的向量,而叉乘运算只针对维度为3的向量

向量的叉乘,即求同时垂直两个向量的向量,即c垂直于a,同时c垂直于b

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;int main()
{Vector3d v(1, 2, 0);Vector3d w(1, 0, 0);cout << "Cross product:\n"<< v.cross(w) << endl;cout << "Area:\n"<< v.cross(w).norm() << endl;
}
// Cross product:
//  0
//  0
// -2
// Area:
// 2

叉乘的几何意义

|c|=|a×b|=|a| |b|sinα (α为a,b向量之间的夹角)
|c| = a,b向量构成的平行四边形的面积 (如下图所示的平行四边形)


5.5.3. 点到线的距离

组合使用,配合增广矩阵,可求点在直线法线上的投影,即2d点到2d线的距离

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;int main()
{Vector2d p_(-1, 2);Vector3d p(p_.x(), p_.y(), 1);Vector2d l_(0, 2);Vector3d lx(l_.x(), 0, 1);Vector3d ly(0, l_.y(), 1);cout << "The normal line is: " << lx.cross(ly).x() << ", " << lx.cross(ly).y() << endl;cout << "The distance to line is: " << p.dot(lx.cross(ly)) / l_.norm() << endl;
}
// The normal line is: -2, 0
// The distance to line is: 1

5.6. reduction运算

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;int main()
{Matrix2d mat;mat << 1, 2,3, 4;cout << "Here is mat.sum():       " << mat.sum() << endl;cout << "Here is mat.prod():      " << mat.prod() << endl;cout << "Here is mat.mean():      " << mat.mean() << endl;cout << "Here is mat.minCoeff():  " << mat.minCoeff() << endl;cout << "Here is mat.maxCoeff():  " << mat.maxCoeff() << endl;cout << "Here is mat.trace():     " << mat.trace() << endl;cout << "\n====\n"<< endl;// 获取minCoeff和maxCoeff元素的索引Matrix3f m = Matrix3f::Random();std::ptrdiff_t i, j;float minOfM = m.minCoeff(&i, &j);cout << "Here is the matrix m:\n"<< m << endl;cout << "Its minimum coefficient (" << minOfM<< ") is at position (" << i << "," << j << ")\n\n";RowVector4i v = RowVector4i::Random();int maxOfV = v.maxCoeff(&i);cout << "Here is the vector v: " << v << endl;cout << "Its maximum coefficient (" << maxOfV<< ") is at position " << i << endl;
}// Here is mat.sum():       10
// Here is mat.prod():      24
// Here is mat.mean():      2.5
// Here is mat.minCoeff():  1
// Here is mat.maxCoeff():  4
// Here is mat.trace():     5// ====// Here is the matrix m:
//  0.680375   0.59688 -0.329554
// -0.211234  0.823295  0.536459
//  0.566198 -0.604897 -0.444451
// Its minimum coefficient (-0.604897) is at position (2,1)// Here is the vector v:  115899597  -48539462  276748203 -290373134
// Its maximum coefficient (276748203) is at position 2

5.7. 分解函数


5.7.1. QR分解

#include <iostream>
#include <Eigen/Dense>using namespace std;
using namespace Eigen;// QR方法解线性方程组
int main()
{Matrix3f A;Vector3f b;A << 1, 2, 3, 4, 5, 6, 7, 8, 10;b << 3, 3, 4;cout << "Here is the matrix A:\n"<< A << endl;cout << "Here is the vector b:\n"<< b << endl;Vector3f x = A.colPivHouseholderQr().solve(b);cout << "The solution is:\n"<< x << endl;cout << "Verify:\n"<< A * x << endl;
}
// Here is the matrix A:
//  1  2  3
//  4  5  6
//  7  8 10
// Here is the vector b:
// 3
// 3
// 4
// The solution is:
//       -2
// 0.999997
//        1
// Verify:
// 3
// 3
// 4

5.7.2. 矩阵求逆

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;// 矩阵求逆
int main()
{Matrix2f A, b;A << 2, -1, -1, 3;b << 1, 2, 3, 1;cout << "Here is the matrix A:\n"<< A << endl;cout << "Here is the right hand side b:\n"<< b << endl;Matrix2f x = A.ldlt().solve(b);cout << "The solution is:\n"<< x << endl;
}// Here is the matrix A:
//  2 -1
// -1  3
// Here is the right hand side b:
// 1 2
// 3 1
// The solution is:
// 1.2 1.4
// 1.4 0.8

5.7.3. 计算数值法求解和真实值的残差

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;// 计算数值法求解和真实值的残差
int main()
{MatrixXd A = MatrixXd::Random(100, 100);MatrixXd b = MatrixXd::Random(100, 50);MatrixXd x = A.fullPivLu().solve(b);double relative_error = (A * x - b).norm() / b.norm(); // norm() is L2 normcout << "The relative error is:\n"<< relative_error << endl;
}

5.7.4. 最小二乘解

#include <iostream>
#include <Eigen/Dense>using namespace Eigen;
using namespace std;// 最小二乘解
int main()
{MatrixXf A = MatrixXf::Random(3, 2);cout << "Here is the matrix A:\n"<< A << endl;VectorXf b = VectorXf::Random(3);cout << "Here is the right hand side b:\n"<< b << endl;cout << "The least-squares solution is:\n"<< A.jacobiSvd(ComputeThinU | ComputeThinV).solve(b) << endl;
}// Here is the matrix A:
//  0.680375   0.59688
// -0.211234  0.823295
//  0.566198 -0.604897
// Here is the right hand side b:
// -0.329554
//  0.536459
// -0.444451
// The least-squares solution is:
// -0.669626
//  0.314253

更多应用请参考 Eigen的接口和应用


谢谢

Eigen 矩阵计算工具相关推荐

  1. [图解tensorflow源码] 入门准备工作附常用的矩阵计算工具[转]

    [图解tensorflow源码] 入门准备工作 附常用的矩阵计算工具[转] Link: https://www.cnblogs.com/yao62995/p/5773142.html tensorfl ...

  2. Eigen: C++开源矩阵计算工具——Eigen的简单用法

    Eigen非常方便矩阵操作,当然它的功能不止如此,由于本人只用到了它的矩阵相关操作,所以这里只给出了它的一些矩阵相关的简单用法,以方便快速入门.矩阵操作在算法研究过程中,非常重要,例如在图像处理中二维 ...

  3. Eigen C++开源矩阵计算工具——Eigen的简单用法

    Eigen非常方便矩阵操作,当然它的功能不止如此,由于本人只用到了它的矩阵相关操作,所以这里只给出了它的一些矩阵相关的简单用法,以方便快速入门.矩阵操作在算法研究过程中,非常重要,例如在图像处理中二维 ...

  4. C++开源矩阵计算工具——Eigen的简单用法(三)

    本节主要涉及Eigen的块操作以及QR分解,Eigen的QR分解非常绕人,搞了很久才搞明白是怎么回事,最后是一个使用Eigen的矩阵操作完成二维高斯拟合求取光点的代码例子,关于二维高斯拟合求取光点的详 ...

  5. C++开源矩阵计算工具——Eigen的简单用法(二)

    本文主要是Eigen中矩阵和向量的算术运算,在Eigen中的这些算术运算重载了C++的+,-,*,所以使用起来非常方便. 1.矩阵的运算 Eigen提供+.-.一元操作符"-".+ ...

  6. C++开源矩阵计算工具——Eigen的简单用法(一)

    Eigen非常方便矩阵操作,当然它的功能不止如此,由于本人只用到了它的矩阵相关操作,所以这里只给出了它的一些矩阵相关的简单用法,以方便快速入门.矩阵操作在算法研究过程中,非常重要,例如在图像处理中二维 ...

  7. C++开源矩阵计算工具——Eigen 在VS2005中的下载、配置与使用

    1.  下载Eigen Eigen的官网下载地址:http://eigen.tuxfamily.org/index.php?title=Main_Page#Download 下载后的文件名为:eige ...

  8. 2020-10-21Eigen: C++开源矩阵计算工具——Eigen的简单用法

    Eigen非常方便矩阵操作,当然它的功能不止如此,由于本人只用到了它的矩阵相关操作,所以这里只给出了它的一些矩阵相关的简单用法,以方便快速入门.矩阵操作在算法研究过程中,非常重要,例如在图像处理中二维 ...

  9. 矩阵基础2-常用的矩阵计算工具

    文章目录 一. R 1.1 R下载和安装 1.2 R简单使用 1.2.1 创建矩阵 1.2.2 矩阵计算 二. Excel 三 matlab 参考: 一. R 1.1 R下载和安装 R软件是免费的,在 ...

最新文章

  1. 使用LeNet对于旋转数字进行识别:合并数字集合
  2. 空投坐标怎么看6_装修时,怎么确认自己买的“瓷砖”是优等品?看“6点”很重要...
  3. 网狐棋牌(六) DataBaseEngine 和 网狐棋牌(七) CEventService
  4. Thread.interrupt()方法理解
  5. 第三章 第一部分 不定积分例题
  6. Mycat高级进阶---Mycat注解
  7. 还被python收智商税?做大数据的朋友告诉我月薪2w的方法
  8. mysql order by
  9. AD16修改规则加宽电源线与地线
  10. eclipse mac常用快捷键
  11. 国内首个SENT 信号解析软件 适配NXP KMA321, melexis MLX90372等SENT信号输出芯片 完美替代PicoScope 解析SENT
  12. java web 项目分模块,javaweb项目模块划分
  13. 微pe工具箱具体分区教程
  14. Maui Shell 来了,开启 Linux 桌面新时代!
  15. ThoughtWorks待遇
  16. R语言notes(1)——行列处理
  17. 77.组合 | 40.组合总和II | 39.组合总和 | 784.字母大小写全排列
  18. WGCLOUD——如何统计用户的日活(dau)、月活(mau)数据
  19. oracle连接失败的原因总结
  20. 三-五功能/半亮/25%亮/全亮/爆闪/SOS_专用应急灯手电筒IC方案

热门文章

  1. 使用nginx作为代理实现动静分离
  2. mac 下更新python
  3. 很久以前的C语言笔记
  4. ipvs学习笔记(二)
  5. Ant部署测试出错(关键字:Ant NoClassDefFoundError xml-apis/jar)
  6. JDK源码(17)-Compiler
  7. HGAT-用于半监督短文本分类的异构图注意力网络
  8. 二分法求正常水深c语言程序,水力學复习.doc
  9. 接口设计的几个注意事项
  10. HBase编程 API入门系列之create(管理端而言)(8)