12#ifndef DUMUX_GMSH_GRID_DATA_HANDLE_HH 
   13#define DUMUX_GMSH_GRID_DATA_HANDLE_HH 
   19#include <dune/common/parallel/communication.hh> 
   20#include <dune/geometry/dimension.hh> 
   21#include <dune/grid/common/partitionset.hh> 
   22#include <dune/grid/common/datahandleif.hh> 
   26#include <dune/grid/uggrid.hh> 
   35template<
class Gr
id, 
class Gr
idFactory, 
class Data>
 
   36struct GmshGridDataHandle : 
public Dune::CommDataHandleIF<GmshGridDataHandle<Grid, GridFactory, Data>, typename Data::value_type>
 
   40    GmshGridDataHandle(
const Grid& grid, 
const GridFactory& gridFactory, Data& elementMarkers, Data& boundaryMarkers, Data& faceMarkers)
 
   41    : gridView_(grid.levelGridView(0))
 
   42    , idSet_(grid.localIdSet())
 
   43    , elementMarkers_(elementMarkers)
 
   44    , boundaryMarkers_(boundaryMarkers)
 
   45    , faceMarkers_(faceMarkers)
 
   47        const auto& indexSet = gridView_.indexSet();
 
   49        for (
const auto& element : elements(gridView_, Dune::Partitions::interior))
 
   50           std::swap(elementMarkers_[gridFactory.insertionIndex(element)], data_[idSet_.id(element)]);
 
   52        for (
const auto& face : entities(gridView_, Dune::Codim<1>()))
 
   53           std::swap(faceMarkers_[indexSet.index(face)], data_[idSet_.id(face)]);
 
 
   58        const auto& indexSet = gridView_.indexSet();
 
   60        elementMarkers_.resize(indexSet.size(0));
 
   61        for (
const auto& element : elements(gridView_))
 
   62           std::swap(elementMarkers_[indexSet.index(element)], data_[idSet_.id(element)]);
 
   64        faceMarkers_.resize(indexSet.size(1));
 
   65        for (
const auto& face : entities(gridView_, Dune::Codim<1>()))
 
   66           std::swap(faceMarkers_[indexSet.index(face)], data_[idSet_.id(face)]);
 
   68        boundaryMarkers_.resize(gridView_.grid().numBoundarySegments(), 0);
 
   69        for (
const auto& element : elements(gridView_.grid().leafGridView()))
 
   71            for (
const auto& intersection : intersections(gridView_.grid().leafGridView(), element))
 
   73                if (intersection.boundary())
 
   75                    const auto marker = faceMarkers_[indexSet.index(element.template subEntity<1>(intersection.indexInInside()))];
 
   76                    boundaryMarkers_[intersection.boundarySegmentIndex()] = marker;
 
 
   82    Dune::CommDataHandleIF<GmshGridDataHandle<Grid, GridFactory, Data>, 
typename Data::value_type>& 
interface()
 
 
   86    { 
return codim == 0 || codim == 1; }
 
 
   92    template<
class EntityType>
 
   93    std::size_t 
size (
const EntityType& e)
 const 
 
   96    template<
class MessageBufferImp, 
class EntityType>
 
   97    void gather (MessageBufferImp& buff, 
const EntityType& e)
 const 
   98    { buff.write(data_[idSet_.id(e)]); }
 
 
  100    template<
class MessageBufferImp, 
class EntityType>
 
  101    void scatter (MessageBufferImp& buff, 
const EntityType& e, std::size_t n)
 
  102    { buff.read(data_[idSet_.id(e)]); }
 
 
  105    using IdSet = 
typename Grid::LocalIdSet;
 
  109    Data& elementMarkers_;
 
  110    Data& boundaryMarkers_;
 
  112    mutable std::map< typename IdSet::IdType, typename Data::value_type> data_;
 
 
  121template<
class Gr
idFactory, 
class Data, 
int dimgr
id>
 
  122struct GmshGridDataHandle<
Dune::UGGrid<dimgrid>, GridFactory, Data>
 
  123: 
public Dune::CommDataHandleIF<GmshGridDataHandle<Dune::UGGrid<dimgrid>, GridFactory, Data>, typename Data::value_type>
 
  125    using Grid = Dune::UGGrid<dimgrid>;
 
  126    using GridView = 
typename Grid::LevelGridView;
 
  128    GmshGridDataHandle(
const Grid& grid, 
const GridFactory& gridFactory, Data& elementMarkers, Data& boundaryMarkers)
 
  129    : gridView_(grid.levelGridView(0))
 
  130    , idSet_(grid.localIdSet())
 
  131    , elementMarkers_(elementMarkers)
 
  132    , boundaryMarkers_(boundaryMarkers)
 
  134        for (
const auto& element : elements(gridView_, Dune::Partitions::interior))
 
  135           std::swap(elementMarkers_[gridFactory.insertionIndex(element)], data_[idSet_.id(element)]);
 
  142        auto bmSizeMin = boundaryMarkers_.size();
 
  143        Dune::MPIHelper::getCommunication().min(&bmSizeMin, 1);
 
  146            auto bmSize = boundaryMarkers_.size();
 
  147            Dune::MPIHelper::getCommunication().broadcast(&bmSize, 1, 0);
 
  148            boundaryMarkers_.resize(bmSize);
 
  149            Dune::MPIHelper::getCommunication().broadcast(&boundaryMarkers_.front(), bmSize, 0);
 
  155        const auto& indexSet = gridView_.indexSet();
 
  156        elementMarkers_.resize(indexSet.size(0));
 
  157        for (
const auto& element : elements(gridView_))
 
  158           std::swap(elementMarkers_[indexSet.index(element)], data_[idSet_.id(element)]);
 
  161    Dune::CommDataHandleIF<GmshGridDataHandle<Grid, GridFactory, Data>, 
typename Data::value_type>& 
interface()
 
  164    bool contains (
int dim, 
int codim)
 const 
  165    { 
return codim == 0 || codim == 1; }
 
  171    template<
class EntityType>
 
  172    std::size_t 
size (
const EntityType& e)
 const 
  175    template<
class MessageBufferImp, 
class EntityType>
 
  176    void gather (MessageBufferImp& buff, 
const EntityType& e)
 const 
  177    { buff.write(data_[idSet_.id(e)]); }
 
  179    template<
class MessageBufferImp, 
class EntityType>
 
  180    void scatter (MessageBufferImp& buff, 
const EntityType& e, std::size_t n)
 
  181    { buff.read(data_[idSet_.id(e)]); }
 
  184    using IdSet = 
typename Grid::LocalIdSet;
 
  188    Data& elementMarkers_;
 
  189    Data& boundaryMarkers_;
 
  190    mutable std::map< typename IdSet::IdType, typename Data::value_type> data_;
 
Definition gridcapabilities.hh:57
Definition common/pdesolver.hh:24
void gather(MessageBufferImp &buff, const EntityType &e) const
Definition gmshgriddatahandle.hh:97
std::size_t size(const EntityType &e) const
Definition gmshgriddatahandle.hh:93
void scatter(MessageBufferImp &buff, const EntityType &e, std::size_t n)
Definition gmshgriddatahandle.hh:101
GmshGridDataHandle(const Grid &grid, const GridFactory &gridFactory, Data &elementMarkers, Data &boundaryMarkers, Data &faceMarkers)
Definition gmshgriddatahandle.hh:40
bool contains(int dim, int codim) const
Definition gmshgriddatahandle.hh:85
bool fixedSize(int dim, int codim) const
returns true if size per entity of given dim and codim is a constant
Definition gmshgriddatahandle.hh:89
~GmshGridDataHandle()
Definition gmshgriddatahandle.hh:56
Dune::CommDataHandleIF< GmshGridDataHandle< Grid, GridFactory, Data >, typename Data::value_type > & interface()
Definition gmshgriddatahandle.hh:82
typename Grid::LevelGridView GridView
Definition gmshgriddatahandle.hh:38