Line data Source code
1 : /*
2 : * Copyright (c) 2011-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 <cassert>
34 : #include "prism_rules.h"
35 : #include "rule_util.h"
36 : #include "grid_object_ids.h"
37 :
38 : namespace ug{
39 : namespace prism_rules{
40 :
41 : /// Output are the vertices of a prism rotated around its vertical axis
42 0 : void RotatePrism(int vrtsOut[NUM_VERTICES], int steps)
43 : {
44 0 : if(steps < 0)
45 0 : steps = (4 - ((-steps) % 4));
46 :
47 0 : for(int i = 0; i < 3; ++i)
48 0 : vrtsOut[(i + steps) % 3] = i;
49 :
50 0 : for(int i = 3; i < 6; ++i)
51 0 : vrtsOut[3 + (i + steps) % 3] = i;
52 0 : }
53 :
54 : // fips the prism upside down and leaves face0 where it is
55 0 : void FlipPrism(int vrtsOut[NUM_VERTICES], const int vrts[NUM_VERTICES])
56 : {
57 0 : vrtsOut[0] = vrts[5];
58 0 : vrtsOut[1] = vrts[4];
59 0 : vrtsOut[2] = vrts[3];
60 0 : vrtsOut[3] = vrts[2];
61 0 : vrtsOut[4] = vrts[1];
62 0 : vrtsOut[5] = vrts[0];
63 0 : }
64 :
65 :
66 0 : int Refine(int* newIndsOut, int* newEdgeVrts, bool& newCenterOut, vector3*, bool* isSnapPoint)
67 : {
68 0 : newCenterOut = false;
69 : // If a refinement rule is not implemented, fillCount will stay at 0.
70 : // Before returning, we'll check, whether fillCount is at 0 and will
71 : // perform refinement with an additional vertex in this case.
72 :
73 : // corner status is used to mark vertices, which are connected to refined edges
74 : // the value tells, how many associated edges are refined.
75 0 : int cornerStatus[6] = {0, 0, 0, 0, 0, 0};
76 :
77 : // here we'll store the index of each edge, which will be refined.
78 : // Size of the valid sub-array is numNewVrts, which is calculated below.
79 : int refEdgeInds[NUM_EDGES];
80 :
81 : // count the number of new vertices and fill newVrtEdgeInds
82 : int numNewVrts = 0;
83 0 : for(int i = 0; i < NUM_EDGES; ++i){
84 0 : if(newEdgeVrts[i]){
85 0 : refEdgeInds[numNewVrts] = i;
86 0 : ++numNewVrts;
87 :
88 : // adjust corner status of associated vertices
89 0 : const int* evi = EDGE_VRT_INDS[i];
90 0 : ++cornerStatus[evi[0]];
91 0 : ++cornerStatus[evi[1]];
92 : }
93 : }
94 :
95 : // snap-point handling
96 : int numSnapPoints = 0;
97 : int snapPointIndex[NUM_VERTICES];
98 0 : if(isSnapPoint){
99 0 : for(int i = 0; i < NUM_VERTICES; ++i){
100 0 : if(isSnapPoint[i]){
101 0 : snapPointIndex[numSnapPoints] = i;
102 0 : ++numSnapPoints;
103 : }
104 : }
105 : }
106 :
107 : // set to true if the algorithm already tried to find a special rule for the
108 : // specified snap-points. Always true if no snap points were provided.
109 0 : bool snapPointsProcessed = (numSnapPoints == 0);
110 :
111 : // the fillCount tells how much data has already been written to newIndsOut.
112 : int fillCount = 0;
113 :
114 : // convenience - indices where new edge-vrts, new face-vrts and new vol-vrts begin.
115 : const int E = NUM_VERTICES;
116 : const int F = NUM_VERTICES + NUM_EDGES;
117 : // const int V = NUM_VERTICES + NUM_EDGES + NUM_FACES;
118 :
119 : // depending on the number of new vertices, we will now apply different
120 : // refinement rules. Further distinction may have to be done.
121 0 : switch(numNewVrts){
122 0 : case 0:
123 : {
124 : // simply put the default prism back to newIndsOut
125 0 : newIndsOut[fillCount++] = GOID_PRISM;
126 0 : newIndsOut[fillCount++] = 0;
127 0 : newIndsOut[fillCount++] = 1;
128 0 : newIndsOut[fillCount++] = 2;
129 0 : newIndsOut[fillCount++] = 3;
130 0 : newIndsOut[fillCount++] = 4;
131 0 : newIndsOut[fillCount++] = 5;
132 0 : }break;
133 :
134 0 : case 1:
135 : {
136 : // we'll iterate over all faces. If a face does not contain the
137 : // refined edge, then we'll create a tetrahedron or pyramid from it.
138 0 : const int refEdge = refEdgeInds[0];
139 0 : const int nVrt = refEdge + E;
140 :
141 : int& fi = fillCount;
142 : int* inds = newIndsOut;
143 :
144 0 : if(numSnapPoints == 2){
145 0 : const int snapEdge = EDGE_FROM_VRTS[snapPointIndex[0]][snapPointIndex[1]];
146 0 : if(snapEdge != -1
147 0 : && !IS_SIDE_EDGE[snapEdge]
148 0 : && IS_SIDE_EDGE[refEdge]
149 0 : && !isSnapPoint[EDGE_VRT_INDS[refEdge][0]]
150 0 : && !isSnapPoint[EDGE_VRT_INDS[refEdge][1]])
151 : {
152 : // find the quad which contains the snapEdge
153 : int snapQuad = -1;
154 0 : for(int iquad = 0; iquad < 3; ++iquad){
155 0 : if(FACE_CONTAINS_EDGE[QUADS[iquad]][snapEdge]){
156 : snapQuad = QUADS[iquad];
157 : break;
158 : }
159 : }
160 :
161 : UG_ASSERT(snapQuad != -1, "The corresponding quad has to exist!");
162 :
163 : int p[NUM_VERTICES];
164 0 : if(IS_BOTTOM_EDGE[snapEdge]){
165 : int tp[NUM_VERTICES];
166 0 : RotatePrism(tp, 3 - snapQuad);
167 0 : FlipPrism(p, tp);
168 : }
169 : else
170 0 : RotatePrism(p, 3 - snapQuad);
171 :
172 :
173 : // UG_LOG("snapEdge: " << snapEdge << "\n");
174 : // UG_LOG("refEdge: " << refEdge << "\n");
175 : // UG_LOG("snapPointIndex[0]: " << snapPointIndex[0] << "\n");
176 : // UG_LOG("snapPointIndex[1]: " << snapPointIndex[1] << "\n");
177 : // UG_LOG("snapQuad: " << snapQuad << "\n");
178 :
179 : // UG_LOG("p: ");
180 : // for(int ivrt = 0; ivrt < NUM_VERTICES; ++ivrt){
181 : // UG_LOG(p[ivrt] << ", ");
182 : // }
183 : // UG_LOG(std::endl);
184 :
185 :
186 :
187 0 : const int e14 = E + EDGE_FROM_VRTS[p[1]][p[4]];
188 :
189 : // we have to create two new prisms
190 : int& fi = fillCount;
191 : int* inds = newIndsOut;
192 :
193 0 : inds[fi++] = GOID_PRISM;
194 0 : inds[fi++] = p[0]; inds[fi++] = p[1]; inds[fi++] = p[2];
195 0 : inds[fi++] = p[3]; inds[fi++] = e14; inds[fi++] = p[5];
196 :
197 0 : inds[fi++] = GOID_TETRAHEDRON;
198 0 : inds[fi++] = p[3]; inds[fi++] = e14; inds[fi++] = p[5];
199 0 : inds[fi++] = p[4];
200 :
201 : snapPointsProcessed = true;
202 : }
203 : }
204 :
205 0 : if(numSnapPoints == 0 || !snapPointsProcessed){
206 0 : for(int i = 0; i < NUM_FACES; ++i){
207 0 : if(!FACE_CONTAINS_EDGE[i][refEdge]){
208 0 : const int* f = FACE_VRT_INDS[i];
209 :
210 : // if f is a tri, then we'll create a tetrahedron, if it is a
211 : // quad, then we'll create a pyramid.
212 0 : if(f[3] == -1){
213 : // create a tetrahedron
214 0 : inds[fi++] = GOID_TETRAHEDRON;
215 0 : inds[fi++] = f[0]; inds[fi++] = f[1];
216 0 : inds[fi++] = f[2]; inds[fi++] = nVrt;
217 : }
218 : else{
219 : // create a prism
220 0 : inds[fi++] = GOID_PYRAMID;
221 0 : inds[fi++] = f[0]; inds[fi++] = f[1];
222 0 : inds[fi++] = f[2]; inds[fi++] = f[3];
223 0 : inds[fi++] = nVrt;
224 : }
225 : }
226 : }
227 : }
228 : }break;
229 :
230 0 : case 2:
231 : {
232 : // three important cases exist
233 0 : const int re0 = refEdgeInds[0];
234 0 : const int re1 = refEdgeInds[1];
235 : // index of the face that contains both edges. This may be -1.
236 0 : const int refFaceInd = FACE_FROM_EDGES[re0][re1];
237 :
238 : // we currently don't handle cases where refFaceInd == -1
239 0 : if(refFaceInd == -1)
240 : break;
241 :
242 0 : if(FACE_VRT_INDS[refFaceInd][3] != -1){
243 : // a side-face (quad) is refined
244 : // rotate the prism, so that the refined face lies on face 3.
245 : int p[NUM_VERTICES];
246 0 : RotatePrism(p, 3 - refFaceInd);
247 :
248 : int& fi = fillCount;
249 : int* inds = newIndsOut;
250 :
251 0 : const int e3 = EDGE_FROM_VRTS[p[0]][p[3]] + E;
252 0 : const int e5 = EDGE_FROM_VRTS[p[2]][p[5]] + E;
253 :
254 : // check whether the cut is horizontal or vertical (ignore others)
255 0 : if(IS_SIDE_EDGE[re0] && IS_SIDE_EDGE[re1]){
256 : // horizontal cut
257 0 : if(numSnapPoints == 1){
258 : // build a pyramid and a prism
259 0 : if(isSnapPoint[p[1]]){
260 0 : inds[fi++] = GOID_PYRAMID;
261 0 : inds[fi++] = p[0]; inds[fi++] = p[2]; inds[fi++] = e5;
262 0 : inds[fi++] = e3; inds[fi++] = p[1];
263 :
264 0 : inds[fi++] = GOID_PRISM;
265 0 : inds[fi++] = e3; inds[fi++] = p[1]; inds[fi++] = e5;
266 0 : inds[fi++] = p[3]; inds[fi++] = p[4]; inds[fi++] = p[5];
267 : snapPointsProcessed = true;
268 : }
269 0 : else if(isSnapPoint[p[4]]){
270 0 : inds[fi++] = GOID_PYRAMID;
271 0 : inds[fi++] = e3; inds[fi++] = e5; inds[fi++] = p[5];
272 0 : inds[fi++] = p[3]; inds[fi++] = p[4];
273 :
274 0 : inds[fi++] = GOID_PRISM;
275 0 : inds[fi++] = p[0]; inds[fi++] = p[1]; inds[fi++] = p[2];
276 0 : inds[fi++] = e3; inds[fi++] = p[4]; inds[fi++] = e5;
277 : snapPointsProcessed = true;
278 : }
279 : }
280 :
281 0 : if(numSnapPoints == 0 || !snapPointsProcessed){
282 : // build two pyramids and a tetrahedron.
283 0 : inds[fi++] = GOID_PYRAMID;
284 0 : inds[fi++] = p[0]; inds[fi++] = p[2]; inds[fi++] = e5;
285 0 : inds[fi++] = e3; inds[fi++] = p[1];
286 :
287 0 : inds[fi++] = GOID_PYRAMID;
288 0 : inds[fi++] = e3; inds[fi++] = e5; inds[fi++] = p[5];
289 0 : inds[fi++] = p[3]; inds[fi++] = p[4];
290 :
291 0 : inds[fi++] = GOID_TETRAHEDRON;
292 0 : inds[fi++] = p[1]; inds[fi++] = e3;
293 0 : inds[fi++] = p[4]; inds[fi++] = e5;
294 : }
295 : }
296 0 : else if(!IS_SIDE_EDGE[re0] && !IS_SIDE_EDGE[re1]){
297 : // vertical cut
298 : // create two prisms
299 0 : const int e2 = EDGE_FROM_VRTS[p[0]][p[2]] + E;
300 0 : const int e8 = EDGE_FROM_VRTS[p[3]][p[5]] + E;
301 :
302 0 : inds[fi++] = GOID_PRISM;
303 0 : inds[fi++] = p[0]; inds[fi++] = p[1]; inds[fi++] = e2;
304 0 : inds[fi++] = p[3]; inds[fi++] = p[4]; inds[fi++] = e8;
305 0 : inds[fi++] = GOID_PRISM;
306 0 : inds[fi++] = p[1]; inds[fi++] = p[2]; inds[fi++] = e2;
307 0 : inds[fi++] = p[4]; inds[fi++] = p[5]; inds[fi++] = e8;
308 : }
309 : }
310 : else{
311 : // refined face is a triangle - either top or bottom.
312 : // some cases still missing...
313 0 : if(numSnapPoints == 2){
314 : // get the non-marked edge of refFace
315 : int nonRefEdge = -1;
316 0 : for(int iedge = 0; iedge < 3; ++iedge){
317 0 : if(!newEdgeVrts[FACE_EDGE_INDS[refFaceInd][iedge]]){
318 : nonRefEdge = FACE_EDGE_INDS[refFaceInd][iedge];
319 : break;
320 : }
321 : }
322 : UG_ASSERT(nonRefEdge == -1, "A non-ref edge has to exist!");
323 :
324 : // get the edge which connects the two snap points
325 0 : int snapEdge = EDGE_FROM_VRTS[snapPointIndex[0]][snapPointIndex[1]];
326 0 : if(snapEdge == -1 || IS_SIDE_EDGE[snapEdge])
327 : break;
328 :
329 0 : int snapQuad = FACE_FROM_EDGES[nonRefEdge][snapEdge];
330 0 : if(snapQuad == -1 || snapQuad == refFaceInd)
331 : break;
332 :
333 :
334 : int p[NUM_VERTICES];
335 0 : if(refFaceInd == BOTTOM_FACE){
336 : int tp[NUM_VERTICES];
337 0 : RotatePrism(tp, 3 - snapQuad);
338 0 : FlipPrism(p, tp);
339 : }
340 : else
341 0 : RotatePrism(p, 3 - snapQuad);
342 :
343 0 : const int e34 = E + EDGE_FROM_VRTS[p[3]][p[4]];
344 0 : const int e45 = E + EDGE_FROM_VRTS[p[4]][p[5]];
345 : // we have to create two new prisms
346 : int& fi = fillCount;
347 : int* inds = newIndsOut;
348 :
349 0 : inds[fi++] = GOID_PRISM;
350 0 : inds[fi++] = p[0]; inds[fi++] = p[3]; inds[fi++] = e34;
351 0 : inds[fi++] = p[2]; inds[fi++] = p[5]; inds[fi++] = e45;
352 :
353 0 : inds[fi++] = GOID_PRISM;
354 0 : inds[fi++] = p[0]; inds[fi++] = p[1]; inds[fi++] = p[2];
355 0 : inds[fi++] = e34; inds[fi++] = p[4]; inds[fi++] = e45;
356 :
357 : snapPointsProcessed = true;
358 : }
359 : }
360 : }break;
361 :
362 0 : case 3:
363 : {
364 : // here we currently only handle a horizontal slice
365 0 : if(IS_SIDE_EDGE[refEdgeInds[0]] && IS_SIDE_EDGE[refEdgeInds[1]]
366 0 : && IS_SIDE_EDGE[refEdgeInds[2]])
367 : {
368 : const int e3 = E + 3;
369 : const int e4 = E + 4;
370 : const int e5 = E + 5;
371 :
372 : // we have to create two new prisms
373 : int& fi = fillCount;
374 : int* inds = newIndsOut;
375 :
376 0 : inds[fi++] = GOID_PRISM;
377 0 : inds[fi++] = 0; inds[fi++] = 1; inds[fi++] = 2;
378 0 : inds[fi++] = e3; inds[fi++] = e4; inds[fi++] = e5;
379 :
380 0 : inds[fi++] = GOID_PRISM;
381 0 : inds[fi++] = e3; inds[fi++] = e4; inds[fi++] = e5;
382 0 : inds[fi++] = 3; inds[fi++] = 4; inds[fi++] = 5;
383 : }
384 : }break;
385 :
386 : case 4:
387 : {
388 : // only a nice, vertical cut is directly supported in the moment
389 : // this means that 2 marked edges are top-edges and two are bottom edges
390 : int numTop = 0;
391 : int numBottom = 0;
392 0 : for(int i = 0; i < 4; ++i){
393 0 : if(IS_TOP_EDGE[refEdgeInds[i]])
394 0 : ++numTop;
395 0 : else if(IS_BOTTOM_EDGE[refEdgeInds[i]])
396 0 : ++numBottom;
397 : }
398 :
399 0 : if((numTop == 2) && (numBottom == 2)){
400 : // find a quadrilateral, whose corner states are all 1.
401 : int freeFaceInd = -1;
402 0 : for(int i = 0; i < NUM_QUADS; ++i){
403 0 : const int* f = FACE_VRT_INDS[QUADS[i]];
404 : bool allOne = true;
405 0 : for(int j = 0; j < 4; ++j){
406 0 : if(cornerStatus[f[j]] != 1){
407 : allOne = false;
408 : break;
409 : }
410 : }
411 :
412 0 : if(allOne){
413 : freeFaceInd = QUADS[i];
414 : break;
415 : }
416 : }
417 :
418 0 : if(freeFaceInd != -1){
419 : // we got the free quad. It is thus clear, that all associated
420 : // edges, which do not lie in the quad itself, will be refined.
421 : // create a hexahedron and a prism.
422 : // First however rotate the prism so that the free quad is at
423 : // index 3.
424 : int p[NUM_VERTICES];
425 0 : RotatePrism(p, 3 - freeFaceInd);
426 :
427 : // some important indices.
428 0 : const int e0 = EDGE_FROM_VRTS[p[0]][p[1]] + E;
429 0 : const int e1 = EDGE_FROM_VRTS[p[1]][p[2]] + E;
430 0 : const int e6 = EDGE_FROM_VRTS[p[3]][p[4]] + E;
431 0 : const int e7 = EDGE_FROM_VRTS[p[4]][p[5]] + E;
432 :
433 : // now create the elements
434 : int& fi = fillCount;
435 : int* inds = newIndsOut;
436 :
437 0 : inds[fi++] = GOID_HEXAHEDRON;
438 0 : inds[fi++] = p[0]; inds[fi++] = e0;
439 0 : inds[fi++] = e1; inds[fi++] = p[2];
440 0 : inds[fi++] = p[3]; inds[fi++] = e6;
441 0 : inds[fi++] = e7; inds[fi++] = p[5];
442 :
443 0 : inds[fi++] = GOID_PRISM;
444 0 : inds[fi++] = e0; inds[fi++] = p[1]; inds[fi++] = e1;
445 0 : inds[fi++] = e6; inds[fi++] = p[4]; inds[fi++] = e7;
446 : }
447 : }
448 : }break;
449 :
450 : case 6:
451 : {
452 : // we currently only support the case, where only top and bottom edges
453 : // are refined, resulting in 4 new prisms.
454 : bool sideRefined = false;
455 0 : for(int i = 0; i < 6; ++i){
456 0 : if(IS_SIDE_EDGE[refEdgeInds[i]]){
457 : sideRefined = true;
458 : break;
459 : }
460 : }
461 :
462 0 : if(!sideRefined){
463 : // create 4 new prisms. First however store the new edge vertices.
464 : const int e0 = E;
465 : const int e1 = E + 1;
466 : const int e2 = E + 2;
467 : const int e6 = E + 6;
468 : const int e7 = E + 7;
469 : const int e8 = E + 8;
470 :
471 : int& fi = fillCount;
472 : int* inds = newIndsOut;
473 :
474 0 : inds[fi++] = GOID_PRISM;
475 0 : inds[fi++] = 0; inds[fi++] = e0; inds[fi++] = e2;
476 0 : inds[fi++] = 3; inds[fi++] = e6; inds[fi++] = e8;
477 :
478 0 : inds[fi++] = GOID_PRISM;
479 0 : inds[fi++] = 1; inds[fi++] = e1; inds[fi++] = e0;
480 0 : inds[fi++] = 4; inds[fi++] = e7; inds[fi++] = e6;
481 :
482 0 : inds[fi++] = GOID_PRISM;
483 0 : inds[fi++] = 2; inds[fi++] = e2; inds[fi++] = e1;
484 0 : inds[fi++] = 5; inds[fi++] = e8; inds[fi++] = e7;
485 :
486 0 : inds[fi++] = GOID_PRISM;
487 0 : inds[fi++] = e0; inds[fi++] = e1; inds[fi++] = e2;
488 0 : inds[fi++] = e6; inds[fi++] = e7; inds[fi++] = e8;
489 : }
490 : }break;
491 :
492 : case 7:
493 : {
494 : // only the case in which two quads are completly refined is regarded.
495 : // In this case corresponding edges of the bottom and top triangle are
496 : // marked/unmarked.
497 : // order bv so, that the one connected to two marked base-edges
498 : // is bv[0].
499 : int bv[3], tv[3];
500 : bool twoQuadsMarked = false;
501 0 : for(int iedge = 0; iedge < 3; ++iedge){
502 0 : if(!newEdgeVrts[iedge]){
503 0 : if(!newEdgeVrts[iedge + 6]){
504 : twoQuadsMarked = true;
505 :
506 0 : bv[0] = (iedge + 2) % 3;
507 0 : bv[1] = (iedge + 3) % 3;
508 0 : bv[2] = (iedge + 4) % 3;
509 :
510 0 : tv[0] = bv[0] + 3;
511 0 : tv[1] = bv[1] + 3;
512 0 : tv[2] = bv[2] + 3;
513 : }
514 : break;
515 : }
516 : }
517 :
518 :
519 : if(twoQuadsMarked){
520 : int& fi = fillCount;
521 : int* inds = newIndsOut;
522 :
523 0 : inds[fi++] = GOID_PRISM;
524 0 : inds[fi++] = bv[0]; inds[fi++] = E + bv[0];
525 0 : inds[fi++] = E + bv[2];
526 :
527 0 : inds[fi++] = E + bv[0] + 3; inds[fi++] = F + QUADS[bv[0]];
528 0 : inds[fi++] = F + QUADS[bv[2]];
529 :
530 :
531 0 : inds[fi++] = GOID_HEXAHEDRON;
532 0 : inds[fi++] = E + bv[0]; inds[fi++] = bv[1];
533 0 : inds[fi++] = bv[2]; inds[fi++] = E + bv[2];
534 :
535 0 : inds[fi++] = F + QUADS[bv[0]]; inds[fi++] = E + bv[1] + 3;
536 0 : inds[fi++] = E + bv[2] + 3; inds[fi++] = F + QUADS[bv[2]];
537 :
538 0 : inds[fi++] = GOID_PRISM;
539 0 : inds[fi++] = E + bv[0] + 3; inds[fi++] = F + QUADS[bv[0]];
540 0 : inds[fi++] = F + QUADS[bv[2]];
541 :
542 0 : inds[fi++] = tv[0]; inds[fi++] = E + 3 + tv[0];
543 0 : inds[fi++] = E + 3 + tv[2];
544 :
545 :
546 0 : inds[fi++] = GOID_HEXAHEDRON;
547 0 : inds[fi++] = F + QUADS[bv[0]]; inds[fi++] = E + bv[1] + 3;
548 0 : inds[fi++] = E + bv[2] + 3; inds[fi++] = F + QUADS[bv[2]];
549 :
550 0 : inds[fi++] = E + 3 + tv[0]; inds[fi++] = tv[1];
551 0 : inds[fi++] = tv[2]; inds[fi++] = E + 3 + tv[2];
552 : }
553 0 : }break;
554 :
555 0 : case 9: // regular refine
556 : {
557 : // we have to create 8 new prisms
558 : int& fi = fillCount;
559 : int* inds = newIndsOut;
560 0 : inds[fi++] = GOID_PRISM;
561 0 : inds[fi++] = 0; inds[fi++] = E; inds[fi++] = E + 2;
562 0 : inds[fi++] = E + 3; inds[fi++] = F + 1; inds[fi++] = F + 3;
563 :
564 0 : inds[fi++] = GOID_PRISM;
565 0 : inds[fi++] = E; inds[fi++] = 1; inds[fi++] = E + 1;
566 0 : inds[fi++] = F + 1; inds[fi++] = E + 4; inds[fi++] = F + 2;
567 :
568 0 : inds[fi++] = GOID_PRISM;
569 0 : inds[fi++] = E + 2; inds[fi++] = E + 1; inds[fi++] = 2;
570 0 : inds[fi++] = F + 3; inds[fi++] = F + 2; inds[fi++] = E + 5;
571 :
572 0 : inds[fi++] = GOID_PRISM;
573 0 : inds[fi++] = E + 1; inds[fi++] = E + 2; inds[fi++] = E;
574 0 : inds[fi++] = F + 2; inds[fi++] = F + 3; inds[fi++] = F + 1;
575 :
576 0 : inds[fi++] = GOID_PRISM;
577 0 : inds[fi++] = E + 3; inds[fi++] = F + 1; inds[fi++] = F + 3;
578 0 : inds[fi++] = 3; inds[fi++] = E + 6; inds[fi++] = E + 8;
579 :
580 0 : inds[fi++] = GOID_PRISM;
581 0 : inds[fi++] = F + 1; inds[fi++] = E + 4; inds[fi++] = F + 2;
582 0 : inds[fi++] = E + 6; inds[fi++] = 4; inds[fi++] = E + 7;
583 :
584 0 : inds[fi++] = GOID_PRISM;
585 0 : inds[fi++] = F + 3; inds[fi++] = F + 2; inds[fi++] = E + 5;
586 0 : inds[fi++] = E + 8; inds[fi++] = E + 7; inds[fi++] = 5;
587 :
588 0 : inds[fi++] = GOID_PRISM;
589 0 : inds[fi++] = F + 2; inds[fi++] = F + 3; inds[fi++] = F + 1;
590 0 : inds[fi++] = E + 7; inds[fi++] = E + 8; inds[fi++] = E + 6;
591 0 : }break;
592 : }
593 :
594 0 : if(!snapPointsProcessed){
595 : UG_LOG("WARNING: Invalid or unsupported snap-point distribution detected. Ignoring snap-points for this element.\n");
596 : }
597 :
598 0 : if(fillCount == 0){
599 : // call recursive refine
600 0 : fillCount = shared_rules::RecursiveRefine(newIndsOut, newEdgeVrts,
601 : FACE_VRT_INDS, FACE_EDGE_INDS, NUM_VERTICES,
602 : NUM_EDGES, NUM_FACES);
603 :
604 : // the rule requires a new center vertex
605 0 : newCenterOut = true;
606 : }
607 :
608 0 : return fillCount;
609 : }
610 :
611 :
612 0 : int CollapseEdge (int* newIndsOut, int v0, int v1)
613 : {
614 : // get the edge which is connected by v0 and v1
615 0 : const int edgeInd = EDGE_FROM_VRTS[v0][v1];
616 :
617 : // if the edge doesn't exist or if one of the triangle-edges is collapsed,
618 : // no volume can be created
619 0 : if((edgeInd < 3) || (edgeInd > 5))
620 : return 0;
621 :
622 : // the index of the opposing edge of the base quadrilateral
623 : // the two non-collapsed vertical edges
624 0 : const int opEdgeInd0 = 3 + (edgeInd + 1) % 3;
625 0 : const int opEdgeInd1 = 3 + (edgeInd + 2) % 3;
626 :
627 : // the base of the pyramid
628 0 : const int baseFaceInd = FACE_FROM_EDGES[opEdgeInd0][opEdgeInd1];
629 :
630 : // the resulting volume is a pyramid
631 : int fi = 0;
632 0 : newIndsOut[fi++] = GOID_PYRAMID;
633 0 : newIndsOut[fi++] = FACE_VRT_INDS[baseFaceInd][0];
634 0 : newIndsOut[fi++] = FACE_VRT_INDS[baseFaceInd][1];
635 0 : newIndsOut[fi++] = FACE_VRT_INDS[baseFaceInd][2];
636 0 : newIndsOut[fi++] = FACE_VRT_INDS[baseFaceInd][3];
637 0 : newIndsOut[fi++] = v0;
638 :
639 0 : return fi;
640 : }
641 :
642 0 : bool IsRegularRefRule(const int edgeMarks)
643 : {
644 : static const int allEdges = 511; // 0b111111111
645 : static const int allHEdges = 455; // 0b111000111
646 : static const int allVEdges = 56; // 0b000111000
647 :
648 0 : return (edgeMarks == allEdges)
649 0 : || (edgeMarks == allHEdges)
650 0 : || (edgeMarks == allVEdges);
651 : }
652 :
653 : }// end of namespace
654 : }// end of namespace
|
