Line data Source code
1 : /*
2 : * Copyright (c) 2010-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 : #include <vector>
34 : #include "boost/container/vector.hpp" // for bool-vectors
35 : #include "regular_refinement.h"
36 : #include "lib_grid/algorithms/geom_obj_util/geom_obj_util.h"
37 : #include "lib_grid/algorithms/selection_util.h"
38 : #include "lib_grid/algorithms/debug_util.h"
39 : #include "lib_grid/callbacks/selection_callbacks.h"
40 :
41 : using namespace std;
42 :
43 : namespace ug
44 : {
45 :
46 0 : bool Refine(Grid& grid, Selector& sel,
47 : RefinementProjector* projector,
48 : bool useSnapPoints)
49 : {
50 : AInt aInt;
51 0 : if(grid.num<Face>() > 0)
52 : grid.attach_to_edges(aInt);
53 0 : if(grid.num<Volume>() > 0)
54 : grid.attach_to_faces(aInt);
55 :
56 0 : bool bSuccess = Refine(grid, sel, aInt, projector, useSnapPoints);
57 :
58 0 : if(grid.num<Face>() > 0)
59 : grid.detach_from_edges(aInt);
60 0 : if(grid.num<Volume>() > 0)
61 : grid.detach_from_faces(aInt);
62 :
63 0 : return bSuccess;
64 : }
65 :
66 0 : static void AdjustSelection(Grid& grid, Selector& sel)
67 : {
68 : // select all edges of selected faces
69 : vector<Edge*> vEdges;
70 0 : for(FaceIterator iter = sel.begin<Face>(); iter != sel.end<Face>(); ++iter){
71 0 : CollectEdges(vEdges, grid, *iter);
72 0 : for(size_t i = 0; i < vEdges.size(); ++i)
73 0 : sel.select(vEdges[i]);
74 : }
75 :
76 : // select all edges of selected volumes
77 : for(VolumeIterator iter = sel.begin<Volume>();
78 0 : iter != sel.end<Volume>(); ++iter)
79 : {
80 0 : CollectEdges(vEdges, grid, *iter);
81 0 : for(size_t i = 0; i < vEdges.size(); ++i)
82 0 : sel.select(vEdges[i]);
83 : }
84 :
85 : // select all faces and volumes which are adjacent to selected edges
86 : vector<Face*> vFaces;
87 : vector<Volume*> vVols;
88 : for(EdgeIterator iter = sel.begin<Edge>();
89 0 : iter != sel.end<Edge>(); ++iter)
90 : {
91 0 : CollectFaces(vFaces, grid, *iter);
92 0 : for(size_t i = 0; i < vFaces.size(); ++i)
93 0 : sel.select(vFaces[i]);
94 :
95 0 : CollectVolumes(vVols, grid, *iter);
96 0 : for(size_t i = 0; i < vVols.size(); ++i)
97 0 : sel.select(vVols[i]);
98 : }
99 0 : }
100 :
101 :
102 0 : bool Refine(Grid& grid, Selector& sel, AInt& aInt,
103 : RefinementProjector* projector,
104 : bool useSnapPoints)
105 : {
106 : // position data is required
107 0 : if(!grid.has_vertex_attachment(aPosition)){
108 : LOG(" WARNING in Refine: aPosition is not attached to the vertices of the grid. Aborting.\n");
109 0 : return false;
110 : }
111 :
112 : // if the refinement-callback is empty, use a linear one.
113 : RefinementProjector defaultProjector;
114 : if(!projector){
115 0 : defaultProjector.set_geometry(make_sp(new Geometry<3, 3>(grid, aPosition)));
116 : projector = &defaultProjector;
117 : }
118 :
119 : // make sure that GRIDOPT_VERTEXCENTRIC_INTERCONNECTION is enabled
120 0 : if(grid.num_edges() && (!grid.option_is_enabled(VRTOPT_STORE_ASSOCIATED_EDGES))){
121 : LOG(" INFO in Refine: autoenabling VRTOPT_STORE_ASSOCIATED_EDGES\n");
122 0 : grid.enable_options(VRTOPT_STORE_ASSOCIATED_EDGES);
123 : }
124 0 : if(grid.num_faces() && (!grid.option_is_enabled(VRTOPT_STORE_ASSOCIATED_FACES))){
125 : LOG(" INFO in Refine: autoenabling VRTOPT_STORE_ASSOCIATED_FACES\n");
126 0 : grid.enable_options(VRTOPT_STORE_ASSOCIATED_FACES);
127 : }
128 0 : if(grid.num_volumes() && (!grid.option_is_enabled(VRTOPT_STORE_ASSOCIATED_VOLUMES))){
129 : LOG(" INFO in Refine: autoenabling VRTOPT_STORE_ASSOCIATED_VOLUMES\n");
130 0 : grid.enable_options(VRTOPT_STORE_ASSOCIATED_VOLUMES);
131 : }
132 :
133 : // make sure that FACEOPT_AUTOGENERATE_EDGES is enabled
134 0 : if(grid.num<Face>() > 0 && (!grid.option_is_enabled(FACEOPT_AUTOGENERATE_EDGES))){
135 : LOG(" INFO in Refine: autoenabling FACEOPT_AUTOGENERATE_EDGES\n");
136 0 : grid.enable_options(FACEOPT_AUTOGENERATE_EDGES);
137 : }
138 :
139 : // if there are volumes, make sure that VOLOPT_AUTOGENERATE_FACES is enabled.
140 0 : if(grid.num<Volume>() > 0 && (!grid.option_is_enabled(VOLOPT_AUTOGENERATE_FACES)))
141 : {
142 : LOG(" INFO in Refine: autoenabling VOLOPT_AUTOGENERATE_FACES\n");
143 0 : grid.enable_options(VOLOPT_AUTOGENERATE_FACES);
144 : }
145 :
146 : // snap-points are recorded here, since we alter the vertex selection and have
147 : // to restore them later on.
148 : vector<Vertex*> snapPoints;
149 0 : if(useSnapPoints)
150 0 : snapPoints.insert(snapPoints.end(), sel.begin<Vertex>(), sel.end<Vertex>());
151 :
152 : // adjust selection
153 0 : AdjustSelection(grid, sel);
154 :
155 : // we will select associated vertices, too, since we have to
156 : // notify the refinement-callback, that they are involved in refinement.
157 0 : sel.clear<Vertex>();
158 0 : SelectAssociatedVertices(sel, sel.begin<Edge>(), sel.end<Edge>(), ISelector::SELECTED);
159 0 : SelectAssociatedVertices(sel, sel.begin<Face>(), sel.end<Face>(), ISelector::SELECTED);
160 0 : SelectAssociatedVertices(sel, sel.begin<Volume>(), sel.end<Volume>(), ISelector::SELECTED);
161 :
162 : // select snap-vertices with a special mark
163 : const ISelector::status_t snapSelVal = ISelector::SELECTED + 1;
164 0 : if(useSnapPoints){
165 0 : for(vector<Vertex*>::iterator iter = snapPoints.begin();
166 0 : iter != snapPoints.end(); ++iter)
167 : {
168 0 : sel.select(*iter, snapSelVal);
169 : }
170 : }
171 :
172 : // aInt has to be attached to the edges of the grid
173 0 : if(sel.num<Face>() > 0 && (!grid.has_edge_attachment(aInt))){
174 : LOG(" WARNING in Refine: aInt is not attached to the edges of the grid. Aborting.\n");
175 : return false;
176 : }
177 :
178 : // if there are selected volumes,
179 : // aInt has to be attached to the faces of the grid
180 0 : if(sel.num<Volume>() && (!grid.has_face_attachment(aInt))){
181 : LOG(" WARNING in Refine: aInt is not attached to the faces of the grid. Aborting.\n");
182 : return false;
183 : }
184 :
185 : // number of edges, faces and volumes that will be refined
186 : const size_t numRefEdges = sel.num<Edge>();
187 : const size_t numRefFaces = sel.num<Face>();
188 : const size_t numRefVols = sel.num<Volume>();
189 :
190 : // we need several arrays.
191 : // this one stores pointers to all edges that shall be refined
192 0 : vector<Edge*> edges(numRefEdges);
193 : // one that stores the new vertex for each edge that shall be refined.
194 0 : vector<RegularVertex*> edgeVrts(numRefEdges);
195 : // one that stores the selected faces
196 0 : vector<Face*> faces(numRefFaces);
197 : // one that stores vertices which are created on faces
198 : vector<Vertex*> faceVrts;
199 0 : if(numRefVols > 0)
200 0 : faceVrts.resize(numRefFaces);
201 : // one that stores selected volumes
202 0 : vector<Volume*> vols(sel.num<Volume>());
203 :
204 : // acces the int-attachment
205 : Grid::EdgeAttachmentAccessor<AInt> aaIntEDGE;
206 : Grid::FaceAttachmentAccessor<AInt> aaIntFACE;
207 :
208 0 : if(numRefFaces > 0)
209 : aaIntEDGE.access(grid, aInt);
210 0 : if(numRefVols > 0)
211 : aaIntFACE.access(grid, aInt);
212 :
213 :
214 : // notify refinement projector about the start of refinement
215 0 : SubGrid<IsSelected> sg(sel.get_grid_objects(), IsSelected(sel));
216 0 : projector->refinement_begins(&sg);
217 :
218 : // store the geometry
219 0 : IGeometry3d& geom = *projector->geometry();
220 :
221 : ////////////////////////////////
222 : // fill the edges- and edgeVrts-array and assign indices to selected edges
223 : {
224 : EdgeIterator edgesEnd = sel.end<Edge>();
225 : int i = 0;
226 : for(EdgeIterator iter = sel.begin<Edge>();
227 0 : iter != edgesEnd; ++iter, ++i)
228 : {
229 : // store the edge
230 0 : edges[i] = *iter;
231 0 : if(numRefFaces > 0)
232 0 : aaIntEDGE[*iter] = i;
233 :
234 : // create the new vertex
235 0 : edgeVrts[i] = *grid.create<RegularVertex>(*iter);
236 : sel.select(edgeVrts[i]);
237 : // calculate new position
238 0 : projector->new_vertex(edgeVrts[i], *iter);
239 : }
240 : }
241 :
242 : ////////////////////////////////
243 : // refine the selected edges
244 : vector<Edge*> newEdges;
245 0 : newEdges.reserve(2);
246 0 : for(size_t i = 0; i < edges.size(); ++i){
247 0 : Edge* e = edges[i];
248 0 : if(e->refine(newEdges, edgeVrts[i])){
249 0 : for(size_t j = 0; j < newEdges.size(); ++j)
250 0 : grid.register_element(newEdges[j], e);
251 : }
252 : else{
253 : LOG(" WARNING in Refine: could not refine edge.\n");
254 : }
255 : }
256 :
257 : ////////////////////////////////
258 : // set up face arrays
259 : {
260 : FaceIterator facesEnd = sel.end<Face>();
261 : int i = 0;
262 : for(FaceIterator iter = sel.begin<Face>();
263 0 : iter != facesEnd; ++iter, ++i)
264 : {
265 : Face* f = *iter;
266 0 : faces[i] = f;
267 0 : if(numRefVols > 0)
268 0 : aaIntFACE[f] = i;
269 : }
270 : }
271 :
272 : ////////////////////////////////
273 : // refine the selected faces
274 : vector<Face*> newFaces;
275 0 : newFaces.reserve(4);
276 : // we need a container that stores the vertex for each edge of a face
277 : // entries will be set to NULL if the associated edge will not be refined
278 : vector<Vertex*> faceEdgeVrts;
279 0 : faceEdgeVrts.reserve(4);
280 :
281 0 : for(size_t i = 0; i < faces.size(); ++i){
282 0 : Face* f = faces[i];
283 : Vertex* newVrt;
284 :
285 : // collect vertices of associated edges
286 : faceEdgeVrts.clear();
287 0 : for(uint j = 0; j < f->num_edges(); ++j){
288 0 : Edge* e = grid.get_edge(f, j);
289 0 : if(sel.is_selected(e))
290 0 : faceEdgeVrts.push_back(edgeVrts[aaIntEDGE[e]]);
291 : else
292 0 : faceEdgeVrts.push_back(NULL);
293 : }
294 :
295 : int snapPointIndex = -1;
296 0 : if(useSnapPoints){
297 0 : Face::ConstVertexArray vrts = f->vertices();
298 0 : const size_t numVrts = f->num_vertices();
299 0 : for(size_t s = 0; s < numVrts; ++s){
300 0 : if(sel.get_selection_status(vrts[s]) == snapSelVal){
301 0 : if(snapPointIndex != -1){
302 : UG_LOG("WARNING: Only one snap-point per face is allowed, "
303 : "but more are present. Ignoring snap-points for this face.\n");
304 : snapPointIndex = -1;
305 : break;
306 : }
307 0 : snapPointIndex = static_cast<int>(s);
308 : }
309 : }
310 : }
311 :
312 0 : if(f->refine(newFaces, &newVrt, &faceEdgeVrts.front(), NULL, NULL, snapPointIndex)){
313 : // if a new vertex was generated, we have to register it
314 0 : if(newVrt){
315 : grid.register_element(newVrt, f);
316 0 : projector->new_vertex(newVrt, f);
317 0 : sel.select(newVrt);
318 : }
319 :
320 : // if volumes are refined too, we have to store the vertex
321 0 : if(numRefVols > 0)
322 0 : faceVrts[i] = newVrt;
323 :
324 : // register the new faces
325 0 : for(size_t j = 0; j < newFaces.size(); ++j)
326 0 : grid.register_element(newFaces[j], f);
327 : }
328 : else{
329 : LOG(" WARNING in Refine: could not refine face.\n");
330 : }
331 : }
332 :
333 : ////////////////////////////////
334 : // set up volume arrays
335 : {
336 : VolumeIterator volsEnd = sel.end<Volume>();
337 : int i = 0;
338 : for(VolumeIterator iter = sel.begin<Volume>();
339 0 : iter != volsEnd; ++iter, ++i)
340 : {
341 : Volume* v = *iter;
342 0 : vols[i] = v;
343 : }
344 : }
345 :
346 : ////////////////////////////////
347 : // refine the selected volumes
348 : vector<Volume*> newVols;
349 0 : newVols.reserve(8);
350 : // we need a container that stores the vertex for each edge of a volume
351 : // entries will be set to NULL if the associated edge will not be refined
352 : vector<Vertex*> volEdgeVrts;
353 0 : volEdgeVrts.reserve(12);
354 : // we need a container that stores the vertex for each face of a volume
355 : // entries will be set to NULL if the associated face will not be refined
356 : vector<Vertex*> volFaceVrts;
357 0 : volFaceVrts.reserve(6);
358 :
359 : // only used for tetrahedron refinement
360 0 : vector<vector3> corners(4, vector3(0, 0, 0));
361 :
362 : // // DEBUG
363 : // UG_LOG("> VOL-REF-BEGIN\n");
364 : // UG_LOG("> DEBUG-ACCESSOR...\n");
365 : // Grid::VertexAttachmentAccessor<APosition> aaPos(grid, aPosition);
366 :
367 : boost::container::vector<bool> isSnapPoint;
368 0 : if(useSnapPoints)
369 : isSnapPoint.reserve(8);
370 :
371 0 : for(size_t i = 0; i < vols.size(); ++i){
372 0 : Volume* v = vols[i];
373 : Vertex* newVrt;
374 :
375 : // collect vertices of associated edges
376 : volEdgeVrts.clear();
377 0 : for(uint j = 0; j < v->num_edges(); ++j){
378 0 : Edge* e = grid.get_edge(v, j);
379 0 : if(sel.is_selected(e))
380 0 : volEdgeVrts.push_back(edgeVrts[aaIntEDGE[e]]);
381 : else
382 0 : volEdgeVrts.push_back(NULL);
383 : }
384 :
385 : // collect vertices of associated faces
386 : volFaceVrts.clear();
387 0 : for(uint j = 0; j < v->num_faces(); ++j){
388 0 : Face* f = grid.get_face(v, j);
389 0 : if(sel.is_selected(f))
390 0 : volFaceVrts.push_back(faceVrts[aaIntFACE[f]]);
391 : else
392 0 : volFaceVrts.push_back(NULL);
393 : }
394 :
395 : // if we're performing tetrahedral refinement, we have to collect
396 : // the corner coordinates, so that the refinement algorithm may choose
397 : // the best interior diagonal.
398 : vector3* pCorners = NULL;
399 0 : if((v->num_vertices() == 4)){
400 0 : for(size_t i = 0; i < 4; ++i){
401 0 : corners[i] = geom.pos(v->vertex(i));
402 : }
403 : pCorners = &corners.front();
404 : }
405 :
406 : bool* pIsSnapPoint = NULL;
407 0 : if(useSnapPoints){
408 0 : Volume::ConstVertexArray vrts = v->vertices();
409 0 : const size_t numVrts = v->num_vertices();
410 : isSnapPoint.clear();
411 :
412 : bool gotOne = false;
413 0 : for(size_t s = 0; s < numVrts; ++s){
414 0 : const bool val = sel.get_selection_status(vrts[s]) == snapSelVal;
415 : isSnapPoint.push_back(val);
416 0 : gotOne = gotOne | val;
417 : }
418 :
419 0 : if(gotOne)
420 : pIsSnapPoint = &isSnapPoint.front();
421 : }
422 :
423 0 : if(v->refine(newVols, &newVrt, &volEdgeVrts.front(),
424 0 : &volFaceVrts.front(), NULL, RegularVertex(), NULL,
425 : pCorners, pIsSnapPoint))
426 : {
427 : // if a new vertex was generated, we have to register it
428 0 : if(newVrt){
429 : grid.register_element(newVrt, v);
430 0 : projector->new_vertex(newVrt, v);
431 0 : sel.select(newVrt);
432 : }
433 :
434 : // register the new volumes
435 0 : for(size_t j = 0; j < newVols.size(); ++j)
436 0 : grid.register_element(newVols[j], v);
437 : }
438 : else{
439 : LOG(" WARNING in Refine: could not refine volume.\n");
440 : }
441 : }
442 : // UG_LOG("> VOL-REF-END\n");
443 :
444 : // erase old volumes
445 : grid.erase(vols.begin(), vols.end());
446 : // erase old faces
447 : grid.erase(faces.begin(), faces.end());
448 : // erase old edges
449 : grid.erase(edges.begin(), edges.end());
450 :
451 : // notify refinement projector about the end of refinement
452 0 : projector->refinement_ends();
453 : return true;
454 0 : }
455 :
456 : }// end of namespace
|
