第 35 章 · Eigen 进阶技巧

最后一章 Eigen讲实战中常用的技巧和模式与 STL 结合、自定义函数、常用工具函数、用Eigen::Map接管外部数据、内存对齐以及几个容易被忽略的坑。35.1 用 vector 管理一批矩阵实战中经常要存一批矩阵/向量用std::vector#includevector#includeEigen/Dense#includeiostreamintmain(){// 存一批 3D 点std::vectorEigen::Vector3dpoints;points.push_back(Eigen::Vector3d(1,0,0));points.push_back(Eigen::Vector3d(0,1,0));points.push_back(Eigen::Vector3d(0,0,1));// 遍历求和Eigen::Vector3d totalEigen::Vector3d::Zero();for(constautop:points){totalp;}std::cout总和 total.transpose()std::endl;// 1 1 1return0;}打印这一行不是可有可无的算出结果却不输出程序跑起来一句话都不说你会以为它坏了。每个完整示例都至少要有像样的输出这是本书写代码的习惯也是你验收自己代码的习惯。35.2 写操作矩阵的函数把常用操作封装成函数参数用const MatrixXd// 计算矩阵所有元素的和doublesumAll(constEigen::MatrixXdm){returnm.sum();}// 把矩阵归一化每列除以自己的范数Eigen::MatrixXdnormalizeColumns(constEigen::MatrixXdm){Eigen::MatrixXd resultm;for(Eigen::Index j0;jm.cols();j){result.col(j).normalize();// 每列归一化}returnresult;}35.3 模板函数处理任意 Eigen 类型想写一个能处理任意类型MatrixXd、Matrix3d、VectorXd的函数用模板templatetypenameDeriveddoublesumAll(constEigen::MatrixBaseDerivedm){returnm.sum();}MatrixBaseDerived是所有矩阵类型的基类模板能接受任何 Eigen 矩阵。35.4 常用工具函数Eigen 内置了很多方便的工具v.normalized();// 返回单位向量不修改原向量v.normalize();// 原地归一化m.cwiseProduct(n);// 逐元素乘等价 .array()*m.cwiseAbs();// 逐元素绝对值m.cwiseMax(other);// 逐元素取最大Eigen::Vector3d vEigen::Vector3d::UnitX();// X 轴单位向量Eigen::Vector3d uEigen::Vector3d::UnitY();// Y 轴Eigen::Vector3d wEigen::Vector3d::UnitZ();// Z 轴35.5 生成单位向量、对角矩阵// 单位轴向量Eigen::Vector3d xEigen::Vector3d::UnitX();// (1,0,0)// 把向量变成对角矩阵Eigen::Vector3dd(1,2,3);Eigen::Matrix3d Dd.asDiagonal();// 对角矩阵 diag(1,2,3)35.6 三维向量的几何操作Eigen::Vector3da(1,0,0),b(0,1,0);doubleanglestd::acos(a.dot(b));// 两向量夹角Eigen::Vector3d normala.cross(b);// 法向量垂直两者doubledist(a-b).norm();// 两点距离35.7 Eigen::Map把外部数据零拷贝接进 Eigen实战中数据常常不在 Eigen 矩阵里而在一段普通内存里C 接口给的double*、传感器缓冲区、std::vectordouble、别的语言传过来的数组。Eigen::Map让你不拷贝地把这块内存当成矩阵用——它只是一个视图改它就等于改原数据。#includeEigen/Dense#includecstdint#includeiostream#includevectorintmain(){// 1) 一段普通 C 数组比如驱动/其他语言给你的缓冲区当矩阵用doubleraw[6]{1,2,3,4,5,6};Eigen::MapEigen::MatrixXdm(raw,2,3);// 默认按「列优先」解读这块内存std::coutMap 出来的 2x3:\nm\n;// 2) 同一块内存换「行优先」解读Eigen::MapEigen::Matrixdouble,2,3,Eigen::RowMajorrm(raw);std::coutRowMajor 视角:\nrm\n;// 3) 改视图就是改原数组全程没有第二份数据m(0,0)100.0;std::cout改 m(0,0) 之后 raw[0] raw[0]\n\n;// 4) 与 std::vector 衔接数据存在 vector 里想让 Eigen 就地算std::vectordoublebuf{1,2,3,4};Eigen::MapEigen::Vector4dv(buf.data());v*2.0;std::coutvector 被就地改成:;for(doublex:buf)std::cout x;std::cout\n\n;// 5) 只取数组中间一段这里内存不满足对齐要求显式写 UnalignedEigen::MapEigen::Vector2d,Eigen::Unalignedpart(raw[1]);std::cout从 raw[1] 起的两个数 part.transpose()\n;// 6) C17 的 new 会遵守对齐要求vector 里直接放 Eigen 类型是安全的std::vectorEigen::Vector4dpts(3);boolaligned16(reinterpret_caststd::uintptr_t(pts.data())%160);std::coutvectorVector4d 的数据地址 16 字节对齐吗 (aligned16?是:否)\n;return0;}存成ch35_map.cpp按附录 C.3 编译运行实际输出-O0与-O2 -DNDEBUG完全一致Map 出来的 2x3: 1 3 5 2 4 6 RowMajor 视角: 1 2 3 4 5 6 改 m(0,0) 之后 raw[0] 100 vector 被就地改成: 2 4 6 8 从 raw[1] 起的两个数 2 3 vectorVector4d 的数据地址 16 字节对齐吗 是注意第 1、2 两段同一块内存读法不同矩阵就不同。raw里的 1 2 3 4 5 6按列优先排就是1 3 5 / 2 4 6按行优先排就是1 2 3 / 4 5 6。对接外部数据时先问清楚对方是哪种排法选错会得到转置了的结果——这正是第 17 章讲的行/列主序在真实工程里的样子。要点Eigen::MapMatrixXd(ptr, rows, cols)动态大小要显式给行列数固定大小可省MapVector4d(ptr)。Map不拥有内存也不拷贝原数据必须活得比 Map 久。想立刻固化一份自己的数据Eigen::MatrixXd copy m;这时才真拷贝。一次算完就丢的场合用 Map 能省掉整块拷贝大数组上差别很明显。35.8 关于「内存对齐」你需要知道的Eigen 用 SIMD 指令一次搬 2、4 甚至 8 个数这类指令要求数据起始地址是 16或 32字节倍数——这就是对齐。三个结论按重要程度排std::vectorEigen::Vector4d这类容器是安全的本书统一 C17 的原因之一。C17 起new会遵守类型的对齐要求上面第 6 段实测打印是。若你在 C11/14 下写这类代码Eigen 老文档会让你给容器加EIGEN_DEFINE_STL_VECTOR_SPECIALIZATION或aligned_allocator——那套已经过时别在 C17 项目里照抄。**把 Eigen 类型作为结构体/类的成员或者放进固定数组编译器一般会自动处理对齐。**真正会出问题的是手工拿一个指针去 Map。Map 一个不保证对齐的地址时显式写Eigen::Unaligned上面第 5 段。默认Map...声称已对齐写错就是未定义行为本书在 MinGW x86-64 上故意用未对齐地址试过MapVector4d实测没有崩、结果也对但这不能当保证——同一段代码换到 ARM 或不同指令集上可能直接段错误。既然多写一个模板参数只要几秒钟别赌。一句话自己 new 的 Eigen 对象不用管对齐拿外来指针做 Map 时才需要想这件事。35.9 几个容易被忽略的坑坑 1MatrixXd::Random()的值范围Random()生成 [-1, 1) 的均匀随机数。想要 0~1 或整数要自己变换Eigen::MatrixXd rEigen::MatrixXd::Random(3,3);// [-1, 1)Eigen::MatrixXd r01(rEigen::MatrixXd::Ones(3,3))/2;// 映射到 [0, 1)坑 2整数矩阵的除法MatrixXi的除法是整数除法截断和 double 矩阵不同Eigen::Matrix2i a;a5,6,7,8;a/2;// 结果是 2 3; 3 4截断不是 2.5 3 ...坑 3临时表达式的生命周期// ❌ 危险auto 捕获临时表达式见第 34 章—— 假设 A 是矩阵、b 和 c 是向量autoxA*bc;// ✅ 用明确类型这里换个名字否则和上面的 x 重名编译不过Eigen::VectorXd x2A*bc;实测这个区别auto x A * b c;之后打印x得51 111此时把A重新赋成全 0再打印x变成了1 1——因为x里存的是待计算的表达式它跟着A变。而Eigen::VectorXd版本一赋值就把结果定下来了。要保留auto就写auto x (A * b c).eval();实测这样推导出的类型就是Eigen::VectorXd。坑 4resize会丢失数据resize到更小尺寸会丢弃超出的数据更大的会保留原有部分但新增部分是垃圾值。35.10 完整示例批量数据处理#includeEigen/Dense#includeiostream#includevectorintmain(){// 模拟一批数据点每个点 3 维std::vectorEigen::Vector3ddata;for(inti0;i5;i){data.push_back(Eigen::Vector3d(i,i*2,i*3));}// 求所有点的平均Eigen::Vector3d meanEigen::Vector3d::Zero();for(constautop:data)meanp;mean/data.size();std::cout平均值 mean.transpose()std::endl;// 求离平均值最近的点doubleminDist1e9;Eigen::Vector3d closest;for(constautop:data){doubled(p-mean).norm();if(dminDist){minDistd;closestp;}}std::cout最近点 closest.transpose()std::endl;return0;}存成ch35.cpp按附录 C.3 编译运行实际输出平均值 2 4 6 最近点 2 4 65 个点是 (0,0,0)、(1,2,3)、(2,4,6)、(3,6,9)、(4,8,12)平均值就是 (2,4,6)它本身恰好在点集里所以离平均值最近的点也是它。把循环里的i * 2改成i * 2 1再跑一次预测一下输出会变成什么再用程序验证——这种改一处、猜结果、跑一遍是学 Eigen 最快的节奏。35.11 小结vectorEigen::...管理一批矩阵。函数参数用const MatrixXd或模板。工具函数normalized、asDiagonal、UnitX/Y/Z、cwiseProduct。Eigen::Map把外部内存零拷贝当矩阵用先确认对方是行主序还是列主序。对齐只在拿外来指针做 Map时才需要考虑不保证对齐就显式写Eigen::Unaligned。注意 Random 范围、整数除法、auto 陷阱、resize 丢数据。 第三部分完成你已经完整掌握了 Eigen 库下一部分我们用两个综合项目把知识串起来并规划你的学习路线。练习题用std::vector存 5 个三维向量求它们的平均值。写一个模板函数求任意 Eigen 矩阵的元素的平均值。用asDiagonal把一个向量变成对角矩阵。求两个三维向量的夹角用acos(dot)。有一行 12 个double的数组用Eigen::Map把它当成 3×4 矩阵打印再换成行主序的 Map 打印一次比较两者差别。从数组中间某个元素开始Map应该加哪个模板参数为什么本书不推荐反正实测没崩的写法总结你学到的所有 Eigen 的坑。