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-2015 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-6;
42 
43 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
44 
46 {
47  if
48  (
49  !empty()
50  && (faceAreas0_.empty() || boundaryMesh().mesh().moving())
51  )
52  {
53  faceAreas0_ = faceAreas();
54  }
55 
56  const cyclicACMIPolyPatch& nbrACMI =
57  refCast<const cyclicACMIPolyPatch>(this->neighbPatch());
58 
59  if
60  (
61  !nbrACMI.empty()
62  && (nbrACMI.faceAreas0().empty() || boundaryMesh().mesh().moving())
63  )
64  {
65  nbrACMI.faceAreas0_ = nbrACMI.faceAreas();
66  }
67 }
68 
69 
71 (
73 ) const
74 {
75  if (owner())
76  {
77  const polyPatch& nonOverlapPatch = this->nonOverlapPatch();
78 
79  initPatchFaceAreas();
80 
81  // Reset patch face areas based on original patch for AMI calculation
82  vectorField::subField Sf = faceAreas();
83  vectorField::subField noSf = nonOverlapPatch.faceAreas();
84 
85  forAll(Sf, faceI)
86  {
87  Sf[faceI] = faceAreas0_[faceI];
88  noSf[faceI] = faceAreas0_[faceI];
89  }
90 
91  // Calculate the AMI using partial face-area-weighted
93  (
95  );
96 
97  srcMask_ =
98  min(scalar(1) - tolerance_, max(tolerance_, AMI().srcWeightsSum()));
99 
100  tgtMask_ =
101  min(scalar(1) - tolerance_, max(tolerance_, AMI().tgtWeightsSum()));
102 
103  forAll(Sf, faceI)
104  {
105  Sf[faceI] *= srcMask_[faceI];
106  noSf[faceI] *= 1.0 - srcMask_[faceI];
107  }
108 
109  setNeighbourFaceAreas();
110 
111  // Set the updated flag
112  updated_ = true;
113  }
114 }
115 
116 
118 {
119  const cyclicACMIPolyPatch& cp =
120  refCast<const cyclicACMIPolyPatch>(this->neighbPatch());
121  const polyPatch& pp = cp.nonOverlapPatch();
122 
123  const vectorField& faceAreas0 = cp.faceAreas0();
124 
125  if (tgtMask_.size() == cp.size())
126  {
128  vectorField::subField noSf = pp.faceAreas();
129 
130  forAll(Sf, faceI)
131  {
132  Sf[faceI] = tgtMask_[faceI]*faceAreas0[faceI];
133  noSf[faceI] = (1.0 - tgtMask_[faceI])*faceAreas0[faceI];
134  }
135  }
136  else
137  {
138  WarningIn("cyclicACMIPolyPatch::setNeighbourFaceAreas() const")
139  << "Target mask size differs to that of the neighbour patch\n"
140  << " May occur when decomposing." << endl;
141  }
142 }
143 
144 
146 {
148 
149  // Initialise the AMI
150  resetAMI();
151 }
152 
153 
155 {
157 }
158 
159 
161 (
162  PstreamBuffers& pBufs,
163  const pointField& p
164 )
165 {
167 
168  // Initialise the AMI
169  resetAMI();
170 }
171 
172 
174 (
175  PstreamBuffers& pBufs,
176  const pointField& p
177 )
178 {
180 }
181 
182 
184 {
186 }
187 
188 
190 {
192 }
193 
194 
196 {
198 }
199 
200 
202 {
203  return srcMask_;
204 }
205 
206 
208 {
209  return tgtMask_;
210 }
211 
212 
213 // * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * * //
214 
216 (
217  const word& name,
218  const label size,
219  const label start,
220  const label index,
221  const polyBoundaryMesh& bm,
222  const word& patchType,
223  const transformType transform
224 )
225 :
226  cyclicAMIPolyPatch(name, size, start, index, bm, patchType, transform),
227  faceAreas0_(),
228  nonOverlapPatchName_(word::null),
229  nonOverlapPatchID_(-1),
230  srcMask_(),
231  tgtMask_(),
232  updated_(false)
233 {
234  AMIRequireMatch_ = false;
235 
236  // Non-overlapping patch might not be valid yet so cannot determine
237  // associated patchID
238 }
239 
240 
242 (
243  const word& name,
244  const dictionary& dict,
245  const label index,
246  const polyBoundaryMesh& bm,
247  const word& patchType
248 )
249 :
250  cyclicAMIPolyPatch(name, dict, index, bm, patchType),
251  faceAreas0_(),
252  nonOverlapPatchName_(dict.lookup("nonOverlapPatch")),
253  nonOverlapPatchID_(-1),
254  srcMask_(),
255  tgtMask_(),
256  updated_(false)
257 {
258  AMIRequireMatch_ = false;
259 
260  if (nonOverlapPatchName_ == name)
261  {
263  (
264  "cyclicACMIPolyPatch::cyclicACMIPolyPatch"
265  "("
266  "const word&, "
267  "const dictionary&, "
268  "const label, "
269  "const polyBoundaryMesh&, "
270  "const word&"
271  ")",
272  dict
273  ) << "Non-overlapping patch name " << nonOverlapPatchName_
274  << " cannot be the same as this patch " << name
275  << exit(FatalIOError);
276  }
277 
278  // Non-overlapping patch might not be valid yet so cannot determine
279  // associated patchID
280 }
281 
282 
284 (
285  const cyclicACMIPolyPatch& pp,
286  const polyBoundaryMesh& bm
287 )
288 :
289  cyclicAMIPolyPatch(pp, bm),
290  faceAreas0_(),
291  nonOverlapPatchName_(pp.nonOverlapPatchName_),
292  nonOverlapPatchID_(-1),
293  srcMask_(),
294  tgtMask_(),
295  updated_(false)
296 {
297  AMIRequireMatch_ = false;
298 
299  // Non-overlapping patch might not be valid yet so cannot determine
300  // associated patchID
301 }
302 
303 
305 (
306  const cyclicACMIPolyPatch& pp,
307  const polyBoundaryMesh& bm,
308  const label index,
309  const label newSize,
310  const label newStart,
311  const word& nbrPatchName,
312  const word& nonOverlapPatchName
313 )
314 :
315  cyclicAMIPolyPatch(pp, bm, index, newSize, newStart, nbrPatchName),
316  faceAreas0_(),
317  nonOverlapPatchName_(nonOverlapPatchName),
318  nonOverlapPatchID_(-1),
319  srcMask_(),
320  tgtMask_(),
321  updated_(false)
322 {
323  AMIRequireMatch_ = false;
324 
325  if (nonOverlapPatchName_ == name())
326  {
328  (
329  "const cyclicACMIPolyPatch& "
330  "const polyBoundaryMesh&, "
331  "const label, "
332  "const label, "
333  "const label, "
334  "const word&, "
335  "const word&"
336  ) << "Non-overlapping patch name " << nonOverlapPatchName_
337  << " cannot be the same as this patch " << name()
338  << exit(FatalError);
339  }
340 
341  // Non-overlapping patch might not be valid yet so cannot determine
342  // associated patchID
343 }
344 
345 
347 (
348  const cyclicACMIPolyPatch& pp,
349  const polyBoundaryMesh& bm,
350  const label index,
351  const labelUList& mapAddressing,
352  const label newStart
353 )
354 :
355  cyclicAMIPolyPatch(pp, bm, index, mapAddressing, newStart),
356  faceAreas0_(),
357  nonOverlapPatchName_(pp.nonOverlapPatchName_),
358  nonOverlapPatchID_(-1),
359  srcMask_(),
360  tgtMask_(),
361  updated_(false)
362 {
363  AMIRequireMatch_ = false;
364 }
365 
366 
367 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
368 
370 {}
371 
372 
373 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
374 
376 {
377  const polyPatch& pp = this->boundaryMesh()[neighbPatchID()];
378  return refCast<const cyclicACMIPolyPatch>(pp);
379 }
380 
381 
383 {
384  if (nonOverlapPatchID_ == -1)
385  {
386  nonOverlapPatchID_ =
387  this->boundaryMesh().findPatchID(nonOverlapPatchName_);
388 
389  if (nonOverlapPatchID_ == -1)
390  {
391  FatalErrorIn("cyclicPolyAMIPatch::neighbPatchID() const")
392  << "Illegal non-overlapping patch name " << nonOverlapPatchName_
393  << nl << "Valid patch names are "
394  << this->boundaryMesh().names()
395  << exit(FatalError);
396  }
397 
398  if (nonOverlapPatchID_ < index())
399  {
400  FatalErrorIn("cyclicPolyAMIPatch::neighbPatchID() const")
401  << "Boundary ordering error: " << type()
402  << " patch must be defined prior to its non-overlapping patch"
403  << nl
404  << type() << " patch: " << name() << ", ID:" << index() << nl
405  << "Non-overlap patch: " << nonOverlapPatchName_
406  << ", ID:" << nonOverlapPatchID_ << nl
407  << exit(FatalError);
408  }
409 
410  const polyPatch& noPp = this->boundaryMesh()[nonOverlapPatchID_];
411 
412  bool ok = true;
413 
414  if (size() == noPp.size())
415  {
416  const scalarField magSf(mag(faceAreas()));
417  const scalarField noMagSf(mag(noPp.faceAreas()));
418 
419  forAll(magSf, faceI)
420  {
421  scalar ratio = mag(magSf[faceI]/(noMagSf[faceI] + ROOTVSMALL));
422 
423  if (ratio - 1 > tolerance_)
424  {
425  ok = false;
426  break;
427  }
428  }
429  }
430  else
431  {
432  ok = false;
433  }
434 
435  if (!ok)
436  {
438  (
439  "Foam::label "
440  "Foam::cyclicACMIPolyPatch::nonOverlapPatchID() const"
441  ) << "Inconsistent ACMI patches " << name() << " and "
442  << noPp.name() << ". Patches should have identical topology"
443  << exit(FatalError);
444  }
445  }
446 
447  return nonOverlapPatchID_;
448 }
449 
450 
452 (
453  const primitivePatch& referPatch,
454  const pointField& thisCtrs,
455  const vectorField& thisAreas,
456  const pointField& thisCc,
457  const pointField& nbrCtrs,
458  const vectorField& nbrAreas,
459  const pointField& nbrCc
460 )
461 {
463  (
464  referPatch,
465  thisCtrs,
466  thisAreas,
467  thisCc,
468  nbrCtrs,
469  nbrAreas,
470  nbrCc
471  );
472 }
473 
474 
476 (
477  PstreamBuffers& pBufs,
478  const primitivePatch& pp
479 ) const
480 {
482 }
483 
484 
486 (
487  PstreamBuffers& pBufs,
488  const primitivePatch& pp,
489  labelList& faceMap,
490  labelList& rotation
491 ) const
492 {
493  return cyclicAMIPolyPatch::order(pBufs, pp, faceMap, rotation);
494 }
495 
496 
498 {
500 
501  os.writeKeyword("nonOverlapPatch") << nonOverlapPatchName_
502  << token::END_STATEMENT << nl;
503 }
504 
505 
506 // ************************************************************************* //
virtual ~cyclicACMIPolyPatch()
Destructor.
virtual void setNeighbourFaceAreas() const
Set neighbour ACMI patch areas.
Pre-declare related SubField type.
Definition: Field.H:61
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
Cyclic patch for Arbitrary Mesh Interface (AMI)
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
const word & name() const
Return name.
const vectorField::subField faceAreas() const
Return face normals.
Definition: polyPatch.C:290
dimensioned< scalar > mag(const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void clearGeom()
Clear geometry.
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
virtual void resetAMI(const AMIPatchToPatchInterpolation::interpolationMethod &AMIMethod=AMIPatchToPatchInterpolation::imFaceAreaWeight) const
Reset the AMI interpolator.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
bool empty() const
Return true if the UList is empty (ie, size() is zero).
Definition: UListI.H:313
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
A class for handling words, derived from string.
Definition: word.H:59
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
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
virtual void movePoints(PstreamBuffers &pBufs, const pointField &)
Correct patches after moving points.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual label nonOverlapPatchID() const
Non-overlapping patch ID.
Namespace for OpenFOAM.
virtual const cyclicACMIPolyPatch & neighbPatch() const
Return a reference to the neighbour patch.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
interpolationMethod
Enumeration specifying interpolation method.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
static const char nl
Definition: Ostream.H:260
const double e
Elementary charge.
Definition: doubleFloat.H:78
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
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.
IOerror FatalIOError
#define WarningIn(functionName)
Report a warning using Foam::Warning.
Foam::polyBoundaryMesh.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual const scalarField & srcMask() const
Return the mask/weighting for the source patch.
Cyclic patch for Arbitrarily Coupled Mesh Interface (ACMI)
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const bMesh & mesh() const
Definition: boundaryMesh.H:202
#define forAll(list, i)
Definition: UList.H:421
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
Macros for easy insertion into run-time selection tables.
virtual void clearGeom()
Clear geometry.
const vectorField & faceAreas0() const
Return access to the original patch face areas.
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
A list of faces which address into the list of points.
virtual void resetAMI(const AMIPatchToPatchInterpolation::interpolationMethod &AMIMethod=AMIPatchToPatchInterpolation::imFaceAreaWeight) const
Reset the AMI interpolator.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
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 initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
static const scalar tolerance_
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const polyPatch & nonOverlapPatch() const
Return a const reference to the non-overlapping patch.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
bool cp(const fileName &src, const fileName &dst)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:756
virtual void movePoints(PstreamBuffers &pBufs, const pointField &)
Correct patches after moving points.
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
virtual void initPatchFaceAreas() const
Initialise patch face areas.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const word null
An empty word.
Definition: word.H:77
defineTypeNameAndDebug(combustionModel, 0)
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
virtual const scalarField & tgtMask() const
Return the mask/weighting for the target patch.