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 {
42  defineTypeNameAndDebug(processorPolyPatch, 0);
43  addToRunTimeSelectionTable(polyPatch, processorPolyPatch, dictionary);
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 
150 (
151  const processorPolyPatch& pp,
152  const polyBoundaryMesh& bm,
153  const label index,
154  const labelUList& mapAddressing,
155  const label newStart
156 )
157 :
158  coupledPolyPatch(pp, bm, index, mapAddressing, newStart),
159  myProcNo_(pp.myProcNo_),
160  neighbProcNo_(pp.neighbProcNo_),
161  neighbFaceCentres_(),
162  neighbFaceAreas_(),
163  neighbFaceCellCentres_()
164 {}
165 
166 
167 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
168 
170 {
171  nbrPointsPtr_.clear();
172  nbrEdgesPtr_.clear();
173 }
174 
175 
176 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
177 
179 (
180  const label myProcNo,
181  const label neighbProcNo
182 )
183 {
184  return
185  "procBoundary"
186  + Foam::name(myProcNo)
187  + "to"
188  + Foam::name(neighbProcNo);
189 }
190 
191 
193 {
194  if (Pstream::parRun())
195  {
196  UOPstream toNeighbProc(neighbProcNo(), pBufs);
197 
198  toNeighbProc
199  << faceCentres()
200  << faceAreas()
201  << faceCellCentres();
202  }
203 }
204 
205 
207 {
208  if (Pstream::parRun())
209  {
210  {
211  UIPstream fromNeighbProc(neighbProcNo(), pBufs);
212 
213  fromNeighbProc
214  >> neighbFaceCentres_
215  >> neighbFaceAreas_
216  >> neighbFaceCellCentres_;
217  }
218 
219  // My normals
220  vectorField faceNormals(size());
221 
222  // Neighbour normals
223  vectorField nbrFaceNormals(neighbFaceAreas_.size());
224 
225  // Face match tolerances
226  scalarField tols = calcFaceTol(*this, points(), faceCentres());
227 
228  // Calculate normals from areas and check
229  forAll(faceNormals, facei)
230  {
231  scalar magSf = magFaceAreas()[facei];
232  scalar nbrMagSf = mag(neighbFaceAreas_[facei]);
233  scalar avSf = (magSf + nbrMagSf)/2.0;
234 
235  // For small face area calculation the results of the area
236  // calculation have been found to only be accurate to ~1e-20
237  if (magSf < small || nbrMagSf < small)
238  {
239  // Undetermined normal. Use dummy normal to force separation
240  // check.
241  faceNormals[facei] = point(1, 0, 0);
242  nbrFaceNormals[facei] = -faceNormals[facei];
243  tols[facei] = great;
244  }
245  else if (mag(magSf - nbrMagSf) > matchTolerance()*sqr(tols[facei]))
246  {
247  const fileName patchOBJName
248  (
249  boundaryMesh().mesh().time().path()/name() + "_faces.obj"
250  );
251 
252  Pout<< "processorPolyPatch::calcGeometry : Writing my "
253  << size() << " faces to " << patchOBJName << endl;
254 
255  writeOBJ(patchOBJName, *this);
256 
257  const fileName centresOBJName
258  (
259  boundaryMesh().mesh().time().path()/name()
260  + "_faceCentresConnections.obj"
261  );
262 
263  Pout<< "processorPolyPatch::calcGeometry :"
264  << " Dumping lines between corresponding face centres to "
265  << centresOBJName.name() << endl;
266 
267  writeOBJ(centresOBJName, neighbFaceCentres_, faceCentres());
268 
270  << "face " << facei << " area does not match neighbour by "
271  << 100*mag(magSf - nbrMagSf)/avSf
272  << "% -- possible face ordering problem." << endl
273  << "patch:" << name()
274  << " my area:" << magSf
275  << " neighbour area:" << nbrMagSf
276  << " matching tolerance:"
277  << matchTolerance()*sqr(tols[facei])
278  << endl
279  << "Mesh face:" << start()+facei
280  << " vertices:"
281  << UIndirectList<point>(points(), operator[](facei))()
282  << endl
283  << "If you are certain your matching is correct"
284  << " you can increase the 'matchTolerance' setting"
285  << " in the patch dictionary in the boundary file."
286  << endl
287  << "Rerun with processor debug flag set for"
288  << " more information." << exit(FatalError);
289  }
290  else
291  {
292  faceNormals[facei] = faceAreas()[facei]/magSf;
293  nbrFaceNormals[facei] = neighbFaceAreas_[facei]/nbrMagSf;
294  }
295  }
296  }
297 }
298 
299 
301 (
302  PstreamBuffers& pBufs,
303  const pointField& p
304 )
305 {
306  polyPatch::movePoints(pBufs, p);
308 }
309 
310 
312 (
313  PstreamBuffers& pBufs,
314  const pointField&
315 )
316 {
318 }
319 
320 
322 {
324 
325  if (Pstream::parRun())
326  {
327  // Express all points as patch face and index in face.
328  labelList pointFace(nPoints());
329  labelList pointIndex(nPoints());
330 
331  for (label patchPointi = 0; patchPointi < nPoints(); patchPointi++)
332  {
333  label facei = pointFaces()[patchPointi][0];
334 
335  pointFace[patchPointi] = facei;
336 
337  const face& f = localFaces()[facei];
338 
339  pointIndex[patchPointi] = findIndex(f, patchPointi);
340  }
341 
342  // Express all edges as patch face and index in face.
343  labelList edgeFace(nEdges());
344  labelList edgeIndex(nEdges());
345 
346  for (label patchEdgeI = 0; patchEdgeI < nEdges(); patchEdgeI++)
347  {
348  label facei = edgeFaces()[patchEdgeI][0];
349 
350  edgeFace[patchEdgeI] = facei;
351 
352  const labelList& fEdges = faceEdges()[facei];
353 
354  edgeIndex[patchEdgeI] = findIndex(fEdges, patchEdgeI);
355  }
356 
357  UOPstream toNeighbProc(neighbProcNo(), pBufs);
358 
359  toNeighbProc
360  << pointFace
361  << pointIndex
362  << edgeFace
363  << edgeIndex;
364  }
365 }
366 
367 
369 {
370  // For completeness
371  polyPatch::topoChange(pBufs);
372 
373  nbrPointsPtr_.clear();
374  nbrEdgesPtr_.clear();
375 
376  if (Pstream::parRun())
377  {
378  labelList nbrPointFace;
379  labelList nbrPointIndex;
380  labelList nbrEdgeFace;
381  labelList nbrEdgeIndex;
382 
383  {
384  // !Note: there is one situation where the opposite points and
385  // do not exactly match and that is while we are distributing
386  // meshes and multiple parts come together from different
387  // processors. This can temporarily create the situation that
388  // we have points which will be merged out later but we still
389  // need the face connectivity to be correct.
390  // So: cannot check here on same points and edges.
391 
392  UIPstream fromNeighbProc(neighbProcNo(), pBufs);
393 
394  fromNeighbProc
395  >> nbrPointFace
396  >> nbrPointIndex
397  >> nbrEdgeFace
398  >> nbrEdgeIndex;
399  }
400 
401  // Convert neighbour faces and indices into face back into
402  // my edges and points.
403 
404  // Convert points.
405  // ~~~~~~~~~~~~~~~
406 
407  nbrPointsPtr_.reset(new labelList(nPoints(), -1));
408  labelList& nbrPoints = nbrPointsPtr_();
409 
410  forAll(nbrPointFace, nbrPointi)
411  {
412  // Find face and index in face on this side.
413  const face& f = localFaces()[nbrPointFace[nbrPointi]];
414 
415  label index = (f.size() - nbrPointIndex[nbrPointi]) % f.size();
416  label patchPointi = f[index];
417 
418  if (nbrPoints[patchPointi] == -1)
419  {
420  // First reference of point
421  nbrPoints[patchPointi] = nbrPointi;
422  }
423  else if (nbrPoints[patchPointi] >= 0)
424  {
425  // Point already visited. Mark as duplicate.
426  nbrPoints[patchPointi] = -2;
427  }
428  }
429 
430  // Reset all duplicate entries to -1.
431  forAll(nbrPoints, patchPointi)
432  {
433  if (nbrPoints[patchPointi] == -2)
434  {
435  nbrPoints[patchPointi] = -1;
436  }
437  }
438 
439  // Convert edges.
440  // ~~~~~~~~~~~~~~
441 
442  nbrEdgesPtr_.reset(new labelList(nEdges(), -1));
443  labelList& nbrEdges = nbrEdgesPtr_();
444 
445  forAll(nbrEdgeFace, nbrEdgeI)
446  {
447  // Find face and index in face on this side.
448  const labelList& f = faceEdges()[nbrEdgeFace[nbrEdgeI]];
449  label index = (f.size() - nbrEdgeIndex[nbrEdgeI] - 1) % f.size();
450  label patchEdgeI = f[index];
451 
452  if (nbrEdges[patchEdgeI] == -1)
453  {
454  // First reference of edge
455  nbrEdges[patchEdgeI] = nbrEdgeI;
456  }
457  else if (nbrEdges[patchEdgeI] >= 0)
458  {
459  // Edge already visited. Mark as duplicate.
460  nbrEdges[patchEdgeI] = -2;
461  }
462  }
463 
464  // Reset all duplicate entries to -1.
465  forAll(nbrEdges, patchEdgeI)
466  {
467  if (nbrEdges[patchEdgeI] == -2)
468  {
469  nbrEdges[patchEdgeI] = -1;
470  }
471  }
472 
473  // Remove any addressing used for shared points/edges calculation
474  // since mostly not needed.
476  }
477 }
478 
479 
481 {
482  if (!nbrPointsPtr_.valid())
483  {
485  << "No extended addressing calculated for patch " << name()
486  << abort(FatalError);
487  }
488  return nbrPointsPtr_();
489 }
490 
491 
493 {
494  if (!nbrEdgesPtr_.valid())
495  {
497  << "No extended addressing calculated for patch " << name()
498  << abort(FatalError);
499  }
500  return nbrEdgesPtr_();
501 }
502 
503 
505 (
506  PstreamBuffers& pBufs,
507  const primitivePatch& pp
508 ) const
509 {
510  if (!Pstream::parRun())
511  {
512  return;
513  }
514 
515  if (owner())
516  {
517  ownToNbrOrderData ownToNbr;
518  autoPtr<ownToNbrDebugOrderData> ownToNbrDebugPtr
519  (
520  coupledPolyPatch::debug
521  ? new ownToNbrDebugOrderData()
522  : nullptr
523  );
524 
526  (
527  ownToNbr,
528  ownToNbrDebugPtr,
529  pp
530  );
531 
532  UOPstream toNeighbour(neighbProcNo(), pBufs);
533  toNeighbour << ownToNbr;
534  if (coupledPolyPatch::debug)
535  {
536  toNeighbour << ownToNbrDebugPtr();
537  }
538  }
539 }
540 
541 
543 (
544  PstreamBuffers& pBufs,
545  const primitivePatch& pp,
546  labelList& faceMap,
547  labelList& rotation
548 ) const
549 {
550  if (!Pstream::parRun())
551  {
552  return false;
553  }
554 
555  ownToNbrOrderData ownToNbr;
556  autoPtr<ownToNbrDebugOrderData> ownToNbrDebugPtr
557  (
558  coupledPolyPatch::debug
559  ? new ownToNbrDebugOrderData()
560  : nullptr
561  );
562 
563  if (!owner())
564  {
565  UIPstream fromOwner(neighbProcNo(), pBufs);
566  fromOwner >> ownToNbr;
567  ownToNbr.transform(transform());
568  if (coupledPolyPatch::debug)
569  {
570  fromOwner >> ownToNbrDebugPtr();
571  ownToNbrDebugPtr->transform(transform());
572  }
573  }
574 
575  return
577  (
578  ownToNbr,
579  ownToNbrDebugPtr,
580  pp,
581  faceMap,
582  rotation
583  );
584 }
585 
586 
588 {
590  writeEntry(os, "myProcNo", myProcNo_);
591  writeEntry(os, "neighbProcNo", neighbProcNo_);
592 }
593 
594 
595 // ************************************************************************* //
virtual bool order(const ownToNbrOrderData &ownToNbr, const autoPtr< ownToNbrDebugOrderData > &ownToNbrDebugPtr, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for the given primitivePatch.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A class for handling file names.
Definition: fileName.H:79
const labelList & nbrPoints() const
Return neighbour point labels. WIP.
const labelList & nbrEdges() const
Return neighbour edge labels. WIP.
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
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
virtual void topoChange(PstreamBuffers &)
Update of the patch topology.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void transform(const transformer &tr)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialise ordering for primitivePatch. Does not.
static word newName(const label myProcNo, const label neighbProcNo)
Return the name of a processorPolyPatch.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
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
void initCalcGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
fvMesh & mesh
virtual void initTopoChange(PstreamBuffers &)
Initialise the update of the patch topology.
Definition: polyPatch.H:115
void movePoints(PstreamBuffers &, const pointField &)
Correct patches after moving points.
Determine correspondence between points. See below.
Macros for easy insertion into run-time selection tables.
Neighbour processor patch.
void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
virtual void initOrder(ownToNbrOrderData &ownToNbr, autoPtr< ownToNbrDebugOrderData > &ownToNbrDebugPtr, const primitivePatch &) const
Initialise ordering for the given primitivePatch. Fills the.
Data to pass from owner.initOrder to nbr.order if debugging.
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:53
A list of faces which address into the list of points.
const pointField & points
A class for handling words, derived from string.
Definition: word.H:59
label nPoints
Data to pass from owner.initOrder to nbr.order.
word name() const
Return file name (part beyond last /)
Definition: fileName.C:195
virtual ~processorPolyPatch()
Destructor.
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
Foam::polyBoundaryMesh.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:54
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
labelList f(nPoints)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void topoChange(PstreamBuffers &)
Update of the patch topology.
Definition: polyPatch.C:69
Template functions to aid in the implementation of demand driven data.
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
virtual void movePoints(const pointField &p)
Correct patches after moving points.
Definition: polyPatch.C:57
vector point
Point is a vector.
Definition: point.H:41
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
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.
A List with indirect addressing.
Definition: fvMatrix.H:106
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:483
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864
virtual void initTopoChange(PstreamBuffers &)
Initialise the update of the patch topology.