#ifndef _RHEO_DISARRAY_H
#define _RHEO_DISARRAY_H
//
// This file is part of Rheolef.
//
// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
//
// Rheolef is free software; you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation; either version 2 of the License, or
// (at your option) any later version.
//
// Rheolef 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 General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with Rheolef; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
// 
// =========================================================================

namespace rheolef {
/**
@linalgclassfile disarray distributed container

Description
===========
This class provides a `std::vector` like
container for a distributed memory model array.

The `disarray` interface is similar to those of the `std::vector`
with the addition of some communication features 
for a distributed memory model.

Example
=======

        int main (int argc, char**argv) {
          environment distributed(argc, argv);
          distributor ownership (100);
          disarray<double> x (ownership, 3.14);
          dout << x << endl;
        }

Notation
========
There are two kind of indexes:

`dis_i`
>        This index refers to the complete array.
>        It is valid on all processes.
`i` 
>	 This index refers to the subset of the array
>        that is owned by the current process.
>        It is valid only on the local process.

Read and write accessors to the 
subset of the array that is owned by the current process
are still denoted by square brackets,
e.g. `value = x[i]` and `x[i] = value`, respectively.

Read and write accessors to
the complete array,
including array subsets that are owned by some others processors,
are introduced with the help of new member functions:

`x.dis_entry(dis_i) = value` 
>	write access at any location
`value = x.dis_at(dis_i)`
>	read access at any location

In order to optimize communications, they should be grouped.

Grouped writes
==============
Loop on any `dis_i` and write your value:

        for (...) { 
          x.dis_entry (dis_i) = value;
        }
 
Finally, perform all communications in one pass:
 
	x.dis_entry_assembly();

After this command, each value is stored in `x`,
at its right place, depending on `dis_i` and its `ownership`.
Note that, when `dis_i` is owned by the current process, 
the value is directly written as `x[i] = value`
and no communication are generated.

Grouped reads
=============
First, define the set of indexes that you want to access:

        std::set<size_t> dis_i_set; 

Then, loop on `dis_i` all these indexes and append it to this set:

        for (...) {
          dis_i_set.insert (dis_i);
        }

Next, perform communications:

        x.set_dis_indexes (dis_i_set);

After this command, each values associated to the `dis_i` index,
and that belongs to the index set, is also available also on the
current processor:
 
        for (...) {
          value = x.dis_at (dis_i);
        }

Note that, when `dis_i` is owned by the current process, 
the value is directly read as `value = x[i]`
and no communication are generated.

Template parameters
===================
Recall that the `std::class`
that takes two template parameters,
a `T` one for the stored elements and
a `A` one for the memory allocator,
the present `disarray` class
take three template parameters:

- `T`: the stored object type
- `M`: the memory model, i.e. `sequential` or `distributed`,
       by default `default_memory_model` with is 
       defined at the @ref configuration_page stage
- `A`: the memory allocator, by default `std::allocator`

Constructor
===========
In the `sequential` case,
the class interface provides a simplified constructor:

	int local_size = 100;
        disarray<double,sequential> x (local_size, init_val);
   
This declaration is similar to those of a `std::vector` one:
no communications are possible.
In order to enable communication, your have to replace the
`local_size` by information on how the array is distributed in memory:

          distributor ownership (100);
          disarray<double> x (ownership, 3.14);

The @ref distributor_4 class does this job.

Implementation
==============
@showfromfile
@snippet disarray.h verbatim_disarray
*/
} // namespace rheolef


/* 
Misc notes:
  "scatter" via "get_dis_entry".

  "gather" via "dis_entry(dis_i) = value"
  or "dis_entry(dis_i) += value". Note that += applies when T=idx_set where
  idx_set is a wrapper class of std::set<size_t> ; the += operator represents the
  union of a set. The operator= is used when T=double or others simple T types
  without algebra. If there is a conflict, i.e. several processes set the dis_i
  index, then the result of operator+= depends upon the order of the process at
  each run and is not deterministic. Such ambiguous behavior is not detected
  yet at run time.

*/

#include "rheolef/communicator.h"
#include "rheolef/distributor.h"
#include "rheolef/diststream.h"
#include "rheolef/heap_allocator.h"
#include "rheolef/msg_util.h"
#ifdef _RHEOLEF_HAVE_MPI
#include "rheolef/mpi_pair_datatype.h"
#endif // _RHEOLEF_HAVE_MPI

namespace rheolef {
/// @brief disarray element output helper
template <class T>
struct _disarray_put_element_type {
  std::ostream& operator() (std::ostream& os, const T& x) const { return os << x; }
};
template <class T>
struct _disarray_put_matlab_type {
  std::ostream& operator() (std::ostream& os, const T& x) const { return os << x << ";"; }
};
/// @brief disarray element input helper
template <class T>
struct _disarray_get_element_type {
  std::istream& operator() (std::istream& is, T& x) const { return is >> x; }
};
} // namespace rheolef
// -------------------------------------------------------------
// the sequential representation
// -------------------------------------------------------------
namespace rheolef {

template <class T, class M, class A> class disarray_rep {};

template <class T, class A>
class disarray_rep<T,sequential,A> : public std::vector<T> { // TODO: vector<T,A> for heap_alloc
public:
    typedef T                                         value_type;
    typedef A                                         allocator_type;
    typedef typename A::difference_type               difference_type;
    typedef std::vector<T>                            base; // TODO: vector<T,A> for heap_alloc
    typedef typename base::size_type                  size_type;
    typedef typename base::iterator                   iterator;
    typedef typename base::const_iterator             const_iterator;
    typedef typename base::const_reference            const_reference;
    typedef typename base::reference                  reference;
    typedef reference	                              dis_reference;
    typedef distributor::communicator_type            communicator_type;
    typedef sequential                                memory_type;

    explicit disarray_rep (const A& alloc = A());
    disarray_rep (const distributor& ownership, const T& init_val = T(), const A& alloc = A());
    void resize   (const distributor& ownership, const T& init_val = T());
    disarray_rep (size_type loc_size = 0,       const T& init_val = T(), const A& alloc = A());
    void resize   (size_type loc_size = 0,       const T& init_val = T());
    disarray_rep (const disarray_rep<T,sequential,A>& x);

    A get_allocator() const { return base::get_allocator(); }
    size_type size() const { return base::size(); }
    iterator begin() { return base::begin(); }
    const_iterator begin() const { return base::begin(); }
    iterator end() { return base::end(); }
    const_iterator end() const { return base::end(); }
    const distributor& ownership() const { return _ownership; }

    reference       operator[] (size_type i)       { return base::operator[] (i); }
    const_reference operator[] (size_type i) const { return base::operator[] (i); }
    const_reference dis_at (size_type dis_i) const { return operator[] (dis_i); }
   
    size_type dis_size () const { return base::size(); }
    size_type first_index () const { return 0; }
    size_type last_index () const { return base::size(); }
    reference dis_entry     (size_type dis_i) { return operator[](dis_i); }
    void get_dis_indexes (std::set<size_type>& ext_idx_set) const;
    void reset_dis_indexes() const {}
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_begin (SetOp = SetOp()) {}
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_end (SetOp = SetOp()) {}
    void repartition (					      // old_numbering for *this
        const disarray_rep<size_type,sequential,A>& partition,	      // old_ownership
        disarray_rep<T,sequential,A>&               new_disarray,	      // new_ownership (created)
        disarray_rep<size_type,sequential,A>&       old_numbering,      // new_ownership
        disarray_rep<size_type,sequential,A>&       new_numbering) const // old_ownership
    {
	error_macro ("not yet");
    }
    template<class A2>
    void reverse_permutation (                                // old_ownership for *this=iold2dis_inew
        disarray_rep<size_type,sequential,A2>& inew2dis_iold) const;    // new_ownership

    idiststream& get_values (idiststream& s);
    odiststream& put_values (odiststream& s) const;
    odiststream& put_matlab (odiststream& s) const;
    template <class GetFunction> idiststream& get_values (idiststream& ips, GetFunction get_element);
    template <class PutFunction> odiststream& put_values (odiststream& ops, PutFunction put_element) const;
    void dump (std::string name) const;
protected:
// data:
    distributor      _ownership;
};
// -------------------------------------------------------------
// the distributed representation
// -------------------------------------------------------------
#ifdef _RHEOLEF_HAVE_MPI
template <class T, class A> struct disarray_dis_reference; /// forward decl

template <class T, class A>
class disarray_rep<T,distributed,A> : public disarray_rep<T,sequential,A> {
public:

// typedefs:

    using base = disarray_rep<T,sequential,A>;
    using size_type       = typename base::size_type;
    using value_type      = typename base::value_type;
    using memory_type     = distributed;
    using iterator        = typename base::iterator;
    using const_iterator  = typename base::const_iterator;
    using reference       = typename base::reference;
    using const_reference = typename base::const_reference;
    using difference_type = typename base::difference_type;
    using communicator_type = distributor::communicator_type;
    using scatter_map_type  = std::map <size_type, T>; // TODO: map<T,A> for heap_alloc
    using dis_reference   = disarray_dis_reference<T,A>;

// allocators:

    disarray_rep (const distributor& ownership, const T&  init_val = T(), const A& alloc = A());
    void resize   (const distributor& ownership, const T&  init_val = T());
    disarray_rep (const disarray_rep<T,distributed,A>& x);

    A get_allocator() const        { return base::get_allocator(); }
    size_type size() const         { return base::size(); }
    const_iterator begin() const   { return base::begin(); }
    const_iterator end() const     { return base::end(); }
    iterator begin()               { return base::begin(); }
    iterator end()                 { return base::end(); }

    const distributor& ownership() const  { return base::_ownership; }
    const mpi::communicator& comm() const { return ownership().comm(); }
    size_type first_index () const        { return ownership().first_index(); }
    size_type last_index () const         { return ownership().last_index(); }
    size_type dis_size () const           { return ownership().dis_size(); }

    dis_reference dis_entry (size_type dis_i) { return dis_reference (*this, dis_i); }

    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_begin (SetOp my_set_op = SetOp());
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_end   (SetOp my_set_op = SetOp());
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly       (SetOp my_set_op = SetOp())
   		{ dis_entry_assembly_begin (my_set_op); dis_entry_assembly_end (my_set_op); }

    template<class Set, class Map>
    void append_dis_entry (const Set& ext_idx_set, Map& ext_idx_map) const;

    template<class Set, class Map>
    void get_dis_entry (const Set& ext_idx_set, Map& ext_idx_map) const {
    	    ext_idx_map.clear();
    	    append_dis_entry (ext_idx_set, ext_idx_map);
	}

    template<class Set>
    void append_dis_indexes (const Set& ext_idx_set) const { append_dis_entry (ext_idx_set, _ext_x); }

    template<class Set>
    void set_dis_indexes    (const Set& ext_idx_set) const { get_dis_entry    (ext_idx_set, _ext_x); }
    void get_dis_indexes (std::set<size_type>& ext_idx_set) const;
    void reset_dis_indexes() const;

    const_reference dis_at (size_type dis_i) const;

    // get all external pairs (dis_i, values):
    const scatter_map_type&    get_dis_map_entries() const { return _ext_x; }

    template<class A2>
    void repartition (					      // old_numbering for *this
        const disarray_rep<size_type,distributed,A2>& partition,	      // old_ownership
        disarray_rep<T,distributed,A>&                new_disarray,	      // new_ownership (created)
        disarray_rep<size_type,distributed,A2>&       old_numbering,     // new_ownership
        disarray_rep<size_type,distributed,A2>&       new_numbering) const; // old_ownership

    template<class A2>
    void permutation_apply (		 		       // old_numbering for *this
        const disarray_rep<size_type,distributed,A2>& new_numbering,      // old_ownership
        disarray_rep<T,distributed,A>&                new_disarray) const;   // new_ownership (already allocated)

    template<class A2>
    void reverse_permutation (                                // old_ownership for *this=iold2dis_inew
        disarray_rep<size_type,distributed,A2>& inew2dis_iold) const;    // new_ownership

    idiststream& get_values (idiststream& s);
    odiststream& put_values (odiststream& s) const;
    odiststream& put_matlab (odiststream& s) const;
    template <class GetFunction> idiststream& get_values (idiststream& ips, GetFunction get_element);
    template <class PutFunction> odiststream& put_values (odiststream& ops, PutFunction put_element) const;
    template <class PutFunction, class A2> odiststream& permuted_put_values (odiststream& ops, const disarray_rep<size_type,distributed,A2>& perm,
		PutFunction put_element) const;
    void dump (std::string name) const;
//protected:
    template<class U, class SetOp>
    void set_dis_entry       (size_type dis_i, const U& val, const SetOp& set_op);
#ifdef TO_CLEAN
    template<class U, class SetOp>
    void set_plus_dis_entry  (size_type dis_i, const U& val, const SetOp& set_op);
    template<class U, class SetOp>
    void set_minus_dis_entry (size_type dis_i, const U& val, const SetOp& set_op);
#endif // TO_CLEAN
// typedefs:
    /** 1) stash: store data before assembly() communications:
      *   select multimap<U> when T=set<U> and map<T> otherwise
      */
    template<class Pair>
    struct remove_const_in_pair {
        typedef Pair type;
    };
    template<class T1, class T2>
    struct remove_const_in_pair<std::pair<T1,T2> > {
        typedef std::pair<typename std::decay<T1>::type,
                          typename std::decay<T2>::type> type;
    };
    template<class U, class IsContainer> struct stash_traits {};
    template<class U>
    struct stash_traits<U,std::false_type> {
        typedef U mapped_type;
        typedef std::map <size_type, U>   map_type; // TODO: map<T,A> for heap_alloc
    };
    template<class U>
    struct stash_traits<U,std::true_type> {
        typedef typename remove_const_in_pair<typename U::value_type>::type mapped_type;
        typedef std::multimap <size_type, mapped_type>   map_type; // TODO: map<T,A> for heap_alloc
    };
    typedef typename details::is_container_of_mpi_datatype<T>::type     is_container;
    typedef typename stash_traits<T,is_container>::mapped_type stash_value;
    typedef typename stash_traits<T,is_container>::map_type    stash_map_type;

    /** 2) message: for communication during assembly_begin(), assembly_end()
      */
    struct message_type {
        std::list<std::pair<size_type,mpi::request> >    waits; // TODO: list<T,A> for heap_alloc
        std::vector<std::pair<size_type,stash_value> >   data;  // TODO: vector<T,A> for heap_alloc
        message_type() : waits(), data() {}
    };
    /** 3) scatter (get_entry): specialized versions for T=container and T=simple type
      */
    template<class Set, class Map>
    void append_dis_entry (const Set& ext_idx_set, Map& ext_idx_map, std::true_type) const;
    template<class Set, class Map>
    void append_dis_entry (const Set& ext_idx_set, Map& ext_idx_map, std::false_type) const;

// data:
    stash_map_type   _stash;		// for assembly msgs:
    message_type     _send;
    message_type     _receive;
    size_type        _receive_max_size;
    mutable scatter_map_type _ext_x;		// for ext values (scatter)
};
// custom reference that interacts with stach
// for ext values on others partitions
template <class T, class A>
struct disarray_dis_reference
{
  using size_type = typename disarray_rep<T,distributed,A>::size_type;

  disarray_dis_reference (disarray_rep<T,distributed,A>& x, size_type dis_i)
    : _x(x), _dis_i(dis_i) {}

#ifndef TO_CLEAN // not usefull: default behavior
  disarray_dis_reference (const disarray_dis_reference<T,A>& r)
    : _x(r._x), _dis_i(r._dis_i) {}

  disarray_dis_reference<T,A>&
  operator= (const disarray_dis_reference<T,A>& r) {
    _x = r._x; _dis_i = r._dis_i; return *this; }
#endif // TO_CLEAN

  disarray_dis_reference<T,A>&
  operator= (const T& value) {
    trace_macro("_x.set_dis_entry (dis_i="<<_dis_i<<", value="<<value<<", op=)");
    _x.set_dis_entry (_dis_i, value, details::generic_set_op());
    return *this;
  }
  template<class U>
  disarray_dis_reference<T,A>&
  operator+= (const U& value) {
    trace_macro("_x.set_dis_entry (dis_i="<<_dis_i<<", value="<<value<<", op+=)");
    _x.set_dis_entry (_dis_i, value, details::generic_set_plus_op());
    return *this;
  }
  template<class U>
  disarray_dis_reference<T,A>&
  operator-= (const U& value) {
    trace_macro("_x.set_dis_entry (dis_i="<<_dis_i<<", value="<<value<<", op-=)");
    _x.set_dis_entry (_dis_i, value, details::generic_set_minus_op());
    return *this;
  }
// data:
protected:
  disarray_rep<T,distributed,A>& _x;
  size_type                      _dis_i;
};

namespace details {
template <class T, class A>
struct is_class_reference<disarray_dis_reference<T,A>> : std::true_type {};
} // namespace details

#endif // _RHEOLEF_HAVE_MPI
// -------------------------------------------------------------
// the basic class with a smart pointer to representation
// the user-level class with memory-model parameter
// -------------------------------------------------------------
template <class T, class M = rheo_default_memory_model, class A = std::allocator<T> >
class disarray {
public:
    typedef M memory_type;
    typedef disarray_rep<T,sequential,A>  	      rep;
    typedef typename rep::base::size_type      size_type;
    typedef typename rep::base::iterator       iterator;
    typedef typename rep::base::const_iterator const_iterator;
};
template <class T, class A>
class disarray<T,sequential,A> : public smart_pointer<disarray_rep<T,sequential,A> > {
public:

// typedefs:

    typedef disarray_rep<T,sequential,A>  	  rep;
    typedef smart_pointer<rep> 		  base;

    typedef sequential 			  memory_type;
    typedef typename rep::size_type 	  size_type;
    typedef typename rep::difference_type difference_type;
    typedef typename rep::value_type 	  value_type;
    typedef typename rep::reference 	  reference;
    typedef typename rep::dis_reference   dis_reference;
    typedef typename rep::iterator 	  iterator;
    typedef typename rep::const_reference const_reference;
    typedef typename rep::const_iterator  const_iterator;

// allocators:


    disarray       (size_type loc_size = 0,       const T& init_val = T(), const A& alloc = A());
    void resize (size_type loc_size = 0,       const T& init_val = T());
    disarray       (const distributor& ownership, const T& init_val = T(), const A& alloc = A());
    void resize (const distributor& ownership, const T& init_val = T());

// local accessors & modifiers:

    A get_allocator() const              { return base::data().get_allocator(); }
    size_type     size () const          { return base::data().size(); }
    size_type dis_size () const          { return base::data().dis_size(); }
    const distributor& ownership() const { return base::data().ownership(); }
    const communicator& comm() const     { return ownership().comm(); }

    reference       operator[] (size_type i)       { return base::data().operator[] (i); }
    const_reference operator[] (size_type i) const { return base::data().operator[] (i); }
    reference       operator() (size_type i)       { return base::data().operator[] (i); }
    const_reference operator() (size_type i) const { return base::data().operator[] (i); }
    const_reference dis_at (size_type dis_i) const { return operator[] (dis_i); }

          iterator begin()       { return base::data().begin(); }
    const_iterator begin() const { return base::data().begin(); }
          iterator end()         { return base::data().end(); }
    const_iterator end() const   { return base::data().end(); }

// global modifiers (for compatibility with distributed interface):

    dis_reference dis_entry (size_type dis_i) { return base::data().dis_entry(dis_i); }
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly (SetOp my_set_op = SetOp()) {}
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_begin (SetOp my_set_op = SetOp()) {}
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_end (SetOp my_set_op = SetOp()) {}

    void dis_entry_assembly_begin() {}
    void dis_entry_assembly_end()   {}
    void dis_entry_assembly()       {}

    void get_dis_indexes (std::set<size_type>& ext_idx_set) const { base::data().get_dis_indexes (ext_idx_set); }
    void reset_dis_indexes() const {}
    template<class Set> void set_dis_indexes    (const Set& ext_idx_set) const {}
    template<class Set> void append_dis_indexes (const Set& ext_idx_set) const {}
    template<class Set, class Map> void append_dis_entry (const Set& ext_idx_set, Map& ext_idx_map) const {}
    template<class Set, class Map> void get_dis_entry    (const Set& ext_idx_set, Map& ext_idx_map) const {}

// apply a partition:
 
    template<class RepSize>
    void repartition (			             // old_numbering for *this
        const RepSize&         partition,	     // old_ownership
        disarray<T,sequential,A>& new_disarray,	     // new_ownership (created)
        RepSize&               old_numbering,	     // new_ownership
        RepSize&               new_numbering) const  // old_ownership
        { return base::data().repartition (partition, new_disarray, old_numbering, new_numbering); }

    template<class RepSize>
    void permutation_apply (                       // old_numbering for *this
        const RepSize&          new_numbering,     // old_ownership
        disarray<T,sequential,A>&  new_disarray) const   // new_ownership (already allocated)
        { return base::data().permutation_apply (new_numbering, new_disarray); }

    void reverse_permutation (                                 // old_ownership for *this=iold2dis_inew
        disarray<size_type,sequential,A>& inew2dis_iold) const   // new_ownership
        { base::data().reverse_permutation (inew2dis_iold.data()); }

// i/o:

    odiststream& put_values (odiststream& ops) const { return base::data().put_values(ops); }
    idiststream& get_values (idiststream& ips)       { return base::data().get_values(ips); }
    template <class GetFunction>
    idiststream& get_values (idiststream& ips, GetFunction get_element)       { return base::data().get_values(ips, get_element); }
    template <class PutFunction>
    odiststream& put_values (odiststream& ops, PutFunction put_element) const { return base::data().put_values(ops, put_element); }
    void dump (std::string name) const { return base::data().dump(name); }
};
template <class T, class A>
inline
disarray<T,sequential,A>::disarray (
    	size_type loc_size,
	const T&  init_val,
        const A&  alloc)
 : base(new_macro(rep(loc_size,init_val,alloc)))
{
}
template <class T, class A>
inline
disarray<T,sequential,A>::disarray (
    	const distributor& ownership,
	const T&           init_val,
        const A&           alloc)
 : base(new_macro(rep(ownership,init_val,alloc)))
{
}
template <class T, class A>
inline
void
disarray<T,sequential,A>::resize (
    	size_type loc_size,
	const T&  init_val)
{
  base::data().resize (loc_size,init_val);
}
template <class T, class A>
inline
void
disarray<T,sequential,A>::resize (
    	const distributor& ownership,
	const T&           init_val)
{
  base::data().resize (ownership,init_val);
}
#ifdef _RHEOLEF_HAVE_MPI
// [verbatim_disarray]
template <class T, class A>
class disarray<T,distributed,A> : public smart_pointer<disarray_rep<T,distributed,A> > {
public:

// typedefs:

    typedef disarray_rep<T,distributed,A> rep;
    typedef smart_pointer<rep> 		  base;

    typedef distributed 		  memory_type;
    typedef typename rep::size_type 	  size_type;
    typedef typename rep::difference_type difference_type;
    typedef typename rep::value_type 	  value_type;
    typedef typename rep::reference 	  reference;
    typedef typename rep::dis_reference   dis_reference;
    typedef typename rep::iterator 	  iterator;
    typedef typename rep::const_reference const_reference;
    typedef typename rep::const_iterator  const_iterator;
    typedef typename rep::scatter_map_type scatter_map_type;

// allocators:

    disarray       (const distributor& ownership = distributor(), const T& init_val = T(), const A& alloc = A());
    void resize (const distributor& ownership = distributor(), const T& init_val = T());

// local accessors & modifiers:

    A get_allocator() const              { return base::data().get_allocator(); }
    size_type     size () const          { return base::data().size(); }
    size_type dis_size () const          { return base::data().dis_size(); }
    const distributor& ownership() const { return base::data().ownership(); }
    const communicator& comm() const     { return base::data().comm(); }

    reference       operator[] (size_type i)       { return base::data().operator[] (i); }
    const_reference operator[] (size_type i) const { return base::data().operator[] (i); }
    reference       operator() (size_type i)       { return base::data().operator[] (i); }
    const_reference operator() (size_type i) const { return base::data().operator[] (i); }

          iterator begin()       { return base::data().begin(); }
    const_iterator begin() const { return base::data().begin(); }
          iterator end()         { return base::data().end(); }
    const_iterator end() const   { return base::data().end(); }

// global accessor:

    template<class Set, class Map>
    void append_dis_entry (const Set& ext_idx_set, Map& ext_idx_map) const { base::data().append_dis_entry (ext_idx_set, ext_idx_map); }

    template<class Set, class Map>
    void get_dis_entry    (const Set& ext_idx_set, Map& ext_idx_map) const { base::data().get_dis_entry (ext_idx_set, ext_idx_map); }

    template<class Set>
    void append_dis_indexes (const Set& ext_idx_set) const { base::data().append_dis_indexes (ext_idx_set); }
    void reset_dis_indexes() const { base::data().reset_dis_indexes(); }
    void get_dis_indexes (std::set<size_type>& ext_idx_set) const { base::data().get_dis_indexes (ext_idx_set); }

    template<class Set>
    void set_dis_indexes    (const Set& ext_idx_set) const { base::data().set_dis_indexes (ext_idx_set); }

    const T& dis_at (size_type dis_i) const { return base::data().dis_at (dis_i); }

    // get all external pairs (dis_i, values):
    const scatter_map_type& get_dis_map_entries() const { return base::data().get_dis_map_entries(); }

// global modifiers (for compatibility with distributed interface):

    dis_reference dis_entry (size_type dis_i) { return base::data().dis_entry(dis_i); }

    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_begin (SetOp my_set_op = SetOp()) { base::data().dis_entry_assembly_begin (my_set_op); }
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly_end   (SetOp my_set_op = SetOp()) { base::data().dis_entry_assembly_end   (my_set_op); }
    template<class SetOp = typename details::default_set_op_traits<T>::type>
    void dis_entry_assembly       (SetOp my_set_op = SetOp()) { base::data().dis_entry_assembly       (my_set_op); }

    void dis_entry_assembly_begin() { base::data().template dis_entry_assembly_begin<typename details::default_set_op_traits<T>::type>(); }
    void dis_entry_assembly_end()   { base::data().template dis_entry_assembly_end<typename details::default_set_op_traits<T>::type>(); }
    void dis_entry_assembly()       { dis_entry_assembly_begin(); dis_entry_assembly_end(); }

// apply a partition:
 
    template<class RepSize>
    void repartition (			            // old_numbering for *this
        const RepSize&        partition,            // old_ownership
        disarray<T,distributed>& new_disarray,            // new_ownership (created)
        RepSize&              old_numbering,        // new_ownership
        RepSize&              new_numbering) const  // old_ownership
        { return base::data().repartition (partition.data(), new_disarray.data(), old_numbering.data(), new_numbering.data()); }

    template<class RepSize>
    void permutation_apply (                       // old_numbering for *this
        const RepSize&          new_numbering,     // old_ownership
        disarray<T,distributed,A>& new_disarray) const   // new_ownership (already allocated)
        { base::data().permutation_apply (new_numbering.data(), new_disarray.data()); }

    void reverse_permutation (                                 // old_ownership for *this=iold2dis_inew
        disarray<size_type,distributed,A>& inew2dis_iold) const   // new_ownership
        { base::data().reverse_permutation (inew2dis_iold.data()); }

// i/o:

    odiststream& put_values (odiststream& ops) const { return base::data().put_values(ops); }
    idiststream& get_values (idiststream& ips)       { return base::data().get_values(ips); }
    void dump (std::string name) const 	    { return base::data().dump(name); }

    template <class GetFunction>
    idiststream& get_values (idiststream& ips, GetFunction get_element)       { return base::data().get_values(ips, get_element); }
    template <class PutFunction>
    odiststream& put_values (odiststream& ops, PutFunction put_element) const { return base::data().put_values(ops, put_element); }
    template <class PutFunction, class A2> odiststream& permuted_put_values (
		odiststream& ops, const disarray<size_type,distributed,A2>& perm, PutFunction put_element) const 
								     { return base::data().permuted_put_values (ops, perm.data(), put_element); }
};
// [verbatim_disarray]
template <class T, class A>
inline
disarray<T,distributed,A>::disarray (
    	const distributor& ownership,
	const T&           init_val,
        const A&           alloc)
 : base(new_macro(rep(ownership,init_val,alloc)))
{
}
template <class T, class A>
inline
void
disarray<T,distributed,A>::resize (
    	const distributor& ownership,
	const T         &  init_val)
{
  base::data().resize (ownership,init_val);
}
#endif // _RHEOLEF_HAVE_MPI

// -------------------------------------------------------------
// i/o with operator<< & >>
// -------------------------------------------------------------
template <class T, class A>
inline
idiststream&
operator >> (idiststream& ips,  disarray<T,sequential,A>& x)
{ 
    return x.get_values(ips); 
}
template <class T, class A>
inline
odiststream&
operator << (odiststream& ops, const disarray<T,sequential,A>& x)
{
    return x.put_values(ops);
}
#ifdef _RHEOLEF_HAVE_MPI
template <class T, class A>
inline
idiststream&
operator >> (idiststream& ips,  disarray<T,distributed,A>& x)
{ 
    return x.get_values(ips); 
}
template <class T, class A>
inline
odiststream&
operator << (odiststream& ops, const disarray<T,distributed,A>& x)
{
    return x.put_values(ops);
}
#endif // _RHEOLEF_HAVE_MPI
} // namespace rheolef

// -------------------------------------------------------------
// not inlined : longer code
// -------------------------------------------------------------
#include "rheolef/disarray_seq.icc"
#include "rheolef/disarray_mpi.icc"
#endif // _RHEO_DISARRAY_H