#include <vector>
#include <math.h>
#include "StlVector.h"

namespace AMDiS {

  double norm(::std::vector<double> *vec) {
    double v = 0;
    for (int i = 0; i < static_cast<int>(vec->size()); i++) {
      v += (*vec)[i] * (*vec)[i];
    }
    return sqrt(v);
  };

  void setValue(::std::vector<double>& x, double value) {  
    for (int i = 0; i < static_cast<int>(x.size()); i++) {  
      x[i] = value;  
    }  
  };   
   
  const ::std::vector<double>& operator*=(::std::vector<double>& x, double scal) {
    for (int i = 0; i < static_cast<int>(x.size()); i++) {
      x[i] *= scal;
    }

    return x;
  };

  const ::std::vector<double>& operator+=(::std::vector<double>& x1,
					  const ::std::vector<double>& x2) {
    for (int i = 0; i < static_cast<int>(x1.size()); i++) {
      x1[i] += x2[i];
    }

    return x1;
  };

  double operator*(::std::vector<double>& x1, ::std::vector<double>& x2) {
    double result = 0.0;

    for (int i = 0; i < static_cast<int>(x1.size()); i++) {
      result += x1[i] * x2[i];
    }

    return result;
  };

  void axpy(double alpha, const ::std::vector<double>& x, ::std::vector<double>& y) {
    for (int i = 0; i < static_cast<int>(y.size()); i++) {
      y[i] += alpha * x[i];
    }
  }

  void xpay(double alpha, const ::std::vector<double>& x, ::std::vector<double>& y) {
    for (int i = 0; i < static_cast<int>(y.size()); i++) {
      y[i] = alpha * y[i] + x[i];
    }
  }

  int logToPhysIndex(::std::vector< ::std::vector<MatEntry> > *matrix, int a, int b)
  {
    int j;

    for (j = 0; j < static_cast<int>((*matrix)[a].size()); j++) 
      if (b == (*matrix)[a][j].col) 
	break;

    return (j == static_cast<int>((*matrix)[a].size())) ? -1 : j;
  }

  void print(::std::vector<double> *vec)
  {
    ::std::cout << "[";

    for (int i = 0; i < static_cast<int>(vec->size()); i++)  {
      ::std::cout << (*vec)[i] << " , ";
    }

    ::std::cout << "]" << ::std::endl;
  }
}