ug4
delaunay_info_impl.h
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2015: G-CSC, Goethe University Frankfurt
3  * Author: Sebastian Reiter
4  *
5  * This file is part of UG4.
6  *
7  * UG4 is free software: you can redistribute it and/or modify it under the
8  * terms of the GNU Lesser General Public License version 3 (as published by the
9  * Free Software Foundation) with the following additional attribution
10  * requirements (according to LGPL/GPL v3 §7):
11  *
12  * (1) The following notice must be displayed in the Appropriate Legal Notices
13  * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14  *
15  * (2) The following notice must be displayed at a prominent place in the
16  * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17  *
18  * (3) The following bibliography is recommended for citation and must be
19  * preserved in all covered files:
20  * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21  * parallel geometric multigrid solver on hierarchically distributed grids.
22  * Computing and visualization in science 16, 4 (2013), 151-164"
23  * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24  * flexible software system for simulating pde based models on high performance
25  * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26  *
27  * This program is distributed in the hope that it will be useful,
28  * but WITHOUT ANY WARRANTY; without even the implied warranty of
29  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30  * GNU Lesser General Public License for more details.
31  */
32 
33 #ifndef __H__UG_delaunay_info_impl
34 #define __H__UG_delaunay_info_impl
35 
37 
38 namespace ug{
39 
40 template <class TAAPos>
41 template <class TIter>
43 init_marks(TIter trisBegin, TIter trisEnd, bool pushFlipCandidates)
44 {
45  using namespace std;
46 
47 // first mark all triangles and associated vertices with mark INNER
48  for(TIter iter = trisBegin; iter != trisEnd; ++iter){
49  Face* f = *iter;
50  if(f->num_vertices() != 3)
51  continue;
52 
53  set_mark(f, INNER);
54  for(size_t i = 0; i < 3; ++i){
55  set_mark(f->vertex(i), INNER);
56  }
57  }
58 
59 // Collect all candidates for flips (only edges with two neighbors, both marked).
61  for(TIter triIter = trisBegin; triIter != trisEnd; ++triIter){
62  Face* t = *triIter;
63  m_grid.associated_elements(edges, t);
64  for(size_t i = 0; i < edges.size(); ++i){
65  Edge* e = edges[i];
66 
67  // treat edges only once
68  if(mark(e) != NONE)
69  continue;
70 
71  // mark inner and outer segments
72  if(m_cbConstrainedEdge(e)){
73  set_mark(e, SEGMENT);
74  set_mark(e->vertex(0), SEGMENT);
75  set_mark(e->vertex(1), SEGMENT);
76  }
77  else{
78  Face* nbrFaces[2];
79  int numNbrs = GetAssociatedFaces(nbrFaces, m_grid, e, 2);
80  // two neighbors, both marked
81  if(numNbrs == 2 && is_inner(nbrFaces[0])
82  && is_inner(nbrFaces[1]))
83  {
84  // the edge is a flip candidate
85  if(pushFlipCandidates)
86  push_candidate(e);
87  }
88  else{
89  // the edge lies on the rim and has to be marked as a segment
90  set_mark(e, SEGMENT);
91  set_mark(e->vertex(0), SEGMENT);
92  set_mark(e->vertex(1), SEGMENT);
93  }
94  }
95  }
96  }
97 
98 //todo: iterate over triangles again, this time examining their corners.
99 // (use grid::mark to only visit each once).
100 // Check for each whether it is connected to two segments with a smaller
101 // angle than Pi/3. If this is the case, it will be marked as DART
102  AAPos aaPos = position_accessor();
103 
104  m_grid.begin_marking();
105  for(TIter triIter = trisBegin; triIter != trisEnd; ++triIter){
106  Face* f = *triIter;
107  Face::ConstVertexArray vrts = f->vertices();
108  const size_t numVrts = f->num_vertices();
109  for(size_t ivrt = 0; ivrt < numVrts; ++ivrt){
110  Vertex* vrt = vrts[ivrt];
111  if(m_grid.is_marked(vrt))
112  continue;
113  m_grid.mark(vrt);
114 
115  if(mark(vrt) != SEGMENT)
116  continue;
117 
118  vector_t vrtPos = aaPos[vrt];
119 
120  bool searching = true;
121  m_grid.associated_elements(edges, vrt);
122  for(size_t iedge = 0; (iedge < edges.size()) && searching; ++iedge){
123  Edge* edge = edges[iedge];
124  if(mark(edge) != SEGMENT)
125  continue;
126 
127  vector_t dir1;
128  VecSubtract(dir1, aaPos[GetConnectedVertex(edge, vrt)], vrtPos);
129 
130  // check angles between all segment-edges. If one is found which
131  // is smaller than PI/3, then the vertex will be marked as DART vertex.
132  for(size_t iotherEdge = iedge + 1; iotherEdge < edges.size(); ++iotherEdge){
133  Edge* otherEdge = edges[iotherEdge];
134  if(mark(otherEdge) != SEGMENT)
135  continue;
136 
137  vector_t dir2;
138  VecSubtract(dir2, aaPos[GetConnectedVertex(otherEdge, vrt)], vrtPos);
139 
140  if(VecAngle(dir1, dir2) < PI / 3. + SMALL){
141  searching = false;
142  set_mark(vrt, DART);
143  break;
144  }
145  }
146  }
147  }
148  }
149  m_grid.end_marking();
150 }
151 
152 
153 template <class TAAPos>
154 template <class TElem>
156 is_inner(TElem* e)
157 {
158  Mark m = mark(e);
159  return (m == INNER) || (m == NEW_INNER);
160 }
161 
162 
163 template <class TAAPos>
164 template <class TElem>
166 is_segment(TElem* e)
167 {
168  Mark m = mark(e);
169  return (m == SEGMENT) || (m == NEW_SEGMENT) || (m == DART) || (m == SHELL);
170 }
171 
172 
173 }// end of namespace
174 
175 #endif //__H__UG_delaunay_info_impl
bool is_inner(TElem *e)
returns true if the specified element is either marked as INNER or as NEW_INNER
Definition: delaunay_info_impl.h:156
void init_marks(TIter trisBegin, TIter trisEnd, bool pushFlipCandidates)
Definition: delaunay_info_impl.h:43
TAAPos::ValueType vector_t
Definition: delaunay_info.h:51
Mark
Definition: delaunay_info.h:54
TAAPos AAPos
Definition: delaunay_info.h:50
bool is_segment(TElem *e)
returns true if the specified element is either marked as SEGMENT, NEW_SEGMENT, DART or SHELL
Definition: delaunay_info_impl.h:166
Base-class for edges.
Definition: grid_base_objects.h:397
virtual Vertex * vertex(size_t index) const
Definition: grid_base_objects.h:366
Faces are 2-dimensional objects.
Definition: grid_base_objects.h:510
virtual ConstVertexArray vertices() const
Definition: grid_base_objects.h:487
virtual Vertex * vertex(size_t index) const
Definition: grid_base_objects.h:486
virtual size_t num_vertices() const
Definition: grid_base_objects.h:488
Vertex *const * ConstVertexArray
Definition: grid_base_objects.h:319
Container which holds an array of pointers.
Definition: pointer_const_array.h:84
size_t size() const
returns the size of the associated array.
Definition: pointer_const_array_impl.hpp:106
Base-class for all vertex-types.
Definition: grid_base_objects.h:231
int GetAssociatedFaces(Face **facesOut, Grid &grid, Edge *e, int maxNumFaces)
writes associated faces of e to facesOut.
Definition: edge_util.cpp:186
Vertex * GetConnectedVertex(Edge *e, Vertex *v)
returns the vertex that is connected to v via e.
Definition: vertex_util.cpp:78
void VecSubtract(vector_t &vOut, const vector_t &v1, const vector_t &v2)
subtracts v2 from v1 and stores the result in a vOut
Definition: math_vector_functions_common_impl.hpp:226
vector_t::value_type VecAngle(const vector_t &v1, const vector_t &v2)
returns the angle between two vectors in radiants
Definition: math_vector_functions_common_impl.hpp:401
Definition: smart_pointer.h:814
the ug namespace
const number PI
Definition: math_constants.h:45
const number SMALL
Definition: math_constants.h:41
@ NONE
Definition: data_import.h:44