dune-foamgrid  2.8-git
dgffoam.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_FOAMGRID_DGFFOAM_HH
4 #define DUNE_FOAMGRID_DGFFOAM_HH
5 
6 //- C++ includes
7 #include <fstream>
8 #include <istream>
9 #include <string>
10 #include <vector>
11 
12 //- dune-common includes
13 #include <dune/common/exceptions.hh>
14 #include <dune/common/fvector.hh>
15 #include <dune/common/parallel/mpihelper.hh>
16 
17 //- dune-grid includes
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>
23 
24 // - dune-foamgrid includes
26 
27 namespace Dune
28 {
29 
30  namespace dgf
31  {
32 
33  // FoamGridParameterBlock
34  // --------------------
35 
37  : public GridParameterBlock
38  {
40  explicit FoamGridParameterBlock ( std::istream &input );
41  };
42 
43  } // namespace dgf
44 
45 
46 
47  template< int dim , int dimworld , class ctype >
48  struct DGFGridInfo< FoamGrid< dim, dimworld, ctype > >
49  {
50  static int refineStepsForHalf ()
51  {
52  return 1;
53  }
54 
55  static double refineWeight ()
56  {
57  return -1.;
58  }
59  };
60 
61 
62 
63  // DGFGridFactory< FoamGrid< dim , dimworld> >
64  // -------------------------------
65 
66  template< int dim , int dimworld , class ctype >
67  struct DGFGridFactory< FoamGrid< dim, dimworld, ctype > >
68  {
72  static const int dimension = dim;
74  static const int dimensionWorld = dimworld;
76  typedef MPIHelper::MPICommunicator MPICommunicatorType;
77 
79  explicit DGFGridFactory ( std::istream &input,
80  MPICommunicatorType comm = MPIHelper::getCommunicator() )
81  : grid_( 0 ),
82  factory_(),
83  dgf_( rank( comm ), size( comm ) )
84  {
85  generate( input );
86  }
87 
89  explicit DGFGridFactory ( const std::string &filename,
90  MPICommunicatorType comm = MPIHelper::getCommunicator() )
91  : grid_( 0 ),
92  factory_(),
93  dgf_( rank( comm ), size( comm ) )
94  {
95  std::ifstream input( filename.c_str() );
96  if ( !input )
97  DUNE_THROW( DGFException, "Error: Macrofile " << filename << " not found" );
98  generate( input );
99  }
100 
102  Grid *grid ()
103  {
104  return grid_;
105  }
106 
108  template< class GG, class II >
109  bool wasInserted ( const Dune::Intersection< GG, II > &intersection ) const
110  {
111  return factory_.wasInserted( intersection );
112  }
113 
115  template< class GG, class II >
116  int boundaryId ( const Dune::Intersection< GG, II > &intersection ) const
117  {
118  return intersection.boundarySegmentIndex();
119  }
120 
122  template< int codim >
123  int numParameters () const
124  {
125  if( codim == 0 )
126  return dgf_.nofelparams;
127  else if( codim == dimension )
128  return dgf_.nofvtxparams;
129  else
130  return 0;
131  }
132 
134  template< class Entity >
135  int numParameters ( const Entity & ) const
136  {
137  return numParameters< Entity::codimension >();
138  }
139 
141  std::vector< double > &parameter ( const typename Grid::template Codim< 0 >::Entity &element )
142  {
143  if( numParameters< 0 >() <= 0 )
144  {
145  DUNE_THROW( InvalidStateException,
146  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
147  }
148  return dgf_.elParams[ factory_.insertionIndex( element ) ];
149  }
150 
152  std::vector< double > &parameter ( const typename Grid::template Codim< dimension >::Entity &vertex )
153  {
154  if( numParameters< dimension >() <= 0 )
155  {
156  DUNE_THROW( InvalidStateException,
157  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
158  }
159  return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
160  }
161 
164  {
165  return dgf_.haveBndParameters;
166  }
167 
169  template< class GG, class II >
170  const DGFBoundaryParameter::type &boundaryParameter ( const Dune::Intersection< GG, II > &intersection ) const
171  {
172  const auto& entity = intersection.inside();
173  const int face = intersection.indexInInside();
174 
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 )
179  {
180  const int k = refElem.subEntity( face, 1, i, dimension );
181  bound[ i ] = factory_.insertionIndex( entity.template subEntity< dimension >( k ) );
182  }
183 
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;
188  else
189  return DGFBoundaryParameter::defaultValue();
190  }
191 
192  private:
193  // create grid
194  void generate ( std::istream &input );
195 
196  // return rank
197  static int rank( MPICommunicatorType MPICOMM )
198  {
199  int rank = 0;
200 #if HAVE_MPI
201  MPI_Comm_rank( MPICOMM, &rank );
202 #endif
203  return rank;
204  }
205 
206  // return size
207  static int size( MPICommunicatorType MPICOMM )
208  {
209  int size = 1;
210 #if HAVE_MPI
211  MPI_Comm_size( MPICOMM, &size );
212 #endif
213  return size;
214  }
215 
216  Grid *grid_;
217  GridFactory< FoamGrid< dim, dimworld, ctype > > factory_;
218  DuneGridFormatParser dgf_;
219  };
220 
221 } // namespace Dune
222 
223 #include "dgffoam.cc"
224 
225 #endif // #ifndef DUNE_FOAMGRID_DGFFOAM_HH
The FoamGrid class.
Definition: dgffoam.cc:6
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:96