processorPolyPatch.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 "processorPolyPatch.H"
28 #include "dictionary.H"
29 #include "SubField.H"
30 #include "demandDrivenData.H"
31 #include "matchPoints.H"
32 #include "OFstream.H"
33 #include "polyMesh.H"
34 #include "Time.H"
35 #include "PstreamBuffers.H"
36 #include "ConstCirculator.H"
37 
38 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
39 
40 namespace Foam
41 {
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
50 (
51  const word& name,
52  const label size,
53  const label start,
54  const label index,
55  const polyBoundaryMesh& bm,
56  const int myProcNo,
57  const int neighbProcNo,
58  const word& patchType
59 )
60 :
61  coupledPolyPatch(name, size, start, index, bm, patchType),
62  myProcNo_(myProcNo),
63  neighbProcNo_(neighbProcNo),
64  neighbFaceCentres_(),
65  neighbFaceAreas_(),
66  neighbFaceCellCentres_()
67 {}
68 
69 
71 (
72  const label size,
73  const label start,
74  const label index,
75  const polyBoundaryMesh& bm,
76  const int myProcNo,
77  const int neighbProcNo,
78  const word& patchType
79 )
80 :
82  (
83  newName(myProcNo, neighbProcNo),
84  size,
85  start,
86  index,
87  bm,
88  patchType
89  ),
90  myProcNo_(myProcNo),
91  neighbProcNo_(neighbProcNo),
92  neighbFaceCentres_(),
93  neighbFaceAreas_(),
94  neighbFaceCellCentres_()
95 {}
96 
97 
99 (
100  const word& name,
101  const dictionary& dict,
102  const label index,
103  const polyBoundaryMesh& bm,
104  const word& patchType
105 )
106 :
107  coupledPolyPatch(name, dict, index, bm, patchType),
108  myProcNo_(dict.lookup<label>("myProcNo")),
109  neighbProcNo_(dict.lookup<label>("neighbProcNo")),
110  neighbFaceCentres_(),
111  neighbFaceAreas_(),
112  neighbFaceCellCentres_()
113 {}
114 
115 
117 (
118  const processorPolyPatch& pp,
119  const polyBoundaryMesh& bm
120 )
121 :
122  coupledPolyPatch(pp, bm),
123  myProcNo_(pp.myProcNo_),
124  neighbProcNo_(pp.neighbProcNo_),
125  neighbFaceCentres_(),
126  neighbFaceAreas_(),
127  neighbFaceCellCentres_()
128 {}
129 
130 
132 (
133  const processorPolyPatch& pp,
134  const polyBoundaryMesh& bm,
135  const label index,
136  const label newSize,
137  const label newStart
138 )
139 :
140  coupledPolyPatch(pp, bm, index, newSize, newStart),
141  myProcNo_(pp.myProcNo_),
142  neighbProcNo_(pp.neighbProcNo_),
143  neighbFaceCentres_(),
144  neighbFaceAreas_(),
145  neighbFaceCellCentres_()
146 {}
147 
148 
149 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
150 
152 {
153  nbrPointsPtr_.clear();
154  nbrEdgesPtr_.clear();
155 }
156 
157 
158 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
159 
161 (
162  const label myProcNo,
163  const label neighbProcNo
164 )
165 {
166  return
167  "procBoundary"
168  + Foam::name(myProcNo)
169  + "to"
170  + Foam::name(neighbProcNo);
171 }
172 
173 
175 {
176  if (Pstream::parRun())
177  {
178  UOPstream toNeighbProc(neighbProcNo(), pBufs);
179 
180  toNeighbProc
181  << faceCentres()
182  << faceAreas()
183  << faceCellCentres();
184  }
185 }
186 
187 
189 {
190  if (Pstream::parRun())
191  {
192  {
193  UIPstream fromNeighbProc(neighbProcNo(), pBufs);
194 
195  fromNeighbProc
196  >> neighbFaceCentres_
197  >> neighbFaceAreas_
198  >> neighbFaceCellCentres_;
199  }
200 
201  // My normals
202  vectorField faceNormals(size());
203 
204  // Neighbour normals
205  vectorField nbrFaceNormals(neighbFaceAreas_.size());
206 
207  // Face match tolerances
208  scalarField tols = calcFaceTol(*this, points(), faceCentres());
209 
210  // Calculate normals from areas and check
211  forAll(faceNormals, facei)
212  {
213  scalar magSf = magFaceAreas()[facei];
214  scalar nbrMagSf = mag(neighbFaceAreas_[facei]);
215  scalar avSf = (magSf + nbrMagSf)/2.0;
216 
217  // For small face area calculation the results of the area
218  // calculation have been found to only be accurate to ~1e-20
219  if (magSf < small || nbrMagSf < small)
220  {
221  // Undetermined normal. Use dummy normal to force separation
222  // check.
223  faceNormals[facei] = point(1, 0, 0);
224  nbrFaceNormals[facei] = -faceNormals[facei];
225  tols[facei] = great;
226  }
227  else if (mag(magSf - nbrMagSf) > matchTolerance()*sqr(tols[facei]))
228  {
229  const fileName patchOBJName
230  (
231  boundaryMesh().mesh().time().path()/name() + "_faces.obj"
232  );
233 
234  Pout<< "processorPolyPatch::calcGeometry : Writing my "
235  << size() << " faces to " << patchOBJName << endl;
236 
237  writeOBJ(patchOBJName, *this);
238 
239  const fileName centresOBJName
240  (
241  boundaryMesh().mesh().time().path()/name()
242  + "_faceCentresConnections.obj"
243  );
244 
245  Pout<< "processorPolyPatch::calcGeometry :"
246  << " Dumping lines between corresponding face centres to "
247  << centresOBJName.name() << endl;
248 
249  writeOBJ(centresOBJName, neighbFaceCentres_, faceCentres());
250 
252  << "face " << facei << " area does not match neighbour by "
253  << 100*mag(magSf - nbrMagSf)/avSf
254  << "% -- possible face ordering problem." << endl
255  << "patch:" << name()
256  << " my area:" << magSf
257  << " neighbour area:" << nbrMagSf
258  << " matching tolerance:"
259  << matchTolerance()*sqr(tols[facei])
260  << endl
261  << "Mesh face:" << start()+facei
262  << " vertices:"
263  << UIndirectList<point>(points(), operator[](facei))()
264  << endl
265  << "If you are certain your matching is correct"
266  << " you can increase the 'matchTolerance' setting"
267  << " in the patch dictionary in the boundary file."
268  << endl
269  << "Rerun with processor debug flag set for"
270  << " more information." << exit(FatalError);
271  }
272  else
273  {
274  faceNormals[facei] = faceAreas()[facei]/magSf;
275  nbrFaceNormals[facei] = neighbFaceAreas_[facei]/nbrMagSf;
276  }
277  }
278  }
279 }
280 
281 
283 (
284  PstreamBuffers& pBufs,
285  const pointField& p
286 )
287 {
288  polyPatch::movePoints(pBufs, p);
290 }
291 
292 
294 (
295  PstreamBuffers& pBufs,
296  const pointField&
297 )
298 {
300 }
301 
302 
304 {
306 
307  if (Pstream::parRun())
308  {
309  // Express all points as patch face and index in face.
310  labelList pointFace(nPoints());
311  labelList pointIndex(nPoints());
312 
313  for (label patchPointi = 0; patchPointi < nPoints(); patchPointi++)
314  {
315  label facei = pointFaces()[patchPointi][0];
316 
317  pointFace[patchPointi] = facei;
318 
319  const face& f = localFaces()[facei];
320 
321  pointIndex[patchPointi] = findIndex(f, patchPointi);
322  }
323 
324  // Express all edges as patch face and index in face.
325  labelList edgeFace(nEdges());
326  labelList edgeIndex(nEdges());
327 
328  for (label patchEdgeI = 0; patchEdgeI < nEdges(); patchEdgeI++)
329  {
330  label facei = edgeFaces()[patchEdgeI][0];
331 
332  edgeFace[patchEdgeI] = facei;
333 
334  const labelList& fEdges = faceEdges()[facei];
335 
336  edgeIndex[patchEdgeI] = findIndex(fEdges, patchEdgeI);
337  }
338 
339  UOPstream toNeighbProc(neighbProcNo(), pBufs);
340 
341  toNeighbProc
342  << pointFace
343  << pointIndex
344  << edgeFace
345  << edgeIndex;
346  }
347 }
348 
349 
351 {
352  // For completeness
353  polyPatch::topoChange(pBufs);
354 
355  nbrPointsPtr_.clear();
356  nbrEdgesPtr_.clear();
357 
358  if (Pstream::parRun())
359  {
360  labelList nbrPointFace;
361  labelList nbrPointIndex;
362  labelList nbrEdgeFace;
363  labelList nbrEdgeIndex;
364 
365  {
366  // !Note: there is one situation where the opposite points and
367  // do not exactly match and that is while we are distributing
368  // meshes and multiple parts come together from different
369  // processors. This can temporarily create the situation that
370  // we have points which will be merged out later but we still
371  // need the face connectivity to be correct.
372  // So: cannot check here on same points and edges.
373 
374  UIPstream fromNeighbProc(neighbProcNo(), pBufs);
375 
376  fromNeighbProc
377  >> nbrPointFace
378  >> nbrPointIndex
379  >> nbrEdgeFace
380  >> nbrEdgeIndex;
381  }
382 
383  // Convert neighbour faces and indices into face back into
384  // my edges and points.
385 
386  // Convert points.
387  // ~~~~~~~~~~~~~~~
388 
389  nbrPointsPtr_.reset(new labelList(nPoints(), -1));
390  labelList& nbrPoints = nbrPointsPtr_();
391 
392  forAll(nbrPointFace, nbrPointi)
393  {
394  // Find face and index in face on this side.
395  const face& f = localFaces()[nbrPointFace[nbrPointi]];
396 
397  label index = (f.size() - nbrPointIndex[nbrPointi]) % f.size();
398  label patchPointi = f[index];
399 
400  if (nbrPoints[patchPointi] == -1)
401  {
402  // First reference of point
403  nbrPoints[patchPointi] = nbrPointi;
404  }
405  else if (nbrPoints[patchPointi] >= 0)
406  {
407  // Point already visited. Mark as duplicate.
408  nbrPoints[patchPointi] = -2;
409  }
410  }
411 
412  // Reset all duplicate entries to -1.
413  forAll(nbrPoints, patchPointi)
414  {
415  if (nbrPoints[patchPointi] == -2)
416  {
417  nbrPoints[patchPointi] = -1;
418  }
419  }
420 
421  // Convert edges.
422  // ~~~~~~~~~~~~~~
423 
424  nbrEdgesPtr_.reset(new labelList(nEdges(), -1));
425  labelList& nbrEdges = nbrEdgesPtr_();
426 
427  forAll(nbrEdgeFace, nbrEdgeI)
428  {
429  // Find face and index in face on this side.
430  const labelList& f = faceEdges()[nbrEdgeFace[nbrEdgeI]];
431  label index = (f.size() - nbrEdgeIndex[nbrEdgeI] - 1) % f.size();
432  label patchEdgeI = f[index];
433 
434  if (nbrEdges[patchEdgeI] == -1)
435  {
436  // First reference of edge
437  nbrEdges[patchEdgeI] = nbrEdgeI;
438  }
439  else if (nbrEdges[patchEdgeI] >= 0)
440  {
441  // Edge already visited. Mark as duplicate.
442  nbrEdges[patchEdgeI] = -2;
443  }
444  }
445 
446  // Reset all duplicate entries to -1.
447  forAll(nbrEdges, patchEdgeI)
448  {
449  if (nbrEdges[patchEdgeI] == -2)
450  {
451  nbrEdges[patchEdgeI] = -1;
452  }
453  }
454 
455  // Remove any addressing used for shared points/edges calculation
456  // since mostly not needed.
458  }
459 }
460 
461 
463 {
464  if (!nbrPointsPtr_.valid())
465  {
467  << "No extended addressing calculated for patch " << name()
468  << abort(FatalError);
469  }
470  return nbrPointsPtr_();
471 }
472 
473 
475 {
476  if (!nbrEdgesPtr_.valid())
477  {
479  << "No extended addressing calculated for patch " << name()
480  << abort(FatalError);
481  }
482  return nbrEdgesPtr_();
483 }
484 
485 
487 (
488  PstreamBuffers& pBufs,
489  const primitivePatch& pp
490 ) const
491 {
492  if (!Pstream::parRun())
493  {
494  return;
495  }
496 
497  if (owner())
498  {
499  ownToNbrOrderData ownToNbr;
500  autoPtr<ownToNbrDebugOrderData> ownToNbrDebugPtr
501  (
502  coupledPolyPatch::debug
503  ? new ownToNbrDebugOrderData()
504  : nullptr
505  );
506 
508  (
509  ownToNbr,
510  ownToNbrDebugPtr,
511  pp
512  );
513 
514  UOPstream toNeighbour(neighbProcNo(), pBufs);
515  toNeighbour << ownToNbr;
516  if (coupledPolyPatch::debug)
517  {
518  toNeighbour << ownToNbrDebugPtr();
519  }
520  }
521 }
522 
523 
525 (
526  PstreamBuffers& pBufs,
527  const primitivePatch& pp,
529  labelList& rotation
530 ) const
531 {
532  if (!Pstream::parRun())
533  {
534  return false;
535  }
536 
537  ownToNbrOrderData ownToNbr;
538  autoPtr<ownToNbrDebugOrderData> ownToNbrDebugPtr
539  (
540  coupledPolyPatch::debug
541  ? new ownToNbrDebugOrderData()
542  : nullptr
543  );
544 
545  if (!owner())
546  {
547  UIPstream fromOwner(neighbProcNo(), pBufs);
548  fromOwner >> ownToNbr;
549  ownToNbr.transform(transform());
550  if (coupledPolyPatch::debug)
551  {
552  fromOwner >> ownToNbrDebugPtr();
553  ownToNbrDebugPtr->transform(transform());
554  }
555  }
556 
557  return
559  (
560  ownToNbr,
561  ownToNbrDebugPtr,
562  pp,
563  faceMap,
564  rotation
565  );
566 }
567 
568 
570 {
572  writeEntry(os, "myProcNo", myProcNo_);
573  writeEntry(os, "neighbProcNo", neighbProcNo_);
574 }
575 
576 
577 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
Buffers for inter-processor communications streams (UOPstream, UIPstream).
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:57
A List with indirect addressing.
Definition: UIndirectList.H:60
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:58
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual bool order(const ownToNbrOrderData &ownToNbr, const autoPtr< ownToNbrDebugOrderData > &ownToNbrDebugPtr, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for the given primitivePatch.
virtual void initOrder(ownToNbrOrderData &ownToNbr, autoPtr< ownToNbrDebugOrderData > &ownToNbrDebugPtr, const primitivePatch &) const
Initialise ordering for the given primitivePatch. Fills the.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:76
A class for handling file names.
Definition: fileName.H:82
word name() const
Return file name (part beyond last /)
Definition: fileName.C:195
Foam::polyBoundaryMesh.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
virtual void movePoints(const pointField &p)
Correct patches after moving points.
Definition: polyPatch.C:57
virtual void initTopoChange(PstreamBuffers &)
Initialise the update of the patch topology.
Definition: polyPatch.H:115
virtual void topoChange(PstreamBuffers &)
Update of the patch topology.
Definition: polyPatch.C:69
Neighbour processor patch.
processorPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const int myProcNo, const int neighbProcNo, const word &patchType=typeName)
Construct from components with specified name.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
virtual ~processorPolyPatch()
Destructor.
void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialise ordering for primitivePatch. Does not.
void initCalcGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
const labelList & nbrPoints() const
Return neighbour point labels. WIP.
const labelList & nbrEdges() const
Return neighbour edge labels. WIP.
static word newName(const label myProcNo, const label neighbProcNo)
Return the name of a processorPolyPatch.
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
void movePoints(PstreamBuffers &, const pointField &)
Correct patches after moving points.
virtual void initTopoChange(PstreamBuffers &)
Initialise the update of the patch topology.
virtual void topoChange(PstreamBuffers &)
Update of the patch topology.
A class for handling words, derived from string.
Definition: word.H:62
Template functions to aid in the implementation of demand driven data.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
const pointField & points
label nPoints
Determine correspondence between points. See below.
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
errorManip< error > abort(error &err)
Definition: errorManip.H:131
vector point
Point is a vector.
Definition: point.H:41
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
dimensioned< scalar > mag(const dimensioned< Type > &)
defineTypeNameAndDebug(combustionModel, 0)
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
error FatalError
labelList f(nPoints)
dictionary dict
volScalarField & p
Data to pass from owner.initOrder to nbr.order if debugging.
Data to pass from owner.initOrder to nbr.order.
void transform(const transformer &tr)