AbstractLinAlgPack_SparseVectorClassDecl.hpp

// @HEADER
// ***********************************************************************
// 
// Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
//                  Copyright (2003) Sandia Corporation
// 
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
// license for use of this work by or on behalf of the U.S. Government.
// 
// This library is free software; you can redistribute it and/or modify
// it under the terms of the GNU Lesser General Public License as
// published by the Free Software Foundation; either version 2.1 of the
// License, or (at your option) any later version.
//  
// This library is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
//  
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
// USA
// Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
// 
// ***********************************************************************
// @HEADER

#ifndef SPARSE_VECTOR_CLASS_DECL_H
#define SPARSE_VECTOR_CLASS_DECL_H

#include <assert.h>

#include <vector>
#include <sstream>

#include "AbstractLinAlgPack_SpVecIndexLookupClass.hpp"

namespace AbstractLinAlgPack {

namespace SparseVectorUtilityPack {

void assert_is_sorted(bool is_sorted);

class DoesNotExistException : public std::logic_error
{public: DoesNotExistException(const std::string& what_arg) : std::logic_error(what_arg) {}};
class NotSortedException : public std::logic_error
{public: NotSortedException(const std::string& what_arg) : std::logic_error(what_arg) {}};
class DuplicateIndexesException : public std::logic_error
{public: DuplicateIndexesException(const std::string& what_arg) : std::logic_error(what_arg) {}};
class OutOfRoomException : public std::logic_error
{public: OutOfRoomException(const std::string& what_arg) : std::logic_error(what_arg) {}};
class UnsizedException : public std::logic_error
{public: UnsizedException(const std::string& what_arg) : std::logic_error(what_arg) {}};
class NoNonZeroElementsException : public std::logic_error
{public: NoNonZeroElementsException(const std::string& what_arg) : std::logic_error(what_arg) {}};

} // end namespace SparseVectorUtilityPack

// /////////////////////////////////////////////////////////////////////////////////////

template<class T_Element>
SparseVectorSlice<T_Element> create_slice(
  const SparseVectorUtilityPack::SpVecIndexLookup<T_Element>& index_lookup
  , size_type size, Range1D rng);


template <class T_Element>
class SparseVectorSlice;

// ///////////////////////////////////////////////////////////////////////
template <class T_Element, class T_Alloc = std::allocator<T_Element> >
class SparseVector {
public:

  typedef T_Alloc                     allocator_type;
  typedef T_Element                   element_type;
  typedef AbstractLinAlgPack::size_type         size_type;
  typedef ptrdiff_t                   difference_type;
  typedef element_type*                 iterator;
  typedef const element_type*               const_iterator;

#if 0 /* defined(_WINDOWS) || defined(_INTEL_CXX) */

  typedef std::reverse_iterator<iterator, element_type
    , element_type&, element_type*, difference_type>  reverse_iterator;

  typedef std::reverse_iterator<const_iterator
    , element_type, const element_type&
    , const element_type*, difference_type>       const_reverse_iterator;

#else

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

#endif

  typedef SparseVectorUtilityPack::DoesNotExistException    DoesNotExistException;
  typedef SparseVectorUtilityPack::NotSortedException     NotSortedException;
  typedef SparseVectorUtilityPack::DuplicateIndexesException  DuplicateIndexesException;
  typedef SparseVectorUtilityPack::OutOfRoomException     OutOfRoomException;
  typedef SparseVectorUtilityPack::UnsizedException     UnsizedException;
  typedef SparseVectorUtilityPack::NoNonZeroElementsException NoNonZeroElementsException;



  SparseVector(const allocator_type& alloc = allocator_type());

  SparseVector(bool assume_sorted, const allocator_type& alloc = allocator_type());

  SparseVector(size_type size, size_type max_nz, difference_type offset = 0
    , bool assume_sorted = false, const allocator_type& alloc = allocator_type());

  SparseVector(const SparseVector<T_Element,T_Alloc>& sp_vec);

  SparseVector( SparseVectorSlice<T_Element> sp_vec_slc
    , const allocator_type& alloc = allocator_type());

  ~SparseVector();


  SparseVector<T_Element,T_Alloc>& operator=(const SparseVector<T_Element,T_Alloc>& sp_vec);

  SparseVector<T_Element,T_Alloc>& operator=(const SparseVectorSlice<T_Element>& sp_vec_slc);


  EOverLap overlap(const SparseVectorSlice<T_Element>& sv) const;


  size_type dim() const;

  size_type nz() const;

  difference_type offset() const;

  bool is_sorted() const;


  iterator begin();

  const_iterator begin() const;

  iterator end();

  const_iterator end() const;

  reverse_iterator rbegin();

  const_reverse_iterator rbegin() const;

  reverse_iterator rend();

  const_reverse_iterator rend() const;
  
  //  end Iterator access to elements

  //  end SparseVectorTemplateInterface


  void resize(size_type size, size_type max_nz, difference_type offset = 0);

  void uninitialized_resize(size_type size, size_type nz, size_type max_nz, difference_type offset = 0);

  size_type max_nz() const;

  void add_element(element_type ele);

  void insert_element(element_type ele);

  void assume_sorted(bool assume_is_sorted);

  void sort();

  void assert_valid_and_sorted() const;



  element_type* lookup_element(size_type i);
  const element_type* lookup_element(size_type i) const;



  operator SparseVectorSlice<T_Element>();

  operator const SparseVectorSlice<T_Element>() const;

  SparseVectorSlice<T_Element> operator()();

  const SparseVectorSlice<T_Element> operator()() const;


  SparseVectorSlice<T_Element> operator()(const Range1D& rng);

  const SparseVectorSlice<T_Element> operator()(const Range1D& rng) const;


  SparseVectorSlice<T_Element> operator()(size_type lbound, size_type ubound);

  const SparseVectorSlice<T_Element> operator()(size_type lbound, size_type ubound) const;


private:

  // /////////////////////////////////////////////////////////////////////////
  // Private types

  typedef SparseVectorUtilityPack::SpVecIndexLookup<element_type> SpVecIndexLookup;

  // /////////////////////////////////////////////////////////////////////////
  // Private data members

  allocator_type      alloc_;     // allocator used to allocate memory
  size_type       size_;      // the number of elements in the full vector
  size_type       max_nz_;    // the amount of storage that has been allocated
//  commented out because of problems with MS Visual C++ 5.0
//  std::vector<element_type, allocator_type> ele_;
  SpVecIndexLookup    index_lookup_;  // Acts as storage for elements and caching of searches.
  bool          assume_sorted_; // true if the client said that you can assume sorted.
  bool          know_is_sorted_; // true if it has been varified that is sorted.

  // //////////////////////////
  // Private member functions

  // Throw a NotSortedException of is_sorted() == false
  void assert_is_sorted() const {
    SparseVectorUtilityPack::assert_is_sorted(is_sorted());
  }

  void assert_space(size_type n) const {
#ifdef LINALGPACK_CHECK_SLICE_SETUP
    if(index_lookup_.nz() + n > max_nz_)
      throw OutOfRoomException("SparseVector<T_Element,T_Alloc>::assert_space():  There is not storage for this many elements");
#endif
  }

  void assert_sized_with_mem_set() const {
    if(!dim())
      throw UnsizedException("SparseVector<...>::assert_sized_with_mem_set() : "
        "Error: The sparse vector is unsized");
    if(!index_lookup_.ele()) {
      throw NoNonZeroElementsException("SparseVector<...>::assert_sized_with_mem_set() : "
        "Error: There is no memory set.");
    }
  }

  SparseVectorSlice<T_Element> get_whole_sp_vec() {
    return SparseVectorSlice<T_Element>(index_lookup_.ele(), index_lookup_.nz()
          , index_lookup_.offset(), size_, is_sorted());
  }

  const SparseVectorSlice<T_Element> get_whole_sp_vec() const {
    return SparseVectorSlice<T_Element>(index_lookup_.ele(), index_lookup_.nz()
          , index_lookup_.offset(), size_, is_sorted());
  }

  SparseVectorSlice<T_Element> get_slice(const Range1D& rng) const {
    assert_is_sorted();
    return create_slice(index_lookup_, size_, rng);
  }

};  // end class SparseVector

// ///////////////////////////////////////////////////////////////////////
template <class T_Element>
class SparseVectorSlice {
public:

  typedef T_Element                   element_type;
  typedef AbstractLinAlgPack::size_type         size_type;
  typedef ptrdiff_t                   difference_type;
  typedef element_type*                 iterator;
  typedef const element_type*               const_iterator;

#if 0 /* defined(_WINDOWS) || defined(_INTEL_CXX) */

  typedef std::reverse_iterator<iterator, element_type
    , element_type&, element_type*, difference_type>  reverse_iterator;

  typedef std::reverse_iterator<const_iterator
    , element_type, const element_type&
    , const element_type*, difference_type>       const_reverse_iterator;

#else

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

#endif

  typedef SparseVectorUtilityPack::DoesNotExistException  DoesNotExistException;
  typedef SparseVectorUtilityPack::NotSortedException   NotSortedException;



  SparseVectorSlice(element_type ele[], size_type nz, difference_type offset, size_type size
    , bool assume_sorted = false);


  void bind(SparseVectorSlice svs);


  EOverLap overlap(const SparseVectorSlice<T_Element>& sv) const;


  size_type dim() const;

  size_type nz() const;

  difference_type offset() const;

  bool is_sorted() const;

  iterator begin();

  const_iterator begin() const;

  iterator end();

  const_iterator end() const;

  reverse_iterator rbegin();

  const_reverse_iterator rbegin() const;

  reverse_iterator rend();

  const_reverse_iterator rend() const;



  element_type* lookup_element(size_type i);
  const element_type* lookup_element(size_type i) const;




  SparseVectorSlice<T_Element>& operator()();

  const SparseVectorSlice<T_Element>& operator()() const;

  SparseVectorSlice* operator&()
  {    return this; }

  const SparseVectorSlice* operator&() const
  {    return this; }


  SparseVectorSlice<T_Element> operator()(const Range1D& rng);

  const SparseVectorSlice<T_Element> operator()(const Range1D& rng) const;


  SparseVectorSlice<T_Element> operator()(size_type lbound, size_type ubound);

  const SparseVectorSlice<T_Element> operator()(size_type lbound, size_type ubound) const;


private:
  // /////////////////////////////////////////////////////////////////////////
  // Private types

  typedef SparseVectorUtilityPack::SpVecIndexLookup<element_type> index_lookup_type;

  // /////////////////////////////////////////////////////////////////////////
  // Private data members

  index_lookup_type   index_lookup_;  // Acts as storage and cacheing
  size_type       size_;      // size of the full vector
  bool          assume_sorted_; // true if the client said that you can assume sorted.


  // /////////////////////////////////////////////////////////////////////////
  // Private member functions

  // Throw a NotSortedException of is_sorted() == false
  void assert_is_sorted() const {
    SparseVectorUtilityPack::assert_is_sorted(is_sorted());
  }

  SparseVectorSlice<T_Element> get_slice(const Range1D& rng) const {
    assert_is_sorted();
    return create_slice(index_lookup_, size_, rng);
  }

  SparseVectorSlice();
  SparseVectorSlice<element_type>& operator=(const SparseVectorSlice<element_type>&);

};  // end class SparseVectorSlice

// ////////////////////////////////////////////////
// Non-member non-public utility functions

namespace SparseVectorUtilityPack {

template< class T_Element >
inline const T_Element* lookup_element( const SpVecIndexLookup<T_Element>& index_lookup
  , typename SpVecIndexLookup<T_Element>::index_type index, bool is_sorted )
{
  size_type poss; 
  return ( ( poss = index_lookup.find_element(index,is_sorted) ) < index_lookup.nz() )
    ? index_lookup.ele() + poss
    : NULL;
}

} // end namespace SparseVectorUtilityPack

// /////////////////////////////////////////////////////////////////////////////////////
// Inline members for SparseVector<>

// constructors

template <class T_Element, class T_Alloc>
inline SparseVector<T_Element,T_Alloc>::SparseVector(const allocator_type& alloc)
  : alloc_(alloc), size_(0), max_nz_(0), assume_sorted_(false), know_is_sorted_(false)
{}

template <class T_Element, class T_Alloc>
inline SparseVector<T_Element,T_Alloc>::SparseVector(bool assume_sorted,const allocator_type& alloc)
  : alloc_(alloc), size_(0), max_nz_(0), assume_sorted_(assume_sorted), know_is_sorted_(false)
{}

template <class T_Element, class T_Alloc>
inline SparseVector<T_Element,T_Alloc>::SparseVector(size_type size, size_type max_nz
    , difference_type offset, bool assume_sorted, const allocator_type& alloc)
  : alloc_(alloc), size_(0), max_nz_(0), assume_sorted_(assume_sorted), know_is_sorted_(false)
{
  resize(size,max_nz,offset);
}

template <class T_Element, class T_Alloc>
inline SparseVector<T_Element,T_Alloc>::~SparseVector() {
  resize(0,0);
}

// SparseVectorTemplateInterface for linear algebra operations

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::size_type SparseVector<T_Element,T_Alloc>::dim() const {
  return size_;
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::size_type SparseVector<T_Element,T_Alloc>::nz() const {
  return index_lookup_.nz();
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::difference_type SparseVector<T_Element,T_Alloc>::offset() const {
  return index_lookup_.offset();
}

template <class T_Element, class T_Alloc>
inline bool SparseVector<T_Element,T_Alloc>::is_sorted() const {
  return nz() <= 1 || assume_sorted_ || know_is_sorted_;
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::iterator SparseVector<T_Element,T_Alloc>::begin() {
  return index_lookup_.nz() ? index_lookup_.ele() : NULL;
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::const_iterator SparseVector<T_Element,T_Alloc>::begin() const {
  return index_lookup_.nz() ? index_lookup_.ele() : NULL;
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::iterator SparseVector<T_Element,T_Alloc>::end() {
  return index_lookup_.nz() ? index_lookup_.ele() + index_lookup_.nz() : NULL;
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::const_iterator SparseVector<T_Element,T_Alloc>::end() const {
  return index_lookup_.nz() ? index_lookup_.ele() + index_lookup_.nz() : NULL;
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::reverse_iterator SparseVector<T_Element,T_Alloc>::rbegin() {
  return reverse_iterator(end());
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::const_reverse_iterator SparseVector<T_Element,T_Alloc>::rbegin() const {
  return const_reverse_iterator(end());
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::reverse_iterator SparseVector<T_Element,T_Alloc>::rend() {
  return reverse_iterator(begin());
}

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::const_reverse_iterator SparseVector<T_Element,T_Alloc>::rend() const {
  return const_reverse_iterator(begin());
}

// Element setup and modification

template <class T_Element, class T_Alloc>
inline typename SparseVector<T_Element,T_Alloc>::size_type SparseVector<T_Element,T_Alloc>::max_nz() const {
  return max_nz_;
}

template <class T_Element, class T_Alloc>
inline void SparseVector<T_Element,T_Alloc>::add_element(element_type ele) {
  assert_space(1);
  assume_sorted_ = know_is_sorted_ = false;
#ifdef _PG_CXX
  new (index_lookup_.ele() + index_lookup_.nz()) element_type;
#else
  alloc_.construct(index_lookup_.ele() + index_lookup_.nz(), ele);
#endif
  index_lookup_.incr_nz();
}

template <class T_Element, class T_Alloc>
inline void SparseVector<T_Element,T_Alloc>::assume_sorted(bool assume_is_sorted) {
  assume_sorted_ = assume_is_sorted;
}

// Lookup an element

template <class T_Element, class T_Alloc>
inline
typename SparseVector<T_Element,T_Alloc>::element_type*
SparseVector<T_Element,T_Alloc>::lookup_element(size_type i)
{
  return const_cast<element_type*>(SparseVectorUtilityPack::lookup_element(index_lookup_,i,assume_sorted_));
}

template <class T_Element, class T_Alloc>
inline
const typename SparseVector<T_Element,T_Alloc>::element_type*
SparseVector<T_Element,T_Alloc>::lookup_element(size_type i) const
{
  return SparseVectorUtilityPack::lookup_element(index_lookup_,i,assume_sorted_);
}

// Creating a slice (subregion) of the sparse vector

template <class T_Element, class T_Alloc>
inline SparseVector<T_Element,T_Alloc>::operator SparseVectorSlice<T_Element>() {
  return get_whole_sp_vec();
}

template <class T_Element, class T_Alloc>
inline SparseVector<T_Element,T_Alloc>::operator const SparseVectorSlice<T_Element>() const {
  return get_whole_sp_vec();
}

template <class T_Element, class T_Alloc>
inline SparseVectorSlice<T_Element> SparseVector<T_Element,T_Alloc>::operator()() {
  return get_whole_sp_vec();
}

template <class T_Element, class T_Alloc>
inline const SparseVectorSlice<T_Element> SparseVector<T_Element,T_Alloc>::operator()() const {
  return get_whole_sp_vec();
}

template <class T_Element, class T_Alloc>
inline SparseVectorSlice<T_Element> SparseVector<T_Element,T_Alloc>::operator()(const Range1D& rng) {
  return get_slice(rng);
}

template <class T_Element, class T_Alloc>
inline const SparseVectorSlice<T_Element> SparseVector<T_Element,T_Alloc>::operator()(const Range1D& rng) const {
  return get_slice(rng);
}

template <class T_Element, class T_Alloc>
inline SparseVectorSlice<T_Element> SparseVector<T_Element,T_Alloc>::operator()(size_type lbound, size_type ubound) {
  return get_slice(Range1D(lbound,ubound));
}

template <class T_Element, class T_Alloc>
inline const SparseVectorSlice<T_Element> SparseVector<T_Element,T_Alloc>::operator()(size_type lbound, size_type ubound) const {
  return get_slice(Range1D(lbound,ubound));
}

// /////////////////////////////////////////////////////////////////////////////////////
// Inline members for SparseVectorSlice<>

// Constuctors

template <class T_Element>
inline SparseVectorSlice<T_Element>::SparseVectorSlice(element_type ele[], size_type nz
    , difference_type offset, size_type size, bool assume_sorted)
  : index_lookup_(ele,nz,offset), size_(size), assume_sorted_(assume_sorted)
{}

template <class T_Element>
inline void SparseVectorSlice<T_Element>::bind(SparseVectorSlice svs)
{
  index_lookup_ = svs.index_lookup_;
  size_     = svs.size_;
  assume_sorted_  = svs.assume_sorted_;
}

// Sparse Vector Templated interface for linear algebra operations

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::size_type SparseVectorSlice<T_Element>::dim() const {
  return size_;
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::size_type SparseVectorSlice<T_Element>::nz() const {
  return index_lookup_.nz();
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::difference_type SparseVectorSlice<T_Element>::offset() const {
  return index_lookup_.offset();
}

template <class T_Element>
inline bool SparseVectorSlice<T_Element>::is_sorted() const {
  return nz() <= 1 || assume_sorted_;
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::iterator SparseVectorSlice<T_Element>::begin() {
  return index_lookup_.ele();
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::const_iterator SparseVectorSlice<T_Element>::begin() const {
  return index_lookup_.ele();
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::iterator SparseVectorSlice<T_Element>::end() {
  return index_lookup_.ele() + index_lookup_.nz();
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::const_iterator SparseVectorSlice<T_Element>::end() const {
  return index_lookup_.ele() +  index_lookup_.nz();
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::reverse_iterator SparseVectorSlice<T_Element>::rbegin() {
  return reverse_iterator(end());
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::const_reverse_iterator SparseVectorSlice<T_Element>::rbegin() const {
  return const_reverse_iterator(end());
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::reverse_iterator SparseVectorSlice<T_Element>::rend() {
  return reverse_iterator(begin());
}

template <class T_Element>
inline typename SparseVectorSlice<T_Element>::const_reverse_iterator SparseVectorSlice<T_Element>::rend() const {
  return const_reverse_iterator(begin());
}

// Lookup an element

template <class T_Element>
inline
typename SparseVectorSlice<T_Element>::element_type*
SparseVectorSlice<T_Element>::lookup_element(size_type i)
{
  return const_cast<element_type*>(SparseVectorUtilityPack::lookup_element(index_lookup_,i,assume_sorted_));
}

template <class T_Element>
inline
const typename SparseVectorSlice<T_Element>::element_type*
SparseVectorSlice<T_Element>::lookup_element(size_type i) const
{
  return SparseVectorUtilityPack::lookup_element(index_lookup_,i,assume_sorted_);
}

// Creating a slice (subregion) of the sparse vector

template <class T_Element>
inline SparseVectorSlice<T_Element>& SparseVectorSlice<T_Element>::operator()() {
  return *this;
}

template <class T_Element>
inline const SparseVectorSlice<T_Element>& SparseVectorSlice<T_Element>::operator()() const {
  return *this;
}

template <class T_Element>
inline SparseVectorSlice<T_Element> SparseVectorSlice<T_Element>::operator()(const Range1D& rng) {
  return get_slice(rng);
}

template <class T_Element>
inline const SparseVectorSlice<T_Element> SparseVectorSlice<T_Element>::operator()(const Range1D& rng) const {
  return get_slice(rng);
}

template <class T_Element>
inline SparseVectorSlice<T_Element> SparseVectorSlice<T_Element>::operator()(size_type lbound, size_type ubound) {
  return get_slice(Range1D(lbound,ubound));
}

template <class T_Element>
inline const SparseVectorSlice<T_Element> SparseVectorSlice<T_Element>::operator()(size_type lbound, size_type ubound) const {
  return get_slice(Range1D(lbound,ubound));
}

} // end namespace AbstractLinAlgPack 

#endif // SPARSE_VECTOR_CLASS_DECL_H