Go to the documentation of this file.
50 cyclicPeriodicAMIPolyPatch,
58 void Foam::cyclicPeriodicAMIPolyPatch::syncTransforms()
const
82 refCast<const coupledPolyPatch>
89 if (
returnReduce((size() && !periodicPatch.size()), orOp<bool>()))
91 if (periodicPatch.separation().size() > 1)
95 <<
" has non-uniform separation vector "
96 << periodicPatch.separation()
97 <<
"This is not allowed inside " <<
type()
98 <<
" patch " <<
name()
102 if (periodicPatch.forwardT().size() > 1)
106 <<
" has non-uniform transformation tensor "
107 << periodicPatch.forwardT()
108 <<
"This is not allowed inside " <<
type()
109 <<
" patch " <<
name()
125 && periodicPatch.parallel()
127 reduce(isParallel, orOp<bool>());
144 if (!periodicPatch.size())
148 if (sep[procI].size())
152 periodicPatch.separation()
156 periodicPatch.collocated()
179 if (!periodicPatch.size())
187 periodicPatch.forwardT()
191 periodicPatch.reverseT()
204 void Foam::cyclicPeriodicAMIPolyPatch::writeOBJ
231 void Foam::cyclicPeriodicAMIPolyPatch::resetAMI()
const
236 const coupledPolyPatch& periodicPatch
238 refCast<const coupledPolyPatch>
240 boundaryMesh()[periodicPatchID()]
249 pointField nbrPoints0(neighbPatch().localPoints());
250 transformPosition(nbrPoints0);
256 if (nTransforms_ > nSectors_/2)
258 nTransforms_ -= nSectors_;
260 else if (nTransforms_ < - nSectors_/2)
262 nTransforms_ += nSectors_;
267 for (label i = 0; i < nTransforms_; ++ i)
269 periodicPatch.transformPosition(thisPoints0);
271 for (label i = 0; i > nTransforms_; -- i)
273 periodicPatch.transformPosition(nbrPoints0);
276 autoPtr<OBJstream> ownStr;
277 autoPtr<OBJstream> neiStr;
280 const Time&
runTime = boundaryMesh().mesh().
time();
285 ownStr.reset(
new OBJstream(dir/
name() + postfix));
286 neiStr.reset(
new OBJstream(dir/neighbPatch().
name() + postfix));
289 <<
"patch:" <<
name()
290 <<
" writing accumulated AMI to " << ownStr().name()
291 <<
" and " << neiStr().name() <<
endl;
303 SubList<face>(localFaces(), size()),
309 SubList<face>(localFaces(), size()),
315 SubList<face>(neighbPatch().localFaces(), neighbPatch().size()),
321 SubList<face>(neighbPatch().localFaces(), neighbPatch().size()),
326 AMIPtr_->setRequireMatch(
false);
327 AMIPtr_->calculate(thisPatch0, nbrPatch0, surfPtr());
331 label nTransformsOld(nTransforms_);
344 scalar srcSum(
gAverage(AMIPtr_->srcWeightsSum()));
345 scalar tgtSum(
gAverage(AMIPtr_->tgtWeightsSum()));
356 scalar srcSumDiff = 0;
359 <<
"patch:" <<
name()
360 <<
" srcSum:" << srcSum
361 <<
" tgtSum:" << tgtSum
370 (1 - srcSum > matchTolerance())
371 || (1 - tgtSum > matchTolerance())
377 periodicPatch.transformPosition(thisPoints);
380 <<
"patch:" <<
name()
381 <<
" moving this side from:"
385 thisPatch.movePoints(thisPoints);
388 <<
"patch:" <<
name()
389 <<
" appending weights with untransformed slave side"
392 AMIPtr_->append(thisPatch, nbrPatch0);
401 periodicPatch.transformPosition(nbrPoints);
404 <<
"patch:" <<
name()
405 <<
" moving neighbour side from:"
409 nbrPatch.movePoints(nbrPoints);
411 AMIPtr_->append(thisPatch0, nbrPatch);
419 const scalar srcSumNew =
gAverage(AMIPtr_->srcWeightsSum());
420 const scalar srcSumDiffNew = srcSumNew - srcSum;
422 if (srcSumDiffNew < srcSumDiff || srcSumDiffNew < SMALL)
426 srcSumDiff = srcSumDiffNew;
430 tgtSum =
gAverage(AMIPtr_->tgtWeightsSum());
437 <<
"patch:" <<
name()
438 <<
" iteration:" << iter
439 <<
" srcSum:" << srcSum
440 <<
" tgtSum:" << tgtSum
447 ownStr.reset(
nullptr);
448 neiStr.reset(
nullptr);
452 nTransforms_ = (nTransforms_ + nTransformsOld)/2;
455 if (iter == maxIter_)
466 <<
"Patches " <<
name() <<
" and " << neighbPatch().name()
467 <<
" do not couple to within a tolerance of "
469 <<
" when transformed according to the periodic patch "
470 << periodicPatch.name() <<
"." <<
nl
471 <<
"The current sum of weights are for owner " <<
name()
472 <<
" : " << srcSum <<
" and for neighbour "
473 << neighbPatch().name() <<
" : " << tgtSum <<
nl
474 <<
"This is only acceptable during post-processing"
475 <<
"; not during running. Improve your mesh or increase"
476 <<
" the 'matchTolerance' setting in the patch specification."
483 mag(srcSum - floor(srcSum + 0.5)) > srcSum*matchTolerance()
484 ||
mag(tgtSum - floor(tgtSum + 0.5)) > tgtSum*matchTolerance()
493 <<
"Patches " <<
name() <<
" and " << neighbPatch().name()
494 <<
" do not overlap an integer number of times when transformed"
495 <<
" according to the periodic patch "
496 << periodicPatch.name() <<
"." <<
nl
497 <<
"The current matchTolerance : " << matchTolerance()
498 <<
", sum of owner weights : " << srcSum
499 <<
", sum of neighbour weights : " << tgtSum
501 <<
"This is only acceptable during post-processing"
502 <<
"; not during running. Improve your mesh or increase"
503 <<
" the 'matchTolerance' setting in the patch specification."
508 const label nFace =
returnReduce(size(), sumOp<label>());
515 srcWghtSum[faceI] =
sum(AMIPtr_->srcWeights()[faceI]);
520 tgtWghtSum[faceI] =
sum(AMIPtr_->tgtWeights()[faceI]);
524 <<
"AMI: Patch " <<
name()
526 <<
" min:" <<
gMin(srcWghtSum)
527 <<
" max:" <<
gMax(srcWghtSum)
530 <<
"AMI: Patch " << neighbPatch().
name()
532 <<
" min:" <<
gMin(tgtWghtSum)
533 <<
" max:" <<
gMax(tgtWghtSum)
549 const word& patchType,
562 faceAreaWeightAMI::typeName
568 AMIPtr_->setRequireMatch(
false);
578 const word& patchType
588 faceAreaWeightAMI::typeName
594 AMIPtr_->setRequireMatch(
false);
605 nTransforms_(pp.nTransforms_),
606 nSectors_(pp.nSectors_),
607 maxIter_(pp.maxIter_)
609 AMIPtr_->setRequireMatch(
false);
619 const label newStart,
620 const word& nbrPatchName
624 nTransforms_(pp.nTransforms_),
625 nSectors_(pp.nSectors_),
626 maxIter_(pp.maxIter_)
628 AMIPtr_->setRequireMatch(
false);
642 nTransforms_(pp.nTransforms_),
643 nSectors_(pp.nSectors_),
644 maxIter_(pp.maxIter_)
646 AMIPtr_->setRequireMatch(
false);
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
#define InfoInFunction
Report an information message using Foam::Info.
virtual bool owner() const
Does this side own the patch?
A class for handling words, derived from Foam::string.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
static constexpr const zero Zero
Global zero (0)
Type gAverage(const FieldField< Field, Type > &f)
word periodicPatchName_
Periodic patch name.
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
virtual const tensorField & forwardT() const
Return face transformation tensor.
static word timeName(const scalar t, const int precision=precision_)
virtual const tensorField & reverseT() const
Return neighbour-cell transformation tensor.
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
static bool master(const label communicator=worldComm)
Am I the master process.
List< bool > boolList
A List of bools.
Ostream & endl(Ostream &os)
Add newline and flush stream.
virtual ~cyclicPeriodicAMIPolyPatch()
Destructor.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
messageStream Info
Information stream (stdout output on master, null elsewhere)
#define DebugInFunction
Report an information message using Foam::Info.
#define SeriousErrorInFunction
Report an error message using Foam::SeriousError.
virtual const fileName & name() const
Return the name of the stream.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
Ostream & indent(Ostream &os)
Indent stream.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
static int myProcNo(const label communicator=worldComm)
Number of this process (starting from masterNo() = 0)
cyclicPeriodicAMIPolyPatch(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 coupled patch) components.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
List< face > faceList
A List of faces.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
fileName globalPath() const
Return global path for the case.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field.
const word & name() const noexcept
The patch name.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Type gMin(const FieldField< Field, Type > &f)
Cyclic patch for periodic Arbitrary Mesh Interface (AMI)
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
label periodicPatchID() const
Periodic patch ID (or -1)
coupledPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType, const transformType transform)
Construct from components.
defineTypeNameAndDebug(combustionModel, 0)
const Time & time() const noexcept
Return time registry.
Type gMax(const FieldField< Field, Type > &f)
static label nProcs(const label communicator=worldComm)
Number of processes in parallel run, and 1 for serial run.
Cyclic patch for Arbitrary Mesh Interface (AMI)