资讯详情

第 29 章 · 求解线性方程组

📅 2026/10/10 9:22:42 | 华诺云谱 👁 阅读
第 29 章 · 求解线性方程组
Eigen 最重要的应用之一解 Ax b。本章讲标准解法、各种求解器怎么选以及为什么不要用逆矩阵。29.1 标准解法solve()解 Ax b 的标准写法Eigen::Matrix2d A;A2,1,1,3;Eigen::Vector2db(5,10);Eigen::Vector2d xA.colPivHouseholderQr().solve(b);std::coutx x.transpose()std::endl;// 1 3A.分解方法().solve(b)是标准模式先对 A 做某种分解再用分解结果求解b29.2 验证解的正确性// 验证 A*x 是否等于 bdoubleerror(A*x-b).norm();std::cout残差 errorstd::endl;// 应接近 0残差residual越接近 0解越准确。29.3 几种求解器怎么选不同矩阵性质用不同求解器求解器适用场景速度colPivHouseholderQr()通用首选稳健适合大多数中llt()对称正定矩阵最快ldlt()对称半正定矩阵快partialPivLu()一般方阵快SVDbdcSvd…()、jacobiSvd…()最小二乘方程多、未知少慢但强⚠️SVD 这一行的选项必须写在尖括号里写法见下面 29.3.1别直接写A.bdcSvd().solve(b)。// 对称正定矩阵用 llt() 最快Eigen::Matrix2d S;S2,1,1,3;Eigen::Vector2d x2S.llt().solve(b);怎么记不确定用哪个 →colPivHouseholderQr()最稳确定对称正定 →llt()最快方程比未知数多超定最小二乘→colPivHouseholderQr()够用且快要 SVD 就bdcSvd…()见下29.3.1 SVD 求解器选项写在尖括号里Eigen 5.0 把 SVD 的计算选项算不算 U、V改成了模板参数。官方CHANGELOG.md在[5.0.0] - 2025-09-30一节写得很明白Runtime SVD options for computing thin/full U/V have been deprecated: use compile-time options instead.也就是运行时传选项这条路已被正式废弃。老写法A.bdcSvd().solve(b)能编译通过但一运行就断言崩溃Assertion failed: computeU() computeV() SVDBase::solve(): Both unitaries U and V are required to be computed (thin unitaries suffice)., file .../Eigen/src/SVD/SVDBase.h, line 322正确写法|是按位或把两个选项合起来#includeEigen/Dense#includeiostreamintmain(){// 3 个方程、2 个未知数典型的超定问题Eigen::MatrixXdA(3,2);A1,1,2,3,4,5;Eigen::VectorXdb(3);b1,2,3;Eigen::VectorXd xA.bdcSvdEigen::ComputeFullU|Eigen::ComputeFullV().solve(b);std::coutx x.transpose()std::endl;std::cout残差 (A*x-b).norm()std::endl;return0;}编译运行编译命令见第 22 章x 0.166667 0.5 残差 0.408248残差不为 0 是正常的——超定问题本来就没有精确解最小二乘给出的是误差平方和最小的那个解。两种 SVD 求解器怎么选解法是矩阵的成员函数A.jacobiSvd选项()、A.bdcSvd选项()注意函数名是小写Svd类名才是大写的JacobiSVD。小矩阵用jacobiSvd大矩阵用bdcSvd更快两者的选项都写在尖括号里实测解一样。Eigen::BDCSVDMatrix3d svd(A, 选项)这种把选项写在圆括号里的老构造函数在 Eigen 5.0 已被标记为deprecated能跑但每次编译都出警告。要写类名形式就写Eigen::BDCSVDMatrix3d, Eigen::ComputeThinU | Eigen::ComputeThinV svd(A);只想解方程、不关心 U/V 具体是什么 → 直接用colPivHouseholderQr().solve(b)它就能给最小二乘解上面那段代码换成它输出完全一样x 0.166667 0.5还更快。第 36 章的拟合项目就是这么写的。29.4 为什么不用逆矩阵数学上 x A⁻¹·b 没错但工程上强烈不推荐// ❌ 不推荐慢且数值不稳定Eigen::Vector2d x_badA.inverse()*b;// ✅ 推荐更快更稳定Eigen::Vector2d x_goodA.colPivHouseholderQr().solve(b);原因慢求逆矩阵本身要更多计算再乘 b 是浪费数值不稳定矩阵接近奇异时逆矩阵误差被放大记住这条铁律解方程用.solve()永远不要A.inverse() * b。29.5 检查矩阵是否可解求解前可以检查std::abs求浮点绝对值需要#include cmath#includecmath// std::abs 对 double 版本在这里// 检查行列式是否为 0if(std::abs(A.determinant())1e-12){std::cout矩阵不可逆可能无解或无穷多解std::endl;}或者用分解的info()检查是否成功第 30 章。29.6 实际例子解一个 3 元方程组3x y - z 1 x 4y z 2 2x - y 5z 3#includeEigen/Dense#includeiostreamintmain(){Eigen::Matrix3d A;A3,1,-1,1,4,1,2,-1,5;Eigen::Vector3db(1,2,3);Eigen::Vector3d xA.colPivHouseholderQr().solve(b);std::cout解 x x.transpose()std::endl;std::cout验证 A*x (A*x).transpose()std::endl;std::cout残差 (A*x-b).norm()std::endl;return0;}存成ch29.cpp按附录 C.3 编译运行实际输出解 x 0.405797 0.275362 0.492754 验证 A*x 1 2 3 残差 6.28037e-16三行连起来读就是 Eigen 解方程的标准验收流程先拿解再乘回原式看是否等于 b最后看残差。验证 A*x打印出1 2 3正是原来的 b残差 6.28037e-16就是 0浮点残渣见第 21 章末尾说明说明解是对的。29.7 小结解方程A.分解().solve(b)。通用选colPivHouseholderQr()对称正定选llt()最小二乘选colPivHouseholderQr()或bdcSvdComputeFullU | ComputeFullV()SVD 的选项必须写在尖括号里见 29.3.1。验证(A*x - b).norm()残差接近 0。铁律用.solve()别用.inverse() * b。下一章深入矩阵分解——solve 背后的原理。练习题解方程组3x y 9; x 2y 8验证结果。解一个 3 元方程组打印残差。用llt()解一个对称正定方程组和colPivHouseholderQr()对比结果。分别用.solve()和.inverse()*b解同一方程说明为什么前者更好。构造一个行列式为 0 的矩阵尝试 solve观察结果。
📝

华诺云谱内容团队

资深建站顾问 · 行业研究员

10年+企业数字化服务经验,专注智能建站、SEO优化与品牌营销,持续输出建站技巧、行业洞察与营销干货,已帮助5000+企业实现数字化增长。

你可能需要的服务

订阅华诺云谱资讯周报

每周一封,精选建站技巧、SEO与营销干货,直达邮箱。已有 8,000+ 企业主订阅,助你少走弯路。

↑