mappedPatchBase.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2022 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 \*---------------------------------------------------------------------------*/
25 
26 #include "mappedPatchBase.H"
28 #include "ListListOps.H"
29 #include "meshSearchMeshObject.H"
30 #include "meshTools.H"
31 #include "OFstream.H"
32 #include "Random.H"
33 #include "treeDataFace.H"
34 #include "treeDataPoint.H"
35 #include "indexedOctree.H"
36 #include "polyMesh.H"
37 #include "polyPatch.H"
38 #include "Time.H"
39 #include "distributionMap.H"
40 #include "SubField.H"
41 #include "triPointRef.H"
42 #include "syncTools.H"
43 #include "treeDataCell.H"
44 #include "DynamicField.H"
45 
46 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47 
48 namespace Foam
49 {
50  defineTypeNameAndDebug(mappedPatchBase, 0);
51 
52  template<>
53  const char* Foam::NamedEnum
54  <
56  6
57  >::names[] =
58  {
59  "nearestCell",
60  "nearestPatchFace",
61  "nearestPatchFaceAMI",
62  "nearestPatchPoint",
63  "nearestFace",
64  "nearestOnlyCell"
65  };
66 
67  template<>
68  const char* Foam::NamedEnum
69  <
71  3
72  >::names[] =
73  {
74  "uniform",
75  "nonuniform",
76  "normal"
77  };
78 }
79 
80 
83 
86 
87 
88 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
89 
91 (
92  const polyPatch& pp
93 ) const
94 {
95  const polyMesh& mesh = pp.boundaryMesh().mesh();
96 
97  // Force construction of min-tet decomp
98  (void)mesh.tetBasePtIs();
99 
100  // Initialise to face-centre
101  tmp<pointField> tfacePoints(new pointField(patch_.size()));
102  pointField& facePoints = tfacePoints.ref();
103 
104  forAll(pp, facei)
105  {
106  facePoints[facei] = facePoint
107  (
108  mesh,
109  pp.start()+facei,
111  ).rawPoint();
112  }
113 
114  return tfacePoints;
115 }
116 
117 
119 (
120  const pointField& facePoints,
121  pointField& samples,
122  labelList& patchFaceProcs,
123  labelList& patchFaces,
124  pointField& patchFc
125 ) const
126 {
127  // Collect all sample points and the faces they come from.
128  {
129  List<pointField> globalFc(Pstream::nProcs());
130  globalFc[Pstream::myProcNo()] = facePoints;
131  Pstream::gatherList(globalFc);
132  Pstream::scatterList(globalFc);
133  // Rework into straight list
134  patchFc = ListListOps::combine<pointField>
135  (
136  globalFc,
138  );
139  }
140 
141  {
142  List<pointField> globalSamples(Pstream::nProcs());
143  globalSamples[Pstream::myProcNo()] = samplePoints(facePoints);
144  Pstream::gatherList(globalSamples);
145  Pstream::scatterList(globalSamples);
146  // Rework into straight list
147  samples = ListListOps::combine<pointField>
148  (
149  globalSamples,
151  );
152  }
153 
154  {
155  labelListList globalFaces(Pstream::nProcs());
156  globalFaces[Pstream::myProcNo()] = identity(patch_.size());
157  // Distribute to all processors
158  Pstream::gatherList(globalFaces);
159  Pstream::scatterList(globalFaces);
160 
161  patchFaces = ListListOps::combine<labelList>
162  (
163  globalFaces,
165  );
166  }
167 
168  {
169  labelList nPerProc(Pstream::nProcs());
170  nPerProc[Pstream::myProcNo()] = patch_.size();
171  Pstream::gatherList(nPerProc);
172  Pstream::scatterList(nPerProc);
173 
174  patchFaceProcs.setSize(patchFaces.size());
175 
176  label sampleI = 0;
177  forAll(nPerProc, proci)
178  {
179  for (label i = 0; i < nPerProc[proci]; i++)
180  {
181  patchFaceProcs[sampleI++] = proci;
182  }
183  }
184  }
185 }
186 
187 
188 // Find the processor/cell containing the samples. Does not account
189 // for samples being found in two processors.
191 (
192  const sampleMode mode,
193  const pointField& samples,
194  labelList& sampleProcs,
195  labelList& sampleIndices,
196  pointField& sampleLocations
197 ) const
198 {
199  // Lookup the correct region
200  const polyMesh& mesh = sampleMesh();
201 
202  // All the info for nearest. Construct to miss
203  List<nearInfo> nearest(samples.size());
204 
205  switch (mode)
206  {
207  case NEARESTCELL:
208  {
209  if (samplePatch_.size() && samplePatch_ != "none")
210  {
212  << "No need to supply a patch name when in "
213  << sampleModeNames_[mode] << " mode." << exit(FatalError);
214  }
215 
216  //- Note: face-diagonal decomposition
217  const indexedOctree<Foam::treeDataCell>& tree = mesh.cellTree();
218 
219  forAll(samples, sampleI)
220  {
221  const point& sample = samples[sampleI];
222 
223  label celli = tree.findInside(sample);
224 
225  if (celli == -1)
226  {
227  nearest[sampleI].second().first() = Foam::sqr(great);
228  nearest[sampleI].second().second() = Pstream::myProcNo();
229  }
230  else
231  {
232  const point& cc = mesh.cellCentres()[celli];
233 
234  nearest[sampleI].first() = pointIndexHit
235  (
236  true,
237  cc,
238  celli
239  );
240  nearest[sampleI].second().first() = magSqr(cc-sample);
241  nearest[sampleI].second().second() = Pstream::myProcNo();
242  }
243  }
244  break;
245  }
246 
247  case NEARESTONLYCELL:
248  {
249  if (samplePatch_.size() && samplePatch_ != "none")
250  {
252  << "No need to supply a patch name when in "
253  << sampleModeNames_[mode] << " mode." << exit(FatalError);
254  }
255 
256  //- Note: face-diagonal decomposition
257  const indexedOctree<Foam::treeDataCell>& tree = mesh.cellTree();
258 
259  forAll(samples, sampleI)
260  {
261  const point& sample = samples[sampleI];
262 
263  nearest[sampleI].first() = tree.findNearest(sample, sqr(great));
264  nearest[sampleI].second().first() = magSqr
265  (
266  nearest[sampleI].first().hitPoint()
267  -sample
268  );
269  nearest[sampleI].second().second() = Pstream::myProcNo();
270  }
271  break;
272  }
273 
274  case NEARESTPATCHFACE:
275  {
276  const polyPatch& pp = samplePolyPatch();
277 
278  if (pp.empty())
279  {
280  forAll(samples, sampleI)
281  {
282  nearest[sampleI].second().first() = Foam::sqr(great);
283  nearest[sampleI].second().second() = Pstream::myProcNo();
284  }
285  }
286  else
287  {
288  // patch faces
289  const labelList patchFaces(identity(pp.size()) + pp.start());
290 
291  treeBoundBox patchBb
292  (
293  treeBoundBox(pp.points(), pp.meshPoints()).extend(1e-4)
294  );
295 
296  indexedOctree<treeDataFace> boundaryTree
297  (
298  treeDataFace // all information needed to search faces
299  (
300  false, // do not cache bb
301  mesh,
302  patchFaces // boundary faces only
303  ),
304  patchBb, // overall search domain
305  8, // maxLevel
306  10, // leafsize
307  3.0 // duplicity
308  );
309 
310  forAll(samples, sampleI)
311  {
312  const point& sample = samples[sampleI];
313 
314  pointIndexHit& nearInfo = nearest[sampleI].first();
315  nearInfo = boundaryTree.findNearest
316  (
317  sample,
318  magSqr(patchBb.span())
319  );
320 
321  if (!nearInfo.hit())
322  {
323  nearest[sampleI].second().first() = Foam::sqr(great);
324  nearest[sampleI].second().second() =
326  }
327  else
328  {
329  point fc(pp[nearInfo.index()].centre(pp.points()));
330  nearInfo.setPoint(fc);
331  nearest[sampleI].second().first() = magSqr(fc-sample);
332  nearest[sampleI].second().second() =
334  }
335  }
336  }
337  break;
338  }
339 
340  case NEARESTPATCHPOINT:
341  {
342  const polyPatch& pp = samplePolyPatch();
343 
344  if (pp.empty())
345  {
346  forAll(samples, sampleI)
347  {
348  nearest[sampleI].second().first() = Foam::sqr(great);
349  nearest[sampleI].second().second() = Pstream::myProcNo();
350  }
351  }
352  else
353  {
354  // patch (local) points
355  treeBoundBox patchBb
356  (
357  treeBoundBox(pp.points(), pp.meshPoints()).extend(1e-4)
358  );
359 
360  indexedOctree<treeDataPoint> boundaryTree
361  (
362  treeDataPoint // all information needed to search faces
363  (
364  mesh.points(),
365  pp.meshPoints() // selection of points to search on
366  ),
367  patchBb, // overall search domain
368  8, // maxLevel
369  10, // leafsize
370  3.0 // duplicity
371  );
372 
373  forAll(samples, sampleI)
374  {
375  const point& sample = samples[sampleI];
376 
377  pointIndexHit& nearInfo = nearest[sampleI].first();
378  nearInfo = boundaryTree.findNearest
379  (
380  sample,
381  magSqr(patchBb.span())
382  );
383 
384  if (!nearInfo.hit())
385  {
386  nearest[sampleI].second().first() = Foam::sqr(great);
387  nearest[sampleI].second().second() =
389  }
390  else
391  {
392  const point& pt = nearInfo.hitPoint();
393 
394  nearest[sampleI].second().first() = magSqr(pt-sample);
395  nearest[sampleI].second().second() =
397  }
398  }
399  }
400  break;
401  }
402 
403  case NEARESTFACE:
404  {
405  if (samplePatch().size() && samplePatch() != "none")
406  {
408  << "No need to supply a patch name when in "
409  << sampleModeNames_[mode] << " mode." << exit(FatalError);
410  }
411 
412  //- Note: face-diagonal decomposition
413  const meshSearchMeshObject& meshSearchEngine =
415 
416  forAll(samples, sampleI)
417  {
418  const point& sample = samples[sampleI];
419 
420  label facei = meshSearchEngine.findNearestFace(sample);
421 
422  if (facei == -1)
423  {
424  nearest[sampleI].second().first() = Foam::sqr(great);
425  nearest[sampleI].second().second() = Pstream::myProcNo();
426  }
427  else
428  {
429  const point& fc = mesh.faceCentres()[facei];
430 
431  nearest[sampleI].first() = pointIndexHit
432  (
433  true,
434  fc,
435  facei
436  );
437  nearest[sampleI].second().first() = magSqr(fc-sample);
438  nearest[sampleI].second().second() = Pstream::myProcNo();
439  }
440  }
441  break;
442  }
443 
444  case NEARESTPATCHFACEAMI:
445  {
446  // nothing to do here
447  return;
448  }
449 
450  default:
451  {
453  << "problem." << abort(FatalError);
454  }
455  }
456 
457 
458  // Find nearest. Combine on master.
461 
462 
463  if (debug)
464  {
466  << "mesh " << sampleRegion() << " : " << endl;
467 
468  forAll(nearest, sampleI)
469  {
470  label proci = nearest[sampleI].second().second();
471  label localI = nearest[sampleI].first().index();
472 
473  Info<< " " << sampleI << " coord:"<< samples[sampleI]
474  << " found on processor:" << proci
475  << " in local cell/face/point:" << localI
476  << " with location:" << nearest[sampleI].first().rawPoint()
477  << endl;
478  }
479  }
480 
481  // Convert back into proc+local index
482  sampleProcs.setSize(samples.size());
483  sampleIndices.setSize(samples.size());
484  sampleLocations.setSize(samples.size());
485 
486  forAll(nearest, sampleI)
487  {
488  if (!nearest[sampleI].first().hit())
489  {
490  sampleProcs[sampleI] = -1;
491  sampleIndices[sampleI] = -1;
492  sampleLocations[sampleI] = vector::max;
493  }
494  else
495  {
496  sampleProcs[sampleI] = nearest[sampleI].second().second();
497  sampleIndices[sampleI] = nearest[sampleI].first().index();
498  sampleLocations[sampleI] = nearest[sampleI].first().hitPoint();
499  }
500  }
501 }
502 
503 
505 {
506  static bool hasWarned = false;
507  if (mapPtr_.valid())
508  {
510  << "Mapping already calculated" << exit(FatalError);
511  }
512 
513  // Get points on face (since cannot use face-centres - might be off
514  // face-diagonal decomposed tets.
515  tmp<pointField> patchPoints(facePoints(patch_));
516 
517  // Get offsetted points
518  const pointField offsettedPoints(samplePoints(patchPoints()));
519 
520  // Do a sanity check - am I sampling my own patch?
521  // This only makes sense for a non-zero offset.
522  bool sampleMyself =
523  (
524  mode_ == NEARESTPATCHFACE
525  && sampleRegion() == patch_.boundaryMesh().mesh().name()
526  && samplePatch() == patch_.name()
527  );
528 
529  // Check offset
530  vectorField d(offsettedPoints-patchPoints());
531  bool coincident = (gAverage(mag(d)) <= rootVSmall);
532 
533  if (sampleMyself && coincident)
534  {
536  << "Invalid offset " << d << endl
537  << "Offset is the vector added to the patch face centres to"
538  << " find the patch face supplying the data." << endl
539  << "Setting it to " << d
540  << " on the same patch, on the same region"
541  << " will find the faces themselves which does not make sense"
542  << " for anything but testing." << endl
543  << "patch_:" << patch_.name() << endl
544  << "sampleRegion_:" << sampleRegion() << endl
545  << "mode_:" << sampleModeNames_[mode_] << endl
546  << "samplePatch_:" << samplePatch() << endl
547  << "offsetMode_:" << offsetModeNames_[offsetMode_] << endl;
548  }
549 
550  // Get global list of all samples and the processor and face they come from.
552  labelList patchFaceProcs;
553  labelList patchFaces;
554  pointField patchFc;
555  collectSamples
556  (
557  patchPoints,
558  samples,
559  patchFaceProcs,
560  patchFaces,
561  patchFc
562  );
563 
564  // Find processor and cell/face samples are in and actual location.
565  labelList sampleProcs;
566  labelList sampleIndices;
567  pointField sampleLocations;
568  findSamples(mode_, samples, sampleProcs, sampleIndices, sampleLocations);
569 
570  // Check for samples that were not found. This will only happen for
571  // NEARESTCELL since finds cell containing a location
572  if (mode_ == NEARESTCELL)
573  {
574  label nNotFound = 0;
575  forAll(sampleProcs, sampleI)
576  {
577  if (sampleProcs[sampleI] == -1)
578  {
579  nNotFound++;
580  }
581  }
582  reduce(nNotFound, sumOp<label>());
583 
584  if (nNotFound > 0)
585  {
586  if (!hasWarned)
587  {
589  << "Did not find " << nNotFound
590  << " out of " << sampleProcs.size() << " total samples."
591  << " Sampling these on owner cell centre instead." << endl
592  << "On patch " << patch_.name()
593  << " on region " << sampleRegion()
594  << " in mode " << sampleModeNames_[mode_] << endl
595  << "with offset mode " << offsetModeNames_[offsetMode_]
596  << ". Suppressing further warnings from " << type() << endl;
597 
598  hasWarned = true;
599  }
600 
601  // Collect the samples that cannot be found
602  DynamicList<label> subMap;
603  DynamicField<point> subSamples;
604 
605  forAll(sampleProcs, sampleI)
606  {
607  if (sampleProcs[sampleI] == -1)
608  {
609  subMap.append(sampleI);
610  subSamples.append(samples[sampleI]);
611  }
612  }
613 
614  // And re-search for pure nearest (should not fail)
615  labelList subSampleProcs;
616  labelList subSampleIndices;
617  pointField subSampleLocations;
618  findSamples
619  (
620  NEARESTONLYCELL,
621  subSamples,
622  subSampleProcs,
623  subSampleIndices,
624  subSampleLocations
625  );
626 
627  // Insert
628  UIndirectList<label>(sampleProcs, subMap) = subSampleProcs;
629  UIndirectList<label>(sampleIndices, subMap) = subSampleIndices;
630  UIndirectList<point>(sampleLocations, subMap) = subSampleLocations;
631  }
632  }
633 
634  // Now we have all the data we need:
635  // - where sample originates from (so destination when mapping):
636  // patchFaces, patchFaceProcs.
637  // - cell/face sample is in (so source when mapping)
638  // sampleIndices, sampleProcs.
639 
640  // forAll(samples, i)
641  // {
642  // Info<< i << " need data in region "
643  // << patch_.boundaryMesh().mesh().name()
644  // << " for proc:" << patchFaceProcs[i]
645  // << " face:" << patchFaces[i]
646  // << " at:" << patchFc[i] << endl
647  // << "Found data in region " << sampleRegion()
648  // << " at proc:" << sampleProcs[i]
649  // << " face:" << sampleIndices[i]
650  // << " at:" << sampleLocations[i]
651  // << nl << endl;
652  // }
653 
654  bool mapSucceeded = true;
655 
656  forAll(samples, i)
657  {
658  if (sampleProcs[i] == -1)
659  {
660  mapSucceeded = false;
661  break;
662  }
663  }
664 
665  if (!mapSucceeded)
666  {
668  << "Mapping failed for " << nl
669  << " patch: " << patch_.name() << nl
670  << " sampleRegion: " << sampleRegion() << nl
671  << " mode: " << sampleModeNames_[mode_] << nl
672  << " samplePatch: " << samplePatch() << nl
673  << " offsetMode: " << offsetModeNames_[offsetMode_]
674  << exit(FatalError);
675  }
676 
677  if (debug && Pstream::master())
678  {
679  OFstream str
680  (
681  patch_.boundaryMesh().mesh().time().path()
682  / patch_.name()
683  + "_mapped.obj"
684  );
685  Pout<< "Dumping mapping as lines from patch faceCentres to"
686  << " sampled cell/faceCentres/points to file " << str.name()
687  << endl;
688 
689  label vertI = 0;
690 
691  forAll(patchFc, i)
692  {
693  meshTools::writeOBJ(str, patchFc[i]);
694  vertI++;
695  meshTools::writeOBJ(str, sampleLocations[i]);
696  vertI++;
697  str << "l " << vertI-1 << ' ' << vertI << nl;
698  }
699  }
700 
701  // Determine schedule.
702  mapPtr_.reset(new distributionMap(sampleProcs, patchFaceProcs));
703 
704  // Rework the schedule from indices into samples to cell data to send,
705  // face data to receive.
706 
707  labelListList& subMap = mapPtr_().subMap();
708  labelListList& constructMap = mapPtr_().constructMap();
709 
710  forAll(subMap, proci)
711  {
712  subMap[proci] = UIndirectList<label>
713  (
714  sampleIndices,
715  subMap[proci]
716  );
717  constructMap[proci] = UIndirectList<label>
718  (
719  patchFaces,
720  constructMap[proci]
721  );
722 
723  // if (debug)
724  //{
725  // Pout<< "To proc:" << proci << " sending values of cells/faces:"
726  // << subMap[proci] << endl;
727  // Pout<< "From proc:" << proci
728  // << " receiving values of patch faces:"
729  // << constructMap[proci] << endl;
730  //}
731  }
732 
733  // Redo constructSize
734  mapPtr_().constructSize() = patch_.size();
735 
736  if (debug)
737  {
738  // Check that all elements get a value.
739  PackedBoolList used(patch_.size());
740  forAll(constructMap, proci)
741  {
742  const labelList& map = constructMap[proci];
743 
744  forAll(map, i)
745  {
746  label facei = map[i];
747 
748  if (used[facei] == 0)
749  {
750  used[facei] = 1;
751  }
752  else
753  {
755  << "On patch " << patch_.name()
756  << " patchface " << facei
757  << " is assigned to more than once."
758  << abort(FatalError);
759  }
760  }
761  }
762  forAll(used, facei)
763  {
764  if (used[facei] == 0)
765  {
767  << "On patch " << patch_.name()
768  << " patchface " << facei
769  << " is never assigned to."
770  << abort(FatalError);
771  }
772  }
773  }
774 }
775 
776 
778 const
779 {
780  const word surfType(surfDict_.lookupOrDefault<word>("type", "none"));
781 
782  if (!surfPtr_.valid() && surfType != "none")
783  {
784  word surfName(surfDict_.lookupOrDefault("name", patch_.name()));
785 
786  const polyMesh& mesh = patch_.boundaryMesh().mesh();
787 
788  surfPtr_ =
790  (
791  surfType,
792  IOobject
793  (
794  surfName,
795  mesh.time().constant(),
797  mesh,
800  ),
801  surfDict_
802  );
803  }
804 
805  return surfPtr_;
806 }
807 
808 
810 {
811  if (AMIPtr_.valid())
812  {
814  << "AMI already calculated" << exit(FatalError);
815  }
816 
817  AMIPtr_.clear();
818 
819  if (debug)
820  {
821  const polyPatch& nbr = samplePolyPatch();
822 
823  pointField nbrPoints(nbr.localPoints());
824 
825  OFstream os(patch_.name() + "_neighbourPatch-org.obj");
826  meshTools::writeOBJ(os, samplePolyPatch().localFaces(), nbrPoints);
827 
828  // transform neighbour patch to local system
829  primitivePatch nbrPatch0
830  (
832  (
833  nbr.localFaces(),
834  nbr.size()
835  ),
836  nbrPoints
837  );
838 
839  OFstream osN(patch_.name() + "_neighbourPatch-trans.obj");
840  meshTools::writeOBJ(osN, nbrPatch0, nbrPoints);
841 
842  OFstream osO(patch_.name() + "_ownerPatch.obj");
843  meshTools::writeOBJ(osO, patch_.localFaces(), patch_.localPoints());
844  }
845 
846  // Construct/apply AMI interpolation to determine addressing and weights
847  AMIPtr_.reset
848  (
849  new AMIInterpolation
850  (
851  patch_,
852  samplePolyPatch(), // nbrPatch0,
853  surfPtr(),
855  true,
856  faceAreaWeightAMI::typeName,
857  -1,
858  AMIReverse_
859  )
860  );
861 }
862 
863 
864 // Hack to read old (List-based) format. See Field.C. The difference
865 // is only that in case of falling back to old format it expects a non-uniform
866 // list instead of a single vector.
868 (
869  const word& keyword,
870  const dictionary& dict,
871  const label size
872 )
873 {
874  tmp<pointField> tfld(new pointField());
875  pointField& fld = tfld.ref();
876 
877  if (size)
878  {
879  ITstream& is = dict.lookup(keyword);
880 
881  // Read first token
882  token firstToken(is);
883 
884  if (firstToken.isWord())
885  {
886  if (firstToken.wordToken() == "uniform")
887  {
888  fld.setSize(size);
889  fld = pTraits<vector>(is);
890  }
891  else if (firstToken.wordToken() == "nonuniform")
892  {
893  is >> static_cast<List<vector>&>(fld);
894  if (fld.size() != size)
895  {
897  (
898  dict
899  ) << "size " << fld.size()
900  << " is not equal to the given value of " << size
901  << exit(FatalIOError);
902  }
903  }
904  else
905  {
907  (
908  dict
909  ) << "expected keyword 'uniform' or 'nonuniform', found "
910  << firstToken.wordToken()
911  << exit(FatalIOError);
912  }
913  }
914  else
915  {
916  if (is.version() == 2.0)
917  {
919  (
920  dict
921  ) << "expected keyword 'uniform' or 'nonuniform', "
922  "assuming List format for backwards compatibility."
923  "Foam version 2.0." << endl;
924 
925  is.putBack(firstToken);
926  is >> static_cast<List<vector>&>(fld);
927  }
928  }
929  }
930  return tfld;
931 }
932 
933 
934 // * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * * * //
935 
937 (
938  const polyPatch& pp
939 )
940 :
941  patch_(pp),
942  sampleRegion_(patch_.boundaryMesh().mesh().name()),
943  mode_(NEARESTPATCHFACE),
944  samplePatch_(""),
945  coupleGroup_(),
946  offsetMode_(UNIFORM),
947  offset_(Zero),
948  offsets_(pp.size(), offset_),
949  distance_(0),
950  sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()),
951  mapPtr_(nullptr),
952  AMIPtr_(nullptr),
953  AMIReverse_(false),
954  surfPtr_(nullptr),
955  surfDict_(fileName("surface"))
956 {}
957 
958 
960 (
961  const polyPatch& pp,
962  const word& sampleRegion,
963  const sampleMode mode,
964  const word& samplePatch,
965  const vectorField& offsets
966 )
967 :
968  patch_(pp),
969  sampleRegion_(sampleRegion),
970  mode_(mode),
971  samplePatch_(samplePatch),
972  coupleGroup_(),
973  offsetMode_(NONUNIFORM),
974  offset_(Zero),
975  offsets_(offsets),
976  distance_(0),
977  sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()),
978  mapPtr_(nullptr),
979  AMIPtr_(nullptr),
980  AMIReverse_(false),
981  surfPtr_(nullptr),
982  surfDict_(fileName("surface"))
983 {}
984 
985 
987 (
988  const polyPatch& pp,
989  const word& sampleRegion,
990  const sampleMode mode,
991  const word& samplePatch,
992  const vector& offset
993 )
994 :
995  patch_(pp),
996  sampleRegion_(sampleRegion),
997  mode_(mode),
998  samplePatch_(samplePatch),
999  coupleGroup_(),
1000  offsetMode_(UNIFORM),
1001  offset_(offset),
1002  offsets_(0),
1003  distance_(0),
1004  sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()),
1005  mapPtr_(nullptr),
1006  AMIPtr_(nullptr),
1007  AMIReverse_(false),
1008  surfPtr_(nullptr),
1009  surfDict_(fileName("surface"))
1010 {}
1011 
1012 
1015  const polyPatch& pp,
1016  const word& sampleRegion,
1017  const sampleMode mode,
1018  const word& samplePatch,
1019  const scalar distance
1020 )
1021 :
1022  patch_(pp),
1023  sampleRegion_(sampleRegion),
1024  mode_(mode),
1025  samplePatch_(samplePatch),
1026  coupleGroup_(),
1027  offsetMode_(NORMAL),
1028  offset_(Zero),
1029  offsets_(0),
1030  distance_(distance),
1031  sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()),
1032  mapPtr_(nullptr),
1033  AMIPtr_(nullptr),
1034  AMIReverse_(false),
1035  surfPtr_(nullptr),
1036  surfDict_(fileName("surface"))
1037 {}
1038 
1039 
1042  const polyPatch& pp,
1043  const dictionary& dict
1044 )
1045 :
1046  patch_(pp),
1047  sampleRegion_(dict.lookupOrDefault<word>("sampleRegion", "")),
1048  mode_(sampleModeNames_.read(dict.lookup("sampleMode"))),
1049  samplePatch_(dict.lookupOrDefault<word>("samplePatch", "")),
1050  coupleGroup_(dict),
1051  offsetMode_(UNIFORM),
1052  offset_(Zero),
1053  offsets_(0),
1054  distance_(0.0),
1055  sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()),
1056  mapPtr_(nullptr),
1057  AMIPtr_(nullptr),
1058  AMIReverse_(dict.lookupOrDefault<bool>("flipNormals", false)),
1059  surfPtr_(nullptr),
1060  surfDict_(dict.subOrEmptyDict("surface"))
1061 {
1062  if (!coupleGroup_.valid())
1063  {
1064  if (sampleRegion_.empty())
1065  {
1066  // If no coupleGroup and no sampleRegion assume local region
1067  sampleRegion_ = patch_.boundaryMesh().mesh().name();
1068  sameRegion_ = true;
1069  }
1070  }
1071 
1072  if (dict.found("offsetMode"))
1073  {
1074  offsetMode_ = offsetModeNames_.read(dict.lookup("offsetMode"));
1075 
1076  switch(offsetMode_)
1077  {
1078  case UNIFORM:
1079  {
1080  offset_ = point(dict.lookup("offset"));
1081  }
1082  break;
1083 
1084  case NONUNIFORM:
1085  {
1086  offsets_ = readListOrField("offsets", dict, patch_.size());
1087  }
1088  break;
1089 
1090  case NORMAL:
1091  {
1092  distance_ = dict.lookup<scalar>("distance");
1093  }
1094  break;
1095  }
1096  }
1097  else if (dict.found("offset"))
1098  {
1099  offsetMode_ = UNIFORM;
1100  offset_ = point(dict.lookup("offset"));
1101  }
1102  else if (dict.found("offsets"))
1103  {
1104  offsetMode_ = NONUNIFORM;
1105  offsets_ = readListOrField("offsets", dict, patch_.size());
1106  }
1107  else if (mode_ != NEARESTPATCHFACE && mode_ != NEARESTPATCHFACEAMI)
1108  {
1110  << "Please supply the offsetMode as one of "
1112  << exit(FatalIOError);
1113  }
1114 }
1115 
1116 
1119  const polyPatch& pp,
1120  const sampleMode mode,
1121  const dictionary& dict
1122 )
1123 :
1124  patch_(pp),
1125  sampleRegion_(dict.lookupOrDefault<word>("sampleRegion", "")),
1126  mode_(mode),
1127  samplePatch_(dict.lookupOrDefault<word>("samplePatch", "")),
1128  coupleGroup_(dict),
1129  offsetMode_(UNIFORM),
1130  offset_(Zero),
1131  offsets_(0),
1132  distance_(0.0),
1133  sameRegion_(sampleRegion_ == patch_.boundaryMesh().mesh().name()),
1134  mapPtr_(nullptr),
1135  AMIPtr_(nullptr),
1136  AMIReverse_(dict.lookupOrDefault<bool>("flipNormals", false)),
1137  surfPtr_(nullptr),
1138  surfDict_(dict.subOrEmptyDict("surface"))
1139 {
1140  if (mode != NEARESTPATCHFACE && mode != NEARESTPATCHFACEAMI)
1141  {
1143  << "Construct from sampleMode and dictionary only applicable for "
1144  << " collocated patches in modes "
1145  << sampleModeNames_[NEARESTPATCHFACE] << ','
1146  << sampleModeNames_[NEARESTPATCHFACEAMI]
1147  << exit(FatalIOError);
1148  }
1149 
1150  if (!coupleGroup_.valid())
1151  {
1152  if (sampleRegion_.empty())
1153  {
1154  // If no coupleGroup and no sampleRegion assume local region
1155  sampleRegion_ = patch_.boundaryMesh().mesh().name();
1156  sameRegion_ = true;
1157  }
1158  }
1159 }
1160 
1161 
1164  const polyPatch& pp,
1165  const mappedPatchBase& mpb
1166 )
1167 :
1168  patch_(pp),
1169  sampleRegion_(mpb.sampleRegion_),
1170  mode_(mpb.mode_),
1171  samplePatch_(mpb.samplePatch_),
1172  coupleGroup_(mpb.coupleGroup_),
1173  offsetMode_(mpb.offsetMode_),
1174  offset_(mpb.offset_),
1175  offsets_(mpb.offsets_),
1176  distance_(mpb.distance_),
1177  sameRegion_(mpb.sameRegion_),
1178  mapPtr_(nullptr),
1179  AMIPtr_(nullptr),
1180  AMIReverse_(mpb.AMIReverse_),
1181  surfPtr_(nullptr),
1182  surfDict_(mpb.surfDict_)
1183 {}
1184 
1185 
1188  const polyPatch& pp,
1189  const mappedPatchBase& mpb,
1190  const labelUList& mapAddressing
1191 )
1192 :
1193  patch_(pp),
1194  sampleRegion_(mpb.sampleRegion_),
1195  mode_(mpb.mode_),
1196  samplePatch_(mpb.samplePatch_),
1197  coupleGroup_(mpb.coupleGroup_),
1198  offsetMode_(mpb.offsetMode_),
1199  offset_(mpb.offset_),
1200  offsets_
1201  (
1202  offsetMode_ == NONUNIFORM
1203  ? vectorField(mpb.offsets_, mapAddressing)
1204  : vectorField(0)
1205  ),
1206  distance_(mpb.distance_),
1207  sameRegion_(mpb.sameRegion_),
1208  mapPtr_(nullptr),
1209  AMIPtr_(nullptr),
1210  AMIReverse_(mpb.AMIReverse_),
1211  surfPtr_(nullptr),
1212  surfDict_(mpb.surfDict_)
1213 {}
1214 
1215 
1216 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
1217 
1219 {
1220  clearOut();
1221 }
1222 
1223 
1225 {
1226  mapPtr_.clear();
1227  AMIPtr_.clear();
1228  surfPtr_.clear();
1229 }
1230 
1231 
1232 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1233 
1235 {
1236  return patch_.boundaryMesh().mesh().time().lookupObject<polyMesh>
1237  (
1238  sampleRegion()
1239  );
1240 }
1241 
1242 
1244 {
1245  const polyMesh& nbrMesh = sampleMesh();
1246 
1247  const label patchi = nbrMesh.boundaryMesh().findPatchID(samplePatch());
1248 
1249  if (patchi == -1)
1250  {
1252  << "Cannot find patch " << samplePatch()
1253  << " in region " << sampleRegion_ << endl
1254  << "Valid patches are " << nbrMesh.boundaryMesh().names()
1255  << exit(FatalError);
1256  }
1257 
1258  return nbrMesh.boundaryMesh()[patchi];
1259 }
1260 
1261 
1264  const pointField& fc
1265 ) const
1266 {
1267  tmp<pointField> tfld(new pointField(fc));
1268  pointField& fld = tfld.ref();
1269 
1270  switch (offsetMode_)
1271  {
1272  case UNIFORM:
1273  {
1274  fld += offset_;
1275  break;
1276  }
1277  case NONUNIFORM:
1278  {
1279  fld += offsets_;
1280  break;
1281  }
1282  case NORMAL:
1283  {
1284  // Get outwards pointing normal
1285  vectorField n(patch_.faceAreas());
1286  n /= mag(n);
1287 
1288  fld += distance_*n;
1289  break;
1290  }
1291  }
1292 
1293  return tfld;
1294 }
1295 
1296 
1298 {
1299  return samplePoints(facePoints(patch_));
1300 }
1301 
1302 
1305  const polyMesh& mesh,
1306  const label facei,
1307  const polyMesh::cellDecomposition decompMode
1308 )
1309 {
1310  const point& fc = mesh.faceCentres()[facei];
1311 
1312  switch (decompMode)
1313  {
1314  case polyMesh::FACE_PLANES:
1316  {
1317  // For both decompositions the face centre is guaranteed to be
1318  // on the face
1319  return pointIndexHit(true, fc, facei);
1320  }
1321  break;
1322 
1324  case polyMesh::CELL_TETS:
1325  {
1326  // Find the intersection of a ray from face centre to cell centre
1327  // Find intersection of (face-centre-decomposition) centre to
1328  // cell-centre with face-diagonal-decomposition triangles.
1329 
1330  const pointField& p = mesh.points();
1331  const face& f = mesh.faces()[facei];
1332 
1333  if (f.size() <= 3)
1334  {
1335  // Return centre of triangle.
1336  return pointIndexHit(true, fc, 0);
1337  }
1338 
1339  label celli = mesh.faceOwner()[facei];
1340  const point& cc = mesh.cellCentres()[celli];
1341  vector d = fc-cc;
1342 
1343  const label fp0 = mesh.tetBasePtIs()[facei];
1344  const point& basePoint = p[f[fp0]];
1345 
1346  label fp = f.fcIndex(fp0);
1347  for (label i = 2; i < f.size(); i++)
1348  {
1349  const point& thisPoint = p[f[fp]];
1350  label nextFp = f.fcIndex(fp);
1351  const point& nextPoint = p[f[nextFp]];
1352 
1353  const triPointRef tri(basePoint, thisPoint, nextPoint);
1354  pointHit hitInfo = tri.intersection
1355  (
1356  cc,
1357  d,
1359  );
1360 
1361  if (hitInfo.hit() && hitInfo.distance() > 0)
1362  {
1363  return pointIndexHit(true, hitInfo.hitPoint(), i-2);
1364  }
1365 
1366  fp = nextFp;
1367  }
1368 
1369  // Fall-back
1370  return pointIndexHit(false, fc, -1);
1371  }
1372  break;
1373 
1374  default:
1375  {
1377  << "problem" << abort(FatalError);
1378  return pointIndexHit();
1379  }
1380  }
1381 }
1382 
1383 
1385 {
1386  writeEntry(os, "sampleMode", sampleModeNames_[mode_]);
1387  if (!sampleRegion_.empty())
1388  {
1389  writeEntry(os, "sampleRegion", sampleRegion_);
1390  }
1391  if (!samplePatch_.empty())
1392  {
1393  writeEntry(os, "samplePatch", samplePatch_);
1394  }
1395  coupleGroup_.write(os);
1396 
1397  if
1398  (
1399  offsetMode_ == UNIFORM
1400  && offset_ == vector::zero
1401  && (mode_ == NEARESTPATCHFACE || mode_ == NEARESTPATCHFACEAMI)
1402  )
1403  {
1404  // Collocated mode. No need to write offset data
1405  }
1406  else
1407  {
1408  writeEntry(os, "offsetMode", offsetModeNames_[offsetMode_]);
1409 
1410  switch (offsetMode_)
1411  {
1412  case UNIFORM:
1413  {
1414  writeEntry(os, "offset", offset_);
1415  break;
1416  }
1417  case NONUNIFORM:
1418  {
1419  writeEntry(os, "offsets", offsets_);
1420  break;
1421  }
1422  case NORMAL:
1423  {
1424  writeEntry(os, "distance", distance_);
1425  break;
1426  }
1427  }
1428 
1429  if (mode_ == NEARESTPATCHFACEAMI)
1430  {
1431  if (AMIReverse_)
1432  {
1433  writeEntry(os, "flipNormals", AMIReverse_);
1434  }
1435 
1436  if (!surfDict_.empty())
1437  {
1438  writeKeyword(os, surfDict_.dictName());
1439  os << surfDict_;
1440  }
1441  }
1442  }
1443 }
1444 
1445 
1446 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:453
word samplePatch_
Patch (if in sampleMode NEARESTPATCH*)
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:82
static wordList words()
The set of names as a list of words.
Definition: NamedEnum.C:105
labelList first(const UList< labelPair > &p)
Definition: patchToPatch.C:38
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:663
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:434
bool isWord() const
Definition: tokenI.H:261
mode_t mode(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file mode.
Definition: POSIX.C:461
A triangle primitive used to calculate face areas and swept volumes.
Definition: triangle.H:57
vector offset_
Offset vector (uniform)
A class for handling file names.
Definition: fileName.H:79
Holds (reference to) pointField. Encapsulation of data needed for octree searches. Used for searching for nearest point. No bounding boxes around points. Only overlaps and calcNearest are implemented, rest makes little sense.
Definition: treeDataPoint.H:59
Ostream & writeKeyword(Foam::Ostream &os, const keyType &kw)
Write the keyword to the Ostream with the current level of indentation.
Definition: keyType.C:155
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
word sampleRegion_
Region to sample.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
MeshObject wrapper around meshSearch(mesh).
tmp< pointField > samplePoints() const
Get the sample points.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
A 2-tuple for storing two objects of different types.
Definition: HashTable.H:65
const polyMesh & sampleMesh() const
Get the region mesh.
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:928
Encapsulation of data needed to search for faces.
Definition: treeDataFace.H:58
bool sameRegion_
Same region.
const Type1 & first() const
Return first.
Definition: Tuple2.H:116
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:297
virtual void write(Ostream &) const
Write as a dictionary.
const word & wordToken() const
Definition: tokenI.H:266
Output to file stream.
Definition: OFstream.H:82
const coupleGroupIdentifier coupleGroup_
PatchGroup (if in sampleMode NEARESTPATCH*)
label findNearestFace(const point &location, const label seedFacei=-1, const bool useTreeSearch=true) const
Definition: meshSearch.C:725
virtual ~mappedPatchBase()
Destructor.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
A token holds items read from Istream.
Definition: token.H:72
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
void putBack(const token &)
Put back token.
Definition: Istream.C:30
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Traits class for primitives.
Definition: pTraits.H:50
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
label findPatchID(const word &patchName) const
Find patch index given a name.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
PointIndexHit< point > pointIndexHit
Definition: pointIndexHit.H:42
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
versionNumber version() const
Return the stream version.
Definition: IOstream.H:396
const Field< PointType > & localPoints() const
Return pointField of points in patch.
T & first()
Return the first element of the list.
Definition: UListI.H:114
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
cellDecomposition
Enumeration defining the decomposition of the cell for.
Definition: polyMesh.H:100
vectorField offsets_
Offset vector (nonuniform)
void setPoint(const Point &p)
Macros for easy insertion into run-time selection tables.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
scalar distance_
Offset distance (normal)
void collectSamples(const pointField &facePoints, pointField &, labelList &patchFaceProcs, labelList &patchFaces, pointField &patchFc) const
Collect single list of samples and originating processor+face.
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1211
scalarField samples(nIntervals, 0)
A list of faces which address into the list of points.
const Point & hitPoint() const
Return hit point.
A List obtained as a section of another List.
Definition: SubList.H:53
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
static const word & geometryDir()
Return the geometry directory name.
A class for handling words, derived from string.
Definition: word.H:59
label size() const
Return the number of elements in the list.
virtual const fileName & name() const
Return the name of the stream.
Definition: IOstream.H:294
const Point & hitPoint() const
Return hit point.
Definition: PointHit.H:126
wordList names() const
Return a list of patch names.
const word & constant() const
Return constant name.
Definition: TimePaths.H:123
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1237
static tmp< pointField > readListOrField(const word &keyword, const dictionary &dict, const label size)
Helper to read field or non-uniform list from dictionary.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
bool hit() const
Is there a hit.
const Field< PointType > & points() const
Return reference to global points.
void calcMapping() const
Calculate mapping.
const polyMesh & mesh() const
Return the mesh reference.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
bool hit() const
Is there a hit.
Definition: PointHit.H:120
Dynamically sized Field.
Definition: DynamicField.H:49
static const zero Zero
Definition: zero.H:97
offsetMode offsetMode_
How to obtain samples.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1224
const vectorField & cellCentres() const
errorManip< error > abort(error &err)
Definition: errorManip.H:131
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
tmp< pointField > facePoints(const polyPatch &) const
Get the points from face-centre-decomposition face centres.
const autoPtr< Foam::searchableSurface > & surfPtr() const
Return a pointer to the AMI projection surface.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
static const char nl
Definition: Ostream.H:260
label findInside(const point &) const
Find shape containing point. Only implemented for certain.
label facePoint(const int facei, const block &block, const label i, const label j)
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
static pointIndexHit facePoint(const polyMesh &, const label facei, const polyMesh::cellDecomposition)
Get a point on the face given a face decomposition method:
static meshSearchMeshObject & New(polyMesh &mesh)
Definition: MeshObject.C:54
labelList f(nPoints)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
offsetMode
How to project face centres.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
static const NamedEnum< sampleMode, 6 > sampleModeNames_
void findSamples(const sampleMode mode, const pointField &, labelList &sampleProcs, labelList &sampleIndices, pointField &sampleLocations) const
Find cells/faces containing samples.
dictionary surfDict_
Dictionary storing projection surface description.
const vectorField & faceCentres() const
void setSize(const label)
Reset size of List.
Definition: List.C:281
const indexedOctree< treeDataCell > & cellTree() const
Return the cell search tree.
Definition: polyMesh.C:962
Class containing processor-to-processor mapping information.
A bit-packed bool list.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
label patchi
void calcAMI() const
Calculate AMI interpolator.
vector point
Point is a vector.
Definition: point.H:41
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define WarningInFunction
Report a warning using Foam::Warning.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:318
Type gAverage(const FieldField< Field, Type > &f)
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
const List< FaceType > & localFaces() const
Return patch faces addressing into local point list.
A List with indirect addressing.
Definition: fvMatrix.H:106
sampleMode
Mesh items to sample.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:306
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
static autoPtr< searchableSurface > New(const word &surfaceType, const IOobject &io, const dictionary &dict)
Return a reference to the selected searchableSurface.
const bool AMIReverse_
Flag to indicate that slave patch should be reversed for AMI.
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
Field< vector > vectorField
Specialisation of Field<T> for vector.
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
static const NamedEnum< offsetMode, 3 > offsetModeNames_
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
pointHit intersection(const point &p, const vector &q, const intersection::algorithm alg, const scalar tol=0.0) const
Fast intersection with a ray.
Definition: triangleI.H:431
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
const polyPatch & samplePolyPatch() const
Get the patch on the region.
label index() const
Return index.
mappedPatchBase(const polyPatch &)
Construct from patch.
DynamicField< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Input token stream.
Definition: ITstream.H:49
dictionary subOrEmptyDict(const word &, const bool mustRead=false) const
Find and return a sub-dictionary as a copy, or.
Definition: dictionary.C:1033
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864
const sampleMode mode_
What to sample.
IOerror FatalIOError
#define InfoInFunction
Report an information message using Foam::Info.