Search code examples
c++xtensor

xarray multiplication via operator * and raw loop give different results


I'm a beginner at xtensor and in the following code, I thought the variables result and sum should be equal but aren't. In this example result == 1000 and sum == 55000. The two variables also hold different results if I compare operations like xt::transpose(x)*A*x and its raw loop implementation (where A has compatible shape with x). How should I use xtensor operations in order to get the same results of raw loops?

#include <stdio>
#include <cmath>
#include <xtensor/xarray.hpp>
#include <xtensor/xbuilder.hpp>
#include <xtensor/xstrided_view.hpp>

int main(int argc, char* argv[]){    
    xt::xarray<double> x = xt::ones<double>({10ul,1ul});
    xt::xarray<double> b = xt::zeros<double>({10ul,1ul});
    
    for(int i = 1; i <= 10; ++i){
        b(i-1,0) = 1000*i;
    }
    double result = (xt::transpose(b)*x)(0);
    double sum = 0;
    for(int i = 0; i < 10; ++i){
        sum += b(i,0)*x(i,0);
    }

    std::cout << result << " " << sum << "\n";
    return 0;
}

Live example: https://godbolt.org/z/GMEeEzn8r (thanks to Daniel Langr in the comments)


Solution

  • There is a small issue in the following line:

        double result = (xt::transpose(b)*x)(0);
    

    Very understandably, you may assume that the multiplication of the vectors gives the sum over the pointwise product, because this is what mathematical expressions would do. However this is not what xtensor does. Fo xtensor, the expression is just the pointwise product of the two vectors, not including summation.

    To get the correct result, change the line to

        double result = xt::sum(xt::transpose(b)*x)();