dune-fem 2.8.0
Loading...
Searching...
No Matches
gridwidth.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_GRIDWIDTH_HH
2#define DUNE_FEM_GRIDWIDTH_HH
3
4//- system includes
5#include <limits>
6
7//- Dune includes
8#include <dune/common/fvector.hh>
9#include <dune/grid/common/capabilities.hh>
14#include <dune/common/binaryfunctions.hh>
15
16namespace Dune
17{
18
19 namespace Fem
20 {
21
22 template <class GridPartImp, class MinMax = Min<double> >
23 class GridWidthProvider;
24
27 {
28 protected:
29 template< class GridPartType, class EntityType,
30 class ElemGeoInfoType, class FaceGeoInfoType >
31 static inline double
32 maxEdgeWidth ( const GridPartType & gridPart,
33 const EntityType& entity,
34 const ElemGeoInfoType &geoInfo,
35 const FaceGeoInfoType &faceGeoInfo )
36 {
37 typedef typename GridPartType :: IntersectionIteratorType IntersectionIteratorType;
38 double faceVol = std::numeric_limits<double>::max() ;
39 int indexInInside = -1;
40 double currVol = std::numeric_limits<double>::min() ;
41 double refFaceVol = std::numeric_limits<double>::min() ;
42
43 const IntersectionIteratorType endit = gridPart.iend( entity );
44 for( IntersectionIteratorType it = gridPart.ibegin( entity );
45 it != endit; ++it)
46 {
47 typedef typename IntersectionIteratorType::Intersection Intersection;
48 const Intersection &intersection = *it;
49
50 typedef typename Intersection::Geometry LocalGeometryType;
51 const LocalGeometryType &interGeo = intersection.geometry();
52
53 const int localIndex = intersection.indexInInside();
54 // calculate face Volume also for non-conforming intersections
55 if( indexInInside != localIndex )
56 {
57 if( indexInInside >= 0 )
58 faceVol = std::min( faceVol, currVol*refFaceVol );
59
60 refFaceVol = faceGeoInfo.referenceVolume( intersection.type() );
61 indexInInside = localIndex;
62 currVol = 0.0;
63 }
64
65 currVol += interGeo.volume();
66 }
67
68 // do check for last set of intersections
69 if( indexInInside >= 0 )
70 faceVol = std::min( faceVol, currVol*refFaceVol );
71
72 const double elemVol = entity.geometry().volume()
73 * geoInfo.referenceVolume( entity.type() );
74 return elemVol / faceVol;
75 }
76
77 public:
78 template <class MinMax>
79 struct MinMaxInit ;
80
81 template < class T >
82 struct MinMaxInit< Min< T > >
83 {
84 static T init ()
85 {
86 return std::numeric_limits< T >::max() ;
87 }
88 };
89
90 template < class T >
91 struct MinMaxInit< Max< T > >
92 {
93 static T init ()
94 {
95 return std::numeric_limits< T >::min() ;
96 }
97 };
98
107 template <class GridPartType, class MinMax>
108 static inline
109 double calcGridWidth (const GridPartType & gridPart,
110 const bool communicate,
111 const MinMax& minmax
112 )
113 {
114 typedef typename GridPartType :: GridType GridType;
115 typedef typename GridPartType :: template Codim<0> :: IteratorType IteratorType;
116
117 // get geo infor for elements
119 GeomInfoType geoInfo( gridPart.indexSet() );
120
121 // get geo infor for faces
122 typedef GeometryInformation< GridType , 1 > FaceGeometryInformationType;
123 FaceGeometryInformationType faceGeoInfo( geoInfo.geomTypes(1) );
124
125 // initialize with big value
126 double width = MinMaxInit< MinMax > :: init() ;
127
128 // cartesian case
129 if( Dune::Capabilities::isCartesian<GridType>::v &&
131 {
132 // here we only need to check one element
133 IteratorType it = gridPart.template begin<0> ();
134 if( it != gridPart.template end<0> () )
135 {
136 width = minmax( maxEdgeWidth(gridPart,
137 *it,
138 geoInfo,
139 faceGeoInfo ),
140 width );
141 }
142 }
143 else
144 {
145 // unstructured case
146 const IteratorType endit = gridPart.template end<0> ();
147 for(IteratorType it = gridPart.template begin<0> ();
148 it != endit; ++it )
149 {
150 width = minmax( maxEdgeWidth(gridPart,
151 *it,
152 geoInfo,
153 faceGeoInfo ),
154 width );
155 }
156 }
157
158 // calculate global minimum
159 if( communicate )
160 {
161 double w = width ;
162 gridPart.comm().template allreduce<MinMax> ( &w, &width, 1 );
163 }
164
165 return width;
166 }
167 template <class GridPartType>
168 static inline
169 double calcGridWidth (const GridPartType & gridPart,
170 const bool communicate = true)
171 {
172 return calcGridWidth ( gridPart, communicate, Min<double>() );
173 }
174 };
175
177 template <class GridType, class MinMax >
179 {
180 typedef GridWidthProvider < GridType, MinMax > ThisType;
181 public:
184
185 protected:
187
189
190 const GridType& grid_;
192
194
195 mutable double width_;
196 mutable int sequence_;
197
199 public:
201 GridWidthProvider(const GridType* grid)
202 : grid_( *grid )
203 , dm_( DofManagerType :: instance( grid_ ))
204 , gridPart_( const_cast<GridType& > (grid_) )
205 , width_(-1.0)
206 , sequence_(-1)
207 {
208 }
209
211 double gridWidth () const
212 {
213 calcWidths();
214 return width_;
215 }
216
217 protected:
218 void calcWidths() const
219 {
220#ifndef NDEBUG
221 // make sure grid width calculation is done on every process
222 int check = (dm_.sequence() != sequence_) ? 1 : 0;
223 int willCalc = gridPart_.comm().min( check );
224 assert( check == willCalc );
225#endif
226
227 if( dm_.sequence() != sequence_ )
228 {
229 // calculate grid width
230 width_ = GridWidth :: calcGridWidth ( gridPart_ , true );
231
232 assert( width_ > 0 );
234 }
235 }
236 };
237
238 } // namespace Fem
239
240} // namespace Dune
241
242#endif
int sequence() const
return number of sequence, if dofmanagers memory was changed by calling some method like resize,...
Definition: dofmanager.hh:978
double min(const Dune::Fem::Double &v, const double p)
Definition: double.hh:953
Definition: bindguard.hh:11
const CollectiveCommunicationType & comm() const
obtain collective communication object
Definition: gridview2gridpart.hh:193
Definition: misc/capabilities.hh:151
utility functions for calculating the maximum grid width
Definition: gridwidth.hh:179
double gridWidth() const
return characteristic grid width
Definition: gridwidth.hh:211
GridPartType gridPart_
Definition: gridwidth.hh:193
GridWidthProvider(const ThisType &)
int sequence_
Definition: gridwidth.hh:196
SingletonList< const GridType *, ThisType > ProviderType
type of singleton provider
Definition: gridwidth.hh:183
void calcWidths() const
Definition: gridwidth.hh:218
DofManager< GridType > DofManagerType
Definition: gridwidth.hh:186
LeafGridPart< GridType > GridPartType
Definition: gridwidth.hh:188
const GridType & grid_
Definition: gridwidth.hh:190
const DofManagerType & dm_
Definition: gridwidth.hh:191
GridWidthProvider(const GridType *grid)
constructor taking grid part
Definition: gridwidth.hh:201
double width_
Definition: gridwidth.hh:195
utility functions for calculating the maximum grid width
Definition: gridwidth.hh:27
static double maxEdgeWidth(const GridPartType &gridPart, const EntityType &entity, const ElemGeoInfoType &geoInfo, const FaceGeoInfoType &faceGeoInfo)
Definition: gridwidth.hh:32
static double calcGridWidth(const GridPartType &gridPart, const bool communicate, const MinMax &minmax)
calculate minimal grid width as with .
Definition: gridwidth.hh:109
static double calcGridWidth(const GridPartType &gridPart, const bool communicate=true)
Definition: gridwidth.hh:169
Definition: gridwidth.hh:79
static T init()
Definition: gridwidth.hh:84
static T init()
Definition: gridwidth.hh:93
ReferenceVolume and local bary center keeper class.
Definition: allgeomtypes.hh:30
default implementation uses method geomTypes of given index set. Used in DiscreteFunctionSpaces.
Definition: allgeomtypes.hh:99
Definition: dofmanager.hh:761
Singleton list for key/object pairs.
Definition: singletonlist.hh:53