DenseVector2dField2d.hpp

Go to the documentation of this file.

/*
 * Copyright(c):
 * Signal Image and Communications (SIC) Department 
 * http://www.sic.sp2mi.univ-poitiers.fr/ 
 *  - University of Poitiers, France http://www.univ-poitiers.fr 
 *  - XLIM  Institute UMR CNRS 7252  http://www.xlim.fr/ 
 * 
 * and 
 * 
 * D2 Fluid, Thermic and Combustion 
 *  - University of Poitiers, France http://www.univ-poitiers.fr 
 *  - PPRIME Institute -  UPR CNRS 3346  http://www.pprime.fr
 *  - ISAE-ENSMA http://www.ensma.fr
 *
 * Contributor(s):
 * The SLIP team,
 * Benoit Tremblais <tremblais_AT_sic.univ-poitiers.fr>,
 * Laurent David <laurent.david_AT_lea.univ-poitiers.fr>, 
 * Ludovic Chatellier <ludovic.chatellier_AT_univ-poitiers.fr>, 
 * Lionel Thomas <lionel.thomas_AT_univ-poitiers.fr>, 
 * Denis Arrivault <arrivault_AT_sic.univ-poitiers.fr>, 
 * Julien Dombre <julien.dombre_AT_univ-poitiers.fr>.
 *
 * Description:
 * The Simple Library of Image Processing (SLIP) is a new image processing 
 * library. It is written in the C++ language following as much as possible 
 * the ISO/ANSI C++ standard. It is consequently compatible with any system  
 * satisfying the ANSI C++ complience. It works on different Unix , Linux ,  
 * Mircrosoft Windows and Mac OS X plateforms. SLIP is a research library that  
 * was created by the Signal, Image and Communications (SIC) departement of  
 * the XLIM, UMR 7252 CNRS Institute in collaboration with the Fluids, Thermic  
 * and Combustion departement of the P', UPR 3346 CNRS Institute of the  
 * University of Poitiers.
 * 
 * The SLIP Library source code has been registered to the APP (French Agency 
 * for the Protection of Programs) by the University of Poitiers and CNRS, 
 * under  registration number IDDN.FR.001.300034.000.S.P.2010.000.21000.

 * http://www.sic.sp2mi.univ-poitiers.fr/slip/
 *
 * This software is governed by the CeCILL-C license under French law and
 * abiding by the rules of distribution of free software.  You can  use, 
 * modify and/ or redistribute the software under the terms of the CeCILL-C 
 * license as circulated by CEA, CNRS and INRIA at the following URL 
 * http://www.cecill.info.
 * As a counterpart to the access to the source code and  rights to copy,
 * modify and redistribute granted by the license, users are provided only
 * with a limited warranty  and the software's author,  the holder of the
 * economic rights,  and the successive licensors  have only  limited
 * liability.
 * 
 * In this respect, the user's attention is drawn to the risks associated
 * with loading,  using,  modifying and/or developing or reproducing the
 * software by the user in light of its specific status of free software,
 * that may mean  that it is complicated to manipulate,  and  that  also
 * therefore means  that it is reserved for developers  and  experienced
 * professionals having in-depth computer knowledge. Users are therefore
 * encouraged to load and test the software's suitability as regards their
 * requirements in conditions enabling the security of their systems and/or
 * data to be ensured and,  more generally, to use and operate it in the 
 * same conditions as regards security. 
 *
 * The fact that you are presently reading this means that you have had
 * knowledge of the CeCILL-C license and that you accept its terms.
 */


#ifndef SLIP_DENSEVECTOR2DFIELD2D_HPP
#define SLIP_DENSEVECTOR2DFIELD2D_HPP

#include <iostream>
#include <fstream>
#include <iterator>
#include <cassert>
#include <numeric>
#include <algorithm>
#include <cmath>
#include <string>
#include <cstddef>
#include "Matrix.hpp"
#include "stride_iterator.hpp"
#include "kstride_iterator.hpp"
#include "iterator2d_box.hpp"
#include "iterator2d_range.hpp"
#include "apply.hpp"
#include "Vector2d.hpp"
#include "GenericMultiComponent2d.hpp"
#include "derivatives.hpp"
#include "arithmetic_op.hpp"
#include "io_tools.hpp"
#include "linear_algebra_traits.hpp"

#include <boost/serialization/access.hpp>
#include <boost/serialization/split_member.hpp>
#include <boost/serialization/version.hpp>
#include <boost/serialization/base_object.hpp>

namespace slip
{

template<typename T>
class stride_iterator;

template<typename T>
class iterator2d_box;

template<typename T>
class iterator2d_range;

template<typename T>
class const_iterator2d_box;

template<typename T>
class const_iterator2d_range;


template <typename T>
class DPoint2d;

template <typename T>
class Point2d;

template <typename T>
class Box2d;

template <typename T>
class DenseVector2dField2d;


template <typename T>
class DenseVector2dField2d:public slip::GenericMultiComponent2d<slip::Vector2d<T> >
{
public :

  typedef slip::Vector2d<T> value_type;
  typedef DenseVector2dField2d<T> self;
  typedef const DenseVector2dField2d<T> const_self;
  typedef GenericMultiComponent2d<slip::Vector2d<T> > base;


  typedef value_type* pointer;
  typedef const value_type* const_pointer;
  typedef value_type& reference;
  typedef const value_type& const_reference;

  typedef ptrdiff_t difference_type;
  typedef std::size_t size_type;

  typedef pointer iterator;
  typedef const_pointer const_iterator;
  
  typedef std::reverse_iterator<iterator> reverse_iterator;
  typedef std::reverse_iterator<const_iterator> const_reverse_iterator;


  typedef pointer row_iterator;
  typedef const_pointer const_row_iterator;
  typedef slip::stride_iterator<pointer> col_iterator;
  typedef slip::stride_iterator<const_pointer> const_col_iterator;

  typedef std::reverse_iterator<iterator> reverse_row_iterator;
  typedef std::reverse_iterator<const_iterator> const_reverse_row_iterator;
  typedef std::reverse_iterator<col_iterator> reverse_col_iterator;
  typedef std::reverse_iterator<const_col_iterator> const_reverse_col_iterator;

  typedef typename slip::GenericMultiComponent2d<slip::Vector2d<T> >::iterator2d iterator2d;
  typedef typename slip::GenericMultiComponent2d<slip::Vector2d<T> >::const_iterator2d const_iterator2d;

  typedef T vector2d_value_type;
  typedef vector2d_value_type* vector2d_pointer;
  typedef const vector2d_value_type* const_vector2d_pointer;
  typedef vector2d_value_type& vector2d_reference;
  typedef const vector2d_value_type const_vector2d_reference;
  typedef slip::kstride_iterator<vector2d_pointer,2> vector2d_iterator;
  typedef slip::kstride_iterator<const_vector2d_pointer,2> const_vector2d_iterator;

  typedef typename slip::lin_alg_traits<value_type>::value_type norm_type;

  static const std::size_t DIM = 2;
public:

  DenseVector2dField2d():
    base()
  {}
 
  DenseVector2dField2d(const size_type height,
                     const size_type width):
    base(height,width)
  {}
  
  DenseVector2dField2d(const size_type height,
                     const size_type width, 
                     const slip::Vector2d<T>& val):
    base(height,width,val)
  {}
 
  DenseVector2dField2d(const size_type height,
                     const size_type width, 
                     const T* val):
    base(height,width,val)
  {}
 
  DenseVector2dField2d(const size_type height,
                                 const size_type width, 
                     const slip::Vector2d<T>* val):
    base(height,width,val)
  {}
 

  template<typename InputIterator>
  DenseVector2dField2d(const size_type height,
                                       const size_type width, 
                                       InputIterator first,
                                       InputIterator last):
    base(height,width,first,last)
  {} 



  
  template<typename InputIterator>
  DenseVector2dField2d(const size_type height,
                                       const size_type width, 
                                       InputIterator first1,
                                       InputIterator last1,
                                       InputIterator first2):
    base(height,width)
  {
    
    std::vector<InputIterator> first_iterators_list(2);
    first_iterators_list[0] = first1;
    first_iterators_list[1] = first2;
    this->fill(first_iterators_list,last1);
  
  } 


  DenseVector2dField2d(const self& rhs):
    base(rhs)
  {}

  ~DenseVector2dField2d()
  {}
  

  std::string name() const;


  
  void write_gnuplot(const std::string& file_path_name);


  void write_tecplot(const std::string& file_path_name,
                                     const std::string& title,
                                     const std::string& zone);


  void read_gnuplot(const std::string& file_path_name);


  void read_tecplot(const std::string& file_path_name);


 

 

  self& operator=(const slip::Vector2d<T>& val);


  self& operator=(const T& val);



  
 
  T& Vx1(const size_type i,
                 const size_type j);

 
 
  const T& Vx1(const size_type i,
                       const size_type j) const;

 
  T& Vx2(const size_type i,
                 const size_type j);
 
  const T& Vx2(const size_type i,
                       const size_type j) const;

  T& u(const size_type i,
                 const size_type j);

 
 
  const T& u(const size_type i,
                       const size_type j) const;

  
  T& v(const size_type i,
       const size_type j);
 
  const T& v(const size_type i,
                       const size_type j) const;


  norm_type norm(const size_type i,
                                 const size_type j) const;

  norm_type angle(const size_type i,
                                  const size_type j) const;
  
 

  self& operator+=(const T& val);
  self& operator-=(const T& val);
  self& operator*=(const T& val);
  self& operator/=(const T& val);

//   self& operator%=(const T& val);
//   self& operator^=(const T& val);
//   self& operator&=(const T& val);
//   self& operator|=(const T& val);
//   self& operator<<=(const T& val);
//   self& operator>>=(const T& val);


   self  operator-() const;
//   self  operator!() const;

  self& operator+=(const slip::Vector2d<T>& val);
  self& operator-=(const slip::Vector2d<T>& val);
  self& operator*=(const slip::Vector2d<T>& val);
  self& operator/=(const slip::Vector2d<T>& val);

  

  self& operator+=(const self& rhs);
  self& operator-=(const self& rhs);
  self& operator*=(const self& rhs);
  self& operator/=(const self& rhs);
    




  template<typename Container2D>
  void norm(Container2D& result)
  {
    //typedef slip::Vector2d<double> vector;
        std::transform(this->begin(),this->end(),result.begin(),std::mem_fun_ref(&value_type::Euclidean_norm));

  }

  template<typename Container2D>
  void angle(Container2D& result)
  {
    
        std::transform(this->begin(),this->end(),result.begin(),std::mem_fun_ref(&slip::Vector2d<T>::angle));

  }

private :
   friend class boost::serialization::access;
 template<class Archive>
    void save(Archive & ar, const unsigned int version) const
    {
      ar & boost::serialization::base_object<slip::GenericMultiComponent2d<slip::Vector2d<T> > >(*this);
    }
   template<class Archive>
    void load(Archive & ar, const unsigned int version)
    {
      ar & boost::serialization::base_object<slip::GenericMultiComponent2d<slip::Vector2d<T> > >(*this);
    }
  BOOST_SERIALIZATION_SPLIT_MEMBER()
 };

  typedef slip::DenseVector2dField2d<double> DenseVector2dField2d_d;
  typedef slip::DenseVector2dField2d<float> DenseVector2dField2d_f;
  typedef slip::DenseVector2dField2d<long> DenseVector2dField2d_l;
  typedef slip::DenseVector2dField2d<unsigned long> DenseVector2dField2d_ul;
  typedef slip::DenseVector2dField2d<short> DenseVector2dField2d_s;
  typedef slip::DenseVector2dField2d<unsigned short> DenseVector2dField2d_us;
  typedef slip::DenseVector2dField2d<int> DenseVector2dField2d_i;
  typedef slip::DenseVector2dField2d<unsigned int> DenseVector2dField2d_ui;
  typedef slip::DenseVector2dField2d<char> DenseVector2dField2d_c;
  typedef slip::DenseVector2dField2d<unsigned char> DenseVector2dField2d_uc;

}//slip::

namespace slip{
 

template<typename T>
DenseVector2dField2d<T> operator+(const DenseVector2dField2d<T>& M1, 
                                                                  const T& val);

  template<typename T>
  DenseVector2dField2d<T> operator+(const T& val, 
                                                                    const DenseVector2dField2d<T>& M1);


template<typename T>
DenseVector2dField2d<T> operator-(const DenseVector2dField2d<T>& M1, 
                                                                  const T& val);

  template<typename T>
  DenseVector2dField2d<T> operator-(const T& val, 
                                                                    const DenseVector2dField2d<T>& M1);



  template<typename T>
  DenseVector2dField2d<T> operator*(const DenseVector2dField2d<T>& M1, 
                                                                    const T& val);

  template<typename T>
  DenseVector2dField2d<T> operator*(const T& val, 
                                                                    const DenseVector2dField2d<T>& M1);
  


  template<typename T>
  DenseVector2dField2d<T> operator/(const DenseVector2dField2d<T>& M1, 
                                                                    const T& val);


template<typename T>
DenseVector2dField2d<T> operator+(const DenseVector2dField2d<T>& M1, 
                                                                  const DenseVector2dField2d<T>& M2);

template<typename T>
DenseVector2dField2d<T> operator+(const DenseVector2dField2d<T>& M1, 
                                                                  const slip::Vector2d<T>& val);

  template<typename T>
  DenseVector2dField2d<T> operator+(const slip::Vector2d<T>& val, 
                                                                    const DenseVector2dField2d<T>& M1);
  

 
template<typename T>
DenseVector2dField2d<T> operator-(const DenseVector2dField2d<T>& M1, 
                                                                  const DenseVector2dField2d<T>& M2);

 
template<typename T>
DenseVector2dField2d<T> operator-(const DenseVector2dField2d<T>& M1, 
                                                                  const slip::Vector2d<T>& val);

  template<typename T>
  DenseVector2dField2d<T> operator-(const slip::Vector2d<T>& val, 
                                                                    const DenseVector2dField2d<T>& M1);
  

  
  template<typename T>
  DenseVector2dField2d<T> operator*(const DenseVector2dField2d<T>& M1, 
                                      const DenseVector2dField2d<T>& M2);


template<typename T>
DenseVector2dField2d<T> operator*(const DenseVector2dField2d<T>& M1, 
                                                                  const slip::Vector2d<T>& val);

  template<typename T>
  DenseVector2dField2d<T> operator*(const slip::Vector2d<T>& val, 
                                                                    const DenseVector2dField2d<T>& M1);
  


  template<typename T>
  DenseVector2dField2d<T> operator/(const DenseVector2dField2d<T>& M1, 
                                      const DenseVector2dField2d<T>& M2);

template<typename T>
DenseVector2dField2d<T> operator/(const DenseVector2dField2d<T>& M1, 
                                                                  const slip::Vector2d<T>& val);

 


}//slip::

namespace slip
{


//   template<typename T>
//   inline
//   DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator=(const DenseVector2dField2d<T> & rhs)
//   {
//     if(this != &rhs)
//       {
//              *matrix_ = *(rhs.matrix_);
//       }
//     return *this;
//   }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator=(const slip::Vector2d<T>& val)
  {
    std::fill_n(this->begin(),this->size(),val);
    return *this;
  }
  
 
 template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator=(const T& val)
  {
    std::fill_n(this->begin(),this->size(),val);
    return *this;
  }


  template<typename T>
  inline
  T& DenseVector2dField2d<T>::Vx1(const typename DenseVector2dField2d::size_type i,
                                                                  const typename DenseVector2dField2d::size_type j)
  {
    return (*this)[i][j][0];
  }

  template<typename T>
  inline
  const T& DenseVector2dField2d<T>::Vx1(const typename DenseVector2dField2d::size_type i,
                                                                                const typename DenseVector2dField2d::size_type j) const
  {
    return (*this)[i][j][0];
  }

 template<typename T>
  inline
  T& DenseVector2dField2d<T>::u(const typename DenseVector2dField2d::size_type i,
                                                                  const typename DenseVector2dField2d::size_type j)
  {
    return this->Vx1(i,j);
  }

  template<typename T>
  inline
  const T& DenseVector2dField2d<T>::u(const typename DenseVector2dField2d::size_type i,
                                                                                const typename DenseVector2dField2d::size_type j) const
  {
    return this->Vx1(i,j);
  }


   template<typename T>
  inline
  T& DenseVector2dField2d<T>::Vx2(const typename DenseVector2dField2d::size_type i,
                                                                  const typename DenseVector2dField2d::size_type j)
  {
    return (*this)[i][j][1];
  }

  template<typename T>
  inline
  const T& DenseVector2dField2d<T>::Vx2(const typename DenseVector2dField2d::size_type i,
                                                                                const typename DenseVector2dField2d::size_type j) const
  {
    return (*this)[i][j][1];
  }

 
  template<typename T>
  inline
  T& DenseVector2dField2d<T>::v(const typename DenseVector2dField2d::size_type i,
                                                                  const typename DenseVector2dField2d::size_type j)
  {
    return this->Vx2(i,j);
  }

  template<typename T>
  inline
  const T& DenseVector2dField2d<T>::v(const typename DenseVector2dField2d::size_type i,
                                                                                const typename DenseVector2dField2d::size_type j) const
  {
    return this->Vx2(i,j);
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator+=(const T&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::plus<slip::Vector2d<T> >(),val));
    return *this;
  }

  template<typename T>
  inline
  typename DenseVector2dField2d<T>::norm_type 
  DenseVector2dField2d<T>::norm(const typename DenseVector2dField2d::size_type i,
                                                                   const typename DenseVector2dField2d::size_type j) const
  {
    return (*this)[i][j].Euclidean_norm();
  }

   template<typename T>
  inline
  typename DenseVector2dField2d<T>::norm_type 
  DenseVector2dField2d<T>::angle(const typename DenseVector2dField2d::size_type i,
                                                                 const typename DenseVector2dField2d::size_type j) const
  {
    return (*this)[i][j].angle();
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator-=(const T&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::minus<slip::Vector2d<T> >(),val));
    return *this;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator*=(const T&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::multiplies<slip::Vector2d<T> >(),val));
    return *this;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator/=(const T&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::divides<slip::Vector2d<T> >(),val));
    return *this;
  }



 template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator+=(const slip::Vector2d<T>&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::plus<slip::Vector2d<T> >(),val));
    return *this;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator-=(const slip::Vector2d<T>&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::minus<slip::Vector2d<T> >(),val));
    return *this;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator*=(const slip::Vector2d<T>&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::multiplies<slip::Vector2d<T> >(),val));
    return *this;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator/=(const slip::Vector2d<T>&val)
  {
     std::transform(this->begin(),this->end(),this->begin(),std::bind2nd(std::divides<slip::Vector2d<T> >(),val));
    return *this;
  }
  

  
  template<typename T>
  inline
  DenseVector2dField2d<T> DenseVector2dField2d<T>::operator-() const
  {
    DenseVector2dField2d<T> tmp(*this);
    std::transform(this->begin(),this->end(),tmp.begin(),std::negate<slip::Vector2d<T> >());
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator+=(const DenseVector2dField2d<T>& rhs)
  {
    assert(this->dim1() == rhs.dim1());
    assert(this->dim2() == rhs.dim2()); 
    std::transform(this->begin(),this->end(),rhs.begin(),this->begin(),std::plus<typename DenseVector2dField2d::value_type>());
    return *this;
  }

   template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator-=(const DenseVector2dField2d<T>& rhs)
  {
    assert(this->dim1() == rhs.dim1());
    assert(this->dim2() == rhs.dim2());
    std::transform(this->begin(),this->end(),rhs.begin(),this->begin(),std::minus<typename DenseVector2dField2d::value_type>());
    return *this;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator*=(const DenseVector2dField2d<T>& rhs)
  {
    assert(this->dim1() == rhs.dim1());
    assert(this->dim2() == rhs.dim2());
    std::transform(this->begin(),this->end(),rhs.begin(),this->begin(),std::multiplies<typename DenseVector2dField2d::value_type>());
    return *this;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T>& DenseVector2dField2d<T>::operator/=(const DenseVector2dField2d<T>& rhs)
  {
    assert(this->dim1() == rhs.dim1());
    assert(this->dim2() == rhs.dim2());
    std::transform(this->begin(),this->end(),rhs.begin(),this->begin(),std::divides<typename DenseVector2dField2d::value_type>());
    return *this;
  }


  template<typename T>
  inline
  std::string 
  DenseVector2dField2d<T>::name() const {return "DenseVector2dField2d";} 
 

  template<typename T>
  inline
  void DenseVector2dField2d<T>::write_gnuplot(const std::string& file_path_name)
  {
    slip::write_gnuplot_vect2d<slip::DenseVector2dField2d<T> >(*this,file_path_name);
  }

  template<typename T>
  inline
  void DenseVector2dField2d<T>::write_tecplot(const std::string& file_path_name,
                                                                                      const std::string& title,
                                                                                      const std::string& zone)
  {
    slip::write_tecplot_vect2d<slip::DenseVector2dField2d<T> >(*this,
                                                                                                                       file_path_name,
                                                                                                                       title,
                                                                                                                       zone);
  }
 template<typename T>
  inline
  void DenseVector2dField2d<T>::read_gnuplot(const std::string& file_path_name)
  {
    slip::read_gnuplot_vect2d<slip::DenseVector2dField2d<T> >(file_path_name,*this);
  }

  template<typename T>
  inline
  void DenseVector2dField2d<T>::read_tecplot(const std::string& file_path_name)

  {
    slip::read_tecplot_vect2d<slip::DenseVector2dField2d<T> >(file_path_name,*this);
  }

  /* @{ */
  template<typename T>
  inline
  DenseVector2dField2d<T> operator+(const DenseVector2dField2d<T>& M1, 
                                                                    const T& val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp += val;
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator+(const T& val,
                                                                    const DenseVector2dField2d<T>& M1)
  {
    return M1 + val;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator-(const DenseVector2dField2d<T>& M1, 
                                                                    const T& val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp -= val;
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator-(const T& val,
                                                                    const DenseVector2dField2d<T>& M1)
  {
    return -(M1 - val);
  }


  template<typename T>
  inline
  DenseVector2dField2d<T> operator*(const DenseVector2dField2d<T>& M1, 
                                                                    const T& val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp *= val;
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator*(const T& val,
                                                                    const DenseVector2dField2d<T>& M1)
  {
    return M1 * val;
  }


  template<typename T>
  inline
  DenseVector2dField2d<T> operator/(const DenseVector2dField2d<T>& M1, 
                                                                    const T& val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp /= val;
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator+(const DenseVector2dField2d<T>& M1, 
                                                                    const DenseVector2dField2d<T>& M2)
  {
    assert(M1.dim1() == M2.dim1());
    assert(M1.dim2() == M2.dim2());
    DenseVector2dField2d<T> tmp(M1.dim1(),M1.dim2());
    std::transform(M1.begin(),M1.end(),M2.begin(),tmp.begin(),std::plus<slip::Vector2d<T> >());
    return tmp;
  }


  template<typename T>
  inline
  DenseVector2dField2d<T> operator+(const DenseVector2dField2d<T>& M1, 
                                                      const slip::Vector2d<T>&val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp += val;
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator+(const slip::Vector2d<T>&val,
                                                      const DenseVector2dField2d<T>& M1)
  {
    return M1 + val;
  }



 template<typename T>
  inline
  DenseVector2dField2d<T> operator-(const DenseVector2dField2d<T>& M1, 
                                                                    const DenseVector2dField2d<T>& M2)
  {
    assert(M1.dim1() == M2.dim1());
    assert(M1.dim2() == M2.dim2());
    DenseVector2dField2d<T> tmp(M1.dim1(),M1.dim2());
    std::transform(M1.begin(),M1.end(),M2.begin(),tmp.begin(),std::minus<slip::Vector2d<T> >());
    return tmp;
  }


  template<typename T>
  inline
  DenseVector2dField2d<T> operator-(const DenseVector2dField2d<T>& M1, 
                                                      const slip::Vector2d<T>&val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp -= val;
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator-(const slip::Vector2d<T>&val,
                                                      const DenseVector2dField2d<T>& M1)
  {
    return -(M1 - val);
  }


  template<typename T>
  inline
  DenseVector2dField2d<T> operator*(const DenseVector2dField2d<T>& M1, 
                                                                    const DenseVector2dField2d<T>& M2)
  {
    assert(M1.dim1() == M2.dim1());
    assert(M1.dim2() == M2.dim2());
    DenseVector2dField2d<T> tmp(M1.dim1(),M1.dim2());
    std::transform(M1.begin(),M1.end(),M2.begin(),tmp.begin(),std::multiplies<slip::Vector2d<T> >());
    return tmp;
  }


  template<typename T>
  inline
  DenseVector2dField2d<T> operator*(const DenseVector2dField2d<T>& M1, 
                                                                    const slip::Vector2d<T>&val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp *= val;
    return tmp;
  }

  template<typename T>
  inline
  DenseVector2dField2d<T> operator*(const slip::Vector2d<T>&val,
                                                                    const DenseVector2dField2d<T>& M1)
  {
    return M1 * val;
  }


 template<typename T>
  inline
  DenseVector2dField2d<T> operator/(const DenseVector2dField2d<T>& M1, 
                                                                    const DenseVector2dField2d<T>& M2)
  {
    assert(M1.dim1() == M2.dim1());
    assert(M1.dim2() == M2.dim2());
    DenseVector2dField2d<T> tmp(M1.dim1(),M1.dim2());
    std::transform(M1.begin(),M1.end(),M2.begin(),tmp.begin(),std::divides<slip::Vector2d<T> >());
    return tmp;
  }


  template<typename T>
  inline
  DenseVector2dField2d<T> operator/(const DenseVector2dField2d<T>& M1, 
                                                                    const slip::Vector2d<T>&val)
  {
    DenseVector2dField2d<T>  tmp = M1;
    tmp /= val;
    return tmp;
  }

/* @} */

}//slip::

#endif //SLIP_DENSEVECTOR2DFIELD2D_HPP