Go to the documentation of this file.
52 Foam::vector Foam::cyclicAMIPolyPatch::findFaceNormalMaxRadius
63 label facei =
findMax(magRadSqr);
66 <<
"Patch: " <<
name() <<
nl
67 <<
" rotFace = " << facei <<
nl
69 <<
" distance = " <<
Foam::sqrt(magRadSqr[facei])
89 <<
" has transform type " << transformTypeNames[
transform()]
90 <<
", neighbour patch " << neighbPatchName()
91 <<
" has transform type "
92 << neighbPatch().transformTypeNames[neighbPatch().transform()]
105 if (rotationAngleDefined_)
107 const tensor T(rotationAxis_*rotationAxis_);
111 0, -rotationAxis_.z(), rotationAxis_.y(),
112 rotationAxis_.z(), 0, -rotationAxis_.x(),
113 -rotationAxis_.y(), rotationAxis_.x(), 0
120 +
sin(rotationAngle_)*S
127 +
sin(-rotationAngle_)*S
132 const vector transformedAreaPos =
gSum(half1Areas & revTPos);
133 const vector transformedAreaNeg =
gSum(half1Areas & revTNeg);
135 const scalar magArea0 =
mag(area0) + ROOTVSMALL;
139 const scalar errorPos =
mag(transformedAreaPos + area0);
140 const scalar errorNeg =
mag(transformedAreaNeg + area0);
142 const scalar normErrorPos = errorPos/magArea0;
143 const scalar normErrorNeg = errorNeg/magArea0;
145 if (errorPos > errorNeg && normErrorNeg < matchTolerance())
148 rotationAngle_ *= -1;
155 const scalar areaError =
min(normErrorPos, normErrorNeg);
157 if (areaError > matchTolerance())
160 <<
"Patch areas are not consistent within "
161 << 100*matchTolerance()
162 <<
" % indicating a possible error in the specified "
163 <<
"angle of rotation" <<
nl
164 <<
" owner patch : " <<
name() <<
nl
165 <<
" neighbour patch : " << neighbPatch().name()
169 <<
" area error : " << 100*areaError <<
" %"
170 <<
" match tolerance : " << matchTolerance()
176 scalar theta =
radToDeg(rotationAngle_);
178 Pout<<
"cyclicAMIPolyPatch::calcTransforms: patch:"
180 <<
" Specified rotation:"
181 <<
" swept angle: " << theta <<
" [deg]"
182 <<
" reverse transform: " << revT
190 if (half0Ctrs.size())
192 n0 = findFaceNormalMaxRadius(half0Ctrs);
194 if (half1Ctrs.size())
196 n1 = -findFaceNormalMaxRadius(half1Ctrs);
199 reduce(n0, maxMagSqrOp<point>());
200 reduce(n1, maxMagSqrOp<point>());
210 (n0 ^ rotationAxis_),
216 (-n1 ^ rotationAxis_),
225 Pout<<
"cyclicAMIPolyPatch::calcTransforms: patch:"
227 <<
" Specified rotation:"
228 <<
" n0:" << n0 <<
" n1:" << n1
229 <<
" swept angle: " << theta <<
" [deg]"
230 <<
" reverse transform: " << revT
246 Pout<<
"cyclicAMIPolyPatch::calcTransforms : patch:" <<
name()
247 <<
" Specified translation : " << separationVector_
267 <<
" Assuming cyclic AMI pairs are colocated" <<
endl;
282 <<
" forwardT = " << forwardT() <<
nl
283 <<
" reverseT = " << reverseT() <<
nl
284 <<
" separation = " << separation() <<
nl
285 <<
" collocated = " << collocated() <<
nl <<
endl;
295 const label periodicID = periodicPatchID();
296 if (periodicID != -1)
301 refCast<const coupledPolyPatch>
310 if (isA<cyclicPolyPatch>(perPp))
313 refCast<const cyclicPolyPatch>(perPp).rotationAxis();
315 refCast<const cyclicPolyPatch>(perPp).rotationCentre();
317 else if (isA<cyclicAMIPolyPatch>(perPp))
320 refCast<const cyclicAMIPolyPatch>(perPp).rotationAxis();
322 refCast<const cyclicAMIPolyPatch>(perPp).rotationCentre();
327 <<
" have unsupported periodicPatch " << perPp.
name()
343 const word surfType(surfDict_.getOrDefault<
word>(
"type",
"none"));
345 if (!surfPtr_ && owner() && surfType !=
"none")
347 word surfName(surfDict_.getOrDefault(
"name",
name()));
398 label patchSize0 = size();
399 label nbrPatchSize0 = nbr.size();
404 if (srcFaceIDs_.size())
406 patchSize0 = srcFaceIDs_.size();
408 if (tgtFaceIDs_.size())
410 nbrPatchSize0 = tgtFaceIDs_.size();
415 transformPosition(nbrPoints);
439 AMIPtr_->upToDate() =
false;
440 AMIPtr_->calculate(patch0, nbrPatch0, surfPtr());
444 AMIPtr_->checkSymmetricWeights(
true);
457 half0Areas[facei] = half0[facei].areaNormal(half0.
points());
464 half1Areas[facei] = half1[facei].areaNormal(half1.
points());
477 <<
"calcTransforms() : patch: " <<
name() <<
nl
478 <<
" forwardT = " << forwardT() <<
nl
479 <<
" reverseT = " << reverseT() <<
nl
480 <<
" separation = " << separation() <<
nl
481 <<
" collocated = " << collocated() <<
nl <<
endl;
490 AMIPtr_->upToDate() =
false;
530 restoreScaledGeometry();
536 AMIPtr_->upToDate() =
false;
596 AMIPtr_->upToDate() =
false;
612 const word& patchType,
614 const word& defaultAMIMethod
622 rotationCentre_(
Zero),
623 rotationAngleDefined_(
false),
625 separationVector_(
Zero),
627 periodicPatchID_(-1),
631 createAMIFaces_(
false),
632 moveFaceCentres_(
false),
650 const word& patchType,
651 const word& defaultAMIMethod
660 rotationCentre_(
Zero),
661 rotationAngleDefined_(
false),
663 separationVector_(
Zero),
665 periodicPatchID_(-1),
678 moveFaceCentres_(
false),
688 <<
"No \"neighbourPatch\" or \"coupleGroup\" provided."
692 if (nbrPatchName_ ==
name)
695 <<
"Neighbour patch name " << nbrPatchName_
696 <<
" cannot be the same as this patch " <<
name
708 rotationAngleDefined_ =
true;
709 rotationAngle_ =
degToRad(rotationAngle_);
713 Info<<
"rotationAngle: " << rotationAngle_ <<
" [rad]"
718 scalar magRot =
mag(rotationAxis_);
722 <<
"Illegal rotationAxis " << rotationAxis_ <<
endl
723 <<
"Please supply a non-zero vector."
726 rotationAxis_ /= magRot;
748 srcFaceIDs_.setSize(
dict.
get<label>(
"srcSize"));
749 tgtFaceIDs_.setSize(
dict.
get<label>(
"tgtSize"));
772 periodicPatchID_(-1),
795 const label newStart,
796 const word& nbrPatchName
800 nbrPatchName_(nbrPatchName),
810 periodicPatchID_(-1),
822 if (nbrPatchName_ ==
name())
825 <<
"Neighbour patch name " << nbrPatchName_
826 <<
" cannot be the same as this patch " <<
name()
855 periodicPatchID_(-1),
879 forAll (addSourceFaces, faceI)
881 const labelList& nbrFaceIs = addSourceFaces[faceI];
885 label nbrFaceI = nbrFaceIs[j];
887 if (nbrFaceI < neighbPatch().size())
904 if (nbrPatchID_ == -1)
906 nbrPatchID_ = this->
boundaryMesh().findPatchID(neighbPatchName());
908 if (nbrPatchID_ == -1)
911 <<
"Illegal neighbourPatch name " << neighbPatchName()
912 <<
nl <<
"Valid patch names are "
919 refCast<const cyclicAMIPolyPatch>
927 <<
"Patch " <<
name()
928 <<
" specifies neighbour patch " << neighbPatchName()
929 <<
nl <<
" but that in return specifies "
946 if (periodicPatchID_ == -1)
948 periodicPatchID_ =
boundaryMesh().findPatchID(periodicPatchName_);
950 if (periodicPatchID_ == -1)
953 <<
"Illegal neighbourPatch name " << periodicPatchName_
954 <<
nl <<
"Valid patch names are "
959 return periodicPatchID_;
966 return index() < neighbPatchID();
973 return refCast<const cyclicAMIPolyPatch>(pp);
982 <<
"AMI interpolator only available to owner patch"
986 if (!AMIPtr_->upToDate())
999 return AMI().applyLowWeightCorrection();
1003 return neighbPatch().AMI().applyLowWeightCorrection();
1022 else if (separated())
1053 forwardT().size() == 1
1067 else if (separated())
1071 separation().size() == 1
1073 : separation()[facei]
1091 reverseT().size() == 1
1105 else if (separated())
1109 separation().size() == 1
1111 : separation()[facei]
1129 reverseT().size() == 1
1186 reverseTransformPosition(prt, facei);
1189 reverseTransformDirection(nrt, facei);
1191 label nbrFacei = -1;
1195 nbrFacei = AMI().tgtPointFace
1206 nbrFacei = neighbPatch().AMI().srcPointFace
1228 if (!nbrPatchName_.empty())
1241 if (rotationAngleDefined_)
1270 if (!surfDict_.empty())
1272 surfDict_.writeEntry(surfDict_.dictName(),
os);
1275 if (createAMIFaces_)
int debug
Static debugging option.
const Field< point_type > & points() const noexcept
Return reference to global points.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
vectorField pointField
pointField is a vectorField.
const bMesh & mesh() const
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.
points setSize(newPointi)
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
virtual void resetAMI() const
Reset the AMI interpolator, use current patch points.
virtual bool owner() const
Does this side own the patch?
A class for handling words, derived from Foam::string.
A class for handling file names.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
word nbrPatchName_
Name of other half.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
static constexpr const zero Zero
Global zero (0)
word periodicPatchName_
Periodic patch name.
point rotationCentre_
Point on axis of rotation for rotational cyclics.
dimensionedScalar sin(const dimensionedScalar &ds)
A List obtained as a section of another List.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
virtual void clearGeom()
Clear geometry.
virtual void reverseTransformPosition(point &l, const label facei) const
Transform a patch-based position from this side to nbr side.
const dictionary surfDict_
Dictionary used during projection surface construction.
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
const autoPtr< searchableSurface > & surfPtr() const
Create and return pointer to the projection surface.
Unit conversion functions.
List< bool > boolList
A List of bools.
const coupleGroupIdentifier coupleGroup_
Optional patchGroup to find neighbPatch.
static const scalar tolerance_
Tolerance used e.g. for area calculations/limits.
Ostream & endl(Ostream &os)
Add newline and flush stream.
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
Type gSum(const FieldField< Field, Type > &f)
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
static autoPtr< AMIInterpolation > New(const word &modelName, const dictionary &dict, const bool reverseTarget=false)
Selector for dictionary.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
virtual void write(Ostream &os) const
Write the polyPatch data as a dictionary.
virtual bool parallel() const
Are the cyclic planes parallel.
virtual void transformPosition(pointField &) const
Transform patch-based positions from nbr side to this side.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
virtual void newInternalProcFaces(label &, label &) const
Return number of new internal of this polyPatch faces.
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
void set(T *p) noexcept
Deprecated(2018-02) Identical to reset().
vector separationVector_
Translation vector.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
autoPtr< AMIPatchToPatchInterpolation > AMIPtr_
AMI interpolation class.
messageStream Info
Information stream (stdout output on master, null elsewhere)
#define DebugInFunction
Report an information message using Foam::Info.
A patch is a list of labels that address the faces in the global face list.
void setSize(const label n)
Alias for resize()
label pointFace(const label facei, const vector &n, point &p) const
virtual void movePoints(PstreamBuffers &pBufs, const pointField &)
Correct patches after moving points.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
constexpr scalar radToDeg(const scalar rad) noexcept
Conversion from radians to degrees.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual void calcTransforms()
Recalculate the transformation tensors.
A cylindrical coordinate system (r-theta-z). The coordinate system angle theta is always in radians.
vector rotationAxis_
Axis of rotation for rotational cyclics.
scalar rotationAngle_
Rotation angle.
virtual label neighbPatchID() const
Neighbour patch ID.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual const cyclicAMIPolyPatch & neighbPatch() const
Return a reference to the neighbour patch.
static bool valid(char c)
Is this character valid for a word?
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
const List< face_type > & localFaces() const
Return patch faces addressing into local point list.
cyclicAMIPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType, const transformType transform=UNKNOWN, const word &defaultAMIMethod=faceAreaWeightAMI::typeName)
Construct from (base coupled patch) components.
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
bool rotationAngleDefined_
Flag to show whether the rotation angle is defined.
errorManip< error > abort(error &err)
Vector< scalar > vector
A scalar version of the templated Vector.
static autoPtr< searchableSurface > New(const word &surfaceType, const IOobject &io, const dictionary &dict)
Return a reference to the selected searchableSurface.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
#define DebugPout
Report an information message using Foam::Pout.
dictionary subOrEmptyDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX, const bool mandatory=false) const
errorManipArg< error, int > exit(error &err, const int errNo=1)
const Field< point_type > & localPoints() const
Return pointField of points in patch.
Output to file stream, using an OSstream.
bool moveFaceCentres_
Move face centres (default = no)
const scalar fraction_
Particle displacement fraction accross AMI.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
const AMIPatchToPatchInterpolation & AMI() const
Return a reference to the AMI interpolator.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
fileName path() const
Return path.
const vectorField::subField faceCentres() const
Return face centres.
bool applyLowWeightCorrection() const
Return true if applying the low weight correction.
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
virtual void clearGeom()
Clear geometry.
label findMax(const ListType &input, label start=0)
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual void reverseTransformDirection(vector &d, const label facei) const
Transform a patch-based direction from this side to.
dimensionedScalar acos(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
const word & neighbPatchName() const
Neighbour patch name.
const dimensionedScalar e
Elementary charge.
static const word null
An empty word.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const Time & time() const
Return the top-level database.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
virtual void movePoints(const Field< point_type > &)
Correct patch after moving points.
const word & name() const noexcept
The patch name.
bool createAMIFaces_
Flag to indicate that new AMI faces will created.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
vector point
Point is a vector.
const word & constant() const
Return constant name.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
autoPtr< coordSystem::cylindrical > cylindricalCS() const
Create a coordinate system from the periodic patch (or nullptr)
label periodicPatchID() const
Periodic patch ID (or -1)
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
Tensor< scalar > tensor
Tensor of scalars, i.e. Tensor<scalar>.
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
A list of faces which address into the list of points.
dimensionedScalar cos(const dimensionedScalar &ds)
Cyclic patch for Arbitrary Mesh Interface (AMI)