Dense matrix and array manipulation » Reductions, visitors and broadcasting module

Norm computations

The (Euclidean a.k.a. ) squared norm of a vector can be obtained squaredNorm(). It is equal to the dot product of the vector by itself, and equivalently to the sum of squared absolute values of its coefficients.

Eigen also provides the norm() method, which returns the square root of squaredNorm().

These operations can also operate on matrices; in that case, a n-by-p matrix is seen as a vector of size (n*p), so for example the norm() method returns the "Frobenius" or "Hilbert-Schmidt" norm. We refrain from speaking of the norm of a matrix because that can mean different things.

If you want other coefficient-wise norms, use the lpNorm<p>() method. The template parameter p can take the special value Infinity if you want the norm, which is the maximum of the absolute values of the coefficients.

The following example demonstrates these methods.

#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;
  
  m << 1,-2,
       -3,4;

  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;
  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;
}

v.squaredNorm() = 5
v.norm() = 2.23607
v.lpNorm<1>() = 3
v.lpNorm<Infinity>() = 2

m.squaredNorm() = 30
m.norm() = 5.47723
m.lpNorm<1>() = 10
m.lpNorm<Infinity>() = 4

Operator norm: The 1-norm and -norm matrix operator norms can easily be computed as follows:

#include <Eigen/Dense>
#include <iostream>

using namespace Eigen;
using namespace std;

int main()
{
  MatrixXf m(2,2);
  m << 1,-2,
       -3,4;

  cout << "1-norm(m)     = " << m.cwiseAbs().colwise().sum().maxCoeff()
       << " == "             << m.colwise().lpNorm<1>().maxCoeff() << endl;

  cout << "infty-norm(m) = " << m.cwiseAbs().rowwise().sum().maxCoeff()
       << " == "             << m.rowwise().lpNorm<1>().maxCoeff() << endl;
}

1-norm(m)     = 6 == 6
infty-norm(m) = 7 == 7

See below for more explanations on the syntax of these expressions.

Boolean reductions

The following reductions operate on boolean values:

• all() returns true if all of the coefficients in a given Matrix or Array evaluate to true .
• any() returns true if at least one of the coefficients in a given Matrix or Array evaluates to true .
• count() returns the number of coefficients in a given Matrix or Array that evaluate to true.

These are typically used in conjunction with the coefficient-wise comparison and equality operators provided by Array. For instance, array > 0 is an Array of the same size as array , with true at those positions where the corresponding coefficient of array is positive. Thus, (array > 0).all() tests whether all coefficients of array are positive. This can be seen in the following example:

#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

User defined reductions

TODO

In the meantime you can have a look at the DenseBase::redux() function.

Combining partial reductions with other operations

It is also possible to use the result of a partial reduction to do further processing. Here is another example that finds the column whose sum of elements is the maximum within a matrix. With column-wise partial reductions this can be coded as:

#include <iostream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;
int main()
{
  MatrixXf mat(2,4);
  mat << 1, 2, 6, 9,
         3, 1, 7, 2;
  
  MatrixXf::Index   maxIndex;
  float maxNorm = mat.colwise().sum().maxCoeff(&maxIndex);
  
  std::cout << "Maximum sum at position " << maxIndex << std::endl;

  std::cout << "The corresponding vector is: " << std::endl;
  std::cout << mat.col( maxIndex ) << std::endl;
  std::cout << "And its sum is is: " << maxNorm << std::endl;
}

Maximum sum at position 2
The corresponding vector is: 
6
7
And its sum is is: 13

The previous example applies the sum() reduction on each column though the colwise() visitor, obtaining a new matrix whose size is 1x4.

Therefore, if

then

The maxCoeff() reduction is finally applied to obtain the column index where the maximum sum is found, which is the column index 2 (third column) in this case.

Combining broadcasting with other operations

Broadcasting can also be combined with other operations, such as Matrix or Array operations, reductions and partial reductions.

Now that broadcasting, reductions and partial reductions have been introduced, we can dive into a more advanced example that finds the nearest neighbour of a vector v within the columns of matrix m. The Euclidean distance will be used in this example, computing the squared Euclidean distance with the partial reduction named squaredNorm():

#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 << m.col(index) << endl;
}

Nearest neighbour is column 0:
1
3

The line that does the job is (m.colwise() - v).colwise().squaredNorm().minCoeff(&index);

We will go step by step to understand what is happening:

• m.colwise() - v is a broadcasting operation, subtracting v from each column in m. The result of this operation is a new matrix whose size is the same as matrix m:

• (m.colwise() - v).colwise().squaredNorm() is a partial reduction, computing the squared norm column-wise. The result of this operation is a row vector where each coefficient is the squared Euclidean distance between each column in m and v:
• Finally, minCoeff(&index) is used to obtain the index of the column in m that is closest to v in terms of Euclidean distance.