cyclicACMIPolyPatch.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 | Copyright (C) 2013-2016 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 (
48 ) const
49 {
50  if (owner())
51  {
52  const polyPatch& nonOverlapPatch = this->nonOverlapPatch();
53 
54  if (debug)
55  {
56  Pout<< "cyclicACMIPolyPatch::resetAMI : recalculating weights"
57  << " for " << name() << " and " << nonOverlapPatch.name()
58  << endl;
59  }
60 
61  if (boundaryMesh().mesh().hasCellCentres())
62  {
63  if (debug)
64  {
65  Pout<< "cyclicACMIPolyPatch::resetAMI : clearing cellCentres"
66  << " for " << name() << " and " << nonOverlapPatch.name()
67  << endl;
68  }
69 
70  //WarningInFunction
71  // << "The mesh already has cellCentres calculated when"
72  // << " resetting ACMI " << name() << "." << endl
73  // << "This is a problem since ACMI adapts the face areas"
74  // << " (to close cells) so this has" << endl
75  // << "to be done before cell centre calculation." << endl
76  // << "This can happen if e.g. the cyclicACMI is after"
77  // << " any processor patches in the boundary." << endl;
78  const_cast<polyMesh&>
79  (
80  boundaryMesh().mesh()
82  }
83 
84 
85  // Trigger re-building of faceAreas
86  (void)boundaryMesh().mesh().faceAreas();
87 
88 
89  // Calculate the AMI using partial face-area-weighted. This leaves
90  // the weights as fractions of local areas (sum(weights) = 1 means
91  // face is fully covered)
93  (
95  );
96 
98  const_cast<AMIPatchToPatchInterpolation&>(this->AMI());
99 
100  srcMask_ =
101  min(scalar(1) - tolerance_, max(tolerance_, AMI.srcWeightsSum()));
102 
103  tgtMask_ =
104  min(scalar(1) - tolerance_, max(tolerance_, AMI.tgtWeightsSum()));
105 
106 
107  // Adapt owner side areas. Note that in uncoupled situations (e.g.
108  // decomposePar) srcMask, tgtMask can be zero size.
109  if (srcMask_.size())
110  {
111  vectorField::subField Sf = faceAreas();
112  vectorField::subField noSf = nonOverlapPatch.faceAreas();
113 
114  forAll(Sf, facei)
115  {
116  Sf[facei] *= srcMask_[facei];
117  noSf[facei] *= 1.0 - srcMask_[facei];
118  }
119  }
120  // Adapt slave side areas
121  if (tgtMask_.size())
122  {
123  const cyclicACMIPolyPatch& cp =
124  refCast<const cyclicACMIPolyPatch>(this->neighbPatch());
125  const polyPatch& pp = cp.nonOverlapPatch();
126 
128  vectorField::subField noSf = pp.faceAreas();
129 
130  forAll(Sf, facei)
131  {
132  Sf[facei] *= tgtMask_[facei];
133  noSf[facei] *= 1.0 - tgtMask_[facei];
134  }
135  }
136 
137  // Re-normalise the weights since the effect of overlap is already
138  // accounted for in the area.
139  {
140  scalarListList& srcWeights = AMI.srcWeights();
141  scalarField& srcWeightsSum = AMI.srcWeightsSum();
142  forAll(srcWeights, i)
143  {
144  scalarList& wghts = srcWeights[i];
145  if (wghts.size())
146  {
147  scalar& sum = srcWeightsSum[i];
148 
149  forAll(wghts, j)
150  {
151  wghts[j] /= sum;
152  }
153  sum = 1.0;
154  }
155  }
156  }
157  {
158  scalarListList& tgtWeights = AMI.tgtWeights();
159  scalarField& tgtWeightsSum = AMI.tgtWeightsSum();
160  forAll(tgtWeights, i)
161  {
162  scalarList& wghts = tgtWeights[i];
163  if (wghts.size())
164  {
165  scalar& sum = tgtWeightsSum[i];
166  forAll(wghts, j)
167  {
168  wghts[j] /= sum;
169  }
170  sum = 1.0;
171  }
172  }
173  }
174 
175  // Set the updated flag
176  updated_ = true;
177  }
178 }
179 
180 
182 {
183  if (debug)
184  {
185  Pout<< "cyclicACMIPolyPatch::initGeometry : " << name() << endl;
186  }
187 
189 
190  // Initialise the AMI
191  resetAMI();
192 }
193 
194 
196 {
197  if (debug)
198  {
199  Pout<< "cyclicACMIPolyPatch::calcGeometry : " << name() << endl;
200  }
202 }
203 
204 
206 (
207  PstreamBuffers& pBufs,
208  const pointField& p
209 )
210 {
211  if (debug)
212  {
213  Pout<< "cyclicACMIPolyPatch::initMovePoints : " << name() << endl;
214  }
216 
217  // Initialise the AMI
218  resetAMI();
219 }
220 
221 
223 (
224  PstreamBuffers& pBufs,
225  const pointField& p
226 )
227 {
228  if (debug)
229  {
230  Pout<< "cyclicACMIPolyPatch::movePoints : " << name() << endl;
231  }
233 }
234 
235 
237 {
238  if (debug)
239  {
240  Pout<< "cyclicACMIPolyPatch::initUpdateMesh : " << name() << endl;
241  }
243 }
244 
245 
247 {
248  if (debug)
249  {
250  Pout<< "cyclicACMIPolyPatch::updateMesh : " << name() << endl;
251  }
253 }
254 
255 
257 {
258  if (debug)
259  {
260  Pout<< "cyclicACMIPolyPatch::clearGeom : " << name() << endl;
261  }
263 }
264 
265 
267 {
268  return srcMask_;
269 }
270 
271 
273 {
274  return tgtMask_;
275 }
276 
277 
278 // * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * //
279 
281 (
282  const word& name,
283  const label size,
284  const label start,
285  const label index,
286  const polyBoundaryMesh& bm,
287  const word& patchType,
288  const transformType transform
289 )
290 :
291  cyclicAMIPolyPatch(name, size, start, index, bm, patchType, transform),
292  nonOverlapPatchName_(word::null),
293  nonOverlapPatchID_(-1),
294  srcMask_(),
295  tgtMask_(),
296  updated_(false)
297 {
298  AMIRequireMatch_ = false;
299 
300  // Non-overlapping patch might not be valid yet so cannot determine
301  // associated patchID
302 }
303 
304 
306 (
307  const word& name,
308  const dictionary& dict,
309  const label index,
310  const polyBoundaryMesh& bm,
311  const word& patchType
312 )
313 :
314  cyclicAMIPolyPatch(name, dict, index, bm, patchType),
315  nonOverlapPatchName_(dict.lookup("nonOverlapPatch")),
316  nonOverlapPatchID_(-1),
317  srcMask_(),
318  tgtMask_(),
319  updated_(false)
320 {
321  AMIRequireMatch_ = false;
322 
323  if (nonOverlapPatchName_ == name)
324  {
326  (
327  dict
328  ) << "Non-overlapping patch name " << nonOverlapPatchName_
329  << " cannot be the same as this patch " << name
330  << exit(FatalIOError);
331  }
332 
333  // Non-overlapping patch might not be valid yet so cannot determine
334  // associated patchID
335 }
336 
337 
339 (
340  const cyclicACMIPolyPatch& pp,
341  const polyBoundaryMesh& bm
342 )
343 :
344  cyclicAMIPolyPatch(pp, bm),
345  nonOverlapPatchName_(pp.nonOverlapPatchName_),
346  nonOverlapPatchID_(-1),
347  srcMask_(),
348  tgtMask_(),
349  updated_(false)
350 {
351  AMIRequireMatch_ = false;
352 
353  // Non-overlapping patch might not be valid yet so cannot determine
354  // associated patchID
355 }
356 
357 
359 (
360  const cyclicACMIPolyPatch& pp,
361  const polyBoundaryMesh& bm,
362  const label index,
363  const label newSize,
364  const label newStart,
365  const word& nbrPatchName,
366  const word& nonOverlapPatchName
367 )
368 :
369  cyclicAMIPolyPatch(pp, bm, index, newSize, newStart, nbrPatchName),
370  nonOverlapPatchName_(nonOverlapPatchName),
371  nonOverlapPatchID_(-1),
372  srcMask_(),
373  tgtMask_(),
374  updated_(false)
375 {
376  AMIRequireMatch_ = false;
377 
378  if (nonOverlapPatchName_ == name())
379  {
381  << "Non-overlapping patch name " << nonOverlapPatchName_
382  << " cannot be the same as this patch " << name()
383  << exit(FatalError);
384  }
385 
386  // Non-overlapping patch might not be valid yet so cannot determine
387  // associated patchID
388 }
389 
390 
392 (
393  const cyclicACMIPolyPatch& pp,
394  const polyBoundaryMesh& bm,
395  const label index,
396  const labelUList& mapAddressing,
397  const label newStart
398 )
399 :
400  cyclicAMIPolyPatch(pp, bm, index, mapAddressing, newStart),
401  nonOverlapPatchName_(pp.nonOverlapPatchName_),
402  nonOverlapPatchID_(-1),
403  srcMask_(),
404  tgtMask_(),
405  updated_(false)
406 {
407  AMIRequireMatch_ = false;
408 }
409 
410 
411 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
412 
414 {}
415 
416 
417 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
418 
420 {
421  const polyPatch& pp = this->boundaryMesh()[neighbPatchID()];
422  return refCast<const cyclicACMIPolyPatch>(pp);
423 }
424 
425 
427 {
428  if (nonOverlapPatchID_ == -1)
429  {
430  nonOverlapPatchID_ =
431  this->boundaryMesh().findPatchID(nonOverlapPatchName_);
432 
433  if (nonOverlapPatchID_ == -1)
434  {
436  << "Illegal non-overlapping patch name " << nonOverlapPatchName_
437  << nl << "Valid patch names are "
438  << this->boundaryMesh().names()
439  << exit(FatalError);
440  }
441 
442  if (nonOverlapPatchID_ < index())
443  {
445  << "Boundary ordering error: " << type()
446  << " patch must be defined prior to its non-overlapping patch"
447  << nl
448  << type() << " patch: " << name() << ", ID:" << index() << nl
449  << "Non-overlap patch: " << nonOverlapPatchName_
450  << ", ID:" << nonOverlapPatchID_ << nl
451  << exit(FatalError);
452  }
453 
454  const polyPatch& noPp = this->boundaryMesh()[nonOverlapPatchID_];
455 
456  bool ok = true;
457 
458  if (size() == noPp.size())
459  {
460  const scalarField magSf(mag(faceAreas()));
461  const scalarField noMagSf(mag(noPp.faceAreas()));
462 
463  forAll(magSf, facei)
464  {
465  scalar ratio = mag(magSf[facei]/(noMagSf[facei] + ROOTVSMALL));
466 
467  if (ratio - 1 > tolerance_)
468  {
469  ok = false;
470  break;
471  }
472  }
473  }
474  else
475  {
476  ok = false;
477  }
478 
479  if (!ok)
480  {
482  << "Inconsistent ACMI patches " << name() << " and "
483  << noPp.name() << ". Patches should have identical topology"
484  << exit(FatalError);
485  }
486  }
487 
488  return nonOverlapPatchID_;
489 }
490 
491 
493 (
494  const primitivePatch& referPatch,
495  const pointField& thisCtrs,
496  const vectorField& thisAreas,
497  const pointField& thisCc,
498  const pointField& nbrCtrs,
499  const vectorField& nbrAreas,
500  const pointField& nbrCc
501 )
502 {
504  (
505  referPatch,
506  thisCtrs,
507  thisAreas,
508  thisCc,
509  nbrCtrs,
510  nbrAreas,
511  nbrCc
512  );
513 }
514 
515 
517 (
518  PstreamBuffers& pBufs,
519  const primitivePatch& pp
520 ) const
521 {
523 }
524 
525 
527 (
528  PstreamBuffers& pBufs,
529  const primitivePatch& pp,
530  labelList& faceMap,
531  labelList& rotation
532 ) const
533 {
534  return cyclicAMIPolyPatch::order(pBufs, pp, faceMap, rotation);
535 }
536 
537 
539 {
541 
542  os.writeKeyword("nonOverlapPatch") << nonOverlapPatchName_
543  << token::END_STATEMENT << nl;
544 }
545 
546 
547 // ************************************************************************* //
const vectorField::subField faceAreas() const
Return face normals.
Definition: polyPatch.C:290
virtual void clearGeom()
Clear geometry.
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
Cyclic patch for Arbitrarily Coupled Mesh Interface (ACMI)
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
const double e
Elementary charge.
Definition: doubleFloat.H:78
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const scalarListList & srcWeights() const
Return const access to source patch weights.
virtual label nonOverlapPatchID() const
Non-overlapping patch ID.
virtual void resetAMI(const AMIPatchToPatchInterpolation::interpolationMethod &AMIMethod=AMIPatchToPatchInterpolation::imFaceAreaWeight) const
Reset the AMI interpolator.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
virtual const cyclicACMIPolyPatch & neighbPatch() const
Return a reference to the neighbour patch.
virtual const scalarField & srcMask() const
Return the mask/weighting for the source patch.
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:253
bool cp(const fileName &src, const fileName &dst)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:620
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
Pre-declare related SubField type.
Definition: Field.H:61
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)
interpolationMethod
Enumeration specifying interpolation method.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
A list of faces which address into the list of points.
virtual void movePoints(PstreamBuffers &pBufs, const pointField &)
Correct patches after moving points.
virtual const scalarField & tgtMask() const
Return the mask/weighting for the target patch.
void clearGeom()
Clear geometry.
dynamicFvMesh & mesh
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
A class for handling words, derived from string.
Definition: word.H:59
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
static const word null
An empty word.
Definition: word.H:77
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)
const word & name() const
Return name.
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
virtual void resetAMI(const AMIPatchToPatchInterpolation::interpolationMethod &AMIMethod=AMIPatchToPatchInterpolation::imFaceAreaWeight) const
Reset the AMI interpolator.
Foam::polyBoundaryMesh.
virtual void clearGeom()
Clear geometry.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:262
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
defineTypeNameAndDebug(combustionModel, 0)
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
const bMesh & mesh() const
Definition: boundaryMesh.H:202
Buffers for inter-processor communications streams (UOPstream, UIPstream).
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
const scalarField & srcWeightsSum() const
Return const access to normalisation factor of source.
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
virtual ~cyclicACMIPolyPatch()
Destructor.
dimensioned< scalar > mag(const dimensioned< Type > &)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const scalarField & tgtWeightsSum() const
Return const access to normalisation factor of target.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const scalarListList & tgtWeights() const
Return const access to target patch weights.
const polyPatch & nonOverlapPatch() const
Return a const reference to the non-overlapping patch.
Namespace for OpenFOAM.
virtual void movePoints(PstreamBuffers &pBufs, const pointField &)
Correct patches after moving points.
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451
IOerror FatalIOError