meshRefinement.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2015-2022 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "meshRefinement.H"
30 #include "volMesh.H"
31 #include "volFields.H"
32 #include "surfaceMesh.H"
33 #include "syncTools.H"
34 #include "Time.H"
35 #include "refinementSurfaces.H"
36 #include "refinementFeatures.H"
37 #include "decompositionMethod.H"
38 #include "regionSplit.H"
39 #include "fvMeshDistribute.H"
40 #include "indirectPrimitivePatch.H"
41 #include "polyTopoChange.H"
42 #include "removeCells.H"
43 #include "mapDistributePolyMesh.H"
44 #include "localPointRegion.H"
45 #include "pointMesh.H"
46 #include "pointFields.H"
47 #include "slipPointPatchFields.H"
48 #include "fixedValuePointPatchFields.H"
49 #include "calculatedPointPatchFields.H"
50 #include "cyclicSlipPointPatchFields.H"
51 #include "processorPointPatch.H"
52 #include "globalIndex.H"
53 #include "meshTools.H"
54 #include "OFstream.H"
55 #include "Random.H"
56 #include "searchableSurfaces.H"
57 #include "treeBoundBox.H"
58 #include "zeroGradientFvPatchFields.H"
59 #include "fvMeshTools.H"
60 #include "motionSmoother.H"
61 #include "faceSet.H"
62 #include "topoDistanceData.H"
63 #include "FaceCellWave.H"
64 #include "PackedBoolList.H"
65 
66 // Leak path
67 #include "shortestPathSet.H"
68 #include "meshSearch.H"
69 
70 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
71 
72 namespace Foam
73 {
74  defineTypeNameAndDebug(meshRefinement, 0);
75 }
76 
77 
78 const Foam::Enum
79 <
80  Foam::meshRefinement::debugType
81 >
82 Foam::meshRefinement::debugTypeNames
83 ({
84  { debugType::MESH, "mesh" },
85  { debugType::OBJINTERSECTIONS, "intersections" },
86  { debugType::FEATURESEEDS, "featureSeeds" },
87  { debugType::ATTRACTION, "attraction" },
88  { debugType::LAYERINFO, "layerInfo" },
89 });
90 
91 
92 //const Foam::Enum
93 //<
94 // Foam::meshRefinement::outputType
95 //>
96 //Foam::meshRefinement::outputTypeNames
97 //({
98 // { outputType::OUTPUTLAYERINFO, "layerInfo" }
99 //});
100 
101 
102 const Foam::Enum
103 <
104  Foam::meshRefinement::writeType
105 >
106 Foam::meshRefinement::writeTypeNames
107 ({
108  { writeType::WRITEMESH, "mesh" },
109  { writeType::NOWRITEREFINEMENT, "noRefinement" },
110  { writeType::WRITELEVELS, "scalarLevels" },
111  { writeType::WRITELAYERSETS, "layerSets" },
112  { writeType::WRITELAYERFIELDS, "layerFields" },
113 });
114 
115 
116 Foam::meshRefinement::writeType Foam::meshRefinement::writeLevel_;
117 
118 //Foam::meshRefinement::outputType Foam::meshRefinement::outputLevel_;
119 
120 // Inside/outside test for polyMesh:.findCell()
121 // 2.4.x : default = polyMesh::FACE_DIAG_TRIS
122 // 1712 : default = polyMesh::CELL_TETS
123 //
124 // - CELL_TETS is better with concave cells, but much slower.
125 // - use faster method (FACE_DIAG_TRIS) here
126 
127 static const Foam::polyMesh::cellDecomposition
128  findCellMode(Foam::polyMesh::FACE_DIAG_TRIS);
129 
130 
131 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
132 
133 Foam::label Foam::meshRefinement::globalFaceCount(const labelList& elems) const
134 {
135  //- Check for duplicates
136  const bitSet isElem(mesh_.nFaces(), elems);
137  if (label(isElem.count()) != elems.size())
138  {
139  FatalErrorInFunction << "Problem Duplicates:"
140  << " isElem:" << isElem.count()
141  << " elems:" << elems.size()
142  << exit(FatalError);
143  }
144 
145  //- Check for same entries on coupled faces
146  {
147  bitSet isElem2(isElem);
148  syncTools::swapFaceList(mesh_, isElem2);
149 
150  for
151  (
152  label facei = mesh_.nInternalFaces();
153  facei < mesh_.nFaces();
154  facei++
155  )
156  {
157  if (isElem2[facei] != isElem[facei])
158  {
159  FatalErrorInFunction
160  << "at face:" << facei
161  << " at:" << mesh_.faceCentres()[facei]
162  << " patch:" << mesh_.boundaryMesh().whichPatch(facei)
163  << " isElem:" << isElem[facei]
164  << " isElem2:" << isElem2[facei]
165  << exit(FatalError);
166  }
167  }
168  }
169 
170  //- Count number of master elements
171  const bitSet isMaster(syncTools::getMasterFaces(mesh_));
172  label count = 0;
173  for (const label i : isElem)
174  {
175  if (isMaster[i])
176  {
177  count++;
178  }
179  }
180  return returnReduce(count, sumOp<label>());
181 }
182 
183 
184 void Foam::meshRefinement::calcNeighbourData
185 (
186  labelList& neiLevel,
187  pointField& neiCc
188 ) const
189 {
190  const labelList& cellLevel = meshCutter_.cellLevel();
191  const pointField& cellCentres = mesh_.cellCentres();
192 
193  const label nBoundaryFaces = mesh_.nBoundaryFaces();
194 
195  if (neiLevel.size() != nBoundaryFaces || neiCc.size() != nBoundaryFaces)
196  {
197  FatalErrorInFunction
198  << nBoundaryFaces << " neiLevel:" << neiLevel.size()
199  << abort(FatalError);
200  }
201 
202  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
203 
204  labelHashSet addedPatchIDSet(meshedPatches());
205 
206  forAll(patches, patchi)
207  {
208  const polyPatch& pp = patches[patchi];
209 
210  const labelUList& faceCells = pp.faceCells();
211  const vectorField::subField faceCentres = pp.faceCentres();
212  const vectorField::subField faceAreas = pp.faceAreas();
213 
214  label bFacei = pp.start()-mesh_.nInternalFaces();
215 
216  if (pp.coupled())
217  {
218  forAll(faceCells, i)
219  {
220  neiLevel[bFacei] = cellLevel[faceCells[i]];
221  neiCc[bFacei] = cellCentres[faceCells[i]];
222  bFacei++;
223  }
224  }
225  else if (addedPatchIDSet.found(patchi))
226  {
227  // Face was introduced from cell-cell intersection. Try to
228  // reconstruct other side cell(centre). Three possibilities:
229  // - cells same size.
230  // - preserved cell smaller. Not handled.
231  // - preserved cell larger.
232  forAll(faceCells, i)
233  {
234  // Extrapolate the face centre.
235  const vector fn = normalised(faceAreas[i]);
236 
237  label own = faceCells[i];
238  label ownLevel = cellLevel[own];
239  label faceLevel = meshCutter_.faceLevel(pp.start()+i);
240  if (faceLevel < 0)
241  {
242  // Due to e.g. face merging no longer a consistent
243  // refinementlevel of face. Assume same as cell
244  faceLevel = ownLevel;
245  }
246 
247  // Normal distance from face centre to cell centre
248  scalar d = ((faceCentres[i] - cellCentres[own]) & fn);
249  if (faceLevel > ownLevel)
250  {
251  // Other cell more refined. Adjust normal distance
252  d *= 0.5;
253  }
254  neiLevel[bFacei] = faceLevel;
255  // Calculate other cell centre by extrapolation
256  neiCc[bFacei] = faceCentres[i] + d*fn;
257  bFacei++;
258  }
259  }
260  else
261  {
262  forAll(faceCells, i)
263  {
264  neiLevel[bFacei] = cellLevel[faceCells[i]];
265  neiCc[bFacei] = faceCentres[i];
266  bFacei++;
267  }
268  }
269  }
270 
271  // Swap coupled boundaries. Apply separation to cc since is coordinate.
272  syncTools::swapBoundaryFacePositions(mesh_, neiCc);
273  syncTools::swapBoundaryFaceList(mesh_, neiLevel);
274 }
275 
276 
277 void Foam::meshRefinement::calcCellCellRays
278 (
279  const pointField& neiCc,
280  const labelList& neiLevel,
281  const labelList& testFaces,
282  pointField& start,
283  pointField& end,
284  labelList& minLevel
285 ) const
286 {
287  const labelList& cellLevel = meshCutter_.cellLevel();
288  const pointField& cellCentres = mesh_.cellCentres();
289 
290 
291  // Mark all non-coupled or coupled+master faces. Leaves only slave of
292  // coupled unset.
293  bitSet isMaster(mesh_.nBoundaryFaces(), true);
294  {
295  for (const polyPatch& pp : mesh_.boundaryMesh())
296  {
297  if (pp.coupled() && !refCast<const coupledPolyPatch>(pp).owner())
298  {
299  isMaster.unset(labelRange(pp.offset(), pp.size()));
300  }
301  }
302  }
303 
304 
305  start.setSize(testFaces.size());
306  end.setSize(testFaces.size());
307  minLevel.setSize(testFaces.size());
308 
309  forAll(testFaces, i)
310  {
311  const label facei = testFaces[i];
312  const label own = mesh_.faceOwner()[facei];
313 
314  if (mesh_.isInternalFace(facei))
315  {
316  const label nei = mesh_.faceNeighbour()[facei];
317 
318  start[i] = cellCentres[own];
319  end[i] = cellCentres[nei];
320  minLevel[i] = min(cellLevel[own], cellLevel[nei]);
321  }
322  else
323  {
324  const label bFacei = facei - mesh_.nInternalFaces();
325 
326  if (isMaster[bFacei])
327  {
328  start[i] = cellCentres[own];
329  end[i] = neiCc[bFacei];
330  }
331  else
332  {
333  // Slave face
334  start[i] = neiCc[bFacei];
335  end[i] = cellCentres[own];
336  }
337  minLevel[i] = min(cellLevel[own], neiLevel[bFacei]);
338  }
339  }
340 
341  // Extend segments a bit
342  {
343  const vectorField smallVec(ROOTSMALL*(end-start));
344  start -= smallVec;
345  end += smallVec;
346  }
347 }
348 
349 
350 void Foam::meshRefinement::updateIntersections(const labelList& changedFaces)
351 {
352  // Stats on edges to test. Count proc faces only once.
353  bitSet isMasterFace(syncTools::getMasterFaces(mesh_));
354 
355  {
356  const label nMasterFaces =
357  returnReduce(isMasterFace.count(), sumOp<label>());
358 
359  label nChangedFaces = 0;
360  forAll(changedFaces, i)
361  {
362  if (isMasterFace.test(changedFaces[i]))
363  {
364  ++nChangedFaces;
365  }
366  }
367  reduce(nChangedFaces, sumOp<label>());
368 
369  if (!dryRun_)
370  {
371  Info<< "Edge intersection testing:" << nl
372  << " Number of edges : " << nMasterFaces << nl
373  << " Number of edges to retest : " << nChangedFaces
374  << endl;
375  }
376  }
377 
378 
379  // Get boundary face centre and level. Coupled aware.
380  labelList neiLevel(mesh_.nBoundaryFaces());
381  pointField neiCc(mesh_.nBoundaryFaces());
382  calcNeighbourData(neiLevel, neiCc);
383 
384  // Collect segments we want to test for
385  pointField start(changedFaces.size());
386  pointField end(changedFaces.size());
387  {
388  labelList minLevel;
389  calcCellCellRays
390  (
391  neiCc,
392  neiLevel,
393  changedFaces,
394  start,
395  end,
396  minLevel
397  );
398  }
399 
400 
401  // Do tests in one go
402  labelList surfaceHit;
403  {
404  labelList surfaceLevel;
405  surfaces_.findHigherIntersection
406  (
407  shells_,
408  start,
409  end,
410  labelList(start.size(), -1), // accept any intersection
411  surfaceHit,
412  surfaceLevel
413  );
414  }
415 
416  // Keep just surface hit
417  forAll(surfaceHit, i)
418  {
419  surfaceIndex_[changedFaces[i]] = surfaceHit[i];
420  }
421 
422  // Make sure both sides have same information. This should be
423  // case in general since same vectors but just to make sure.
424  syncTools::syncFaceList(mesh_, surfaceIndex_, maxEqOp<label>());
425 
426  label nTotHits = returnReduce(countHits(), sumOp<label>());
427 
428  if (!dryRun_)
429  {
430  Info<< " Number of intersected edges : " << nTotHits << endl;
431  }
432 
433  // Set files to same time as mesh
434  setInstance(mesh_.facesInstance());
435 }
436 
437 
438 void Foam::meshRefinement::nearestFace
439 (
440  const labelUList& startFaces,
441  const bitSet& isBlockedFace,
442 
443  autoPtr<mapDistribute>& mapPtr,
444  labelList& faceToStart,
445  const label nIter
446 ) const
447 {
448  // From startFaces walk out (but not through isBlockedFace). Returns
449  // faceToStart which is the index into startFaces (or rather distributed
450  // version of it). E.g.
451  // pointField startFc(mesh.faceCentres(), startFaces);
452  // mapPtr().distribute(startFc);
453  // forAll(faceToStart, facei)
454  // {
455  // const label sloti = faceToStart[facei];
456  // if (sloti != -1)
457  // {
458  // Pout<< "face:" << mesh.faceCentres()[facei]
459  // << " nearest:" << startFc[sloti]
460  // << endl;
461  // }
462  // }
463 
464 
465  // Consecutive global numbering for start elements
466  const globalIndex globalStart(startFaces.size());
467 
468  // Field on cells and faces.
469  List<topoDistanceData<label>> cellData(mesh_.nCells());
470  List<topoDistanceData<label>> faceData(mesh_.nFaces());
471 
472  // Mark blocked faces to there not visited
473  for (const label facei : isBlockedFace)
474  {
475  faceData[facei] = topoDistanceData<label>(0, -1);
476  }
477 
478  List<topoDistanceData<label>> startData(startFaces.size());
479  forAll(startFaces, i)
480  {
481  const label facei = startFaces[i];
482  if (isBlockedFace[facei])
483  {
484  FatalErrorInFunction << "Start face:" << facei
485  << " at:" << mesh_.faceCentres()[facei]
486  << " is also blocked" << exit(FatalError);
487  }
488  startData[i] = topoDistanceData<label>(0, globalStart.toGlobal(i));
489  }
490 
491 
492  // Propagate information inwards
493  FaceCellWave<topoDistanceData<label>> deltaCalc
494  (
495  mesh_,
496  startFaces,
497  startData,
498  faceData,
499  cellData,
500  0
501  );
502  deltaCalc.iterate(nIter);
503 
504  // And extract
505 
506  faceToStart.setSize(mesh_.nFaces());
507  faceToStart = -1;
508  bool haveWarned = false;
509  forAll(faceData, facei)
510  {
511  if (!faceData[facei].valid(deltaCalc.data()))
512  {
513  if (!haveWarned)
514  {
515  WarningInFunction
516  << "Did not visit some faces, e.g. face " << facei
517  << " at " << mesh_.faceCentres()[facei]
518  << endl;
519  haveWarned = true;
520  }
521  }
522  else
523  {
524  faceToStart[facei] = faceData[facei].data();
525  }
526  }
527 
528  // Compact
529  List<Map<label>> compactMap;
530  mapPtr.reset(new mapDistribute(globalStart, faceToStart, compactMap));
531 }
532 
533 
534 void Foam::meshRefinement::nearestPatch
535 (
536  const labelList& adaptPatchIDs,
537  labelList& nearestPatch,
538  labelList& nearestZone
539 ) const
540 {
541  // Determine nearest patch for all mesh faces. Used when removing cells
542  // to give some reasonable patch to exposed faces.
543 
544  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
545 
546  nearestZone.setSize(mesh_.nFaces(), -1);
547 
548  if (adaptPatchIDs.size())
549  {
550  nearestPatch.setSize(mesh_.nFaces(), adaptPatchIDs[0]);
551 
552 
553  // Get per-face the zone or -1
554  labelList faceToZone(mesh_.nFaces(), -1);
555  {
556  for (const faceZone& fz : mesh_.faceZones())
557  {
558  UIndirectList<label>(faceToZone, fz) = fz.index();
559  }
560  }
561 
562 
563  // Count number of faces in adaptPatchIDs
564  label nFaces = 0;
565  forAll(adaptPatchIDs, i)
566  {
567  const polyPatch& pp = patches[adaptPatchIDs[i]];
568  nFaces += pp.size();
569  }
570 
571  // Field on cells and faces.
572  List<topoDistanceData<labelPair>> cellData(mesh_.nCells());
573  List<topoDistanceData<labelPair>> faceData(mesh_.nFaces());
574 
575  // Start of changes
576  labelList patchFaces(nFaces);
577  List<topoDistanceData<labelPair>> patchData(nFaces);
578  nFaces = 0;
579  forAll(adaptPatchIDs, i)
580  {
581  label patchi = adaptPatchIDs[i];
582  const polyPatch& pp = patches[patchi];
583 
584  forAll(pp, i)
585  {
586  patchFaces[nFaces] = pp.start()+i;
587  patchData[nFaces] = topoDistanceData<labelPair>
588  (
589  0,
590  labelPair
591  (
592  patchi,
593  faceToZone[pp.start()+i]
594  )
595  );
596  nFaces++;
597  }
598  }
599 
600  // Propagate information inwards
601  FaceCellWave<topoDistanceData<labelPair>> deltaCalc
602  (
603  mesh_,
604  patchFaces,
605  patchData,
606  faceData,
607  cellData,
608  mesh_.globalData().nTotalCells()+1
609  );
610 
611  // And extract
612 
613  bool haveWarned = false;
614  forAll(faceData, facei)
615  {
616  if (!faceData[facei].valid(deltaCalc.data()))
617  {
618  if (!haveWarned)
619  {
620  WarningInFunction
621  << "Did not visit some faces, e.g. face " << facei
622  << " at " << mesh_.faceCentres()[facei] << endl
623  << "Assigning these faces to patch "
624  << adaptPatchIDs[0]
625  << endl;
626  haveWarned = true;
627  }
628  }
629  else
630  {
631  const labelPair& data = faceData[facei].data();
632  nearestPatch[facei] = data.first();
633  nearestZone[facei] = data.second();
634  }
635  }
636  }
637  else
638  {
639  // Use patch 0
640  nearestPatch.setSize(mesh_.nFaces(), 0);
641  }
642 }
643 
644 
645 Foam::labelList Foam::meshRefinement::nearestPatch
646 (
647  const labelList& adaptPatchIDs
648 ) const
649 {
650  labelList nearestAdaptPatch;
651  labelList nearestAdaptZone;
652  nearestPatch(adaptPatchIDs, nearestAdaptPatch, nearestAdaptZone);
653  return nearestAdaptPatch;
654 }
655 
656 
657 void Foam::meshRefinement::nearestIntersection
658 (
659  const labelList& surfacesToTest,
660  const labelList& testFaces,
661 
662  labelList& surface1,
663  List<pointIndexHit>& hit1,
664  labelList& region1,
665  labelList& surface2,
666  List<pointIndexHit>& hit2,
667  labelList& region2
668 ) const
669 {
670  // Swap neighbouring cell centres and cell level
671  labelList neiLevel(mesh_.nBoundaryFaces());
672  pointField neiCc(mesh_.nBoundaryFaces());
673  calcNeighbourData(neiLevel, neiCc);
674 
675 
676  // Collect segments
677 
678  pointField start(testFaces.size());
679  pointField end(testFaces.size());
680  labelList minLevel(testFaces.size());
681 
682  calcCellCellRays
683  (
684  neiCc,
685  neiLevel,
686  testFaces,
687  start,
688  end,
689  minLevel
690  );
691 
692  // Do tests in one go
693 
694  surfaces_.findNearestIntersection
695  (
696  surfacesToTest,
697  start,
698  end,
699 
700  surface1,
701  hit1,
702  region1,
703  surface2,
704  hit2,
705  region2
706  );
707 }
708 
709 
710 Foam::labelList Foam::meshRefinement::nearestIntersection
711 (
712  const labelList& surfacesToTest,
713  const label defaultRegion
714 ) const
715 {
716  // Determine nearest intersection for all mesh faces. Used when removing
717  // cells to give some reasonable patch to exposed faces. Use this
718  // function instead of nearestPatch if you don't have patches yet.
719 
720 
721  // All faces with any intersection
722  const labelList testFaces(intersectedFaces());
723 
724  // Find intersection (if any)
725  labelList surface1;
726  List<pointIndexHit> hit1;
727  labelList region1;
728  labelList surface2;
729  List<pointIndexHit> hit2;
730  labelList region2;
731  nearestIntersection
732  (
733  surfacesToTest,
734  testFaces,
735 
736  surface1,
737  hit1,
738  region1,
739  surface2,
740  hit2,
741  region2
742  );
743 
744 
745  labelList nearestRegion(mesh_.nFaces(), defaultRegion);
746 
747  // Field on cells and faces.
748  List<topoDistanceData<label>> cellData(mesh_.nCells());
749  List<topoDistanceData<label>> faceData(mesh_.nFaces());
750 
751  // Start walking from all intersected faces
752  DynamicList<label> patchFaces(testFaces.size());
753  DynamicList<topoDistanceData<label>> patchData(testFaces.size());
754  forAll(testFaces, i)
755  {
756  label facei = testFaces[i];
757  if (surface1[i] != -1)
758  {
759  patchFaces.append(facei);
760  label regioni = surfaces_.globalRegion(surface1[i], region1[i]);
761  patchData.append(topoDistanceData<label>(0, regioni));
762  }
763  else if (surface2[i] != -1)
764  {
765  patchFaces.append(facei);
766  label regioni = surfaces_.globalRegion(surface2[i], region2[i]);
767  patchData.append(topoDistanceData<label>(0, regioni));
768  }
769  }
770 
771  // Propagate information inwards
772  FaceCellWave<topoDistanceData<label>> deltaCalc
773  (
774  mesh_,
775  patchFaces,
776  patchData,
777  faceData,
778  cellData,
779  mesh_.globalData().nTotalCells()+1
780  );
781 
782  // And extract
783 
784  bool haveWarned = false;
785  forAll(faceData, facei)
786  {
787  if (!faceData[facei].valid(deltaCalc.data()))
788  {
789  if (!haveWarned)
790  {
791  WarningInFunction
792  << "Did not visit some faces, e.g. face " << facei
793  << " at " << mesh_.faceCentres()[facei] << endl
794  << "Assigning these faces to global region "
795  << defaultRegion << endl;
796  haveWarned = true;
797  }
798  }
799  else
800  {
801  nearestRegion[facei] = faceData[facei].data();
802  }
803  }
804 
805  return nearestRegion;
806 }
807 
808 
809 void Foam::meshRefinement::testSyncPointList
810 (
811  const string& msg,
812  const polyMesh& mesh,
813  const List<scalar>& fld
814 )
815 {
816  if (fld.size() != mesh.nPoints())
817  {
818  FatalErrorInFunction
819  << msg << endl
820  << "fld size:" << fld.size() << " mesh points:" << mesh.nPoints()
821  << abort(FatalError);
822  }
823 
824  Pout<< "Checking field " << msg << endl;
825  scalarField minFld(fld);
826  syncTools::syncPointList
827  (
828  mesh,
829  minFld,
830  minEqOp<scalar>(),
831  GREAT
832  );
833  scalarField maxFld(fld);
834  syncTools::syncPointList
835  (
836  mesh,
837  maxFld,
838  maxEqOp<scalar>(),
839  -GREAT
840  );
841  forAll(minFld, pointi)
842  {
843  const scalar& minVal = minFld[pointi];
844  const scalar& maxVal = maxFld[pointi];
845  if (mag(minVal-maxVal) > SMALL)
846  {
847  Pout<< msg << " at:" << mesh.points()[pointi] << nl
848  << " minFld:" << minVal << nl
849  << " maxFld:" << maxVal << nl
850  << endl;
851  }
852  }
853 }
854 
855 
856 void Foam::meshRefinement::testSyncPointList
857 (
858  const string& msg,
859  const polyMesh& mesh,
860  const List<point>& fld
861 )
862 {
863  if (fld.size() != mesh.nPoints())
864  {
865  FatalErrorInFunction
866  << msg << endl
867  << "fld size:" << fld.size() << " mesh points:" << mesh.nPoints()
868  << abort(FatalError);
869  }
870 
871  Pout<< "Checking field " << msg << endl;
872  pointField minFld(fld);
873  syncTools::syncPointList
874  (
875  mesh,
876  minFld,
877  minMagSqrEqOp<point>(),
878  point(GREAT, GREAT, GREAT)
879  );
880  pointField maxFld(fld);
881  syncTools::syncPointList
882  (
883  mesh,
884  maxFld,
885  maxMagSqrEqOp<point>(),
886  vector::zero
887  );
888  forAll(minFld, pointi)
889  {
890  const point& minVal = minFld[pointi];
891  const point& maxVal = maxFld[pointi];
892  if (mag(minVal-maxVal) > SMALL)
893  {
894  Pout<< msg << " at:" << mesh.points()[pointi] << nl
895  << " minFld:" << minVal << nl
896  << " maxFld:" << maxVal << nl
897  << endl;
898  }
899  }
900 }
901 
902 
903 void Foam::meshRefinement::checkData()
904 {
905  Pout<< "meshRefinement::checkData() : Checking refinement structure."
906  << endl;
907  meshCutter_.checkMesh();
908 
909  Pout<< "meshRefinement::checkData() : Checking refinement levels."
910  << endl;
911  meshCutter_.checkRefinementLevels(1, labelList(0));
912 
913 
914  const label nBnd = mesh_.nBoundaryFaces();
915 
916  Pout<< "meshRefinement::checkData() : Checking synchronization."
917  << endl;
918 
919  // Check face centres
920  {
921  // Boundary face centres
922  pointField::subList boundaryFc
923  (
924  mesh_.faceCentres(),
925  nBnd,
926  mesh_.nInternalFaces()
927  );
928 
929  // Get neighbouring face centres
930  pointField neiBoundaryFc(boundaryFc);
931  syncTools::syncBoundaryFacePositions
932  (
933  mesh_,
934  neiBoundaryFc,
935  eqOp<point>()
936  );
937 
938  // Compare
939  testSyncBoundaryFaceList
940  (
941  mergeDistance_,
942  "testing faceCentres : ",
943  boundaryFc,
944  neiBoundaryFc
945  );
946  }
947 
948  // Check meshRefinement
949  const labelList& surfIndex = surfaceIndex();
950 
951 
952  {
953  // Get boundary face centre and level. Coupled aware.
954  labelList neiLevel(nBnd);
955  pointField neiCc(nBnd);
956  calcNeighbourData(neiLevel, neiCc);
957 
958  // Collect segments we want to test for
959  pointField start(mesh_.nFaces());
960  pointField end(mesh_.nFaces());
961 
962  forAll(start, facei)
963  {
964  start[facei] = mesh_.cellCentres()[mesh_.faceOwner()[facei]];
965 
966  if (mesh_.isInternalFace(facei))
967  {
968  end[facei] = mesh_.cellCentres()[mesh_.faceNeighbour()[facei]];
969  }
970  else
971  {
972  end[facei] = neiCc[facei-mesh_.nInternalFaces()];
973  }
974  }
975 
976  // Extend segments a bit
977  {
978  const vectorField smallVec(ROOTSMALL*(end-start));
979  start -= smallVec;
980  end += smallVec;
981  }
982 
983 
984  // Do tests in one go
985  labelList surfaceHit;
986  {
987  labelList surfaceLevel;
988  surfaces_.findHigherIntersection
989  (
990  shells_,
991  start,
992  end,
993  labelList(start.size(), -1), // accept any intersection
994  surfaceHit,
995  surfaceLevel
996  );
997  }
998  // Get the coupled hit
999  labelList neiHit
1000  (
1001  SubList<label>
1002  (
1003  surfaceHit,
1004  nBnd,
1005  mesh_.nInternalFaces()
1006  )
1007  );
1008  syncTools::swapBoundaryFaceList(mesh_, neiHit);
1009 
1010  // Check
1011  forAll(surfaceHit, facei)
1012  {
1013  if (surfIndex[facei] != surfaceHit[facei])
1014  {
1015  if (mesh_.isInternalFace(facei))
1016  {
1017  WarningInFunction
1018  << "Internal face:" << facei
1019  << " fc:" << mesh_.faceCentres()[facei]
1020  << " cached surfaceIndex_:" << surfIndex[facei]
1021  << " current:" << surfaceHit[facei]
1022  << " ownCc:"
1023  << mesh_.cellCentres()[mesh_.faceOwner()[facei]]
1024  << " neiCc:"
1025  << mesh_.cellCentres()[mesh_.faceNeighbour()[facei]]
1026  << endl;
1027  }
1028  else if
1029  (
1030  surfIndex[facei]
1031  != neiHit[facei-mesh_.nInternalFaces()]
1032  )
1033  {
1034  WarningInFunction
1035  << "Boundary face:" << facei
1036  << " fc:" << mesh_.faceCentres()[facei]
1037  << " cached surfaceIndex_:" << surfIndex[facei]
1038  << " current:" << surfaceHit[facei]
1039  << " ownCc:"
1040  << mesh_.cellCentres()[mesh_.faceOwner()[facei]]
1041  << " neiCc:"
1042  << neiCc[facei-mesh_.nInternalFaces()]
1043  << " end:" << end[facei]
1044  << " ownLevel:"
1045  << meshCutter_.cellLevel()[mesh_.faceOwner()[facei]]
1046  << " faceLevel:"
1047  << meshCutter_.faceLevel(facei)
1048  << endl;
1049  }
1050  }
1051  }
1052  }
1053  {
1054  labelList::subList boundarySurface
1055  (
1056  surfaceIndex_,
1057  mesh_.nBoundaryFaces(),
1058  mesh_.nInternalFaces()
1059  );
1060 
1061  labelList neiBoundarySurface(boundarySurface);
1062  syncTools::swapBoundaryFaceList
1063  (
1064  mesh_,
1065  neiBoundarySurface
1066  );
1067 
1068  // Compare
1069  testSyncBoundaryFaceList
1070  (
1071  0, // tolerance
1072  "testing surfaceIndex() : ",
1073  boundarySurface,
1074  neiBoundarySurface
1075  );
1076  }
1077 
1078 
1079  // Find duplicate faces
1080  Pout<< "meshRefinement::checkData() : Counting duplicate faces."
1081  << endl;
1082 
1083  labelList duplicateFace
1084  (
1085  localPointRegion::findDuplicateFaces
1086  (
1087  mesh_,
1088  identity(mesh_.nBoundaryFaces(), mesh_.nInternalFaces())
1089  )
1090  );
1091 
1092  // Count
1093  {
1094  label nDup = 0;
1095 
1096  forAll(duplicateFace, i)
1097  {
1098  if (duplicateFace[i] != -1)
1099  {
1100  nDup++;
1101  }
1102  }
1103  nDup /= 2; // will have counted both faces of duplicate
1104  Pout<< "meshRefinement::checkData() : Found " << nDup
1105  << " duplicate pairs of faces." << endl;
1106  }
1107 }
1108 
1109 
1110 void Foam::meshRefinement::setInstance(const fileName& inst)
1111 {
1112  meshCutter_.setInstance(inst);
1113  surfaceIndex_.instance() = inst;
1114 }
1115 
1116 
1117 Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::doRemoveCells
1118 (
1119  const labelList& cellsToRemove,
1120  const labelList& exposedFaces,
1121  const labelList& exposedPatchIDs,
1122  removeCells& cellRemover
1123 )
1124 {
1125  polyTopoChange meshMod(mesh_);
1126 
1127  // Arbitrary: put exposed faces into last patch.
1128  cellRemover.setRefinement
1129  (
1130  cellsToRemove,
1131  exposedFaces,
1132  exposedPatchIDs,
1133  meshMod
1134  );
1135 
1136  // Remove any unnecessary fields
1137  mesh_.clearOut();
1138  mesh_.moving(false);
1139 
1140  // Change the mesh (no inflation)
1141  autoPtr<mapPolyMesh> mapPtr = meshMod.changeMesh(mesh_, false, true);
1142  mapPolyMesh& map = *mapPtr;
1143 
1144  // Update fields
1145  mesh_.updateMesh(map);
1146 
1147  // Move mesh (since morphing might not do this)
1148  if (map.hasMotionPoints())
1149  {
1150  mesh_.movePoints(map.preMotionPoints());
1151  }
1152  else
1153  {
1154  // Delete mesh volumes. No other way to do this?
1155  mesh_.clearOut();
1156  }
1157 
1158  // Reset the instance for if in overwrite mode
1159  mesh_.setInstance(timeName());
1160  setInstance(mesh_.facesInstance());
1161 
1162  // Update local mesh data
1163  cellRemover.updateMesh(map);
1164 
1165  // Update intersections. Recalculate intersections for exposed faces.
1166  labelList newExposedFaces = renumber
1167  (
1168  map.reverseFaceMap(),
1169  exposedFaces
1170  );
1171 
1172  //Pout<< "removeCells : updating intersections for "
1173  // << newExposedFaces.size() << " newly exposed faces." << endl;
1174 
1175  updateMesh(map, newExposedFaces);
1176 
1177  return mapPtr;
1178 }
1179 
1180 
1181 void Foam::meshRefinement::doSplitFaces
1182 (
1183  const labelList& splitFaces,
1184  const labelPairList& splits,
1185  //const List<Pair<point>>& splitPoints,
1186  polyTopoChange& meshMod
1187 ) const
1188 {
1189  forAll(splitFaces, i)
1190  {
1191  label facei = splitFaces[i];
1192  const face& f = mesh_.faces()[facei];
1193 
1194  // Split as start and end index in face
1195  const labelPair& split = splits[i];
1196 
1197  label nVerts = split[1]-split[0];
1198  if (nVerts < 0)
1199  {
1200  nVerts += f.size();
1201  }
1202  nVerts += 1;
1203 
1204 
1205  // Split into f0, f1
1206  face f0(nVerts);
1207 
1208  label fp = split[0];
1209  forAll(f0, i)
1210  {
1211  f0[i] = f[fp];
1212  fp = f.fcIndex(fp);
1213  }
1214 
1215  face f1(f.size()-f0.size()+2);
1216  fp = split[1];
1217  forAll(f1, i)
1218  {
1219  f1[i] = f[fp];
1220  fp = f.fcIndex(fp);
1221  }
1222 
1223 
1224  // Determine face properties
1225  label own = mesh_.faceOwner()[facei];
1226  label nei = -1;
1227  label patchi = -1;
1228  if (facei >= mesh_.nInternalFaces())
1229  {
1230  patchi = mesh_.boundaryMesh().whichPatch(facei);
1231  }
1232  else
1233  {
1234  nei = mesh_.faceNeighbour()[facei];
1235  }
1236 
1237  label zonei = mesh_.faceZones().whichZone(facei);
1238  bool zoneFlip = false;
1239  if (zonei != -1)
1240  {
1241  const faceZone& fz = mesh_.faceZones()[zonei];
1242  zoneFlip = fz.flipMap()[fz.whichFace(facei)];
1243  }
1244 
1245 
1246  if (debug)
1247  {
1248  Pout<< "face:" << facei << " verts:" << f
1249  << " split into f0:" << f0
1250  << " f1:" << f1 << endl;
1251  }
1252 
1253  // Change/add faces
1254  meshMod.modifyFace
1255  (
1256  f0, // modified face
1257  facei, // label of face
1258  own, // owner
1259  nei, // neighbour
1260  false, // face flip
1261  patchi, // patch for face
1262  zonei, // zone for face
1263  zoneFlip // face flip in zone
1264  );
1265 
1266  meshMod.addFace
1267  (
1268  f1, // modified face
1269  own, // owner
1270  nei, // neighbour
1271  -1, // master point
1272  -1, // master edge
1273  facei, // master face
1274  false, // face flip
1275  patchi, // patch for face
1276  zonei, // zone for face
1277  zoneFlip // face flip in zone
1278  );
1279 
1280 
1282  //meshMod.modifyPoint
1283  //(
1284  // f[split[0]],
1285  // splitPoints[i][0],
1286  // -1,
1287  // true
1288  //);
1289  //meshMod.modifyPoint
1290  //(
1291  // f[split[1]],
1292  // splitPoints[i][1],
1293  // -1,
1294  // true
1295  //);
1296  }
1297 }
1298 
1299 
1300 Foam::label Foam::meshRefinement::splitFacesUndo
1301 (
1302  const labelList& splitFaces,
1303  const labelPairList& splits,
1304  const dictionary& motionDict,
1305 
1306  labelList& duplicateFace,
1307  List<labelPair>& baffles
1308 )
1309 {
1310  label nSplit = returnReduce(splitFaces.size(), sumOp<label>());
1311  Info<< nl
1312  << "Splitting " << nSplit
1313  << " faces across diagonals" << "." << nl << endl;
1314 
1315  // To undo: store original faces
1316  faceList originalFaces(UIndirectList<face>(mesh_.faces(), splitFaces));
1317  labelPairList facePairs(splitFaces.size(), labelPair(-1, -1));
1318 
1319 
1320  {
1321  polyTopoChange meshMod(mesh_);
1322  meshMod.setCapacity
1323  (
1324  meshMod.points().size(),
1325  meshMod.faces().size()+splitFaces.size(),
1326  mesh_.nCells()
1327  );
1328 
1329  // Insert the mesh changes
1330  doSplitFaces(splitFaces, splits, meshMod);
1331 
1332  // Remove any unnecessary fields
1333  mesh_.clearOut();
1334  mesh_.moving(false);
1335 
1336  // Change the mesh (no inflation)
1337  autoPtr<mapPolyMesh> mapPtr = meshMod.changeMesh(mesh_, false, true);
1338  mapPolyMesh& map = *mapPtr;
1339 
1340  // Update fields
1341  mesh_.updateMesh(map);
1342 
1343  // Move mesh (since morphing might not do this)
1344  if (map.hasMotionPoints())
1345  {
1346  mesh_.movePoints(map.preMotionPoints());
1347  }
1348  else
1349  {
1350  mesh_.clearOut();
1351  }
1352 
1353  // Reset the instance for if in overwrite mode
1354  mesh_.setInstance(timeName());
1355  setInstance(mesh_.facesInstance());
1356 
1357 
1358  // Update local mesh data
1359  // ~~~~~~~~~~~~~~~~~~~~~~
1360 
1361  forAll(originalFaces, i)
1362  {
1363  inplaceRenumber(map.reversePointMap(), originalFaces[i]);
1364  }
1365 
1366  {
1367  Map<label> splitFaceToIndex(2*splitFaces.size());
1368  forAll(splitFaces, i)
1369  {
1370  splitFaceToIndex.insert(splitFaces[i], i);
1371  }
1372 
1373  forAll(map.faceMap(), facei)
1374  {
1375  label oldFacei = map.faceMap()[facei];
1376 
1377  const auto oldFaceFnd = splitFaceToIndex.cfind(oldFacei);
1378 
1379  if (oldFaceFnd.found())
1380  {
1381  labelPair& twoFaces = facePairs[oldFaceFnd.val()];
1382  if (twoFaces[0] == -1)
1383  {
1384  twoFaces[0] = facei;
1385  }
1386  else if (twoFaces[1] == -1)
1387  {
1388  twoFaces[1] = facei;
1389  }
1390  else
1391  {
1392  FatalErrorInFunction
1393  << "problem: twoFaces:" << twoFaces
1394  << exit(FatalError);
1395  }
1396  }
1397  }
1398  }
1399 
1400 
1401  // Update baffle data
1402  // ~~~~~~~~~~~~~~~~~~
1403 
1404  if (duplicateFace.size())
1405  {
1406  meshRefinement::updateList
1407  (
1408  map.faceMap(),
1409  label(-1),
1410  duplicateFace
1411  );
1412  }
1413 
1414  const labelList& oldToNewFaces = map.reverseFaceMap();
1415  forAll(baffles, i)
1416  {
1417  labelPair& baffle = baffles[i];
1418  baffle.first() = oldToNewFaces[baffle.first()];
1419  baffle.second() = oldToNewFaces[baffle.second()];
1420 
1421  if (baffle.first() == -1 || baffle.second() == -1)
1422  {
1423  FatalErrorInFunction
1424  << "Removed baffle : faces:" << baffle
1425  << exit(FatalError);
1426  }
1427  }
1428 
1429 
1430  // Update intersections
1431  // ~~~~~~~~~~~~~~~~~~~~
1432 
1433  {
1434  DynamicList<label> changedFaces(facePairs.size());
1435  forAll(facePairs, i)
1436  {
1437  changedFaces.append(facePairs[i].first());
1438  changedFaces.append(facePairs[i].second());
1439  }
1440 
1441  // Update intersections on changed faces
1442  updateMesh(map, growFaceCellFace(changedFaces));
1443  }
1444  }
1445 
1446 
1447 
1448  // Undo loop
1449  // ~~~~~~~~~
1450  // Maintains two bits of information:
1451  // facePairs : two faces originating from the same face
1452  // originalFaces : original face in current vertices
1453 
1454 
1455  for (label iteration = 0; iteration < 100; iteration++)
1456  {
1457  Info<< nl
1458  << "Undo iteration " << iteration << nl
1459  << "----------------" << endl;
1460 
1461 
1462  // Check mesh for errors
1463  // ~~~~~~~~~~~~~~~~~~~~~
1464 
1465  faceSet errorFaces(mesh_, "errorFaces", mesh_.nBoundaryFaces());
1466  bool hasErrors = motionSmoother::checkMesh
1467  (
1468  false, // report
1469  mesh_,
1470  motionDict,
1471  errorFaces,
1472  dryRun_
1473  );
1474  if (!hasErrors)
1475  {
1476  break;
1477  }
1478 
1479  // Extend faces
1480  {
1481  const labelList grownFaces(growFaceCellFace(errorFaces));
1482  errorFaces.clear();
1483  errorFaces.insert(grownFaces);
1484  }
1485 
1486 
1487  // Merge face pairs
1488  // ~~~~~~~~~~~~~~~~
1489  // (if one of the faces is in the errorFaces set)
1490 
1491  polyTopoChange meshMod(mesh_);
1492 
1493  // Indices (in facePairs) of merged faces
1494  labelHashSet mergedIndices(facePairs.size());
1495  forAll(facePairs, index)
1496  {
1497  const labelPair& twoFaces = facePairs[index];
1498 
1499  if
1500  (
1501  errorFaces.found(twoFaces.first())
1502  || errorFaces.found(twoFaces.second())
1503  )
1504  {
1505  const face& originalFace = originalFaces[index];
1506 
1507 
1508  // Determine face properties
1509  label own = mesh_.faceOwner()[twoFaces[0]];
1510  label nei = -1;
1511  label patchi = -1;
1512  if (twoFaces[0] >= mesh_.nInternalFaces())
1513  {
1514  patchi = mesh_.boundaryMesh().whichPatch(twoFaces[0]);
1515  }
1516  else
1517  {
1518  nei = mesh_.faceNeighbour()[twoFaces[0]];
1519  }
1520 
1521  label zonei = mesh_.faceZones().whichZone(twoFaces[0]);
1522  bool zoneFlip = false;
1523  if (zonei != -1)
1524  {
1525  const faceZone& fz = mesh_.faceZones()[zonei];
1526  zoneFlip = fz.flipMap()[fz.whichFace(twoFaces[0])];
1527  }
1528 
1529  // Modify first face
1530  meshMod.modifyFace
1531  (
1532  originalFace, // modified face
1533  twoFaces[0], // label of face
1534  own, // owner
1535  nei, // neighbour
1536  false, // face flip
1537  patchi, // patch for face
1538  zonei, // zone for face
1539  zoneFlip // face flip in zone
1540  );
1541  // Merge second face into first
1542  meshMod.removeFace(twoFaces[1], twoFaces[0]);
1543 
1544  mergedIndices.insert(index);
1545  }
1546  }
1547 
1548  label n = returnReduce(mergedIndices.size(), sumOp<label>());
1549 
1550  Info<< "Detected " << n
1551  << " split faces that will be restored to their original faces."
1552  << nl << endl;
1553 
1554  if (n == 0)
1555  {
1556  // Nothing to be restored
1557  break;
1558  }
1559 
1560  nSplit -= n;
1561 
1562 
1563  // Remove any unnecessary fields
1564  mesh_.clearOut();
1565  mesh_.moving(false);
1566 
1567  // Change the mesh (no inflation)
1568  autoPtr<mapPolyMesh> mapPtr = meshMod.changeMesh(mesh_, false, true);
1569  mapPolyMesh& map = *mapPtr;
1570 
1571  // Update fields
1572  mesh_.updateMesh(map);
1573 
1574  // Move mesh (since morphing might not do this)
1575  if (map.hasMotionPoints())
1576  {
1577  mesh_.movePoints(map.preMotionPoints());
1578  }
1579  else
1580  {
1581  mesh_.clearOut();
1582  }
1583 
1584  // Reset the instance for if in overwrite mode
1585  mesh_.setInstance(timeName());
1586  setInstance(mesh_.facesInstance());
1587 
1588 
1589  // Update local mesh data
1590  // ~~~~~~~~~~~~~~~~~~~~~~
1591 
1592  {
1593  const labelList& oldToNewFaces = map.reverseFaceMap();
1594  const labelList& oldToNewPoints = map.reversePointMap();
1595 
1596  // Compact out merged faces
1597  DynamicList<label> changedFaces(mergedIndices.size());
1598 
1599  label newIndex = 0;
1600  forAll(facePairs, index)
1601  {
1602  const labelPair& oldSplit = facePairs[index];
1603  label new0 = oldToNewFaces[oldSplit[0]];
1604  label new1 = oldToNewFaces[oldSplit[1]];
1605 
1606  if (!mergedIndices.found(index))
1607  {
1608  // Faces still split
1609  if (new0 < 0 || new1 < 0)
1610  {
1611  FatalErrorInFunction
1612  << "Problem: oldFaces:" << oldSplit
1613  << " newFaces:" << labelPair(new0, new1)
1614  << exit(FatalError);
1615  }
1616 
1617  facePairs[newIndex] = labelPair(new0, new1);
1618  originalFaces[newIndex] = renumber
1619  (
1620  oldToNewPoints,
1621  originalFaces[index]
1622  );
1623  newIndex++;
1624  }
1625  else
1626  {
1627  // Merged face. Only new0 kept.
1628  if (new0 < 0 || new1 == -1)
1629  {
1630  FatalErrorInFunction
1631  << "Problem: oldFaces:" << oldSplit
1632  << " newFace:" << labelPair(new0, new1)
1633  << exit(FatalError);
1634  }
1635  changedFaces.append(new0);
1636  }
1637  }
1638 
1639  facePairs.setSize(newIndex);
1640  originalFaces.setSize(newIndex);
1641 
1642 
1643  // Update intersections
1644  updateMesh(map, growFaceCellFace(changedFaces));
1645  }
1646 
1647  // Update baffle data
1648  // ~~~~~~~~~~~~~~~~~~
1649  {
1650  if (duplicateFace.size())
1651  {
1652  meshRefinement::updateList
1653  (
1654  map.faceMap(),
1655  label(-1),
1656  duplicateFace
1657  );
1658  }
1659 
1660  const labelList& reverseFaceMap = map.reverseFaceMap();
1661  forAll(baffles, i)
1662  {
1663  labelPair& baffle = baffles[i];
1664  baffle.first() = reverseFaceMap[baffle.first()];
1665  baffle.second() = reverseFaceMap[baffle.second()];
1666 
1667  if (baffle.first() == -1 || baffle.second() == -1)
1668  {
1669  FatalErrorInFunction
1670  << "Removed baffle : faces:" << baffle
1671  << exit(FatalError);
1672  }
1673  }
1674  }
1675 
1676  }
1677 
1678  return nSplit;
1679 }
1680 
1681 
1682 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
1683 
1684 Foam::meshRefinement::meshRefinement
1685 (
1686  fvMesh& mesh,
1687  const scalar mergeDistance,
1688  const bool overwrite,
1689  const refinementSurfaces& surfaces,
1690  const refinementFeatures& features,
1691  const shellSurfaces& shells,
1692  const shellSurfaces& limitShells,
1693  const labelUList& checkFaces,
1694  const bool dryRun
1695 )
1696 :
1697  mesh_(mesh),
1698  mergeDistance_(mergeDistance),
1699  overwrite_(overwrite),
1700  oldInstance_(mesh.pointsInstance()),
1701  surfaces_(surfaces),
1702  features_(features),
1703  shells_(shells),
1704  limitShells_(limitShells),
1705  dryRun_(dryRun),
1706  meshCutter_
1707  (
1708  mesh,
1709  false // do not try to read history.
1710  ),
1711  surfaceIndex_
1712  (
1713  IOobject
1714  (
1715  "surfaceIndex",
1716  mesh_.facesInstance(),
1717  fvMesh::meshSubDir,
1718  mesh_,
1719  IOobject::NO_READ,
1720  IOobject::NO_WRITE,
1721  false
1722  ),
1723  labelList(mesh_.nFaces(), -1)
1724  ),
1725  userFaceData_(0)
1726 {
1727  // recalculate intersections for all faces
1728  updateIntersections(checkFaces);
1729 }
1730 
1731 
1732 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1733 
1734 const Foam::labelList& Foam::meshRefinement::surfaceIndex() const
1735 {
1736  if (surfaceIndex_.size() != mesh_.nFaces())
1737  {
1738  const_cast<meshRefinement&>(*this).updateIntersections
1739  (
1740  identity(mesh_.nFaces())
1741  );
1742  }
1743  return surfaceIndex_;
1744 }
1745 
1746 
1747 Foam::labelList& Foam::meshRefinement::surfaceIndex()
1748 {
1749  if (surfaceIndex_.size() != mesh_.nFaces())
1750  {
1751  updateIntersections(identity(mesh_.nFaces()));
1752  }
1753  return surfaceIndex_;
1754 }
1755 
1756 
1757 Foam::label Foam::meshRefinement::countHits() const
1758 {
1759  // Stats on edges to test. Count proc faces only once.
1760  bitSet isMasterFace(syncTools::getMasterFaces(mesh_));
1761 
1762  label nHits = 0;
1763 
1764  const labelList& surfIndex = surfaceIndex();
1765 
1766  forAll(surfIndex, facei)
1767  {
1768  if (surfIndex[facei] >= 0 && isMasterFace.test(facei))
1769  {
1770  ++nHits;
1771  }
1772  }
1773  return nHits;
1774 }
1775 
1776 
1777 Foam::autoPtr<Foam::mapDistributePolyMesh> Foam::meshRefinement::balance
1778 (
1779  const bool keepZoneFaces,
1780  const bool keepBaffles,
1781  const scalarField& cellWeights,
1782  decompositionMethod& decomposer,
1783  fvMeshDistribute& distributor
1784 )
1785 {
1786  autoPtr<mapDistributePolyMesh> map;
1787 
1788  if (Pstream::parRun())
1789  {
1790  // Wanted distribution
1791  labelList distribution;
1792 
1793 
1794  // Faces where owner and neighbour are not 'connected' so can
1795  // go to different processors.
1796  boolList blockedFace;
1797  label nUnblocked = 0;
1798 
1799  // Faces that move as block onto single processor
1800  PtrList<labelList> specifiedProcessorFaces;
1801  labelList specifiedProcessor;
1802 
1803  // Pairs of baffles
1804  List<labelPair> couples;
1805 
1806  // Constraints from decomposeParDict
1807  decomposer.setConstraints
1808  (
1809  mesh_,
1810  blockedFace,
1811  specifiedProcessorFaces,
1812  specifiedProcessor,
1813  couples
1814  );
1815 
1816 
1817  if (keepZoneFaces || keepBaffles)
1818  {
1819  if (keepZoneFaces)
1820  {
1821  // Determine decomposition to keep/move surface zones
1822  // on one processor. The reason is that snapping will make these
1823  // into baffles, move and convert them back so if they were
1824  // proc boundaries after baffling&moving the points might be no
1825  // longer synchronised so recoupling will fail. To prevent this
1826  // keep owner&neighbour of such a surface zone on the same
1827  // processor.
1828 
1829  const faceZoneMesh& fZones = mesh_.faceZones();
1830  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
1831 
1832  // Get faces whose owner and neighbour should stay together,
1833  // i.e. they are not 'blocked'.
1834 
1835  forAll(fZones, zonei)
1836  {
1837  const faceZone& fZone = fZones[zonei];
1838 
1839  forAll(fZone, i)
1840  {
1841  label facei = fZone[i];
1842  if (blockedFace[facei])
1843  {
1844  if
1845  (
1846  mesh_.isInternalFace(facei)
1847  || pbm[pbm.whichPatch(facei)].coupled()
1848  )
1849  {
1850  blockedFace[facei] = false;
1851  nUnblocked++;
1852  }
1853  }
1854  }
1855  }
1856 
1857 
1858  // If the faceZones are not synchronised the blockedFace
1859  // might not be synchronised. If you are sure the faceZones
1860  // are synchronised remove below check.
1861  syncTools::syncFaceList
1862  (
1863  mesh_,
1864  blockedFace,
1865  andEqOp<bool>() // combine operator
1866  );
1867  }
1868  reduce(nUnblocked, sumOp<label>());
1869 
1870  if (keepZoneFaces)
1871  {
1872  Info<< "Found " << nUnblocked
1873  << " zoned faces to keep together." << endl;
1874  }
1875 
1876 
1877  // Extend un-blockedFaces with any cyclics
1878  {
1879  boolList separatedCoupledFace(mesh_.nFaces(), false);
1880  selectSeparatedCoupledFaces(separatedCoupledFace);
1881 
1882  label nSeparated = 0;
1883  forAll(separatedCoupledFace, facei)
1884  {
1885  if (separatedCoupledFace[facei])
1886  {
1887  if (blockedFace[facei])
1888  {
1889  blockedFace[facei] = false;
1890  nSeparated++;
1891  }
1892  }
1893  }
1894  reduce(nSeparated, sumOp<label>());
1895  Info<< "Found " << nSeparated
1896  << " separated coupled faces to keep together." << endl;
1897 
1898  nUnblocked += nSeparated;
1899  }
1900 
1901 
1902  if (keepBaffles)
1903  {
1904  const label nBnd = mesh_.nBoundaryFaces();
1905 
1906  labelList coupledFace(mesh_.nFaces(), -1);
1907  {
1908  // Get boundary baffles that need to stay together
1909  List<labelPair> allCouples
1910  (
1911  localPointRegion::findDuplicateFacePairs(mesh_)
1912  );
1913 
1914  // Merge with any couples from
1915  // decompositionMethod::setConstraints
1916  forAll(couples, i)
1917  {
1918  const labelPair& baffle = couples[i];
1919  coupledFace[baffle.first()] = baffle.second();
1920  coupledFace[baffle.second()] = baffle.first();
1921  }
1922  forAll(allCouples, i)
1923  {
1924  const labelPair& baffle = allCouples[i];
1925  coupledFace[baffle.first()] = baffle.second();
1926  coupledFace[baffle.second()] = baffle.first();
1927  }
1928  }
1929 
1930  couples.setSize(nBnd);
1931  label nCpl = 0;
1932  forAll(coupledFace, facei)
1933  {
1934  if (coupledFace[facei] != -1 && facei < coupledFace[facei])
1935  {
1936  couples[nCpl++] = labelPair(facei, coupledFace[facei]);
1937  }
1938  }
1939  couples.setSize(nCpl);
1940  }
1941  label nCouples = returnReduce(couples.size(), sumOp<label>());
1942 
1943  if (keepBaffles)
1944  {
1945  Info<< "Found " << nCouples << " baffles to keep together."
1946  << endl;
1947  }
1948  }
1949 
1950 
1951  // Make sure blockedFace not set on couples
1952  forAll(couples, i)
1953  {
1954  const labelPair& baffle = couples[i];
1955  blockedFace[baffle.first()] = false;
1956  blockedFace[baffle.second()] = false;
1957  }
1958 
1959  distribution = decomposer.decompose
1960  (
1961  mesh_,
1962  cellWeights,
1963  blockedFace,
1964  specifiedProcessorFaces,
1965  specifiedProcessor,
1966  couples // explicit connections
1967  );
1968 
1969  if (debug)
1970  {
1971  labelList nProcCells(distributor.countCells(distribution));
1972  Pout<< "Wanted distribution:" << nProcCells << endl;
1973 
1974  Pstream::listCombineReduce(nProcCells, plusEqOp<label>());
1975 
1976  Pout<< "Wanted resulting decomposition:" << endl;
1977  forAll(nProcCells, proci)
1978  {
1979  Pout<< " " << proci << '\t' << nProcCells[proci] << endl;
1980  }
1981  Pout<< endl;
1982  }
1983 
1984  // Do actual sending/receiving of mesh
1985  map = distributor.distribute(distribution);
1986 
1987  // Update numbering of meshRefiner
1988  distribute(map());
1989 
1990  // Set correct instance (for if overwrite)
1991  mesh_.setInstance(timeName());
1992  setInstance(mesh_.facesInstance());
1993 
1994 
1995  if (debug && keepZoneFaces)
1996  {
1997  const faceZoneMesh& fZones = mesh_.faceZones();
1998  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
1999 
2000  // Check that faceZone faces are indeed internal
2001  forAll(fZones, zonei)
2002  {
2003  const faceZone& fZone = fZones[zonei];
2004 
2005  forAll(fZone, i)
2006  {
2007  label facei = fZone[i];
2008  label patchi = pbm.whichPatch(facei);
2009 
2010  if (patchi >= 0 && pbm[patchi].coupled())
2011  {
2012  WarningInFunction
2013  << "Face at " << mesh_.faceCentres()[facei]
2014  << " on zone " << fZone.name()
2015  << " is on coupled patch " << pbm[patchi].name()
2016  << endl;
2017  }
2018  }
2019  }
2020  }
2021  }
2022 
2023  return map;
2024 }
2025 
2026 
2027 Foam::labelList Foam::meshRefinement::intersectedFaces() const
2028 {
2029  label nBoundaryFaces = 0;
2030 
2031  const labelList& surfIndex = surfaceIndex();
2032 
2033  forAll(surfIndex, facei)
2034  {
2035  if (surfIndex[facei] != -1)
2036  {
2037  nBoundaryFaces++;
2038  }
2039  }
2040 
2041  labelList surfaceFaces(nBoundaryFaces);
2042  nBoundaryFaces = 0;
2043 
2044  forAll(surfIndex, facei)
2045  {
2046  if (surfIndex[facei] != -1)
2047  {
2048  surfaceFaces[nBoundaryFaces++] = facei;
2049  }
2050  }
2051  return surfaceFaces;
2052 }
2053 
2054 
2055 Foam::labelList Foam::meshRefinement::intersectedPoints() const
2056 {
2057  const faceList& faces = mesh_.faces();
2058 
2059  // Mark all points on faces that will become baffles
2060  bitSet isBoundaryPoint(mesh_.nPoints());
2061  label nBoundaryPoints = 0;
2062 
2063  const labelList& surfIndex = surfaceIndex();
2064 
2065  forAll(surfIndex, facei)
2066  {
2067  if (surfIndex[facei] != -1)
2068  {
2069  const face& f = faces[facei];
2070 
2071  forAll(f, fp)
2072  {
2073  if (isBoundaryPoint.set(f[fp]))
2074  {
2075  nBoundaryPoints++;
2076  }
2077  }
2078  }
2079  }
2080 
2082  //labelList adaptPatchIDs(meshedPatches());
2083  //forAll(adaptPatchIDs, i)
2084  //{
2085  // label patchi = adaptPatchIDs[i];
2086  //
2087  // if (patchi != -1)
2088  // {
2089  // const polyPatch& pp = mesh_.boundaryMesh()[patchi];
2090  //
2091  // label facei = pp.start();
2092  //
2093  // forAll(pp, i)
2094  // {
2095  // const face& f = faces[facei];
2096  //
2097  // forAll(f, fp)
2098  // {
2099  // if (isBoundaryPoint.set(f[fp]))
2100  // nBoundaryPoints++;
2101  // }
2102  // }
2103  // facei++;
2104  // }
2105  // }
2106  //}
2107 
2108 
2109  // Pack
2110  labelList boundaryPoints(nBoundaryPoints);
2111  nBoundaryPoints = 0;
2112  forAll(isBoundaryPoint, pointi)
2113  {
2114  if (isBoundaryPoint.test(pointi))
2115  {
2116  boundaryPoints[nBoundaryPoints++] = pointi;
2117  }
2118  }
2119 
2120  return boundaryPoints;
2121 }
2122 
2123 
2124 //- Create patch from set of patches
2125 Foam::autoPtr<Foam::indirectPrimitivePatch> Foam::meshRefinement::makePatch
2126 (
2127  const polyMesh& mesh,
2128  const labelList& patchIDs
2129 )
2130 {
2131  const polyBoundaryMesh& patches = mesh.boundaryMesh();
2132 
2133  // Count faces.
2134  label nFaces = 0;
2135 
2136  forAll(patchIDs, i)
2137  {
2138  const polyPatch& pp = patches[patchIDs[i]];
2139 
2140  nFaces += pp.size();
2141  }
2142 
2143  // Collect faces.
2144  labelList addressing(nFaces);
2145  nFaces = 0;
2146 
2147  forAll(patchIDs, i)
2148  {
2149  const polyPatch& pp = patches[patchIDs[i]];
2150 
2151  label meshFacei = pp.start();
2152 
2153  forAll(pp, i)
2154  {
2155  addressing[nFaces++] = meshFacei++;
2156  }
2157  }
2158 
2159  return autoPtr<indirectPrimitivePatch>::New
2160  (
2161  IndirectList<face>(mesh.faces(), addressing),
2162  mesh.points()
2163  );
2164 }
2165 
2166 
2167 Foam::tmp<Foam::pointVectorField> Foam::meshRefinement::makeDisplacementField
2168 (
2169  const pointMesh& pMesh,
2170  const labelList& adaptPatchIDs
2171 )
2172 {
2173  const polyMesh& mesh = pMesh();
2174 
2175  // Construct displacement field.
2176  const pointBoundaryMesh& pointPatches = pMesh.boundary();
2177 
2178  wordList patchFieldTypes
2179  (
2180  pointPatches.size(),
2181  slipPointPatchVectorField::typeName
2182  );
2183 
2184  forAll(adaptPatchIDs, i)
2185  {
2186  patchFieldTypes[adaptPatchIDs[i]] =
2187  fixedValuePointPatchVectorField::typeName;
2188  }
2189 
2190  forAll(pointPatches, patchi)
2191  {
2192  if (isA<processorPointPatch>(pointPatches[patchi]))
2193  {
2194  patchFieldTypes[patchi] = calculatedPointPatchVectorField::typeName;
2195  }
2196  else if (isA<cyclicPointPatch>(pointPatches[patchi]))
2197  {
2198  patchFieldTypes[patchi] = cyclicSlipPointPatchVectorField::typeName;
2199  }
2200  }
2201 
2202  // Note: time().timeName() instead of meshRefinement::timeName() since
2203  // postprocessable field.
2204  return tmp<pointVectorField>::New
2205  (
2206  IOobject
2207  (
2208  "pointDisplacement",
2209  mesh.time().timeName(),
2210  mesh,
2211  IOobject::NO_READ,
2212  IOobject::AUTO_WRITE
2213  ),
2214  pMesh,
2215  dimensionedVector(dimLength, Zero),
2216  patchFieldTypes
2217  );
2218 }
2219 
2220 
2221 void Foam::meshRefinement::checkCoupledFaceZones(const polyMesh& mesh)
2222 {
2223  const faceZoneMesh& fZones = mesh.faceZones();
2224 
2225  // Check any zones are present anywhere and in same order
2226 
2227  {
2228  List<wordList> zoneNames(Pstream::nProcs());
2229  zoneNames[Pstream::myProcNo()] = fZones.names();
2230  Pstream::allGatherList(zoneNames);
2231 
2232  // All have same data now. Check.
2233  forAll(zoneNames, proci)
2234  {
2235  if
2236  (
2237  proci != Pstream::myProcNo()
2238  && (zoneNames[proci] != zoneNames[Pstream::myProcNo()])
2239  )
2240  {
2241  FatalErrorInFunction
2242  << "faceZones are not synchronised on processors." << nl
2243  << "Processor " << proci << " has faceZones "
2244  << zoneNames[proci] << nl
2245  << "Processor " << Pstream::myProcNo()
2246  << " has faceZones "
2247  << zoneNames[Pstream::myProcNo()] << nl
2248  << exit(FatalError);
2249  }
2250  }
2251  }
2252 
2253  // Check that coupled faces are present on both sides.
2254 
2255  labelList faceToZone(mesh.nBoundaryFaces(), -1);
2256 
2257  forAll(fZones, zonei)
2258  {
2259  const faceZone& fZone = fZones[zonei];
2260 
2261  forAll(fZone, i)
2262  {
2263  label bFacei = fZone[i]-mesh.nInternalFaces();
2264 
2265  if (bFacei >= 0)
2266  {
2267  if (faceToZone[bFacei] == -1)
2268  {
2269  faceToZone[bFacei] = zonei;
2270  }
2271  else if (faceToZone[bFacei] == zonei)
2272  {
2273  FatalErrorInFunction
2274  << "Face " << fZone[i] << " in zone "
2275  << fZone.name()
2276  << " is twice in zone!"
2277  << abort(FatalError);
2278  }
2279  else
2280  {
2281  FatalErrorInFunction
2282  << "Face " << fZone[i] << " in zone "
2283  << fZone.name()
2284  << " is also in zone "
2285  << fZones[faceToZone[bFacei]].name()
2286  << abort(FatalError);
2287  }
2288  }
2289  }
2290  }
2291 
2292  labelList neiFaceToZone(faceToZone);
2293  syncTools::swapBoundaryFaceList(mesh, neiFaceToZone);
2294 
2295  forAll(faceToZone, i)
2296  {
2297  if (faceToZone[i] != neiFaceToZone[i])
2298  {
2299  FatalErrorInFunction
2300  << "Face " << mesh.nInternalFaces()+i
2301  << " is in zone " << faceToZone[i]
2302  << ", its coupled face is in zone " << neiFaceToZone[i]
2303  << abort(FatalError);
2304  }
2305  }
2306 }
2307 
2308 
2309 void Foam::meshRefinement::calculateEdgeWeights
2310 (
2311  const polyMesh& mesh,
2312  const bitSet& isMasterEdge,
2313  const labelList& meshPoints,
2314  const edgeList& edges,
2315  scalarField& edgeWeights,
2316  scalarField& invSumWeight
2317 )
2318 {
2319  const pointField& pts = mesh.points();
2320 
2321  // Calculate edgeWeights and inverse sum of edge weights
2322  edgeWeights.setSize(isMasterEdge.size());
2323  invSumWeight.setSize(meshPoints.size());
2324 
2325  forAll(edges, edgei)
2326  {
2327  const edge& e = edges[edgei];
2328  scalar eMag = max
2329  (
2330  SMALL,
2331  mag
2332  (
2333  pts[meshPoints[e[1]]]
2334  - pts[meshPoints[e[0]]]
2335  )
2336  );
2337  edgeWeights[edgei] = 1.0/eMag;
2338  }
2339 
2340  // Sum per point all edge weights
2341  weightedSum
2342  (
2343  mesh,
2344  isMasterEdge,
2345  meshPoints,
2346  edges,
2347  edgeWeights,
2348  scalarField(meshPoints.size(), 1.0), // data
2349  invSumWeight
2350  );
2351 
2352  // Inplace invert
2353  forAll(invSumWeight, pointi)
2354  {
2355  scalar w = invSumWeight[pointi];
2356 
2357  if (w > 0.0)
2358  {
2359  invSumWeight[pointi] = 1.0/w;
2360  }
2361  }
2362 }
2363 
2364 
2365 Foam::label Foam::meshRefinement::appendPatch
2366 (
2367  fvMesh& mesh,
2368  const label insertPatchi,
2369  const word& patchName,
2370  const dictionary& patchDict
2371 )
2372 {
2373  // Clear local fields and e.g. polyMesh parallelInfo.
2374  mesh.clearOut();
2375 
2376  polyBoundaryMesh& polyPatches =
2377  const_cast<polyBoundaryMesh&>(mesh.boundaryMesh());
2378  fvBoundaryMesh& fvPatches = const_cast<fvBoundaryMesh&>(mesh.boundary());
2379 
2380  label patchi = polyPatches.size();
2381 
2382  // Add polyPatch at the end
2383  polyPatches.setSize(patchi+1);
2384  polyPatches.set
2385  (
2386  patchi,
2387  polyPatch::New
2388  (
2389  patchName,
2390  patchDict,
2391  insertPatchi,
2392  polyPatches
2393  )
2394  );
2395  fvPatches.setSize(patchi+1);
2396  fvPatches.set
2397  (
2398  patchi,
2399  fvPatch::New
2400  (
2401  polyPatches[patchi], // point to newly added polyPatch
2402  mesh.boundary()
2403  )
2404  );
2405 
2406  addPatchFields<volScalarField>
2407  (
2408  mesh,
2409  calculatedFvPatchField<scalar>::typeName
2410  );
2411  addPatchFields<volVectorField>
2412  (
2413  mesh,
2414  calculatedFvPatchField<vector>::typeName
2415  );
2416  addPatchFields<volSphericalTensorField>
2417  (
2418  mesh,
2419  calculatedFvPatchField<sphericalTensor>::typeName
2420  );
2421  addPatchFields<volSymmTensorField>
2422  (
2423  mesh,
2424  calculatedFvPatchField<symmTensor>::typeName
2425  );
2426  addPatchFields<volTensorField>
2427  (
2428  mesh,
2429  calculatedFvPatchField<tensor>::typeName
2430  );
2431 
2432  // Surface fields
2433 
2434  addPatchFields<surfaceScalarField>
2435  (
2436  mesh,
2437  calculatedFvPatchField<scalar>::typeName
2438  );
2439  addPatchFields<surfaceVectorField>
2440  (
2441  mesh,
2442  calculatedFvPatchField<vector>::typeName
2443  );
2444  addPatchFields<surfaceSphericalTensorField>
2445  (
2446  mesh,
2447  calculatedFvPatchField<sphericalTensor>::typeName
2448  );
2449  addPatchFields<surfaceSymmTensorField>
2450  (
2451  mesh,
2452  calculatedFvPatchField<symmTensor>::typeName
2453  );
2454  addPatchFields<surfaceTensorField>
2455  (
2456  mesh,
2457  calculatedFvPatchField<tensor>::typeName
2458  );
2459  return patchi;
2460 }
2461 
2462 
2463 Foam::label Foam::meshRefinement::addPatch
2464 (
2465  fvMesh& mesh,
2466  const word& patchName,
2467  const dictionary& patchInfo
2468 )
2469 {
2470  polyBoundaryMesh& polyPatches =
2471  const_cast<polyBoundaryMesh&>(mesh.boundaryMesh());
2472  fvBoundaryMesh& fvPatches = const_cast<fvBoundaryMesh&>(mesh.boundary());
2473 
2474  const label patchi = polyPatches.findPatchID(patchName);
2475  if (patchi != -1)
2476  {
2477  // Already there
2478  return patchi;
2479  }
2480 
2481 
2482  label insertPatchi = polyPatches.size();
2483  label startFacei = mesh.nFaces();
2484 
2485  forAll(polyPatches, patchi)
2486  {
2487  const polyPatch& pp = polyPatches[patchi];
2488 
2489  if (isA<processorPolyPatch>(pp))
2490  {
2491  insertPatchi = patchi;
2492  startFacei = pp.start();
2493  break;
2494  }
2495  }
2496 
2497  dictionary patchDict(patchInfo);
2498  patchDict.set("nFaces", 0);
2499  patchDict.set("startFace", startFacei);
2500 
2501  // Below is all quite a hack. Feel free to change once there is a better
2502  // mechanism to insert and reorder patches.
2503 
2504  label addedPatchi = appendPatch(mesh, insertPatchi, patchName, patchDict);
2505 
2506 
2507  // Create reordering list
2508  // patches before insert position stay as is
2509  labelList oldToNew(addedPatchi+1);
2510  for (label i = 0; i < insertPatchi; i++)
2511  {
2512  oldToNew[i] = i;
2513  }
2514  // patches after insert position move one up
2515  for (label i = insertPatchi; i < addedPatchi; i++)
2516  {
2517  oldToNew[i] = i+1;
2518  }
2519  // appended patch gets moved to insert position
2520  oldToNew[addedPatchi] = insertPatchi;
2521 
2522  // Shuffle into place
2523  polyPatches.reorder(oldToNew, true);
2524  fvPatches.reorder(oldToNew);
2525 
2526  reorderPatchFields<volScalarField>(mesh, oldToNew);
2527  reorderPatchFields<volVectorField>(mesh, oldToNew);
2528  reorderPatchFields<volSphericalTensorField>(mesh, oldToNew);
2529  reorderPatchFields<volSymmTensorField>(mesh, oldToNew);
2530  reorderPatchFields<volTensorField>(mesh, oldToNew);
2531  reorderPatchFields<surfaceScalarField>(mesh, oldToNew);
2532  reorderPatchFields<surfaceVectorField>(mesh, oldToNew);
2533  reorderPatchFields<surfaceSphericalTensorField>(mesh, oldToNew);
2534  reorderPatchFields<surfaceSymmTensorField>(mesh, oldToNew);
2535  reorderPatchFields<surfaceTensorField>(mesh, oldToNew);
2536 
2537  return insertPatchi;
2538 }
2539 
2540 
2541 Foam::label Foam::meshRefinement::addMeshedPatch
2542 (
2543  const word& name,
2544  const dictionary& patchInfo
2545 )
2546 {
2547  label meshedi = meshedPatches_.find(name);
2548 
2549  if (meshedi != -1)
2550  {
2551  // Already there. Get corresponding polypatch
2552  return mesh_.boundaryMesh().findPatchID(name);
2553  }
2554  else
2555  {
2556  // Add patch
2557  label patchi = addPatch(mesh_, name, patchInfo);
2558 
2559 // dictionary patchDict(patchInfo);
2560 // patchDict.set("nFaces", 0);
2561 // patchDict.set("startFace", 0);
2562 // autoPtr<polyPatch> ppPtr
2563 // (
2564 // polyPatch::New
2565 // (
2566 // name,
2567 // patchDict,
2568 // 0,
2569 // mesh_.boundaryMesh()
2570 // )
2571 // );
2572 // label patchi = fvMeshTools::addPatch
2573 // (
2574 // mesh_,
2575 // ppPtr(),
2576 // dictionary(), // optional field values
2577 // calculatedFvPatchField<scalar>::typeName,
2578 // true
2579 // );
2580 
2581  // Store
2582  meshedPatches_.append(name);
2583 
2584  // Clear patch based addressing
2585  faceToCoupledPatch_.clear();
2586 
2587  return patchi;
2588  }
2589 }
2590 
2591 
2592 Foam::labelList Foam::meshRefinement::meshedPatches() const
2593 {
2594  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
2595 
2596  DynamicList<label> patchIDs(meshedPatches_.size());
2597  forAll(meshedPatches_, i)
2598  {
2599  label patchi = patches.findPatchID(meshedPatches_[i]);
2600 
2601  if (patchi == -1)
2602  {
2603  WarningInFunction
2604  << "Problem : did not find patch " << meshedPatches_[i]
2605  << endl << "Valid patches are " << patches.names()
2606  << endl; //abort(FatalError);
2607  }
2608  else if (!polyPatch::constraintType(patches[patchi].type()))
2609  {
2610  patchIDs.append(patchi);
2611  }
2612  }
2613 
2614  return patchIDs;
2615 }
2616 
2617 
2618 Foam::label Foam::meshRefinement::addFaceZone
2619 (
2620  const word& fzName,
2621  const word& masterPatch,
2622  const word& slavePatch,
2623  const surfaceZonesInfo::faceZoneType& fzType
2624 )
2625 {
2626  label zonei = surfaceZonesInfo::addFaceZone
2627  (
2628  fzName, //name
2629  labelList(0), //addressing
2630  boolList(0), //flipmap
2631  mesh_
2632  );
2633 
2634  faceZoneToMasterPatch_.insert(fzName, masterPatch);
2635  faceZoneToSlavePatch_.insert(fzName, slavePatch);
2636  faceZoneToType_.insert(fzName, fzType);
2637 
2638  return zonei;
2639 }
2640 
2641 
2642 bool Foam::meshRefinement::getFaceZoneInfo
2643 (
2644  const word& fzName,
2645  label& masterPatchID,
2646  label& slavePatchID,
2647  surfaceZonesInfo::faceZoneType& fzType
2648 ) const
2649 {
2650  const polyBoundaryMesh& pbm = mesh_.boundaryMesh();
2651 
2652  if (!faceZoneToMasterPatch_.found(fzName))
2653  {
2654  return false;
2655  }
2656  else
2657  {
2658  const word& masterName = faceZoneToMasterPatch_[fzName];
2659  masterPatchID = pbm.findPatchID(masterName);
2660 
2661  const word& slaveName = faceZoneToSlavePatch_[fzName];
2662  slavePatchID = pbm.findPatchID(slaveName);
2663 
2664  fzType = faceZoneToType_[fzName];
2665 
2666  return true;
2667  }
2668 }
2669 
2670 
2671 Foam::label Foam::meshRefinement::addPointZone(const word& name)
2672 {
2673  pointZoneMesh& pointZones = mesh_.pointZones();
2674 
2675  label zoneI = pointZones.findZoneID(name);
2676 
2677  if (zoneI == -1)
2678  {
2679  zoneI = pointZones.size();
2680  pointZones.clearAddressing();
2681  pointZones.setSize(zoneI+1);
2682  pointZones.set
2683  (
2684  zoneI,
2685  new pointZone
2686  (
2687  name, // name
2688  labelList(0), // addressing
2689  zoneI, // index
2690  pointZones // pointZoneMesh
2691  )
2692  );
2693  }
2694  return zoneI;
2695 }
2696 
2697 
2698 void Foam::meshRefinement::selectSeparatedCoupledFaces(boolList& selected) const
2699 {
2700  for (const polyPatch& pp : mesh_.boundaryMesh())
2701  {
2702  // Check all coupled. Avoid using .coupled() so we also pick up AMI.
2703  const auto* cpp = isA<coupledPolyPatch>(pp);
2704 
2705  if (cpp && (cpp->separated() || !cpp->parallel()))
2706  {
2707  SubList<bool>(selected, pp.size(), pp.start()) = true;
2708  }
2709  }
2710 }
2711 
2712 
2713 Foam::labelList Foam::meshRefinement::countEdgeFaces
2714 (
2715  const uindirectPrimitivePatch& pp
2716 ) const
2717 {
2718  // Count number of faces per edge. Parallel consistent.
2719 
2720  const labelListList& edgeFaces = pp.edgeFaces();
2721  labelList nEdgeFaces(edgeFaces.size());
2722  forAll(edgeFaces, edgei)
2723  {
2724  nEdgeFaces[edgei] = edgeFaces[edgei].size();
2725  }
2726 
2727  // Sync across processor patches
2728  if (Pstream::parRun())
2729  {
2730  const globalMeshData& globalData = mesh_.globalData();
2731  const mapDistribute& map = globalData.globalEdgeSlavesMap();
2732  const indirectPrimitivePatch& cpp = globalData.coupledPatch();
2733 
2734  // Match pp edges to coupled edges
2735  labelList patchEdges;
2736  labelList coupledEdges;
2737  PackedBoolList sameEdgeOrientation;
2738  PatchTools::matchEdges
2739  (
2740  pp,
2741  cpp,
2742  patchEdges,
2743  coupledEdges,
2744  sameEdgeOrientation
2745  );
2746 
2747  // Convert patch-edge data into cpp-edge data
2748  labelList coupledNEdgeFaces(map.constructSize(), Zero);
2749  UIndirectList<label>(coupledNEdgeFaces, coupledEdges) =
2750  UIndirectList<label>(nEdgeFaces, patchEdges);
2751 
2752  // Synchronise
2753  globalData.syncData
2754  (
2755  coupledNEdgeFaces,
2756  globalData.globalEdgeSlaves(),
2757  globalData.globalEdgeTransformedSlaves(),
2758  map,
2759  plusEqOp<label>()
2760  );
2761 
2762  // Convert back from cpp-edge to patch-edge
2763  UIndirectList<label>(nEdgeFaces, patchEdges) =
2764  UIndirectList<label>(coupledNEdgeFaces, coupledEdges);
2765  }
2766  return nEdgeFaces;
2767 }
2768 
2769 
2770 Foam::label Foam::meshRefinement::findCell
2771 (
2772  const polyMesh& mesh,
2773  const vector& perturbVec,
2774  const point& p
2775 )
2776 {
2777  // Force calculation of base points (needs to be synchronised)
2778  (void)mesh.tetBasePtIs();
2779 
2780  label celli = mesh.findCell(p, findCellMode);
2781 
2782  if (returnReduceAnd(celli < 0))
2783  {
2784  // See if we can perturb a bit
2785  celli = mesh.findCell(p+perturbVec, findCellMode);
2786  }
2787  return celli;
2788 }
2789 
2790 
2791 Foam::label Foam::meshRefinement::findRegion
2792 (
2793  const polyMesh& mesh,
2794  const labelList& cellToRegion,
2795  const vector& perturbVec,
2796  const point& p
2797 )
2798 {
2799  label regioni = -1;
2800 
2801  // Force calculation of base points (needs to be synchronised)
2802  (void)mesh.tetBasePtIs();
2803 
2804  label celli = mesh.findCell(p, findCellMode);
2805  if (celli != -1)
2806  {
2807  regioni = cellToRegion[celli];
2808  }
2809  reduce(regioni, maxOp<label>());
2810 
2811  if (regioni == -1)
2812  {
2813  // See if we can perturb a bit
2814  celli = mesh.findCell(p+perturbVec, findCellMode);
2815  if (celli != -1)
2816  {
2817  regioni = cellToRegion[celli];
2818  }
2819  reduce(regioni, maxOp<label>());
2820  }
2821  return regioni;
2822 }
2823 
2824 
2825 Foam::fileName Foam::meshRefinement::writeLeakPath
2826 (
2827  const polyMesh& mesh,
2828  const pointField& locationsInMesh,
2829  const pointField& locationsOutsideMesh,
2830  const boolList& blockedFace,
2831  coordSetWriter& writer
2832 )
2833 {
2834  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
2835 
2836  fileName outputDir
2837  (
2838  mesh.time().globalPath()
2839  / functionObject::outputPrefix
2840  / mesh.pointsInstance()
2841  );
2842  outputDir.clean(); // Remove unneeded ".."
2843 
2844  // Write the leak path
2845 
2846  meshSearch searchEngine(mesh);
2847  shortestPathSet leakPath
2848  (
2849  "leakPath",
2850  mesh,
2851  searchEngine,
2852  coordSet::coordFormatNames[coordSet::coordFormat::DISTANCE],
2853  false, //true,
2854  50, // tbd. Number of iterations
2855  pbm.groupPatchIDs()["wall"],
2856  locationsInMesh,
2857  locationsOutsideMesh,
2858  blockedFace
2859  );
2860 
2861  // Split leak path according to segment. Note: segment index
2862  // is global (= index in locationsInsideMesh)
2863  List<pointList> segmentPoints;
2864  List<scalarList> segmentDist;
2865  {
2866  label nSegments = 0;
2867  if (leakPath.segments().size())
2868  {
2869  nSegments = max(leakPath.segments())+1;
2870  }
2871  reduce(nSegments, maxOp<label>());
2872 
2873  labelList nElemsPerSegment(nSegments, Zero);
2874  for (label segmenti : leakPath.segments())
2875  {
2876  nElemsPerSegment[segmenti]++;
2877  }
2878  segmentPoints.setSize(nElemsPerSegment.size());
2879  segmentDist.setSize(nElemsPerSegment.size());
2880  forAll(nElemsPerSegment, i)
2881  {
2882  segmentPoints[i].setSize(nElemsPerSegment[i]);
2883  segmentDist[i].setSize(nElemsPerSegment[i]);
2884  }
2885  nElemsPerSegment = 0;
2886 
2887  forAll(leakPath, elemi)
2888  {
2889  label segmenti = leakPath.segments()[elemi];
2890  pointList& points = segmentPoints[segmenti];
2891  scalarList& dist = segmentDist[segmenti];
2892  label& n = nElemsPerSegment[segmenti];
2893 
2894  points[n] = leakPath[elemi];
2895  dist[n] = leakPath.distance()[elemi];
2896  n++;
2897  }
2898  }
2899 
2900  PtrList<coordSet> allLeakPaths(segmentPoints.size());
2901  forAll(allLeakPaths, segmenti)
2902  {
2903  // Collect data from all processors
2904  List<pointList> gatheredPts(Pstream::nProcs());
2905  gatheredPts[Pstream::myProcNo()] = std::move(segmentPoints[segmenti]);
2906  Pstream::gatherList(gatheredPts);
2907 
2908  List<scalarList> gatheredDist(Pstream::nProcs());
2909  gatheredDist[Pstream::myProcNo()] = std::move(segmentDist[segmenti]);
2910  Pstream::gatherList(gatheredDist);
2911 
2912  // Combine processor lists into one big list.
2913  pointList allPts
2914  (
2915  ListListOps::combine<pointList>
2916  (
2917  gatheredPts, accessOp<pointList>()
2918  )
2919  );
2920  scalarList allDist
2921  (
2922  ListListOps::combine<scalarList>
2923  (
2924  gatheredDist, accessOp<scalarList>()
2925  )
2926  );
2927 
2928  // Sort according to distance
2929  labelList indexSet(Foam::sortedOrder(allDist));
2930 
2931  allLeakPaths.set
2932  (
2933  segmenti,
2934  new coordSet
2935  (
2936  leakPath.name(),
2937  leakPath.axis(),
2938  pointList(allPts, indexSet),
2939  //scalarList(allDist, indexSet)
2940  scalarList(allPts.size(), scalar(segmenti))
2941  )
2942  );
2943  }
2944 
2945  fileName fName;
2946  if (Pstream::master())
2947  {
2948  List<scalarField> allLeakData(allLeakPaths.size());
2949  forAll(allLeakPaths, segmenti)
2950  {
2951  allLeakData[segmenti] = allLeakPaths[segmenti].distance();
2952  }
2953 
2954  writer.nFields(1);
2955 
2956  writer.open
2957  (
2958  allLeakPaths,
2959  (outputDir / allLeakPaths[0].name())
2960  );
2961 
2962  fName = writer.write("leakPath", allLeakData);
2963 
2964  // Force writing to finish before FatalError exit
2965  writer.close(true);
2966  }
2967 
2968  // Probably do not need to broadcast name (only written on master anyhow)
2969  Pstream::broadcast(fName);
2970 
2971  return fName;
2972 }
2973 
2974 
2975 // Modify cellRegion to be consistent with locationsInMesh.
2976 // - all regions not in locationsInMesh are set to -1
2977 // - check that all regions inside locationsOutsideMesh are set to -1
2978 Foam::label Foam::meshRefinement::findRegions
2979 (
2980  const polyMesh& mesh,
2981  const vector& perturbVec,
2982  const pointField& locationsInMesh,
2983  const pointField& locationsOutsideMesh,
2984  const label nRegions,
2985  labelList& cellRegion,
2986  const boolList& blockedFace,
2987  // Leak-path
2988  const bool exitIfLeakPath,
2989  const refPtr<coordSetWriter>& leakPathFormatter
2990 )
2991 {
2992  bitSet insideCell(mesh.nCells());
2993 
2994  // Mark all cells reachable from locationsInMesh
2995  labelList insideRegions(locationsInMesh.size());
2996  forAll(insideRegions, i)
2997  {
2998  // Find the region containing the point
2999  label regioni = findRegion
3000  (
3001  mesh,
3002  cellRegion,
3003  perturbVec,
3004  locationsInMesh[i]
3005  );
3006 
3007  insideRegions[i] = regioni;
3008 
3009  // Mark all cells in the region as being inside
3010  forAll(cellRegion, celli)
3011  {
3012  if (cellRegion[celli] == regioni)
3013  {
3014  insideCell.set(celli);
3015  }
3016  }
3017  }
3018 
3019 
3020 
3021  // Check that all the locations outside the
3022  // mesh do not conflict with those inside
3023  forAll(locationsOutsideMesh, i)
3024  {
3025  // Find the region containing the point,
3026  // and the corresponding inside region index
3027 
3028  label indexi;
3029  label regioni = findRegion
3030  (
3031  mesh,
3032  cellRegion,
3033  perturbVec,
3034  locationsOutsideMesh[i]
3035  );
3036 
3037  if (regioni == -1 && (indexi = insideRegions.find(regioni)) != -1)
3038  {
3039  if (leakPathFormatter)
3040  {
3041  const fileName fName
3042  (
3043  writeLeakPath
3044  (
3045  mesh,
3046  locationsInMesh,
3047  locationsOutsideMesh,
3048  blockedFace,
3049  leakPathFormatter.constCast()
3050  )
3051  );
3052  Info<< "Dumped leak path to " << fName << endl;
3053  }
3054 
3055  auto& err =
3056  (
3057  exitIfLeakPath
3058  ? FatalErrorInFunction
3059  : WarningInFunction
3060  );
3061 
3062  err << "Location in mesh " << locationsInMesh[indexi]
3063  << " is inside same mesh region " << regioni
3064  << " as one of the locations outside mesh "
3065  << locationsOutsideMesh << endl;
3066 
3067  if (exitIfLeakPath)
3068  {
3069  FatalError << exit(FatalError);
3070  }
3071  }
3072  }
3073 
3074 
3075  label nRemove = 0;
3076 
3077  // Now update cellRegion to -1 for unreachable cells
3078  forAll(insideCell, celli)
3079  {
3080  if (!insideCell.test(celli))
3081  {
3082  cellRegion[celli] = -1;
3083  ++nRemove;
3084  }
3085  else if (cellRegion[celli] == -1)
3086  {
3087  ++nRemove;
3088  }
3089  }
3090 
3091  return nRemove;
3092 }
3093 
3094 
3095 Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::splitMeshRegions
3096 (
3097  const labelList& globalToMasterPatch,
3098  const labelList& globalToSlavePatch,
3099  const pointField& locationsInMesh,
3100  const pointField& locationsOutsideMesh,
3101  const bool exitIfLeakPath,
3102  const refPtr<coordSetWriter>& leakPathFormatter
3103 )
3104 {
3105  // Force calculation of face decomposition (used in findCell)
3106  (void)mesh_.tetBasePtIs();
3107 
3108  // Determine connected regions. regionSplit is the labelList with the
3109  // region per cell.
3110 
3111  boolList blockedFace(mesh_.nFaces(), false);
3112  selectSeparatedCoupledFaces(blockedFace);
3113 
3114  regionSplit cellRegion(mesh_, blockedFace);
3115 
3116  label nRemove = findRegions
3117  (
3118  mesh_,
3119  vector::uniform(mergeDistance_), // perturbVec
3120  locationsInMesh,
3121  locationsOutsideMesh,
3122  cellRegion.nRegions(),
3123  cellRegion,
3124  blockedFace,
3125  // Leak-path
3126  exitIfLeakPath,
3127  leakPathFormatter
3128  );
3129 
3130  // Subset
3131  // ~~~~~~
3132 
3133  // Get cells to remove
3134  DynamicList<label> cellsToRemove(nRemove);
3135  forAll(cellRegion, celli)
3136  {
3137  if (cellRegion[celli] == -1)
3138  {
3139  cellsToRemove.append(celli);
3140  }
3141  }
3142  cellsToRemove.shrink();
3143 
3144  label nTotCellsToRemove = returnReduce
3145  (
3146  cellsToRemove.size(),
3147  sumOp<label>()
3148  );
3149 
3150 
3151  autoPtr<mapPolyMesh> mapPtr;
3152  if (nTotCellsToRemove > 0)
3153  {
3154  label nCellsToKeep =
3155  mesh_.globalData().nTotalCells()
3156  - nTotCellsToRemove;
3157 
3158  Info<< "Keeping all cells containing points " << locationsInMesh << endl
3159  << "Selected for keeping : "
3160  << nCellsToKeep
3161  << " cells." << endl;
3162 
3163 
3164  // Remove cells
3165  removeCells cellRemover(mesh_);
3166 
3167  labelList exposedFaces(cellRemover.getExposedFaces(cellsToRemove));
3168  labelList exposedPatch;
3169 
3170  label nExposedFaces = returnReduce(exposedFaces.size(), sumOp<label>());
3171  if (nExposedFaces)
3172  {
3173  // FatalErrorInFunction
3174  // << "Removing non-reachable cells should only expose"
3175  // << " boundary faces" << nl
3176  // << "ExposedFaces:" << exposedFaces << abort(FatalError);
3177 
3178  // Patch for exposed faces for lack of anything sensible.
3179  label defaultPatch = 0;
3180  if (globalToMasterPatch.size())
3181  {
3182  defaultPatch = globalToMasterPatch[0];
3183  }
3184 
3185  WarningInFunction
3186  << "Removing non-reachable cells exposes "
3187  << nExposedFaces << " internal or coupled faces." << endl
3188  << " These get put into patch " << defaultPatch << endl;
3189  exposedPatch.setSize(exposedFaces.size(), defaultPatch);
3190  }
3191 
3192  mapPtr = doRemoveCells
3193  (
3194  cellsToRemove,
3195  exposedFaces,
3196  exposedPatch,
3197  cellRemover
3198  );
3199  }
3200  return mapPtr;
3201 }
3202 
3203 
3204 void Foam::meshRefinement::distribute(const mapDistributePolyMesh& map)
3205 {
3206  // mesh_ already distributed; distribute my member data
3207 
3208  // surfaceQueries_ ok.
3209 
3210  // refinement
3211  meshCutter_.distribute(map);
3212 
3213  // surfaceIndex is face data.
3214  map.distributeFaceData(surfaceIndex_);
3215 
3216  // faceToPatch (baffles that were on coupled faces) is not maintained
3217  // (since baffling also disconnects points)
3218  faceToCoupledPatch_.clear();
3219 
3220  // maintainedFaces are indices of faces.
3221  forAll(userFaceData_, i)
3222  {
3223  map.distributeFaceData(userFaceData_[i].second());
3224  }
3225 
3226  // Redistribute surface and any fields on it.
3227  {
3228  Random rndGen(653213);
3229 
3230  // Get local mesh bounding box. Single box for now.
3231  List<treeBoundBox> meshBb
3232  (
3233  1,
3234  treeBoundBox(mesh_.points()).extend(rndGen, 1e-4)
3235  );
3236 
3237  // Distribute all geometry (so refinementSurfaces and shellSurfaces)
3238  searchableSurfaces& geometry =
3239  const_cast<searchableSurfaces&>(surfaces_.geometry());
3240 
3241  forAll(geometry, i)
3242  {
3243  autoPtr<mapDistribute> faceMap;
3244  autoPtr<mapDistribute> pointMap;
3245  geometry[i].distribute
3246  (
3247  meshBb,
3248  false, // do not keep outside triangles
3249  faceMap,
3250  pointMap
3251  );
3252 
3253  if (faceMap)
3254  {
3255  // (ab)use the instance() to signal current modification time
3256  geometry[i].instance() = geometry[i].time().timeName();
3257  }
3258 
3259  faceMap.clear();
3260  pointMap.clear();
3261  }
3262  }
3263 }
3264 
3265 
3266 void Foam::meshRefinement::updateMesh
3267 (
3268  const mapPolyMesh& map,
3269  const labelList& changedFaces
3270 )
3271 {
3272  Map<label> dummyMap(0);
3273 
3274  updateMesh(map, changedFaces, dummyMap, dummyMap, dummyMap);
3275 }
3276 
3277 
3278 void Foam::meshRefinement::storeData
3279 (
3280  const labelList& pointsToStore,
3281  const labelList& facesToStore,
3282  const labelList& cellsToStore
3283 )
3284 {
3285  // For now only meshCutter has storable/retrievable data.
3286  meshCutter_.storeData
3287  (
3288  pointsToStore,
3289  facesToStore,
3290  cellsToStore
3291  );
3292 }
3293 
3294 
3295 void Foam::meshRefinement::updateMesh
3296 (
3297  const mapPolyMesh& map,
3298  const labelList& changedFaces,
3299  const Map<label>& pointsToRestore,
3300  const Map<label>& facesToRestore,
3301  const Map<label>& cellsToRestore
3302 )
3303 {
3304  // For now only meshCutter has storable/retrievable data.
3305 
3306  // Update numbering of cells/vertices.
3307  meshCutter_.updateMesh
3308  (
3309  map,
3310  pointsToRestore,
3311  facesToRestore,
3312  cellsToRestore
3313  );
3314 
3315  // Update surfaceIndex
3316  updateList(map.faceMap(), label(-1), surfaceIndex_);
3317 
3318  // Update faceToCoupledPatch_
3319  {
3320  Map<label> newFaceToPatch(faceToCoupledPatch_.size());
3321  forAllConstIters(faceToCoupledPatch_, iter)
3322  {
3323  const label newFacei = map.reverseFaceMap()[iter.key()];
3324 
3325  if (newFacei >= 0)
3326  {
3327  newFaceToPatch.insert(newFacei, iter.val());
3328  }
3329  }
3330  faceToCoupledPatch_.transfer(newFaceToPatch);
3331  }
3332 
3333 
3334  // Update cached intersection information
3335  updateIntersections(changedFaces);
3336 
3337  // Update maintained faces
3338  forAll(userFaceData_, i)
3339  {
3340  labelList& data = userFaceData_[i].second();
3341 
3342  if (userFaceData_[i].first() == KEEPALL)
3343  {
3344  // extend list with face-from-face data
3345  updateList(map.faceMap(), label(-1), data);
3346  }
3347  else if (userFaceData_[i].first() == MASTERONLY)
3348  {
3349  // keep master only
3350  labelList newFaceData(map.faceMap().size(), -1);
3351 
3352  forAll(newFaceData, facei)
3353  {
3354  label oldFacei = map.faceMap()[facei];
3355 
3356  if (oldFacei >= 0 && map.reverseFaceMap()[oldFacei] == facei)
3357  {
3358  newFaceData[facei] = data[oldFacei];
3359  }
3360  }
3361  data.transfer(newFaceData);
3362  }
3363  else
3364  {
3365  // remove any face that has been refined i.e. referenced more than
3366  // once.
3367 
3368  // 1. Determine all old faces that get referenced more than once.
3369  // These get marked with -1 in reverseFaceMap
3370  labelList reverseFaceMap(map.reverseFaceMap());
3371 
3372  forAll(map.faceMap(), facei)
3373  {
3374  label oldFacei = map.faceMap()[facei];
3375 
3376  if (oldFacei >= 0)
3377  {
3378  if (reverseFaceMap[oldFacei] != facei)
3379  {
3380  // facei is slave face. Mark old face.
3381  reverseFaceMap[oldFacei] = -1;
3382  }
3383  }
3384  }
3385 
3386  // 2. Map only faces with intact reverseFaceMap
3387  labelList newFaceData(map.faceMap().size(), -1);
3388  forAll(newFaceData, facei)
3389  {
3390  label oldFacei = map.faceMap()[facei];
3391 
3392  if (oldFacei >= 0)
3393  {
3394  if (reverseFaceMap[oldFacei] == facei)
3395  {
3396  newFaceData[facei] = data[oldFacei];
3397  }
3398  }
3399  }
3400  data.transfer(newFaceData);
3401  }
3402  }
3403 }
3404 
3405 
3406 bool Foam::meshRefinement::write() const
3407 {
3408  bool writeOk = mesh_.write();
3409 
3410  // Make sure that any distributed surfaces (so ones which probably have
3411  // been changed) get written as well.
3412  // Note: should ideally have some 'modified' flag to say whether it
3413  // has been changed or not.
3414  searchableSurfaces& geometry =
3415  const_cast<searchableSurfaces&>(surfaces_.geometry());
3416 
3417  forAll(geometry, i)
3418  {
3419  searchableSurface& s = geometry[i];
3420 
3421  // Check if instance() of surface is not constant or system.
3422  // Is good hint that surface is distributed.
3423  if
3424  (
3425  s.instance() != s.time().system()
3426  && s.instance() != s.time().caseSystem()
3427  && s.instance() != s.time().constant()
3428  && s.instance() != s.time().caseConstant()
3429  )
3430  {
3431  // Make sure it gets written to current time, not constant.
3432  s.instance() = s.time().timeName();
3433  writeOk = writeOk && s.write();
3434  }
3435  }
3436 
3437  return writeOk;
3438 }
3439 
3440 
3441 Foam::bitSet Foam::meshRefinement::getMasterPoints
3442 (
3443  const polyMesh& mesh,
3444  const labelList& meshPoints
3445 )
3446 {
3447  const globalIndex globalPoints(meshPoints.size());
3448 
3449  labelList myPoints
3450  (
3451  identity(globalPoints.localSize(), globalPoints.localStart())
3452  );
3453 
3454  syncTools::syncPointList
3455  (
3456  mesh,
3457  meshPoints,
3458  myPoints,
3459  minEqOp<label>(),
3460  labelMax
3461  );
3462 
3463 
3464  bitSet isPatchMasterPoint(meshPoints.size());
3465  forAll(meshPoints, pointi)
3466  {
3467  if (myPoints[pointi] == globalPoints.toGlobal(pointi))
3468  {
3469  isPatchMasterPoint.set(pointi);
3470  }
3471  }
3472 
3473  return isPatchMasterPoint;
3474 }
3475 
3476 
3477 Foam::bitSet Foam::meshRefinement::getMasterEdges
3478 (
3479  const polyMesh& mesh,
3480  const labelList& meshEdges
3481 )
3482 {
3483  const globalIndex globalEdges(meshEdges.size());
3484 
3485  labelList myEdges
3486  (
3487  identity(globalEdges.localSize(), globalEdges.localStart())
3488  );
3489 
3490  syncTools::syncEdgeList
3491  (
3492  mesh,
3493  meshEdges,
3494  myEdges,
3495  minEqOp<label>(),
3496  labelMax
3497  );
3498 
3499 
3500  bitSet isMasterEdge(meshEdges.size());
3501  forAll(meshEdges, edgei)
3502  {
3503  if (myEdges[edgei] == globalEdges.toGlobal(edgei))
3504  {
3505  isMasterEdge.set(edgei);
3506  }
3507  }
3508 
3509  return isMasterEdge;
3510 }
3511 
3512 
3513 void Foam::meshRefinement::printMeshInfo(const bool debug, const string& msg)
3514 const
3515 {
3516  const globalMeshData& pData = mesh_.globalData();
3517 
3518  if (debug)
3519  {
3520  Pout<< msg.c_str()
3521  << " : cells(local):" << mesh_.nCells()
3522  << " faces(local):" << mesh_.nFaces()
3523  << " points(local):" << mesh_.nPoints()
3524  << endl;
3525  }
3526 
3527  {
3528  bitSet isMasterFace(syncTools::getMasterFaces(mesh_));
3529  label nMasterFaces = isMasterFace.count();
3530 
3531  bitSet isMeshMasterPoint(syncTools::getMasterPoints(mesh_));
3532  label nMasterPoints = isMeshMasterPoint.count();
3533 
3534  Info<< msg.c_str()
3535  << " : cells:" << pData.nTotalCells()
3536  << " faces:" << returnReduce(nMasterFaces, sumOp<label>())
3537  << " points:" << returnReduce(nMasterPoints, sumOp<label>())
3538  << endl;
3539  }
3540 
3541 
3542  //if (debug)
3543  {
3544  const labelList& cellLevel = meshCutter_.cellLevel();
3545 
3546  labelList nCells(gMax(cellLevel)+1, Zero);
3547 
3548  forAll(cellLevel, celli)
3549  {
3550  nCells[cellLevel[celli]]++;
3551  }
3552 
3553  Pstream::listCombineGather(nCells, plusEqOp<label>());
3554 
3556  if (Pstream::master())
3557  {
3558  Info<< "Cells per refinement level:" << endl;
3559  forAll(nCells, leveli)
3560  {
3561  Info<< " " << leveli << '\t' << nCells[leveli]
3562  << endl;
3563  }
3564  }
3565  }
3566 }
3567 
3568 
3569 Foam::word Foam::meshRefinement::timeName() const
3570 {
3571  if (overwrite_ && mesh_.time().timeIndex() == 0)
3572  {
3573  return oldInstance_;
3574  }
3575 
3576  return mesh_.time().timeName();
3577 }
3578 
3579 
3580 void Foam::meshRefinement::dumpRefinementLevel() const
3581 {
3582  // Note: use time().timeName(), not meshRefinement::timeName()
3583  // so as to dump the fields to 0, not to constant.
3584  {
3585  volScalarField volRefLevel
3586  (
3587  IOobject
3588  (
3589  "cellLevel",
3590  mesh_.time().timeName(),
3591  mesh_,
3592  IOobject::NO_READ,
3593  IOobject::AUTO_WRITE,
3594  false
3595  ),
3596  mesh_,
3597  dimensionedScalar(dimless, Zero),
3598  zeroGradientFvPatchScalarField::typeName
3599  );
3600 
3601  const labelList& cellLevel = meshCutter_.cellLevel();
3602 
3603  forAll(volRefLevel, celli)
3604  {
3605  volRefLevel[celli] = cellLevel[celli];
3606  }
3607 
3608  volRefLevel.write();
3609  }
3610 
3611  // Dump pointLevel
3612  {
3613  const pointMesh& pMesh = pointMesh::New(mesh_);
3614 
3615  pointScalarField pointRefLevel
3616  (
3617  IOobject
3618  (
3619  "pointLevel",
3620  mesh_.time().timeName(),
3621  mesh_,
3622  IOobject::NO_READ,
3623  IOobject::NO_WRITE,
3624  false
3625  ),
3626  pMesh,
3627  dimensionedScalar(dimless, Zero)
3628  );
3629 
3630  const labelList& pointLevel = meshCutter_.pointLevel();
3631 
3632  forAll(pointRefLevel, pointi)
3633  {
3634  pointRefLevel[pointi] = pointLevel[pointi];
3635  }
3636 
3637  pointRefLevel.write();
3638  }
3639 }
3640 
3641 
3642 void Foam::meshRefinement::dumpIntersections(const fileName& prefix) const
3643 {
3644  {
3645  OFstream str(prefix + "_edges.obj");
3646  label verti = 0;
3647  Pout<< "meshRefinement::dumpIntersections :"
3648  << " Writing cellcentre-cellcentre intersections to file "
3649  << str.name() << endl;
3650 
3651 
3652  // Redo all intersections
3653  // ~~~~~~~~~~~~~~~~~~~~~~
3654 
3655  // Get boundary face centre and level. Coupled aware.
3656  labelList neiLevel(mesh_.nBoundaryFaces());
3657  pointField neiCc(mesh_.nBoundaryFaces());
3658  calcNeighbourData(neiLevel, neiCc);
3659 
3660  labelList intersectionFaces(intersectedFaces());
3661 
3662  // Collect segments we want to test for
3663  pointField start(intersectionFaces.size());
3664  pointField end(intersectionFaces.size());
3665  {
3666  labelList minLevel;
3667  calcCellCellRays
3668  (
3669  neiCc,
3670  labelList(neiCc.size(), -1),
3671  intersectionFaces,
3672  start,
3673  end,
3674  minLevel
3675  );
3676  }
3677 
3678 
3679  // Do tests in one go
3680  labelList surfaceHit;
3681  List<pointIndexHit> surfaceHitInfo;
3682  surfaces_.findAnyIntersection
3683  (
3684  start,
3685  end,
3686  surfaceHit,
3687  surfaceHitInfo
3688  );
3689 
3690  forAll(intersectionFaces, i)
3691  {
3692  if (surfaceHit[i] != -1)
3693  {
3694  meshTools::writeOBJ(str, start[i]);
3695  verti++;
3696  meshTools::writeOBJ(str, surfaceHitInfo[i].point());
3697  verti++;
3698  meshTools::writeOBJ(str, end[i]);
3699  verti++;
3700  str << "l " << verti-2 << ' ' << verti-1 << nl
3701  << "l " << verti-1 << ' ' << verti << nl;
3702  }
3703  }
3704  }
3705 
3706  Pout<< endl;
3707 }
3708 
3709 
3710 void Foam::meshRefinement::write
3711 (
3712  const debugType debugFlags,
3713  const writeType writeFlags,
3714  const fileName& prefix
3715 ) const
3716 {
3717  if (writeFlags & WRITEMESH)
3718  {
3719  write();
3720  }
3721 
3722  if (writeFlags && !(writeFlags & NOWRITEREFINEMENT))
3723  {
3724  meshCutter_.write();
3725 
3726  // Force calculation before writing
3727  (void)surfaceIndex();
3728  surfaceIndex_.write();
3729  }
3730 
3731  if (writeFlags & WRITELEVELS)
3732  {
3733  dumpRefinementLevel();
3734  }
3735 
3736  if ((debugFlags & OBJINTERSECTIONS) && prefix.size())
3737  {
3738  dumpIntersections(prefix);
3739  }
3740 }
3741 
3742 
3743 void Foam::meshRefinement::removeFiles(const polyMesh& mesh)
3744 {
3745  IOobject io
3746  (
3747  "dummy",
3748  mesh.facesInstance(),
3749  mesh.meshSubDir,
3750  mesh
3751  );
3752  fileName setsDir(io.path());
3753 
3754  if (topoSet::debug) DebugVar(setsDir);
3755 
3756  // Remove local files
3757  if (exists(setsDir/"surfaceIndex"))
3758  {
3759  rm(setsDir/"surfaceIndex");
3760  }
3761 
3762  // Remove other files
3763  hexRef8::removeFiles(mesh);
3764 }
3765 
3766 
3767 Foam::meshRefinement::writeType Foam::meshRefinement::writeLevel()
3768 {
3769  return writeLevel_;
3770 }
3771 
3772 
3773 void Foam::meshRefinement::writeLevel(const writeType flags)
3774 {
3775  writeLevel_ = flags;
3776 }
3777 
3778 
3779 //Foam::meshRefinement::outputType Foam::meshRefinement::outputLevel()
3780 //{
3781 // return outputLevel_;
3782 //}
3783 //
3784 //
3785 //void Foam::meshRefinement::outputLevel(const outputType flags)
3786 //{
3787 // outputLevel_ = flags;
3788 //}
3789 
3790 
3791 const Foam::dictionary& Foam::meshRefinement::subDict
3792 (
3793  const dictionary& dict,
3794  const word& keyword,
3795  const bool noExit,
3796  enum keyType::option matchOpt
3797 )
3798 {
3799  const auto finder(dict.csearch(keyword, matchOpt));
3800 
3801  if (!finder.good())
3802  {
3803  auto& err = FatalIOErrorInFunction(dict);
3804 
3805  err << "Entry '" << keyword << "' not found in dictionary "
3806  << dict.relativeName() << nl;
3807 
3808  if (noExit)
3809  {
3810  return dictionary::null;
3811  }
3812  else
3813  {
3814  err << exit(FatalIOError);
3815  }
3816  }
3817 
3818  return finder.dict();
3819 }
3820 
3821 
3822 Foam::ITstream& Foam::meshRefinement::lookup
3823 (
3824  const dictionary& dict,
3825  const word& keyword,
3826  const bool noExit,
3827  enum keyType::option matchOpt
3828 )
3829 {
3830  const entry* eptr = dict.findEntry(keyword, matchOpt);
3831 
3832  if (!eptr)
3833  {
3834  auto& err = FatalIOErrorInFunction(dict);
3835 
3836  err << "Entry '" << keyword << "' not found in dictionary "
3837  << dict.relativeName() << nl;
3838 
3839  if (noExit)
3840  {
3841  // Fake entry
3842  return dict.first()->stream();
3843  }
3844  else
3845  {
3846  err << exit(FatalIOError);
3847  }
3848  }
3849 
3850  return eptr->stream();
3851 }
3852 
3853 
3854 // ************************************************************************* //
Foam::polyBoundaryMesh::findPatchID
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
Definition: polyBoundaryMesh.C:737
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:584
Foam::regionSplit
This class separates the mesh into distinct unconnected regions, each of which is then given a label ...
Definition: regionSplit.H:136
Foam::motionSmootherAlgo::checkMesh
static bool checkMesh(const bool report, const polyMesh &mesh, const dictionary &dict, labelHashSet &wrongFaces, const bool dryRun=false)
Check mesh with mesh settings in dict. Collects incorrect faces.
Definition: motionSmootherAlgoCheck.C:497
Foam::meshSearch
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
Definition: meshSearch.H:56
Foam::fvMesh::clearOut
void clearOut()
Clear all geometry and addressing.
Definition: fvMesh.C:232
Foam::Pair::first
const T & first() const noexcept
Access the first element.
Definition: Pair.H:136
Foam::meshRefinement::printMeshInfo
void printMeshInfo(const bool, const string &) const
Print some mesh stats.
Definition: meshRefinement.C:3514
Foam::meshRefinement::addFaceZone
label addFaceZone(const word &fzName, const word &masterPatch, const word &slavePatch, const surfaceZonesInfo::faceZoneType &fzType)
Add/lookup faceZone and update information. Return index of.
Definition: meshRefinement.C:2620
dict
dictionary dict
Definition: searchingEngine.H:11
Foam::UList::size
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:118
Foam::removeCells::getExposedFaces
labelList getExposedFaces(const bitSet &removedCell) const
Get labels of faces exposed after cells removal.
Definition: removeCells.C:77
Foam::UList::unset
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type unset(const label i)
Unset the bool entry at specified position, always false for out-of-range access. ...
Definition: UList.H:741
Foam::bitSet::count
unsigned int count(const bool on=true) const
Count number of bits set.
Definition: bitSetI.H:493
Foam::PackedList::size
label size() const noexcept
Number of entries.
Definition: PackedListI.H:370
Foam::meshRefinement::subDict
static const dictionary & subDict(const dictionary &dict, const word &keyword, const bool noExit, enum keyType::option matchOpt=keyType::REGEX)
Wrapper around dictionary::subDict which does not exit.
Definition: meshRefinement.C:3793
Foam::ZoneMesh::clearAddressing
void clearAddressing()
Clear addressing.
Definition: ZoneMesh.C:703
Foam::mapPolyMesh::reversePointMap
const labelList & reversePointMap() const
Reverse point map.
Definition: mapPolyMesh.H:579
Foam::eqOp
Definition: ops.H:67
Foam::bitSet::set
void set(const bitSet &bitset)
Set specified bits from another bitset.
Definition: bitSetI.H:583
Foam::fileName
A class for handling file names.
Definition: fileName.H:71
Foam::faceSet
A list of face labels.
Definition: faceSet.H:47
Foam::polyTopoChange::faces
const DynamicList< face > & faces() const
Definition: polyTopoChange.H:555
Foam::meshRefinement::addPointZone
label addPointZone(const word &name)
Add pointZone if does not exist. Return index of zone.
Definition: meshRefinement.C:2672
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:125
Foam::pointBoundaryMesh
Foam::pointBoundaryMesh.
Definition: pointBoundaryMesh.H:52
Foam::polyMesh::facesInstance
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:853
Foam::dictionary::transfer
void transfer(dictionary &dict)
Transfer the contents of the argument and annul the argument.
Definition: dictionary.C:865
Foam::face
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:68
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::syncTools::getMasterFaces
static bitSet getMasterFaces(const polyMesh &mesh)
Get per face whether it is uncoupled or a master of a coupled set of faces.
Definition: syncTools.C:119
Foam::meshRefinement::countHits
label countHits() const
Count number of intersections (local)
Definition: meshRefinement.C:1756
Foam::List::transfer
void transfer(List< T > &list)
Transfer the contents of the argument List into this list and annul the argument list.
Definition: List.C:439
Foam::primitiveMesh::nPoints
label nPoints() const noexcept
Number of mesh points.
Definition: primitiveMeshI.H:30
Foam::FatalError
error FatalError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL ERROR&#39; header text and sta...
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:120
Foam::meshRefinement::removeFiles
static void removeFiles(const polyMesh &)
Helper: remove all relevant files from mesh instance.
Definition: meshRefinement.C:3744
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:578
Foam::List::append
void append(const T &val)
Append an element at the end of the list.
Definition: List.H:491
Foam::sortedOrder
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
Foam::polyMesh::tetBasePtIs
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:893
Foam::meshRefinement::getFaceZoneInfo
bool getFaceZoneInfo(const word &fzName, label &masterPatchID, label &slavePatchID, surfaceZonesInfo::faceZoneType &fzType) const
Lookup faceZone information. Return false if no information.
Definition: meshRefinement.C:2644
Foam::meshRefinement::findRegions
static label findRegions(const polyMesh &, const vector &perturbVec, const pointField &locationsInMesh, const pointField &locationsOutsideMesh, const label nRegions, labelList &cellRegion, const boolList &blockedFace, const bool exitIfLeakPath, const refPtr< coordSetWriter > &leakPathFormatter)
Find regions points are in.
Definition: meshRefinement.C:2980
Foam::IOobject::name
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:150
Foam::hexRef8::removeFiles
static void removeFiles(const polyMesh &)
Helper: remove all relevant files from mesh instance.
Definition: hexRef8.C:5711
Foam::globalMeshData::globalEdgeSlavesMap
const mapDistribute & globalEdgeSlavesMap() const
Definition: globalMeshData.C:2207
Foam::MeshObject< polyMesh, UpdateableMeshObject, pointMesh >::New
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
Definition: MeshObject.C:41
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:40
Foam::meshRefinement::findRegion
static label findRegion(const polyMesh &, const labelList &cellRegion, const vector &perturbVec, const point &p)
Find region point is in. Uses optional perturbation to re-test.
Definition: meshRefinement.C:2793
Foam::globalMeshData
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
Definition: globalMeshData.H:103
Foam::polyPatch::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:309
Foam::renumber
IntListType renumber(const labelUList &oldToNew, const IntListType &input)
Renumber the values (not the indices) of a list.
Definition: ListOpsTemplates.C:26
Foam::OFstream
Output to file stream, using an OSstream.
Definition: OFstream.H:49
Foam::nl
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:49
Foam::meshRefinement::doRemoveCells
autoPtr< mapPolyMesh > doRemoveCells(const labelList &cellsToRemove, const labelList &exposedFaces, const labelList &exposedPatchIDs, removeCells &cellRemover)
Remove cells. Put exposedFaces into exposedPatchIDs.
Definition: meshRefinement.C:1117
Foam::meshRefinement::lookup
static ITstream & lookup(const dictionary &dict, const word &keyword, const bool noExit, enum keyType::option matchOpt=keyType::REGEX)
Wrapper around dictionary::lookup which does not exit.
Definition: meshRefinement.C:3824
Foam::meshRefinement::splitFacesUndo
label splitFacesUndo(const labelList &splitFaces, const labelPairList &splits, const dictionary &motionDict, labelList &duplicateFace, List< labelPair > &baffles)
Split faces along diagonal. Maintain mesh quality. Return.
Definition: meshRefinement.C:1300
Foam::syncTools::syncBoundaryFacePositions
static void syncBoundaryFacePositions(const polyMesh &mesh, UList< point > &positions, const CombineOp &cop)
Synchronize locations on boundary faces only.
Definition: syncTools.H:412
Foam::PatchTools::matchEdges
static void matchEdges(const PrimitivePatch< FaceList1, PointField1 > &p1, const PrimitivePatch< FaceList2, PointField2 > &p2, labelList &p1EdgeLabels, labelList &p2EdgeLabels, bitSet &sameOrientation)
Find corresponding edges on patches sharing the same points.
Definition: PatchToolsMatch.C:70
Foam::UList::data
T * data() noexcept
Return pointer to the underlying array serving as data storage.
Definition: UListI.H:230
Foam::removeCells
Given list of cells to remove, insert all the topology changes.
Definition: removeCells.H:59
coupled
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
rndGen
Random rndGen
Definition: createFields.H:23
Foam::polyMesh::meshSubDir
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:402
Foam::decompositionMethod::decompose
virtual labelList decompose(const pointField &points, const scalarField &pointWeights) const
Return the wanted processor number for every coordinate.
Definition: decompositionMethod.C:1312
Foam::polyBoundaryMesh::reorder
void reorder(const labelUList &oldToNew, const bool validBoundary)
Reorders patches. Ordering does not have to be done in.
Definition: polyBoundaryMesh.C:1196
Foam::dimensionedVector
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Definition: dimensionedVector.H:46
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:487
Foam::decompositionMethod::setConstraints
void setConstraints(const polyMesh &mesh, boolList &blockedFace, PtrList< labelList > &specifiedProcessorFaces, labelList &specifiedProcessor, List< labelPair > &explicitConnections) const
Helper: extract constraints:
Definition: decompositionMethod.C:1204
Foam::UPstream::parRun
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:639
Foam::UPtrList::reorder
void reorder(const labelUList &oldToNew, const bool check=false)
Reorder elements. Reordering must be unique (ie, shuffle).
Definition: UPtrList.C:62
Foam::autoPtr::clear
void clear() noexcept
Same as reset(nullptr)
Definition: autoPtr.H:255
Foam::mapDistributeBase::constructSize
label constructSize() const noexcept
Constructed data size.
Definition: mapDistributeBase.H:568
Foam::meshRefinement::intersectedFaces
labelList intersectedFaces() const
Get faces with intersection.
Definition: meshRefinement.C:2026
Foam::meshRefinement::writeLevel
static writeType writeLevel()
Get/set write level.
Definition: meshRefinement.C:3768
Foam::syncTools::syncFaceList
static void syncFaceList(const polyMesh &mesh, UList< T > &faceValues, const CombineOp &cop)
Synchronize values on all mesh faces.
Definition: syncTools.H:432
Foam::mapDistributePolyMesh
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Definition: mapDistributePolyMesh.H:62
Foam::meshTools::writeOBJ
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of a point.
Definition: meshTools.C:196
Foam::faceZone::whichFace
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
Definition: faceZone.C:365
Foam::meshRefinement::getMasterEdges
static bitSet getMasterEdges(const polyMesh &mesh, const labelList &meshEdges)
Determine master edge for subset of edges. If coupled.
Definition: meshRefinement.C:3479
Foam::searchableSurface
Base class of (analytical or triangulated) surface. Encapsulates all the search routines. WIP.
Definition: searchableSurface.H:65
Foam::IOobjectOption::NO_WRITE
Ignore writing from objectRegistry::writeObject()
Definition: IOobjectOption.H:81
Foam::globalMeshData::globalEdgeSlaves
const labelListList & globalEdgeSlaves() const
Definition: globalMeshData.C:2176
Foam::normalised
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Definition: quaternionI.H:673
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:182
Foam::UPstream::myProcNo
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:688
Foam::meshRefinement::calculateEdgeWeights
static void calculateEdgeWeights(const polyMesh &mesh, const bitSet &isMasterEdge, const labelList &meshPoints, const edgeList &edges, scalarField &edgeWeights, scalarField &invSumWeight)
Helper: calculate edge weights (1/length)
Definition: meshRefinement.C:2311
Foam::pointList
List< point > pointList
A List of points.
Definition: pointList.H:61
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:361
Foam::primitiveMesh::nFaces
label nFaces() const noexcept
Number of mesh faces.
Definition: primitiveMeshI.H:83
Foam::meshRefinement::intersectedPoints
labelList intersectedPoints() const
Get points on surfaces with intersection and boundary faces.
Definition: meshRefinement.C:2054
Foam::polyTopoChange::changeMesh
autoPtr< mapPolyMesh > changeMesh(polyMesh &mesh, const labelUList &patchMap, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
Definition: polyTopoChange.C:2999
Foam::returnReduce
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Perform reduction on a copy, using specified binary operation.
Definition: PstreamReduceOps.H:377
Foam::polyMesh::cellDecomposition
cellDecomposition
Enumeration defining the decomposition of the cell for.
Definition: polyMesh.H:102
Foam::shortestPathSet
Finds shortest path (in terms of cell centres) to walk on mesh from any point in insidePoints to any ...
Definition: shortestPathSet.H:111
Foam::refPtr
A class for managing references or pointers (no reference counting)
Definition: HashPtrTable.H:49
Foam::mapPolyMesh::hasMotionPoints
bool hasMotionPoints() const
Has valid preMotionPoints?
Definition: mapPolyMesh.H:764
Foam::mapPolyMesh
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:157
Foam::labelUList
UList< label > labelUList
A UList of labels.
Definition: UList.H:80
Foam::faceMap
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
Definition: blockMeshMergeTopological.C:87
Foam::Pstream::broadcast
static void broadcast(Type &value, const label comm=UPstream::worldComm)
Broadcast content (contiguous or non-contiguous) to all processes in communicator.
Definition: PstreamBroadcast.C:28
Foam::localPointRegion::findDuplicateFacePairs
static labelPairList findDuplicateFacePairs(const polyMesh &)
Helper routine to find all baffles (two boundary faces.
Definition: localPointRegion.C:618
Foam::IOobject::path
fileName path() const
The complete path.
Definition: IOobject.C:465
Foam::refinementSurfaces
Container for data on surfaces used for surface-driven refinement. Contains all the data about the le...
Definition: refinementSurfaces.H:59
Foam::meshRefinement::addPatch
static label addPatch(fvMesh &, const word &name, const dictionary &)
Helper:add patch to mesh. Update all registered fields.
Definition: meshRefinement.C:2465
Foam::polyMesh::points
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1066
Foam::pointMesh
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:45
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:413
Foam::globalMeshData::syncData
static void syncData(List< Type > &elems, const labelListList &slaves, const labelListList &transformedSlaves, const mapDistribute &slavesMap, const globalIndexAndTransform &, const CombineOp &cop, const TransformOp &top)
Helper: synchronise data with transforms.
Definition: globalMeshDataTemplates.C:30
Foam::UList::fcIndex
label fcIndex(const label i) const noexcept
The forward circular index. The next index in the list which returns to the first at the end of the l...
Definition: UListI.H:53
Foam::refPtr::constCast
T & constCast() const
Return non-const reference to the object or to the contents of a (non-null) managed pointer...
Definition: refPtr.H:277
Foam::autoPtr::reset
void reset(T *p=nullptr) noexcept
Delete managed object and set to new given pointer.
Definition: autoPtrI.H:37
Foam::meshRefinement::writeTypeNames
static const Enum< writeType > writeTypeNames
Definition: meshRefinement.H:122
Foam::labelHashSet
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
Definition: HashSet.H:85
Foam::refinementFeatures
Encapsulates queries for features.
Definition: refinementFeatures.H:49
Foam::mapDistributePolyMesh::distributeFaceData
void distributeFaceData(List< T > &values) const
Distribute list of face data.
Definition: mapDistributePolyMesh.H:332
Foam::meshRefinement::appendPatch
static label appendPatch(fvMesh &, const label insertPatchi, const word &, const dictionary &)
Helper:append patch to end of mesh.
Definition: meshRefinement.C:2367
Foam::Pstream::allGatherList
static void allGatherList(List< T > &values, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Gather data, but keep individual values separate. Uses linear/tree communication. ...
Definition: PstreamGatherList.C:340
Foam::meshRefinement::doSplitFaces
void doSplitFaces(const labelList &splitFaces, const labelPairList &splits, polyTopoChange &meshMod) const
Split faces into two.
Definition: meshRefinement.C:1181
Foam::meshRefinement::distribute
void distribute(const mapDistributePolyMesh &)
Update local numbering for mesh redistribution.
Definition: meshRefinement.C:3205
Foam::meshRefinement::updateList
static void updateList(const labelList &newToOld, const T &nullValue, List< T > &elems)
Helper: reorder list according to map.
Definition: meshRefinementTemplates.C:31
timeName
word timeName
Definition: getTimeIndex.H:3
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: POSIX.C:752
Foam::vtk::write
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Definition: foamVtkOutputTemplates.C:29
Foam::PrimitivePatch
A list of faces which address into the list of points.
Definition: PrimitivePatch.H:75
Foam::BitOps::count
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of &#39;true&#39; entries.
Definition: BitOps.H:73
Foam::globalIndex
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
Foam::SubList
A List obtained as a section of another List.
Definition: SubList.H:50
Foam::UPstream::nProcs
static label nProcs(const label communicator=worldComm)
Number of ranks in parallel run (for given communicator) is 1 for serial run.
Definition: UPstream.H:656
Foam::fvMeshDistribute
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
Definition: fvMeshDistribute.H:66
Foam::returnReduceAnd
bool returnReduceAnd(const bool value, const label comm=UPstream::worldComm)
Perform logical (and) MPI Allreduce on a copy. Uses UPstream::reduceAnd.
Definition: PstreamReduceOps.H:396
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:38
Foam::Pstream::gatherList
static void gatherList(const List< commsStruct > &comms, List< T > &values, const int tag, const label comm)
Gather data, but keep individual values separate. Uses the specified communication schedule...
Definition: PstreamGatherList.C:41
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
Foam::polyMesh::pointsInstance
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:847
io
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:289
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:52
Foam::polyTopoChange::addFace
label addFace(const face &f, const label own, const label nei, const label masterPointID, const label masterEdgeID, const label masterFaceID, const bool flipFaceFlux, const label patchID, const label zoneID, const bool zoneFlip)
Add face to cells. Return new face label.
Definition: polyTopoChange.C:2766
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::ZoneMesh::findZoneID
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
Definition: ZoneMesh.C:513
Foam::edge
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:59
Foam::identity
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
Definition: labelList.C:31
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:63
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:46
Foam::syncTools::getMasterPoints
static bitSet getMasterPoints(const polyMesh &mesh)
Get per point whether it is uncoupled or a master of a coupled set of points.
Definition: syncTools.C:61
Foam::polyTopoChange::points
const DynamicList< point > & points() const
Points. Shrunk after constructing mesh (or calling of compact())
Definition: polyTopoChange.H:550
Foam::polyBoundaryMesh::names
wordList names() const
Return a list of patch names.
Definition: polyBoundaryMesh.C:521
Foam::UPtrList::size
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:99
Foam::scalarList
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:61
Foam::calculatedFvPatchField
This boundary condition is not designed to be evaluated; it is assmued that the value is assigned via...
Definition: calculatedFvPatchField.H:62
Foam::coordSetWriter
Base class for writing coordSet(s) and tracks with fields.
Definition: coordSetWriter.H:81
Foam::meshRefinement::checkCoupledFaceZones
static void checkCoupledFaceZones(const polyMesh &)
Helper function: check that face zones are synced.
Definition: meshRefinement.C:2222
Foam::dictionary::null
static const dictionary null
An empty dictionary, which is also the parent for all dictionaries.
Definition: dictionary.H:465
Foam::searchableSurfaces
Container for searchableSurfaces. The collection is specified as a dictionary. For example...
Definition: searchableSurfaces.H:89
Foam::tmp::New
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
Definition: tmp.H:203
Foam::exists
bool exists(const fileName &name, const bool checkGzip=true, const bool followLink=true)
Does the name exist (as DIRECTORY or FILE) in the file system?
Definition: POSIX.C:788
Foam::meshRefinement::surfaceIndex
const labelList & surfaceIndex() const
Per start-end edge the index of the surface hit.
Definition: meshRefinement.C:1733
Foam::PrimitivePatch::edgeFaces
const labelListList & edgeFaces() const
Return edge-face addressing.
Definition: PrimitivePatch.C:255
Foam::decompositionMethod
Abstract base class for domain decomposition.
Definition: decompositionMethod.H:47
Foam::primitiveMesh::nInternalFaces
label nInternalFaces() const noexcept
Number of internal faces.
Definition: primitiveMeshI.H:71
Foam::shellSurfaces
Encapsulates queries for volume refinement (&#39;refine all cells within shell&#39;).
Definition: shellSurfaces.H:53
Foam::meshRefinement::balance
autoPtr< mapDistributePolyMesh > balance(const bool keepZoneFaces, const bool keepBaffles, const scalarField &cellWeights, decompositionMethod &decomposer, fvMeshDistribute &distributor)
Redecompose according to cell count.
Definition: meshRefinement.C:1777
Foam::vector
Vector< scalar > vector
Definition: vector.H:57
Foam::polyMesh::faces
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1091
Foam::distribution
Accumulating histogram of values. Specified bin resolution automatic generation of bins...
Definition: distribution.H:57
Foam::polyBoundaryMesh::groupPatchIDs
const HashTable< labelList > & groupPatchIDs() const
The patch indices per patch group.
Definition: polyBoundaryMesh.C:447
Foam::removeCells::updateMesh
void updateMesh(const mapPolyMesh &)
Force recalculation of locally stored data on topological change.
Definition: removeCells.H:149
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:26
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:139
Foam::PtrList::set
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
Definition: PtrList.H:163
Foam::syncTools::syncPointList
static void syncPointList(const polyMesh &mesh, List< T > &pointValues, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh points.
Definition: syncToolsTemplates.C:695
Foam::bitSet::test
bool test(const label pos) const
Test value at specified position, never auto-vivify entries.
Definition: bitSetI.H:514
Foam::Random
Random number generator.
Definition: Random.H:55
Foam::removeCells::setRefinement
void setRefinement(const bitSet &removedCell, const labelUList &facesToExpose, const labelUList &patchIDs, polyTopoChange &) const
Play commands into polyTopoChange to remove cells.
Definition: removeCells.C:181
Foam::meshRefinement::timeName
word timeName() const
Replacement for Time::timeName() that returns oldInstance (if overwrite_)
Definition: meshRefinement.C:3570
Foam::syncTools::swapFaceList
static void swapFaceList(const polyMesh &mesh, UList< T > &faceValues)
Swap coupled face values. Uses eqOp.
Definition: syncTools.H:522
Foam::polyBoundaryMesh
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO...
Definition: polyBoundaryMesh.H:58
Foam::polyTopoChange::setCapacity
void setCapacity(const label nPoints, const label nFaces, const label nCells)
Explicitly pre-size the dynamic storage for expected mesh size for if construct-without-mesh.
Definition: polyTopoChange.C:2468
Foam::vector
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Return time name of given scalar time formatted with the given precision.
Definition: Time.C:760
Foam::meshRefinement::checkData
void checkData()
Debugging: check that all faces still obey start()>end()
Definition: meshRefinement.C:902
Foam::ensightOutput::debug
int debug
Static debugging option.
Foam::labelPair
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:50
writer
vtk::internalMeshWriter writer(topoMesh, topoCells, vtk::formatType::INLINE_ASCII, runTime.path()/"blockTopology")
Foam::gMax
Type gMax(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:616
stdFoam::end
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
Definition: stdFoam.H:193
Foam::polyMesh::faceZones
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
Definition: polyMesh.H:646
Foam::vtk::fileWriter::open
virtual bool open(const fileName &file, bool parallel=UPstream::parRun())
Open file for writing (creates parent directory).
Definition: foamVtkFileWriter.C:311
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::data
Database for solution data, solver performance and other reduced data.
Definition: data.H:51
Foam::fileName::clean
static bool clean(std::string &str)
Cleanup filename string, possibly applies other transformations such as changing the path separator e...
Definition: fileName.C:192
Foam::meshRefinement::addMeshedPatch
label addMeshedPatch(const word &name, const dictionary &)
Add patch originating from meshing. Update meshedPatches_.
Definition: meshRefinement.C:2543
Foam::meshRefinement::meshedPatches
labelList meshedPatches() const
Get patchIDs for patches added in addMeshedPatch.
Definition: meshRefinement.C:2593
Foam::fvMeshDistribute::distribute
autoPtr< mapDistributePolyMesh > distribute(const labelList &dist)
Send cells to neighbours according to distribution.
Definition: fvMeshDistribute.C:1891
f
labelList f(nPoints)
Foam::PtrList::setSize
void setSize(const label newLen)
Same as resize()
Definition: PtrList.H:183
Foam::globalMeshData::nTotalCells
label nTotalCells() const noexcept
Return total number of cells in decomposed mesh.
Definition: globalMeshData.H:464
Foam::OFstream::name
virtual const fileName & name() const
Read/write access to the name of the stream.
Definition: OSstream.H:128
Foam::meshRefinement::makePatch
static autoPtr< indirectPrimitivePatch > makePatch(const polyMesh &, const labelList &)
Create patch from set of patches.
Definition: meshRefinement.C:2127
Foam::dictionary::relativeName
fileName relativeName(const bool caseTag=false) const
The dictionary name relative to the case.
Definition: dictionary.C:179
Foam::mapPolyMesh::faceMap
const labelList & faceMap() const
Old face map.
Definition: mapPolyMesh.H:500
Foam::meshRefinement::updateMesh
void updateMesh(const mapPolyMesh &, const labelList &changedFaces)
Update for external change to mesh. changedFaces are in new mesh.
Definition: meshRefinement.C:3268
Foam::meshRefinement::updateIntersections
void updateIntersections(const labelList &changedFaces)
Find any intersection of surface. Store in surfaceIndex_.
Definition: meshRefinement.C:349
findCellMode
static const Foam::polyMesh::cellDecomposition findCellMode(Foam::polyMesh::FACE_DIAG_TRIS)
Foam::meshRefinement::splitMeshRegions
autoPtr< mapPolyMesh > splitMeshRegions(const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const pointField &locationsInMesh, const pointField &locationsOutsideMesh, const bool exitIfLeakPath, const refPtr< coordSetWriter > &leakPathFormatter)
Split mesh. Keep part containing point. Return empty map if.
Definition: meshRefinement.C:3097
fld
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputLagrangian.H:23
Foam::meshRefinement
Helper class which maintains intersections of (changing) mesh with (static) surfaces.
Definition: meshRefinement.H:82
Foam::coordSet::coordFormat::DISTANCE
Use additional distance field for (scalar) axis.
split
static bool split(const std::string &line, std::string &key, std::string &val)
Definition: cpuInfo.C:32
Foam::globalMeshData::coupledPatch
const indirectPrimitivePatch & coupledPatch() const
Return patch of all coupled faces.
Definition: globalMeshData.C:2008
Foam::meshRefinement::setInstance
void setInstance(const fileName &)
Set instance of all local IOobjects.
Definition: meshRefinement.C:1109
Foam::globalPoints
Calculates points shared by more than two processor patches or cyclic patches.
Definition: globalPoints.H:98
Foam::polyTopoChange::modifyFace
void modifyFace(const face &f, const label facei, const label own, const label nei, const bool flipFaceFlux, const label patchID, const label zoneID, const bool zoneFlip)
Modify vertices or cell of face.
Definition: polyTopoChange.C:2830
Foam::inplaceRenumber
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
Definition: ListOpsTemplates.C:50
Foam::mapDistribute
Class containing processor-to-processor mapping information.
Definition: mapDistribute.H:156
Foam::pointZone
A subset of mesh points.
Definition: pointZone.H:61
Foam::polyPatch::constraintType
static bool constraintType(const word &patchType)
Return true if the given type is a constraint type.
Definition: polyPatch.C:270
Foam::meshRefinement::storeData
void storeData(const labelList &pointsToStore, const labelList &facesToStore, const labelList &cellsToStore)
Signal points/face/cells for which to store data.
Definition: meshRefinement.C:3280
Foam::point
vector point
Point is a vector.
Definition: point.H:37
Foam::bitSet
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition: bitSet.H:59
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:400
Foam::Enum
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: IOstreamOption.H:55
Foam::meshRefinement::testSyncPointList
static void testSyncPointList(const string &msg, const polyMesh &mesh, const List< scalar > &fld)
Definition: meshRefinement.C:809
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
FatalIOErrorInFunction
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:607
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:36
Foam::primitiveMesh::nCells
label nCells() const noexcept
Number of mesh cells.
Definition: primitiveMeshI.H:89
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
Definition: List.H:55
Foam::fvBoundaryMesh
Foam::fvBoundaryMesh.
Definition: fvBoundaryMesh.H:53
Foam::meshRefinement::getMasterPoints
static bitSet getMasterPoints(const polyMesh &mesh, const labelList &meshPoints)
Determine master point for subset of points. If coupled.
Definition: meshRefinement.C:3443
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:79
Foam::UIndirectList
A List with indirect addressing. Like IndirectList but does not store addressing. ...
Definition: faMatrix.H:57
Foam::polyTopoChange
Direct mesh changes based on v1.3 polyTopoChange syntax.
Definition: polyTopoChange.H:95
Foam::vtk::fileWriter::close
void close()
End the file contents and close the file after writing.
Definition: foamVtkFileWriter.C:369
Foam::polyBoundaryMesh::whichPatch
label whichPatch(const label faceIndex) const
Return patch index for a given mesh face index.
Definition: polyBoundaryMesh.H:343
Foam::UPstream::master
static bool master(const label communicator=worldComm)
Am I the master rank.
Definition: UPstream.H:672
Foam::polyPatch::start
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:441
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:59
Foam::zoneIdentifier::name
const word & name() const noexcept
The zone name.
Definition: zoneIdentifier.H:145
Foam::IOobjectOption::NO_READ
Nothing to be read.
Definition: IOobjectOption.H:63
Foam::keyType::option
option
Enumeration for the data type and search/match modes (bitmask)
Definition: keyType.H:79
Foam::IOobjectOption::AUTO_WRITE
Automatically write from objectRegistry::writeObject()
Definition: IOobjectOption.H:84
Foam::functionObject::outputPrefix
static word outputPrefix
Directory prefix.
Definition: functionObject.H:390
Foam::syncTools::swapBoundaryFacePositions
static void swapBoundaryFacePositions(const polyMesh &mesh, UList< point > &positions)
Swap coupled positions. Uses eqOp.
Definition: syncTools.H:503
Foam::topoDistanceData
For use with FaceCellWave. Determines topological distance to starting faces. Templated on passive tr...
Definition: topoDistanceData.H:49
Foam::regionSplit::nRegions
label nRegions() const
Return total number of regions.
Definition: regionSplit.H:328
Foam::dictionary::findEntry
const entry * findEntry(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find an entry (const access) with the given keyword.
Definition: dictionaryI.H:80
Foam::reduce
void reduce(const List< UPstream::commsStruct > &comms, T &value, const BinaryOp &bop, const int tag, const label comm)
Reduce inplace (cf. MPI Allreduce) using specified communication schedule.
Definition: PstreamReduceOps.H:49
Foam::polyMesh::findCell
label findCell(const point &p, const cellDecomposition=CELL_TETS) const
Find cell enclosing this location and return index.
Definition: polyMesh.C:1505
Foam::meshRefinement::debugTypeNames
static const Enum< debugType > debugTypeNames
Definition: meshRefinement.H:100
Foam::faceZone::flipMap
const boolList & flipMap() const noexcept
Return face flip map.
Definition: faceZone.H:325
Foam::localPointRegion::findDuplicateFaces
static labelList findDuplicateFaces(const primitiveMesh &, const labelList &)
Helper routine to find baffles (two boundary faces using the.
Definition: localPointRegion.C:537
Foam::meshRefinement::selectSeparatedCoupledFaces
void selectSeparatedCoupledFaces(boolList &) const
Select coupled faces that are not collocated.
Definition: meshRefinement.C:2699
Foam::treeBoundBox
Standard boundBox with extra functionality for use in octree.
Definition: treeBoundBox.H:90
Foam::ZoneMesh::names
wordList names() const
A list of the zone names.
Definition: ZoneMesh.C:292
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::labelMax
constexpr label labelMax
Definition: label.H:55
Foam::Field< vector >::subField
SubField< vector > subField
Declare type of subField.
Definition: Field.H:128
Foam::meshRefinement::writeType
writeType
Enumeration for what to write. Used as a bit-pattern.
Definition: meshRefinement.H:113
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::vectorField
Field< vector > vectorField
Specialisation of Field<T> for vector.
Definition: primitiveFieldsFwd.H:48
Foam::polyPatch::New
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return pointer to a new patch created on freestore from components.
Definition: polyPatchNew.C:28
Foam::surfaceZonesInfo::faceZoneType
faceZoneType
What to do with faceZone faces.
Definition: surfaceZonesInfo.H:88
Foam::meshRefinement::debugType
debugType
Enumeration for what to debug. Used as a bit-pattern.
Definition: meshRefinement.H:91
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
Definition: HashPtrTable.H:48
DebugVar
#define DebugVar(var)
Report a variable name and value.
Definition: messageStream.H:523
Foam::FixedList::data
T * data() noexcept
Return pointer to the underlying array serving as data storage.
Definition: FixedListI.H:169
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:73
Foam::dictionary::set
entry * set(entry *entryPtr)
Assign a new entry, overwriting any existing entry.
Definition: dictionary.C:777
Foam::Pstream::listCombineGather
static void listCombineGather(const List< commsStruct > &comms, List< T > &values, const CombineOp &cop, const int tag, const label comm)
Definition: PstreamCombineGather.C:301
Foam::fvPatch::New
static autoPtr< fvPatch > New(const polyPatch &, const fvBoundaryMesh &)
Return a pointer to a new patch created on freestore from polyPatch.
Definition: fvPatchNew.C:28
Foam::faceZone
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:60
Foam::syncTools::syncEdgeList
static void syncEdgeList(const polyMesh &mesh, List< T > &edgeValues, const CombineOp &cop, const T &nullValue, const TransformOp &top, const FlipOp &fop)
Synchronize values on all mesh edges.
Definition: syncToolsTemplates.C:774
Foam::meshRefinement::write
bool write() const
Write mesh and all data.
Definition: meshRefinement.C:3407
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:62
Foam::dictionary::csearch
const_searcher csearch(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search dictionary for given keyword.
Definition: dictionarySearch.C:258
Foam::meshRefinement::makeDisplacementField
static tmp< pointVectorField > makeDisplacementField(const pointMesh &pMesh, const labelList &adaptPatchIDs)
Helper function to make a pointVectorField with correct.
Definition: meshRefinement.C:2169
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::tmp
A class for managing temporary objects.
Definition: HashPtrTable.H:50
Foam::UILList::first
T * first()
The first entry in the list.
Definition: UILList.H:544
Foam::autoPtr::New
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
Definition: autoPtr.H:178
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:69
Foam::Pair::second
const T & second() const noexcept
Access the second element.
Definition: Pair.H:146
s
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputSpray.H:25
Foam::meshRefinement::dumpIntersections
void dumpIntersections(const fileName &prefix) const
Debug: Write intersection information to OBJ format.
Definition: meshRefinement.C:3643
Foam::entry::stream
virtual ITstream & stream() const =0
Return token stream, if entry is a primitive entry.
Foam::primitiveMesh::nBoundaryFaces
label nBoundaryFaces() const noexcept
Number of boundary faces (== nFaces - nInternalFaces)
Definition: primitiveMeshI.H:77
Foam::pointMesh::boundary
const pointBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: pointMesh.H:130
Foam::fvMeshDistribute::countCells
static labelList countCells(const labelList &)
Helper function: count cells per processor in wanted distribution.
Definition: fvMeshDistribute.C:1867
Foam::polyTopoChange::removeFace
void removeFace(const label facei, const label mergeFacei)
Remove/merge face.
Definition: polyTopoChange.C:2867
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported...
Definition: IOobject.H:166
Foam::boolList
List< bool > boolList
A List of bools.
Definition: List.H:60
Foam::IndirectList
A List with indirect addressing.
Definition: IndirectList.H:60
Foam::meshRefinement::dumpRefinementLevel
void dumpRefinementLevel() const
Write refinement level as volScalarFields for postprocessing.
Definition: meshRefinement.C:3581
Foam::meshRefinement::countEdgeFaces
labelList countEdgeFaces(const uindirectPrimitivePatch &pp) const
Count number of faces per patch edge. Parallel consistent.
Definition: meshRefinement.C:2715
Foam::Pout
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
meshBb
List< treeBoundBox > meshBb(1, treeBoundBox(coarseMesh.points()).extend(rndGen, 1e-3))
Foam::ITstream
An input stream of tokens.
Definition: ITstream.H:48
Foam::surfaceZonesInfo::addFaceZone
static label addFaceZone(const word &name, const labelList &addressing, const boolList &flipMap, polyMesh &mesh)
Definition: surfaceZonesInfo.C:502
Foam::mapPolyMesh::reverseFaceMap
const labelList & reverseFaceMap() const
Reverse face map.
Definition: mapPolyMesh.H:617
Foam::mapPolyMesh::preMotionPoints
const pointField & preMotionPoints() const
Pre-motion point positions.
Definition: mapPolyMesh.H:756
Foam::meshRefinement::findCell
static label findCell(const polyMesh &, const vector &perturbVec, const point &p)
Find cell point is in. Uses optional perturbation to re-test.
Definition: meshRefinement.C:2772
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:26
forAllConstIters
forAllConstIters(mixture.phases(), phase)
Definition: pEqn.H:28
Foam::entry
A keyword and a list of tokens is an &#39;entry&#39;.
Definition: entry.H:63
Foam::syncTools::swapBoundaryFaceList
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &faceValues)
Swap coupled boundary face values. Uses eqOp.
Definition: syncTools.H:485
Foam::TimePaths::globalPath
fileName globalPath() const
Return global path for the case.
Definition: TimePathsI.H:73
Foam::Pstream::listCombineReduce
static void listCombineReduce(List< T > &values, const CombineOp &cop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
After completion all processors have the same data.
Definition: PstreamCombineGather.C:530
Foam::coordSet::coordFormatNames
static const Enum< coordFormat > coordFormatNames
String representation of coordFormat enum.
Definition: coordSet.H:74
Foam::globalMeshData::globalEdgeTransformedSlaves
const labelListList & globalEdgeTransformedSlaves() const
Definition: globalMeshData.C:2186
Foam::rm
bool rm(const fileName &file)
Remove a file (or its gz equivalent), returning true if successful.
Definition: POSIX.C:1357
pts
const pointField & pts
Definition: searchingEngine_CGAL.H:64
Foam::fvMesh::boundary
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:705
Foam::FatalIOError
IOerror FatalIOError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL IO ERROR&#39; header text and ...
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:157

Copyright © 2011-2018  OpenFOAM | OPENFOAM® is a registered trademark of OpenCFD Ltd.

Generated by   doxygen 1.8.14