cyclicACMIPolyPatch.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) 2013-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 "cyclicACMIPolyPatch.H"
27 #include "SubField.H"
28 #include "Time.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(cyclicACMIPolyPatch, 0);
36 
37  addToRunTimeSelectionTable(polyPatch, cyclicACMIPolyPatch, word);
38  addToRunTimeSelectionTable(polyPatch, cyclicACMIPolyPatch, dictionary);
39 }
40 
41 const Foam::scalar Foam::cyclicACMIPolyPatch::tolerance_ = 1e-10;
42 
43 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
44 
46 {
47  if (owner())
48  {
49  const polyPatch& nonOverlapPatch = this->nonOverlapPatch();
50 
51  if (debug)
52  {
53  Pout<< "cyclicACMIPolyPatch::resetAMI : recalculating weights"
54  << " for " << name() << " and " << nonOverlapPatch.name()
55  << endl;
56  }
57 
58  if (boundaryMesh().mesh().hasCellCentres())
59  {
60  if (debug)
61  {
62  Pout<< "cyclicACMIPolyPatch::resetAMI : clearing cellCentres"
63  << " for " << name() << " and " << nonOverlapPatch.name()
64  << endl;
65  }
66 
67  const_cast<polyMesh&>
68  (
69  boundaryMesh().mesh()
71  }
72 
73  // Trigger re-building of faceAreas
74  (void)boundaryMesh().mesh().faceAreas();
75 
76  // Calculate the AMI using partial face-area-weighted. This leaves
77  // the weights as fractions of local areas (sum(weights) = 1 means
78  // face is fully covered)
80 
81  AMIInterpolation& AMI = this->AMIs_[0];
82 
83  srcMask_ =
84  min(scalar(1) - tolerance_, max(tolerance_, AMI.srcWeightsSum()));
85 
86  tgtMask_ =
87  min(scalar(1) - tolerance_, max(tolerance_, AMI.tgtWeightsSum()));
88 
89 
90  // Adapt owner side areas. Note that in uncoupled situations (e.g.
91  // decomposePar) srcMask, tgtMask can be zero size.
92  if (srcMask_.size())
93  {
94  vectorField::subField Sf = faceAreas();
95  vectorField::subField noSf = nonOverlapPatch.faceAreas();
96 
97  forAll(Sf, facei)
98  {
99  Sf[facei] *= srcMask_[facei];
100  noSf[facei] *= 1.0 - srcMask_[facei];
101  }
102  }
103  // Adapt slave side areas
104  if (tgtMask_.size())
105  {
106  const cyclicACMIPolyPatch& cp =
107  refCast<const cyclicACMIPolyPatch>(this->neighbPatch());
108  const polyPatch& pp = cp.nonOverlapPatch();
109 
111  vectorField::subField noSf = pp.faceAreas();
112 
113  forAll(Sf, facei)
114  {
115  Sf[facei] *= tgtMask_[facei];
116  noSf[facei] *= 1.0 - tgtMask_[facei];
117  }
118  }
119 
120  // Re-normalise the weights since the effect of overlap is already
121  // accounted for in the area.
122  {
123  scalarListList& srcWeights = AMI.srcWeights();
124  scalarField& srcWeightsSum = AMI.srcWeightsSum();
125  forAll(srcWeights, i)
126  {
127  scalarList& wghts = srcWeights[i];
128  if (wghts.size())
129  {
130  scalar& sum = srcWeightsSum[i];
131 
132  forAll(wghts, j)
133  {
134  wghts[j] /= sum;
135  }
136  sum = 1.0;
137  }
138  }
139  }
140  {
141  scalarListList& tgtWeights = AMI.tgtWeights();
142  scalarField& tgtWeightsSum = AMI.tgtWeightsSum();
143  forAll(tgtWeights, i)
144  {
145  scalarList& wghts = tgtWeights[i];
146  if (wghts.size())
147  {
148  scalar& sum = tgtWeightsSum[i];
149  forAll(wghts, j)
150  {
151  wghts[j] /= sum;
152  }
153  sum = 1.0;
154  }
155  }
156  }
157 
158  // Set the updated flag
159  updated_ = true;
160  }
161 }
162 
163 
165 {
167 
168  // Initialise the AMI
169  resetAMI();
170 }
171 
172 
174 (
175  PstreamBuffers& pBufs,
176  const pointField& p
177 )
178 {
180 
181  // Initialise the AMI
182  resetAMI();
183 }
184 
185 
187 {
188  return srcMask_;
189 }
190 
191 
193 {
194  return tgtMask_;
195 }
196 
197 
198 // * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * //
199 
201 (
202  const word& name,
203  const label size,
204  const label start,
205  const label index,
206  const polyBoundaryMesh& bm,
207  const word& patchType,
208  const transformType transform
209 )
210 :
212  (
213  name,
214  size,
215  start,
216  index,
217  bm,
218  patchType,
219  transform,
220  false,
222  ),
223  nonOverlapPatchName_(word::null),
224  nonOverlapPatchID_(-1),
225  srcMask_(),
226  tgtMask_(),
227  updated_(false)
228 {
229  // Non-overlapping patch might not be valid yet so cannot determine
230  // associated patchID
231 }
232 
233 
235 (
236  const word& name,
237  const dictionary& dict,
238  const label index,
239  const polyBoundaryMesh& bm,
240  const word& patchType
241 )
242 :
244  (
245  name,
246  dict,
247  index,
248  bm,
249  patchType,
250  false,
252  ),
253  nonOverlapPatchName_(dict.lookup("nonOverlapPatch")),
254  nonOverlapPatchID_(-1),
255  srcMask_(),
256  tgtMask_(),
257  updated_(false)
258 {
259  if (nonOverlapPatchName_ == name)
260  {
262  (
263  dict
264  ) << "Non-overlapping patch name " << nonOverlapPatchName_
265  << " cannot be the same as this patch " << name
266  << exit(FatalIOError);
267  }
268 
269  // Non-overlapping patch might not be valid yet so cannot determine
270  // associated patchID
271 }
272 
273 
275 (
276  const cyclicACMIPolyPatch& pp,
277  const polyBoundaryMesh& bm
278 )
279 :
280  cyclicAMIPolyPatch(pp, bm),
281  nonOverlapPatchName_(pp.nonOverlapPatchName_),
282  nonOverlapPatchID_(-1),
283  srcMask_(),
284  tgtMask_(),
285  updated_(false)
286 {
287  // Non-overlapping patch might not be valid yet so cannot determine
288  // associated patchID
289 }
290 
291 
293 (
294  const cyclicACMIPolyPatch& pp,
295  const polyBoundaryMesh& bm,
296  const label index,
297  const label newSize,
298  const label newStart,
299  const word& nbrPatchName,
300  const word& nonOverlapPatchName
301 )
302 :
303  cyclicAMIPolyPatch(pp, bm, index, newSize, newStart, nbrPatchName),
304  nonOverlapPatchName_(nonOverlapPatchName),
305  nonOverlapPatchID_(-1),
306  srcMask_(),
307  tgtMask_(),
308  updated_(false)
309 {
310  if (nonOverlapPatchName_ == name())
311  {
313  << "Non-overlapping patch name " << nonOverlapPatchName_
314  << " cannot be the same as this patch " << name()
315  << exit(FatalError);
316  }
317 
318  // Non-overlapping patch might not be valid yet so cannot determine
319  // associated patchID
320 }
321 
322 
324 (
325  const cyclicACMIPolyPatch& pp,
326  const polyBoundaryMesh& bm,
327  const label index,
328  const labelUList& mapAddressing,
329  const label newStart
330 )
331 :
332  cyclicAMIPolyPatch(pp, bm, index, mapAddressing, newStart),
333  nonOverlapPatchName_(pp.nonOverlapPatchName_),
334  nonOverlapPatchID_(-1),
335  srcMask_(),
336  tgtMask_(),
337  updated_(false)
338 {
339  // Non-overlapping patch might not be valid yet so cannot determine
340  // associated patchID
341 }
342 
343 
344 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
345 
347 {}
348 
349 
350 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
351 
353 {
354  const polyPatch& pp = this->boundaryMesh()[neighbPatchID()];
355  return refCast<const cyclicACMIPolyPatch>(pp);
356 }
357 
358 
360 {
361  if (nonOverlapPatchID_ == -1)
362  {
363  nonOverlapPatchID_ =
364  this->boundaryMesh().findPatchID(nonOverlapPatchName_);
365 
366  if (nonOverlapPatchID_ == -1)
367  {
369  << "Illegal non-overlapping patch name " << nonOverlapPatchName_
370  << nl << "Valid patch names are "
371  << this->boundaryMesh().names()
372  << exit(FatalError);
373  }
374 
375  if (nonOverlapPatchID_ < index())
376  {
378  << "Boundary ordering error: " << type()
379  << " patch must be defined prior to its non-overlapping patch"
380  << nl
381  << type() << " patch: " << name() << ", ID:" << index() << nl
382  << "Non-overlap patch: " << nonOverlapPatchName_
383  << ", ID:" << nonOverlapPatchID_ << nl
384  << exit(FatalError);
385  }
386 
387  const polyPatch& noPp = this->boundaryMesh()[nonOverlapPatchID_];
388 
389  bool ok = true;
390 
391  if (size() == noPp.size())
392  {
393  const scalarField magSf(mag(faceAreas()));
394  const scalarField noMagSf(mag(noPp.faceAreas()));
395 
396  forAll(magSf, facei)
397  {
398  scalar ratio = mag(magSf[facei]/(noMagSf[facei] + rootVSmall));
399 
400  if (ratio - 1 > tolerance_)
401  {
402  ok = false;
403  break;
404  }
405  }
406  }
407  else
408  {
409  ok = false;
410  }
411 
412  if (!ok)
413  {
415  << "Inconsistent ACMI patches " << name() << " and "
416  << noPp.name() << ". Patches should have identical topology"
417  << exit(FatalError);
418  }
419  }
420 
421  return nonOverlapPatchID_;
422 }
423 
424 
426 (
427  const primitivePatch& referPatch,
428  const pointField& thisCtrs,
429  const vectorField& thisAreas,
430  const pointField& thisCc,
431  const pointField& nbrCtrs,
432  const vectorField& nbrAreas,
433  const pointField& nbrCc
434 )
435 {
437  (
438  referPatch,
439  thisCtrs,
440  thisAreas,
441  thisCc,
442  nbrCtrs,
443  nbrAreas,
444  nbrCc
445  );
446 }
447 
448 
450 (
451  PstreamBuffers& pBufs,
452  const primitivePatch& pp
453 ) const
454 {
456 }
457 
458 
460 (
461  PstreamBuffers& pBufs,
462  const primitivePatch& pp,
463  labelList& faceMap,
464  labelList& rotation
465 ) const
466 {
467  return cyclicAMIPolyPatch::order(pBufs, pp, faceMap, rotation);
468 }
469 
470 
472 {
474  writeEntry(os, "nonOverlapPatch", nonOverlapPatchName_);
475 }
476 
477 
478 // ************************************************************************* //
virtual void resetAMI() const
Reset the AMI interpolator.
Cyclic patch for Arbitrarily Coupled Mesh Interface (ACMI)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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.
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
virtual const cyclicACMIPolyPatch & neighbPatch() const
Return a reference to the neighbour patch.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
const scalarListList & tgtWeights() const
Return const access to target patch weights.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:753
virtual const scalarField & tgtMask() const
Return the mask/weighting for the target patch.
Pre-declare related SubField type.
Definition: Field.H:60
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
Macros for easy insertion into run-time selection tables.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
const polyPatch & nonOverlapPatch() const
Return a const reference to the non-overlapping patch.
A list of faces which address into the list of points.
void clearGeom()
Clear geometry.
dynamicFvMesh & mesh
virtual const scalarField & srcMask() const
Return the mask/weighting for the source patch.
A class for handling words, derived from string.
Definition: word.H:59
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
static const word null
An empty word.
Definition: word.H:77
const scalarField & tgtWeightsSum() const
Return const access to normalisation factor of target.
cyclicACMIPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType, const transformType transform=UNKNOWN)
Construct from (base couped patch) components.
Cyclic patch for Arbitrary Mesh Interface (AMI)
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:53
const vectorField::subField faceAreas() const
Return face normals.
Definition: polyPatch.C:290
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const scalarListList & srcWeights() const
Return const access to source patch weights.
static const char nl
Definition: Ostream.H:265
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
defineTypeNameAndDebug(combustionModel, 0)
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const scalarField & srcWeightsSum() const
Return const access to normalisation factor of source.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
virtual void resetAMI() const
Reset the AMI interpolator.
virtual void calcGeometry(const primitivePatch &referPatch, const pointField &thisCtrs, const vectorField &thisAreas, const pointField &thisCc, const pointField &nbrCtrs, const vectorField &nbrAreas, const pointField &nbrCc)
Calculate the patch geometry.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
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...
virtual ~cyclicACMIPolyPatch()
Destructor.
dimensioned< scalar > mag(const dimensioned< Type > &)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const bMesh & mesh() const
Definition: boundaryMesh.H:199
mesh Sf()
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
virtual label nonOverlapPatchID() const
Non-overlapping patch ID.
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
IOerror FatalIOError