closure-test.cc (flens/examples/closure-test.cc)

      1
      2
      3
      4
      5
      6
      7
      8
      9
     10
     11
     12
     13
     14
     15
     16
     17
     18
     19
     20
     21
     22
     23
     24
     25
     26
     27
     28
     29
     30
     31
     32
     33
     34
     35
     36
     37
     38
     39
     40
     41
     42
     43
     44
     45
     46
     47
     48
     49
     50
     51
     52
     53
     54
     55
     56
     57
     58
     59
     60
     61
     62
     63
     64
     65
     66
     67
     68
     69
     70
     71
     72
     73
     74
     75
     76
     77

#include <flens/flens.h>

template <typename INT, typename T, typename VX, typename VY, typename VZ>
void
myComputationalRoutine(INT n,
                       const T &alpha,
                       const VX *x, INT incX,
                       const VY *y, INT incY,
                       VZ  *z, INT incZ)
{
    for (INT i=0, iX=0, iY=0, iZ=0; i<n; ++i, iX+=incX, iY+=incY, iZ+=incZ) {
        z[iZ] = alpha*(x[iX] + y[iY]);
    }
}

namespace flens { namespace blas {

//
//  z = alpha*(x+y)
//

template <typename T, typename VX, typename VY>
using MyCrazyVectorClosure =

    VectorClosure<OpMult,
                  ScalarValue<T>,
                  VectorClosure<OpAdd,
                                DenseVector<VX>,
                                DenseVector<VY>
                                >
                 >;

template <typename T, typename VX, typename VY, typename VZ>
void
copy(const MyCrazyVectorClosure<T, VX, VY> &expr, DenseVector<VZ> &z)
{
    const auto alpha = expr.left().value();
    const auto x     = expr.right().left();
    const auto y     = expr.right().right();

    ASSERT(x.length()==z.length());
    ASSERT(y.length()==z.length());

    std::cerr << "calling myComputationalRoutine" << std::endl;

    myComputationalRoutine(z.length(), alpha,
                           x.data(), x.stride(),
                           y.data(), y.stride(),
                           z.data(), z.stride());
}

} } // namespace blas, flens

#include <flens/flens.cxx>

using namespace flens;
using namespace std;

int
main()
{
    typedef DenseVector<Array<double> >   RealDenseVector;

    RealDenseVector      x(5), y(5), z(5);

    x = 1;
    y = 2;
    z = 0;

    z = 3*(x +y);

    cout << "z = " << z << endl;
}