meshToMesh.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) 2012-2019 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 "meshToMesh.H"
27 #include "Time.H"
28 #include "globalIndex.H"
29 #include "meshToMeshMethod.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(meshToMesh, 0);
36 
37  template<>
38  const char* Foam::NamedEnum
39  <
41  3
42  >::names[] =
43  {
44  "direct",
45  "mapNearest",
46  "cellVolumeWeight"
47  };
48 
51 }
52 
53 
54 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
55 
56 template<>
57 void Foam::meshToMesh::mapAndOpSrcToTgt
58 (
59  const AMIInterpolation& AMI,
60  const Field<scalar>& srcField,
61  Field<scalar>& tgtField,
62  const plusEqOp<scalar>& cop
63 ) const
64 {}
65 
66 
67 template<>
68 void Foam::meshToMesh::mapAndOpSrcToTgt
69 (
70  const AMIInterpolation& AMI,
71  const Field<vector>& srcField,
72  Field<vector>& tgtField,
73  const plusEqOp<vector>& cop
74 ) const
75 {}
76 
77 
78 template<>
79 void Foam::meshToMesh::mapAndOpSrcToTgt
80 (
81  const AMIInterpolation& AMI,
82  const Field<sphericalTensor>& srcField,
83  Field<sphericalTensor>& tgtField,
84  const plusEqOp<sphericalTensor>& cop
85 ) const
86 {}
87 
88 
89 template<>
90 void Foam::meshToMesh::mapAndOpSrcToTgt
91 (
92  const AMIInterpolation& AMI,
93  const Field<symmTensor>& srcField,
94  Field<symmTensor>& tgtField,
95  const plusEqOp<symmTensor>& cop
96 ) const
97 {}
98 
99 
100 template<>
101 void Foam::meshToMesh::mapAndOpSrcToTgt
102 (
103  const AMIInterpolation& AMI,
104  const Field<tensor>& srcField,
105  Field<tensor>& tgtField,
106  const plusEqOp<tensor>& cop
107 ) const
108 {}
109 
110 
111 template<>
112 void Foam::meshToMesh::mapAndOpTgtToSrc
113 (
114  const AMIInterpolation& AMI,
115  Field<scalar>& srcField,
116  const Field<scalar>& tgtField,
117  const plusEqOp<scalar>& cop
118 ) const
119 {}
120 
121 
122 template<>
123 void Foam::meshToMesh::mapAndOpTgtToSrc
124 (
125  const AMIInterpolation& AMI,
126  Field<vector>& srcField,
127  const Field<vector>& tgtField,
128  const plusEqOp<vector>& cop
129 ) const
130 {}
131 
132 
133 template<>
134 void Foam::meshToMesh::mapAndOpTgtToSrc
135 (
136  const AMIInterpolation& AMI,
137  Field<sphericalTensor>& srcField,
138  const Field<sphericalTensor>& tgtField,
139  const plusEqOp<sphericalTensor>& cop
140 ) const
141 {}
142 
143 
144 template<>
145 void Foam::meshToMesh::mapAndOpTgtToSrc
146 (
147  const AMIInterpolation& AMI,
148  Field<symmTensor>& srcField,
149  const Field<symmTensor>& tgtField,
150  const plusEqOp<symmTensor>& cop
151 ) const
152 {}
153 
154 
155 template<>
156 void Foam::meshToMesh::mapAndOpTgtToSrc
157 (
158  const AMIInterpolation& AMI,
159  Field<tensor>& srcField,
160  const Field<tensor>& tgtField,
161  const plusEqOp<tensor>& cop
162 ) const
163 {}
164 
165 
166 Foam::labelList Foam::meshToMesh::maskCells
167 (
168  const polyMesh& src,
169  const polyMesh& tgt
170 ) const
171 {
172  boundBox intersectBb
173  (
174  max(src.bounds().min(), tgt.bounds().min()),
175  min(src.bounds().max(), tgt.bounds().max())
176  );
177 
178  intersectBb.inflate(0.01);
179 
180  const cellList& srcCells = src.cells();
181  const faceList& srcFaces = src.faces();
182  const pointField& srcPts = src.points();
183 
185  forAll(srcCells, srcI)
186  {
187  boundBox cellBb(srcCells[srcI].points(srcFaces, srcPts), false);
188  if (intersectBb.overlaps(cellBb))
189  {
190  cells.append(srcI);
191  }
192  }
193 
194  if (debug)
195  {
196  Pout<< "participating source mesh cells: " << cells.size() << endl;
197  }
198 
199  return move(cells);
200 }
201 
202 
203 void Foam::meshToMesh::normaliseWeights
204 (
205  const word& descriptor,
206  const labelListList& addr,
207  scalarListList& wght
208 ) const
209 {
210  const label nCell = returnReduce(wght.size(), sumOp<label>());
211 
212  if (nCell > 0)
213  {
214  forAll(wght, celli)
215  {
216  scalarList& w = wght[celli];
217  scalar s = sum(w);
218 
219  forAll(w, i)
220  {
221  // note: normalise by s instead of cell volume since
222  // 1-to-1 methods duplicate contributions in parallel
223  w[i] /= s;
224  }
225  }
226  }
227 }
228 
229 
230 void Foam::meshToMesh::calcAddressing
231 (
232  const word& methodName,
233  const polyMesh& src,
234  const polyMesh& tgt
235 )
236 {
237  autoPtr<meshToMeshMethod> methodPtr
238  (
240  (
241  methodName,
242  src,
243  tgt
244  )
245  );
246 
247  methodPtr->calculate
248  (
249  srcToTgtCellAddr_,
250  srcToTgtCellWght_,
251  tgtToSrcCellAddr_,
252  tgtToSrcCellWght_
253  );
254 
255  V_ = methodPtr->V();
256 
257  if (debug)
258  {
259  methodPtr->writeConnectivity(src, tgt, srcToTgtCellAddr_);
260  }
261 }
262 
263 
264 void Foam::meshToMesh::calculate(const word& methodName)
265 {
266  Info<< "Creating mesh-to-mesh addressing for " << srcRegion_.name()
267  << " and " << tgtRegion_.name() << " regions using "
268  << methodName << endl;
269 
270  singleMeshProc_ = calcDistribution(srcRegion_, tgtRegion_);
271 
272  if (singleMeshProc_ == -1)
273  {
274  // create global indexing for src and tgt meshes
275  globalIndex globalSrcCells(srcRegion_.nCells());
276  globalIndex globalTgtCells(tgtRegion_.nCells());
277 
278  // Create processor map of overlapping cells. This map gets
279  // (possibly remote) cells from the tgt mesh such that they (together)
280  // cover all of the src mesh
281  autoPtr<mapDistribute> mapPtr = calcProcMap(srcRegion_, tgtRegion_);
282  const mapDistribute& map = mapPtr();
283 
284  pointField newTgtPoints;
285  faceList newTgtFaces;
286  labelList newTgtFaceOwners;
287  labelList newTgtFaceNeighbours;
288  labelList newTgtCellIDs;
289 
290  distributeAndMergeCells
291  (
292  map,
293  tgtRegion_,
294  globalTgtCells,
295  newTgtPoints,
296  newTgtFaces,
297  newTgtFaceOwners,
298  newTgtFaceNeighbours,
299  newTgtCellIDs
300  );
301 
302 
303  // create a new target mesh
304  polyMesh newTgt
305  (
306  IOobject
307  (
308  "newTgt." + Foam::name(Pstream::myProcNo()),
309  tgtRegion_.time().timeName(),
310  tgtRegion_.time(),
312  ),
313  move(newTgtPoints),
314  move(newTgtFaces),
315  move(newTgtFaceOwners),
316  move(newTgtFaceNeighbours),
317  false // no parallel comms
318  );
319 
320  // create some dummy patch info
322  patches[0] = new polyPatch
323  (
324  "defaultFaces",
325  newTgt.nFaces() - newTgt.nInternalFaces(),
326  newTgt.nInternalFaces(),
327  0,
328  newTgt.boundaryMesh(),
329  word::null
330  );
331 
332  newTgt.addPatches(patches);
333 
334  // force calculation of tet-base points used for point-in-cell
335  (void)newTgt.tetBasePtIs();
336 
337  // force construction of cell tree
338 // (void)newTgt.cellTree();
339 
340  if (debug)
341  {
342  Pout<< "Created newTgt mesh:" << nl
343  << " old cells = " << tgtRegion_.nCells()
344  << ", new cells = " << newTgt.nCells() << nl
345  << " old faces = " << tgtRegion_.nFaces()
346  << ", new faces = " << newTgt.nFaces() << endl;
347 
348  if (debug > 1)
349  {
350  Pout<< "Writing newTgt mesh: " << newTgt.name() << endl;
351  newTgt.write();
352  }
353  }
354 
355  calcAddressing(methodName, srcRegion_, newTgt);
356 
357  // per source cell the target cell address in newTgt mesh
358  forAll(srcToTgtCellAddr_, i)
359  {
360  labelList& addressing = srcToTgtCellAddr_[i];
361  forAll(addressing, addrI)
362  {
363  addressing[addrI] = newTgtCellIDs[addressing[addrI]];
364  }
365  }
366 
367  // convert target addresses in newTgtMesh into global cell numbering
368  forAll(tgtToSrcCellAddr_, i)
369  {
370  labelList& addressing = tgtToSrcCellAddr_[i];
371  forAll(addressing, addrI)
372  {
373  addressing[addrI] = globalSrcCells.toGlobal(addressing[addrI]);
374  }
375  }
376 
377  // set up as a reverse distribute
379  (
381  List<labelPair>(),
382  tgtRegion_.nCells(),
383  map.constructMap(),
384  false,
385  map.subMap(),
386  false,
387  tgtToSrcCellAddr_,
389  flipOp(),
390  labelList()
391  );
392 
393  // set up as a reverse distribute
395  (
397  List<labelPair>(),
398  tgtRegion_.nCells(),
399  map.constructMap(),
400  false,
401  map.subMap(),
402  false,
403  tgtToSrcCellWght_,
405  flipOp(),
406  scalarList()
407  );
408 
409  // weights normalisation
410  normaliseWeights
411  (
412  "source",
413  srcToTgtCellAddr_,
414  srcToTgtCellWght_
415  );
416 
417  normaliseWeights
418  (
419  "target",
420  tgtToSrcCellAddr_,
421  tgtToSrcCellWght_
422  );
423 
424  // cache maps and reset addresses
425  List<Map<label>> cMap;
426  srcMapPtr_.reset
427  (
428  new mapDistribute(globalSrcCells, tgtToSrcCellAddr_, cMap)
429  );
430  tgtMapPtr_.reset
431  (
432  new mapDistribute(globalTgtCells, srcToTgtCellAddr_, cMap)
433  );
434 
435  // collect volume intersection contributions
436  reduce(V_, sumOp<scalar>());
437  }
438  else
439  {
440  calcAddressing(methodName, srcRegion_, tgtRegion_);
441 
442  normaliseWeights
443  (
444  "source",
445  srcToTgtCellAddr_,
446  srcToTgtCellWght_
447  );
448 
449  normaliseWeights
450  (
451  "target",
452  tgtToSrcCellAddr_,
453  tgtToSrcCellWght_
454  );
455  }
456 
457  Info<< " Overlap volume: " << V_ << endl;
458 }
459 
460 
463 {
464  switch (method)
465  {
466  case imDirect:
467  {
469  break;
470  }
471  case imMapNearest:
472  {
474  break;
475  }
476  case imCellVolumeWeight:
477  {
479  break;
480  }
481  default:
482  {
484  << "Unhandled enumeration " << method
485  << abort(FatalError);
486  }
487  }
488 
490 }
491 
492 
493 void Foam::meshToMesh::calculatePatchAMIs(const word& AMIMethodName)
494 {
495  if (!patchAMIs_.empty())
496  {
498  << "patch AMI already calculated"
499  << exit(FatalError);
500  }
501 
502  patchAMIs_.setSize(srcPatchID_.size());
503 
504  forAll(srcPatchID_, i)
505  {
506  label srcPatchi = srcPatchID_[i];
507  label tgtPatchi = tgtPatchID_[i];
508 
509  const polyPatch& srcPP = srcRegion_.boundaryMesh()[srcPatchi];
510  const polyPatch& tgtPP = tgtRegion_.boundaryMesh()[tgtPatchi];
511 
512  Info<< "Creating AMI between source patch " << srcPP.name()
513  << " and target patch " << tgtPP.name()
514  << " using " << AMIMethodName
515  << endl;
516 
517  Info<< incrIndent;
518 
519  patchAMIs_.set
520  (
521  i,
522  new AMIInterpolation
523  (
524  srcPP,
525  tgtPP,
527  false,
528  AMIMethodName,
529  -1,
530  true // flip target patch since patch normals are aligned
531  )
532  );
533 
534  Info<< decrIndent;
535  }
536 }
537 
538 
539 void Foam::meshToMesh::constructNoCuttingPatches
540 (
541  const word& methodName,
542  const word& AMIMethodName,
543  const bool interpAllPatches
544 )
545 {
546  if (interpAllPatches)
547  {
548  const polyBoundaryMesh& srcBM = srcRegion_.boundaryMesh();
549  const polyBoundaryMesh& tgtBM = tgtRegion_.boundaryMesh();
550 
551  DynamicList<label> srcPatchID(srcBM.size());
552  DynamicList<label> tgtPatchID(tgtBM.size());
553  forAll(srcBM, patchi)
554  {
555  const polyPatch& pp = srcBM[patchi];
556 
557  // We want to map all the global patches, including constraint
558  // patches (since they might have mappable properties, e.g.
559  // jumpCyclic). We'll fix the value afterwards.
560  if (!isA<processorPolyPatch>(pp))
561  {
562  srcPatchID.append(pp.index());
563 
564  label tgtPatchi = tgtBM.findPatchID(pp.name());
565 
566  if (tgtPatchi != -1)
567  {
568  tgtPatchID.append(tgtPatchi);
569  }
570  else
571  {
573  << "Source patch " << pp.name()
574  << " not found in target mesh. "
575  << "Available target patches are " << tgtBM.names()
576  << exit(FatalError);
577  }
578  }
579  }
580 
581  srcPatchID_.transfer(srcPatchID);
582  tgtPatchID_.transfer(tgtPatchID);
583  }
584 
585  // calculate volume addressing and weights
586  calculate(methodName);
587 
588  // calculate patch addressing and weights
589  calculatePatchAMIs(AMIMethodName);
590 }
591 
592 
593 void Foam::meshToMesh::constructFromCuttingPatches
594 (
595  const word& methodName,
596  const word& AMIMethodName,
597  const HashTable<word>& patchMap,
598  const wordList& cuttingPatches
599 )
600 {
601  srcPatchID_.setSize(patchMap.size());
602  tgtPatchID_.setSize(patchMap.size());
603 
604  label i = 0;
605  forAllConstIter(HashTable<word>, patchMap, iter)
606  {
607  const word& tgtPatchName = iter.key();
608  const word& srcPatchName = iter();
609 
610  const polyPatch& srcPatch = srcRegion_.boundaryMesh()[srcPatchName];
611  const polyPatch& tgtPatch = tgtRegion_.boundaryMesh()[tgtPatchName];
612 
613  srcPatchID_[i] = srcPatch.index();
614  tgtPatchID_[i] = tgtPatch.index();
615  i++;
616  }
617 
618  // calculate volume addressing and weights
619  calculate(methodName);
620 
621  // calculate patch addressing and weights
622  calculatePatchAMIs(AMIMethodName);
623 
624  // set IDs of cutting patches on target mesh
625  cuttingPatches_.setSize(cuttingPatches.size());
626  forAll(cuttingPatches_, i)
627  {
628  const word& patchName = cuttingPatches[i];
629  cuttingPatches_[i] = tgtRegion_.boundaryMesh().findPatchID(patchName);
630  }
631 }
632 
633 
634 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
635 
637 (
638  const polyMesh& src,
639  const polyMesh& tgt,
640  const interpolationMethod& method,
641  bool interpAllPatches
642 )
643 :
644  srcRegion_(src),
645  tgtRegion_(tgt),
646  srcPatchID_(),
647  tgtPatchID_(),
648  patchAMIs_(),
649  cuttingPatches_(),
650  srcToTgtCellAddr_(),
651  tgtToSrcCellAddr_(),
652  srcToTgtCellWght_(),
653  tgtToSrcCellWght_(),
654  V_(0.0),
655  singleMeshProc_(-1),
656  srcMapPtr_(nullptr),
657  tgtMapPtr_(nullptr)
658 {
659  constructNoCuttingPatches
660  (
661  interpolationMethodNames_[method],
663  (
664  interpolationMethodAMI(method)
665  ),
666  interpAllPatches
667  );
668 }
669 
670 
672 (
673  const polyMesh& src,
674  const polyMesh& tgt,
675  const word& methodName,
676  const word& AMIMethodName,
677  bool interpAllPatches
678 )
679 :
680  srcRegion_(src),
681  tgtRegion_(tgt),
682  srcPatchID_(),
683  tgtPatchID_(),
684  patchAMIs_(),
685  cuttingPatches_(),
686  srcToTgtCellAddr_(),
687  tgtToSrcCellAddr_(),
688  srcToTgtCellWght_(),
689  tgtToSrcCellWght_(),
690  V_(0.0),
691  singleMeshProc_(-1),
692  srcMapPtr_(nullptr),
693  tgtMapPtr_(nullptr)
694 {
695  constructNoCuttingPatches(methodName, AMIMethodName, interpAllPatches);
696 }
697 
698 
700 (
701  const polyMesh& src,
702  const polyMesh& tgt,
703  const interpolationMethod& method,
704  const HashTable<word>& patchMap,
705  const wordList& cuttingPatches
706 )
707 :
708  srcRegion_(src),
709  tgtRegion_(tgt),
710  srcPatchID_(),
711  tgtPatchID_(),
712  patchAMIs_(),
713  cuttingPatches_(),
714  srcToTgtCellAddr_(),
715  tgtToSrcCellAddr_(),
716  srcToTgtCellWght_(),
717  tgtToSrcCellWght_(),
718  V_(0.0),
719  singleMeshProc_(-1),
720  srcMapPtr_(nullptr),
721  tgtMapPtr_(nullptr)
722 {
723  constructFromCuttingPatches
724  (
725  interpolationMethodNames_[method],
727  (
728  interpolationMethodAMI(method)
729  ),
730  patchMap,
731  cuttingPatches
732  );
733 }
734 
735 
737 (
738  const polyMesh& src,
739  const polyMesh& tgt,
740  const word& methodName, // internal mapping
741  const word& AMIMethodName, // boundary mapping
742  const HashTable<word>& patchMap,
743  const wordList& cuttingPatches
744 )
745 :
746  srcRegion_(src),
747  tgtRegion_(tgt),
748  srcPatchID_(),
749  tgtPatchID_(),
750  patchAMIs_(),
751  cuttingPatches_(),
752  srcToTgtCellAddr_(),
753  tgtToSrcCellAddr_(),
754  srcToTgtCellWght_(),
755  tgtToSrcCellWght_(),
756  V_(0.0),
757  singleMeshProc_(-1),
758  srcMapPtr_(nullptr),
759  tgtMapPtr_(nullptr)
760 {
761  constructFromCuttingPatches
762  (
763  methodName,
764  AMIMethodName,
765  patchMap,
766  cuttingPatches
767  );
768 }
769 
770 
771 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
772 
774 {}
775 
776 
777 // ************************************************************************* //
static word interpolationMethodToWord(const interpolationMethod &method)
Convert interpolationMethod to word representation.
Class containing functor to negate primitives. Dummy for all other types.
Definition: flipOp.H:50
bool overlaps(const boundBox &) const
Overlaps/touches boundingBox?
Definition: boundBoxI.H:120
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:82
void append(T *)
Append an element at the end of the list.
Definition: PtrListI.H:39
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual Ostream & write(const char)=0
Write character.
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
const word & name() const
Return name.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:278
const labelListList & subMap() const
From subsetted data back to original data.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Plus op for FixedList<scalar>
static const NamedEnum< interpolationMethod, 3 > interpolationMethodNames_
Definition: meshToMesh.H:75
label findPatchID(const word &patchName) const
Find patch index given a name.
scalar V() const
Return const access to the overlap volume.
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:58
const cellList & cells() const
patches[0]
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
void writeConnectivity(const polyMesh &mesh1, const polyMesh &mesh2, const labelListList &mesh1ToMesh2Addr) const
Write the connectivity (debugging)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
interpolationMethod
Enumeration specifying interpolation method.
virtual void calculate(labelListList &srcToTgtAddr, scalarListList &srcToTgtWght, labelListList &tgtToTgtAddr, scalarListList &tgtToTgtWght)=0
Calculate addressing and weights.
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
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){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const cellShapeList & cells
const pointField & points
A class for handling words, derived from string.
Definition: word.H:59
wordList names() const
Return a list of patch names.
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
static void distribute(const Pstream::commsTypes commsType, const List< labelPair > &schedule, const label constructSize, const labelListList &subMap, const bool subHasFlip, const labelListList &constructMap, const bool constructHasFlip, List< T > &, const negateOp &negOp, const int tag=UPstream::msgType())
Distribute data. Note:schedule only used for.
static const word null
An empty word.
Definition: word.H:77
List< label > labelList
A List of labels.
Definition: labelList.H:56
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
An STL-conforming hash table.
Definition: HashTable.H:61
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Foam::polyBoundaryMesh.
void inflate(const scalar s)
Inflate box by factor*mag(span) in all dimensions.
Definition: boundBox.C:210
static AMIInterpolation::interpolationMethod interpolationMethodAMI(const interpolationMethod method)
Conversion between mesh and patch interpolation methods.
Definition: meshToMesh.C:462
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:235
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
const point & max() const
Maximum describing the bounding box.
Definition: boundBoxI.H:60
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:440
meshToMesh(const polyMesh &src, const polyMesh &tgt, const interpolationMethod &method, const bool interpAllPatches=true)
Construct from source and target meshes.
Definition: meshToMesh.C:637
label patchi
Class containing processor-to-processor mapping information.
virtual ~meshToMesh()
Destructor.
Definition: meshToMesh.C:773
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
label index() const
Return the index of this patch in the boundaryMesh.
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
messageStream Info
const labelListList & constructMap() const
From subsetted data to new reconstructed data.
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:52
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const point & min() const
Minimum describing the bounding box.
Definition: boundBoxI.H:54
interpolationMethod
Enumeration specifying interpolation method.
Definition: meshToMesh.H:67
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:228
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
static autoPtr< meshToMeshMethod > New(const word &methodName, const polyMesh &src, const polyMesh &tgt)
Selector.
Namespace for OpenFOAM.