Eigen  5.0.1
 
Loading...
Searching...
No Matches
Expression templates in Eigen

Arithmetic operators and most methods in Eigen do not compute anything by themselves. Instead they return a small expression object that merely describes the computation to be performed and stores references to the operands. The actual computation happens later, typically when the whole expression is assigned to a destination. This technique is known as expression templates, and understanding its user-visible consequences helps to write both correct and fast Eigen code.

How it works

Consider:

VectorXf a(50), b(50), c(50), d(50);
a = 3*b + 4*c + 5*d;
Matrix< float, Dynamic, 1 > VectorXf
DynamicĂ—1 vector of type float.
Definition Matrix.h:488

The right-hand side does not create any temporary vector. The subexpression 3*b is an object of type CwiseBinaryOp<scalar_product_op, ..., ...> that just stores the scalar and a reference to b; the additions nest similar objects around it. The complete type of the right-hand side encodes the whole expression tree at compile time. Only when this expression is assigned to a does Eigen generate a single evaluation loop, morally equivalent to:

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

The benefits are that the arrays are traversed only once, no temporary objects are created, and the single fused loop can be unrolled and vectorized (see Vectorization). Consequently, you should not be afraid of writing large expressions: doing so gives Eigen more opportunities for optimization. A step-by-step tour of the machinery behind this example is given in What happens inside Eigen, on a simple example, and the involved class hierarchy is described in The class hierarchy.

Not every subexpression is evaluated lazily: matrix products, for instance, are evaluated into temporaries in most contexts, both for performance and for correctness in the face of aliasing. The rules governing which subexpressions are evaluated where are documented in Lazy Evaluation and Aliasing.

What it means in practice

Because an "expression" is not a matrix, a few things deserve attention:

  • Do not use auto to capture an expression unless you know exactly what you are doing: the deduced type is the expression type, not a plain matrix, so the computation re-runs every time the variable is used, and the expression may hold dangling references to destroyed temporaries. See The auto keyword.
  • Assignment evaluates coefficient-wise into the destination, so if the destination also appears on the right-hand side, coefficients may be read after they have been overwritten. This is the aliasing problem, explained in Aliasing. Use .eval(), an xxxInPlace() method, or noalias() as appropriate.
  • Functions should accept expressions, not just matrices. Take parameters of type const MatrixBase<Derived>& (or ArrayBase, DenseBase, EigenBase, or Ref) so that arbitrary expressions can be passed without evaluating them into a temporary; see Writing Functions Taking Eigen Types as Parameters.
  • Materializing an expression is done with .eval(), which returns a plain matrix or array (and is a no-op when the expression already is one), or by assigning the expression to a plain object.
  • The two branches of the ternary operator ?: must have a common type, which two different expression types usually do not; use if/else or evaluate the branches. See Ternary operator.

Extending the expression system

Coefficient-wise custom operations rarely need a new expression type: a functor passed to unaryExpr(), binaryExpr(), or NullaryExpr() is usually enough (see Writing custom functors and Matrix manipulation via nullary-expressions). When a genuinely new expression type is required, Adding a new expression type walks through a complete example.