dune-grid  2.4.1-rc2
dgfug.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_GRID_IO_FILE_DGFPARSER_DGFUG_HH
4 #define DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_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
19 #include <dune/grid/uggrid.hh>
20 
21 //- local includes
22 #include "dgfparser.hh"
23 #include "blocks/gridparameter.hh"
24 
25 
26 namespace Dune
27 {
28 
29  namespace dgf
30  {
31 
32  // UGGridParameterBlock
33  // --------------------
34 
36  : public GridParameterBlock
37  {
39  explicit UGGridParameterBlock ( std::istream &input );
40 
42  bool noClosure () const { return noClosure_; }
44  bool noCopy () const { return noCopy_; }
46  size_t heapSize () const { return heapSize_; }
47 
48  protected:
49  bool noClosure_; // no closure for UGGrid
50  bool noCopy_; // no copies for UGGrid
51  size_t heapSize_; // heap size for UGGrid
52  };
53 
54  } // namespace dgf
55 
56 
57 
58 #if ENABLE_UG
59  template< int dim >
60  struct DGFGridInfo< UGGrid< dim > >
61  {
62  static int refineStepsForHalf ()
63  {
64  return 1;
65  }
66 
67  static double refineWeight ()
68  {
69  return -1.;
70  }
71  };
72 
73 
74 
75  // DGFGridFactory< UGGrid< dim > >
76  // -------------------------------
77 
78  template< int dim >
79  struct DGFGridFactory< UGGrid< dim > >
80  {
82  typedef UGGrid< dim > Grid;
84  static const int dimension = dim;
86  typedef MPIHelper::MPICommunicator MPICommunicatorType;
87 
89  explicit DGFGridFactory ( std::istream &input,
90  MPICommunicatorType comm = MPIHelper::getCommunicator() )
91  : grid_( 0 ),
92  factory_(),
93  dgf_( rank( comm ), size( comm ) )
94  {
95  generate( input );
96  }
97 
99  explicit DGFGridFactory ( const std::string &filename,
100  MPICommunicatorType comm = MPIHelper::getCommunicator() )
101  : grid_( 0 ),
102  factory_(),
103  dgf_( rank( comm ), size( comm ) )
104  {
105  std::ifstream input( filename.c_str() );
106  if ( !input )
107  DUNE_THROW( DGFException, "Error: Macrofile " << filename << " not found" );
108  generate( input );
109  }
110 
112  Grid *grid ()
113  {
114  return grid_;
115  }
116 
118  template< class GG, class II >
119  bool wasInserted ( const Dune::Intersection< GG, II > &intersection ) const
120  {
121  return factory_.wasInserted( intersection );
122  }
123 
125  template< class GG, class II >
126  int boundaryId ( const Dune::Intersection< GG, II > &intersection ) const
127  {
128  return intersection.boundarySegmentIndex();
129  }
130 
132  template< int codim >
133  int numParameters () const
134  {
135  if( codim == 0 )
136  return dgf_.nofelparams;
137  else if( codim == dimension )
138  return dgf_.nofvtxparams;
139  else
140  return 0;
141  }
142 
144  template< class Entity >
145  int numParameters ( const Entity & ) const
146  {
147  return numParameters< Entity::codimension >();
148  }
149 
151  std::vector< double > &parameter ( const typename Grid::template Codim< 0 >::Entity &element )
152  {
153  if( numParameters< 0 >() <= 0 )
154  {
155  DUNE_THROW( InvalidStateException,
156  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
157  }
158  return dgf_.elParams[ factory_.insertionIndex( element ) ];
159  }
160 
162  std::vector< double > &parameter ( const typename Grid::template Codim< dimension >::Entity &vertex )
163  {
164  if( numParameters< dimension >() <= 0 )
165  {
166  DUNE_THROW( InvalidStateException,
167  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
168  }
169  return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
170  }
171 
173  bool haveBoundaryParameters () const
174  {
175  return dgf_.haveBndParameters;
176  }
177 
179  template< class GG, class II >
180  const DGFBoundaryParameter::type &boundaryParameter ( const Dune::Intersection< GG, II > &intersection ) const
181  {
183  typename Intersection::Entity entity = intersection.inside();
184  const int face = intersection.indexInInside();
185 
186  const ReferenceElement< double, dimension > &refElem
187  = ReferenceElements< double, dimension >::general( entity.type() );
188  int corners = refElem.size( face, 1, dimension );
189  std::vector< unsigned int > bound( corners );
190  for( int i = 0; i < corners; ++i )
191  {
192  const int k = refElem.subEntity( face, 1, i, dimension );
193  bound[ i ] = factory_.insertionIndex( entity.template subEntity< dimension >( k ) );
194  }
195 
196  DuneGridFormatParser::facemap_t::key_type key( bound, false );
197  const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
198  if( pos != dgf_.facemap.end() )
199  return dgf_.facemap.find( key )->second.second;
200  else
202  }
203 
204  private:
205  // create grid
206  void generate ( std::istream &input );
207 
208  // return rank
209  static int rank( MPICommunicatorType MPICOMM )
210  {
211  int rank = 0;
212 #if HAVE_MPI
213  MPI_Comm_rank( MPICOMM, &rank );
214 #endif
215  return rank;
216  }
217 
218  // return size
219  static int size( MPICommunicatorType MPICOMM )
220  {
221  int size = 1;
222 #if HAVE_MPI
223  MPI_Comm_size( MPICOMM, &size );
224 #endif
225  return size;
226  }
227 
228  Grid *grid_;
229  GridFactory< UGGrid< dim > > factory_;
231  };
232 #endif // #if ENABLE_UG
233 
234 } // namespace Dune
235 
236 #endif // #ifndef DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
bool noCopy_
Definition: dgfug.hh:50
bool noClosure() const
returns true if no closure should be used for UGGrid
Definition: dgfug.hh:42
Intersection of a mesh entities of codimension 0 ("elements") with a "neighboring" element or with th...
Definition: albertagrid/dgfparser.hh:26
Definition: dgfug.hh:35
Front-end for the grid manager of the finite element toolbox UG.
Definition: uggrid.hh:207
static const type & defaultValue()
default constructor
Definition: parser.hh:26
Wrapper class for entities.
Definition: common/entity.hh:61
size_t heapSize_
Definition: dgfug.hh:51
size_t boundarySegmentIndex() const
index of the boundary segment within the macro grid
Definition: common/intersection.hh:260
Provide a generic factory class for unstructured grids.
Definition: common/gridfactory.hh:263
std::string type
type of additional boundary parameters
Definition: parser.hh:23
Entity inside() const
return EntityPointer to the Entity on the inside of this intersection. That is the Entity where we st...
Definition: common/intersection.hh:286
bool noClosure_
Definition: dgfug.hh:49
UGGridParameterBlock(std::istream &input)
constructor taking istream
Definition: dgfug.cc:16
The UGGrid class.
size_t heapSize() const
returns heap size used on construction of the grid
Definition: dgfug.hh:46
Grid abstract base classThis class is the base class for all grid implementations. Although no virtual functions are used we call it abstract since its methods do not contain an implementation but forward to the methods of the derived class via the Barton-Nackman trick.
Definition: common/grid.hh:388
Definition: common.hh:179
Include standard header files.
Definition: agrid.hh:59
The DuneGridFormatParser class: reads a DGF file and stores build information in vector structures us...
Definition: parser.hh:44
int indexInInside() const
Local index of codim 1 entity in the inside() entity where intersection is contained in...
Definition: common/intersection.hh:393
Definition: agrid.hh:66
GridImp::template Codim< 0 >::Entity Entity
Type of entity that this Intersection belongs to.
Definition: common/intersection.hh:185
Common Grid parametersFor each grid implementation there is a set of parameters that can be passed vi...
Definition: gridparameter.hh:31
Some simple static information for a given GridType.
Definition: io/file/dgfparser/dgfparser.hh:54
bool noCopy() const
returns true if no copies are made for UGGrid elements
Definition: dgfug.hh:44
exception class for IO errors in the DGF parser
Definition: dgfexception.hh:12