5#ifndef DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
6#define DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
15#include <dune/common/exceptions.hh>
16#include <dune/common/fvector.hh>
17#include <dune/common/parallel/mpihelper.hh>
95 dgf_( rank( comm ), size( comm ) )
105 dgf_( rank( comm ), size( comm ) )
107 std::ifstream input( filename.c_str() );
109 DUNE_THROW(
DGFException,
"Error: Macrofile " << filename <<
" not found" );
120 template<
class GG,
class II >
123 return factory_.wasInserted( intersection );
127 template<
class GG,
class II >
134 template<
int codim >
138 return dgf_.nofelparams;
140 return dgf_.nofvtxparams;
146 template<
class Entity >
149 return numParameters< Entity::codimension >();
153 std::vector< double > &
parameter (
const typename Grid::template Codim< 0 >::Entity &element )
155 if( numParameters< 0 >() <= 0 )
157 DUNE_THROW( InvalidStateException,
158 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
160 return dgf_.elParams[ factory_.insertionIndex( element ) ];
164 std::vector< double > &
parameter (
const typename Grid::template Codim< dimension >::Entity &
vertex )
166 if( numParameters< dimension >() <= 0 )
168 DUNE_THROW( InvalidStateException,
169 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
171 return dgf_.vtxParams[ factory_.insertionIndex(
vertex ) ];
177 return dgf_.haveBndParameters;
181 template<
class GG,
class II >
188 auto refElem = referenceElement< double, dimension >( entity.type() );
189 int corners = refElem.size( face, 1,
dimension );
190 std::vector< unsigned int > bound( corners );
191 for(
int i = 0; i < corners; ++i )
193 const int k = refElem.subEntity( face, 1, i,
dimension );
194 bound[ i ] = factory_.insertionIndex( entity.template subEntity< dimension >( k ) );
197 DuneGridFormatParser::facemap_t::key_type key( bound,
false );
198 const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
199 if( pos != dgf_.facemap.end() )
200 return dgf_.facemap.find( key )->second.second;
207 void generate ( std::istream &input );
214 MPI_Comm_rank( MPICOMM, &rank );
224 MPI_Comm_size( MPICOMM, &size );
230 GridFactory< UGGrid< dim > > factory_;
231 DuneGridFormatParser dgf_;
Include standard header files.
Definition: agrid.hh:60
@ vertex
Definition: common.hh:133
Definition: dgfgridfactory.hh:38
MPIHelper::MPICommunicator MPICommunicatorType
Definition: dgfgridfactory.hh:41
G Grid
Definition: dgfgridfactory.hh:39
static const int dimension
Definition: dgfgridfactory.hh:40
Intersection of a mesh entity of codimension 0 ("element") with a "neighboring" element or with the d...
Definition: common/intersection.hh:164
int indexInInside() const
Local index of codim 1 entity in the inside() entity where intersection is contained in.
Definition: common/intersection.hh:346
size_t boundarySegmentIndex() const
index of the boundary segment within the macro grid
Definition: common/intersection.hh:236
Entity inside() const
return Entity on the inside of this intersection. That is the Entity where we started this.
Definition: common/intersection.hh:250
GridImp::template Codim< 0 >::Entity Entity
Type of entity that this Intersection belongs to.
Definition: common/intersection.hh:192
Wrapper class for entities.
Definition: common/entity.hh:66
Common Grid parameters.
Definition: gridparameter.hh:35
exception class for IO errors in the DGF parser
Definition: dgfexception.hh:16
Some simple static information for a given GridType.
Definition: io/file/dgfparser/dgfparser.hh:56
size_t heapSize() const
returns heap size used on construction of the grid
Definition: dgfug.hh:48
bool noCopy_
Definition: dgfug.hh:52
UGGridParameterBlock(std::istream &input)
constructor taking istream
Definition: dgfug.cc:20
size_t heapSize_
Definition: dgfug.hh:53
bool noClosure_
Definition: dgfug.hh:51
bool noCopy() const
returns true if no copies are made for UGGrid elements
Definition: dgfug.hh:46
bool noClosure() const
returns true if no closure should be used for UGGrid
Definition: dgfug.hh:44
static double refineWeight()
Definition: dgfug.hh:69
static int refineStepsForHalf()
Definition: dgfug.hh:64
const DGFBoundaryParameter::type & boundaryParameter(const Dune::Intersection< GG, II > &intersection) const
return invalid value
Definition: dgfug.hh:182
std::vector< double > & parameter(const typename Grid::template Codim< dimension >::Entity &vertex)
return parameter for vertex
Definition: dgfug.hh:164
DGFGridFactory(const std::string &filename, MPICommunicatorType comm=MPIHelper::getCommunicator())
constructor taking filename
Definition: dgfug.hh:101
MPIHelper::MPICommunicator MPICommunicatorType
MPI communicator type.
Definition: dgfug.hh:88
int numParameters(const Entity &) const
return number of parameters
Definition: dgfug.hh:147
int boundaryId(const Dune::Intersection< GG, II > &intersection) const
will return boundary segment index
Definition: dgfug.hh:128
UGGrid< dim > Grid
grid type
Definition: dgfug.hh:84
bool wasInserted(const Dune::Intersection< GG, II > &intersection) const
please doc me
Definition: dgfug.hh:121
DGFGridFactory(std::istream &input, MPICommunicatorType comm=MPIHelper::getCommunicator())
constructor taking istream
Definition: dgfug.hh:91
int numParameters() const
return number of parameters
Definition: dgfug.hh:135
bool haveBoundaryParameters() const
UGGrid does not support boundary parameters.
Definition: dgfug.hh:175
std::vector< double > & parameter(const typename Grid::template Codim< 0 >::Entity &element)
return parameter for codim 0 entity
Definition: dgfug.hh:153
Grid * grid()
return grid
Definition: dgfug.hh:114
static const type & defaultValue()
default constructor
Definition: parser.hh:28
std::string type
type of additional boundary parameters
Definition: parser.hh:25
Front-end for the grid manager of the finite element toolbox UG3.
Definition: uggrid.hh:207