#include <iostream>
#include <valarray>
#include <boost/timer.hpp>
#include <boost/numeric/ublas/vector.hpp>

const size_t size = 1024;
const int iterations = 1000000;
typedef double RealType;

template<int iterations, class vector_t>
void arithmetic_benchmark(const vector_t& v)
{
  {
    vector_t tmp = v;
    for ( int i = 0; i < iterations; ++i)
    {
      tmp += v;
      tmp -= v;
    }
  }
}

template<int iterations, class vector_t>
void arithmetic_benchmark3(const vector_t& v)
{
  {
    vector_t tmp = v;
    for ( int i = 0; i < iterations; ++i)
    {
      noalias( tmp ) += v;
      noalias( tmp ) -= v;
    }
  }
}

template<int iterations, class vector_t>
void arithmetic_benchmark2(const vector_t& v)
{
  {
    vector_t tmp = v;
    for ( int i = 0; i < iterations; ++i)
    {
      tmp.plus_assign( v );
      tmp.minus_assign( v );
    }
  }
}

int main()
{
  {
    std::valarray<RealType> v(size);
    v = 1.0;
    boost::timer timer;
    arithmetic_benchmark<iterations>(v);
    std::cout << timer.elapsed() << std::endl;
  }

  {
    boost::numeric::ublas::vector<RealType> v(size);
    std::fill(v.begin(), v.end(), 1.0);
    boost::timer timer;
    arithmetic_benchmark3<iterations>(v);
    std::cout << timer.elapsed() << std::endl;
  }

  {
    boost::numeric::ublas::vector<RealType> v(size);
    std::fill(v.begin(), v.end(), 1.0);
    boost::timer timer;
    arithmetic_benchmark2<iterations>(v);
    std::cout << timer.elapsed() << std::endl;
  }

  return 0;
}