0

I cannot seem to understand why passing directly a function as argument to ublas::element_prod() produces a wrong result.

If I run the following code:

#include <cmath>
#include <boost/numeric/ublas/vector.hpp>
#include <boost/numeric/ublas/io.hpp>

namespace ublas = boost::numeric::ublas;

ublas::vector<double>
vector_ln(ublas::vector<double> x) {
    for(auto it = x.begin(); it != x.end(); ++it) {
        (*it) = log(*it);
    }
    return x;
}

int main(void) {
    ublas::vector<double> x(2, 2.0);
    ublas::vector<double> y(2, 10.0);

    std::cout << "x = " << x << std::endl;
    std::cout << "y = " << y << std::endl;

    auto tmp = vector_ln(y);
    auto ret1 = ublas::element_prod(
            x,
            tmp);
    std::cout << ret1 << std::endl;

    std::cout << "x = " << x << std::endl;
    std::cout << "y = " << y << std::endl;

    auto ret2 = ublas::element_prod(
            x,
            vector_ln(y));
    std::cout << ret2 << std::endl;
}

I get the following output:

x = [2](2,2)
y = [2](10,10)
[2](4.60517,4.60517)
x = [2](2,2)
y = [2](10,10)
[2](0,4.60517)

Can anybody enlighten me on why the second coding style produces a wrong result, with no compile error?

1 Answers1

1

The problem is that ublas uses Expression templates, where the result of a number of operations is a temporary which simply keeps a pointer/reference to its inputs, and is evaluated when assigning to a variable. This is done to reduce unnecessary computations and copies. See https://en.wikipedia.org/wiki/Expression_templates

However, with the introduction of C++11, it creates a dangerous interaction with the use of auto, as this saves a copy of the expression template, and not the result. This expression template has dangling references to the temporary returned by vector_ln(y), and causes the problem you are seeing.

Since the primary problem is the interaction with auto, the solution is to save it to the correct ublas matrix type as the result of element_prod(). It only worked in the first case because none of the stored references were temporaries.

Dave S
  • 20,507
  • 3
  • 48
  • 68
  • Thank you Dave! Indeed, by declaring ret2's type explicitly the problem is fixed. Does this mean that using uBLAS in a C++11 code-base is generally not a good idea? – Michele Petracca Aug 24 '19 at 14:25
  • It's not bad, but you have to understand you really shouldn't use 'auto' with its operations. – Dave S Aug 24 '19 at 21:50