dune-grid  2.2.0
dgfug.hh
Go to the documentation of this file.
1 #ifndef DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
2 #define DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
3 
4 //- C++ includes
5 #include <fstream>
6 #include <istream>
7 #include <string>
8 #include <vector>
9 
10 //- dune-common includes
11 #include <dune/common/exceptions.hh>
12 #include <dune/common/fvector.hh>
13 #include <dune/common/mpihelper.hh>
14 
15 //- dune-grid includes
17 #include <dune/grid/uggrid.hh>
18 
19 //- local includes
20 #include "dgfparser.hh"
21 #include "blocks/gridparameter.hh"
22 
23 
24 namespace Dune
25 {
26 
27  namespace dgf
28  {
29 
30  // UGGridParameterBlock
31  // --------------------
32 
34  : public GridParameterBlock
35  {
37  explicit UGGridParameterBlock ( std::istream &input );
38 
40  bool noClosure () const { return noClosure_; }
42  bool noCopy () const { return noCopy_; }
44  size_t heapSize () const { return heapSize_; }
45 
46  protected:
47  bool noClosure_; // no closure for UGGrid
48  bool noCopy_; // no copies for UGGrid
49  size_t heapSize_; // heap size for UGGrid
50  };
51 
52  } // namespace dgf
53 
54 
55 
56 #if ENABLE_UG
57  template< int dim >
58  struct DGFGridInfo< UGGrid< dim > >
59  {
60  static int refineStepsForHalf ()
61  {
62  return 1;
63  }
64 
65  static double refineWeight ()
66  {
67  return -1.;
68  }
69  };
70 
71 
72 
73  // DGFGridFactory< UGGrid< dim > >
74  // -------------------------------
75 
76  template< int dim >
77  struct DGFGridFactory< UGGrid< dim > >
78  {
80  typedef UGGrid< dim > Grid;
82  static const int dimension = dim;
84  typedef MPIHelper::MPICommunicator MPICommunicatorType;
85 
87  explicit DGFGridFactory ( std::istream &input,
88  MPICommunicatorType comm = MPIHelper::getCommunicator() )
89  : grid_( 0 ),
90  factory_(),
91  dgf_( rank( comm ), size( comm ) )
92  {
93  generate( input );
94  }
95 
97  explicit DGFGridFactory ( const std::string &filename,
98  MPICommunicatorType comm = MPIHelper::getCommunicator() )
99  : grid_( 0 ),
100  factory_(),
101  dgf_( rank( comm ), size( comm ) )
102  {
103  std::ifstream input( filename.c_str() );
104  if ( !input )
105  DUNE_THROW( DGFException, "Error: Macrofile " << filename << " not found" );
106  generate( input );
107  }
108 
110  Grid *grid ()
111  {
112  return grid_;
113  }
114 
116  template< class GG, template< class > class II >
117  bool wasInserted ( const Dune::Intersection< GG, II > &intersection ) const
118  {
119  return factory_.wasInserted( intersection );
120  }
121 
123  template < class GG, template< class > class II >
124  int boundaryId ( const Dune::Intersection< GG, II > &intersection ) const
125  {
126  return intersection.boundarySegmentIndex();
127  }
128 
130  template< int codim >
131  int numParameters () const
132  {
133  if( codim == 0 )
134  return dgf_.nofelparams;
135  else if( codim == dimension )
136  return dgf_.nofvtxparams;
137  else
138  return 0;
139  }
140 
142  template< class Entity >
143  int numParameters ( const Entity & ) const
144  {
145  return numParameters< Entity::codimension >();
146  }
147 
149  std::vector< double > &parameter ( const typename Grid::template Codim< 0 >::Entity &element )
150  {
151  if( numParameters< 0 >() <= 0 )
152  {
153  DUNE_THROW( InvalidStateException,
154  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
155  }
156  return dgf_.elParams[ factory_.insertionIndex( element ) ];
157  }
158 
160  std::vector< double > &parameter ( const typename Grid::template Codim< dimension >::Entity &vertex )
161  {
162  if( numParameters< dimension >() <= 0 )
163  {
164  DUNE_THROW( InvalidStateException,
165  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
166  }
167  return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
168  }
169 
171  bool haveBoundaryParameters () const
172  {
173  return dgf_.haveBndParameters;
174  }
175 
177  template < class GG, template< class > class II >
179  {
180  typedef Dune::Intersection< GG, II > Intersection;
181  typename Intersection::EntityPointer inside = intersection.inside();
182  const typename Intersection::Entity &entity = *inside;
183  const int face = intersection.indexInInside();
184 
185  const GenericReferenceElement< double, dimension > &refElem
186  = GenericReferenceElements< double, dimension >::general( entity.type() );
187  int corners = refElem.size( face, 1, dimension );
188  std::vector< unsigned int > bound( corners );
189  for( int i = 0; i < corners; ++i )
190  {
191  const int k = refElem.subEntity( face, 1, i, dimension );
192  bound[ i ] = factory_.insertionIndex( *entity.template subEntity< dimension >( k ) );
193  }
194 
195  DuneGridFormatParser::facemap_t::key_type key( bound, false );
196  const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
197  if( pos != dgf_.facemap.end() )
198  return dgf_.facemap.find( key )->second.second;
199  else
201  }
202 
203  private:
204  // create grid
205  void generate ( std::istream &input );
206 
207  // return rank
208  static int rank( MPICommunicatorType MPICOMM )
209  {
210  int rank = 0;
211 #if HAVE_MPI
212  MPI_Comm_rank( MPICOMM, &rank );
213 #endif
214  return rank;
215  }
216 
217  // return size
218  static int size( MPICommunicatorType MPICOMM )
219  {
220  int size = 1;
221 #if HAVE_MPI
222  MPI_Comm_size( MPICOMM, &size );
223 #endif
224  return size;
225  }
226 
227  Grid *grid_;
228  GridFactory< UGGrid< dim > > factory_;
229  DuneGridFormatParser dgf_;
230  };
231 #endif // #if ENABLE_UG
232 
233 } // namespace Dune
234 
235 #endif // #ifndef DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH