"Fossies" - the Fresh Open Source Software Archive

Member "dune-istl-2.7.1/dune/istl/owneroverlapcopy.hh" (26 Nov 2020, 24046 Bytes) of package /linux/misc/dune/dune-istl-2.7.1.tar.gz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style: standard) with prefixed line numbers and code folding option. Alternatively you can here view or download the uninterpreted source code file. For more information about "owneroverlapcopy.hh" see the Fossies "Dox" file reference documentation and the last Fossies "Diffs" side-by-side code changes report: 2.6.0_vs_2.7.0.

    1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
    2 // vi: set et ts=4 sw=2 sts=2:
    3 #ifndef DUNE_ISTL_OWNEROVERLAPCOPY_HH
    4 #define DUNE_ISTL_OWNEROVERLAPCOPY_HH
    5 
    6 #include <new>
    7 #include <iostream>
    8 #include <vector>
    9 #include <list>
   10 #include <map>
   11 #include <set>
   12 #include <tuple>
   13 
   14 #include "cmath"
   15 
   16 // MPI header
   17 #if HAVE_MPI
   18 #include <mpi.h>
   19 #endif
   20 
   21 #include <dune/common/enumset.hh>
   22 
   23 #if HAVE_MPI
   24 #include <dune/common/parallel/indexset.hh>
   25 #include <dune/common/parallel/communicator.hh>
   26 #include <dune/common/parallel/remoteindices.hh>
   27 #include <dune/common/parallel/mpicommunication.hh>
   28 #endif
   29 
   30 #include "solvercategory.hh"
   31 #include "istlexception.hh"
   32 #include <dune/common/parallel/communication.hh>
   33 #include <dune/istl/matrixmarket.hh>
   34 
   35 template<int dim, template<class,class> class Comm>
   36 void testRedistributed(int s);
   37 
   38 
   39 namespace Dune {
   40 
   41   /**
   42      @addtogroup ISTL_Comm
   43      @{
   44    */
   45 
   46   /**
   47    * @file
   48    * @brief Classes providing communication interfaces for
   49    * overlapping Schwarz methods.
   50    * @author Peter Bastian
   51    */
   52 
   53   /**
   54    * @brief Attribute set for overlapping Schwarz.
   55    */
   56   struct OwnerOverlapCopyAttributeSet
   57   {
   58     enum AttributeSet {
   59       owner=1, overlap=2, copy=3
   60     };
   61   };
   62 
   63   /**
   64    * @brief Information about the index distribution.
   65    *
   66    * This class contains information about indices local to
   67    * the process together with information about on which
   68    * processes those indices are also present together with the
   69    * attribute they have there.
   70    *
   71    * This information might be used to set up an IndexSet together with
   72    * an RemoteIndices object needed for the ISTL communication classes.
   73    */
   74   template <class G, class L>
   75   class IndexInfoFromGrid
   76   {
   77   public:
   78     /** @brief The type of the global index. */
   79     typedef G GlobalIdType;
   80 
   81     /** @brief The type of the local index. */
   82     typedef L LocalIdType;
   83 
   84     /**
   85      * @brief A triple describing a local index.
   86      *
   87      * The triple consists of the global index and the local
   88      * index and an attribute
   89      */
   90     typedef std::tuple<GlobalIdType,LocalIdType,int> IndexTripel;
   91     /**
   92      * @brief A triple describing a remote index.
   93      *
   94      * The triple consists of a process number and the global index and
   95      * the attribute of the index at the remote process.
   96      */
   97     typedef std::tuple<int,GlobalIdType,int> RemoteIndexTripel;
   98 
   99     /**
  100      * @brief Add a new index triple to the set of local indices.
  101      *
  102      * @param x The index triple.
  103      */
  104     void addLocalIndex (const IndexTripel& x)
  105     {
  106       if (std::get<2>(x)!=OwnerOverlapCopyAttributeSet::owner &&
  107           std::get<2>(x)!=OwnerOverlapCopyAttributeSet::overlap &&
  108           std::get<2>(x)!=OwnerOverlapCopyAttributeSet::copy)
  109         DUNE_THROW(ISTLError,"OwnerOverlapCopyCommunication: global index not in index set");
  110       localindices.insert(x);
  111     }
  112 
  113     /**
  114      * @brief Add a new remote index triple to the set of remote indices.
  115      *
  116      * @param x The index triple to add.
  117      */
  118     void addRemoteIndex (const RemoteIndexTripel& x)
  119     {
  120       if (std::get<2>(x)!=OwnerOverlapCopyAttributeSet::owner &&
  121           std::get<2>(x)!=OwnerOverlapCopyAttributeSet::overlap &&
  122           std::get<2>(x)!=OwnerOverlapCopyAttributeSet::copy)
  123         DUNE_THROW(ISTLError,"OwnerOverlapCopyCommunication: global index not in index set");
  124       remoteindices.insert(x);
  125     }
  126 
  127     /**
  128      * @brief Get the set of indices local to the process.
  129      * @return The set of local indices.
  130      */
  131     const std::set<IndexTripel>& localIndices () const
  132     {
  133       return localindices;
  134     }
  135 
  136     /**
  137      * @brief Get the set of remote indices.
  138      * @return the set of remote indices.
  139      */
  140     const std::set<RemoteIndexTripel>& remoteIndices () const
  141     {
  142       return remoteindices;
  143     }
  144 
  145     /**
  146      * @brief Remove all indices from the sets.
  147      */
  148     void clear ()
  149     {
  150       localindices.clear();
  151       remoteindices.clear();
  152     }
  153 
  154   private:
  155     /** @brief The set of local indices. */
  156     std::set<IndexTripel> localindices;
  157     /** @brief The set of remote indices. */
  158     std::set<RemoteIndexTripel> remoteindices;
  159   };
  160 
  161 
  162 #if HAVE_MPI
  163 
  164   /**
  165    * @brief A class setting up standard communication for a two-valued
  166    * attribute set with owner/overlap/copy semantics.
  167    *
  168    * set up communication from known distribution with owner/overlap/copy semantics
  169    */
  170   template <class GlobalIdType, class LocalIdType=int>
  171   class OwnerOverlapCopyCommunication
  172   {
  173     template<typename M, typename G, typename L>
  174     friend void loadMatrixMarket(M&,
  175                                  const std::string&,
  176                                  OwnerOverlapCopyCommunication<G,L>&,
  177                                  bool);
  178     // used types
  179     typedef typename IndexInfoFromGrid<GlobalIdType,LocalIdType>::IndexTripel IndexTripel;
  180     typedef typename IndexInfoFromGrid<GlobalIdType,LocalIdType>::RemoteIndexTripel RemoteIndexTripel;
  181     typedef typename std::set<IndexTripel>::const_iterator localindex_iterator;
  182     typedef typename std::set<RemoteIndexTripel>::const_iterator remoteindex_iterator;
  183     typedef typename OwnerOverlapCopyAttributeSet::AttributeSet AttributeSet;
  184     typedef Dune::ParallelLocalIndex<AttributeSet> LI;
  185   public:
  186     typedef Dune::ParallelIndexSet<GlobalIdType,LI,512> PIS;
  187     typedef Dune::RemoteIndices<PIS> RI;
  188     typedef Dune::RemoteIndexListModifier<PIS,typename RI::Allocator,false> RILM;
  189     typedef typename RI::RemoteIndex RX;
  190     typedef Dune::BufferedCommunicator BC;
  191     typedef Dune::Interface IF;
  192     typedef EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::owner> OwnerSet;
  193     typedef EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::copy> CopySet;
  194     typedef Combine<EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::owner>,EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::overlap>,AttributeSet> OwnerOverlapSet;
  195     typedef Dune::AllSet<AttributeSet> AllSet;
  196   protected:
  197 
  198 
  199     /** \brief gather/scatter callback for communcation */
  200     template<typename T>
  201     struct CopyGatherScatter
  202     {
  203       typedef typename CommPolicy<T>::IndexedType V;
  204 
  205       static V gather(const T& a, std::size_t i)
  206       {
  207         return a[i];
  208       }
  209 
  210       static void scatter(T& a, V v, std::size_t i)
  211       {
  212         a[i] = v;
  213       }
  214     };
  215     template<typename T>
  216     struct AddGatherScatter
  217     {
  218       typedef typename CommPolicy<T>::IndexedType V;
  219 
  220       static V gather(const T& a, std::size_t i)
  221       {
  222         return a[i];
  223       }
  224 
  225       static void scatter(T& a, V v, std::size_t i)
  226       {
  227         a[i] += v;
  228       }
  229     };
  230 
  231     void buildOwnerOverlapToAllInterface () const
  232     {
  233       if (OwnerOverlapToAllInterfaceBuilt)
  234         OwnerOverlapToAllInterface.free();
  235       typedef Combine<EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::owner>,EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::overlap>,AttributeSet> OwnerOverlapSet;
  236       typedef Combine<OwnerOverlapSet,EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::copy>,AttributeSet> AllSet;
  237       OwnerOverlapSet sourceFlags;
  238       AllSet destFlags;
  239       OwnerOverlapToAllInterface.build(ri,sourceFlags,destFlags);
  240       OwnerOverlapToAllInterfaceBuilt = true;
  241     }
  242 
  243     void buildOwnerToAllInterface () const
  244     {
  245       if (OwnerToAllInterfaceBuilt)
  246         OwnerToAllInterface.free();
  247       OwnerSet sourceFlags;
  248       AllSet destFlags;
  249       OwnerToAllInterface.build(ri,sourceFlags,destFlags);
  250       OwnerToAllInterfaceBuilt = true;
  251     }
  252 
  253     void buildOwnerCopyToAllInterface () const
  254     {
  255       if (OwnerCopyToAllInterfaceBuilt)
  256         OwnerCopyToAllInterface.free();
  257       typedef Combine<EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::owner>,EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::copy>,AttributeSet> OwnerCopySet;
  258       typedef Combine<OwnerCopySet,EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::overlap>,AttributeSet> AllSet;
  259       OwnerCopySet sourceFlags;
  260       AllSet destFlags;
  261       OwnerCopyToAllInterface.build(ri,sourceFlags,destFlags);
  262       OwnerCopyToAllInterfaceBuilt = true;
  263     }
  264 
  265     void buildOwnerCopyToOwnerCopyInterface () const
  266     {
  267       if (OwnerCopyToOwnerCopyInterfaceBuilt)
  268         OwnerCopyToOwnerCopyInterface.free();
  269       typedef Combine<EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::owner>,EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::copy>,AttributeSet> OwnerCopySet;
  270       OwnerCopySet sourceFlags;
  271       OwnerCopySet destFlags;
  272       OwnerCopyToOwnerCopyInterface.build(ri,sourceFlags,destFlags);
  273       OwnerCopyToOwnerCopyInterfaceBuilt = true;
  274     }
  275 
  276     void buildCopyToAllInterface () const
  277     {
  278       if (CopyToAllInterfaceBuilt)
  279         CopyToAllInterface.free();
  280       CopySet sourceFlags;
  281       AllSet destFlags;
  282       CopyToAllInterface.build(ri,sourceFlags,destFlags);
  283       CopyToAllInterfaceBuilt = true;
  284     }
  285 
  286   public:
  287 
  288     /**
  289      * @brief Set right Solver Category (default is overlapping).
  290      */
  291     void
  292     DUNE_DEPRECATED_MSG("The solver category can only be set in the constructor. This method is deprecated and will be removed after Dune 2.7")
  293     setSolverCategory (SolverCategory::Category set) {
  294       category_ = set;
  295     }
  296 
  297     SolverCategory::Category
  298     DUNE_DEPRECATED_MSG("This method is deprecated and will be removed after Dune 2.7, use category() instead.")
  299     getSolverCategory () const {
  300       return category_;
  301     }
  302 
  303     /**
  304      * @brief Get Solver Category.
  305      * @return The Solver Category.
  306      */
  307     SolverCategory::Category category () const {
  308       return category_;
  309     }
  310 
  311     const CollectiveCommunication<MPI_Comm>& communicator() const
  312     {
  313       return cc;
  314     }
  315 
  316     /**
  317      * @brief Communicate values from owner data points to all other data points.
  318      *
  319      * @brief source The data to send from.
  320      * @brief dest The data to send to.
  321      */
  322     template<class T>
  323     void copyOwnerToAll (const T& source, T& dest) const
  324     {
  325       if (!OwnerToAllInterfaceBuilt)
  326         buildOwnerToAllInterface ();
  327       BC communicator;
  328       communicator.template build<T>(OwnerToAllInterface);
  329       communicator.template forward<CopyGatherScatter<T> >(source,dest);
  330       communicator.free();
  331     }
  332 
  333     /**
  334      * @brief Communicate values from copy data points to all other data points.
  335      *
  336      * @brief source The data to send from.
  337      * @brief dest The data to send to.
  338      */
  339     template<class T>
  340     void copyCopyToAll (const T& source, T& dest) const
  341     {
  342       if (!CopyToAllInterfaceBuilt)
  343         buildCopyToAllInterface ();
  344       BC communicator;
  345       communicator.template build<T>(CopyToAllInterface);
  346       communicator.template forward<CopyGatherScatter<T> >(source,dest);
  347       communicator.free();
  348     }
  349 
  350     /**
  351      * @brief Communicate values from owner data points to all other data points and add them to those values.
  352      *
  353      * @brief source The data to send from.
  354      * @brief dest The data to add them communicated values to.
  355      */
  356     template<class T>
  357     void addOwnerOverlapToAll (const T& source, T& dest) const
  358     {
  359       if (!OwnerOverlapToAllInterfaceBuilt)
  360         buildOwnerOverlapToAllInterface ();
  361       BC communicator;
  362       communicator.template build<T>(OwnerOverlapToAllInterface);
  363       communicator.template forward<AddGatherScatter<T> >(source,dest);
  364       communicator.free();
  365     }
  366 
  367     /**
  368      * @brief Communicate values from owner and copy data points to all other data points and add them to those values.
  369      *
  370      * @brief source The data to send from.
  371      * @brief dest The data to add them communicated values to.
  372      */
  373     template<class T>
  374     void addOwnerCopyToAll (const T& source, T& dest) const
  375     {
  376       if (!OwnerCopyToAllInterfaceBuilt)
  377         buildOwnerCopyToAllInterface ();
  378       BC communicator;
  379       communicator.template build<T>(OwnerCopyToAllInterface);
  380       communicator.template forward<AddGatherScatter<T> >(source,dest);
  381       communicator.free();
  382     }
  383 
  384     /**
  385      * @brief Communicate values from owner and copy data points to owner and copy data points and add them to those values.
  386      *
  387      * @brief source The data to send from.
  388      * @brief dest The data to add the communicated values to.
  389      */
  390     template<class T>
  391     void addOwnerCopyToOwnerCopy (const T& source, T& dest) const
  392     {
  393       if (!OwnerCopyToOwnerCopyInterfaceBuilt)
  394         buildOwnerCopyToOwnerCopyInterface ();
  395       BC communicator;
  396       communicator.template build<T>(OwnerCopyToOwnerCopyInterface);
  397       communicator.template forward<AddGatherScatter<T> >(source,dest);
  398       communicator.free();
  399     }
  400 
  401 
  402     /**
  403      * @brief Compute a global dot product of two vectors.
  404      *
  405      * @param x The first vector of the product.
  406      * @param y The second vector of the product.
  407      * @param result Reference to store the result in.
  408      */
  409     template<class T1, class T2>
  410     void dot (const T1& x, const T1& y, T2& result) const
  411     {
  412       // set up mask vector
  413       if (mask.size()!=static_cast<typename std::vector<double>::size_type>(x.size()))
  414       {
  415         mask.resize(x.size());
  416         for (typename std::vector<double>::size_type i=0; i<mask.size(); i++)
  417           mask[i] = 1;
  418         for (typename PIS::const_iterator i=pis.begin(); i!=pis.end(); ++i)
  419           if (i->local().attribute()!=OwnerOverlapCopyAttributeSet::owner)
  420             mask[i->local().local()] = 0;
  421       }
  422       result = T2(0.0);
  423 
  424       for (typename T1::size_type i=0; i<x.size(); i++)
  425         result += x[i]*(y[i])*mask[i];
  426       result = cc.sum(result);
  427     }
  428 
  429     /**
  430      * @brief Compute the global Euclidean norm of a vector.
  431      *
  432      * @param x The vector to compute the norm of.
  433      * @return The global Euclidean norm of that vector.
  434      */
  435     template<class T1>
  436     typename FieldTraits<typename T1::field_type>::real_type norm (const T1& x) const
  437     {
  438       using real_type = typename FieldTraits<typename T1::field_type>::real_type;
  439 
  440       // set up mask vector
  441       if (mask.size()!=static_cast<typename std::vector<double>::size_type>(x.size()))
  442       {
  443         mask.resize(x.size());
  444         for (typename std::vector<double>::size_type i=0; i<mask.size(); i++)
  445           mask[i] = 1;
  446         for (typename PIS::const_iterator i=pis.begin(); i!=pis.end(); ++i)
  447           if (i->local().attribute()!=OwnerOverlapCopyAttributeSet::owner)
  448             mask[i->local().local()] = 0;
  449       }
  450       auto result = real_type(0.0);
  451       for (typename T1::size_type i=0; i<x.size(); i++)
  452         result += Impl::asVector(x[i]).two_norm2()*mask[i];
  453       return sqrt(cc.sum(result));
  454     }
  455 
  456     typedef Dune::EnumItem<AttributeSet,OwnerOverlapCopyAttributeSet::copy> CopyFlags;
  457 
  458     /** @brief The type of the parallel index set. */
  459     typedef Dune::ParallelIndexSet<GlobalIdType,LI,512> ParallelIndexSet;
  460 
  461     /** @brief The type of the remote indices. */
  462     typedef Dune::RemoteIndices<PIS> RemoteIndices;
  463 
  464     /**
  465      * @brief The type of the reverse lookup of indices. */
  466     typedef Dune::GlobalLookupIndexSet<ParallelIndexSet> GlobalLookupIndexSet;
  467 
  468     /**
  469      * @brief Get the underlying parallel index set.
  470      * @return The underlying parallel index set.
  471      */
  472     const ParallelIndexSet& indexSet() const
  473     {
  474       return pis;
  475     }
  476 
  477     /**
  478      * @brief Get the underlying remote indices.
  479      * @return The underlying remote indices.
  480      */
  481     const RemoteIndices& remoteIndices() const
  482     {
  483       return ri;
  484     }
  485 
  486     /**
  487      * @brief Get the underlying parallel index set.
  488      * @return The underlying parallel index set.
  489      */
  490     ParallelIndexSet& indexSet()
  491     {
  492       return pis;
  493     }
  494 
  495 
  496     /**
  497      * @brief Get the underlying remote indices.
  498      * @return The underlying remote indices.
  499      */
  500     RemoteIndices& remoteIndices()
  501     {
  502       return ri;
  503     }
  504 
  505     void buildGlobalLookup()
  506     {
  507       if(globalLookup_) {
  508         if(pis.seqNo()==oldseqNo)
  509           // Nothing changed!
  510           return;
  511         delete globalLookup_;
  512       }
  513 
  514       globalLookup_ = new GlobalLookupIndexSet(pis);
  515       oldseqNo = pis.seqNo();
  516     }
  517 
  518     void buildGlobalLookup(std::size_t size)
  519     {
  520       if(globalLookup_) {
  521         if(pis.seqNo()==oldseqNo)
  522           // Nothing changed!
  523           return;
  524         delete globalLookup_;
  525       }
  526       globalLookup_ = new GlobalLookupIndexSet(pis, size);
  527       oldseqNo = pis.seqNo();
  528     }
  529 
  530     void freeGlobalLookup()
  531     {
  532       delete globalLookup_;
  533       globalLookup_=0;
  534     }
  535 
  536     const GlobalLookupIndexSet& globalLookup() const
  537     {
  538       assert(globalLookup_ != 0);
  539       return *globalLookup_;
  540     }
  541 
  542     /**
  543      * @brief Set vector to zero at copy dofs
  544      *
  545      * @param x The vector to project.
  546      */
  547     template<class T1>
  548     void project (T1& x) const
  549     {
  550       for (typename PIS::const_iterator i=pis.begin(); i!=pis.end(); ++i)
  551         if (i->local().attribute()==OwnerOverlapCopyAttributeSet::copy)
  552           x[i->local().local()] = 0;
  553     }
  554 
  555     /**
  556      * @brief Construct the communication without any indices.
  557      *
  558      * The local index set and the remote indices have to be set up
  559      * later on.
  560      * @param comm_ The MPI Communicator to use, e. g. MPI_COMM_WORLD
  561      * @param cat_ The Solver category, default is overlapping
  562      * @param freecomm_ Whether to free the communicator comm_ in the destructor, default is false
  563      */
  564     OwnerOverlapCopyCommunication (MPI_Comm comm_,
  565                                    SolverCategory::Category cat_ = SolverCategory::overlapping,
  566                                    bool freecomm_ = false)
  567       : comm(comm_), cc(comm_), pis(), ri(pis,pis,comm_),
  568         OwnerToAllInterfaceBuilt(false), OwnerOverlapToAllInterfaceBuilt(false),
  569         OwnerCopyToAllInterfaceBuilt(false), OwnerCopyToOwnerCopyInterfaceBuilt(false),
  570         CopyToAllInterfaceBuilt(false), globalLookup_(0), category_(cat_),
  571         freecomm(freecomm_)
  572     {}
  573 
  574     /**
  575      * @brief Construct the communication without any indices using MPI_COMM_WORLD.
  576      *
  577      * The local index set and the remote indices have to be set up
  578      * later on.
  579      * @param cat_ The Solver category, default is overlapping
  580        is false
  581      */
  582     OwnerOverlapCopyCommunication (SolverCategory::Category cat_ = SolverCategory::overlapping)
  583       : comm(MPI_COMM_WORLD), cc(MPI_COMM_WORLD), pis(), ri(pis,pis,MPI_COMM_WORLD),
  584         OwnerToAllInterfaceBuilt(false), OwnerOverlapToAllInterfaceBuilt(false),
  585         OwnerCopyToAllInterfaceBuilt(false), OwnerCopyToOwnerCopyInterfaceBuilt(false),
  586         CopyToAllInterfaceBuilt(false), globalLookup_(0), category_(cat_), freecomm(false)
  587     {}
  588 
  589     /**
  590      * @brief Constructor
  591      * @param indexinfo The set of IndexTripels describing the local and remote indices.
  592      * @param comm_ The communicator to use in the communication.
  593      * @param cat_ The Solver category, default is overlapping
  594      * @param freecomm_ Whether to free the communicator comm_ in the destructor, default is false
  595      */
  596     OwnerOverlapCopyCommunication (const IndexInfoFromGrid<GlobalIdType, LocalIdType>& indexinfo,
  597                                    MPI_Comm comm_,
  598                                    SolverCategory::Category cat_ = SolverCategory::overlapping,
  599                                    bool freecomm_ = false)
  600       : comm(comm_), cc(comm_), OwnerToAllInterfaceBuilt(false),
  601         OwnerOverlapToAllInterfaceBuilt(false), OwnerCopyToAllInterfaceBuilt(false),
  602         OwnerCopyToOwnerCopyInterfaceBuilt(false), CopyToAllInterfaceBuilt(false),
  603         globalLookup_(0), category_(cat_), freecomm(freecomm_)
  604     {
  605       // set up an ISTL index set
  606       pis.beginResize();
  607       for (localindex_iterator i=indexinfo.localIndices().begin(); i!=indexinfo.localIndices().end(); ++i)
  608       {
  609         if (std::get<2>(*i)==OwnerOverlapCopyAttributeSet::owner)
  610           pis.add(std::get<0>(*i),LI(std::get<1>(*i),OwnerOverlapCopyAttributeSet::owner,true));
  611         if (std::get<2>(*i)==OwnerOverlapCopyAttributeSet::overlap)
  612           pis.add(std::get<0>(*i),LI(std::get<1>(*i),OwnerOverlapCopyAttributeSet::overlap,true));
  613         if (std::get<2>(*i)==OwnerOverlapCopyAttributeSet::copy)
  614           pis.add(std::get<0>(*i),LI(std::get<1>(*i),OwnerOverlapCopyAttributeSet::copy,true));
  615         //                std::cout << cc.rank() << ": adding index " << std::get<0>(*i) << " " << std::get<1>(*i) << " " << std::get<2>(*i) << std::endl;
  616       }
  617       pis.endResize();
  618 
  619       // build remote indices WITHOUT communication
  620       //          std::cout << cc.rank() << ": build remote indices" << std::endl;
  621       ri.setIndexSets(pis,pis,cc);
  622       if (indexinfo.remoteIndices().size()>0)
  623       {
  624         remoteindex_iterator i=indexinfo.remoteIndices().begin();
  625         int p = std::get<0>(*i);
  626         RILM modifier = ri.template getModifier<false,true>(p);
  627         typename PIS::const_iterator pi=pis.begin();
  628         for ( ; i!=indexinfo.remoteIndices().end(); ++i)
  629         {
  630           // handle processor change
  631           if (p!=std::get<0>(*i))
  632           {
  633             p = std::get<0>(*i);
  634             modifier = ri.template getModifier<false,true>(p);
  635             pi=pis.begin();
  636           }
  637 
  638           // position to correct entry in parallel index set
  639           while (pi->global()!=std::get<1>(*i) && pi!=pis.end())
  640             ++pi;
  641           if (pi==pis.end())
  642             DUNE_THROW(ISTLError,"OwnerOverlapCopyCommunication: global index not in index set");
  643 
  644           // insert entry
  645           //                      std::cout << cc.rank() << ": adding remote index " << std::get<0>(*i) << " " << std::get<1>(*i) << " " << std::get<2>(*i) << std::endl;
  646           if (std::get<2>(*i)==OwnerOverlapCopyAttributeSet::owner)
  647             modifier.insert(RX(OwnerOverlapCopyAttributeSet::owner,&(*pi)));
  648           if (std::get<2>(*i)==OwnerOverlapCopyAttributeSet::overlap)
  649             modifier.insert(RX(OwnerOverlapCopyAttributeSet::overlap,&(*pi)));
  650           if (std::get<2>(*i)==OwnerOverlapCopyAttributeSet::copy)
  651             modifier.insert(RX(OwnerOverlapCopyAttributeSet::copy,&(*pi)));
  652         }
  653       }else{
  654         // Force remote indices to be synced!
  655         ri.template getModifier<false,true>(0);
  656       }
  657     }
  658 
  659     // destructor: free memory in some objects
  660     ~OwnerOverlapCopyCommunication ()
  661     {
  662       ri.free();
  663       if (OwnerToAllInterfaceBuilt) OwnerToAllInterface.free();
  664       if (OwnerOverlapToAllInterfaceBuilt) OwnerOverlapToAllInterface.free();
  665       if (OwnerCopyToAllInterfaceBuilt) OwnerCopyToAllInterface.free();
  666       if (OwnerCopyToOwnerCopyInterfaceBuilt) OwnerCopyToOwnerCopyInterface.free();
  667       if (CopyToAllInterfaceBuilt) CopyToAllInterface.free();
  668       if (globalLookup_) delete globalLookup_;
  669       if (freecomm==true)
  670         if(comm!=MPI_COMM_NULL)
  671         {
  672 #ifdef MPI_2
  673           // If it is possible to query whether MPI_Finalize
  674           // was called, only free the communicator before
  675           // calling MPI_Finalize.
  676           int wasFinalized = 0;
  677           MPI_Finalized( &wasFinalized );
  678           if(!wasFinalized)
  679 #endif
  680           MPI_Comm_free(&comm);
  681         }
  682     }
  683 
  684   private:
  685     OwnerOverlapCopyCommunication (const OwnerOverlapCopyCommunication&)
  686     {}
  687     MPI_Comm comm;
  688     CollectiveCommunication<MPI_Comm> cc;
  689     PIS pis;
  690     RI ri;
  691     mutable IF OwnerToAllInterface;
  692     mutable bool OwnerToAllInterfaceBuilt;
  693     mutable IF OwnerOverlapToAllInterface;
  694     mutable bool OwnerOverlapToAllInterfaceBuilt;
  695     mutable IF OwnerCopyToAllInterface;
  696     mutable bool OwnerCopyToAllInterfaceBuilt;
  697     mutable IF OwnerCopyToOwnerCopyInterface;
  698     mutable bool OwnerCopyToOwnerCopyInterfaceBuilt;
  699     mutable IF CopyToAllInterface;
  700     mutable bool CopyToAllInterfaceBuilt;
  701     mutable std::vector<double> mask;
  702     int oldseqNo;
  703     GlobalLookupIndexSet* globalLookup_;
  704     // re-intruduce const qualifier once deprecated setCategory got removed after Dune 2.7
  705     /* const */ SolverCategory::Category category_;
  706     bool freecomm;
  707   };
  708 
  709 #endif
  710 
  711 
  712   /** @} end documentation */
  713 
  714 } // end namespace
  715 
  716 #endif