轉載於:https://www.cnblogs.com/houkai/p/6348044.html

矩陣和向量的運算

提供一些概述和細節：關於矩陣、向量以及標量的運算。

介紹

Eigen提供了matrix/vector的運算操作，既包括過載了c++的算術運算子+/-/*，也引入了一些特殊的運算比如點乘dot、叉乘cross等。

對於Matrix類（matrix和vectors）這些操作只支援線性代數運算，比如：matrix1*matrix2表示矩陣的乘機，vetor+scalar是不允許的。如果你想執行非線性代數操作，請看下一篇（暫時放下）。

加減

左右兩側變數具有相同的尺寸（行和列），並且元素型別相同（Eigen不自動轉化型別）操作包括：

• 二元運算 + 如a+b
• 二元運算 - 如a-b
• 一元運算 - 如-a
• 複合運算 += 如a+=b
• 複合運算 -= 如a-=b

#include <iostream>
#include <Eigen/Dense>
using namespace Eigen;
int main()
{
  Matrix2d a;
  a << 1, 2,
       3, 4;
  MatrixXd b(2,2);
  b << 2, 3,
       1, 4;
  std::cout << "a + b =\n" << a + b << std::endl;
  std::cout << "a - b =\n" << a - b << std::endl;
  std::cout << "Doing a += b;" << std::endl;
  a += b;
  std::cout << "Now a =\n" << a << std::endl;
  Vector3d v(1,2,3);
  Vector3d w(1,0,0);
  std::cout << "-v + w - v =\n" << -v + w - v << std::endl;
}

輸出：

a + b =
3 5
4 8
a - b =
-1 -1
 2  0
Doing a += b;
Now a =
3 5
4 8
-v + w - v =
-1
-4
-6

標量乘法和除法

乘/除標量是非常簡單的，如下：

• 二元運算 * 如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;
  Vector3d v(1,2,3);
  std::cout << "a * 2.5 =\n" << a * 2.5 << std::endl;
  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

表示式模板

這裡簡單介紹，在高階主題中會詳細解釋。在Eigen中，線性運算比如+不會對變數自身做任何操作，會返回一個“表示式物件”來描述被執行的計算。當整個表示式被評估完（一般是遇到=號），實際的操作才執行。

這樣做主要是為了優化，比如

VectorXf a(50), b(50), c(50), d(50);
...
a = 3*b + 4*c + 5*d;

Eigen會編譯這段程式碼最終遍歷一次即可運算完成。

for(int i = 0; i < 50; ++i)
  a[i] = 3*b[i] + 4*c[i] + 5*d[i];

因此，我們不必要擔心大的線性表示式的運算效率。

轉置和共軛

 表示transpose轉置

 表示conjugate共軛

 表示adjoint(共軛轉置) 伴隨矩陣

MatrixXcf a = MatrixXcf::Random(2,2);
cout << "Here is the matrix a\n" << a << endl;
cout << "Here is the matrix a^T\n" << a.transpose() << endl;
cout << "Here is the conjugate of a\n" << a.conjugate() << endl;
cout << "Here is the matrix a^*\n" << a.adjoint() << endl;

輸出

Here is the matrix a
 (-0.211,0.68) (-0.605,0.823)
 (0.597,0.566)  (0.536,-0.33)
Here is the matrix a^T
 (-0.211,0.68)  (0.597,0.566)
(-0.605,0.823)  (0.536,-0.33)
Here is the conjugate of a
 (-0.211,-0.68) (-0.605,-0.823)
 (0.597,-0.566)    (0.536,0.33)
Here is the matrix a^*
 (-0.211,-0.68)  (0.597,-0.566)
(-0.605,-0.823)    (0.536,0.33)

對於實數矩陣，conjugate不執行任何操作，adjoint等價於transpose。

transpose和adjoint會簡單的返回一個代理物件並不對本省做轉置。如果執行 b=a.transpose() ，a不變，轉置結果被賦值給b。如果執行 a=a.transpose() Eigen在轉置結束之前結果會開始寫入a，所以a的最終結果不一定等於a的轉置。

Matrix2i a; a << 1, 2, 3, 4;
cout << "Here is the matrix a:\n" << a << endl;
a = a.transpose(); // !!! do NOT do this !!!
cout << "and the result of the aliasing effect:\n" << a << endl;

Here is the matrix a:
1 2
3 4
and the result of the aliasing effect:
1 2
2 4

這被稱為“別名問題”。在debug模式，當assertions開啟的情況加，這種常見陷阱可以被自動檢測到。

對 a=a.transpose() 這種操作，可以執行in-palce轉置。類似還有adjointInPlace。

MatrixXf a(2,3); a << 1, 2, 3, 4, 5, 6;
cout << "Here is the initial matrix a:\n" << a << endl;
a.transposeInPlace();
cout << "and after being transposed:\n" << a << endl;

Here is the initial matrix a:
1 2 3
4 5 6
and after being transposed:
1 4
2 5
3 6

矩陣-矩陣的乘法和矩陣-向量的乘法

向量也是一種矩陣，實質都是矩陣-矩陣的乘法。

• 二元運算 *如a*b
• 複合運算 *=如a*=b

#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

m=m*m並不會導致別名問題，Eigen在這裡做了特殊處理，引入了臨時變數。實質將編譯為：

tmp = m*m
m = tmp

如果你確定矩陣乘法是安全的（並沒有別名問題），你可以使用noalias()函式來避免臨時變數 c.noalias() += a*b 。

點運算和叉運算

dot()執行點積，cross()執行叉積，點運算得到1*1的矩陣。當然，點運算也可以用u.adjoint()*v來代替。

#include <iostream>
#include <Eigen/Dense>
using namespace Eigen;
using namespace std;
int main()
{
  Vector3d v(1,2,3);
  Vector3d w(0,1,2);
  cout << "Dot product: " << v.dot(w) << endl;
  double dp = v.adjoint()*w; // automatic conversion of the inner product to a scalar
  cout << "Dot product via a matrix product: " << dp << endl;
  cout << "Cross product:\n" << v.cross(w) << endl;
}

輸出

Dot product: 8
Dot product via a matrix product: 8
Cross product:
 1
-2
 1

注意：點積只對三維vector有效。對於複數，Eigen的點積是第一個變數共軛和第二個變數的線性積。

基礎的歸約操作

Eigen提供了而一些歸約函式：sum()、prod()、maxCoeff()和minCoeff()，他們對所有元素進行操作。

#include <iostream>
#include <Eigen/Dense>
using namespace std;
int main()
{
  Eigen::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;
}

輸出

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

trace表示矩陣的跡，對角元素的和等價於 a.diagonal().sum() 。

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 the matrix m:
  0.68  0.597  -0.33
-0.211  0.823  0.536
 0.566 -0.605 -0.444
Its minimum coefficient (-0.605) is at position (2,1)

Here is the vector v:  1  0  3 -3
Its maximum coefficient (3) is at position 2

操作的有效性

Eigen會檢測執行操作的有效性，在編譯階段Eigen會檢測它們，錯誤資訊是繁冗的，但錯誤資訊會大寫字母突出，比如:

Matrix3f m;
Vector4f v;
v = m*v;      // Compile-time error: YOU_MIXED_MATRICES_OF_DIFFERENT_SIZES

當然動態尺寸的錯誤要在執行時發現，如果在debug模式，assertions會觸發後，程式將崩潰。

MatrixXf m(3,3);
VectorXf v(4);
v = m * v; // Run-time assertion failure here: "invalid matrix product"