dune-grid-glue 2.5-git
codim0extractor.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/*
4 * Filename: codim0extractor.hh
5 * Version: 1.0
6 * Created on: Jun 23, 2009
7 * Author: Oliver Sander, Christian Engwer
8 * ---------------------------------
9 * Project: dune-grid-glue
10 * Description: base class for grid extractors extracting surface grids
11 *
12 */
18#ifndef DUNE_GRIDGLUE_EXTRACTORS_CODIM0EXTRACTOR_HH
19#define DUNE_GRIDGLUE_EXTRACTORS_CODIM0EXTRACTOR_HH
20
21#include <deque>
22#include <functional>
23
24#include <dune/common/deprecated.hh>
25#include <dune/grid/common/mcmgmapper.hh>
26
27#include "extractor.hh"
28#include "extractorpredicate.hh"
29
30namespace Dune {
31
32 namespace GridGlue {
33
34template<typename GV>
35class Codim0Extractor : public Extractor<GV,0>
36{
37
38public:
39
40 /* E X P O R T E D T Y P E S A N D C O N S T A N T S */
41 using Extractor<GV,0>::codim;
42 typedef typename Extractor<GV,0>::ctype ctype;
43 using Extractor<GV,0>::dim;
44 using Extractor<GV,0>::dimworld;
46
47 typedef typename GV::Traits::template Codim<dim>::Entity Vertex;
48 typedef typename GV::Traits::template Codim<0>::Entity Element;
49 typedef std::function<bool(const Element&, unsigned int subentity)> Predicate;
50
51 static const Dune::PartitionIteratorType PType = Dune::Interior_Partition;
52 typedef typename GV::Traits::template Codim<0>::template Partition<PType>::Iterator ElementIter;
53
54 // import typedefs from base class
60
66 DUNE_DEPRECATED_MSG("Please use a std::function<bool(const Element&, unsigned int)> in favor of the ExtractorPredicate.")
67 Codim0Extractor(const GV& gv, const ExtractorPredicate<GV,0>& descr)
68 : Extractor<GV,0>(gv), positiveNormalDirection_(false)
69 {
70 std::cout << "This is Codim0Extractor on a <"
71 << GV::dimension << "," << GV::dimensionworld << "> grid!" << std::endl;
72 const auto predicate = [&](const Element& element, unsigned int subentity) -> bool {
73 return descr.contains(element, subentity);
74 };
75 update(predicate);
76 }
77
83 Codim0Extractor(const GV& gv, const Predicate& predicate)
84 : Extractor<GV,0>(gv), positiveNormalDirection_(false)
85 {
86 std::cout << "This is Codim0Extractor on a <"
87 << GV::dimension << "," << GV::dimensionworld << "> grid!" << std::endl;
88 update(predicate);
89 }
90
92 const bool & positiveNormalDirection() const { return positiveNormalDirection_; }
93
94protected:
96private:
97 void update(const Predicate& predicate);
98};
99
100
101template<typename GV>
102void Codim0Extractor<GV>::update(const Predicate& predicate)
103{
104 // In this first pass iterate over all entities of codim 0.
105 // Get its corner vertices, find resp. store them together with their associated index,
106 // and remember the indices of the corners.
107
108 // free everything there is in this object
109 this->clear();
110
111 // several counter for consecutive indexing are needed
112 size_t element_index = 0;
113 size_t vertex_index = 0;
114
115 // a temporary container where newly acquired face
116 // information can be stored at first
117 std::deque<SubEntityInfo> temp_faces;
118
119 // iterate over all codim 0 elements on the grid
120 for (ElementIter elit = this->gv_.template begin<0, PType>();
121 elit != this->gv_.template end<0, PType>(); ++elit)
122 {
123 const auto& elmt = *elit;
124 const auto geometry = elmt.geometry();
125 IndexType eindex = this->cellMapper_.index(elmt);
126
127 // only do sth. if this element is "interesting"
128 // implicit cast is done automatically
129 if (predicate(elmt,0))
130 {
131 // add an entry to the element info map, the index will be set properly later
132 this->elmtInfo_[eindex] = new ElementInfo(element_index, elmt, 1);
133
134 unsigned int numCorners = elmt.subEntities(dim);
135 unsigned int vertex_indices[numCorners]; // index in global vector
136 unsigned int vertex_numbers[numCorners]; // index in parent entity
137
138 // try for each of the faces vertices whether it is already inserted or not
139 for (unsigned int i = 0; i < numCorners; ++i)
140 {
141 vertex_numbers[i] = i;
142
143 // get the vertex pointer and the index from the index set
144 const Vertex vertex = elit->template subEntity<dim>(vertex_numbers[i]);
145 IndexType vindex = this->gv_.indexSet().template index<dim>(vertex);
146
147 // if the vertex is not yet inserted in the vertex info map
148 // it is a new one -> it will be inserted now!
149 typename VertexInfoMap::iterator vimit = this->vtxInfo_.find(vindex);
150 if (vimit == this->vtxInfo_.end())
151 {
152 // insert into the map
153 this->vtxInfo_[vindex] = new VertexInfo(vertex_index, vertex);
154 // remember this vertex' index
155 vertex_indices[i] = vertex_index;
156 // increase the current index
157 vertex_index++;
158 }
159 else
160 {
161 // only remember the vertex' index
162 vertex_indices[i] = vimit->second->idx;
163 }
164 }
165
166 // flip cell if necessary
167 {
168 switch (int(dim))
169 {
170 case 0 :
171 break;
172 case 1 :
173 {
174 // The following test only works if the zero-th coordinate is the
175 // one that defines the orientation. A sufficient condition for
176 // this is dimworld == 1
177 /* assert(dimworld==1); */
178 bool elementNormalDirection =
179 (geometry.corner(1)[0] < geometry.corner(0)[0]);
180 if ( positiveNormalDirection_ != elementNormalDirection )
181 {
182 std::swap(vertex_indices[0], vertex_indices[1]);
183 std::swap(vertex_numbers[0], vertex_numbers[1]);
184 }
185 break;
186 }
187 case 2 :
188 {
189 Dune::FieldVector<ctype, dimworld>
190 v0 = geometry.corner(1),
191 v1 = geometry.corner(2);
192 v0 -= geometry.corner(0);
193 v1 -= geometry.corner(0);
194 ctype normal_sign = v0[0]*v1[1] - v0[1]*v1[0];
195 bool elementNormalDirection = (normal_sign < 0);
196 if ( positiveNormalDirection_ != elementNormalDirection )
197 {
198 std::cout << "swap\n";
199 if (elmt.type().isCube())
200 {
201 for (int i = 0; i < (1<<dim); i+=2)
202 {
203 // swap i and i+1
204 std::swap(vertex_indices[i], vertex_indices[i+1]);
205 std::swap(vertex_numbers[i], vertex_numbers[i+1]);
206 }
207 } else if (elmt.type().isSimplex()) {
208 std::swap(vertex_indices[0], vertex_indices[1]);
209 std::swap(vertex_numbers[0], vertex_numbers[1]);
210 } else {
211 DUNE_THROW(Dune::Exception, "Unexpected Geometrytype");
212 }
213 }
214 break;
215 }
216 }
217 }
218
219 // add a new face to the temporary collection
220 temp_faces.push_back(SubEntityInfo(eindex,0,elmt.type()));
221 element_index++;
222 for (unsigned int i=0; i<numCorners; i++) {
223 temp_faces.back().corners[i].idx = vertex_indices[i];
224 // remember the vertices' numbers in parent element's vertices
225 temp_faces.back().corners[i].num = vertex_numbers[i];
226 }
227
228 }
229 } // end loop over elements
230
231 // allocate the array for the face specific information...
232 this->subEntities_.resize(element_index);
233 // ...and fill in the data from the temporary containers
234 copy(temp_faces.begin(), temp_faces.end(), this->subEntities_.begin());
235
236 // now first write the array with the coordinates...
237 this->coords_.resize(this->vtxInfo_.size());
238 typename VertexInfoMap::const_iterator it1 = this->vtxInfo_.begin();
239 for (; it1 != this->vtxInfo_.end(); ++it1)
240 {
241 // get a pointer to the associated info object
242 CoordinateInfo* current = &this->coords_[it1->second->idx];
243 // store this coordinates index // NEEDED?
244 current->index = it1->second->idx;
245 // store the vertex' index for the index2vertex mapping
246 current->vtxindex = it1->first;
247 // store the vertex' coordinates under the associated index
248 // in the coordinates array
249 const auto vtx = this->grid().entity(it1->second->p);
250 current->coord = vtx.geometry().corner(0);
251 }
252
253}
254
255} // namespace GridGlue
256
257} // namespace Dune
258
259#endif // DUNE_GRIDGLUE_EXTRACTORS_CODIM0EXTRACTOR_HH
extractor base class
Base class for predicates selecting the part of a grid to be extracted.
Definition gridglue.hh:30
Definition codim0extractor.hh:36
bool & positiveNormalDirection()
Definition codim0extractor.hh:91
Extractor< GV, 0 >::IndexType IndexType
Definition codim0extractor.hh:45
const bool & positiveNormalDirection() const
Definition codim0extractor.hh:92
GV::Traits::template Codim< dim >::Entity Vertex
Definition codim0extractor.hh:47
Extractor< GV, 0 >::CoordinateInfo CoordinateInfo
Definition codim0extractor.hh:58
static const Dune::PartitionIteratorType PType
Definition codim0extractor.hh:51
Extractor< GV, 0 >::VertexInfo VertexInfo
Definition codim0extractor.hh:57
Extractor< GV, 0 >::ctype ctype
Definition codim0extractor.hh:42
GV::Traits::template Codim< 0 >::template Partition< PType >::Iterator ElementIter
Definition codim0extractor.hh:52
bool positiveNormalDirection_
Definition codim0extractor.hh:95
Extractor< GV, 0 >::VertexInfoMap VertexInfoMap
Definition codim0extractor.hh:59
std::function< bool(const Element &, unsigned int subentity)> Predicate
Definition codim0extractor.hh:49
Codim0Extractor(const GV &gv, const Predicate &predicate)
Constructor.
Definition codim0extractor.hh:83
Extractor< GV, 0 >::ElementInfo ElementInfo
Definition codim0extractor.hh:56
GV::Traits::template Codim< 0 >::Entity Element
Definition codim0extractor.hh:48
Extractor< GV, 0 >::SubEntityInfo SubEntityInfo
Definition codim0extractor.hh:55
Provides codimension-independent methods for grid extraction.
Definition extractor.hh:43
Element element(unsigned int index) const
gets the parent element for a given face index, throws an exception if index not valid
Definition extractor.hh:375
@ dim
Definition extractor.hh:48
int IndexType
Definition extractor.hh:74
@ dimworld
Definition extractor.hh:47
GV::Grid::ctype ctype
Definition extractor.hh:59
std::map< IndexType, VertexInfo * > VertexInfoMap
Definition extractor.hh:184
@ codim
Definition extractor.hh:49
simple struct holding a vertex pointer and an index
Definition extractor.hh:117
simple struct holding an element seed and an index
Definition extractor.hh:129
Holds some information about an element's subentity involved in a coupling.
Definition extractor.hh:148
Base class for subentity-selecting predicates.
Definition extractorpredicate.hh:31