3 #ifndef DUNE_FOAMGRID_DGFFOAM_HH
4 #define DUNE_FOAMGRID_DGFFOAM_HH
13 #include <dune/common/exceptions.hh>
14 #include <dune/common/fvector.hh>
15 #include <dune/common/parallel/mpihelper.hh>
18 #include <dune/grid/common/intersection.hh>
19 #include <dune/grid/io/file/dgfparser.hh>
20 #include <dune/grid/io/file/dgfparser/dgfparser.hh>
21 #include <dune/grid/io/file/dgfparser/parser.hh>
22 #include <dune/grid/io/file/dgfparser/blocks/gridparameter.hh>
37 :
public GridParameterBlock
47 template<
int dim ,
int dimworld ,
class ctype >
48 struct DGFGridInfo<
FoamGrid< dim, dimworld, ctype > >
66 template<
int dim ,
int dimworld ,
class ctype >
67 struct DGFGridFactory<
FoamGrid< dim, dimworld, ctype > >
72 static const int dimension = dim;
74 static const int dimensionWorld = dimworld;
83 dgf_( rank( comm ), size( comm ) )
93 dgf_( rank( comm ), size( comm ) )
95 std::ifstream input( filename.c_str() );
97 DUNE_THROW( DGFException,
"Error: Macrofile " << filename <<
" not found" );
108 template<
class GG,
class II >
109 bool wasInserted (
const Dune::Intersection< GG, II > &intersection )
const
111 return factory_.wasInserted( intersection );
115 template<
class GG,
class II >
116 int boundaryId (
const Dune::Intersection< GG, II > &intersection )
const
118 return intersection.boundarySegmentIndex();
122 template<
int codim >
126 return dgf_.nofelparams;
127 else if( codim == dimension )
128 return dgf_.nofvtxparams;
134 template<
class Entity >
137 return numParameters< Entity::codimension >();
141 std::vector< double > &
parameter (
const typename Grid::template Codim< 0 >::Entity &element )
143 if( numParameters< 0 >() <= 0 )
145 DUNE_THROW( InvalidStateException,
146 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
148 return dgf_.elParams[ factory_.insertionIndex( element ) ];
152 std::vector< double > &
parameter (
const typename Grid::template Codim< dimension >::Entity &vertex )
154 if( numParameters< dimension >() <= 0 )
156 DUNE_THROW( InvalidStateException,
157 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
159 return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
165 return dgf_.haveBndParameters;
169 template<
class GG,
class II >
170 const DGFBoundaryParameter::type &
boundaryParameter (
const Dune::Intersection< GG, II > &intersection )
const
172 const auto& entity = intersection.inside();
173 const int face = intersection.indexInInside();
175 const auto refElem = ReferenceElements< double, dimension >::general( entity.type() );
176 int corners = refElem.size( face, 1, dimension );
177 std::vector< unsigned int > bound( corners );
178 for(
int i = 0; i < corners; ++i )
180 const int k = refElem.subEntity( face, 1, i, dimension );
181 bound[ i ] = factory_.insertionIndex( entity.template subEntity< dimension >( k ) );
184 DuneGridFormatParser::facemap_t::key_type key( bound,
false );
185 const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
186 if( pos != dgf_.facemap.end() )
187 return dgf_.facemap.find( key )->second.second;
189 return DGFBoundaryParameter::defaultValue();
194 void generate ( std::istream &input );
197 static int rank( MPICommunicatorType MPICOMM )
201 MPI_Comm_rank( MPICOMM, &rank );
207 static int size( MPICommunicatorType MPICOMM )
211 MPI_Comm_size( MPICOMM, &size );
217 GridFactory< FoamGrid< dim, dimworld, ctype > > factory_;
218 DuneGridFormatParser dgf_;
Definition: dgffoam.hh:38
FoamGridParameterBlock(std::istream &input)
constructor taking istream
static int refineStepsForHalf()
Definition: dgffoam.hh:50
static double refineWeight()
Definition: dgffoam.hh:55
std::vector< double > & parameter(const typename Grid::template Codim< dimension >::Entity &vertex)
return parameter for vertex
Definition: dgffoam.hh:152
bool haveBoundaryParameters() const
FoamGrid does not support boundary parameters.
Definition: dgffoam.hh:163
DGFGridFactory(const std::string &filename, MPICommunicatorType comm=MPIHelper::getCommunicator())
constructor taking filename
Definition: dgffoam.hh:89
DGFGridFactory(std::istream &input, MPICommunicatorType comm=MPIHelper::getCommunicator())
constructor taking istream
Definition: dgffoam.hh:79
Grid * grid()
return grid
Definition: dgffoam.hh:102
int boundaryId(const Dune::Intersection< GG, II > &intersection) const
will return boundary segment index
Definition: dgffoam.hh:116
const DGFBoundaryParameter::type & boundaryParameter(const Dune::Intersection< GG, II > &intersection) const
return invalid value
Definition: dgffoam.hh:170
MPIHelper::MPICommunicator MPICommunicatorType
MPI communicator type.
Definition: dgffoam.hh:76
bool wasInserted(const Dune::Intersection< GG, II > &intersection) const
please doc me
Definition: dgffoam.hh:109
std::vector< double > & parameter(const typename Grid::template Codim< 0 >::Entity &element)
return parameter for codim 0 entity
Definition: dgffoam.hh:141
FoamGrid< dim, dimworld, ctype > Grid
grid type
Definition: dgffoam.hh:70
int numParameters(const Entity &) const
return number of parameters
Definition: dgffoam.hh:135
int numParameters() const
return number of parameters
Definition: dgffoam.hh:123
An implementation of the Dune grid interface: a 1- or 2-dimensional simplicial grid in an n-dimension...
Definition: foamgrid.hh:92