Line data Source code
1 : /*
2 : * DiamondsEstablish3D.cpp
3 : *
4 : * Created on: 08.12.2025
5 : * Author: Markus Knodel
6 : */
7 :
8 : #include <lib_grid/algorithms/extrusion/DiamondsEstablish3D.h>
9 : #include <typeinfo>
10 :
11 : namespace ug
12 : {
13 :
14 : namespace arte
15 : {
16 :
17 : namespace diamonds
18 : {
19 :
20 0 : DiamondsEstablish3D::DiamondsEstablish3D( Grid & grid,
21 : SubsetHandler & sh,
22 : DiamondsEstablish3D::VecVolManifVrtxCombi const & vecVolManifVrtxC
23 0 : )
24 : :
25 0 : m_grid(grid),
26 0 : m_sh(sh),
27 : m_aaPos(Grid::VertexAttachmentAccessor<APosition>()),
28 0 : m_vecVolManifVrtxCombiToShrink4Diams(vecVolManifVrtxC),
29 : m_vecVolElmFac5(VecVolumeElementFaceQuintuplet()),
30 : m_vecElems2BQuenched(VecElems2BQuenched()),
31 : m_disappearingVols(std::vector<Volume*>()),
32 : m_disappearingFacs(std::vector<Face*>()),
33 : m_vecElemGroupVrtx2BQuenched(VecElemGroupVrtx2BQnchd4D()),
34 : m_attElmGrpVrtx2BQnchd(AttElemGrpVrtx2BQuenchd()),
35 : m_attAccsElmGrpVrtx2BQnchd(Grid::VertexAttachmentAccessor<AttElemGrpVrtx2BQuenchd>()),
36 : m_attMarkTwinFace(ABool()),
37 : m_attAccsFacIsTwinFac(Grid::FaceAttachmentAccessor<ABool>()),
38 : m_attMarkShiftFace(ABool()),
39 : m_attAccsFacIsShiftFac(Grid::FaceAttachmentAccessor<ABool>()),
40 : m_attMarkShiftTriangleFace(ABool()),
41 : m_attAccsFacIsShiftTriangleFac(Grid::FaceAttachmentAccessor<ABool>()),
42 : m_attMarkShiftQuadriliteralFace(ABool()),
43 : m_attAccsFacIsShiftQuadriliteralFac(Grid::FaceAttachmentAccessor<ABool>()),
44 : m_attMarkVolGetsShrinked(ABool()),
45 : m_attAccsVolGetsShrinked(Grid::VolumeAttachmentAccessor<ABool>()),
46 : m_attVrtVec(AttVrtVec()),
47 : m_attAccsVrtVecVol(Grid::VolumeAttachmentAccessor<AttVrtVec>()),
48 : // m_attVrtVecFace(AttVrtVec()),
49 : // m_attAccsVrtVecFace(Grid::VolumeAttachmentAccessor<AttVrtVec>()),
50 : m_attMarkVrtxIsCenterVrtx(ABool()),
51 : m_attAccsVrtxIsCenterVrtx(Grid::VertexAttachmentAccessor<ABool>()),
52 : m_attMarkVrtxIsShiftVrtx(ABool()),
53 : m_attAccsVrtxIsShiftVrtx(Grid::VertexAttachmentAccessor<ABool>()),
54 : m_attMarkEdgeIsShiftEdge(ABool()),
55 : m_attAccsEdgeIsShiftEdge(Grid::EdgeAttachmentAccessor<ABool>()),
56 : m_attInfoVecSudosTouchingVrtx(AttVecInt()),
57 : m_attAccsInfoVecSudosTouchingVrtx(Grid::VertexAttachmentAccessor<AttVecInt>()),
58 : m_attVrtxIsMidPtOfShiftVrtx(ABool()),
59 : m_attAccsVrtxIsMidPtOfShiftVrtcs(Grid::VertexAttachmentAccessor<ABool>()),
60 : m_vecCombiShiftVrtxMidVrtx(VecCombiShiftVrtxMidVrtx()),
61 : m_attMidPtVrtxOfShiftVrtx(AVertex()),
62 : m_attAccsMidPtVrtxOfShiftVrtx(Grid::VertexAttachmentAccessor<AVertex>()),
63 : m_attCenterVrtxOfShiftVrtx(AVertex()),
64 : m_attAccsCenterVrtxOfShiftVrtx(Grid::VertexAttachmentAccessor<AVertex>()),
65 : m_vecCombiNewVolsTwoCross(VecCombiNewVolsProps()),
66 : m_vecCombiNewVolsThreeCross(VecCombiNewVolsProps()),
67 : m_sudosTable(VecIndxVec()),
68 : m_vecCombiCntrVrtxSudo(VecCombiCntrVrtxSudo()),
69 : m_attEdgeCanBeRemoved(ABool()),
70 : m_attAccsEdgeCanBeRemoved(Grid::EdgeAttachmentAccessor<ABool>()),
71 : m_attNewSudoOfVrtx(AInt()),
72 : m_attAccsNewSudoOfVrtx(Grid::VertexAttachmentAccessor<AInt>()),
73 : m_faces2BDeletedAtLastStep(std::vector<Face*>()),
74 0 : m_edges2BDeletedAtLastStep(std::vector<Edge*>())
75 : {
76 : // // Notloesung, nicht in die erste Initialisierung vor geschweifter Klammer, da copy constructor privat
77 0 : m_sel = Selector();
78 0 : }
79 :
80 : //bool DiamondsEstablish3D::createTheDiamonds()
81 : //{
82 : // IndexType sudosVols = m_sh.num_subsets();
83 : //
84 : // IndexType sudosEdges = sudosVols + 1;
85 : //
86 : // IndexType sudosFaces = sudosVols + 2;
87 : //
88 : // for( auto & vmvcd : m_vecVolManifVrtxCombiToShrink4Diams )
89 : // {
90 : // Volume* vol;
91 : // vmvcd.spuckVol(vol);
92 : //
93 : // m_sh.assign_subset(vol, sudosVols);
94 : //
95 : // IndexType numLowdimElmsFnd = vmvcd.computeTheLowdimElm(m_grid);
96 : //
97 : // if( numLowdimElmsFnd != 1 )
98 : // {
99 : // UG_LOG("number of lowdim elems found strange " << numLowdimElmsFnd << std::endl);
100 : // UG_THROW("number of lowdim elems found strange " << numLowdimElmsFnd << std::endl);
101 : // }
102 : //
103 : // Edge* edge;
104 : // vmvcd.spuckLowdimElem( edge );
105 : //
106 : // if( edge == nullptr )
107 : // {
108 : // UG_LOG("Edge nicht gefunden " << std::endl);
109 : // UG_THROW("Edge nicht gefunden " << std::endl);
110 : // }
111 : // m_sh.assign_subset( edge, sudosEdges);
112 : //
113 : // Face * fac;
114 : // vmvcd.spuckManif(fac);
115 : // m_sh.assign_subset( fac, sudosFaces );
116 : // }
117 : //
118 : // UG_LOG("Established diamonds" << std::endl);
119 : //
120 : //
121 : // return true;
122 : //}
123 :
124 0 : DiamondsEstablish3D::~DiamondsEstablish3D()
125 : {
126 : // Auto-generated destructor stub
127 0 : }
128 :
129 0 : bool DiamondsEstablish3D::initialize()
130 : {
131 0 : if(!m_grid.has_vertex_attachment(aPosition) )
132 : {
133 : UG_LOG("Error in ExpandFractures Arte 3D: Missing position attachment");
134 0 : return false;
135 : }
136 :
137 0 : m_aaPos = Grid::VertexAttachmentAccessor<APosition>(m_grid, aPosition);
138 :
139 0 : return true;
140 : }
141 :
142 : /////////////////////////////////////////////////////////////////
143 :
144 :
145 0 : bool DiamondsEstablish3D::createTheDiamonds()
146 : {
147 0 : if( ! initialize())
148 : {
149 : UG_LOG("initialization diamonds did not work " << std::endl);
150 0 : return false;
151 : }
152 :
153 0 : if( ! setSelector() )
154 : return false;
155 :
156 : UG_LOG("selektiert" << std::endl);
157 :
158 0 : if( ! figureOutTheEdges() )
159 : {
160 : UG_LOG("Edges not found " << std::endl);
161 0 : return false;
162 : }
163 :
164 : UG_LOG("edges figured out " << std::endl);
165 :
166 0 : if( ! findRegions2BShrinked())
167 : {
168 : UG_LOG("Regions to be shrinked not found " << std::endl);
169 0 : return false;
170 : }
171 :
172 : UG_LOG("regions to be shrinked found " << std::endl);
173 :
174 0 : if( ! establishElems2BeQuenched() )
175 : {
176 : UG_LOG("quenching elements not establishable " << std::endl);
177 0 : return false;
178 : }
179 :
180 : UG_LOG("quenching elements established " << std::endl);
181 :
182 0 : if( ! sortElems2BQuenched())
183 : {
184 : UG_LOG("sorting quench impossible" << std::endl);
185 0 : return false;
186 : }
187 :
188 : UG_LOG("sorted elements finished " << std::endl);
189 :
190 0 : if( ! attachMarkers())
191 : {
192 : UG_LOG("Markers do not want D" << std::endl);
193 0 : return false;
194 : }
195 :
196 : UG_LOG("markers attached" << std::endl);
197 :
198 0 : if( ! assignBasicAtts())
199 : {
200 : UG_LOG("unassignale basic att" << std::endl);
201 0 : return false;
202 : }
203 :
204 : UG_LOG("basic attcs assigned" << std::endl);
205 :
206 0 : if( ! trafoCollectedInfo2Attachments())
207 : {
208 : UG_LOG("untrafoable infos" << std::endl);
209 0 : return false;
210 : }
211 :
212 : UG_LOG("trafo collected infos " << std::endl);
213 :
214 : // disable selection inheritance to avoid infinite recursion.
215 0 : m_sel.enable_selection_inheritance(false);
216 :
217 0 : if( ! createConditionForNewVrtcs())
218 : {
219 : UG_LOG("conditions for new vertices do not work " << std::endl);
220 0 : return false;
221 : }
222 :
223 : UG_LOG("conditions for new vertices created" << std::endl);
224 :
225 0 : if( ! distributeInfosForShrinkingVols())
226 : {
227 : UG_LOG("info distribution did not work " << std::endl);
228 0 : return false;
229 : }
230 :
231 : UG_LOG("info distributed" << std::endl);
232 :
233 0 : if( ! determineShiftFaces())
234 : {
235 : UG_LOG("shift faces not determinable " << std::endl);
236 0 : return false;
237 : }
238 :
239 : UG_LOG("detect removable edges " << std::endl);
240 :
241 0 : if( ! detectRemovableEdges())
242 : {
243 : UG_LOG("removable edges not detected " << std::endl);
244 0 : return false;
245 : }
246 :
247 : UG_LOG("shrink volumes " << std::endl);
248 :
249 0 : if( ! shrinkVolumes())
250 : {
251 : UG_LOG("shrinken schief gegangen " << std::endl);
252 0 : return false;
253 : }
254 :
255 : UG_LOG("volumes shrinked " << std::endl);
256 :
257 0 : if( ! postprocessNewDiamVols())
258 : {
259 : UG_LOG("postprocessing diam vols not working" << std::endl);
260 0 : return false;
261 : }
262 :
263 : UG_LOG("diam vols postprocessed" << std::endl);
264 :
265 0 : if( ! detachMarkers())
266 : {
267 : UG_LOG("Markers do not get detouched D" << std::endl);
268 0 : return false;
269 : }
270 :
271 : UG_LOG("markers detached" << std::endl);
272 :
273 : UG_LOG("noch nicht soweit: Established diamonds" << std::endl);
274 :
275 :
276 :
277 0 : return true;
278 : }
279 :
280 : /////////////////////////////////////////////////////////////////
281 :
282 :
283 0 : bool DiamondsEstablish3D::setSelector()
284 : {
285 :
286 0 : m_sel.assign_grid(m_grid);
287 :
288 0 : m_sel.enable_autoselection(false);
289 0 : m_sel.enable_selection_inheritance(true); //required for select and mark, disabled later
290 0 : m_sel.enable_strict_inheritance(false);
291 :
292 :
293 0 : return true;
294 : }
295 :
296 : /////////////////////////////////////////////////////////////////
297 :
298 0 : bool DiamondsEstablish3D::figureOutTheEdges()
299 : {
300 :
301 : // compute the edges connecting the old and shift vertex
302 :
303 0 : for( auto & vmvcd : m_vecVolManifVrtxCombiToShrink4Diams )
304 : {
305 : // IndexType numLowdimElmsFnd = vmvcd.checkIntegrity(m_grid);
306 : // if( numLowdimElmsFnd != 1 )
307 0 : if( ! vmvcd.checkIntegrity(m_grid) )
308 : {
309 : UG_LOG("number of lowdim elems found strange " << std::endl);
310 : return false;
311 : }
312 : }
313 :
314 : UG_LOG("Edges found" << std::endl);
315 :
316 0 : return true;
317 : }
318 :
319 : ///////////////////////////////////////////////////////////////////////////////////////
320 :
321 :
322 : ///////////////////////////////////////////////////////////////////////////////
323 :
324 : ///////////////////////////////////////////////////////////////////////////////////////
325 :
326 0 : bool DiamondsEstablish3D::findRegions2BShrinked()
327 : {
328 : UG_LOG("want to find regions to be shrinked" << std::endl);
329 :
330 :
331 0 : VecVolManifVrtxCombi vecVolManifVrtxCopy = m_vecVolManifVrtxCombiToShrink4Diams;
332 :
333 : int d_out = 0;
334 :
335 0 : for( typename VecVolManifVrtxCombi::iterator itVMVOuter = vecVolManifVrtxCopy.begin();
336 0 : itVMVOuter != vecVolManifVrtxCopy.end(); )
337 : {
338 :
339 : // UG_LOG("out round " << d_out << std::endl);
340 :
341 : bool partnerFound = false;
342 :
343 0 : VolManifVrtxCombi outer = *itVMVOuter;
344 :
345 : VrtxPair oldAndShiftVrtxOuter;
346 : outer.spuckOldAndShiftVrtx( oldAndShiftVrtxOuter );
347 :
348 : Vertex * oldVrtxOuter = oldAndShiftVrtxOuter.first;
349 : Vertex * shiftVrtxOuter = oldAndShiftVrtxOuter.second;
350 :
351 : Face * faceOuter;
352 : outer.spuckManif(faceOuter);
353 :
354 : IndexType sudoOuter = outer.spuckSudo();
355 :
356 : int d_in = 0;
357 :
358 0 : for( typename VecVolManifVrtxCombi::iterator itVMVInner = itVMVOuter + 1;
359 0 : itVMVInner != vecVolManifVrtxCopy.end(); )
360 : {
361 : // UG_LOG("start in " << d_in << std::endl);
362 :
363 0 : VolManifVrtxCombi inner = *itVMVInner;
364 :
365 : VrtxPair oldAndShiftVrtxInner;
366 : inner.spuckOldAndShiftVrtx( oldAndShiftVrtxInner );
367 :
368 : Vertex * oldVrtxInner = oldAndShiftVrtxInner.first;
369 :
370 0 : if( oldVrtxInner == oldVrtxOuter )
371 : {
372 : Vertex * shiftVrtxInner = oldAndShiftVrtxInner.second;
373 :
374 : // if connected by a face, twin to outer found, partner may exist or not
375 : // check if the face is identical, else is from the partner twin
376 :
377 : Face * faceInner;
378 : inner.spuckManif(faceInner);
379 :
380 0 : if( faceInner == faceOuter )
381 : {
382 : // twin found
383 :
384 : PairVolFacVrtxCmb prVolFacVrtxC( outer, inner );
385 : VolumeElementFaceQuintuplet vef5;
386 :
387 0 : trafoVolFacVrtxCombiPair2FullLowDimManifQuintuplet( prVolFacVrtxC, vef5 );
388 :
389 0 : if( ! vef5.checkIntegrity() )
390 : {
391 : UG_LOG("strange twin produced" << std::endl);
392 0 : return false;
393 : }
394 :
395 0 : m_vecVolElmFac5.push_back(vef5);
396 :
397 : itVMVInner = vecVolManifVrtxCopy.erase(itVMVInner);
398 :
399 : partnerFound = true;
400 :
401 : }
402 : else
403 : {
404 : itVMVInner++;
405 : }
406 :
407 : if( partnerFound )
408 : {
409 : break;
410 : }
411 : }
412 : else
413 : {
414 : itVMVInner++;
415 : }
416 : if( partnerFound )
417 : {
418 : break;
419 : }
420 :
421 :
422 : // UG_LOG("end in " << d_in << std::endl);
423 : d_in++;
424 : }
425 :
426 : // UG_LOG("first round finished " << std::endl);
427 0 : if( ! partnerFound )
428 : {
429 : UG_LOG("no partner found " << std::endl);
430 : return false;
431 : }
432 : else
433 : {
434 : itVMVOuter = vecVolManifVrtxCopy.erase(itVMVOuter);
435 : }
436 :
437 : // UG_LOG("end out " << d_out << std::endl);
438 : d_out++;
439 : }
440 :
441 : UG_LOG("end of search reached diams " << std::endl);
442 :
443 : // for( auto & v: m_vecVolElmFac5 )
444 : // {
445 : // IndexType sudoNum = m_sh.num_subsets();
446 : //
447 : // IndexType sudoVols = sudoNum;
448 : // IndexType sudoFacs = sudoNum+1;
449 : // IndexType sudoEdgs = sudoNum+2;
450 : // IndexType sudoVrtx = sudoNum+3;
451 : //
452 : // Volume * vol1;
453 : // Volume * vol2;
454 : // Face * fac;
455 : // Edge * edg1;
456 : // Edge * edg2;
457 : // Vertex * vrtC;
458 : // Vertex * vrtE1;
459 : // Vertex * vrtE2;
460 : //
461 : // std::pair<VolumeElementTwin,VolumeElementTwin> pvv;
462 : //
463 : // v.spuckPairFullLowDimTwin(pvv);
464 : //
465 : // pvv.first.spuckFullDimElem(vol1);
466 : // pvv.second.spuckFullDimElem(vol2);
467 : //
468 : // v.spuckManifElem( fac );
469 : //
470 : // pvv.first.spuckLowDimElem(edg1);
471 : // pvv.second.spuckLowDimElem(edg2);
472 : //
473 : // v.spuckCenterVertex(vrtC);
474 : //
475 : // VrtxPair vp;
476 : // v.spuckShiftVrtcs(vp);
477 : //
478 : // vrtE1 = vp.first;
479 : // vrtE2 = vp.second;
480 : //
481 : // m_sh.assign_subset(vol1, sudoVols);
482 : // m_sh.assign_subset(vol2, sudoVols);
483 : // m_sh.assign_subset(fac, sudoFacs);
484 : // m_sh.assign_subset(edg1, sudoEdgs);
485 : // m_sh.assign_subset(edg2, sudoEdgs);
486 : // m_sh.assign_subset(vrtC, sudoVrtx);
487 : // m_sh.assign_subset(vrtE1, sudoVrtx);
488 : // m_sh.assign_subset(vrtE2, sudoVrtx);
489 : //
490 : // }
491 :
492 :
493 : UG_LOG("found regions to be shrinked" << std::endl);
494 :
495 :
496 : return true;
497 0 : }
498 :
499 :
500 : ////////////////////////////////////////////////////////////////////////////
501 :
502 0 : bool DiamondsEstablish3D::trafoVolFacVrtxCombiPair2FullLowDimManifQuintuplet(
503 : PairVolFacVrtxCmb & prVolFacVrtxC, VolumeElementFaceQuintuplet & vef5 )
504 : {
505 :
506 : VolManifVrtxCombi & mvcOne = prVolFacVrtxC.first;
507 : VolManifVrtxCombi & mvcTwo = prVolFacVrtxC.second;
508 :
509 : IndexType sudoOne = mvcOne.spuckSudo();
510 : IndexType sudoTwo = mvcTwo.spuckSudo();
511 :
512 0 : if( sudoOne != sudoTwo )
513 : {
514 : UG_LOG("sudos differ " << std::endl);
515 0 : return false;
516 : }
517 :
518 : IndexType sudo = sudoOne;
519 :
520 : Face * faceOne;
521 : Face * faceTwo;
522 :
523 : mvcOne.spuckManif(faceOne);
524 : mvcTwo.spuckManif(faceTwo);
525 :
526 0 : if( faceOne != faceTwo )
527 : {
528 : UG_LOG("faces differ " << std::endl);
529 0 : return false;
530 : }
531 :
532 : Face * connectingFace = faceOne;
533 :
534 : VrtxPair oldAndShiftVrtxOne;
535 : mvcOne.spuckOldAndShiftVrtx( oldAndShiftVrtxOne );
536 :
537 : VrtxPair oldAndShiftVrtxTwo;
538 : mvcTwo.spuckOldAndShiftVrtx( oldAndShiftVrtxTwo );
539 :
540 : Vertex * oldVrtxOne = oldAndShiftVrtxOne.first;
541 : Vertex * oldVrtxTwo = oldAndShiftVrtxTwo.first;
542 :
543 0 : if( oldVrtxOne != oldVrtxTwo )
544 : {
545 : UG_LOG("center vertices not identical " << std::endl);
546 0 : return false;
547 : }
548 :
549 : Vertex * oldVrtx = oldVrtxOne;
550 :
551 : Vertex * shiftVrtxOne = oldAndShiftVrtxOne.second;
552 : Vertex * shiftVrtxTwo = oldAndShiftVrtxTwo.second;
553 :
554 0 : if( shiftVrtxOne == shiftVrtxTwo )
555 : {
556 : UG_LOG("shift vertices coincide but should not " << std::endl);
557 0 : return false;
558 : }
559 :
560 : Volume * volOne;
561 : Volume * volTwo;
562 :
563 : mvcOne.spuckFulldimElem( volOne );
564 : mvcTwo.spuckFulldimElem( volTwo );
565 :
566 0 : if( volOne == volTwo )
567 : {
568 : UG_LOG("volumes coincide but should not " << std::endl);
569 0 : return false;
570 : }
571 :
572 : Edge * edgeOne;
573 : Edge * edgeTwo;
574 :
575 : mvcOne.spuckLowDimElem( edgeOne );
576 : mvcTwo.spuckLowDimElem( edgeTwo );
577 0 : if( edgeOne == edgeTwo )
578 : {
579 : UG_LOG("edges coincide but should not " << std::endl);
580 0 : return false;
581 : }
582 :
583 : VolumeElementTwin volElTwinOne( volOne, edgeOne, sudo );
584 : VolumeElementTwin volElTwinTwo( volTwo, edgeTwo, sudo );
585 :
586 0 : if( ! volElTwinOne.checkIntegrity() || ! volElTwinTwo.checkIntegrity() )
587 : {
588 : UG_LOG("twins of vol edge not integer " << std::endl);
589 0 : return false;
590 : }
591 :
592 : std::pair<VolumeElementTwin,VolumeElementTwin> volElTwinPair( volElTwinOne, volElTwinTwo );
593 :
594 0 : vef5 = VolumeElementFaceQuintuplet( volElTwinPair, connectingFace );
595 :
596 0 : if( ! vef5.checkIntegrity() )
597 : {
598 : UG_LOG( "quitent not integer " << std::endl );
599 0 : return false;
600 : }
601 :
602 : return true;
603 : }
604 :
605 :
606 : ////////////////////////////////////////////////////////////////////////////
607 :
608 0 : bool DiamondsEstablish3D::establishElems2BeQuenched()
609 : {
610 0 : VecVolumeElementFaceQuintuplet vvef5 = m_vecVolElmFac5;
611 :
612 0 : for( typename VecVolumeElementFaceQuintuplet::iterator itOuter = vvef5.begin();
613 0 : itOuter != vvef5.end();
614 : )
615 : {
616 0 : VolumeElementFaceQuintuplet vef5Outer = *itOuter;
617 :
618 : VecVolumeElementFaceQuintuplet vfld5ThisEdgePr;
619 :
620 0 : vfld5ThisEdgePr.push_back(vef5Outer);
621 :
622 : // now collect all quintuplets belonging to the same edge pair
623 :
624 : EdgePair edgPrOuter;
625 :
626 : vef5Outer.spuckPairLowDimElem(edgPrOuter);
627 :
628 0 : for( typename VecVolumeElementFaceQuintuplet::iterator itInner = vvef5.begin() + 1;
629 0 : itInner != vvef5.end();
630 : )
631 : {
632 0 : VolumeElementFaceQuintuplet vef5Inner = *itInner;
633 :
634 : EdgePair edgPrInner;
635 : vef5Inner.spuckPairLowDimElem(edgPrInner);
636 :
637 : if( edgPrInner == edgPrOuter )
638 : {
639 0 : vfld5ThisEdgePr.push_back(vef5Inner);
640 :
641 : itInner = vvef5.erase(itInner);
642 : }
643 : else // if( edgPrInner != edgPrOuter )
644 : {
645 : EdgePair testEdges = edgPrInner;
646 : std::swap( testEdges.first, testEdges.second );
647 :
648 : if( testEdges == edgPrOuter )
649 : {
650 0 : if( ! vef5Inner.swapEntries() )
651 : {
652 : UG_LOG("swapping not worked " << std::endl);
653 : return false;
654 : }
655 :
656 0 : if( ! vef5Inner.checkIntegrity())
657 : {
658 : UG_LOG("not integer any more after sapping test" << std::endl);
659 : return false;
660 : }
661 :
662 : EdgePair testAgainEdgPr;
663 : vef5Inner.spuckPairLowDimElem(testAgainEdgPr);
664 :
665 : if( testAgainEdgPr == edgPrOuter )
666 : {
667 0 : vfld5ThisEdgePr.push_back(vef5Inner);
668 :
669 : itInner = vvef5.erase(itInner);
670 : }
671 : else
672 : {
673 : UG_LOG("should be correct but is not " << std::endl);
674 : return false;
675 : }
676 : }
677 : else
678 : {
679 : itInner++;
680 : }
681 : }
682 : }
683 :
684 : Elems2BQuenched elem2BQuenched( vfld5ThisEdgePr );
685 :
686 0 : if( ! elem2BQuenched.checkIntegrity())
687 : {
688 : UG_LOG("an elem to be quenched not integer" << std::endl);
689 : return false;
690 : }
691 :
692 0 : if( ! assignMidPointOfShiftVrtcs(elem2BQuenched))
693 : {
694 : UG_LOG("mid point not assignable" << std::endl);
695 : return false;
696 : }
697 :
698 0 : m_vecElems2BQuenched.push_back(elem2BQuenched);
699 :
700 : itOuter = vvef5.erase(itOuter);
701 0 : }
702 :
703 : // int d_q = 10;
704 :
705 : return true;
706 0 : }
707 :
708 : //////////////////////////////////////////////////////////////////////////////7
709 :
710 0 : void DiamondsEstablish3D::debugE2bQ(Elems2BQuenched & e2bq)
711 : {
712 : // for( Elems2BQuenched & e2bq : m_vecElems2BQuenched )
713 : {
714 0 : IndexType sudoNum = m_sh.num_subsets();
715 :
716 : IndexType sudoVols = sudoNum;
717 0 : IndexType sudoFacs = sudoNum+1;
718 0 : IndexType sudoEdgs = sudoNum+2;
719 0 : IndexType sudoVrtx = sudoNum+3;
720 :
721 : VecVolumeElementFaceQuintuplet vve5;
722 :
723 : EdgePair edgP;
724 : e2bq.spuckPairLowDimElem(edgP);
725 :
726 : e2bq.spuckVecFullLowDimManifQuintuplet(vve5);
727 :
728 : Vertex * centerVrtx;
729 : e2bq.spuckCenterVertex(centerVrtx);
730 :
731 : vector3 vertexLocation;
732 :
733 0 : if( centerVrtx != nullptr )
734 : vertexLocation = m_aaPos[centerVrtx];
735 : else
736 : {
737 : UG_LOG("center null " << std::endl);
738 : }
739 :
740 : number x = vertexLocation[0];
741 : number y = vertexLocation[1];
742 : number z = vertexLocation[2];
743 :
744 0 : number distSq = ( x - 0.5 )*( x - 0.5 ) + ( y - 0.5 )*( y - 0.5 ) + ( z - 0.5 )*( z - 0.5 );
745 :
746 0 : if( distSq < 0.02 )
747 : {
748 :
749 0 : for( VolumeElementFaceQuintuplet & v: vve5 )
750 : {
751 : Volume * vol1;
752 : Volume * vol2;
753 : Face * fac;
754 : Edge * edg1;
755 : Edge * edg2;
756 : Vertex * vrtC;
757 : Vertex * vrtE1;
758 : Vertex * vrtE2;
759 :
760 : std::pair<VolumeElementTwin,VolumeElementTwin> pvv;
761 :
762 : v.spuckPairFullLowDimTwin(pvv);
763 :
764 : pvv.first.spuckFullDimElem(vol1);
765 : pvv.second.spuckFullDimElem(vol2);
766 :
767 : v.spuckManifElem( fac );
768 :
769 : pvv.first.spuckLowDimElem(edg1);
770 : pvv.second.spuckLowDimElem(edg2);
771 :
772 : v.spuckCenterVertex(vrtC);
773 :
774 : VrtxPair vp;
775 : v.spuckShiftVrtcs(vp);
776 :
777 : vrtE1 = vp.first;
778 : vrtE2 = vp.second;
779 :
780 :
781 : Volume * newVol1;
782 : Volume * newVol2;
783 :
784 0 : newVol1 = *m_grid.create<Prism>(
785 0 : PrismDescriptor( vol1->vertex(0),
786 0 : vol1->vertex(1),
787 0 : vol1->vertex(2),
788 0 : vol1->vertex(3),
789 0 : vol1->vertex(4),
790 0 : vol1->vertex(5)
791 : )
792 : );
793 :
794 0 : newVol2 = *m_grid.create<Prism>(
795 0 : PrismDescriptor( vol2->vertex(0),
796 0 : vol2->vertex(1),
797 0 : vol2->vertex(2),
798 0 : vol2->vertex(3),
799 0 : vol2->vertex(4),
800 0 : vol2->vertex(5)
801 : )
802 : );
803 :
804 :
805 0 : m_sh.assign_subset(newVol1, sudoVols);
806 0 : m_sh.assign_subset(newVol2, sudoVols);
807 :
808 : Face * newFac;
809 :
810 0 : newFac = *m_grid.create<Triangle>(TriangleDescriptor( fac->vertex(0), fac->vertex(1), fac->vertex(2) ));
811 :
812 0 : m_sh.assign_subset(newFac, sudoFacs);
813 :
814 :
815 : // m_sh.assign_subset(vol1, sudoVols);
816 : // m_sh.assign_subset(vol2, sudoVols);
817 :
818 : // m_sh.assign_subset(fac, sudoFacs);
819 :
820 :
821 : // m_sh.assign_subset(fac, sudoFacs);
822 : // m_sh.assign_subset(edg1, sudoEdgs);
823 : // m_sh.assign_subset(edg2, sudoEdgs);
824 : // m_sh.assign_subset(vrtC, sudoVrtx);
825 : // m_sh.assign_subset(vrtE1, sudoVrtx);
826 : // m_sh.assign_subset(vrtE2, sudoVrtx);
827 :
828 : Edge * edgOne;
829 : Edge * edgTwo;
830 :
831 0 : edgOne = *m_grid.create<RegularEdge>(EdgeDescriptor( edgP.first->vertex(0), edgP.first->vertex(1) ));
832 0 : edgTwo = *m_grid.create<RegularEdge>(EdgeDescriptor( edgP.second->vertex(0), edgP.second->vertex(1) ));
833 :
834 0 : m_sh.assign_subset(edgOne,sudoEdgs);
835 0 : m_sh.assign_subset(edgTwo,sudoEdgs);
836 :
837 : Vertex * centerV;
838 :
839 : e2bq.spuckCenterVertex(centerV);
840 :
841 0 : Vertex * newCenter = *m_grid.create<RegularVertex>();
842 :
843 : m_aaPos[newCenter] = m_aaPos[centerV];
844 :
845 0 : m_sh.assign_subset(newCenter, sudoVrtx);
846 :
847 : }
848 :
849 :
850 :
851 :
852 : // if( d_q == 20 )
853 : // return true;
854 : //
855 : // d_q++;
856 : //
857 : }
858 0 : }
859 0 : }
860 :
861 : ////////////////////////////////////////////////////////////////////////////
862 :
863 0 : bool DiamondsEstablish3D::sortElems2BQuenched()
864 : {
865 0 : VecElems2BQuenched ve2bq = m_vecElems2BQuenched;
866 :
867 0 : for( typename VecElems2BQuenched::iterator itOuter = ve2bq.begin();
868 0 : itOuter != ve2bq.end();
869 : )
870 : {
871 : Elems2BQuenched e2bqOuter = *itOuter;
872 :
873 : Vertex * centerVrtxOuter;
874 : VecElems2BQuenched vecElems2Q4ThisVrtx;
875 : e2bqOuter.spuckOrigCenterVertex( centerVrtxOuter );
876 :
877 0 : vecElems2Q4ThisVrtx.push_back(e2bqOuter);
878 :
879 0 : for( typename VecElems2BQuenched::iterator itInner = ve2bq.begin() + 1;
880 0 : itInner != ve2bq.end();
881 : )
882 : {
883 : Elems2BQuenched e2bqInner = *itInner;
884 :
885 : Vertex * centerVrtxInner;
886 : e2bqInner.spuckOrigCenterVertex( centerVrtxInner );
887 :
888 0 : if( centerVrtxOuter == centerVrtxInner )
889 : {
890 0 : vecElems2Q4ThisVrtx.push_back(e2bqInner);
891 : itInner = ve2bq.erase(itInner);
892 : }
893 : else
894 : {
895 : itInner++;
896 : }
897 : }
898 :
899 : ElemGroupVrtx2BQuenched4Diams egv2b( vecElems2Q4ThisVrtx );
900 :
901 0 : if( ! egv2b.checkIntegrity() )
902 : {
903 : UG_LOG("elem group not integer for vertex" << std::endl);
904 : return false;
905 : }
906 :
907 0 : m_vecElemGroupVrtx2BQuenched.push_back( egv2b );
908 :
909 : itOuter = ve2bq.erase(itOuter);
910 :
911 0 : }
912 :
913 : return true;
914 0 : }
915 :
916 : ////////////////////////////////////////////////////////////////////////////
917 :
918 0 : bool DiamondsEstablish3D::attachMarkers()
919 : {
920 : // m_aAdjMarkerVFP = AttVertFracProp();
921 : //
922 : //// support::VertexFracturePropertiesVol<IndexType> vfp0; // false, 0 );
923 : // VertxFracPropts vfp0; // false, 0 );
924 : // // default value: no boundary fracture, no fractures crossing
925 : //
926 : // m_grid.attach_to_vertices_dv( m_aAdjMarkerVFP, vfp0 );
927 : // m_aaMarkVrtVFP = Grid::VertexAttachmentAccessor<AttVertFracProp> ( m_grid, m_aAdjMarkerVFP );
928 :
929 : m_attElmGrpVrtx2BQnchd = AttElemGrpVrtx2BQuenchd();
930 :
931 : ElemGroupVrtx2BQuenched4Diams egv2bQEmpty;
932 :
933 0 : m_grid.attach_to_vertices_dv( m_attElmGrpVrtx2BQnchd, egv2bQEmpty );
934 :
935 0 : m_attAccsElmGrpVrtx2BQnchd = Grid::VertexAttachmentAccessor<AttElemGrpVrtx2BQuenchd>( m_grid, m_attElmGrpVrtx2BQnchd);
936 :
937 : m_attMarkTwinFace = ABool(); // used to know if an face is twin face
938 :
939 0 : m_grid.attach_to_faces_dv( m_attMarkTwinFace, false );
940 :
941 0 : m_attAccsFacIsTwinFac = Grid::FaceAttachmentAccessor<ABool>( m_grid, m_attMarkTwinFace );
942 :
943 : m_attMarkShiftFace = ABool(); // used to know if an face is shift face
944 :
945 0 : m_grid.attach_to_faces_dv( m_attMarkShiftFace, false );
946 :
947 0 : m_attAccsFacIsShiftFac = Grid::FaceAttachmentAccessor<ABool>( m_grid, m_attMarkShiftFace );
948 :
949 : m_attMarkShiftTriangleFace = ABool(); // used to know if an face is shift face
950 :
951 0 : m_grid.attach_to_faces_dv(m_attMarkShiftTriangleFace, false);
952 :
953 0 : m_attAccsFacIsShiftTriangleFac = Grid::FaceAttachmentAccessor<ABool>(m_grid, m_attMarkShiftTriangleFace);
954 :
955 : m_attMarkShiftQuadriliteralFace = ABool(); // used to know if an face is shift face
956 :
957 0 : m_grid.attach_to_faces_dv(m_attMarkShiftQuadriliteralFace, false);
958 :
959 0 : m_attAccsFacIsShiftQuadriliteralFac = Grid::FaceAttachmentAccessor<ABool>(m_grid, m_attMarkShiftQuadriliteralFace);
960 :
961 : m_attMarkVolGetsShrinked = ABool();
962 :
963 0 : m_grid.attach_to_volumes_dv( m_attMarkVolGetsShrinked, false );
964 :
965 0 : m_attAccsVolGetsShrinked = Grid::VolumeAttachmentAccessor<ABool>( m_grid, m_attMarkVolGetsShrinked );
966 :
967 : m_attVrtVec = AttVrtVec();
968 :
969 0 : m_grid.attach_to_volumes_dv(m_attVrtVec, std::vector<Vertex*>());
970 :
971 : // m_attVrtVecFace = AttVrtVec();
972 : //
973 : // m_grid.attach_to_faces(m_attVrtVecFace, std::vector<Vertex*>());
974 : //
975 : // m_attAccsVrtVecFace = Grid::VolumeAttachmentAccessor<AttVrtVec>(m_grid, m_attVrtVecFace);
976 :
977 0 : m_attAccsVrtVecVol = Grid::VolumeAttachmentAccessor<AttVrtVec>(m_grid, m_attVrtVec);
978 :
979 : m_attMarkVrtxIsCenterVrtx = ABool();
980 :
981 0 : m_grid.attach_to_vertices_dv(m_attMarkVrtxIsCenterVrtx,false);
982 :
983 0 : m_attAccsVrtxIsCenterVrtx = Grid::VertexAttachmentAccessor<ABool>( m_grid, m_attMarkVrtxIsCenterVrtx);
984 :
985 : m_attMarkVrtxIsShiftVrtx = ABool();
986 :
987 0 : m_grid.attach_to_vertices_dv(m_attMarkVrtxIsShiftVrtx,false);
988 :
989 0 : m_attAccsVrtxIsShiftVrtx = Grid::VertexAttachmentAccessor<ABool>( m_grid, m_attMarkVrtxIsShiftVrtx);
990 :
991 : m_attMarkEdgeIsShiftEdge = ABool();
992 :
993 0 : m_grid.attach_to_edges_dv(m_attMarkEdgeIsShiftEdge,false);
994 :
995 0 : m_attAccsEdgeIsShiftEdge = Grid::EdgeAttachmentAccessor<ABool>( m_grid, m_attMarkEdgeIsShiftEdge );
996 :
997 : m_attInfoVecSudosTouchingVrtx = AttVecInt();
998 :
999 0 : m_grid.attach_to_vertices_dv(m_attInfoVecSudosTouchingVrtx, std::vector<IndexType>());
1000 :
1001 0 : m_attAccsInfoVecSudosTouchingVrtx = Grid::VertexAttachmentAccessor<AttVecInt>( m_grid, m_attInfoVecSudosTouchingVrtx);
1002 :
1003 : m_attVrtxIsMidPtOfShiftVrtx = ABool();
1004 :
1005 0 : m_grid.attach_to_vertices_dv(m_attVrtxIsMidPtOfShiftVrtx,false);
1006 :
1007 0 : m_attAccsVrtxIsMidPtOfShiftVrtcs = Grid::VertexAttachmentAccessor<ABool>(m_grid, m_attVrtxIsMidPtOfShiftVrtx);
1008 :
1009 : m_attMidPtVrtxOfShiftVrtx = AVertex();
1010 :
1011 0 : m_grid.attach_to_vertices_dv(m_attMidPtVrtxOfShiftVrtx, nullptr);
1012 :
1013 0 : m_attAccsMidPtVrtxOfShiftVrtx = Grid::VertexAttachmentAccessor<AVertex>( m_grid, m_attMidPtVrtxOfShiftVrtx );
1014 :
1015 : m_attCenterVrtxOfShiftVrtx = AVertex();
1016 :
1017 0 : m_grid.attach_to_vertices_dv( m_attCenterVrtxOfShiftVrtx, nullptr);
1018 :
1019 0 : m_attAccsCenterVrtxOfShiftVrtx = Grid::VertexAttachmentAccessor<AVertex>( m_grid, m_attCenterVrtxOfShiftVrtx );
1020 :
1021 : m_attEdgeCanBeRemoved = ABool();
1022 :
1023 0 : m_grid.attach_to_edges_dv(m_attEdgeCanBeRemoved, false);
1024 :
1025 0 : m_attAccsEdgeCanBeRemoved = Grid::EdgeAttachmentAccessor<ABool>( m_grid, m_attEdgeCanBeRemoved);
1026 :
1027 : m_attNewSudoOfVrtx = AInt();
1028 :
1029 0 : m_grid.attach_to_vertices_dv( m_attNewSudoOfVrtx, -1 );
1030 :
1031 0 : m_attAccsNewSudoOfVrtx = Grid::VertexAttachmentAccessor<AInt>( m_grid, m_attNewSudoOfVrtx );
1032 :
1033 0 : return true;
1034 : }
1035 :
1036 : ////////////////////////////////////////////////////////////////////////////
1037 :
1038 0 : bool DiamondsEstablish3D::detachMarkers()
1039 : {
1040 0 : m_grid.detach_from_vertices(m_attElmGrpVrtx2BQnchd);
1041 :
1042 0 : m_grid.detach_from_faces(m_attMarkTwinFace);
1043 0 : m_grid.detach_from_faces(m_attMarkShiftFace);
1044 :
1045 0 : m_grid.detach_from_faces(m_attMarkShiftTriangleFace);
1046 0 : m_grid.detach_from_faces(m_attMarkShiftQuadriliteralFace);
1047 :
1048 0 : m_grid.detach_from_volumes(m_attMarkVolGetsShrinked);
1049 :
1050 0 : m_grid.detach_from_volumes( m_attVrtVec );
1051 :
1052 : // m_grid.detach_from_volumes(m_attVrtVecFace);
1053 :
1054 0 : m_grid.detach_from_vertices(m_attMarkVrtxIsCenterVrtx);
1055 :
1056 0 : m_grid.detach_from_vertices(m_attMarkVrtxIsShiftVrtx);
1057 :
1058 0 : m_grid.detach_from_edges(m_attMarkEdgeIsShiftEdge);
1059 :
1060 0 : m_grid.detach_from_vertices(m_attInfoVecSudosTouchingVrtx);
1061 :
1062 0 : m_grid.detach_from_vertices(m_attVrtxIsMidPtOfShiftVrtx);
1063 :
1064 0 : m_grid.detach_from_vertices(m_attMidPtVrtxOfShiftVrtx);
1065 :
1066 0 : m_grid.detach_from_vertices(m_attCenterVrtxOfShiftVrtx);
1067 :
1068 0 : m_grid.detach_from_edges(m_attEdgeCanBeRemoved);
1069 :
1070 0 : m_grid.detach_from_vertices(m_attNewSudoOfVrtx);
1071 :
1072 0 : return true;
1073 : }
1074 :
1075 : ////////////////////////////////////////////////////////////////////////////
1076 :
1077 0 : bool DiamondsEstablish3D::assignBasicAtts()
1078 : {
1079 :
1080 0 : for( ElemGroupVrtx2BQuenched4Diams egv2bq : m_vecElemGroupVrtx2BQuenched )
1081 : {
1082 : Vertex * centerVrtx;
1083 : egv2bq.spuckOrigCenterVertex(centerVrtx);
1084 :
1085 0 : m_sel.select(centerVrtx);
1086 :
1087 0 : m_attAccsElmGrpVrtx2BQnchd[centerVrtx] = egv2bq;
1088 0 : m_attAccsVrtxIsCenterVrtx[centerVrtx] = true;
1089 : IndxVec sudoList = egv2bq.spuckSudoList();
1090 0 : m_attAccsInfoVecSudosTouchingVrtx[centerVrtx] = sudoList;
1091 :
1092 0 : m_sel.select( m_grid.associated_edges_begin(centerVrtx), m_grid.associated_edges_end(centerVrtx) );
1093 0 : m_sel.select( m_grid.associated_faces_begin(centerVrtx), m_grid.associated_faces_end(centerVrtx) );
1094 0 : m_sel.select( m_grid.associated_volumes_begin(centerVrtx), m_grid.associated_volumes_end(centerVrtx) );
1095 :
1096 0 : bool sudoAlreadyKnown = generateNewDiamSudos(centerVrtx, sudoList);
1097 :
1098 0 : }
1099 :
1100 0 : return true;
1101 : }
1102 :
1103 : ////////////////////////////////////////////////////////////////////////////
1104 :
1105 0 : bool DiamondsEstablish3D::trafoCollectedInfo2Attachments()
1106 : {
1107 : //for(VertexIterator iterV = m_sel.vertices_begin(); iterV != m_sel.vertices_end(); iterV++)
1108 0 : for( ElemGroupVrtx2BQuenched4Diams & egv2bQ : m_vecElemGroupVrtx2BQuenched )
1109 : {
1110 : Vertex * centerV;
1111 :
1112 : egv2bQ.spuckOrigCenterVertex(centerV);
1113 :
1114 : VecElems2BQuenched ve2bq;
1115 : egv2bQ.spuckVecElems2BQuenched4Diams(ve2bq);
1116 :
1117 0 : for( Elems2BQuenched & e2bq : ve2bq )
1118 : {
1119 : VecVolumeElementFaceQuintuplet vvef5;
1120 :
1121 : e2bq.spuckVecFullLowDimManifQuintuplet(vvef5);
1122 :
1123 0 : for( VolumeElementFaceQuintuplet & vef5 : vvef5 )
1124 : {
1125 0 : if( ! trafoQuintupleInfo2Attachments(vef5))
1126 : {
1127 : UG_LOG("Quintuples do not get to attachments" << std::endl);
1128 : return false;
1129 : }
1130 : }
1131 : // debugE2bQ(e2bq);
1132 :
1133 : Vertex * midPtVrtx;
1134 : e2bq.spuckMidPointOfShiftVrtcs(midPtVrtx);
1135 : m_attAccsVrtxIsMidPtOfShiftVrtcs[midPtVrtx] = true;
1136 :
1137 0 : m_sel.select(midPtVrtx);
1138 :
1139 : EdgePair ep;
1140 : e2bq.spuckPairLowDimElem(ep);
1141 : m_sel.select(ep.first);
1142 : m_sel.select(ep.second);
1143 :
1144 : VrtxPair vp;
1145 : e2bq.spuckShiftVrtcs(vp);
1146 : Vertex * shiftVrtxOne = vp.first;
1147 : Vertex * shiftVrtxTwo = vp.second;
1148 : m_sel.select(shiftVrtxOne);
1149 : m_sel.select(shiftVrtxTwo);
1150 :
1151 0 : m_attAccsMidPtVrtxOfShiftVrtx[shiftVrtxOne] = midPtVrtx;
1152 0 : m_attAccsMidPtVrtxOfShiftVrtx[shiftVrtxTwo] = midPtVrtx;
1153 :
1154 0 : m_attAccsCenterVrtxOfShiftVrtx[shiftVrtxOne] = centerV;
1155 0 : m_attAccsCenterVrtxOfShiftVrtx[shiftVrtxTwo] = centerV;
1156 :
1157 : // TODO FIXME hier die sudo Liste auch noch an den Center Vertex hängen als attachment
1158 :
1159 0 : }
1160 0 : }
1161 :
1162 :
1163 :
1164 : return true;
1165 : }
1166 :
1167 :
1168 : //////////////////////////////////////////////////////////////////////////////////
1169 :
1170 0 : bool DiamondsEstablish3D::trafoQuintupleInfo2Attachments(VolumeElementFaceQuintuplet & vef5 )
1171 : {
1172 : VrtxPair vp;
1173 : vef5.spuckShiftVrtcs(vp);
1174 :
1175 : m_attAccsVrtxIsShiftVrtx[vp.first] = true;
1176 : m_attAccsVrtxIsShiftVrtx[vp.second] = true;
1177 :
1178 : PairVolumeEdgeTwin pvet;
1179 : vef5.spuckPairFullLowDimTwin(pvet);
1180 :
1181 : VolumeEdgeTwin vetOne = pvet.first;
1182 : VolumeEdgeTwin vetTwo = pvet.second;
1183 :
1184 : Volume * volOne;
1185 : Volume * volTwo;
1186 : Edge * edgOne;
1187 : Edge * edgTwo;
1188 :
1189 : vetOne.spuckFullDimElem(volOne);
1190 : vetOne.spuckLowDimElem(edgOne);
1191 :
1192 : vetTwo.spuckFullDimElem(volTwo);
1193 : vetTwo.spuckLowDimElem(edgTwo);
1194 :
1195 : m_attAccsVolGetsShrinked[volOne] = true;
1196 : m_attAccsVolGetsShrinked[volTwo] = true;
1197 :
1198 : m_attAccsEdgeIsShiftEdge[edgOne] = true;
1199 : m_attAccsEdgeIsShiftEdge[edgTwo] = true;
1200 :
1201 0 : m_sel.select(volOne);
1202 : m_sel.select(volTwo);
1203 :
1204 : m_sel.select(edgOne);
1205 : m_sel.select(edgTwo);
1206 :
1207 : Face * fac;
1208 :
1209 : vef5.spuckManifElem(fac);
1210 :
1211 : m_attAccsFacIsTwinFac[fac] = true;
1212 :
1213 : m_sel.select(fac);
1214 :
1215 0 : return true;
1216 : }
1217 :
1218 : //////////////////////////////////////////////////////////////////////////////////
1219 :
1220 0 : bool DiamondsEstablish3D::assignMidPointOfShiftVrtcs(Elems2BQuenched & e2bq)
1221 : {
1222 : VrtxPair shiftVrtcs;
1223 :
1224 : e2bq.spuckShiftVrtcs(shiftVrtcs);
1225 :
1226 : Vertex * midPtVrtx;
1227 :
1228 0 : if( !createMidVrtx(shiftVrtcs, midPtVrtx))
1229 : {
1230 : UG_LOG("mid vertex already created " << std::endl);
1231 : // return false;
1232 : }
1233 :
1234 : e2bq.assignMidPointOfShiftVrtcs(midPtVrtx);
1235 :
1236 0 : return true;
1237 : }
1238 :
1239 : //////////////////////////////////////////////////////////////////////////////////
1240 :
1241 0 : bool DiamondsEstablish3D::createMidVrtx( VrtxPair const & shiftVrtcs, Vertex * & midPtVrtx )
1242 : {
1243 :
1244 : bool pairAlreadyCombined = false;
1245 :
1246 0 : for( CombiShiftVrtxMidVrtx & csvmv : m_vecCombiShiftVrtxMidVrtx )
1247 : {
1248 : VrtxPair prNoSwap;
1249 :
1250 : csvmv.spuckShiftVrtxPair(prNoSwap);
1251 :
1252 : if( prNoSwap == shiftVrtcs )
1253 : {
1254 : pairAlreadyCombined = true;
1255 : csvmv.spuckSinglVrtx(midPtVrtx);
1256 :
1257 : break;
1258 : }
1259 : else
1260 : {
1261 : std::swap( prNoSwap.first, prNoSwap.second );
1262 :
1263 : if( prNoSwap == shiftVrtcs )
1264 : {
1265 : pairAlreadyCombined = true;
1266 : csvmv.spuckSinglVrtx(midPtVrtx);
1267 :
1268 : break;
1269 : }
1270 :
1271 : if( pairAlreadyCombined )
1272 : break;
1273 :
1274 : }
1275 : }
1276 :
1277 : if( ! pairAlreadyCombined )
1278 : {
1279 :
1280 : Vertex * vrtxOne;
1281 : Vertex * vrtxTwo;
1282 :
1283 0 : vrtxOne = shiftVrtcs.first;
1284 0 : vrtxTwo = shiftVrtcs.second;
1285 :
1286 :
1287 : vector3 posOne = m_aaPos[vrtxOne];
1288 : vector3 posTwo = m_aaPos[vrtxTwo];
1289 :
1290 : vector3 sum;
1291 :
1292 : VecAdd(sum,posOne,posTwo);
1293 :
1294 : vector3 scal;
1295 :
1296 : VecScale(scal,sum,0.5);
1297 :
1298 0 : midPtVrtx = *m_grid.create<RegularVertex>();
1299 :
1300 : m_aaPos[midPtVrtx] = scal;
1301 :
1302 : CombiShiftVrtxMidVrtx csvmvNew( shiftVrtcs, midPtVrtx );
1303 :
1304 0 : m_vecCombiShiftVrtxMidVrtx.push_back(csvmvNew);
1305 : }
1306 :
1307 0 : return (! pairAlreadyCombined);
1308 :
1309 : }
1310 :
1311 : //////////////////////////////////////////////////////////////////////////////////
1312 :
1313 0 : bool DiamondsEstablish3D::createConditionForNewVrtcs()
1314 : {
1315 :
1316 : // iterate over all surrounding volumes to enable volume changes, this loop taken from SR but shortened
1317 0 : for(VolumeIterator iterSurrVol = m_sel.volumes_begin(); iterSurrVol != m_sel.volumes_end(); iterSurrVol++ )
1318 : {
1319 : Volume * sv = *iterSurrVol;
1320 :
1321 : std::vector<Vertex*> & newVrts = m_attAccsVrtVecVol[sv];
1322 0 : newVrts.resize(sv->num_vertices());
1323 :
1324 0 : for(size_t iVrt = 0; iVrt < sv->num_vertices(); iVrt++ )
1325 : {
1326 0 : newVrts[iVrt] = nullptr;
1327 : }
1328 : // erstmal so tun, als ob keine neuen Vertizes erzeugt werden an den alten Vertizes
1329 : }
1330 :
1331 : // for( FaceIterator iterFac = m_sel.faces_begin(); iterFac != m_sel.faces_end(); iterFac++)
1332 : // {
1333 : // Face * fac = *iterFac;
1334 : //
1335 : // std::vector<Vertex*> & newVrts = m_attAccsVrtVecFace[fac];
1336 : //
1337 : // newVrts.resize(fac->num_vertices());
1338 : //
1339 : // for(size_t iVrt = 0; iVrt < sv->num_vertices(); iVrt++ )
1340 : // {
1341 : // newVrts[iVrt] = nullptr;
1342 : // }
1343 : //
1344 : // }
1345 :
1346 0 : return true;
1347 : }
1348 :
1349 :
1350 : //////////////////////////////////////////////////////////////////////////////////
1351 :
1352 0 : bool DiamondsEstablish3D::distributeInfosForShrinkingVols()
1353 : {
1354 0 : for( ElemGroupVrtx2BQuenched4Diams & egv2bq : m_vecElemGroupVrtx2BQuenched )
1355 : {
1356 : VecElems2BQuenched ve2bq;
1357 :
1358 : egv2bq.spuckVecElems2BQuenched4Diams(ve2bq);
1359 :
1360 0 : for( Elems2BQuenched & e2bq : ve2bq )
1361 : {
1362 : VecVolumeElementFaceQuintuplet vvef5;
1363 :
1364 : e2bq.spuckVecFullLowDimManifQuintuplet(vvef5);
1365 :
1366 : Vertex * centerVrtx;
1367 :
1368 : egv2bq.spuckOrigCenterVertex(centerVrtx);
1369 :
1370 : Vertex * newMidVrtx;
1371 :
1372 : e2bq.spuckMidPointOfShiftVrtcs(newMidVrtx);
1373 :
1374 0 : for( VolumeElementFaceQuintuplet & vef5 : vvef5 )
1375 : {
1376 : PairVolumeEdgeTwin pvet;
1377 :
1378 : vef5.spuckPairFullLowDimTwin(pvet);
1379 :
1380 : VolumeEdgeTwin vetOne, vetTwo;
1381 :
1382 : vetOne = pvet.first;
1383 : vetTwo = pvet.second;
1384 :
1385 : Volume * volOne;
1386 : Volume * volTwo;
1387 :
1388 : vetOne.spuckFullDimElem(volOne);
1389 : vetTwo.spuckFullDimElem(volTwo);
1390 :
1391 0 : if( !teachMidVrtx2Vol(volOne, centerVrtx, newMidVrtx) || ! teachMidVrtx2Vol(volTwo,centerVrtx,newMidVrtx) )
1392 : {
1393 : UG_LOG("not taughtable" << std::endl);
1394 0 : return false;
1395 : }
1396 : }
1397 0 : }
1398 0 : }
1399 :
1400 : // TODO FIXME FACES detektieren, die verschoben werden müssen!!!
1401 :
1402 : return true;
1403 : }
1404 :
1405 :
1406 : //////////////////////////////////////////////////////////////////////////////////
1407 :
1408 0 : bool DiamondsEstablish3D::teachMidVrtx2Vol( Volume * const & vol, Vertex * const & origVrtx, Vertex * const & midVrtx )
1409 : {
1410 : IndexType taught = 0;
1411 :
1412 0 : std::vector<Vertex*> & newVrts4Fac = m_attAccsVrtVecVol[ vol ];
1413 :
1414 0 : for( IndexType indVrt = 0; indVrt < (vol)->num_vertices(); indVrt++ )
1415 : {
1416 0 : Vertex* volVrt = (vol)->vertex(indVrt);
1417 :
1418 0 : if( origVrtx == volVrt )
1419 : {
1420 0 : newVrts4Fac[ indVrt ] = midVrtx;
1421 0 : taught++;
1422 : }
1423 : }
1424 :
1425 0 : if( taught != 1 )
1426 : {
1427 : UG_LOG("strange problem with the volume " << taught << std::endl);
1428 0 : m_sh.assign_subset(vol, m_sh.num_subsets());
1429 0 : return false;
1430 : }
1431 :
1432 : return true;
1433 : }
1434 :
1435 : //////////////////////////////////////////////////////////////////////////////////
1436 :
1437 0 : bool DiamondsEstablish3D::shrinkVolumes()
1438 : {
1439 : // holds local side vertex indices
1440 : std::vector<size_t> locVrtInds;
1441 :
1442 : // first we create new edges from selected ones which are connected to
1443 : // inner vertices. This allows to preserve old subsets.
1444 : // Since we have to make sure that we use the right vertices,
1445 : // we have to iterate over the selected volumes and perform all actions on the edges
1446 : // of those volumes.
1447 :
1448 : int volNum = 0;
1449 :
1450 0 : for(VolumeIterator iter_sv = m_sel.volumes_begin(); iter_sv != m_sel.volumes_end(); ++iter_sv)
1451 : {
1452 : volNum++;
1453 :
1454 : Volume* sv = *iter_sv;
1455 :
1456 0 : UG_LOG("Diamond entering volume to create new elems " << CalculateCenter(sv, m_aaPos) << std::endl);
1457 :
1458 : // check for each edge whether it has to be copied.
1459 0 : for(size_t i_edge = 0; i_edge < sv->num_edges(); ++i_edge)
1460 : {
1461 0 : Edge* e = m_grid.get_edge(sv, i_edge);
1462 :
1463 0 : if(m_sel.is_selected(e))
1464 : {
1465 : // check the associated vertices through the volumes aaVrtVecVol attachment.
1466 : // If at least one has an associated new vertex and if no edge between the
1467 : // new vertices already exists, we'll create the new edge.
1468 : size_t ind0, ind1;
1469 0 : sv->get_vertex_indices_of_edge(ind0, ind1, i_edge);
1470 0 : Vertex* nv0 = (m_attAccsVrtVecVol[sv])[ind0];
1471 0 : Vertex* nv1 = (m_attAccsVrtVecVol[sv])[ind1];
1472 :
1473 0 : if(nv0 || nv1)
1474 : {
1475 : // if one vertex has no associated new one, then we use the vertex itself
1476 0 : if(!nv0)
1477 0 : nv0 = sv->vertex(ind0);
1478 0 : if(!nv1)
1479 0 : nv1 = sv->vertex(ind1);
1480 :
1481 : // create the new edge if it not already exists.
1482 0 : if( ! m_grid.get_edge(nv0, nv1))
1483 0 : m_grid.create_by_cloning(e, EdgeDescriptor(nv0, nv1), e);
1484 : }
1485 : }
1486 : }
1487 : }
1488 :
1489 :
1490 :
1491 : UG_LOG("Vol enter clone finished " << std::endl);
1492 :
1493 : // now we create new faces from selected ones which are connected to
1494 : // inner vertices. This allows to preserve old subsets.
1495 : // Since we have to make sure that we use the right vertices,
1496 : // we have to iterate over the selected volumes and perform all actions on the side-faces
1497 : // of those volumes.
1498 :
1499 0 : FaceDescriptor fd;
1500 :
1501 : int volNum2 = 0;
1502 :
1503 0 : for(VolumeIterator iter_sv = m_sel.volumes_begin(); iter_sv != m_sel.volumes_end(); ++iter_sv)
1504 : {
1505 : volNum2++;
1506 :
1507 : Volume* sv = *iter_sv;
1508 : // check for each face whether it has to be copied.
1509 :
1510 0 : UG_LOG("Diamond Face descriptor for vol " << CalculateCenter(sv, m_aaPos) << std::endl);
1511 :
1512 0 : for(size_t i_face = 0; i_face < sv->num_faces(); ++i_face)
1513 : {
1514 0 : Face* sf = m_grid.get_face(sv, i_face);
1515 :
1516 0 : if( m_sel.is_selected(sf))
1517 : {
1518 : // check the associated vertices through the volumes aaVrtVecVol attachment.
1519 : // If no face between the new vertices already exists, we'll create the new face.
1520 0 : sv->get_vertex_indices_of_face(locVrtInds, i_face);
1521 0 : fd.set_num_vertices(sf->num_vertices());
1522 :
1523 0 : for(size_t i = 0; i < sf->num_vertices(); ++i)
1524 : {
1525 0 : Vertex* nVrt = (m_attAccsVrtVecVol[sv])[locVrtInds[i]];
1526 :
1527 0 : if(nVrt)
1528 0 : fd.set_vertex(i, nVrt);
1529 : else
1530 0 : fd.set_vertex(i, sv->vertex(locVrtInds[i]));
1531 : }
1532 :
1533 : // if the new face does not already exist, we'll create it
1534 0 : if(!m_grid.get_face(fd))
1535 0 : m_grid.create_by_cloning(sf, fd, sf);
1536 : }
1537 : }
1538 : }
1539 :
1540 : UG_LOG("Face descriptor left" << std::endl);
1541 :
1542 : // Expand all faces.
1543 : // Since volumes are replaced on the fly, we have to take care with the iterator.
1544 : // record all new volumes in a vector. This will help to adjust positions later on.
1545 :
1546 : std::vector<Volume*> newFractureVolumes;
1547 : std::vector<IndexType> subsOfNewVolumes;
1548 :
1549 : // create alternative volumes where there are ending crossing clefts
1550 :
1551 0 : VolumeDescriptor vd;
1552 :
1553 : // beim Volumen fängt das Abfangen an, es muss abgefragt werden, ob das Volumen aus
1554 : // einem Segment ist, wo bisher falsch expandiert wird
1555 : // erst beim Face ab zu fangen ist viel zu spät TODO FIXME
1556 : // nicht zu vergessen, die Edges ordentlich sortiert zu sammeln, die zwei endende Vertizes enthalten,
1557 : // und dann zu splitten, wenn alle Attachements entfernt sind, dann Neustart anfordern mit neuer Geometrie
1558 :
1559 : int sudoNewVols = m_sh.num_subsets();
1560 :
1561 : // std::vector<Volume*> newVolsTwoCut;
1562 : // std::vector<Volume*> newVolsThreeCut;
1563 :
1564 :
1565 0 : for(VolumeIterator iter_sv = m_sel.volumes_begin(); iter_sv != m_sel.volumes_end();)
1566 : {
1567 :
1568 0 : Volume* sv = *iter_sv;
1569 : ++iter_sv;
1570 :
1571 :
1572 0 : UG_LOG("Diamond Volume new creation try at " << CalculateCenter(sv, m_aaPos) << std::endl);
1573 :
1574 :
1575 : // now expand the fracture faces of sv to volumes.
1576 0 : for(IndexType iSideVol = 0; iSideVol < sv->num_sides(); iSideVol++)
1577 : {
1578 : // get the local vertex indices of the side of the volume
1579 0 : sv->get_vertex_indices_of_face(locVrtInds, iSideVol);
1580 :
1581 0 : Face* sideFace = m_grid.get_side(sv, iSideVol);
1582 :
1583 0 : Volume * shiftVol = nullptr;
1584 :
1585 0 : if( sideFace )
1586 : {
1587 :
1588 0 : if( m_attAccsFacIsShiftFac[sideFace] )
1589 : {
1590 :
1591 0 : if( m_attAccsFacIsShiftQuadriliteralFac[sideFace] )
1592 : {
1593 :
1594 : std::vector<Vertex*> centerVrtcs;
1595 : std::vector<Vertex*> shiftVrtcs;
1596 : std::vector<Vertex*> midPtVrtcs;
1597 :
1598 :
1599 0 : for( IndexType vrtxInd = 0; vrtxInd < sideFace->num_vertices(); vrtxInd++ )
1600 : {
1601 0 : Vertex * sideVrtx = sideFace->vertex(vrtxInd);
1602 :
1603 0 : if( m_attAccsVrtxIsCenterVrtx[sideVrtx] )
1604 : {
1605 0 : centerVrtcs.push_back(sideVrtx);
1606 : }
1607 0 : else if( m_attAccsVrtxIsShiftVrtx[sideVrtx] )
1608 : {
1609 0 : shiftVrtcs.push_back(sideVrtx);
1610 : }
1611 : else
1612 : {
1613 0 : UG_LOG("not center not shift but shift face? "
1614 : << CalculateCenter( sideFace, m_aaPos ) << std::endl );
1615 0 : return false;
1616 : }
1617 : }
1618 :
1619 0 : if( centerVrtcs.size() != 2 || shiftVrtcs.size() != 2 )
1620 : {
1621 0 : UG_LOG("strange face behaviour for shift face of fac and vol "
1622 : << CalculateCenter( sideFace, m_aaPos )
1623 : << " ---- "
1624 : << CalculateCenter( sv, m_aaPos )
1625 : << std::endl);
1626 0 : return false;
1627 : }
1628 :
1629 : std::swap( shiftVrtcs[0], shiftVrtcs[1] );
1630 :
1631 0 : if( ! findShiftFaceVertices( sv, centerVrtcs, midPtVrtcs ))
1632 : {
1633 0 : UG_LOG("vertices of shift face strange "
1634 : << CalculateCenter( sideFace, m_aaPos ) << std::endl);
1635 0 : return false;
1636 : }
1637 :
1638 0 : UG_LOG("Diamentenerzeugung fuer " << CalculateCenter( sideFace, m_aaPos ) << std::endl);
1639 :
1640 : // centers.push_back(centerVrtcs);
1641 : // shifts.push_back(shiftVrtcs);
1642 : // midPts.push_back(midPtVrtcs);
1643 :
1644 : // // sehr sonderbar, wieso folgende Reihenfolge notwendig, damit es klappt:
1645 : //
1646 : // shiftVol = *m_grid.create<Prism>(
1647 : // PrismDescriptor( centerVrtcs[0], midPtVrtcs[0], shiftVrtcs[1],
1648 : // centerVrtcs[1], midPtVrtcs[1], shiftVrtcs[0]
1649 : // )
1650 : // );
1651 : // // sehr sonderbar, wieso folgende Reihenfolge notwendig, damit es klappt:
1652 :
1653 0 : shiftVol = *m_grid.create<Prism>(
1654 0 : PrismDescriptor( centerVrtcs[0], midPtVrtcs[0], shiftVrtcs[0],
1655 : centerVrtcs[1], midPtVrtcs[1], shiftVrtcs[1]
1656 : )
1657 : );
1658 :
1659 : int isThreeCross = -1;
1660 : IndexType foundThreeCross = 0;
1661 :
1662 : VrtxIndxCombi vrtxSudosCombiCenters;
1663 :
1664 0 : for( int i = 0; i < 2; i++ )
1665 : {
1666 0 : Vertex * cntrVrtx = centerVrtcs[i];
1667 :
1668 0 : IndxVec sudoVec = m_attAccsInfoVecSudosTouchingVrtx[cntrVrtx];
1669 :
1670 0 : IndexType sudoNums = sudoVec.size();
1671 :
1672 : VrtxIndxPair centerVrtxSudos( cntrVrtx, sudoVec );
1673 :
1674 0 : vrtxSudosCombiCenters.push_back(centerVrtxSudos);
1675 :
1676 0 : if( sudoNums < 2 )
1677 : {
1678 : UG_LOG("no cross but new prism " << std::endl);
1679 : return false;
1680 : }
1681 0 : else if( sudoNums == 2 )
1682 : {
1683 : ;
1684 : }
1685 0 : else if( sudoNums == 3 )
1686 : {
1687 : isThreeCross = i;
1688 0 : foundThreeCross++;
1689 : }
1690 : else
1691 : {
1692 : UG_LOG("strange sudo number " << std::endl);
1693 : return false;
1694 : }
1695 0 : }
1696 :
1697 0 : if( foundThreeCross > 0 )
1698 : {
1699 0 : if( foundThreeCross != 1 )
1700 : {
1701 : UG_LOG("three cross vertices too much " << std::endl);
1702 0 : return false;
1703 : }
1704 :
1705 0 : CombiNewVolsProps cnvp( shiftVol, sideFace, vrtxSudosCombiCenters, shiftVrtcs, midPtVrtcs, isThreeCross );
1706 :
1707 0 : m_vecCombiNewVolsThreeCross.push_back(cnvp);
1708 0 : }
1709 : else
1710 : {
1711 0 : CombiNewVolsProps cnvp( shiftVol, sideFace, vrtxSudosCombiCenters, shiftVrtcs, midPtVrtcs );
1712 :
1713 0 : m_vecCombiNewVolsTwoCross.push_back(cnvp);
1714 0 : }
1715 :
1716 : //m_attAccsInfoVecSudosTouchingVrtx
1717 :
1718 : // newVols.push_back(shiftVol);
1719 : // numFacs++;
1720 :
1721 : // if( numFacs == 1 )
1722 : // return true;
1723 : // Volume * shiftVol2 = *m_grid.create<Prism>(
1724 : // PrismDescriptor(
1725 : // shiftVrtcs[1], midPtVrtcs[1], centerVrtcs[1],
1726 : // shiftVrtcs[0], midPtVrtcs[0], centerVrtcs[0]
1727 : // )
1728 : // );
1729 :
1730 : // m_sh.assign_subset(shiftVol,sudoNewVols);
1731 : // m_sh.assign_subset(shiftVol2,m_sh.num_subsets()+2);
1732 :
1733 0 : UG_LOG("Diamentenerzeugung geschafft " << CalculateCenter( sideFace, m_aaPos ) << std::endl);
1734 :
1735 : // return true;
1736 :
1737 0 : }
1738 :
1739 : // if( shiftVol && shiftVol2 )
1740 : // {
1741 : // m_sh.assign_subset(shiftVol, sudoNewVols);
1742 : // m_sh.assign_subset(shiftVol2, sudoNewVols+1);
1743 : //
1744 : // m_sh.assign_subset( centerVrtcs[0], m_sh.num_subsets() );
1745 : // m_sh.assign_subset( centerVrtcs[1], m_sh.num_subsets() );
1746 : // m_sh.assign_subset( midPtVrtcs[0], m_sh.num_subsets() );
1747 : // m_sh.assign_subset( midPtVrtcs[1], m_sh.num_subsets() );
1748 : // m_sh.assign_subset( shiftVrtcs[0], m_sh.num_subsets() );
1749 : // m_sh.assign_subset( shiftVrtcs[1], m_sh.num_subsets() );
1750 : //
1751 : //
1752 : // return true;
1753 : // }
1754 : }
1755 : }
1756 :
1757 : // if( shiftVol )
1758 : // {
1759 : // m_sh.assign_subset(shiftVol,m_sh.num_subsets());
1760 : // }
1761 :
1762 : }
1763 :
1764 : // now set up a new volume descriptor and replace the volume.
1765 0 : if(vd.num_vertices() != sv->num_vertices())
1766 0 : vd.set_num_vertices(sv->num_vertices());
1767 :
1768 0 : for(size_t i_vrt = 0; i_vrt < sv->num_vertices(); ++i_vrt)
1769 : {
1770 0 : if( (m_attAccsVrtVecVol[sv])[i_vrt] )
1771 0 : vd.set_vertex(i_vrt, (m_attAccsVrtVecVol[sv])[i_vrt]);
1772 : else
1773 0 : vd.set_vertex(i_vrt, sv->vertex(i_vrt));
1774 : }
1775 :
1776 0 : m_grid.create_by_cloning(sv, vd, sv);
1777 0 : m_grid.erase(sv);
1778 :
1779 : }
1780 :
1781 : // for( Volume * v : newVols )
1782 : // {
1783 : // m_sh.assign_subset(v,m_sh.num_subsets());
1784 : // }
1785 :
1786 : UG_LOG("Volumes erzeugt " << std::endl);
1787 :
1788 0 : for( EdgeIterator itEdg = m_sel.begin<Edge>(); itEdg != m_sel.end<Edge>(); )
1789 : {
1790 : Edge * edg = *itEdg;
1791 : itEdg++;
1792 :
1793 0 : if( m_attAccsEdgeCanBeRemoved[edg] )
1794 : {
1795 0 : m_grid.erase(edg);
1796 : }
1797 : }
1798 :
1799 : UG_LOG("unnoetige Ecken entfernt " << std::endl);
1800 :
1801 0 : for(FaceIterator iter = m_sel.begin<Face>(); iter != m_sel.end<Face>();)
1802 : {
1803 : Face* fac = *iter;
1804 : ++iter;
1805 :
1806 0 : if( ! m_attAccsFacIsShiftFac[fac] )
1807 0 : m_grid.erase(fac);
1808 : }
1809 :
1810 : UG_LOG("Gesichter entfernt " << std::endl);
1811 :
1812 :
1813 : // int suse = m_sh.num_subsets();
1814 : //
1815 : // for( int i = 0; i < centers.size(); i++ )
1816 : // {
1817 : // std::vector<Vertex*> centerVrtcs = centers[i];
1818 : // std::vector<Vertex*> shiftVrtcs = shifts[i];
1819 : // std::vector<Vertex*> midPtVrtcs = midPts[i];
1820 : //
1821 : // UG_LOG("versuche neues Volumen nr " << i << std::endl);
1822 : //
1823 : // Volume * shiftVol = *m_grid.create<Prism>(
1824 : // PrismDescriptor( centerVrtcs[0], shiftVrtcs[0], midPtVrtcs[0],
1825 : // centerVrtcs[1], shiftVrtcs[1], midPtVrtcs[1]
1826 : // )
1827 : // );
1828 : //
1829 : // m_sh.assign_subset(shiftVol, suse);
1830 : //
1831 : // UG_LOG("neue Volumen " << i << CalculateCenter(shiftVol,m_aaPos) << std::endl);
1832 : //
1833 : // }
1834 :
1835 : UG_LOG("neue Volumen erzeugt " << std::endl);
1836 :
1837 : return true;
1838 0 : }
1839 :
1840 : //////////////////////////////////////////////////////////////////////////////////
1841 :
1842 0 : bool DiamondsEstablish3D::determineShiftFaces()
1843 : {
1844 :
1845 0 : for(FaceIterator iterF = m_sel.begin<Face>(); iterF != m_sel.end<Face>();)
1846 : {
1847 : Face* fac = *iterF;
1848 :
1849 : IndexType numShiftVrcs = 0;
1850 : IndexType numCentrVrtcs = 0;
1851 :
1852 0 : if( typeid(*fac) == typeid(Triangle) )
1853 : {
1854 0 : UG_LOG("we have a triangle " << CalculateCenter( fac, m_aaPos ) << " -> " << typeid(*fac).name() << std::endl);
1855 : }
1856 0 : else if( typeid(*fac) == typeid(Quadrilateral) )
1857 : {
1858 0 : UG_LOG("we have a quadriliteral " << CalculateCenter( fac, m_aaPos ) << " -> " << typeid(*fac).name() << std::endl);
1859 : }
1860 : else
1861 : {
1862 0 : UG_LOG("we have whatever at " << CalculateCenter( fac, m_aaPos ) << " -> " << typeid(*fac).name() << std::endl);
1863 : }
1864 :
1865 0 : for( IndexType iVrt = 0; iVrt < fac->num_vertices(); iVrt++ )
1866 : {
1867 0 : Vertex * vrt = fac->vertex(iVrt);
1868 :
1869 0 : if( m_attAccsVrtxIsCenterVrtx[vrt] )
1870 : {
1871 0 : numCentrVrtcs++;
1872 : }
1873 0 : else if( m_attAccsVrtxIsShiftVrtx[vrt] )
1874 : {
1875 0 : numShiftVrcs++;
1876 : }
1877 : }
1878 :
1879 0 : if( numShiftVrcs == 2 && numCentrVrtcs == 2 )
1880 : {
1881 : m_attAccsFacIsShiftFac[fac] = true;
1882 : m_attAccsFacIsShiftQuadriliteralFac[fac] = true;
1883 :
1884 0 : if( typeid(*fac) != typeid(Quadrilateral) )
1885 : {
1886 0 : UG_LOG("4 vertices but not Quadriliteral? " << CalculateCenter( fac, m_aaPos ) << " -> " << typeid(*fac).name() << std::endl);
1887 : }
1888 :
1889 : }
1890 0 : else if( numShiftVrcs == 2 && numCentrVrtcs == 1 )
1891 : {
1892 : m_attAccsFacIsShiftFac[fac] = true;
1893 : m_attAccsFacIsShiftTriangleFac[fac] = true;
1894 0 : UG_LOG("was ist das für ein Typ " << CalculateCenter( fac, m_aaPos ) << " -> " << typeid(*fac).name() << std::endl );
1895 : }
1896 :
1897 0 : if( numCentrVrtcs > 2 || numShiftVrcs > 2 )
1898 : {
1899 : UG_LOG("too much shift or center vertices in one face " << std::endl);
1900 0 : UG_LOG("weia was ist das für ein Typ " << CalculateCenter( fac, m_aaPos ) << " -> " << typeid(*fac).name() << std::endl );
1901 : return false;
1902 : }
1903 :
1904 : iterF++;
1905 : }
1906 :
1907 : return true;
1908 : }
1909 :
1910 : //////////////////////////////////////////////////////////////////////////////////
1911 :
1912 : //bool DiamondsEstablish3D::findShiftFaceVertices(
1913 : // std::vector<Vertex*> & centerVrtcs,
1914 : // std::vector<Vertex*> & shiftVrtcs,
1915 : // std::vector<Vertex*> & midPtVrtcs
1916 : //)
1917 : //{
1918 : //
1919 : // // figure out to which center vertex which shift vertex belongs
1920 : //
1921 : // if( centerVrtcs.size() != 2 || shiftVrtcs.size() != 2 )
1922 : // {
1923 : // UG_LOG("strange face behaviour for shift face" << std::endl);
1924 : //// << CalculateCenter( fac, m_aaPos ) << std::endl);
1925 : // return false;
1926 : // }
1927 : //
1928 : // if( ! checkAttsOfShiftFaceVrtcs( centerVrtcs, shiftVrtcs) )
1929 : // {
1930 : // std::swap( shiftVrtcs[0], shiftVrtcs[1] );
1931 : //
1932 : // if( ! checkAttsOfShiftFaceVrtcs( centerVrtcs, shiftVrtcs) )
1933 : // {
1934 : // UG_LOG("shift and swap vertices do not agree independent of ordering" << std::endl);
1935 : // return false;
1936 : // }
1937 : // }
1938 : //
1939 : // for( Vertex * & spv : shiftVrtcs )
1940 : // {
1941 : // midPtVrtcs.push_back( m_attAccsMidPtVrtxOfShiftVrtx[spv] );
1942 : // }
1943 : //
1944 : // return true;
1945 : //}
1946 :
1947 : ////////////////////////////////////////////////////////////////////////////////////
1948 :
1949 0 : bool DiamondsEstablish3D::checkAttsOfShiftFaceVrtcs( std::vector<Vertex*> const & centerVrtcs, std::vector<Vertex*> const & shiftVrtcs )
1950 : {
1951 : Vertex * const & centerVrtxOne = centerVrtcs[0];
1952 : Vertex * const & centerVrtxTwo = centerVrtcs[1];
1953 :
1954 : Vertex * const & shiftVrtxOne = shiftVrtcs[0];
1955 : Vertex * const & shiftVrtxTwo = shiftVrtcs[1];
1956 :
1957 0 : if( m_attAccsCenterVrtxOfShiftVrtx[shiftVrtxOne] == centerVrtxOne )
1958 : {
1959 0 : if( m_attAccsCenterVrtxOfShiftVrtx[shiftVrtxTwo] == centerVrtxTwo )
1960 : {
1961 : return true;
1962 : }
1963 : else
1964 : {
1965 : // UG_LOG("shift and center vertex do not fit together " << std::endl);
1966 0 : return false;
1967 : }
1968 : }
1969 :
1970 : return false;
1971 : }
1972 :
1973 : ////////////////////////////////////////////////////////////////////////////////////
1974 :
1975 0 : bool DiamondsEstablish3D::findShiftFaceVertices( Volume * & vol,
1976 : std::vector<Vertex*> & centerVrtcs,
1977 : std::vector<Vertex*> & midPtVrtcs
1978 : )
1979 : {
1980 0 : midPtVrtcs = std::vector<Vertex*>(2);
1981 :
1982 0 : std::vector<Vertex*> & newVrts4Fac = m_attAccsVrtVecVol[ vol ];
1983 :
1984 0 : for( IndexType indVrt = 0; indVrt < (vol)->num_vertices(); indVrt++ )
1985 : {
1986 0 : Vertex* volVrt = (vol)->vertex(indVrt);
1987 :
1988 0 : for( int i = 0; i < 2; i++ )
1989 : {
1990 0 : Vertex * ceVe = centerVrtcs[i];
1991 :
1992 0 : if( ceVe == volVrt )
1993 : {
1994 0 : midPtVrtcs[i] = newVrts4Fac[indVrt];
1995 : }
1996 : }
1997 : }
1998 :
1999 :
2000 0 : return true;
2001 : }
2002 :
2003 : ////////////////////////////////////////////////////////////////////////////////////
2004 :
2005 0 : bool DiamondsEstablish3D::postprocessNewDiamVols()
2006 : {
2007 0 : for( CombiNewVolsProps & cnvp : m_vecCombiNewVolsTwoCross )
2008 : {
2009 : VrtxIndxCombi viCombiCenter;
2010 : cnvp.spuckCenterVrtcsSudos(viCombiCenter);
2011 : Volume * vol;
2012 : cnvp.spuckFulldimElem(vol);
2013 :
2014 : VrtxIndxPair vrtxIndxPr( viCombiCenter[0] );
2015 :
2016 : Vertex * center = vrtxIndxPr.first;
2017 0 : IndxVec sudoVec = vrtxIndxPr.second;
2018 :
2019 0 : for( CombiCntrVrtxSudo & ccvs: m_vecCombiCntrVrtxSudo )
2020 : {
2021 0 : if( ccvs.spuckSudoVec() == sudoVec )
2022 : {
2023 : // m_sh.assign_subset(vol, ccvs.spuckNewSudo());
2024 0 : if( ! assignSudoOfNewVols2VolAndSubElems(vol,ccvs.spuckNewSudo() ) )
2025 : {
2026 : UG_LOG("sudos not assignable" << std::endl);
2027 : return false;
2028 : }
2029 : break;
2030 : }
2031 : }
2032 :
2033 0 : }
2034 :
2035 0 : for( CombiNewVolsProps & cnvp : m_vecCombiNewVolsThreeCross )
2036 : {
2037 : VrtxIndxCombi viCombiCenter;
2038 : cnvp.spuckCenterVrtcsSudos(viCombiCenter);
2039 : Volume * vol;
2040 : cnvp.spuckFulldimElem(vol);
2041 :
2042 0 : IndexType i3c = cnvp.spuckThreeCrossIndex();
2043 :
2044 0 : VrtxIndxPair vrtxIndxPr( viCombiCenter[i3c] );
2045 :
2046 : Vertex * center = vrtxIndxPr.first;
2047 0 : IndxVec sudoVec = vrtxIndxPr.second;
2048 :
2049 0 : for( CombiCntrVrtxSudo & ccvs: m_vecCombiCntrVrtxSudo )
2050 : {
2051 0 : if( ccvs.spuckSudoVec() == sudoVec )
2052 : {
2053 : // m_sh.assign_subset(vol, ccvs.spuckNewSudo());
2054 0 : if( ! assignSudoOfNewVols2VolAndSubElems(vol,ccvs.spuckNewSudo() ) )
2055 : {
2056 : UG_LOG("sudos not assignable" << std::endl);
2057 : return false;
2058 : }
2059 : break;
2060 : }
2061 : }
2062 :
2063 0 : }
2064 :
2065 : UG_LOG("try to split diam vols " << std::endl);
2066 :
2067 0 : for( CombiNewVolsProps & cnvp : m_vecCombiNewVolsThreeCross )
2068 : {
2069 0 : if( ! splitThreeCrossLargeDiams(cnvp) )
2070 : {
2071 : UG_LOG("splitting not possible");
2072 : return false;
2073 : }
2074 : }
2075 :
2076 0 : for( Face * fac : m_faces2BDeletedAtLastStep )
2077 : {
2078 0 : if( ! fac )
2079 : {
2080 : UG_LOG("want to delete null fac " << std::endl);
2081 : return false;
2082 : }
2083 :
2084 0 : m_grid.erase(fac);
2085 : }
2086 :
2087 0 : for( Edge * edg : m_edges2BDeletedAtLastStep )
2088 : {
2089 0 : if( ! edg )
2090 : {
2091 : UG_LOG("want to delete edge unexisting at end " << std::endl);
2092 : return false;
2093 : }
2094 :
2095 0 : m_grid.erase(edg);
2096 : }
2097 :
2098 :
2099 : UG_LOG("all diam vols splitted " << std::endl);
2100 :
2101 0 : return true;
2102 : }
2103 :
2104 : ////////////////////////////////////////////////////////////////////////////////////
2105 :
2106 0 : bool DiamondsEstablish3D::generateNewDiamSudos(Vertex * & centerV, IndxVec sudoList )
2107 : {
2108 0 : bool sudoCombiUnKnown = addElem(m_sudosTable, sudoList);
2109 :
2110 0 : if( ! sudoCombiUnKnown )
2111 : {
2112 0 : for( CombiCntrVrtxSudo & ccvs : m_vecCombiCntrVrtxSudo )
2113 : {
2114 0 : if( ccvs.spuckSudoVec() == sudoList )
2115 : {
2116 : ccvs.schluckVertex(centerV);
2117 0 : m_sh.assign_subset(centerV, ccvs.spuckNewSudo());
2118 0 : m_attAccsNewSudoOfVrtx[centerV] = ccvs.spuckNewSudo();
2119 0 : break;
2120 : }
2121 : }
2122 : }
2123 : else
2124 : {
2125 0 : IndexType newSudo = m_sh.num_subsets();
2126 0 : m_sh.assign_subset(centerV,newSudo);
2127 0 : m_attAccsNewSudoOfVrtx[centerV] = newSudo;
2128 : CombiCntrVrtxSudo ccvs( sudoList, newSudo );
2129 0 : m_vecCombiCntrVrtxSudo.push_back(ccvs);
2130 : UG_LOG("NNNNNNNNNNNNNNNNNNNN" << std::endl);
2131 : UG_LOG("creating new sudo of " << std::endl);
2132 0 : for( auto & i : sudoList )
2133 : {
2134 0 : UG_LOG("list part " << i << std::endl);
2135 : }
2136 : UG_LOG("DDDDDDDDDDDDDDDDDDDD" << std::endl);
2137 :
2138 : }
2139 :
2140 0 : return ( sudoCombiUnKnown );
2141 : }
2142 :
2143 : ////////////////////////////////////////////////////////////////////////////////////
2144 :
2145 0 : bool DiamondsEstablish3D::assignSudoOfNewVols2VolAndSubElems(Volume * & vol, IndexType sudo)
2146 : {
2147 0 : if( ! vol )
2148 : {
2149 : UG_LOG("vol to assign sudo null " << std::endl);
2150 0 : return false;
2151 : }
2152 :
2153 0 : m_sh.assign_subset(vol, sudo);
2154 :
2155 0 : for(IndexType iFace = 0; iFace < vol->num_faces(); ++iFace)
2156 : {
2157 0 : Face * fac = m_grid.get_face(vol, iFace);
2158 :
2159 0 : if( ! fac )
2160 : {
2161 : UG_LOG("face null sudo " << std::endl);
2162 0 : return false;
2163 : }
2164 :
2165 0 : m_sh.assign_subset( fac, sudo );
2166 :
2167 : }
2168 :
2169 0 : for(IndexType iEdge = 0; iEdge < vol->num_edges(); ++iEdge)
2170 : {
2171 0 : Edge* edg = m_grid.get_edge(vol, iEdge);
2172 :
2173 0 : if( ! edg )
2174 : {
2175 : UG_LOG("edge null sudo " << std::endl);
2176 0 : return false;
2177 : }
2178 :
2179 0 : m_sh.assign_subset( edg, sudo );
2180 :
2181 : }
2182 :
2183 0 : for( IndexType iVrt = 0; iVrt < vol->num_vertices(); iVrt++ )
2184 : {
2185 0 : Vertex * vrt = vol->vertex(iVrt);
2186 :
2187 0 : if( !vrt )
2188 : {
2189 : UG_LOG("vertex null sudo " << std::endl);
2190 0 : return false;
2191 : }
2192 :
2193 0 : m_sh.assign_subset( vrt, sudo );
2194 : }
2195 :
2196 : return true;
2197 : }
2198 :
2199 :
2200 : ///////////////////////////////////////////////////////////////////////////////////
2201 :
2202 0 : bool DiamondsEstablish3D::detectRemovableEdges()
2203 : {
2204 :
2205 0 : for( FaceIterator itFac = m_sel.begin<Face>(); itFac != m_sel.end<Face>(); itFac++)
2206 : {
2207 : Face * fac = *itFac;
2208 :
2209 : if( ! fac )
2210 : {
2211 : UG_LOG("face does not exist " << std::endl);
2212 : return false;
2213 : }
2214 :
2215 : UG_LOG("looking into twin facs try " << std::endl);
2216 :
2217 0 : if( m_attAccsFacIsTwinFac[fac] )
2218 : {
2219 : UG_LOG("inside twin facs try " << std::endl);
2220 :
2221 0 : for( IndexType iEdg = 0; iEdg < fac->num_edges(); iEdg++)
2222 : {
2223 : UG_LOG("inside edges of twins try " << std::endl);
2224 :
2225 0 : Edge * edg = m_grid.get_edge(fac,iEdg);
2226 :
2227 0 : if( ! edg )
2228 : {
2229 : UG_LOG("edge does not exist " << std::endl);
2230 0 : return false;
2231 : }
2232 :
2233 : IndexType numAssoCenterVrcs = 0;
2234 :
2235 0 : for( IndexType iV = 0; iV < 2; iV++ )
2236 : {
2237 0 : Vertex * vrt = edg->vertex(iV);
2238 :
2239 0 : if( ! vrt )
2240 : {
2241 : UG_LOG("vertex does not exist " << std::endl);
2242 0 : return false;
2243 : }
2244 :
2245 0 : if( m_attAccsVrtxIsCenterVrtx[vrt] )
2246 : {
2247 0 : numAssoCenterVrcs++;
2248 : }
2249 : }
2250 :
2251 : UG_LOG("number of associatec center vertices is " << numAssoCenterVrcs << std::endl);
2252 :
2253 0 : if( numAssoCenterVrcs == 1 )
2254 : {
2255 : m_attAccsEdgeCanBeRemoved[edg] = true;
2256 : }
2257 :
2258 : }
2259 : }
2260 : }
2261 :
2262 0 : return true;
2263 : }
2264 :
2265 : ///////////////////////////////////////////////////////////////////////////////////
2266 :
2267 0 : bool DiamondsEstablish3D::splitThreeCrossLargeDiams( CombiNewVolsProps & combiNewVolsProps )
2268 : {
2269 : VrtxIndxCombi viCombiCenter;
2270 : combiNewVolsProps.spuckCenterVrtcsSudos(viCombiCenter);
2271 : Volume * vol;
2272 : combiNewVolsProps.spuckFulldimElem(vol);
2273 :
2274 0 : IndexType ind3CrossSide = combiNewVolsProps.spuckThreeCrossIndex();
2275 :
2276 : IndexType indOtherSide = ( ind3CrossSide + 1 ) % 2;
2277 :
2278 0 : VrtxIndxPair diamCenterVrtxPr( viCombiCenter[ind3CrossSide] );
2279 :
2280 0 : Vertex * centerDiam = diamCenterVrtxPr.first;
2281 :
2282 0 : VrtxIndxPair diamOtherVrtxPr( viCombiCenter[ indOtherSide ] );
2283 :
2284 0 : Vertex * centerOther = diamOtherVrtxPr.first;
2285 :
2286 : std::vector<Vertex*> midPtVrtcs;
2287 : combiNewVolsProps.spuckMidPtVrtcs(midPtVrtcs);
2288 :
2289 0 : Vertex * midPtDiam = midPtVrtcs[ind3CrossSide];
2290 0 : Vertex * midPtOther = midPtVrtcs[indOtherSide];
2291 :
2292 : vector3 centerDiamVec3 = m_aaPos[centerDiam];
2293 : vector3 midPtDiamVec3 = m_aaPos[midPtDiam];
2294 : vector3 centerOtherVec3 = m_aaPos[centerOther];
2295 :
2296 : vector3 cutPtVec3;
2297 :
2298 0 : DropAPerpendicular(cutPtVec3, midPtDiamVec3, centerDiamVec3, centerOtherVec3);
2299 :
2300 0 : Vertex * cutVrtx = *m_grid.create<RegularVertex>();
2301 :
2302 : m_aaPos[cutVrtx] = cutPtVec3;
2303 :
2304 : combiNewVolsProps.schluckNewSplitVrtx(cutVrtx);
2305 :
2306 : std::vector<Vertex*> shiftVrtcs;
2307 : combiNewVolsProps.spuckShiftVrtcs(shiftVrtcs);
2308 :
2309 0 : Vertex * shiftVrtxDiam = shiftVrtcs[ind3CrossSide];
2310 0 : Vertex * shiftVrtxOther = shiftVrtcs[ indOtherSide];
2311 :
2312 0 : Volume * splitVolPrism = *m_grid.create<Prism>(
2313 0 : PrismDescriptor( cutVrtx, midPtDiam, shiftVrtxDiam,
2314 : centerOther, midPtOther, shiftVrtxOther
2315 : )
2316 0 : );
2317 :
2318 0 : assignSudoOfNewVols2VolAndSubElems(splitVolPrism, m_attAccsNewSudoOfVrtx[centerOther]);
2319 : // m_sh.assign_subset( splitVolPrism, m_attAccsNewSudoOfVrtx[centerOther] );
2320 :
2321 0 : Volume * splitVolTetra = *m_grid.create<Tetrahedron>(
2322 0 : TetrahedronDescriptor( shiftVrtxDiam, midPtDiam, cutVrtx, centerDiam ) );
2323 :
2324 : // m_sh.assign_subset( splitVolTetra, m_attAccsNewSudoOfVrtx[centerDiam] );
2325 0 : assignSudoOfNewVols2VolAndSubElems(splitVolTetra, m_attAccsNewSudoOfVrtx[centerDiam]);
2326 :
2327 : // figure out the faces and the edges which later should be erased
2328 :
2329 : IndexType numEdgesFound = 0;
2330 :
2331 : Edge* centersEdg = nullptr;
2332 :
2333 0 : for( IndexType iEdg = 0; iEdg < vol->num_edges(); iEdg++)
2334 : {
2335 0 : Edge * edg = m_grid.get_edge(vol, iEdg);
2336 :
2337 0 : if( EdgeContains(edg, centerDiam ) && EdgeContains(edg, centerOther))
2338 : {
2339 : centersEdg = edg;
2340 0 : numEdgesFound++;
2341 : }
2342 : }
2343 :
2344 0 : if( numEdgesFound != 1 )
2345 : {
2346 : UG_LOG("edges found at end " << numEdgesFound << std::endl);
2347 : return false;
2348 : }
2349 :
2350 0 : addElem(m_edges2BDeletedAtLastStep,centersEdg);
2351 :
2352 : std::vector<Face*> faces2Del;
2353 :
2354 0 : for( IndexType iFac = 0; iFac < vol->num_faces(); iFac++)
2355 : {
2356 0 : Face * fac = m_grid.get_face(vol, iFac);
2357 :
2358 0 : if( FaceContains(fac, centersEdg) )
2359 : {
2360 0 : faces2Del.push_back(fac);
2361 : }
2362 : }
2363 :
2364 0 : if( faces2Del.size() != 2 )
2365 : {
2366 : UG_LOG("strange number of faces to del at fin " << faces2Del.size() << std::endl);
2367 : return false;
2368 : }
2369 :
2370 0 : for( Face * f : faces2Del )
2371 : {
2372 0 : addElem(m_faces2BDeletedAtLastStep, f);
2373 : }
2374 :
2375 0 : m_grid.erase(vol);
2376 :
2377 :
2378 : return true;
2379 0 : }
2380 :
2381 :
2382 : ///////////////////////////////////////////////////////////////////////////////////
2383 :
2384 :
2385 : } /* namespace diamonds */
2386 :
2387 : } /* namespace arte */
2388 :
2389 : } /* namespace ug */
|
