4. Block operations
1) Using block operations
Block of size(p, q), starting at (i, j)
dynamic-size block expression: matrix.block(i, j, p, q);
fixed-size block expression: matrix.block<p, q>(i, j);
int main()
{
Eigen::MatrixXf m(,);
m << , , , ,
, , , ,
,,,,
,,,;
cout << "Block in the middle" << endl;
cout << m.block<,>(,) << endl << endl;
for (int i = ; i <= ; ++i)
{
cout << "Block of size " << i << "x" << i << endl;
cout << m.block(,,i,i) << endl << endl;
}
}
Output:
Block in the middle Block of size 1x1 Block of size 2x2 Block of size 3x3
.block() 的返回值可以是左值, 也可以是右值.
2) Columns and rows
matrix.row(i), matrix.col(j)
int main()
{
Eigen::MatrixXf m(,);
m << ,,,
,,,
,,;
cout << "Here is the matrix m:" << endl << m << endl;
cout << "2nd Row: " << m.row() << endl;
m.col() += * m.col();
cout << "After adding 3 times the first column into the third column, the matrix m is:\n";
cout << m << endl;
}
3) Corner-related operations
太复杂, 没兴趣.
4) Block operations for vectors
同上, 用时再查.
5. Advanced initialization
1) The comma initializer
RowVectorXd vec1();
vec1 << , , ; RowVectorXd vec2();
vec2 << , , , ; RowVectorXd joined();
joined << vec1, vec2; MatrixXf matA(, );
matA << , , , ; MatrixXf matB(, );
matB << matA, matA/, matA/, matA; Matrix3f m;
m.row() << , , ;
m.block(,,,) << , , , ;
m.col().tail() << , ;
2) Special matrices and arrays
Zero(), Constant(value), Constant(rows, cols, value), Random(), Identity(), LinSpaced(size, low, high),
Identity() only for matrix;
LinSpaced only for vectors and one-dimensional arrays; high也包含在数组中. 默认都是列向量.
std::cout << "A fixed-size array:\n";
Array33f a1 = Array33f::Zero(); std::cout << "A one-dimensional dynamic-size array:\n";
ArrayXf a2 = ArrayXf::Zero(); std::cout << "A two-dimensional dynamic-size array:\n";ArrayXXf a3 = ArrayXXf::Zero(, );
ArrayXXf table(, );
table.col() = ArrayXf::LinSpaced(, , );
table.col() = M_PI / * table.col();
table.col() = table.col().sin();
table.col() = table.col().cos(); Output:
Degrees Radians Sine Cosine 0.175 0.174 0.985
0.349 0.342 0.94
0.524 0.5 0.866
0.698 0.643 0.766
0.873 0.766 0.643
1.05 0.866 0.5
1.22 0.94 0.342
1.4 0.985 0.174
1.57 -4.37e-08
类似的, 有 setZero(), MatrixBase::setIdentity(), DenseBase::setLinSpaced()
3) Usage as temporary objects
int main()
{
MatrixXd m = MatrixXd::Random(,);
m = (m + MatrixXd::Constant(,,1.2)) * ;
cout << "m =" << endl << m << endl;
VectorXd v();
v << , , ;
cout << "m * v =" << endl << m * v << endl;
}
comma-initializer
MatrixXf mat = MatrixXf::Random(, );
std::cout << mat << std::endl << std::endl;
mat = (MatrixXf(,) << , , , ).finished() * mat;
注意到, .finished() 是必须的.
6. Reductions, visitors and broadcasting
1) Reducations
sum(), prod(), mean(), minCoeff(), maxCoeff(), trace()
2) Norm computations
应该用不到.
3) Boolean reductions
all(), any(), count()
int main()
{
ArrayXXf a(,); a << ,,
,;
cout << "(a > 0).all() = " << (a > ).all() << endl;
cout << "(a > 0).any() = " << (a > ).any() << endl;
cout << "(a > 0).count() = " << (a > ).count() << endl;
cout << endl;
cout << "(a > 2).all() = " << (a > ).all() << endl;
cout << "(a > 2).any() = " << (a > ).any() << endl;
cout << "(a > 2).count() = " << (a > ).count() << endl;
} Output:
(a > ).all() =
(a > ).any() =
(a > ).count() = (a > ).all() =
(a > ).any() =
(a > ).count() =
4) Visitors
maxCoeff(&x, &y), minCoeff(&x, &y)
x, y 即为行索引和列索引, MatrixXf::Index, 注意, 不是 int
//get location of maximum
MatrixXf::Index maxRow, maxCol;
float max = m.maxCoeff(&maxRow, &maxCol); //get location of minimum
MatrixXf::Index minRow, minCol;
float min = m.minCoeff(&minRow, &minCol);
5) Partial reductions
colwise(), rowwise()
int main() {
Eigen::MatrixXf mat(, );
mat << , , , ,
, , , ;
cout << "Column's maximum: " << endl
<< mat.colwise().maxCoeff() << endl;
cout << "Row's maximum: " << endl
<< mat.rowwise().maxCoeff() << endl;
} Output:
Column's maximum: Row's maximum:
combining partial reductions with other operations
int main() {
MatrixXf mat(, );
mat << , , , ,
, , , ;
MatrixXf::Index maxIndex;
float maxNorm = mat.colwise().sum().maxCoeff(&maxIndex); cout << "Maximum sum at position: " << maxIndex << endl;
cout << "The corresponding vector is: " << endl;
cout << mat.col(maxIndex) << endl;
cout << "And its sum is: " << maxNorm << endl;
} Output:
Maximum sum at position
The corresponding vector is: And its sum is is:
6) Broadcasting
int main() {
MatrixXf mat(, );
VectorXf v(); mat << , , , ,
, , , ;
v << , ;
w << , , , ; mat.colwise() += v; cout << "Broadcasting result: " << endl;
cout << mat << endl; mat.rowwise() += w.transpose();
cout << "Broadcasting result: " << endl;
cout << mat << endl;
} Output:
Broadcasting result:
1 2 6 9
4 2 8 3
Broadcasting result:
7) Combining broadcasting with other operations
略.
7. Interfacing with raw buffers: the Map class
1) Map types and declaring Map variables
RowsAtCompileTime 和 ColsAtCompileTime 是生成矩阵的行数和列数
Map<MatrixXf> mf(pf, rows, columns);
pf 指向一个定义 coefficient array 的 region 的起始, 这里注意 pf 是一个指针, 指向的是一个 coefficient, rows 和 columns 是该 mf 的行和列, 默认全部.
这里有一个问题, Map 模板参数指定的 RowsAtCompileTime 和 ColsAtCompileTime 与括号中的 rows 和 columns 有什么关系?
a. Matrix
Dynamic
Map<Matrix<int, Dynamic, Dynamic> > 即 Map<MatrixXi>, 此时生成 matrix 的行和列由 rows 和 columns 指定
这里还要注意一下, pf 的属性, 如果 pf 是指向一个内存区域, 内存区域的数字怎么获取, 要看这个区域存储的是什么.
存储 array, 假设 array 长度为 15, 要转换成3行4列矩阵, 那么取前12个数
存储 matrix, 注意 matrix 默认是按列存储的, 所以会取按列排布的前12个数, 也就是前4列
fixed-size
Map<Matrix<int, 3, 4> > , 此时生成 Matrix 的行和列已定, rows 和 columns 要与 3, 4 保持一致, 或者省略否则报错.
b. Vector
Dynamic
Map<Matrix<int, 1, Dynamic> > or Map<Matrix<int, Dynamic, 1> >
Vector 有一个特别地方, rows 和 columns 可以只写一个, 比如 Map<VectorXi> m2map(p, m2.size()), 也可以全写, 但是维度要与 Map 中的一致.
fixed-size
注意维度一致.
Map<const Vector4i> mi(pi);
声明一个 fixed-size read-only vector.
Map 还有几个可选的模板参数
Map<typename MatrixType, int MapOptions, typename StrideType>
MapOptions: Aligned or Unaligned, 默认 Unaligned
StrideType: 决定 array 的排布, 使用 Stride class
int array[];
for(int i = ; i < ; ++i) array[i] = i;
cout << "Column-major:\n" << Map<Matrix<int,,> >(array) << endl;
cout << "Row-major:\n" << Map<Matrix<int,,,RowMajor> >(array) << endl;
cout << "Row-major using stride:\n" <<
Map<Matrix<int,,>, Unaligned, Stride<,> >(array) << endl; Output: Column-major: Row-major: Row-major using stride:
2) Using Map variables
注意: a Map type is not identical to its Dense equivalent.
typedef Matrix<float,,Dynamic> MatrixType;
typedef Map<MatrixType> MapType;
typedef Map<const MatrixType> MapTypeConst; // a read-only map
const int n_dims = ; MatrixType m1(n_dims), m2(n_dims);
m1.setRandom();
m2.setRandom();
float *p = &m2(); // get the address storing the data for m2
MapType m2map(p,m2.size()); // m2map shares data with m2
MapTypeConst m2mapconst(p,m2.size()); // a read-only accessor for m2
cout << "m1: " << m1 << endl;
cout << "m2: " << m2 << endl;
cout << "Squared euclidean distance: " << (m1-m2).squaredNorm() << endl;
cout << "Squared euclidean distance, using map: " <<
(m1-m2map).squaredNorm() << endl;
m2map() = ; // this will change m2, since they share the same array
cout << "Updated m2: " << m2 << endl;
cout << "m2 coefficient 2, constant accessor: " << m2mapconst() << endl;
/* m2mapconst(2) = 5; */ // this yields a compile-time error Output: m1: 0.68 -0.211 0.566 0.597 0.823
m2: -0.605 -0.33 0.536 -0.444 0.108
Squared euclidean distance: 3.26
Squared euclidean distance, using map: 3.26
Updated m2: -0.605 -0.33 0.536 0.108
m2 coefficient , constant accessor: 0.536
这里注意一下, m2map, m2mapconst 与 m2 共享内存空间, 所以对 m2map 的修改, 会同步到 m2 和 m2mapconst.
但是, 如果如果 *p 指向的是 array 中的地址, 那么 map 转化出的 matrix 不会与 array 共享内存空间.
3) Changing the mapped array
用 new 可以对声明后的 Map object 进行修改
int data[] = {,,,,,,,,};
Map<RowVectorXi> v(data,);
cout << "The mapped vector v is: " << v << "\n";
new (&v) Map<RowVectorXi>(data+,);
cout << "Now v is: " << v << "\n"; Output:
The mapped vector v is:
Now v is:
所以也可以声明未知 array 的 map object
Map<Matrix3f> A(NULL); // don't try to use this matrix yet!
VectorXf b(n_matrices);
for (int i = ; i < n_matrices; i++)
{
new (&A) Map<Matrix3f>(get_matrix_pointer(i));
b(i) = A.trace();
}
8. Reshape and Slicing
1) Reshape
MatrixXf M1(,); // Column-major storage
M1 << , , ,
, , ,
, , ;
Map<RowVectorXf> v1(M1.data(), M1.size());
cout << "v1:" << endl << v1 << endl;
Matrix<float,Dynamic,Dynamic,RowMajor> M2(M1);
Map<RowVectorXf> v2(M2.data(), M2.size());
cout << "v2:" << endl << v2 << endl; Output:
v1: v2:
MatrixXf M1(,); // Column-major storage
M1 << , , , , , ,
, , , , , ;
Map<MatrixXf> M2(M1.data(), ,);
cout << "M2:" << endl << M2 << endl; Output:
M2:
注意这里, 存储方式的不同, 会导致输出不同.
2) Slicing
RowVectorXf v = RowVectorXf::LinSpaced(,,);
cout << "Input:" << endl << v << endl;
Map<RowVectorXf,,InnerStride<> > v2(v.data(), v.size()/);
cout << "Even:" << v2 << endl; Output:
Input: Even:
MatrixXf M1 = MatrixXf::Random(,);
cout << "Column major input:" << endl << M1 << "\n";
Map<MatrixXf,,OuterStride<> > M2(M1.data(), M1.rows(), (M1.cols()+)/, OuterStride<>(M1.outerStride()*));
cout << "1 column over 3:" << endl << M2 << "\n";
typedef Matrix<float,Dynamic,Dynamic,RowMajor> RowMajorMatrixXf;
RowMajorMatrixXf M3(M1);
cout << "Row major input:" << endl << M3 << "\n";
Map<RowMajorMatrixXf,,Stride<Dynamic,> > M4(M3.data(), M3.rows(), (M3.cols()+)/,
Stride<Dynamic,>(M3.outerStride(),));
cout << "1 column over 3:" << endl << M4 << "\n"; Output:
Column major input:
0.68 0.597 -0.33 0.108 -0.27 0.832 -0.717 -0.514
-0.211 0.823 0.536 -0.0452 0.0268 0.271 0.214 -0.726
0.566 -0.605 -0.444 0.258 0.904 0.435 -0.967 0.608
column over :
0.68 0.108 -0.717
-0.211 -0.0452 0.214
0.566 0.258 -0.967
Row major input:
0.68 0.597 -0.33 0.108 -0.27 0.832 -0.717 -0.514
-0.211 0.823 0.536 -0.0452 0.0268 0.271 0.214 -0.726
0.566 -0.605 -0.444 0.258 0.904 0.435 -0.967 0.608
column over :
0.68 0.108 -0.717
-0.211 -0.0452 0.214
0.566 0.258 -0.967
注意 M1.data(), 返回一个指针.
9. Aliasing
m = m.transpose() 会造成 aliasing
a) m = m.transpose().eval();
b) m = m.transposeInPlace();
component-wise operations is safe.
mat = 2 * mat;
mat = mat - MatrixXf::Identity(2, 2);
matB = matA * matA; // Simple bu not quite as efficient
matB.noalias() = matA * matA; // Use this
10. Storage orders
默认 column-major, 也可以选择 row-major
11. Quick reference guide
Eigen/LU Inverse, determinant, LU decompositions with solver (FullPivLU, PartialPivLU)
Eigen/SVD SVD decompositions with least-squares solver (JacobiSVD, BDCSVD)
Eigen/QR decomposition with solver (HouseholderQR, ColPivHouseholderQR, FullPivHouseholderQR)
...
Finished.