12#ifndef DUMUX_IO_GRID_DATA_HH 
   13#define DUMUX_IO_GRID_DATA_HH 
   19#include <dune/common/exceptions.hh> 
   20#include <dune/common/parallel/communication.hh> 
   21#include <dune/common/parallel/mpihelper.hh> 
   22#include <dune/grid/common/gridfactory.hh> 
   23#include <dune/grid/io/file/dgfparser/parser.hh> 
   24#include <dune/grid/io/file/dgfparser/gridptr.hh> 
   29#include <dune/grid/uggrid.hh> 
   40struct isUG : 
public std::false_type {};
 
   44struct isUG<
Dune::UGGrid<dim>> : 
public std::true_type {};
 
   56    using Intersection = 
typename Grid::LeafIntersection;
 
   57    using Element = 
typename Grid::template Codim<0>::Entity;
 
   58    using Vertex = 
typename Grid::template Codim<Grid::dimension>::Entity;
 
   62    enum DataSourceType { dgf, gmsh, vtk };
 
   66    GridData(std::shared_ptr<Grid> grid, std::shared_ptr<Dune::GridFactory<Grid>> factory,
 
   67             std::vector<int>&& elementMarkers, std::vector<int>&& boundaryMarkers, std::vector<int>&& faceMarkers = std::vector<int>{})
 
   69    , gridFactory_(factory)
 
   70    , elementMarkers_(elementMarkers)
 
   71    , boundaryMarkers_(boundaryMarkers)
 
   72    , faceMarkers_(faceMarkers)
 
   73    , dataSourceType_(DataSourceType::gmsh)
 
 
   79    , dataSourceType_(DataSourceType::dgf)
 
 
   83    GridData(std::shared_ptr<Grid> grid, std::shared_ptr<Dune::GridFactory<Grid>> factory,
 
   86    , gridFactory_(factory)
 
   88    , pointData_(pointData)
 
   89    , dataSourceType_(DataSourceType::vtk)
 
 
  102    const std::vector<double>& 
parameters(
const Vertex& vertex)
 const 
  104        if (dataSourceType_ == DataSourceType::dgf)
 
  106            if (vertex.level() != 0)
 
  107                DUNE_THROW(Dune::IOError, 
"You can only obtain parameters for level 0 vertices!");
 
  109            return dgfGrid_.parameters(vertex);
 
  112            DUNE_THROW(Dune::InvalidStateException, 
"The parameters method is only available if the grid was constructed with a DGF file.");
 
 
  118    const std::vector<double>& 
parameters(
const Element& element)
 const 
  120        if (dataSourceType_ == DataSourceType::dgf)
 
  122            if (element.hasFather())
 
  124                auto level0Element = element;
 
  125                while(level0Element.hasFather())
 
  126                    level0Element = level0Element.father();
 
  128                return dgfGrid_.parameters(level0Element);
 
  132                return dgfGrid_.parameters(element);
 
  136            DUNE_THROW(Dune::InvalidStateException, 
"The parameters method is only available if the grid was constructed with a DGF file.");
 
 
  142    template <
class Gr
idImp, 
class IntersectionImp>
 
  143    const Dune::DGFBoundaryParameter::type& 
parameters(
const Dune::Intersection<GridImp, IntersectionImp>& intersection)
 const 
  145        if (dataSourceType_ == DataSourceType::dgf)
 
  146            return dgfGrid_.parameters(intersection);
 
  148            DUNE_THROW(Dune::InvalidStateException, 
"The parameters method is only available if the grid was constructed with a DGF file.");
 
 
  165        if (dataSourceType_ != DataSourceType::gmsh)
 
  166            DUNE_THROW(Dune::InvalidStateException, 
"Domain markers are only available for gmsh grids.");
 
  167        if (boundarySegmentIndex >= boundaryMarkers_.size())
 
  168            DUNE_THROW(Dune::RangeError, 
"Boundary segment index "<< boundarySegmentIndex << 
" bigger than number of boundary segments in grid.\n" 
  169                                         "Make sure to call this function only for boundaries that were defined as physical entities in gmsh.");
 
  170        return boundaryMarkers_[boundarySegmentIndex];
 
 
  185    { 
return gridFactory_->wasInserted(intersection); }
 
 
  194        if (dataSourceType_ != DataSourceType::gmsh)
 
  195            DUNE_THROW(Dune::InvalidStateException, 
"Domain markers are only available for gmsh grids.");
 
  198        auto level0element = element;
 
  199        while (level0element.hasFather())
 
  200            level0element = level0element.father();
 
  205            return elementMarkers_[gridPtr_->levelGridView(0).indexSet().index(level0element)];
 
  207            return elementMarkers_[gridFactory_->insertionIndex(level0element)];
 
 
  214    template<bool ug = Detail::isUG<Grid>::value, 
typename std::enable_if_t<!ug, int> = 0>
 
  217        return GmshDataHandle(*gridPtr_, *gridFactory_, elementMarkers_, boundaryMarkers_, faceMarkers_);
 
 
  224    template<bool ug = Detail::isUG<Grid>::value, 
typename std::enable_if_t<ug, int> = 0>
 
  227        return GmshDataHandle(*gridPtr_, *gridFactory_, elementMarkers_, boundaryMarkers_);
 
 
  242    double getParameter(
const Element& element, 
const std::string& fieldName)
 const 
 
  252    template<
class T, std::
size_t size>
 
  253    std::array<T, size> 
getArrayParameter(
const Element& element, 
const std::string& fieldName)
 const 
  255        if (dataSourceType_ != DataSourceType::vtk)
 
  256            DUNE_THROW(Dune::InvalidStateException, 
"This access function is only available for data from VTK files.");
 
  258        if (cellData_.count(fieldName) == 0)
 
  259            DUNE_THROW(Dune::IOError, 
"No field with the name " << fieldName << 
" found in cell data.");
 
  262        auto level0element = element;
 
  263        while (level0element.hasFather())
 
  264            level0element = level0element.father();
 
  267        const auto index = [&]{
 
  268            if (gridPtr_->comm().size() > 1)
 
  269                return gridPtr_->levelGridView(0).indexSet().index(level0element);
 
  271                return gridFactory_->insertionIndex(level0element);
 
  274        std::array<T, size> param;
 
  275        const auto& data = cellData_.at(fieldName);
 
  277        if (
const std::size_t nc = data.size()/gridPtr_->levelGridView(0).size(0); size != nc)
 
  278            DUNE_THROW(Dune::IOError, 
"Array size " << size << 
" does not match number of components of field " << nc);
 
  280        for (std::size_t i = 0; i < size; ++i)
 
  281            param[i] = 
static_cast<T
>(data[i + size*index]);
 
 
  291    double getParameter(
const Vertex& vertex, 
const std::string& fieldName)
 const 
 
  302    template<
class T, std::
size_t size>
 
  305        if (dataSourceType_ != DataSourceType::vtk)
 
  306            DUNE_THROW(Dune::InvalidStateException, 
"This access function is only available for data from VTK files.");
 
  308        if (vertex.level() != 0)
 
  309            DUNE_THROW(Dune::IOError, 
"You can only obtain parameters for level 0 vertices!");
 
  311        if (pointData_.count(fieldName) == 0)
 
  312            DUNE_THROW(Dune::IOError, 
"No field with the name " << fieldName << 
" found in point data");
 
  315        const auto index = [&]{
 
  316            if (gridPtr_->comm().size() > 1)
 
  317                return gridPtr_->levelGridView(0).indexSet().index(vertex);
 
  319                return gridFactory_->insertionIndex(vertex);
 
  322        std::array<T, size> param;
 
  323        const auto& data = pointData_.at(fieldName);
 
  325        if (
const std::size_t nc = data.size()/gridPtr_->levelGridView(0).size(Grid::dimension); size != nc)
 
  326            DUNE_THROW(Dune::IOError, 
"Array size " << size << 
" does not match number of components of field " << nc);
 
  328        for (std::size_t i = 0; i < size; ++i)
 
  329            param[i] = 
static_cast<T
>(data[i + size*index]);
 
 
  338        return VtkDataHandle(*gridPtr_, *gridFactory_, cellData_, pointData_);
 
 
  345    std::shared_ptr<Grid> gridPtr_;
 
  346    std::shared_ptr<Dune::GridFactory<Grid>> gridFactory_;
 
  352    std::vector<int> elementMarkers_;
 
  353    std::vector<int> boundaryMarkers_;
 
  354    std::vector<int> faceMarkers_;
 
  362    Dune::GridPtr<Grid> dgfGrid_;
 
  366    DataSourceType dataSourceType_;
 
 
GridData(std::shared_ptr< Grid > grid, std::shared_ptr< Dune::GridFactory< Grid > > factory, std::vector< int > &&elementMarkers, std::vector< int > &&boundaryMarkers, std::vector< int > &&faceMarkers=std::vector< int >{})
constructor for gmsh grid data
Definition griddata.hh:66
double getParameter(const Vertex &vertex, const std::string &fieldName) const
Call the parameters function of the DGF grid pointer if available for vertex data.
Definition griddata.hh:291
int getBoundaryDomainMarker(const Intersection &intersection) const
Return the boundary domain marker (Gmsh physical entity number) of an intersection Only available whe...
Definition griddata.hh:178
double getParameter(const Element &element, const std::string &fieldName) const
Get an element parameter.
Definition griddata.hh:242
GridData(std::shared_ptr< Grid > grid, std::shared_ptr< Dune::GridFactory< Grid > > factory, VTKReader::Data &&cellData, VTKReader::Data &&pointData)
constructor for VTK grid data
Definition griddata.hh:83
const Dune::DGFBoundaryParameter::type & parameters(const Dune::Intersection< GridImp, IntersectionImp > &intersection) const
Call the parameters function of the DGF grid pointer if available.
Definition griddata.hh:143
GridData(Dune::GridPtr< Grid > grid)
constructor for dgf grid data
Definition griddata.hh:77
VtkDataHandle createVtkDataHandle()
Create a data handle for communication of the data in parallel simulations.
Definition griddata.hh:336
const std::vector< double > & parameters(const Vertex &vertex) const
Call the parameters function of the DGF grid pointer if available for vertex data.
Definition griddata.hh:102
std::array< T, size > getArrayParameter(const Vertex &vertex, const std::string &fieldName) const
Call the parameters function of the DGF grid pointer if available for vertex data.
Definition griddata.hh:303
int getBoundaryDomainMarker(int boundarySegmentIndex) const
Return the boundary domain marker (Gmsh physical entity number) of an intersection Only available whe...
Definition griddata.hh:163
GmshDataHandle createGmshDataHandle()
Create a data handle for communication of the data in parallel simulations.
Definition griddata.hh:215
std::array< T, size > getArrayParameter(const Element &element, const std::string &fieldName) const
Get an element parameter that is an array.
Definition griddata.hh:253
int getElementDomainMarker(const Element &element) const
Return the element domain marker (Gmsh physical entity number) of an element. Only available when usi...
Definition griddata.hh:192
const std::vector< double > & parameters(const Element &element) const
Call the parameters function of the DGF grid pointer if available for element data.
Definition griddata.hh:118
bool wasInserted(const Intersection &intersection) const
Returns true if an intersection was inserted during grid creation.
Definition griddata.hh:184
std::unordered_map< std::string, std::vector< double > > Data
the cell / point data type for point data read from a grid file
Definition vtkreader.hh:277
A data handle for commucating grid data for gmsh grids.
Distance implementation details.
Definition cvfelocalresidual.hh:25
Definition common/pdesolver.hh:24
Definition griddata.hh:40
A data handle for commucating grid data for gmsh grids.
Definition gmshgriddatahandle.hh:37
A data handle for communicating grid data for VTK grids.
Definition vtkgriddatahandle.hh:36
A data handle for communicating grid data for VTK grids.
A vtk file reader using tinyxml2 as xml backend.