dune-grid  2.9.0
dgfgridfactory.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_DGF_GRIDFACTORY_HH
6 #define DUNE_DGF_GRIDFACTORY_HH
7 
8 #include <iostream>
9 #include <string>
10 #include <vector>
11 #include <map>
12 #include <assert.h>
13 
14 #include <dune/common/parallel/mpihelper.hh>
17 
20 
21 
22 namespace Dune
23 {
24 
25  // External Forward Declarations
26  // -----------------------------
27 
28  template < class GridImp, class IntersectionImp >
29  class Intersection;
30 
31 
32 
33  // DGFGridFactory
34  // --------------
35 
36  template < class G >
38  {
39  typedef G Grid;
40  const static int dimension = Grid::dimension;
41  typedef MPIHelper::MPICommunicator MPICommunicatorType;
42 
43  private:
44  typedef typename Grid::template Codim< 0 >::Entity Element;
45 
46  typedef typename Grid::template Codim< dimension >::Entity Vertex;
47 
48  public:
49 
50  explicit DGFGridFactory ( const std::string &filename,
51  MPICommunicatorType comm = MPIHelper::getCommunicator() )
52  : macroGrid_( filename.c_str(), comm )
53  {
54  grid_ = macroGrid_.template createGrid< Grid >();
55 
56  if( macroGrid_.nofelparams > 0 )
57  {
58  const size_t nofElements = macroGrid_.elements.size();
59  for( size_t i = 0; i < nofElements; ++i )
60  {
61  std::vector< double > coord;
62 
63  DomainType p(0);
64  const size_t nofCorners = macroGrid_.elements[i].size();
65  for (size_t k=0; k<nofCorners; ++k)
66  for (int j=0; j<DomainType::dimension; ++j)
67  p[j]+=macroGrid_.vtx[macroGrid_.elements[i][k]][j];
68  p/=double(nofCorners);
69 
70  elInsertOrder_.insert( std::make_pair( p, i ) );
71  }
72  }
73 
74  if( macroGrid_.nofvtxparams > 0 )
75  {
76  const size_t nofVertices = macroGrid_.vtx.size();
77  for( size_t i = 0; i < nofVertices; ++i )
78  {
79  std::vector< double > coord;
80 
81  DomainType p;
82  for( int k = 0; k < DomainType::dimension; ++k )
83  p[ k ] = macroGrid_.vtx[i][k];
84 
85  vtxInsertOrder_.insert( std::make_pair( p, i ) );
86  }
87  }
88  }
89 
91  {
92  return grid_;
93  }
94 
95  template <class Intersection>
96  bool wasInserted(const Intersection &intersection) const
97  {
98  return intersection.boundary();
99  }
100 
101  template <class Intersection>
102  int boundaryId(const Intersection &intersection) const
103  {
104  return (intersection.boundary()) ? int(intersection.indexInInside()+1) : int(0);
105  }
106 
107  template< int codim >
108  int numParameters () const
109  {
110  if( codim == 0 )
111  return macroGrid_.nofelparams;
112  else if( codim == dimension )
113  return macroGrid_.nofvtxparams;
114  else
115  return 0;
116  }
117 
118  template < class Entity >
119  int numParameters ( const Entity & ) const
120  {
121  return numParameters< Entity::codimension >();
122  }
123 
124  std::vector<double>& parameter(const Element &element)
125  {
126  const typename Element::Geometry &geo = element.geometry();
127  DomainType coord( geo.corner( 0 ) );
128  for( int i = 1; i < geo.corners(); ++i )
129  coord += geo.corner( i );
130  coord /= double( geo.corners() );
131 
132  InsertOrderIterator it = elInsertOrder_.find( coord );
133  if( it != elInsertOrder_.end() )
134  return macroGrid_.elParams[ it->second ];
135  assert(0);
136  return emptyParam;
137  }
138 
139  std::vector<double>& parameter(const Vertex &vertex)
140  {
141  const typename Vertex::Geometry &geo = vertex.geometry();
142  DomainType coord( geo.corner( 0 ) );
143 
144  InsertOrderIterator it = vtxInsertOrder_.find( coord );
145  if( it != vtxInsertOrder_.end() )
146  return macroGrid_.vtxParams[ it->second ];
147  return emptyParam;
148  }
149 
150  // return true if boundary parameters found
152  {
153  return false;
154  }
155 
156  template< class GG, class II >
157  const typename DGFBoundaryParameter::type &
158  boundaryParameter ( const Intersection< GG, II > & intersection ) const
159  {
161  }
162 
163  private:
164  typedef FieldVector<typename Grid::ctype,Grid::dimensionworld> DomainType;
165  struct Compare
166  {
167  bool operator() ( const DomainType &a, const DomainType &b ) const
168  {
169  // returns true, if a < b; c[i] < -eps;
170  const DomainType c = a - b;
171  const double eps = 1e-8;
172 
173  for( int i = 0; i < DomainType::dimension; ++i )
174  {
175  if( c[ i ] <= -eps )
176  return true;
177  if( c[ i ] >= eps )
178  return false;
179  }
180  return false;
181  }
182  };
183  typedef std::map< DomainType, size_t, Compare > InsertOrderMap;
184  typedef typename InsertOrderMap::const_iterator InsertOrderIterator;
185 
186  MacroGrid macroGrid_;
187  Grid *grid_;
188  InsertOrderMap elInsertOrder_;
189  InsertOrderMap vtxInsertOrder_;
190  std::vector<double> emptyParam;
191  };
192 
193 } // end namespace Dune
194 
195 #endif
Include standard header files.
Definition: agrid.hh:60
@ vertex
Definition: common.hh:133
Definition: dgfgridfactory.hh:38
int numParameters() const
Definition: dgfgridfactory.hh:108
const DGFBoundaryParameter::type & boundaryParameter(const Intersection< GG, II > &intersection) const
Definition: dgfgridfactory.hh:158
int boundaryId(const Intersection &intersection) const
Definition: dgfgridfactory.hh:102
MPIHelper::MPICommunicator MPICommunicatorType
Definition: dgfgridfactory.hh:41
int numParameters(const Entity &) const
Definition: dgfgridfactory.hh:119
bool wasInserted(const Intersection &intersection) const
Definition: dgfgridfactory.hh:96
DGFGridFactory(const std::string &filename, MPICommunicatorType comm=MPIHelper::getCommunicator())
Definition: dgfgridfactory.hh:50
std::vector< double > & parameter(const Vertex &vertex)
Definition: dgfgridfactory.hh:139
std::vector< double > & parameter(const Element &element)
Definition: dgfgridfactory.hh:124
G Grid
Definition: dgfgridfactory.hh:39
Grid * grid()
Definition: dgfgridfactory.hh:90
static const int dimension
Definition: dgfgridfactory.hh:40
bool haveBoundaryParameters() const
Definition: dgfgridfactory.hh:151
Intersection of a mesh entity of codimension 0 ("element") with a "neighboring" element or with the d...
Definition: common/intersection.hh:164
bool boundary() const
Return true if intersection is with interior or exterior boundary (see the cases above)
Definition: common/intersection.hh:216
int indexInInside() const
Local index of codim 1 entity in the inside() entity where intersection is contained in.
Definition: common/intersection.hh:346
Wrapper class for entities.
Definition: common/entity.hh:66
constexpr static int dimension
The dimension of the grid.
Definition: common/grid.hh:387
static const type & defaultValue()
default constructor
Definition: parser.hh:28
std::string type
type of additional boundary parameters
Definition: parser.hh:25
int nofvtxparams
Definition: parser.hh:163
std::vector< std::vector< double > > vtxParams
Definition: parser.hh:165
int nofelparams
Definition: parser.hh:163
std::vector< std::vector< double > > elParams
Definition: parser.hh:165
std::vector< std::vector< double > > vtx
Definition: parser.hh:125
std ::vector< std ::vector< unsigned int > > elements
Definition: parser.hh:134