//  NumVec Template         copyright 1996, Information Disciplines, Inc.
//  ---------------

#ifndef NUMVEC
#define NUMVEC const NumVec

// This class supports one dimensional arrays with numeric elements and fixed
// bounds.  It's purpose is not to compete with STL's Vector class, but rather 
// to serve as a base class for Matrix and other higher-dimensionality classes
// and as a more understandable container example for teaching.

// Arithmetic operators on an element-by-element basis are supported by 
// inline functions, in order to avoid instantiating operations that are 
// unsupported by the element type. 

template<class T>
  class NumVec {
protected:
  sizet  low_bound, extent;             //  Dimension information
  T     *data;                          //  Internal representation


//  Private methods used by constructors and assignment operator
//  ------------------------------------------------------------
  void allocate ()   {data = new T [extent];  assert(data);}
  void copy  (NUMVEC& rs);      //  Copy element by element (Bounds
                                //   needn't be identical, just extent.
  void copy  (T    rs)   ;      //  Set all elements to rs


public:

//  Constructors and destructor
//  ---------------------------
  NumVec(SIZET lb, SIZET hb)            //  Specified lower and upper bound
     : low_bound(lb), extent(hb-lb+1) 
    {assert(lb <= hb); allocate();}

  explicit NumVec(SIZET dim)            //  0-origin by default, like C
     : low_bound(0) , extent(dim)       //  (explicit to avoid unintended
    {assert(dim > 0); allocate();}      //   allocation from scalar operand.) 

  NumVec(NUMVEC<T>& x)                  //  Copy constructor
     : low_bound(x.low_bound), extent(x.extent)
    {allocate(); copy(x);}

  ~NumVec() {delete [] data;}       //  Not virtual -- clients won't delete 
                                    //      derived objects through a pointer

//  Accessors
//  ---------
  sizet  lbound() const {return low_bound;}
  sizet  hbound() const {return low_bound+extent-1;}
  sizet  dim   () const {return extent;}


//  Subscript operator
//  ------------------
   T&    operator[](INT rs) 
   {assert(rs>=low_bound && rs <= low_bound+extent); return data[rs-low_bound];}
   T     operator[](INT rs) const
   {assert(rs>=low_bound && rs <= low_bound+extent); return data[rs-low_bound];}


//  Member operators and other functions
//  ------------------------------------
    NumVec& operator=(NUMVEC& rs)           //  Array assignment
    {copy(rs); return *this;}               //    assignment (deep copy)
                                            //  (Note that bounds needn't
                                            //     be identical, just extent)

    NumVec& operator=(const T  rs)           //  Set all elements to a
    {copy(rs); return *this;}                //    scalar value

    NumVec&  reverse()                       //  Note:  in-place reverse
    {for (sizet ctr = 0; ctr < extent/2; ctr++)
      {swap(*(data+ctr), *(data+extent-ctr-1));}
     return *this;
    }

    ostream& put(ostream& s) const;          //  Output-stream operator 


//  Arithmetic operators
//  --------------------

   NumVec  operator-()                         //  Unary minus operator
   {NumVec result(lbound(),hbound()); 
    for (sizet ctr=0;  ctr < extent; ++ctr)
       *(result.data + ctr) = *(data + ctr);
    return result; 
   }

  NumVec& operator += (const NumVec<T>& rs)    //  Compound add
     {assert(dim()==rs.dim()); 
      for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) += *(rs.data+ctr);           
	  return *this;}	                                                       
  NumVec& operator += (const T rs) 
     {for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) += rs;           
	  return *this;}	                                                       
  
  NumVec& operator -= (const NumVec<T>& rs)    //  Compound subtract
     {assert(dim()==rs.dim()); 
      for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) -= *(rs.data+ctr);           
	  return *this;}	                                                       
  NumVec& operator -= (const T rs) 
     {for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) -= rs;           
	  return *this;}	                                                       
 
NumVec& operator *= (const NumVec<T>& rs)       //  Compound multiply
     {assert(dim()==rs.dim()); 
      for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) *= *(rs.data+ctr);           
	  return *this;}	                                                       
NumVec& operator *= (const T rs) 
     {for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) *= rs;           
	  return *this;}	                                                       

NumVec& operator /= (const NumVec<T>& rs)       //  Compound divide
     {assert(dim()==rs.dim()); 
      for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) /= *(rs.data+ctr);           
	  return *this;}	                                                       
NumVec& operator /= (const T rs) 
     {for (sizet ctr = 0; ctr < extent; ++ctr) 
           *(data+ctr) /= rs;           
	  return *this;}	                                                       


//  Reduction functions 
//  -------------------

   T sum() const
   {T  result = 0;                          //  (0 must be convertible to T)
    for (sizet ctr = 0; ctr < extent; ctr++)
       result += *(data+ctr);
    return result;
   }

  T prod () const                           
  {T  result = 1;                           //  (1 must be convertible to T)
   for (sizet ctr = 0; ctr < dim(); ctr++)
      result *= *(data+ctr);
   return result;
  }

};                  // *************  End of class definition
//  The following functions should be moved to a separate implementation file.

template<class T>
  void NumVec<T>::copy(const NumVec<T>& rs) 
  {assert(dim() == rs.dim());
   for (sizet  count = 0; count < dim(); count++)
      *(data+count) = *(rs.data+count);
   return;
  }

template<class T>                   //  
  void NumVec<T>::copy(T rs) 
  {for (sizet count  = 0; count < dim(); count++)
     *(data+count) = rs;
   return;
  }

template<class T>                   //  
  ostream& NumVec<T>::put(ostream& s) const
  {sizet  ctr = lbound();
   s << data[ctr++];
   
   while (ctr <= hbound())
     s << ", " << data[ctr++];
   return s;
  }


//  Non-member functions
//  --------------------

//  Arithmetic operations

template<class T> 
   NumVec<T> operator+(const NumVec<T> ls, const NumVec<T> rs) 
      {return NumVec<T>(ls)+= rs;}
template<class T> 
   NumVec<T> operator+(const NumVec<T> ls, const T rs) 
      {return NumVec<T>(ls)+= rs;}
template<class T> 
   NumVec<T> operator+(const T ls, const NumVec<T> rs) {return rs + ls;}

template<class T> 
   NumVec<T> operator-(const NumVec<T> ls, const NumVec<T> rs) 
      {return NumVec<T>(ls)-= rs;}
template<class T> 
   NumVec<T> operator-(const NumVec<T> ls, const T rs) 
      {return NumVec<T>(ls)-= rs;}

template<class T> 
   NumVec<T> operator*(const NumVec<T> ls, const NumVec<T> rs)
      {return NumVec<T>(ls)*= rs;}
template<class T>
   NumVec<T> operator*(const NumVec<T> ls, const T rs)
      {return NumVec<T>(ls)*=rs;}
template<class T>
   NumVec<T> operator*(const T ls, const NumVec<T> rs) {return rs * ls;}


template<class T> 
   NumVec<T> operator/(const NumVec<T> ls, const NumVec<T> rs) 
      {return NumVec<T>(ls)/= rs;}
template<class T> 
   NumVec<T> operator/(const NumVec<T> ls, const T rs) 
      {return NumVec<T>(ls)/= rs;}




template<class T>       
  ostream& operator<< (ostream& ls, NUMVEC<T>& rs){return rs.put(ls);}

#endif