compressed1d.h

MTL C++ Numeric Library

// -*- c++ -*-
//
// Copyright 1997, 1998, 1999 University of Notre Dame.
// Authors: Andrew Lumsdaine, Jeremy G. Siek, Lie-Quan Lee
//
// This file is part of the Matrix Template Library
//
// You should have received a copy of the License Agreement for the
// Matrix Template Library along with the software;  see the
// file LICENSE.  If not, contact Office of Research, University of Notre
// Dame, Notre Dame, IN  46556.
//
// Permission to modify the code and to distribute modified code is
// granted, provided the text of this NOTICE is retained, a notice that
// the code was modified is included with the above COPYRIGHT NOTICE and
// with the COPYRIGHT NOTICE in the LICENSE file, and that the LICENSE
// file is distributed with the modified code.
//
// LICENSOR MAKES NO REPRESENTATIONS OR WARRANTIES, EXPRESS OR IMPLIED.
// By way of example, but not limitation, Licensor MAKES NO
// REPRESENTATIONS OR WARRANTIES OF MERCHANTABILITY OR FITNESS FOR ANY
// PARTICULAR PURPOSE OR THAT THE USE OF THE LICENSED SOFTWARE COMPONENTS
// OR DOCUMENTATION WILL NOT INFRINGE ANY PATENTS, COPYRIGHTS, TRADEMARKS
// OR OTHER RIGHTS.
//
//
//===========================================================================


#ifndef MTL_COMPRESSED1D_H
#define MTL_COMPRESSED1D_H


#include <vector>
#include <algorithm>

#include "mtl/mtl_config.h"
#include "mtl/entry.h"
#include "mtl/compressed_iter.h"
#include "mtl/reverse_iter.h"
#include "mtl/matrix_traits.h"
#include "mtl/refcnt_ptr.h"
#include "mtl/scaled1D.h"
#include "mtl/light1D.h"
#include "mtl/dense1D.h" /* for VC++ workaround -JGS */

namespace mtl {

using std::lower_bound;
using std::find;

//: Compressed Sparse Vector
//
// The compressed1D Vector is a sparse vector implemented with
// a pair of parallel arrays. One array is for the element values,
// and the other array is for the indices of those elements.
// The elements are ordered by their index as they are inserted
// into the compressed1D. compressed1D's can be used to build
// matrices with the array storage format, and they can also
// be used on their own.
//
// <pre>
// [ (1.2, 3), (4.6, 5), (1.0, 10), (3.7, 32) ] A Sparse Vector
//
// [ 1.2, 4.6, 1.0, 3.7 ]  Value Array
// [  3,   5,  10,  32 ]   Index Array
// </pre>
//
// <p>The <tt>compressed1D::iterator</tt> dereferences (<tt>*i</tt>) to
// return the element value. One can access the index of that element
// with the <tt>i.index()</tt> function. One can also access
// the array of indices through the <tt>nz_struct()</tt> method.
//
// <p>One particularly useful fact is that one can perform
// scatters and gathers of sparse elements by using the
// <tt>mtl::copy(x,y)</tt> function with a sparse and a
// dense vector.
//
//!category: containers
//!component: type
//!models: Vector
//!definition: compressed1D.h
//!tparam: T - the element type
//!tparam: SizeType - the type for the stored indices - int
//!tparam: IND_OFFSET - To handle indexing from 0 or 1 - index_from_zero
//!example: gather_scatter.cc, array2D.cc, sparse_copy.cc
template <class T, class SizeType = int, int IND_OFFSET = index_from_zero>
class compressed1D {

  typedef std::vector<T> values_vec;
  typedef std::vector<SizeType> indices_vec;
  typedef indices_vec indices_t;
  typedef values_vec  values_t;

  typedef refcnt_ptr< values_vec > values_ptr;
  typedef refcnt_ptr< indices_vec > indices_ptr;

  typedef compressed1D<T,SizeType,IND_OFFSET> self;

  typedef typename indices_vec::iterator index_iterator;
  typedef typename values_vec::iterator value_iterator;
  typedef typename indices_vec::const_iterator index_const_iterator;
  typedef typename values_vec::const_iterator value_const_iterator;

  friend class elt_ref<self>;
  friend class const_elt_ref<self>;

public:
  /**@name Type Definitions */

  enum { N = 0 };

  //: This is a sparse vector  
  typedef sparse_tag sparsity;
  //: This is a 1D container
  typedef oned_tag dimension;
  //: Scaled type of this vector
  typedef scaled1D< self > scaled_type;
  //: Element type
  typedef typename values_t::value_type value_type;
#if defined(_MSVCPP_)
  typedef value_type* pointer;
#else
  //: A pointer to the element type
  typedef typename values_t::pointer pointer;
#endif
  //: Unsigned integral type for dimensions and indices
  typedef SizeType size_type;
  //: Integral type for differences in iterators
  typedef typename values_t::difference_type difference_type;
  //: Reference to the value type
  typedef elt_ref<self> reference;
  //: Const reference to the value type
  typedef const_elt_ref<self> const_reference;
  //: Iterator type
  typedef compressed_iter<0,values_vec, indices_vec,IND_OFFSET> iterator;
  //: Const iterator type
  typedef compressed_iter<1,values_vec, indices_vec,IND_OFFSET> const_iterator;
  //: Reverse iterator type
  typedef reverse_iter<iterator> reverse_iterator;
  //: The const reverse iterator type
  typedef reverse_iter<const_iterator> const_reverse_iterator;
  //: Reference to the index array
  typedef const indices_vec& IndexArrayRef; 
  //: The type for the index array
  typedef const indices_t IndexArray; 

  //: The type for the subrange vector
  typedef self subrange_type; /* JGS need to think about this */

  

  /**@name Constructors */
  
  //: Default Constructor
  inline compressed1D()
    : values(new values_t(0)), 
      indices(new indices_t(0)), 
      size_(0) { }
      
  //: Length N Constructor
  inline compressed1D(size_type n)
    : values(new values_t(0)), 
      indices(new indices_t(0)), 
      size_(n) { }

  //: Copy Constructor
  inline compressed1D(const self& x)
    : values(x.values), indices(x.indices), size_(x.size_) { }

  //: Index Array Constructor
  template <class IndexIter>
  inline compressed1D(IndexIter first, IndexIter last,  size_type n)
    : values(new values_t(n)),
      indices(new indices_t(n)),
      size_(n)
  { 
    std::copy(first, last, indices->begin());
  }

  //: Assignment Operator
  inline self& operator=(const self& x) {
    values = x.values;
    indices = x.indices;
    size_ = x.size_;
    return *this;
  }

  /**@name Access Methods */
  

  /**@name Iterator Access Methods */
  
  //: Return an iterator pointing to the beginning of the vector
  //!wheredef: Container
  inline iterator begin() {
    return iterator(values->begin(), indices->begin(), 0);
  }
  //: Return an iterator pointing past the end of the vector
  //!wheredef: Container
  inline iterator end() {
    return iterator(values->begin(), indices->begin(), indices->size());
  }
  //: Return a const iterator pointing to the begining of the vector
  //!wheredef: Container
  inline const_iterator begin() const {
    return const_iterator(values->begin(), indices->begin(), 0); 
  }
  //: Return a const iterator pointing past the end of the vector
  //!wheredef: Container
  inline const_iterator end() const { 
    return const_iterator(values->begin(), indices->begin(), indices->size());
  }
  //: Return a reverse iterator pointing to the last element of the vector
  //!wheredef: Reversible Container
  inline reverse_iterator rbegin() {
    return reverse_iterator(end());
  }
  //: Return a reverse iterator pointing past the end of the vector
  //!wheredef: Reversible Container
  inline reverse_iterator rend() {
    return reverse_iterator(begin());
  }
  //: Return a const reverse iterator pointing to the last element of the vector
  //!wheredef: Reversible Container
  inline const_reverse_iterator rbegin() const {
    return const_reverse_iterator(end());
  }
  //: Return a const reverse iterator pointing past the end of the vector
  //!wheredef: Reversible Container
  inline const_reverse_iterator rend() const { 
    return const_reverse_iterator(begin());
  }
  

  /**@name Element Access Methods */
  
  //: Access the element with index i
  inline reference operator[](size_type i) MTL_THROW_ASSERTION {
    MTL_ASSERT(i < size(), "compressed1D::operator[]");
    return reference(*this, i);
  }
  //: Access the element with index i
  inline const_reference operator[](size_type i) const MTL_THROW_ASSERTION {
    MTL_ASSERT(i < size(), "compressed1D::operator[]");
    return const_reference(*this, i);
  }
  
  //for BEAM only
  value_type& get_ref(size_type i) {
    iterator iter = find(i);
    if ( iter != end() ) { 
      if ( iter.index() != i )
	iter = insert(iter, i, T(0));
    } else 
      iter = insert(iter, i, T(0));

    return *iter;
  }
  

  /*
   * list operations
   */

  /* OneD required insert method: */
  //: Insert val into the vector at index i
  inline iterator insert(size_type i, const T& val) MTL_THROW_ASSERTION {
    MTL_ASSERT(i < size(), "compressed1D::insert");
    index_iterator index_iter = lower_bound(indices->begin(),
                                            indices->end(), 
                                            i - IND_OFFSET);/* F to C */
    index_iter = indices->insert(index_iter, i - IND_OFFSET); /* F to C */
    size_type n = index_iter - indices->begin();
    value_iterator val_iter = values->insert(values->begin() + n, val);
    return iterator(index_iter, val_iter);
  }

  //: Push back, warning: must insert elements ordered by index.
  // This method does not change the size() of the vector.
  inline void push_back(size_type i,  const T& val) MTL_THROW_ASSERTION {
    MTL_ASSERT(i < size(), "compressed1D::push_back");
    values->push_back(val);
    indices->push_back(i);
  }
  
  //: Erase the vector
  inline void clear() { indices->clear(); values->clear(); }
  //: The size of the vector (including non-zeroes)
  inline size_type size() const { return size_; }
  //: The number of non-zero elements (the number stored)
  inline size_type nnz() const { return values->size(); }
  //: Resize the vector to size n
  inline void resize(size_type n) { 
    size_ = n;
  }
  //: Reserve storage for n non-zero elements
  inline void reserve(size_type n) {
    values->reserve(n);
    indices->reserve(n);
  }
  //: Returns the array of indices
  inline IndexArrayRef nz_struct() const {
    return *indices;
  }

  //: Returns the array of indices
  inline IndexArrayRef nz_struct() {
    return *indices;
  }

protected:

  inline iterator find(size_type i) {
    index_iterator iter = lower_bound(indices->begin(), indices->end(), 
                                      i - IND_OFFSET); /* F to C */
    size_type n = iter - indices->begin();
    return iterator(values->begin(), indices->begin(), n);
  }

  inline const_iterator find(size_type i) const {
    index_const_iterator iter = lower_bound(indices->begin(), indices->end(),
                                    i - IND_OFFSET); /* F to C */
    size_type n = iter - indices->begin();
    return const_iterator(values->begin(), indices->begin(), n);
  }

  inline iterator insert(iterator iter, size_type i, T v) {
    index_iterator ind = indices->insert(iter.index_iter(), 
                                 i - IND_OFFSET); /* F to C */
    values->insert(iter.value_iter(), v);
    size_type n = ind - indices->begin();
    return iterator(values->begin(), indices->begin(), n);
  }

  
protected:
  values_ptr values;
  indices_ptr indices;
  size_type size_;
};

} /* namespace mtl */

#endif




#if 0

  class elt_inserter {
  public:
    inline elt_inserter(const elt_inserter& x) : val(x.val), ind(x.ind) { }
    inline elt_inserter(T& v, size_type& i) : val(v), ind(i) { }
    void operator=(const twod_elt<T,size_type>& elt) {
      val = elt.val();
      ind = elt.minor();
    }
  protected:
    T& val;
    size_type& ind;
  };

  class insert_iterator {
    typedef values_vec::iterator val_iter_t;
    typedef indices_vec::iterator ind_iter_t;
    typedef insert_iterator self;
  public:
    inline insert_iterator(const val_iter_t& v, const ind_iter_t& i)
      : val_iter(v), ind_iter(i) { }
    inline elt_inserter operator * () { 
      return elt_inserter(*val_iter, *ind_iter); 
    }
    inline self& operator ++ () { ++val_iter; ++ind_iter; return *this; }
    inline self& operator ++ (int) { 
      self t = *this; ++val_iter; ++ind_iter; return t;
    }
  protected:
    val_iter_t val_iter;
    ind_iter_t ind_iter;
  };

  inline insert_iterator inserter() { 
    return insert_iterator(values->begin(), indices->begin());
  }
#endif