简介
本文介绍一下Dense Matrix的3中操作: reduction, visitor, broadcasting.
归约计算reduction.
Eigen的归约计算是这样的一类计算,它是对矩阵或者数组进行的计算,并返回一个标量。可能其中最常用的就是sum()
计算,它返回的是系数的和。
比如,前面曾经用到的示例:
#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
范式计算 (Norm comptation)
以下列出常用的向量范数和矩阵范数的定义,以及对应的MATLAB的函数:
1、向量范数
- 1-范数:,即向量元素绝对值之和,matlab调用函数norm(x, 1) 。
- 2-范数:,Euclid范数(欧几里得范数,常用计算向量长度),即向量元素绝对值的平方和再开方,matlab调用函数norm(x, 2)。
- ∞-范数:,即所有向量元素绝对值中的最大值,matlab调用函数norm(x, inf)。
- -∞-范数:,即所有向量元素绝对值中的最小值,matlab调用函数norm(x, -inf)。
- p-范数:,即向量元素绝对值的p次方和的1/p次幂,matlab调用函数norm(x, p)。
2、矩阵范数
- 1-范数:, 列和范数,即所有矩阵列向量绝对值之和的最大值,matlab调用函数norm(A, 1)。
- 2-范数:,谱范数,即A’A矩阵的最大特征值的开平方。matlab调用函数norm(x, 2)。
- ∞-范数:,行和范数,即所有矩阵行向量绝对值之和的最大值,matlab调用函数norm(A, inf)。
- F-范数:,Frobenius范数,即矩阵元素绝对值的平方和再开平方,matlab调用函数norm(A, ’fro‘)。
在Eigen中,提供了计算范式的一些方法。
- squaredNorm 向量的平方范式,2-范式,
ℓ
2
\ell^2
ℓ2的范式,可以使用方法
squaredNorm()
,即所有系数的平方之和,数学上是一个向量与自身的点乘。 - norm 即称欧拉范式,squaredNorm结果再开平方,得到平方根。
当这些计算作用于矩阵时, 一个 nXp的矩阵,可以当做一个大小为(n*p)的向量,比如"Frobenius" or "Hilbert-Schmidt" norm
。
我们避免提矩阵的 ℓ 2 \ell^2 ℓ2的范式(欧拉平方范式)运算,这样会导致歧义。
Eigen也额外提供了一些对系数进行欧拉范式操作方法,如
ℓ
p
\ell^p
ℓp范式,可以使用lpNorm<p>()
方法。如果要计算
ℓ
∞
\ell^\infty
ℓ∞ 范式,模板参数p取特殊值Infinity
即可,这得到系数绝对值的最大值。
示例:
//matrix_norm1.cpp
#include <Eigen/Dense>
#include <iostream>
using namespace std;
using namespace Eigen;
int main()
{
VectorXf v(2);
MatrixXf m(2,2), n(2,2);
v << -1,
2;
cout << "vector v: "<<endl<<v<<endl;
cout << "v.squaredNorm() = " << v.squaredNorm() << endl;
cout << "v.norm() = " << v.norm() << endl;
cout << "v.lpNorm<1>() = " << v.lpNorm<1>() << endl;
cout << "v.lpNorm<Infinity>() = " << v.lpNorm<Infinity>() << endl;
cout << endl;
m << 1,-2,
-3,4;
cout << "matrix m: "<<endl<<m<<endl;
cout << "m.squaredNorm() = " << m.squaredNorm() << endl;
cout << "m.norm() = " << m.norm() << endl;
cout << "m.lpNorm<1>() = " << m.lpNorm<1>() << endl;
cout << "m.lpNorm<Infinity>() = " << m.lpNorm<Infinity>() << endl;
}
执行结果:
$ g++ -I /usr/local/include/eigen3 matrix_norm1.cpp -o matrix_norm1
$ ./matrix_norm1
vector v:
-1
2
v.squaredNorm() = 5
v.norm() = 2.23607
v.lpNorm<1>() = 3
v.lpNorm<Infinity>() = 2
matrix m:
1 -2
-3 4
m.squaredNorm() = 30
m.norm() = 5.47723
m.lpNorm<1>() = 10
m.lpNorm<Infinity>() = 4
boolean 归约
Eigen提供了一些针对Boolean运算的归约操作。
- all() : 如果系数全部满足Boolean计算结果,则结果为true;否则为false。
- any() : 如果系数有任何满足Boolean计算的结果,则结果为true;否则为false。
- count() : 返回满足条件的系数的数。
给一个示例:
#include <Eigen/Dense>
#include <iostream>
using namespace std;
using namespace Eigen;
int main()
{
ArrayXXf a(2,2);
a << 1,2,
3,4;
cout << "(a > 0).all() = " << (a > 0).all() << endl;
cout << "(a > 0).any() = " << (a > 0).any() << endl;
cout << "(a > 0).count() = " << (a > 0).count() << endl;
cout << endl;
cout << "(a > 2).all() = " << (a > 2).all() << endl;
cout << "(a > 2).any() = " << (a > 2).any() << endl;
cout << "(a > 2).count() = " << (a > 2).count() << endl;
}
结果:
(a > 0).all() = 1
(a > 0).any() = 1
(a > 0).count() = 4
(a > 2).all() = 0
(a > 2).any() = 1
(a > 2).count() = 2
访问子 visitor
访问子visitor用于访问矩阵或数组内的系数及系数的位置。最简单的是前面曾经使用到的maxCoeff(&x,&y) , minCoeff(&x,&y)
。
示例请看下面的Partial reduction中的例子程序。
Partial reductions 局部归约
Eigen提供了针对矩阵的列或者行的归约操作,它们对矩阵的每行/每列进行操作,并返回结果向量。
对应的函数为colwise() , rowwise()
。
Note colwise()
返回的是一个行向量,而rowwise()
返回的是一个列向量。
示例:
//matrix_norm2.cpp
#include <iostream>
#include <Eigen/Dense>
using namespace std;
int main()
{
Eigen::MatrixXf mat(2,4);
mat << 1, 2, 6, 9,
3, 1, 7, 2;
std::cout << "Column's maximum: " << std::endl
<< mat.colwise().maxCoeff() << std::endl;
std::cout << "Row's maximum: " << std::endl
<< mat.rowwise().maxCoeff() << std::endl;
}
执行:
$ g++ -I /usr/local/include/eigen3 matrix_norm2.cpp -o matrix_norm2
$ ./matrix_norm2
Column's maximum:
3 2 7 9
Row's maximum:
9
7
broadcasting
广播broadcasting比Partial reduction多了一个操作。其把一个向量赋值多次,结合成一个矩阵,可以和以及矩阵进行运算。
下面是一个简单的示例程序:
//matrix_broadcast.cpp
#include <iostream>
#include <Eigen/Dense>
using namespace std;
int main()
{
Eigen::MatrixXf mat(2,4);
Eigen::VectorXf v(2);
mat << 1, 2, 6, 9,
3, 1, 7, 2;
v << 0,
1;
//add v to each column of m
mat.colwise() += v;
std::cout << "Broadcasting result1: " << std::endl;
std::cout << mat << std::endl;
// -------------
Eigen::MatrixXf mm(2,2);
mm<<1,2,3,4;
v << 1,1;
mm.rowwise() += v.transpose();
std::cout << "Broadcasting result2: " << std::endl;
std::cout << mm << std::endl;
}
执行的情况类似:
[
1
2
6
9
3
1
7
2
]
+
[
0
0
0
0
1
1
1
1
]
=
[
1
2
6
9
4
2
8
3
]
.
\begin{bmatrix} 1 & 2 & 6 & 9 \\ 3 & 1 & 7 & 2 \end{bmatrix} + \begin{bmatrix} 0 & 0 & 0 & 0 \\ 1 & 1 & 1 & 1 \end{bmatrix} = \begin{bmatrix} 1 & 2 & 6 & 9 \\ 4 & 2 & 8 & 3 \end{bmatrix}.
[13216792]+[01010101]=[14226893].
[ 1 2 3 4 ] + [ 1 1 1 1 ] = [ 2 3 4 5 ] . \begin{bmatrix} 1 & 2 \\ 3 & 4 \end{bmatrix} + \begin{bmatrix} 1 & 1 \\ 1 & 1 \end{bmatrix} = \begin{bmatrix} 2 & 3 \\ 4 & 5 \end{bmatrix}. [1324]+[1111]=[2435].
执行结果:
$ g++ -I /usr/local/include/eigen3 matrix_broadcast.cpp -o matrix_broadcast
$
$ ./matrix_broadcast
Broadcasting result1:
1 2 6 9
4 2 8 3
Broadcasting result2:
2 3
4 5
综合: 结合broadcasting和其他的操作及运算
//matrix_broadcast2.cpp
#include <iostream>
#include <Eigen/Dense>
using namespace std;
using namespace Eigen;
int main()
{
Eigen::MatrixXf m(2,4);
Eigen::VectorXf v(2);
m << 1, 23, 6, 9,
3, 11, 7, 2;
v << 2,
3;
MatrixXf::Index index;
// find nearest neighbour
(m.colwise() - v).colwise().squaredNorm().minCoeff(&index);
cout << "Nearest neighbour is column " << index << ":" << endl;
cout << "Nearest column vector: " << endl << m.col(index) << endl;
}
-
计算
m.colwise() - v
:
[ 1 23 6 9 3 11 7 2 ] − [ 2 2 2 2 3 3 3 3 ] = [ − 1 21 4 7 0 8 4 − 1 ] \begin{bmatrix} 1 & 23 & 6 & 9 \\ 3 & 11 & 7 & 2 \end{bmatrix} - \begin{bmatrix} 2 & 2 & 2 & 2 \\ 3 & 3 & 3 & 3 \end{bmatrix} = \begin{bmatrix} -1 & 21 & 4 & 7 \\ 0 & 8 & 4 & -1 \end{bmatrix} [1323116792]−[23232323]=[−10218447−1] -
计算
(m.colwise() - v).colwise().squaredNorm()
:
(m.colwise() - v).colwise().squaredNorm() = [ \begin{bmatrix} 1 & 505 & 32 & 50 \end{bmatrix} ]
- minCoeff(&index) 用来获得在m矩阵内,最小的系数位置索引,该列向量与向量v具有最短欧拉距离。
执行结果:
$ g++ -I /usr/local/include/eigen3 matrix_broadcast2.cpp -o matrix_broadcast2
$ ./matrix_broadcast2
Nearest neighbour is column 0:
Nearest column vector:
1
3
[ 4 3 13 11 4 3 13 11 4 3 13 11 ] \begin{bmatrix} 4 & 3 & 13 & 11 \\ 4 & 3 & 13 & 11\\4 & 3 & 13 & 11 \end{bmatrix} ⎣⎡444333131313111111⎦⎤
[ 4 3 13 11 ] \begin{bmatrix} 4 & 3 & 13 & 11 \end{bmatrix} [431311]