createPatch.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 | Copyright (C) 2011-2016 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Application
25  createPatch
26 
27 Description
28  Utility to create patches out of selected boundary faces. Faces come either
29  from existing patches or from a faceSet.
30 
31  More specifically it:
32  - creates new patches (from selected boundary faces). Synchronise faces
33  on coupled patches.
34  - synchronises points on coupled boundaries
35  - remove patches with 0 faces in them
36 
37 \*---------------------------------------------------------------------------*/
38 
39 #include "cyclicPolyPatch.H"
40 #include "syncTools.H"
41 #include "argList.H"
42 #include "polyMesh.H"
43 #include "Time.H"
44 #include "SortableList.H"
45 #include "OFstream.H"
46 #include "meshTools.H"
47 #include "faceSet.H"
48 #include "IOPtrList.H"
49 #include "polyTopoChange.H"
50 #include "polyModifyFace.H"
51 #include "wordReList.H"
52 
53 using namespace Foam;
54 
55 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
56 
57 namespace Foam
58 {
60 }
61 
62 void changePatchID
63 (
64  const polyMesh& mesh,
65  const label faceID,
66  const label patchID,
67  polyTopoChange& meshMod
68 )
69 {
70  const label zoneID = mesh.faceZones().whichZone(faceID);
71 
72  bool zoneFlip = false;
73 
74  if (zoneID >= 0)
75  {
76  const faceZone& fZone = mesh.faceZones()[zoneID];
77 
78  zoneFlip = fZone.flipMap()[fZone.whichFace(faceID)];
79  }
80 
81  meshMod.setAction
82  (
84  (
85  mesh.faces()[faceID], // face
86  faceID, // face ID
87  mesh.faceOwner()[faceID], // owner
88  -1, // neighbour
89  false, // flip flux
90  patchID, // patch ID
91  false, // remove from zone
92  zoneID, // zone ID
93  zoneFlip // zone flip
94  )
95  );
96 }
97 
98 
99 // Filter out the empty patches.
100 void filterPatches(polyMesh& mesh, const HashSet<word>& addedPatchNames)
101 {
102  const polyBoundaryMesh& patches = mesh.boundaryMesh();
103 
104  // Patches to keep
105  DynamicList<polyPatch*> allPatches(patches.size());
106 
107  label nOldPatches = returnReduce(patches.size(), sumOp<label>());
108 
109  // Copy old patches.
110  forAll(patches, patchi)
111  {
112  const polyPatch& pp = patches[patchi];
113 
114  // Note: reduce possible since non-proc patches guaranteed in same order
115  if (!isA<processorPolyPatch>(pp))
116  {
117 
118  // Add if
119  // - non zero size
120  // - or added from the createPatchDict
121  // - or cyclic (since referred to by other cyclic half or
122  // proccyclic)
123 
124  if
125  (
126  addedPatchNames.found(pp.name())
127  || returnReduce(pp.size(), sumOp<label>()) > 0
128  || isA<coupledPolyPatch>(pp)
129  )
130  {
131  allPatches.append
132  (
133  pp.clone
134  (
135  patches,
136  allPatches.size(),
137  pp.size(),
138  pp.start()
139  ).ptr()
140  );
141  }
142  else
143  {
144  Info<< "Removing zero-sized patch " << pp.name()
145  << " type " << pp.type()
146  << " at position " << patchi << endl;
147  }
148  }
149  }
150  // Copy non-empty processor patches
151  forAll(patches, patchi)
152  {
153  const polyPatch& pp = patches[patchi];
154 
155  if (isA<processorPolyPatch>(pp))
156  {
157  if (pp.size())
158  {
159  allPatches.append
160  (
161  pp.clone
162  (
163  patches,
164  allPatches.size(),
165  pp.size(),
166  pp.start()
167  ).ptr()
168  );
169  }
170  else
171  {
172  Info<< "Removing empty processor patch " << pp.name()
173  << " at position " << patchi << endl;
174  }
175  }
176  }
177 
178  label nAllPatches = returnReduce(allPatches.size(), sumOp<label>());
179  if (nAllPatches != nOldPatches)
180  {
181  Info<< "Removing patches." << endl;
182  allPatches.shrink();
183  mesh.removeBoundary();
184  mesh.addPatches(allPatches);
185  }
186  else
187  {
188  Info<< "No patches removed." << endl;
189  forAll(allPatches, i)
190  {
191  delete allPatches[i];
192  }
193  }
194 }
195 
196 
197 // Dump for all patches the current match
198 void dumpCyclicMatch(const fileName& prefix, const polyMesh& mesh)
199 {
200  const polyBoundaryMesh& patches = mesh.boundaryMesh();
201 
202  forAll(patches, patchi)
203  {
204  if
205  (
206  isA<cyclicPolyPatch>(patches[patchi])
207  && refCast<const cyclicPolyPatch>(patches[patchi]).owner()
208  )
209  {
210  const cyclicPolyPatch& cycPatch =
211  refCast<const cyclicPolyPatch>(patches[patchi]);
212 
213  // Dump patches
214  {
215  OFstream str(prefix+cycPatch.name()+".obj");
216  Pout<< "Dumping " << cycPatch.name()
217  << " faces to " << str.name() << endl;
219  (
220  str,
221  cycPatch,
222  cycPatch.points()
223  );
224  }
225 
226  const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
227  {
228  OFstream str(prefix+nbrPatch.name()+".obj");
229  Pout<< "Dumping " << nbrPatch.name()
230  << " faces to " << str.name() << endl;
232  (
233  str,
234  nbrPatch,
235  nbrPatch.points()
236  );
237  }
238 
239 
240  // Lines between corresponding face centres
241  OFstream str(prefix+cycPatch.name()+nbrPatch.name()+"_match.obj");
242  label vertI = 0;
243 
244  Pout<< "Dumping cyclic match as lines between face centres to "
245  << str.name() << endl;
246 
247  forAll(cycPatch, facei)
248  {
249  const point& fc0 = mesh.faceCentres()[cycPatch.start()+facei];
250  meshTools::writeOBJ(str, fc0);
251  vertI++;
252  const point& fc1 = mesh.faceCentres()[nbrPatch.start()+facei];
253  meshTools::writeOBJ(str, fc1);
254  vertI++;
255 
256  str<< "l " << vertI-1 << ' ' << vertI << nl;
257  }
258  }
259  }
260 }
261 
262 
263 void separateList
264 (
265  const vectorField& separation,
266  UList<vector>& field
267 )
268 {
269  if (separation.size() == 1)
270  {
271  // Single value for all.
272 
273  forAll(field, i)
274  {
275  field[i] += separation[0];
276  }
277  }
278  else if (separation.size() == field.size())
279  {
280  forAll(field, i)
281  {
282  field[i] += separation[i];
283  }
284  }
285  else
286  {
288  << "Sizes of field and transformation not equal. field:"
289  << field.size() << " transformation:" << separation.size()
290  << abort(FatalError);
291  }
292 }
293 
294 
295 // Synchronise points on both sides of coupled boundaries.
296 template<class CombineOp>
297 void syncPoints
298 (
299  const polyMesh& mesh,
300  pointField& points,
301  const CombineOp& cop,
302  const point& nullValue
303 )
304 {
305  if (points.size() != mesh.nPoints())
306  {
308  << "Number of values " << points.size()
309  << " is not equal to the number of points in the mesh "
310  << mesh.nPoints() << abort(FatalError);
311  }
312 
313  const polyBoundaryMesh& patches = mesh.boundaryMesh();
314 
315  // Is there any coupled patch with transformation?
316  bool hasTransformation = false;
317 
318  if (Pstream::parRun())
319  {
320  // Send
321 
322  forAll(patches, patchi)
323  {
324  const polyPatch& pp = patches[patchi];
325 
326  if
327  (
328  isA<processorPolyPatch>(pp)
329  && pp.nPoints() > 0
330  && refCast<const processorPolyPatch>(pp).owner()
331  )
332  {
333  const processorPolyPatch& procPatch =
334  refCast<const processorPolyPatch>(pp);
335 
336  // Get data per patchPoint in neighbouring point numbers.
337  pointField patchInfo(procPatch.nPoints(), nullValue);
338 
339  const labelList& meshPts = procPatch.meshPoints();
340  const labelList& nbrPts = procPatch.neighbPoints();
341 
342  forAll(nbrPts, pointi)
343  {
344  label nbrPointi = nbrPts[pointi];
345  if (nbrPointi >= 0 && nbrPointi < patchInfo.size())
346  {
347  patchInfo[nbrPointi] = points[meshPts[pointi]];
348  }
349  }
350 
351  OPstream toNbr(Pstream::blocking, procPatch.neighbProcNo());
352  toNbr << patchInfo;
353  }
354  }
355 
356 
357  // Receive and set.
358 
359  forAll(patches, patchi)
360  {
361  const polyPatch& pp = patches[patchi];
362 
363  if
364  (
365  isA<processorPolyPatch>(pp)
366  && pp.nPoints() > 0
367  && !refCast<const processorPolyPatch>(pp).owner()
368  )
369  {
370  const processorPolyPatch& procPatch =
371  refCast<const processorPolyPatch>(pp);
372 
373  pointField nbrPatchInfo(procPatch.nPoints());
374  {
375  // We do not know the number of points on the other side
376  // so cannot use Pstream::read.
377  IPstream fromNbr
378  (
380  procPatch.neighbProcNo()
381  );
382  fromNbr >> nbrPatchInfo;
383  }
384  // Null any value which is not on neighbouring processor
385  nbrPatchInfo.setSize(procPatch.nPoints(), nullValue);
386 
387  if (!procPatch.parallel())
388  {
389  hasTransformation = true;
390  transformList(procPatch.forwardT(), nbrPatchInfo);
391  }
392  else if (procPatch.separated())
393  {
394  hasTransformation = true;
395  separateList(-procPatch.separation(), nbrPatchInfo);
396  }
397 
398  const labelList& meshPts = procPatch.meshPoints();
399 
400  forAll(meshPts, pointi)
401  {
402  label meshPointi = meshPts[pointi];
403  points[meshPointi] = nbrPatchInfo[pointi];
404  }
405  }
406  }
407  }
408 
409  // Do the cyclics.
410  forAll(patches, patchi)
411  {
412  const polyPatch& pp = patches[patchi];
413 
414  if
415  (
416  isA<cyclicPolyPatch>(pp)
417  && refCast<const cyclicPolyPatch>(pp).owner()
418  )
419  {
420  const cyclicPolyPatch& cycPatch =
421  refCast<const cyclicPolyPatch>(pp);
422 
423  const edgeList& coupledPoints = cycPatch.coupledPoints();
424  const labelList& meshPts = cycPatch.meshPoints();
425  const cyclicPolyPatch& nbrPatch = cycPatch.neighbPatch();
426  const labelList& nbrMeshPts = nbrPatch.meshPoints();
427 
428  pointField half0Values(coupledPoints.size());
429 
430  forAll(coupledPoints, i)
431  {
432  const edge& e = coupledPoints[i];
433  label point0 = meshPts[e[0]];
434  half0Values[i] = points[point0];
435  }
436 
437  if (!cycPatch.parallel())
438  {
439  hasTransformation = true;
440  transformList(cycPatch.reverseT(), half0Values);
441  }
442  else if (cycPatch.separated())
443  {
444  hasTransformation = true;
445  separateList(cycPatch.separation(), half0Values);
446  }
447 
448  forAll(coupledPoints, i)
449  {
450  const edge& e = coupledPoints[i];
451  label point1 = nbrMeshPts[e[1]];
452  points[point1] = half0Values[i];
453  }
454  }
455  }
456 
457  //- Note: hasTransformation is only used for warning messages so
458  // reduction not strictly nessecary.
459  //reduce(hasTransformation, orOp<bool>());
460 
461  // Synchronize multiple shared points.
462  const globalMeshData& pd = mesh.globalData();
463 
464  if (pd.nGlobalPoints() > 0)
465  {
466  if (hasTransformation)
467  {
469  << "There are decomposed cyclics in this mesh with"
470  << " transformations." << endl
471  << "This is not supported. The result will be incorrect"
472  << endl;
473  }
474 
475 
476  // Values on shared points.
477  pointField sharedPts(pd.nGlobalPoints(), nullValue);
478 
479  forAll(pd.sharedPointLabels(), i)
480  {
481  label meshPointi = pd.sharedPointLabels()[i];
482  // Fill my entries in the shared points
483  sharedPts[pd.sharedPointAddr()[i]] = points[meshPointi];
484  }
485 
486  // Combine on master.
487  Pstream::listCombineGather(sharedPts, cop);
488  Pstream::listCombineScatter(sharedPts);
489 
490  // Now we will all have the same information. Merge it back with
491  // my local information.
492  forAll(pd.sharedPointLabels(), i)
493  {
494  label meshPointi = pd.sharedPointLabels()[i];
495  points[meshPointi] = sharedPts[pd.sharedPointAddr()[i]];
496  }
497  }
498 }
499 
500 
501 
502 int main(int argc, char *argv[])
503 {
504  #include "addOverwriteOption.H"
505  #include "addRegionOption.H"
506  #include "addDictOption.H"
507  #include "setRootCase.H"
508  #include "createTime.H"
509  runTime.functionObjects().off();
510 
511  Foam::word meshRegionName = polyMesh::defaultRegion;
512  args.optionReadIfPresent("region", meshRegionName);
513 
514  const bool overwrite = args.optionFound("overwrite");
515 
516  #include "createNamedPolyMesh.H"
517 
518  const word oldInstance = mesh.pointsInstance();
519 
520  const word dictName("createPatchDict");
521  #include "setSystemMeshDictionaryIO.H"
522 
523  Info<< "Reading " << dictName << nl << endl;
524 
526 
527  // Whether to synchronise points
528  const Switch pointSync(dict.lookup("pointSync"));
529 
530 
531  const polyBoundaryMesh& patches = mesh.boundaryMesh();
532 
533  // If running parallel check same patches everywhere
534  patches.checkParallelSync(true);
535 
536 
537  dumpCyclicMatch("initial_", mesh);
538 
539  // Read patch construct info from dictionary
540  PtrList<dictionary> patchSources(dict.lookup("patches"));
541 
542  HashSet<word> addedPatchNames;
543  forAll(patchSources, addedI)
544  {
545  const dictionary& dict = patchSources[addedI];
546  addedPatchNames.insert(dict.lookup("name"));
547  }
548 
549 
550  // 1. Add all new patches
551  // ~~~~~~~~~~~~~~~~~~~~~~
552 
553  if (patchSources.size())
554  {
555  // Old and new patches.
556  DynamicList<polyPatch*> allPatches(patches.size()+patchSources.size());
557 
558  label startFacei = mesh.nInternalFaces();
559 
560  // Copy old patches.
561  forAll(patches, patchi)
562  {
563  const polyPatch& pp = patches[patchi];
564 
565  if (!isA<processorPolyPatch>(pp))
566  {
567  allPatches.append
568  (
569  pp.clone
570  (
571  patches,
572  patchi,
573  pp.size(),
574  startFacei
575  ).ptr()
576  );
577  startFacei += pp.size();
578  }
579  }
580 
581  forAll(patchSources, addedI)
582  {
583  const dictionary& dict = patchSources[addedI];
584 
585  word patchName(dict.lookup("name"));
586 
587  label destPatchi = patches.findPatchID(patchName);
588 
589  if (destPatchi == -1)
590  {
591  dictionary patchDict(dict.subDict("patchInfo"));
592 
593  destPatchi = allPatches.size();
594 
595  Info<< "Adding new patch " << patchName
596  << " as patch " << destPatchi
597  << " from " << patchDict << endl;
598 
599  patchDict.set("nFaces", 0);
600  patchDict.set("startFace", startFacei);
601 
602  // Add an empty patch.
603  allPatches.append
604  (
606  (
607  patchName,
608  patchDict,
609  destPatchi,
610  patches
611  ).ptr()
612  );
613  }
614  else
615  {
616  Info<< "Patch '" << patchName << "' already exists. Only "
617  << "moving patch faces - type will remain the same" << endl;
618  }
619  }
620 
621  // Copy old patches.
622  forAll(patches, patchi)
623  {
624  const polyPatch& pp = patches[patchi];
625 
626  if (isA<processorPolyPatch>(pp))
627  {
628  allPatches.append
629  (
630  pp.clone
631  (
632  patches,
633  patchi,
634  pp.size(),
635  startFacei
636  ).ptr()
637  );
638  startFacei += pp.size();
639  }
640  }
641 
642  allPatches.shrink();
643  mesh.removeBoundary();
644  mesh.addPatches(allPatches);
645 
646  Info<< endl;
647  }
648 
649 
650 
651  // 2. Repatch faces
652  // ~~~~~~~~~~~~~~~~
653 
654  polyTopoChange meshMod(mesh);
655 
656 
657  forAll(patchSources, addedI)
658  {
659  const dictionary& dict = patchSources[addedI];
660 
661  const word patchName(dict.lookup("name"));
662  label destPatchi = patches.findPatchID(patchName);
663 
664  if (destPatchi == -1)
665  {
667  << "patch " << patchName << " not added. Problem."
668  << abort(FatalError);
669  }
670 
671  const word sourceType(dict.lookup("constructFrom"));
672 
673  if (sourceType == "patches")
674  {
675  labelHashSet patchSources
676  (
677  patches.patchSet
678  (
679  wordReList(dict.lookup("patches"))
680  )
681  );
682 
683  // Repatch faces of the patches.
684  forAllConstIter(labelHashSet, patchSources, iter)
685  {
686  const polyPatch& pp = patches[iter.key()];
687 
688  Info<< "Moving faces from patch " << pp.name()
689  << " to patch " << destPatchi << endl;
690 
691  forAll(pp, i)
692  {
693  changePatchID
694  (
695  mesh,
696  pp.start() + i,
697  destPatchi,
698  meshMod
699  );
700  }
701  }
702  }
703  else if (sourceType == "set")
704  {
705  const word setName(dict.lookup("set"));
706 
707  faceSet faces(mesh, setName);
708 
709  Info<< "Read " << returnReduce(faces.size(), sumOp<label>())
710  << " faces from faceSet " << faces.name() << endl;
711 
712  // Sort (since faceSet contains faces in arbitrary order)
713  labelList faceLabels(faces.toc());
714 
715  SortableList<label> patchFaces(faceLabels);
716 
717  forAll(patchFaces, i)
718  {
719  label facei = patchFaces[i];
720 
721  if (mesh.isInternalFace(facei))
722  {
724  << "Face " << facei << " specified in set "
725  << faces.name()
726  << " is not an external face of the mesh." << endl
727  << "This application can only repatch existing boundary"
728  << " faces." << exit(FatalError);
729  }
730 
731  changePatchID
732  (
733  mesh,
734  facei,
735  destPatchi,
736  meshMod
737  );
738  }
739  }
740  else
741  {
743  << "Invalid source type " << sourceType << endl
744  << "Valid source types are 'patches' 'set'" << exit(FatalError);
745  }
746  }
747  Info<< endl;
748 
749 
750  // Change mesh, use inflation to reforce calculation of transformation
751  // tensors.
752  Info<< "Doing topology modification to order faces." << nl << endl;
753  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, true);
754  mesh.movePoints(map().preMotionPoints());
755 
756  dumpCyclicMatch("coupled_", mesh);
757 
758  // Synchronise points.
759  if (!pointSync)
760  {
761  Info<< "Not synchronising points." << nl << endl;
762  }
763  else
764  {
765  Info<< "Synchronising points." << nl << endl;
766 
767  // This is a bit tricky. Both normal and position might be out and
768  // current separation also includes the normal
769  // ( separation_ = (nf&(Cr - Cf))*nf ).
770 
771  // For cyclic patches:
772  // - for separated ones use user specified offset vector
773 
774  forAll(mesh.boundaryMesh(), patchi)
775  {
776  const polyPatch& pp = mesh.boundaryMesh()[patchi];
777 
778  if (pp.size() && isA<coupledPolyPatch>(pp))
779  {
780  const coupledPolyPatch& cpp =
781  refCast<const coupledPolyPatch>(pp);
782 
783  if (cpp.separated())
784  {
785  Info<< "On coupled patch " << pp.name()
786  << " separation[0] was "
787  << cpp.separation()[0] << endl;
788 
789  if (isA<cyclicPolyPatch>(pp) && pp.size())
790  {
791  const cyclicPolyPatch& cycpp =
792  refCast<const cyclicPolyPatch>(pp);
793 
795  {
796  // Force to wanted separation
797  Info<< "On cyclic translation patch " << pp.name()
798  << " forcing uniform separation of "
799  << cycpp.separationVector() << endl;
800  const_cast<vectorField&>(cpp.separation()) =
801  pointField(1, cycpp.separationVector());
802  }
803  else
804  {
805  const cyclicPolyPatch& nbr = cycpp.neighbPatch();
806  const_cast<vectorField&>(cpp.separation()) =
807  pointField
808  (
809  1,
810  nbr[0].centre(mesh.points())
811  - cycpp[0].centre(mesh.points())
812  );
813  }
814  }
815  Info<< "On coupled patch " << pp.name()
816  << " forcing uniform separation of "
817  << cpp.separation() << endl;
818  }
819  else if (!cpp.parallel())
820  {
821  Info<< "On coupled patch " << pp.name()
822  << " forcing uniform rotation of "
823  << cpp.forwardT()[0] << endl;
824 
825  const_cast<tensorField&>
826  (
827  cpp.forwardT()
828  ).setSize(1);
829  const_cast<tensorField&>
830  (
831  cpp.reverseT()
832  ).setSize(1);
833 
834  Info<< "On coupled patch " << pp.name()
835  << " forcing uniform rotation of "
836  << cpp.forwardT() << endl;
837  }
838  }
839  }
840 
841  Info<< "Synchronising points." << endl;
842 
843  pointField newPoints(mesh.points());
844 
845  syncPoints
846  (
847  mesh,
848  newPoints,
850  point(GREAT, GREAT, GREAT)
851  );
852 
853  scalarField diff(mag(newPoints-mesh.points()));
854  Info<< "Points changed by average:" << gAverage(diff)
855  << " max:" << gMax(diff) << nl << endl;
856 
857  mesh.movePoints(newPoints);
858  }
859 
860  // 3. Remove zeros-sized patches
861  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
862 
863  Info<< "Removing patches with no faces in them." << nl<< endl;
864  filterPatches(mesh, addedPatchNames);
865 
866 
867  dumpCyclicMatch("final_", mesh);
868 
869 
870  // Set the precision of the points data to 10
872 
873  if (!overwrite)
874  {
875  runTime++;
876  }
877  else
878  {
879  mesh.setInstance(oldInstance);
880  }
881 
882  // Write resulting mesh
883  Info<< "Writing repatched mesh to " << runTime.timeName() << nl << endl;
884  mesh.write();
885 
886  Info<< "End\n" << endl;
887 
888  return 0;
889 }
890 
891 
892 // ************************************************************************* //
const cyclicPolyPatch & neighbPatch() const
label nPoints() const
Return number of points supporting patch faces.
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
Definition: triad.C:394
A HashTable with keys but without contents.
Definition: HashSet.H:59
void append(T *)
Append an element at the end of the list.
Definition: PtrListI.H:39
dictionary dict
const labelList & sharedPointAddr() const
Return addressing into the complete globally shared points.
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
Definition: polyMesh.C:1048
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A class for handling file names.
Definition: fileName.H:69
A list of face labels.
Definition: faceSet.H:48
bool checkParallelSync(const bool report=false) const
Check whether all procs have all patches and in same order. Return.
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
Definition: polyPatch.H:219
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const boolList & flipMap() const
Return face flip map.
Definition: faceZone.H:253
Class describing modification of a face.
virtual const tensorField & reverseT() const
Return neighbour-cell transformation tensor.
const double e
Elementary charge.
Definition: doubleFloat.H:78
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual const vectorField & separation() const
If the planes are separated the separation vector.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
A list that is sorted upon construction or when explicitly requested with the sort() method...
Definition: List.H:68
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
Output to file stream.
Definition: OFstream.H:81
int neighbProcNo() const
Return neigbour processor number.
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:461
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
Definition: argListI.H:198
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
const Field< PointType > & points() const
Return reference to global points.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
const edgeList & coupledPoints() const
Return connected points (from patch local to neighbour patch local)
const vectorField & faceCentres() const
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
label size() const
Return number of elements in list.
Definition: DLListBaseI.H:77
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
virtual transformType transform() const
Type of transform.
void transformList(const tensor &, UList< T > &)
Apply transformation to list. Either single transformation tensor.
Definition: transformList.C:49
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:979
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
Neighbour processor patch.
Input inter-processor communications stream.
Definition: IPstream.H:50
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
Definition: faceZone.C:305
points setSize(newPointi)
virtual bool separated() const
Are the planes separated.
word name() const
Return file name (part beyond last /)
Definition: fileName.C:179
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
A class for handling words, derived from string.
Definition: word.H:59
Cyclic plane patch.
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:88
virtual bool parallel() const
Are the cyclic planes parallel.
const word & name() const
Return name.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
defineTemplateTypeNameAndDebug(IOPtrList< ensightPart >, 0)
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
bool found(const Key &) const
Return true if hashedEntry is found in table.
Foam::polyBoundaryMesh.
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1140
void removeBoundary()
Remove boundary patches.
Definition: polyMeshClear.C:36
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
static const char nl
Definition: Ostream.H:262
Type gMax(const FieldField< Field, Type > &f)
Output inter-processor communications stream.
Definition: OPstream.H:50
IOobject dictIO(dictName, runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE)
void setInstance(const fileName &)
Set the instance for mesh files.
Definition: polyMeshIO.C:32
void addPatches(const List< polyPatch * > &, const bool validBoundary=true)
Add boundary patches.
Definition: polyMesh.C:878
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
A PtrList of objects of type <T> with automated input and output.
Definition: IOPtrList.H:50
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:393
label patchi
vector point
Point is a vector.
Definition: point.H:41
#define WarningInFunction
Report a warning using Foam::Warning.
word dictName("noiseDict")
Type gAverage(const FieldField< Field, Type > &f)
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:772
label whichZone(const label objectIndex) const
Given a global object index, return the zone it is in.
Definition: ZoneMesh.C:221
Direct mesh changes based on v1.3 polyTopoChange syntax.
virtual bool write() const
Write using setting from DB.
const labelList & neighbPoints() const
Return neighbour point labels. WIP.
label nPoints() const
messageStream Info
const vector & separationVector() const
Translation vector for translational cyclics.
virtual const tensorField & forwardT() const
Return face transformation tensor.
dimensioned< scalar > mag(const dimensioned< Type > &)
const faceZoneMesh & faceZones() const
Return face zone mesh.
Definition: polyMesh.H:463
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:64
autoPtr< mapPolyMesh > changeMesh(polyMesh &mesh, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
const labelList & sharedPointLabels() const
Return indices of local points that are globally shared.
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1017
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1004
Foam::argList args(argc, argv)
label nGlobalPoints() const
Return number of globally shared points.
label nInternalFaces() const
Namespace for OpenFOAM.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451