3 #ifndef DUNE_ALBERTA_DGFPARSER_HH
4 #define DUNE_ALBERTA_DGFPARSER_HH
25 template<
class Gr
idImp,
class IntersectionImp >
33 template<
int dim,
int dimworld >
34 struct DGFGridFactory< AlbertaGrid< dim, dimworld > >
36 typedef AlbertaGrid<dim,dimworld>
Grid;
39 typedef typename Grid::template Codim<0>::Entity Element;
40 typedef typename Grid::template Codim<dimension>::Entity Vertex;
53 template<
class Intersection >
54 bool wasInserted (
const Intersection &intersection )
const
56 return factory_.wasInserted( intersection );
59 template<
class Intersection >
60 int boundaryId (
const Intersection &intersection )
const
62 return intersection.impl().boundaryId();
68 return dgf_.haveBndParameters;
71 template <
class GG,
class II >
77 const int face = intersection.indexInInside();
79 auto refElem = referenceElement< double, dimension >( entity.type() );
80 int corners = refElem.size( face, 1,
dimension );
81 std :: vector< unsigned int > bound( corners );
82 for(
int i=0; i < corners; ++i )
84 const int k = refElem.subEntity( face, 1, i,
dimension );
85 bound[ i ] = factory_.insertionIndex( entity.template subEntity< dimension >( k ) );
88 DuneGridFormatParser::facemap_t::key_type key( bound,
false );
89 const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
90 if( pos != dgf_.facemap.end() )
91 return dgf_.facemap.find( key )->second.second;
100 return dgf_.nofelparams;
102 return dgf_.nofvtxparams;
107 std::vector< double > &
parameter (
const Element &element )
109 if( numParameters< 0 >() <= 0 )
111 DUNE_THROW( InvalidStateException,
112 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
114 return dgf_.elParams[ factory_.insertionIndex( element ) ];
119 if( numParameters< dimension >() <= 0 )
121 DUNE_THROW( InvalidStateException,
122 "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
124 return dgf_.vtxParams[ factory_.insertionIndex(
vertex ) ];
128 bool generate( std::istream &input );
131 GridFactory factory_;
132 DuneGridFormatParser dgf_;
140 template<
int dim,
int dimworld >
141 struct DGFGridInfo< AlbertaGrid< dim, dimworld > >
159 template<
int dim,
int dimworld >
160 inline DGFGridFactory< AlbertaGrid< dim, dimworld > >
161 ::DGFGridFactory ( std::istream &input, MPICommunicatorType comm )
167 DUNE_THROW(DGFException,
"Error resetting input stream." );
172 template<
int dim,
int dimworld >
173 inline DGFGridFactory< AlbertaGrid< dim, dimworld > >
174 ::DGFGridFactory (
const std::string &filename, MPICommunicatorType comm )
177 std::ifstream input( filename.c_str() );
179 DUNE_THROW( DGFException,
"Macrofile " << filename <<
" not found." );
180 if( !generate( input ) )
181 grid_ =
new AlbertaGrid< dim, dimworld >( filename.c_str() );
187 #endif // #if HAVE_ALBERTA
189 #endif // #ifndef DUNE_ALBERTA_DGFPARSER_HH