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-2020 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 
77  // Calculate the AMI using partial face-area-weighted. This leaves
78  // the weights as fractions of local areas (sum(weights) = 1 means
79  // face is fully covered)
81 
82  AMIInterpolation& AMI = this->AMIs_[0];
83 
84  srcMask_ =
85  min(scalar(1) - tolerance_, max(tolerance_, AMI.srcWeightsSum()));
86 
87  tgtMask_ =
88  min(scalar(1) - tolerance_, max(tolerance_, AMI.tgtWeightsSum()));
89 
90 
91  // Adapt owner side areas. Note that in uncoupled situations (e.g.
92  // decomposePar) srcMask, tgtMask can be zero size.
93  if (srcMask_.size())
94  {
95  vectorField::subField Sf = faceAreas();
96  scalarField::subField magSf = magFaceAreas();
97  vectorField::subField noSf = nonOverlapPatch.faceAreas();
98  scalarField::subField noMagSf = nonOverlapPatch.magFaceAreas();
99 
100  forAll(Sf, facei)
101  {
102  Sf[facei] *= srcMask_[facei];
103  magSf[facei] = mag(Sf[facei]);
104  noSf[facei] *= 1.0 - srcMask_[facei];
105  noMagSf[facei] = mag(noSf[facei]);
106  }
107  }
108  // Adapt slave side areas
109  if (tgtMask_.size())
110  {
111  const cyclicACMIPolyPatch& cp =
112  refCast<const cyclicACMIPolyPatch>(this->nbrPatch());
113  const polyPatch& pp = cp.nonOverlapPatch();
114 
116  scalarField::subField magSf = cp.magFaceAreas();
117  vectorField::subField noSf = pp.faceAreas();
118  scalarField::subField noMagSf = pp.magFaceAreas();
119 
120  forAll(Sf, facei)
121  {
122  Sf[facei] *= tgtMask_[facei];
123  magSf[facei] = mag(Sf[facei]);
124  noSf[facei] *= 1.0 - tgtMask_[facei];
125  noMagSf[facei] = mag(noSf[facei]);
126  }
127  }
128 
129  // Re-normalise the weights since the effect of overlap is already
130  // accounted for in the area.
131  {
132  scalarListList& srcWeights = AMI.srcWeights();
133  scalarField& srcWeightsSum = AMI.srcWeightsSum();
134  forAll(srcWeights, i)
135  {
136  scalarList& wghts = srcWeights[i];
137  if (wghts.size())
138  {
139  scalar& sum = srcWeightsSum[i];
140 
141  forAll(wghts, j)
142  {
143  wghts[j] /= sum;
144  }
145  sum = 1.0;
146  }
147  }
148  }
149  {
150  scalarListList& tgtWeights = AMI.tgtWeights();
151  scalarField& tgtWeightsSum = AMI.tgtWeightsSum();
152  forAll(tgtWeights, i)
153  {
154  scalarList& wghts = tgtWeights[i];
155  if (wghts.size())
156  {
157  scalar& sum = tgtWeightsSum[i];
158  forAll(wghts, j)
159  {
160  wghts[j] /= sum;
161  }
162  sum = 1.0;
163  }
164  }
165  }
166 
167  // Set the updated flag
168  updated_ = true;
169  }
170 }
171 
172 
174 {
176 
177  // Initialise the AMI
178  resetAMI();
179 }
180 
181 
183 {
184  static_cast<cyclicTransform&>(*this) =
186  (
187  name(),
188  faceAreas(),
189  *this,
190  nbrPatchName(),
191  nbrPatch(),
192  matchTolerance()
193  );
194 }
195 
196 
198 (
199  PstreamBuffers& pBufs,
200  const pointField& p
201 )
202 {
204 
205  // Initialise the AMI
206  resetAMI();
207 }
208 
209 
211 {
212  return srcMask_;
213 }
214 
215 
217 {
218  return tgtMask_;
219 }
220 
221 
222 // * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * //
223 
225 (
226  const word& name,
227  const label size,
228  const label start,
229  const label index,
230  const polyBoundaryMesh& bm,
231  const word& patchType
232 )
233 :
235  (
236  name,
237  size,
238  start,
239  index,
240  bm,
241  patchType,
242  false,
244  ),
245  nonOverlapPatchName_(word::null),
246  nonOverlapPatchID_(-1),
247  srcMask_(),
248  tgtMask_(),
249  updated_(false)
250 {
251  // Non-overlapping patch might not be valid yet so cannot determine
252  // associated patchID
253 }
254 
255 
257 (
258  const word& name,
259  const dictionary& dict,
260  const label index,
261  const polyBoundaryMesh& bm,
262  const word& patchType
263 )
264 :
266  (
267  name,
268  dict,
269  index,
270  bm,
271  patchType,
272  false,
274  ),
275  nonOverlapPatchName_(dict.lookup("nonOverlapPatch")),
276  nonOverlapPatchID_(-1),
277  srcMask_(),
278  tgtMask_(),
279  updated_(false)
280 {
281  if (nonOverlapPatchName_ == name)
282  {
284  (
285  dict
286  ) << "Non-overlapping patch name " << nonOverlapPatchName_
287  << " cannot be the same as this patch " << name
288  << exit(FatalIOError);
289  }
290 
291  // Non-overlapping patch might not be valid yet so cannot determine
292  // associated patchID
293 }
294 
295 
297 (
298  const cyclicACMIPolyPatch& pp,
299  const polyBoundaryMesh& bm
300 )
301 :
302  cyclicAMIPolyPatch(pp, bm),
303  nonOverlapPatchName_(pp.nonOverlapPatchName_),
304  nonOverlapPatchID_(-1),
305  srcMask_(),
306  tgtMask_(),
307  updated_(false)
308 {
309  // Non-overlapping patch might not be valid yet so cannot determine
310  // associated patchID
311 }
312 
313 
315 (
316  const cyclicACMIPolyPatch& pp,
317  const polyBoundaryMesh& bm,
318  const label index,
319  const label newSize,
320  const label newStart,
321  const word& nbrPatchName,
322  const word& nonOverlapPatchName
323 )
324 :
325  cyclicAMIPolyPatch(pp, bm, index, newSize, newStart, nbrPatchName),
326  nonOverlapPatchName_(nonOverlapPatchName),
327  nonOverlapPatchID_(-1),
328  srcMask_(),
329  tgtMask_(),
330  updated_(false)
331 {
332  if (nonOverlapPatchName_ == name())
333  {
335  << "Non-overlapping patch name " << nonOverlapPatchName_
336  << " cannot be the same as this patch " << name()
337  << exit(FatalError);
338  }
339 
340  // Non-overlapping patch might not be valid yet so cannot determine
341  // associated patchID
342 }
343 
344 
346 (
347  const cyclicACMIPolyPatch& pp,
348  const polyBoundaryMesh& bm,
349  const label index,
350  const labelUList& mapAddressing,
351  const label newStart
352 )
353 :
354  cyclicAMIPolyPatch(pp, bm, index, mapAddressing, newStart),
355  nonOverlapPatchName_(pp.nonOverlapPatchName_),
356  nonOverlapPatchID_(-1),
357  srcMask_(),
358  tgtMask_(),
359  updated_(false)
360 {
361  // Non-overlapping patch might not be valid yet so cannot determine
362  // associated patchID
363 }
364 
365 
366 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
367 
369 {}
370 
371 
372 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
373 
375 {
376  const polyPatch& pp = this->boundaryMesh()[nbrPatchID()];
377  return refCast<const cyclicACMIPolyPatch>(pp);
378 }
379 
380 
382 {
383  if (nonOverlapPatchID_ == -1)
384  {
385  nonOverlapPatchID_ =
386  this->boundaryMesh().findPatchID(nonOverlapPatchName_);
387 
388  if (nonOverlapPatchID_ == -1)
389  {
391  << "Illegal non-overlapping patch name " << nonOverlapPatchName_
392  << nl << "Valid patch names are "
393  << this->boundaryMesh().names()
394  << exit(FatalError);
395  }
396 
397  if (nonOverlapPatchID_ < index())
398  {
400  << "Boundary ordering error: " << type()
401  << " patch must be defined prior to its non-overlapping patch"
402  << nl
403  << type() << " patch: " << name() << ", ID:" << index() << nl
404  << "Non-overlap patch: " << nonOverlapPatchName_
405  << ", ID:" << nonOverlapPatchID_ << nl
406  << exit(FatalError);
407  }
408 
409  const polyPatch& noPp = this->boundaryMesh()[nonOverlapPatchID_];
410 
411  bool ok = true;
412 
413  if (size() == noPp.size())
414  {
415  const scalarField magSf(magFaceAreas());
416  const scalarField noMagSf(noPp.magFaceAreas());
417 
418  forAll(magSf, facei)
419  {
420  scalar ratio = mag(magSf[facei]/(noMagSf[facei] + rootVSmall));
421 
422  if (ratio - 1 > tolerance_)
423  {
424  ok = false;
425  break;
426  }
427  }
428  }
429  else
430  {
431  ok = false;
432  }
433 
434  if (!ok)
435  {
437  << "Inconsistent ACMI patches " << name() << " and "
438  << noPp.name() << ". Patches should have identical topology"
439  << exit(FatalError);
440  }
441  }
442 
443  return nonOverlapPatchID_;
444 }
445 
446 
448 (
449  PstreamBuffers& pBufs,
450  const primitivePatch& pp
451 ) const
452 {
454 }
455 
456 
458 (
459  PstreamBuffers& pBufs,
460  const primitivePatch& pp,
461  labelList& faceMap,
462  labelList& rotation
463 ) const
464 {
465  return cyclicAMIPolyPatch::order(pBufs, pp, faceMap, rotation);
466 }
467 
468 
470 {
472  writeEntry(os, "nonOverlapPatch", nonOverlapPatchName_);
473 }
474 
475 
476 // ************************************************************************* //
virtual void resetAMI() const
Reset the AMI interpolator.
virtual void initCalcGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual void initCalcGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
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.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
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:323
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:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
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.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
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)
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.
Cyclic patch for Arbitrary Mesh Interface (AMI)
cyclicACMIPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType)
Construct from (base couped patch) components.
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
const vectorField::subField faceAreas() const
Return face areas.
Definition: polyPatch.C:303
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const scalarListList & srcWeights() const
Return const access to source patch weights.
static const char nl
Definition: Ostream.H:260
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
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 const cyclicACMIPolyPatch & nbrPatch() const
Return a reference to the neighbour patch.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialise ordering for primitivePatch. Does not.
const scalarField::subField magFaceAreas() const
Return face area magnitudes.
Definition: polyPatch.C:309
virtual void resetAMI() const
Reset the AMI interpolator.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:335
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:105
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
Cyclic plane tranformation.
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
Initialise 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:844
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
IOerror FatalIOError