processorFaPatch.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2016-2017 Wikki Ltd
9  Copyright (C) 2019 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "processorFaPatch.H"
31 #include "IPstream.H"
32 #include "OPstream.H"
33 #include "transformField.H"
34 #include "faBoundaryMesh.H"
35 #include "faMesh.H"
36 #include "globalMeshData.H"
37 #include "demandDrivenData.H"
38 
39 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43  defineTypeNameAndDebug(processorFaPatch, 0);
44  addToRunTimeSelectionTable(faPatch, processorFaPatch, dictionary);
45 }
46 
47 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
48 
50 {
51  deleteDemandDrivenData(neighbPointsPtr_);
52  deleteDemandDrivenData(nonGlobalPatchPointsPtr_);
53 }
54 
55 
56 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
57 
58 Foam::label Foam::processorFaPatch::comm() const
59 {
60  return boundaryMesh().mesh().comm();
61 }
62 
63 
65 {
66  return Pstream::msgType();
67 }
68 
69 
71 {
72  // If it is not running parallel or there are no global points
73  // create a 1->1 map
74 
75  // Can not use faGlobalMeshData at this point yet
76 
77  if
78  (
80  || !boundaryMesh().mesh()().globalData().nGlobalPoints()
81  )
82  {
83  nonGlobalPatchPointsPtr_ = new labelList(identity(nPoints()));
84  }
85  else
86  {
87  // Get reference to shared points
88  const labelList& sharedPoints =
89  boundaryMesh().mesh()().globalData().sharedPointLabels();
90 
91  nonGlobalPatchPointsPtr_ = new labelList(nPoints());
92  labelList& ngpp = *nonGlobalPatchPointsPtr_;
93 
94  const labelList& faMeshPatchPoints = pointLabels();
95 
96  const labelList& meshPoints =
97  boundaryMesh().mesh().patch().meshPoints();
98 
99  label noFiltPoints = 0;
100 
101  forAll(faMeshPatchPoints, pointI)
102  {
103  label curP = meshPoints[faMeshPatchPoints[pointI]];
104 
105  bool found = false;
106 
107  forAll(sharedPoints, sharedI)
108  {
109  if (sharedPoints[sharedI] == curP)
110  {
111  found = true;
112  break;
113  }
114  }
115 
116  if (!found)
117  {
118  ngpp[noFiltPoints] = pointI;
119  ++noFiltPoints;
120  }
121  }
122 
123  ngpp.setSize(noFiltPoints);
124 
125 
126 // // Get reference to shared points
127 // const labelList& sharedPoints =
128 // boundaryMesh().mesh().globalData().sharedPointLabels();
129 
130 // nonGlobalPatchPointsPtr_ = new labelList(nPoints());
131 // labelList& ngpp = *nonGlobalPatchPointsPtr_;
132 
133 // const labelList& patchPoints = pointLabels();
134 
135 // label noFiltPoints = 0;
136 
137 // forAll(patchPoints, pointI)
138 // {
139 // label curP = patchPoints[pointI];
140 
141 // bool found = false;
142 
143 // forAll(sharedPoints, pI)
144 // {
145 // if (sharedPoints[pI] == curP)
146 // {
147 // found = true;
148 // break;
149 // }
150 // }
151 
152 // if (!found)
153 // {
154 // ngpp[noFiltPoints] = pointI;
155 // noFiltPoints++;
156 // }
157 // }
158 
159 // ngpp.setSize(noFiltPoints);
160  }
161 }
162 
163 
165 {
166  if (Pstream::parRun())
167  {
168  OPstream toNeighbProc
169  (
171  neighbProcNo(),
172  3*(sizeof(label) + size()*sizeof(vector))
173  );
174 
175  toNeighbProc
176  << edgeCentres()
177  << edgeLengths()
178  << edgeFaceCentres();
179  }
180 }
181 
182 
184 {
185  if (Pstream::parRun())
186  {
187  {
188  IPstream fromNeighbProc
189  (
191  neighbProcNo(),
192  3*(sizeof(label) + size()*sizeof(vector))
193  );
194  fromNeighbProc
195  >> neighbEdgeCentres_
196  >> neighbEdgeLengths_
197  >> neighbEdgeFaceCentres_;
198  }
199 
200  const scalarField& magEl = magEdgeLengths();
201 
202  forAll(magEl, edgei)
203  {
204  scalar nmagEl = mag(neighbEdgeLengths_[edgei]);
205  scalar avEl = (magEl[edgei] + nmagEl)/2.0;
206 
207  if (mag(magEl[edgei] - nmagEl)/avEl > 1e-6)
208  {
210  << "edge " << edgei
211  << " length does not match neighbour by "
212  << 100*mag(magEl[edgei] - nmagEl)/avEl
213  << "% -- possible edge ordering problem"
214  << exit(FatalError);
215  }
216  }
217 
218  calcTransformTensors
219  (
220  edgeCentres(),
221  neighbEdgeCentres_,
222  edgeNormals(),
223  neighbEdgeLengths_/mag(neighbEdgeLengths_)
224  );
225  }
226 }
227 
228 
230 {
232  initGeometry();
233 }
234 
235 
237 {
238  calcGeometry();
239 }
240 
241 
243 {
244  // For completeness
246 
247  deleteDemandDrivenData(neighbPointsPtr_);
248 
249  if (Pstream::parRun())
250  {
251  // Express all points as patch edge and index in edge.
252  labelList patchEdge(nPoints());
253  labelList indexInEdge(nPoints());
254 
255  const edgeList::subList patchEdges =
256  patchSlice(boundaryMesh().mesh().edges());
257 
258  const labelListList& ptEdges = pointEdges();
259 
260  for (label patchPointI = 0; patchPointI < nPoints(); ++patchPointI)
261  {
262  label edgeI = ptEdges[patchPointI][0];
263 
264  patchEdge[patchPointI] = edgeI;
265 
266  const edge& e = patchEdges[edgeI];
267 
268  indexInEdge[patchPointI] = e.find(pointLabels()[patchPointI]);
269  }
270 
271  OPstream toNeighbProc
272  (
274  neighbProcNo(),
275  2*sizeof(label) + 2*nPoints()*sizeof(label)
276  );
277 
278  toNeighbProc
279  << patchEdge
280  << indexInEdge;
281  }
282 }
283 
284 
286 {
287  // For completeness
289 
290  if (Pstream::parRun())
291  {
292  labelList nbrPatchEdge(nPoints());
293  labelList nbrIndexInEdge(nPoints());
294 
295  {
296  // Note cannot predict exact size since edgeList not (yet) sent as
297  // binary entity but as List of edges.
298  IPstream fromNeighbProc
299  (
301  neighbProcNo()
302  );
303 
304  fromNeighbProc
305  >> nbrPatchEdge
306  >> nbrIndexInEdge;
307  }
308 
309  if (nbrPatchEdge.size() == nPoints())
310  {
311  // Convert neighbour edges and indices into face back into
312  // my edges and points.
313  neighbPointsPtr_ = new labelList(nPoints());
314  labelList& neighbPoints = *neighbPointsPtr_;
315 
316  const edgeList::subList patchEdges =
317  patchSlice(boundaryMesh().mesh().edges());
318 
319  forAll(nbrPatchEdge, nbrPointI)
320  {
321  // Find edge and index in edge on this side.
322  const edge& e = patchEdges[nbrPatchEdge[nbrPointI]];
323 
324  const label index = 1 - nbrIndexInEdge[nbrPointI];
325 
326  const label patchPointI = pointLabels().find(e[index]);
327 
328  neighbPoints[patchPointI] = nbrPointI;
329  }
330  }
331  else
332  {
333  // Differing number of points. Probably patch includes
334  // part of a cyclic.
335  neighbPointsPtr_ = nullptr;
336  }
337  }
338 }
339 
340 
342 {
343  if (!neighbPointsPtr_)
344  {
345  // Was probably created from cyclic patch and hence the
346  // number of edges or points might differ on both
347  // sides of the processor patch since one side might have
348  // it merged with another bit of geometry
349 
351  << "No extended addressing calculated for patch " << name()
352  << nl
353  << "This can happen if the number of points on both"
354  << " sides of the two coupled patches differ." << nl
355  << "This happens if the processorPatch was constructed from"
356  << " part of a cyclic patch."
357  << abort(FatalError);
358  }
359 
360  return *neighbPointsPtr_;
361 }
362 
363 
365 {
366  if (Pstream::parRun())
367  {
368  // The face normals point in the opposite direction on the other side
369  scalarField neighbEdgeCentresCn
370  (
371  (
372  neighbEdgeLengths()
373  /mag(neighbEdgeLengths())
374  )
375  & (
376  neighbEdgeCentres()
377  - neighbEdgeFaceCentres()
378  )
379  );
380 
381  w = neighbEdgeCentresCn/
382  (
383  (edgeNormals() & faPatch::delta())
384  + neighbEdgeCentresCn
385  );
386  }
387  else
388  {
389  w = 1.0;
390  }
391 }
392 
393 
395 {
396  if (Pstream::parRun())
397  {
398  dc = (1.0 - weights())/(edgeNormals() & faPatch::delta());
399  }
400  else
401  {
402  dc = 1.0/(edgeNormals() & faPatch::delta());
403  }
404 }
405 
406 
408 {
409  if (Pstream::parRun())
410  {
411  // To the transformation if necessary
412  if (parallel())
413  {
414  return
416  - (
417  neighbEdgeCentres()
418  - neighbEdgeFaceCentres()
419  );
420  }
421  else
422  {
423  return
425  - transform
426  (
427  forwardT(),
428  (
429  neighbEdgeCentres()
430  - neighbEdgeFaceCentres()
431  )
432  );
433  }
434  }
435  else
436  {
437  return faPatch::delta();
438  }
439 }
440 
441 
443 {
444  if (!nonGlobalPatchPointsPtr_)
445  {
446  makeNonGlobalPatchPoints();
447  }
448 
449  return *nonGlobalPatchPointsPtr_;
450 }
451 
452 
454 (
455  const labelUList& internalData
456 ) const
457 {
458  return patchInternalField(internalData);
459 }
460 
461 
463 (
464  const labelUList& internalData,
465  const labelUList& edgeFaces
466 ) const
467 {
468  return patchInternalField(internalData, edgeFaces);
469 }
470 
471 
473 (
474  const Pstream::commsTypes commsType,
475  const labelUList& interfaceData
476 ) const
477 {
478  send(commsType, interfaceData);
479 }
480 
481 
483 (
484  const Pstream::commsTypes commsType,
485  const labelUList&
486 ) const
487 {
488  return receive<label>(commsType, this->size());
489 }
490 
491 
493 (
494  const Pstream::commsTypes commsType,
495  const labelUList& iF
496 ) const
497 {
498  send(commsType, patchInternalField(iF)());
499 }
500 
501 
503 (
504  const Pstream::commsTypes commsType,
505  const labelUList&
506 ) const
507 {
508  return receive<label>(commsType, this->size());
509 }
510 
511 
513 (
514  const Pstream::commsTypes commsType,
515  const labelUList&,
516  const labelUList&
517 ) const
518 {
519  return receive<label>(commsType, this->size());
520 }
521 
522 
524 {
526  os.writeEntry("myProcNo", myProcNo_);
527  os.writeEntry("neighbProcNo", neighbProcNo_);
528 }
529 
530 
531 // ************************************************************************* //
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
Foam::processorFaPatch::comm
virtual label comm() const
Return communicator used for communication.
Definition: processorFaPatch.C:58
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:67
Foam::boundaryMesh::mesh
const bMesh & mesh() const
Definition: boundaryMesh.H:206
Foam::processorFaPatch::movePoints
void movePoints(const pointField &)
Correct patches after moving points.
Definition: processorFaPatch.C:236
Foam::processorFaPatch::updateMesh
virtual void updateMesh()
Update of the patch topology.
Definition: processorFaPatch.C:285
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::processorFaPatch::initGeometry
void initGeometry()
Initialise the calculation of the patch geometry.
Definition: processorFaPatch.C:164
Foam::processorFaPatch::write
virtual void write(Ostream &os) const
Write the patch data as a dictionary.
Definition: processorFaPatch.C:523
Foam::processorFaPatch::initTransfer
virtual void initTransfer(const Pstream::commsTypes commsType, const labelUList &interfaceData) const
Initialise interface data transfer.
Definition: processorFaPatch.C:473
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::faPatch::movePoints
virtual void movePoints(const pointField &)
Correct patch after moving points.
Definition: faPatch.C:463
Foam::OPstream
Output inter-processor communications stream.
Definition: OPstream.H:53
globalMeshData.H
demandDrivenData.H
Template functions to aid in the implementation of demand driven data.
Foam::SubList
A List obtained as a section of another List.
Definition: SubList.H:54
processorFaPatch.H
Foam::edge
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:63
OPstream.H
faMesh.H
Foam::processorFaPatch::transfer
virtual tmp< labelField > transfer(const Pstream::commsTypes commsType, const labelUList &interfaceData) const
Transfer and return neighbour field.
Definition: processorFaPatch.C:483
Foam::transform
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Definition: dimensionSet.C:521
Foam::processorFaPatch::initInternalFieldTransfer
virtual void initInternalFieldTransfer(const Pstream::commsTypes commsType, const labelUList &internalData) const
Initialise neighbour field transfer.
Definition: processorFaPatch.C:493
Foam::faPatch::delta
virtual tmp< vectorField > delta() const
Return cell-centre to face-centre vector.
Definition: faPatch.C:423
Foam::faPatch::initUpdateMesh
virtual void initUpdateMesh()
Initialise the update of the patch topology.
Definition: faPatch.H:133
transformField.H
Spatial transformation functions for primitive fields.
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
nPoints
label nPoints
Definition: gmvOutputHeader.H:2
Foam::deleteDemandDrivenData
void deleteDemandDrivenData(DataPtr &dataPtr)
Definition: demandDrivenData.H:42
IPstream.H
Foam::Field< scalar >
Foam::processorFaPatch::neighbPoints
const labelList & neighbPoints() const
Return neighbour point labels. This is for my local point the.
Definition: processorFaPatch.C:341
Foam::processorFaPatch::initUpdateMesh
virtual void initUpdateMesh()
Initialise the update of the patch topology.
Definition: processorFaPatch.C:242
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:222
Foam::FatalError
error FatalError
faBoundaryMesh.H
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::processorFaPatch::makeNonGlobalPatchPoints
void makeNonGlobalPatchPoints() const
Find non-globa patch points.
Definition: processorFaPatch.C:70
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::faPatch::updateMesh
virtual void updateMesh()
Update of the patch topology.
Definition: faPatch.H:137
Foam::processorFaPatch::makeDeltaCoeffs
void makeDeltaCoeffs(scalarField &) const
Make patch face - neighbour cell distances.
Definition: processorFaPatch.C:394
found
bool found
Definition: TABSMDCalcMethod2.H:32
Foam::UPstream::commsTypes
commsTypes
Types of communications.
Definition: UPstream.H:69
Foam::processorFaPatch::nonGlobalPatchPoints
const labelList & nonGlobalPatchPoints() const
Return the set of labels of the processor patch points which are.
Definition: processorFaPatch.C:442
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::processorFaPatch::~processorFaPatch
virtual ~processorFaPatch()
Destructor.
Definition: processorFaPatch.C:49
Foam::UPstream::msgType
static int & msgType() noexcept
Message tag of standard messages.
Definition: UPstream.H:540
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::faPatch::write
virtual void write(Ostream &) const
Write.
Definition: faPatch.C:481
Foam::Vector< scalar >
Foam::processorFaPatch::interfaceInternalField
virtual tmp< labelField > interfaceInternalField(const labelUList &internalData) const
Definition: processorFaPatch.C:454
Foam::UPstream::parRun
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:433
Foam::List< label >
Foam::processorFaPatch::makeWeights
void makeWeights(scalarField &) const
Make patch weighting factors.
Definition: processorFaPatch.C:364
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::processorFaPatch::internalFieldTransfer
virtual tmp< labelField > internalFieldTransfer(const Pstream::commsTypes commsType, const labelUList &internalData) const
Return neighbour field.
Definition: processorFaPatch.C:503
Foam::UList< label >
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
Foam::identity
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
Definition: labelList.C:38
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::boundaryMesh
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
Definition: boundaryMesh.H:62
Foam::IPstream
Input inter-processor communications stream.
Definition: IPstream.H:53
Foam::processorFaPatch::delta
virtual tmp< vectorField > delta() const
Return delta (P to N) vectors across coupled patch.
Definition: processorFaPatch.C:407
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::processorFaPatch::calcGeometry
void calcGeometry()
Calculate the patch geometry.
Definition: processorFaPatch.C:183
Foam::PrimitivePatch::meshPoints
const labelList & meshPoints() const
Return labelList of mesh points in patch.
Definition: PrimitivePatch.C:330
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::processorFaPatch::initMovePoints
void initMovePoints(const pointField &)
Initialise the patches for moving points.
Definition: processorFaPatch.C:229
pointLabels
labelList pointLabels(nPoints, -1)
Foam::processorFaPatch::tag
virtual int tag() const
Return message tag to use for communication.
Definition: processorFaPatch.C:64
Foam::UPstream::commsTypes::blocking