#include <blitz/array.h>
#include <iostream>

int main (int argc, char *argv[])
{
    using namespace blitz;

    typedef Array<double,1> Vec;
    typedef Array<double,2> Mat;

    Mat m(3,3);
    m = 0;
    m(0,0) = 1;
    m(1,1) = 2;
    m(2,2) = 3;

    Mat v(3,1);
    v(0,0) = 3;
    v(1,0) = 6; 
    v(2,0) = 9;

    Vec v2(3);
    v2(0) = 3;
    v2(1) = 6;
    v2(2) = 9;

    Mat w(3,1);
    Vec w2(3);

    w = sum(m(tensor::i,tensor::k)*v(tensor::k,tensor::j),tensor::k);

    w2 = sum(m(tensor::i,tensor::j)*v2(tensor::j),tensor::j);

    cerr << "m: " << m;
    cerr << "v: " << v;
    cerr << "w: " << w;

    cerr << "v2: " << v2 << endl;
    cerr << "w2: " << w2 << endl;

    return 0;
} // end of main()
