dune-grid  2.9.0
dgfug.hh
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (C) DUNE Project contributors, see file LICENSE.md in module root
2 // SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
3 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
4 // vi: set et ts=4 sw=2 sts=2:
5 #ifndef DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
6 #define DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
7 
8 //- C++ includes
9 #include <fstream>
10 #include <istream>
11 #include <string>
12 #include <vector>
13 
14 //- dune-common includes
15 #include <dune/common/exceptions.hh>
16 #include <dune/common/fvector.hh>
17 #include <dune/common/parallel/mpihelper.hh>
18 
19 //- dune-grid includes
21 #include <dune/grid/uggrid.hh>
22 
23 //- local includes
24 #include "dgfparser.hh"
25 #include "blocks/gridparameter.hh"
26 
27 
28 namespace Dune
29 {
30 
31  namespace dgf
32  {
33 
34  // UGGridParameterBlock
35  // --------------------
36 
38  : public GridParameterBlock
39  {
41  explicit UGGridParameterBlock ( std::istream &input );
42 
44  bool noClosure () const { return noClosure_; }
46  bool noCopy () const { return noCopy_; }
48  size_t heapSize () const { return heapSize_; }
49 
50  protected:
51  bool noClosure_; // no closure for UGGrid
52  bool noCopy_; // no copies for UGGrid
53  size_t heapSize_; // heap size for UGGrid
54  };
55 
56  } // namespace dgf
57 
58 
59 
60 #if HAVE_DUNE_UGGRID
61  template< int dim >
62  struct DGFGridInfo< UGGrid< dim > >
63  {
64  static int refineStepsForHalf ()
65  {
66  return 1;
67  }
68 
69  static double refineWeight ()
70  {
71  return -1.;
72  }
73  };
74 
75 
76 
77  // DGFGridFactory< UGGrid< dim > >
78  // -------------------------------
79 
80  template< int dim >
81  struct DGFGridFactory< UGGrid< dim > >
82  {
86  static const int dimension = dim;
88  typedef MPIHelper::MPICommunicator MPICommunicatorType;
89 
91  explicit DGFGridFactory ( std::istream &input,
92  MPICommunicatorType comm = MPIHelper::getCommunicator() )
93  : grid_( 0 ),
94  factory_(),
95  dgf_( rank( comm ), size( comm ) )
96  {
97  generate( input );
98  }
99 
101  explicit DGFGridFactory ( const std::string &filename,
102  MPICommunicatorType comm = MPIHelper::getCommunicator() )
103  : grid_( 0 ),
104  factory_(),
105  dgf_( rank( comm ), size( comm ) )
106  {
107  std::ifstream input( filename.c_str() );
108  if ( !input )
109  DUNE_THROW( DGFException, "Error: Macrofile " << filename << " not found" );
110  generate( input );
111  }
112 
114  Grid *grid ()
115  {
116  return grid_;
117  }
118 
120  template< class GG, class II >
121  bool wasInserted ( const Dune::Intersection< GG, II > &intersection ) const
122  {
123  return factory_.wasInserted( intersection );
124  }
125 
127  template< class GG, class II >
128  int boundaryId ( const Dune::Intersection< GG, II > &intersection ) const
129  {
130  return intersection.boundarySegmentIndex();
131  }
132 
134  template< int codim >
135  int numParameters () const
136  {
137  if( codim == 0 )
138  return dgf_.nofelparams;
139  else if( codim == dimension )
140  return dgf_.nofvtxparams;
141  else
142  return 0;
143  }
144 
146  template< class Entity >
147  int numParameters ( const Entity & ) const
148  {
149  return numParameters< Entity::codimension >();
150  }
151 
153  std::vector< double > &parameter ( const typename Grid::template Codim< 0 >::Entity &element )
154  {
155  if( numParameters< 0 >() <= 0 )
156  {
157  DUNE_THROW( InvalidStateException,
158  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
159  }
160  return dgf_.elParams[ factory_.insertionIndex( element ) ];
161  }
162 
164  std::vector< double > &parameter ( const typename Grid::template Codim< dimension >::Entity &vertex )
165  {
166  if( numParameters< dimension >() <= 0 )
167  {
168  DUNE_THROW( InvalidStateException,
169  "Calling DGFGridFactory::parameter is only allowed if there are parameters." );
170  }
171  return dgf_.vtxParams[ factory_.insertionIndex( vertex ) ];
172  }
173 
176  {
177  return dgf_.haveBndParameters;
178  }
179 
181  template< class GG, class II >
183  {
185  typename Intersection::Entity entity = intersection.inside();
186  const int face = intersection.indexInInside();
187 
188  auto refElem = referenceElement< double, dimension >( entity.type() );
189  int corners = refElem.size( face, 1, dimension );
190  std::vector< unsigned int > bound( corners );
191  for( int i = 0; i < corners; ++i )
192  {
193  const int k = refElem.subEntity( face, 1, i, dimension );
194  bound[ i ] = factory_.insertionIndex( entity.template subEntity< dimension >( k ) );
195  }
196 
197  DuneGridFormatParser::facemap_t::key_type key( bound, false );
198  const DuneGridFormatParser::facemap_t::const_iterator pos = dgf_.facemap.find( key );
199  if( pos != dgf_.facemap.end() )
200  return dgf_.facemap.find( key )->second.second;
201  else
203  }
204 
205  private:
206  // create grid
207  void generate ( std::istream &input );
208 
209  // return rank
210  static int rank( MPICommunicatorType MPICOMM )
211  {
212  int rank = 0;
213 #if HAVE_MPI
214  MPI_Comm_rank( MPICOMM, &rank );
215 #endif
216  return rank;
217  }
218 
219  // return size
220  static int size( MPICommunicatorType MPICOMM )
221  {
222  int size = 1;
223 #if HAVE_MPI
224  MPI_Comm_size( MPICOMM, &size );
225 #endif
226  return size;
227  }
228 
229  Grid *grid_;
230  GridFactory< UGGrid< dim > > factory_;
231  DuneGridFormatParser dgf_;
232  };
233 #endif // #if HAVE_DUNE_UGGRID
234 
235 } // namespace Dune
236 
237 #endif // #ifndef DUNE_GRID_IO_FILE_DGFPARSER_DGFUG_HH
The UGGrid class.
Include standard header files.
Definition: agrid.hh:60
@ vertex
Definition: common.hh:133
Definition: dgfgridfactory.hh:38
MPIHelper::MPICommunicator MPICommunicatorType
Definition: dgfgridfactory.hh:41
G Grid
Definition: dgfgridfactory.hh:39
static const int dimension
Definition: dgfgridfactory.hh:40
Intersection of a mesh entity of codimension 0 ("element") with a "neighboring" element or with the d...
Definition: common/intersection.hh:164
int indexInInside() const
Local index of codim 1 entity in the inside() entity where intersection is contained in.
Definition: common/intersection.hh:346
size_t boundarySegmentIndex() const
index of the boundary segment within the macro grid
Definition: common/intersection.hh:236
Entity inside() const
return Entity on the inside of this intersection. That is the Entity where we started this.
Definition: common/intersection.hh:250
GridImp::template Codim< 0 >::Entity Entity
Type of entity that this Intersection belongs to.
Definition: common/intersection.hh:192
Wrapper class for entities.
Definition: common/entity.hh:66
Common Grid parameters.
Definition: gridparameter.hh:35
exception class for IO errors in the DGF parser
Definition: dgfexception.hh:16
Some simple static information for a given GridType.
Definition: io/file/dgfparser/dgfparser.hh:56
Definition: dgfug.hh:39
size_t heapSize() const
returns heap size used on construction of the grid
Definition: dgfug.hh:48
bool noCopy_
Definition: dgfug.hh:52
UGGridParameterBlock(std::istream &input)
constructor taking istream
Definition: dgfug.cc:20
size_t heapSize_
Definition: dgfug.hh:53
bool noClosure_
Definition: dgfug.hh:51
bool noCopy() const
returns true if no copies are made for UGGrid elements
Definition: dgfug.hh:46
bool noClosure() const
returns true if no closure should be used for UGGrid
Definition: dgfug.hh:44
static double refineWeight()
Definition: dgfug.hh:69
static int refineStepsForHalf()
Definition: dgfug.hh:64
std::vector< double > & parameter(const typename Grid::template Codim< 0 >::Entity &element)
return parameter for codim 0 entity
Definition: dgfug.hh:153
DGFGridFactory(const std::string &filename, MPICommunicatorType comm=MPIHelper::getCommunicator())
constructor taking filename
Definition: dgfug.hh:101
MPIHelper::MPICommunicator MPICommunicatorType
MPI communicator type.
Definition: dgfug.hh:88
int numParameters(const Entity &) const
return number of parameters
Definition: dgfug.hh:147
std::vector< double > & parameter(const typename Grid::template Codim< dimension >::Entity &vertex)
return parameter for vertex
Definition: dgfug.hh:164
Grid * grid()
return grid
Definition: dgfug.hh:114
const DGFBoundaryParameter::type & boundaryParameter(const Dune::Intersection< GG, II > &intersection) const
return invalid value
Definition: dgfug.hh:182
int boundaryId(const Dune::Intersection< GG, II > &intersection) const
will return boundary segment index
Definition: dgfug.hh:128
UGGrid< dim > Grid
grid type
Definition: dgfug.hh:84
bool wasInserted(const Dune::Intersection< GG, II > &intersection) const
please doc me
Definition: dgfug.hh:121
DGFGridFactory(std::istream &input, MPICommunicatorType comm=MPIHelper::getCommunicator())
constructor taking istream
Definition: dgfug.hh:91
int numParameters() const
return number of parameters
Definition: dgfug.hh:135
bool haveBoundaryParameters() const
UGGrid does not support boundary parameters.
Definition: dgfug.hh:175
static const type & defaultValue()
default constructor
Definition: parser.hh:28
std::string type
type of additional boundary parameters
Definition: parser.hh:25
Front-end for the grid manager of the finite element toolbox UG3.
Definition: uggrid.hh:207