1. 程式人生 > >Eigen之矩陣、向量、標量的操作運算

Eigen之矩陣、向量、標量的操作運算

介紹

Eigen是通過中過載C++操作運算子如+、-、*或通過dot()、cross()等來實現矩陣/向量的操作運算

加法和減法

  • binary operator + as in a+b
  • binary operator - as in a-b
  • unary operator - as in -a
  • compound operator += as in a+=b
  • compound operator -= as in 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; }

Output:

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

標量的乘除

  • binary operator * as in matrix*scalar
  • binary operator * as in scalar*matrix
  • binary operator / as in matrix/scalar
  • compound operator = as in matrix=scalar
  • compound operator /= as in 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;
}

Output:

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

A note about expression templates

這是我們在這個頁面上解釋的一個高階主題,但是現在只提到它是很有用的。在Eigen,算術運算子,比如運算子+不執行任何計算,他們只是返回一個描述計算的“表示式物件”。實際的計算髮生在稍後,當整個表示式被求值時,通常在操作符=。雖然這聽起來可能很重,但是任何現代的優化編譯器都能夠優化抽象,結果是完美的優化程式碼。例如,當你這樣做的時候:

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

Eigen將其編譯為一個for迴圈,這樣陣列就只遍歷一次。簡化(例如忽略SIMD優化),這個迴圈是這樣的:

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

因此,你不應該害怕使用Eigen的相對大的算術表示式:這隻給了Eigen更多的機會來優化。

矩陣的轉置、共軛及伴隨

例:

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(),則結果值會賦值給b,但是如果進行的運算為a = a.transpose(),則Eigen會在轉置完成之前將值賦值給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
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

矩陣與矩陣、矩陣與向量間的乘法

  • binary operator * as in a*b
  • compound operator = as in a=b (this multiplies on the right: a*=b
    is equivalent to a = 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

提示:對於矩陣的乘法不會出現上述的混疊問題,因為Eigen將乘法作為一個特例

tmp = m*m;
m = tmp;

點乘和叉乘

Eigen中通過dot()和cross()來實現點乘和叉乘
例:

#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

其他基礎操作

Eigen還提供了一些函式實現一些其他基本操作,如sum()進行求和,prod()實現矩陣內值的積,maxCoeff()求矩陣內最大值等。
例:

#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

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