Eigen学习笔记(14)-原位矩阵分解-卡核

原文：Eigen官网-Inplace matrix decompositions

从Eigen3.3开始，LU、Cholesky和QR分解可以就地操作，即直接在给定的输入矩阵内操作。当处理大型矩阵或可用内存非常有限（嵌入式系统）时，此功能特别有用。

为此，必须使用Ref<>矩阵类型实例化相应的分解类，并且必须使用输入矩阵作为参数构造分解对象。作为一个例子，让我们考虑一个局部旋转的就地LU分解。

#include <iostream>
#include <Eigen/Dense>
using namespace std;
using namespace Eigen;
int main()
{
  MatrixXd A(2,2);
  A << 2, -1, 1, 3;
  cout << "Here is the input matrix A before decomposition:\\n" << A << endl;
  //声明inplace LU 对象 lu,并输出此时矩阵A中的内容
  //这里，对象lu 计算L和U因子并将其存储在由矩阵A所在的内存中。因此，A的系数在因子分解期间被破坏，并被L和U因子替换，因为可以验证：
  PartialPivLU<Ref<MatrixXd> > lu(A);
  cout << "Here is the input matrix A after decomposition:\\n" << A << endl;
  
  //输出lu对象的L和U因子，可以发现，结果和矩阵A的结果一样
  cout << "Here is the matrix storing the L and U factors:\\n" << lu.matrixLU() << endl;
  
  //然后，可以像往常一样使用lu对象，例如解决Ax=b问题：
  //由于原始矩阵A的内容已经丢失，我们不得不声明一个新的矩阵A0来验证结果。
  MatrixXd A0(2,2); A0 << 2, -1, 1, 3;
  VectorXd b(2);    b << 1, 2;
  VectorXd x = lu.solve(b);
  cout << "Residual: " << (A0 * x - b).norm() << endl;
  // 输出结果：Residual: 0
  
  //由于内存在A和lu之间共享，修改矩阵A将使lu无效。通过修改A的内容并再次尝试解决初始问题，可以很容易地验证这一点：
  A << 3, 4, -2, 1;
  x = lu.solve(b);
  cout << "Residual: " << (A0 * x - b).norm() << endl;
  // 输出结果：Residual: 15.8114
  
  //如果要用修改后的A更新因子分解，必须像往常一样调用compute方法：
  A0 = A; // save A
  lu.compute(A);
  x = lu.solve(b);
  cout << "Residual: " << (A0 * x - b).norm() << endl;
  // 输出结果：Residual: 0
 
  //请注意，调用compute不会更改lu对象引用的内存。因此，如果使用与A不同的另一个矩阵A1调用计算方法，则不会修改A1的内容。这仍然是用于存储矩阵A1的L和U因子的A的内容。这一点很容易验证如下：
  MatrixXd A1(2,2);
  A1 << 5,-2,3,4;
  lu.compute(A1);
  cout << "Here is the input matrix A1 after decomposition:\\n" << A1 << endl;
  // 输出结果： Here is the input matrix A1 after decomposition:
 // 5 -2
 // 3  4
 //The matrix A1 is unchanged, and one can thus solve A1*x=b, and directly check the residual without any copy of A1:
  x = lu.solve(b);
  cout << "Residual: " << (A1 * x - b).norm() << endl;
  // 输出结果： Residual: 2.48253e-16

结果如下：

Here is the input matrix A before decomposition:
 2 -1
 1  3
 
 Here is the input matrix A after decomposition:
  2  -1
0.5 3.5

Here is the matrix storing the L and U factors:
  2  -1
0.5 3.5

如下是支持inplace机制的矩阵分解：