dense2d.h

很好用的库

  /* Access Methods */
  

  /* Iterator Access Methods */
  
  //: Return an iterator pointing to the first 1D container
  inline iterator begin() {
    return iterator(data(), ld_, 0, starts, offset);
  }
  //: Return an iterator pointing past the end of the 2D container
  inline iterator end() {
    return iterator(data(), ld_, offset.twod_length(), starts, offset);
  }
  //: Return a const iterator pointing to the first 1D container
  inline const_iterator begin() const {
    return const_iterator(data(), ld_, 0, starts, offset);
  }
  //: Return a const iterator pointing past the end of the 2D container
  inline const_iterator end() const {
    return const_iterator(data(), ld_, offset.twod_length(), 
                          starts, offset);
  }

  /* reverse iterators */

  //: Return a reverse iterator pointing to the last 1D container
  inline reverse_iterator rbegin() {
    return reverse_iterator(end());
  }
  //: Return a reverse iterator pointing past the start of the 2D container
  inline reverse_iterator rend() {
    return reverse_iterator(begin());
  }
  //: Return a const reverse iterator pointing to the last 1D container
  inline const_reverse_iterator rbegin() const {
    return const_reverse_iterator(end());
  }
  //: Return a const reverse iterator pointing past the start of the 2D container
  inline const_reverse_iterator rend() const {
    return const_reverse_iterator(begin());
  }


  /* Element Access Methods */
  //: Return a reference to the (i,j) element, where (i,j) is in the 2D coordinate system  
  inline const elt_type& operator()(size_type i, size_type j) const {
    return *(data() + offset.elt(i, j));
  }
  //: Return a const reference to the (i,j) element, where (i,j) is in the 2D coordinate system  
  inline elt_type& operator()(size_type i, size_type j) {
    return *(data() + offset.elt(i, j));
  }
  
  /* Size Methods */

  //: Number of non-zeroes
  inline size_type nnz() const { return offset.major() * offset.minor(); }

  //: Capacity
  inline size_type capacity() const { return offset.major() * offset.minor(); }
  
  //: Major axis size
  inline size_type major() const { return offset.major(); }

  //: Minor axis size
  inline size_type minor() const { return offset.minor(); }

  //: Leading Dimension
  inline size_type ld() const { return ld_; }
 
  //: Memory Access
  inline const elt_type* data() const { return &(*data_)[0]; }
  inline elt_type* data() { return &(*data_)[0]; }

  /* obsolete
  inline const elt_type* get_contiguous() const { return data_->data(); }
  inline elt_type* get_contiguous() { return data_->data(); }

  inline void set_contiguous(elt_type*) { }
  */

  /* Vector Access Methods */
  
  //: OneD Access
  inline OneD operator[](size_type i) const {
    typedef OneD* oned_ptr;
    typedef InnerOneD* inner_oned_ptr;
    return __bracket<is_strided>()((elt_type*)data() + offset.oned_offset(i),
                                   offset.oned_length(i), 
                                   starts.first, ld_,
                                   oned_ptr(), 
                                   inner_oned_ptr());
  }



  /* All the submatrix stuff is in matrix_implementation for now

  inline submatrix_type sub_matrix(size_type m_start, size_type m_finish,
                              size_type n_start, size_type n_finish) {
    return submatrix_type(data_->data() + m_start * ld_ + n_start,
                     dim_type(m_finish - m_start, n_finish - n_start), ld_);
  }
  inline submatrix_type sub_matrix(size_type m_start, size_type n_start,
                              size_type m, size_type n) {
    return submatrix_type(data_->data() + m_start * ld_ + n_start, 
                     dim_type(m, n), ld_);
  }

  inline submatrix_type section(size_type m_start, size_type n_start,
                           size_type m, size_type n) {
    return submatrix_type(data_->data() + m_start * ld_ + n_start, 
                     dim_type(m, n), ld_, dim_type(m_start, n_start));
  }
  typedef range<size_type> Range;
  inline generic_dense2D operator()(Range m, Range n) {
    return generic_dense2D(data_->data() + m.start * ld_ + n.start, 
                m.finish - m.start, n.finish - n.start, ld_);    
  }

  inline OneD::subrange_type operator()(size_type i, Range n) {
    return operator[i](n);
  }
  typedef strided1D< InnerOneD > MinorVector;
  inline MinorVector minor_vector(size_type i) const {
    InnerOneD vec((elt_type*)data_->data() + i, 
                  offset.major() * ld_, starts.first);
    return MinorVector(vec, ld_);
  }
  inline MinorVector::subrange_type operator()(Range m, size_type j) {
    return minor_vector(j)(m);
  }

  */

  /*JGS friend not working for transpose constructor  
  protected:
  */
  size_type ld_;/* JGS redundant */
  rep_ptr data_;
  pair_type starts;
  Offset offset;
};


template <class T, class OffsetGen, int M, int N>
struct gen_dense2D;

/* why didn't I use std::vector here?
 or perhaps I should just use plain old memory here?
 will that work with the reference counting in terms
 of deallocating?
 */

//: Dense2D Storage Type
//
// Inherits from generic_dense2D. The class "owns" its data.
//
//!category: containers
//!component: type
//!tparam: T - the element type
//!tparam: OffsetGen - the Offset class generator
//!tparam: MM - For static sized matrix, the major dimension
//!tparam: NN - For static sized matrix, the minor dimension
//!models: TwoDStorage

template<class T, class OffsetGen, int MM = 0, int NN = 0>
class dense2D
 : public generic_dense2D< std::vector<T> ,
               refcnt_ptr< std::vector<T> >, OffsetGen, MM, NN >
/* : public generic_dense2D< bare_bones_array<T> ,
               refcnt_ptr< bare_bones_array<T> >, OffsetGen, MM, NN >
*/
{
public:
  typedef generic_dense2D< std::vector<T> ,
               refcnt_ptr< std::vector<T> >, OffsetGen, MM, NN> super;
  /*  typedef generic_dense2D< bare_bones_array<T> ,
                 refcnt_ptr< bare_bones_array<T> >, OffsetGen, MM, NN> super;
  */
  typedef typename super::Offset Offset;
  //: Pair type for dimension
  typedef typename Offset::dim_type dim_type;
  //: Pair type for bandwidth
  typedef typename Offset::band_type band_type;
  typedef typename super::reptype reptype;
  typedef typename super::rep_ptr rep_ptr;
  //: Unsigned integral type for dimensions and indices
  typedef typename super::size_type size_type;
  //: The transpose type
  typedef dense2D<T, typename OffsetGen::transpose_type,
                  MM, NN> transpose_type;
// VC++ doesn't like this
  //friend class transpose_type;
  //: This has internal storage
  typedef internal_tag storage_loc;

  //: Default Constructor
  inline dense2D() { }

  //: Constructor from Dimension Pair
  inline dense2D(dim_type dim)
    : super(new reptype(Offset::size(dim.first(), dim.second(), 0, 0)),
            dim.first(), 
            dim.second(),
            dim.second()) { }

  //: Constructor from Dimension and Bandwidth Pairs
  inline dense2D(dim_type dim, band_type bw)
    : super(new reptype(Offset::size(dim.first(),dim.second(),
                                     bw.first(), bw.second())), 
            dim.first(), 
            dim.second(), 
            dim.second(),
            bw) { }

  //: Copy Constructor
  inline dense2D(const dense2D& x)
    : super(x) { }

  //: Assignment Operator
  inline dense2D& operator=(const dense2D& x) {
    super::operator=(x);
    return *this;
  }

  //: Transpose Constructor
  inline dense2D(const transpose_type& x, do_transpose t, do_transpose)
    : super(x, t, t) { }

#if !defined(_MSVCPP_)
  // JGS, use actual stream types
  //: Matrix Stream Constructor
  template <class MatrixStream, class Orien>
  inline dense2D(MatrixStream& s, Orien)
    : super(new reptype(Offset::size(Orien::map(dim_type(s.nrows(), 
                                                         s.ncols())).first(),
                                     Orien::map(dim_type(s.nrows(), 
                                                         s.ncols())).second(),
                                     0, 0)),
            s,
            Orien()) { }

  //: Matrix Stream Constructor with bandwidth
  template <class MatrixStream, class Orien>
  inline dense2D(MatrixStream& s, Orien, band_type bw)
    : super(new reptype(Offset::size(Orien::map(dim_type(s.nrows(), 
                                                         s.ncols())).first(),
                                     Orien::map(dim_type(s.nrows(), 
                                                         s.ncols())).second(),
                                     bw.first(), bw.second())),
            s,
            Orien(), 
            bw) { }
#endif

#if 0
  // deprecated
  template <class SubMatrix>
  struct partitioned {
    typedef dense2D<SubMatrix, OffsetGen> type;
    typedef gen_dense2D<SubMatrix, OffsetGen> generator;
  };
#endif
  
#if 1 // This makes no sense. dense2D can not be a "view"
  //: banded view constructor
  template <class TwoD>
  inline dense2D(const TwoD& x, band_type bw, banded_tag)
    : super(x.data_, x, bw, banded_tag()) { }
#endif

  //: Destructor
  inline ~dense2D() { }

  inline void resize(size_type m, size_type n) {
    rep_ptr newdata = new reptype(Offset::size(m, n, 0, 0));
    size_type i, j;
    size_type M = MTL_MIN(m, offset.major());
    size_type N = MTL_MIN(n, offset.minor());
    for (i = 0; i < M; ++i)
      for (j = 0; j < N; ++j)
	(*newdata)[i * n + j] = (*this)(i,j);
    for (; i < m; ++i)
      for (; j < n; ++j)
      (*newdata)[i * n + j] = T();

    data_ = newdata;
    ld_ = n;
    offset.dim = dim_type(m, n);
    offset.ld = n;
  }

};

template <class T, class OffsetGen, int M, int N>
struct gen_external2D; 

#ifndef MTL_DISABLE_BLOCKING
template <class Block, class OffsetGen, int M, int N>
struct gen_block2D;
#endif

//: blah
//!noindex:
template <class T, class OffsetGen, int M, int N>
struct gen_dense2D {
  typedef gen_dense2D<T, typename OffsetGen::transpose_type,N,M> transpose_type;
  typedef gen_external2D<T, OffsetGen,M,N> submatrix_type;

#ifndef MTL_DISABLE_BLOCKING
  template <class Block>
  struct blocked_view {
    typedef gen_block2D<Block, OffsetGen, M, N> type;
  };
#endif

  typedef gen_dense2D<T, typename OffsetGen::banded_view_type,M,N>
           banded_view_type;

  typedef dense2D<T, OffsetGen, M, N> type;
};


//: External2D Storage Type
//
// Inherits from generic_dense2D. The class does not "own" its data.
//
//!category: containers
//!component: type
//!tparam: T - the element type
//!tparam: OffsetGen - the Offset class generator
//!tparam: MM - For static sized matrix, the major dimension
//!tparam: NN - For static sized matrix, the minor dimension
//!models: TwoDStorage
//
template <class T, class OffsetGen, int MM = 0, int NN = 0>
class external2D
 : public generic_dense2D< external_vec<T,NN>,
                           external_vec<T,NN>*, OffsetGen, MM, NN >
{
  typedef generic_dense2D< external_vec<T,NN>,
                           external_vec<T,NN>*, OffsetGen, MM, NN > super;
public:
  external_vec<T,NN> rep;
  typedef dimension<T> dyn_dim;
  typedef typename super::Offset Offset;
  //: Pair type for dimension
  typedef typename Offset::dim_type dim_type;
  //: Pair type for bandwidth
  typedef typename Offset::band_type band_type;

  typedef typename super::reptype reptype;

  //: Unsigned integral type for dimensions and indices
  typedef typename super::size_type size_type;
  //: Type for the transpose
  typedef external2D<T, typename OffsetGen::transpose_type, 
                     MM, NN> transpose_type;
// VC++ doesn't like this
  //friend class transpose_type;
  //: This has external storage
  typedef external_tag storage_loc;

  //: Default Constructor
  inline external2D() { }

  //: Construct from pointer and dimensions
  inline external2D(T* data, dim_type dim)
    : super(&rep, dim.first(), dim.second(), dim.second()), 
      rep(data, dim.first() * dim.second()) 
  { }

  //: Construct from pointer, dimensions, and leading dimension
  inline external2D(T* data, dim_type dim, size_type ld)
    : super(&rep, dim.first(), dim.second(), ld),
      rep(data, dim.first() * ld)
  { }

  //: non-zero indices in upper left corner
  inline external2D(T* data, dim_type dim, size_type ld, 
		    dyn_dim s, char)
    : super(&rep, dim.first(), dim.second(), ld, s, char()),
      rep(data, dim.first() * ld)
  { }

  //: Constructor with bandwith
  inline external2D(T* data, dim_type dim, band_type bw)
    : super(&rep, dim.first(), dim.second(), dim.second(), bw),
      rep(data, dim.first() * dim.second())
  { }
  //: Constructor with leading dimension and bandwith
  inline external2D(T* data, dim_type dim, size_type ld, band_type bw)
    : super(&rep, dim.first(), dim.second(), ld, bw),
      rep(data, dim.first() * ld)
  { }

  //: Copy Constructor
  inline external2D(const external2D& x)
    :  super(&rep, x), rep(x.rep)
  { }

  //: Assignment Operator
  inline external2D& operator=(const external2D& x) {
    rep = x.rep;
    super::operator=(x);
    data_ = &rep;
    return *this;
  }

  //: Transpose Constructor
  inline external2D(const transpose_type& x, do_transpose t, do_transpose)
    : super(x, t, t), rep(x.rep) { }

  /* JGS This conflicts with the external2D(T* data, dim_type dim,
     band_type bw) constructor, and I am not sure this is really
     needed anyway.

  //: Matrix Stream Constructor
  template <class MatrixStream, class Orien>
  inline external2D(T* data, MatrixStream& s, Orien o)
    : rep(data, s.nrows() * s.ncols()), super(&rep, s, o) { }
  template <class MatrixStream, class Orien>
  inline external2D(T* data, MatrixStream& s, Orien o, 
                    band_type bw)
    : super(&rep, s, o, bw),
      rep(data, s.nrows() * s.ncols()) 
  { }  
  */

  //: banded view constructor
  template <class TwoD>
  inline external2D(const TwoD& x, band_type bw, banded_tag) 
    : super(&rep, x, bw, banded_tag()), rep((T*)x.data(), x.major() * x.ld()) { }


  inline ~external2D() { }
#if 0
  // deprecated
  template <class SubMatrix>
  struct partitioned {
    typedef dense2D<SubMatrix, OffsetGen> type;
    typedef gen_dense2D<SubMatrix, OffsetGen> generator;
  };
#endif
};

//: blah
//!noindex:
template <class T, class OffsetGen, int M, int N>
struct gen_external2D {
  typedef gen_external2D<T, typename OffsetGen::transpose_type, N, M> transpose_type;
  typedef gen_external2D<T, OffsetGen,M,N> submatrix_type;
  typedef gen_external2D<T, typename OffsetGen::banded_view_type,M,N>
           banded_view_type;

#ifndef MTL_DISABLE_BLOCKING
  template <class Block>
  struct blocked_view {
    typedef gen_block2D<Block, OffsetGen, M, N> type;
  };
#endif

  typedef external2D<T, OffsetGen, M, N> type;
};



} /* namespace mtl */


#endif /* MTL_DENSE2D_H */