synaps/linalg/VectDse.h

Go to the documentation of this file.

/*********************************************************************
*   This file is part of the source code of the SYNAPS environment.  *
*   Author(s): B. Mourrain, GALAAD, INRIA                            *
**********************************************************************/
#ifndef synaps_linalg_VectDse_H
#define synaps_linalg_VectDse_H
//----------------------------------------------------------------------
#include <sstream>
#include <vector>
#include <synaps/init.h>
#include <synaps/base/shared_object.h>
#include <synaps/arithm/texp.h>
#include <synaps/base/Range.h>
#include <synaps/base/seq1d.h>
#include <synaps/arithm/let.h>
#include <synaps/arithm/REF.h>
#include <synaps/linalg/rep1d.h>
#include <synaps/linalg/VECTOR.m>
#include <synaps/linalg/MATRIX.m>
//----------------------------------------------------------------------
__BEGIN_NAMESPACE_SYNAPS
//----------------------------------------------------------------------
template <class C, class R = linalg::rep1d<C> >
struct VectDse   {

  shared_object<R>  data;
  R &       rep()       {return (*data);}
  const R & rep() const {return (*data);}

  typedef typename R::size_type          size_type;
  typedef C                              value_type;
  typedef typename R::iterator           iterator;
  typedef typename R::const_iterator     const_iterator;
  typedef typename R::reverse_iterator   reverse_iterator;
  typedef typename R::const_reverse_iterator const_reverse_iterator;
  typedef R                              rep_type;
  typedef VectDse<C,R>                   self_t;
  
  operator arithm::template_expression<self_t>(){ return arithm::template_expression<self_t>(*this); };
  
  VectDse(): data() {}
  VectDse(const R & r):data(r){}
  VectDse(R* r):data(r){}
  VectDse(const shared_object<R> & r):data(r){}
  VectDse(size_type s) {rep().resize(s);}
  VectDse(size_type s,const type::AsSize & mth) {rep().resize(s);}
  VectDse(size_type s, value_type* t) {rep()=R(t,t+s);}
  VectDse(iterator b,iterator e) {rep()=R(b,e);}
  VectDse(size_type m, const char * str);

  //  inline operator generic () { return generic(data.body); }

  // Copy constructor
  VectDse(const self_t & r): data(r.data) {}
  template<class X>   VectDse( const arithm::template_expression<X> & x );
  template<class C1,class R1> VectDse( const VectDse<C1,R1> & v );
  template<class X>   self_t & operator=(const X & exp);
  
  iterator               begin()        {return rep().begin();}
  iterator               end()          {return rep().end();}
  const_iterator         begin()  const {return rep().begin();}
  const_iterator         end()    const {return rep().end();}
  reverse_iterator       rbegin()       {return rep().rbegin();}
  reverse_iterator       rend()         {return rep().rend();}
  const_reverse_iterator rbegin() const {return rep().rbegin();}
  const_reverse_iterator rend  () const {return rep().rbegin();}

  self_t & operator =(const self_t & V) { using let::assign; assign(data,V.data); return *this; }
  self_t & operator+=(const self_t & v)     { add(*this,v); return *this; };
  self_t & operator-=(const self_t & v)     { sub(*this,v); return *this; };
  self_t & operator*=(const value_type & c) { mul(*this,c); return *this; };
  self_t & operator/=(const value_type & c) { div(*this,c); return *this; };
  template<class X> 
  self_t & operator+=(const X & x) { *this = *this + x; return *this; };
  template<class X>
  self_t & operator-=(const X & x) { *this = *this - x; return *this; };
  template<class X>
  self_t & operator*=(const X & x) { *this = *this * x; return *this; };
  template<class X>
  self_t & operator/=(const X & x) { *this = *this / x; return *this; };
  
  value_type & operator[] (size_type i)       {return rep()[i];}
  const value_type & operator[] (size_type i) const {return rep()[i];}

  bool operator== (const self_t & v) 
  { 
    using VECTOR::equal;
    return equal(this->rep(),v.rep());
  };
  
  template<class S>
  self_t & operator=(const VAL<S> & V) 
  { 
    using let::assign;
    //using let::copy; 
    assign(*this,V.rep()); return *this;
  }
  
  bool operator!= (const self_t & v)  { return !(*this==v); }

  VectDse<C,seq1d<value_type,iterator> > operator[](const Range & I);
  
  size_type size() const {return rep().size();}
  self_t resize(int n)  {rep().resize(n); return *this;}

};

//======================================================================
/*{
 The @value_type@ is the type of the coefficients. The @size_type@ is the
 type of the indices used to access the elements of the vectors.  The
 @iterator@ is the type of the forward iterators on the vectors.  The
 internal representation type @R@ is @rep_type@.
}*/

namespace let
{
  template<class Ca,class Ra,class Cb, class Rb> inline
  void assign( VectDse<Ca,Ra> & a, const VectDse<Cb,Rb> & b )
  {
    using VECTOR::assign;
    assign(a.rep(),b.rep());
  };
  //----------------------------------------------------------------------
  template<class Ca, class Ra,class Cb, class  Rb> inline
  void copy(VectDse<Ca,Ra> & v, VectDse<Cb,Rb>* t)
  {
    using VECTOR::assign;
    assign(v.rep(), (*t).rep());
  }
  
  template<class Ca, class Ra, class Cb, class Rb> inline
  void copy( VectDse<Ca,Ra> & a, const VectDse<Cb,Rb> & b )
  {
    using VECTOR::assign; assign(a.rep(),b.rep());
  };

  template<class Ca, class Ra, class TEXP> inline
  void copy( VectDse<Ca,Ra> & a, const TEXP & b )
  {
    using let::assign; assign(a,b);
  }
}


//======================================================================
// INPUT OUTPUT
//======================================================================
template<class C, class R>
std::ostream & operator<<(std::ostream & os, const VectDse<C,R> & V)
{
  using namespace VECTOR;
  return print(os,V.rep());
}
//----------------------------------------------------------------------
template<class C, class R>  inline
std::istream & operator>>(std::istream & is, VectDse<C,R> & V)
{
  using namespace VECTOR;
  return read(is,V);
}

//--------------------------------------------------------------------
template <class C, class R> 
VectDse<C,R>::VectDse(size_type m, const char * str)
{
  rep()= R(m);
  std::istringstream ins(str);
  for(size_type i=0; i< m && ins.good(); i++)   ins >> rep()[i];
}

//----------------------------------------------------------------------
template <class C, class R> 
VectDse<C,seq1d<typename VectDse<C,R>::value_type,typename VectDse<C,R>::iterator> > 
   VectDse<C,R>::operator[](const Range & I) 
{
  using namespace VECTOR;
  typedef typename VectDse<C,R>::value_type value_type;
  typedef typename VectDse<C,R>::iterator   iterator;
  typedef seq1d<value_type,iterator>     rep_type;
      
  typename VectDse<C,R>::size_type n=I.i2-I.i1+1;
  return VectDse<C,rep_type>(rep_type(n,create(iterator(),rep(),I.i1),
                                      create(iterator(),rep(),I.i2+1)));
}
//----------------------------------------------------------------------

template<class C, class R> 
void add( VectDse<C,R> & r, const VectDse<C,R> & a )
{
  using VECTOR::add;
  add(r.rep(),a.rep());
};

template<class C, class R>
void add( VectDse<C,R> & r, const VectDse<C,R> & a, const VectDse<C,R> & b )
{
  using VECTOR::add;
  add(r.rep(),a.rep(),b.rep());
};

template<class C, class R>
void sub( VectDse<C,R> & r, const VectDse<C,R> & a )
{
  using VECTOR::sub;
  sub(r.rep(),a.rep());
};

template<class C, class R> 
void sub( VectDse<C,R> & r, const VectDse<C,R> & a, const VectDse<C,R> & b )
{
  using VECTOR::sub;
  sub(r.rep(),a.rep(),b.rep());
};

template<class C, class R> 
void mul( C & r, const VectDse<C,R> & a, const VectDse<C,R> & b )
{
  using VECTOR::innerprod;
  r = innerprod(a.rep(),b.rep());
};

template<class C, class R> 
void mul( VectDse<C,R> & r, const VectDse<C,R> & a, const C & b )
{
  using VECTOR::mul_ext;
  mul_ext(r.rep(),a.rep(),b);
};

template<class C, class R>
void mul( VectDse<C,R> & r, const C & b, const VectDse<C,R> & a )
{
  using VECTOR::mul_ext;
  mul_ext(r.rep(),a.rep(),b);
};

template<class C, class R>
void mul( VectDse<C,R> & r, const C & c )
{
  using VECTOR::mul_ext;
  mul_ext(r.rep(),c);
};

template<class C, class R>
void div( VectDse<C,R> & r, const C & c )
{
  using VECTOR::div_ext;
  div_ext(r.rep(),c);
};

template<class C, class R>
void div( VectDse<C,R> & r, const VectDse<C,R> & a, const C & c )
{
  using VECTOR::div_ext;
  div_ext(r.rep(),a.rep(),c);
};

template<class C, class R> 
void neg( VectDse<C,R> & r, const VectDse<C,R> & a )
{
  using VECTOR::neg;
  neg( r.rep(), a.rep()) ;
};

namespace arithm 
{
  template<class C, class R> struct structureof< VectDse<C,R> > { typedef structure::vector T; };
  
#define template_arg_list template<class Ca, class Ra, class Cb, class Rb>
#define parm0 VectDse<Ca,Ra>
#define parm1 VectDse<Cb,Rb>
  declare_binary_operator(template_arg_list,parm0,parm1,_add_,operator+);
  declare_binary_operator(template_arg_list,parm0,parm1,_sub_,operator-);
  declare_binary_operator(template_arg_list,parm0,parm1,_mul_,operator*);
#undef  template_arg_list
#define template_arg_list template<class C, class R>
#undef  parm0
#define parm0 VectDse<C,R>
#undef  parm1
#define parm1 typename VectDse<C,R>::value_type
  declare_binary_operator(template_arg_list,parm0,parm1,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm1,parm0,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm0,parm1,_div_,operator/);
#undef  parm1
#undef  template_arg_list
#define template_arg_list template<class C, class R, class K>
#define parm1 typename K::integer
  declare_binary_operator(template_arg_list,parm0,parm1,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm1,parm0,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm0,parm1,_div_,operator/);
#undef  parm1
#define parm1 typename K::rational
  declare_binary_operator(template_arg_list,parm0,parm1,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm1,parm0,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm0,parm1,_div_,operator/);
#undef parm1
#define parm1 typename K::floating
  declare_binary_operator(template_arg_list,parm0,parm1,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm1,parm0,_mul_,operator*);
  declare_binary_operator(template_arg_list,parm0,parm1,_div_,operator/);
#undef template_arg_list
#undef parm0
#undef parm1
#define template_arg_list template<class C, class R>
#define parm0 VectDse<C,R>
  declare_unary_operator(template_arg_list,parm0,_neg_,operator-);
#undef template_arg_list
#undef parm0
};


template<class C, class R>
template<class X> inline
VectDse<C,R>::VectDse( const arithm::template_expression<X> & x )
{
  using let::assign;
  assign(*this,x);
};

template<class C, class R>
template<class X> inline
VectDse<C,R> & VectDse<C,R>::operator= ( const X & x )
{
  using let::assign;  assign(*this,x);
};

template<class C, class R> 
template<class C1,class R1>  VectDse<C,R>::VectDse( const VectDse<C1,R1> & v )
{
  let::assign(*this,v);
};


namespace type
{
  template<class Y, class C> struct instanceof< VectDse<C>, Y  >  { typedef VectDse<Y> T; };
};
__END_NAMESPACE_SYNAPS
//====================================================================
#endif //synaps_linalg_VectDse_H