How to integrate a library that uses expression templates?

Question

I would like to use the Eigen matrix library as the linear algebra engine in my program. Eigen uses expression templates to implement lazy evaluation and to simplify loops and calculations.

For example:

#include<Eigen/Core>

int main()
{
  int size = 40;
  // VectorXf is a vector of floats, with dynamic size.
  Eigen::VectorXf u(size), v(size), w(size), z(size);
  u = 2*v + w + 0.2*z;
}

Since Eigen uses expression templates, code like

u = 2*v + w + 0.2*z;

In the above mentioned sample reduces to a single loop of length 10 (not 40, the floats are put into regiser by chunks of 4) without creating a temporary. How cool is that?

But if I integrate the library like this:

class UsingEigen
{
    public:  
        UsingEigen(const Eigen::VectorXf& data):
            data_(data)
        {}

        UsingEigen operator + (const UsingEigen& adee)const
        {
            return UsingEigen(data_ + adee.data_);
        }

        ...
    private:
        Eigen::VectorXf data_;
}

Then the expressions like:

UsingEigen a, b, c, d;
a = b + c + d;

cannot take advantage of the way Eigen is implemented. And this is not the last of it. There are many other examples, where expression templates are used in Eigen.

The easy solution would be not to define the operators by myself, make data_ public and just write expressions like:

UsingEigen a, b, c, d;
a.data_ = b.data_ + c.data_ + d.data_;

This breaks encapsulation, but it preserves the efficiency of Eigen.

Other way could be to make my own operators, but let them return expression templates. But since I am a beginner in C++, I do not know if this is the right way to go.

I am sorry if the question is too general in nature. I am a beginner and have noone to ask. Up until now I was using std::vector<float> everywhere, but now I need to use matrices also. To switch from std::vector<float> to Eigen in my entire project is a big step and I am afraid of making a wrong call right in the start. Any advice is welcomed!

Answer 1

Why would exposing data_ break encapsulation? Encapsulation means hiding the implementation details and only exposing the interface. If your wrapper class UsingEigen does not add any behavior or state to the native Eigen library, the interface does not change. In this case, you should drop this wrapper altogether and write your program using the Eigen data structures.

Exposing a matrix or a vector is not breaking encapsulation: only exposing the implementation of the matrix or vector would do that. The Eigen library exposes the arithmetic operators but not their implementation.

With expression template libraries, the most common way for users to extend the library functionality is by adding behavior, not adding by adding state. And for adding behavior you do not need to write wrapper classes: you can also add non-member functions that are implemented in terms of the Eigen class member functions. See this column "How Non-Member Functions Improve Encapsulation" by Scott Meyers.

As for your concern that the transformation of your current program to a version that explicitly uses the Eigen functionality: you can perform the change step-by-step, changing small parts of your program each time, making sure your unit tests (you do have unit tests, don't you?) do not break as you go along.

Answer 2

In my opinion, this looks more of a object oriented design problem rather than a library usage problem. Whatever you read from the books are the right recommendations. i.e, do not expose member variables and shield the upper layers from the nuances of the 3rd party layer usage.

What you could look forward is right abstractions of mathematical functions that can be implemented using this library internally. i.e, you could expose a library of your own with high level functions than elementary vector and matrix operations. In this way you can utilize the peculiarities of the interactions among the library objects and at the same time you don't have to expose your member variables to upper layers.

For e.g you could abstract away my higher level APIs like computing the distance from a point to a plane, distance between two planes, computing the new coordinates of a point w.r.t another coordinate system using the transformation matrices etc. To implement these methods internally you can utilize the library objects. You can restrict to not to have any of the library classes used in the API signatures to avoid dependency for the upper layers on this library.

Upper layers of your program should be higher in the level of abstraction and need not bother about the elementary implementation details such as how the calculation of the distance from a point to the plane is implemented etc. Also, they even need not know if this lower layer is implemented using this library or something else. They would just use the interfaces of your library.

Answer 3

Set up a class template to hold general Eigen expressions and make UsingEigen a special instance of it:

template<typename expr_t>
class UsingEigenExpr
{
    UsingEigen(expr_t const& expr) : expr(expr) {}
    expr_t expr;
    operator UsingEigenExpr<Eigen::VectorXf>() const
    {
        return {expr};
    }
};
using UsingEigen = UsingEigenExpr<Eigen::VectorXf>;

Then overload any required function, e.g. as

template<typename expr1_t, typename expr2_t, typename function_t>
auto binary_op(UsingEigenExpr<expr1_t> const& x, UsingEigenExpr<expr2_t> const& y, function_t function)
{
    return UsingEigenExpr<decltype(function(std::declval<expr1_t>(),std::declval<expr2_t>()))>(function(x.expr,y.expr));
}

template<typename expr1_t, typename expr2_t>
auto operator+(UsingEigenExpr<expr1_t> const& x, UsingEigenExpr<expr2_t> const& y)
{
    return binary_op(x,y,[](auto const& x, auto const& y) {return x+y;});
}

and so on for other binary operators like operator-, for unary operators, and more generally for all the other stuff you want to use. Further, you could add some other member functions to UsingEigenExpr, e.g. size(), norm(), etc.

Use it as

UsingEigen b, c, d;
auto a = b + c + d;

to store the expression, or

UsingEigen b, c, d;
UsingEigen a = b + c + d;

to directly evaluate it.

Although this approach works, in the end you find yourself duplicating all the required functionality, so use it carefully.

Answer 4

I don't understand all your question I will try to answer you most of them. In this sentence:

 UsingEigen operator + (const UsingEigen& adee)const
    {
        return UsingEigen(data_ + adee.data_);
    }

You have an overload operator ( sorry I don't know if this is the correct way to write in English), for this reason you can write:

a = b + c + d;

instead of:

a.data_ = b.data_ + c.data_ + d.data_;

You won't have any problem, cost of your program will be the same. In addition you will have encapsulation and efficiency.

On the other way if you want define your own operator you can do it like the template do it. You can find information on the web searching "overload operator" but is similar to this:

UsingEigen operator + (const UsingEigen& adee)const
    {
        return UsingEigen(data_ + adee.data_);
    }

Instead of "+" you can put the operator and do the operations you need.

If you want to create a matrix it is simple. You only need to create a array of array or vector of vector.

I think is something like this:

std::vector<vector<float>>

I am not sure but it is easy, on the other hand you can use a simple matrix on this way:

float YourMatrix [size][size];

I hope it could help you. I don't understand all your question if you need something more add me on google+ and I will try to help you.

Sorry for my English, I hope you can understand all and it helps you.