Go to the documentation of this file.
51 Foam::vector Foam::cyclicAMIPolyPatch::findFaceNormalMaxRadius
66 Info<<
"findFaceMaxRadius(const pointField&) : patch: " <<
name() <<
nl
67 <<
" rotFace = " << facei <<
nl
69 <<
" distance = " <<
Foam::sqrt(magRadSqr[facei])
90 <<
" has transform type " << transformTypeNames[
transform()]
91 <<
", neighbour patch " << neighbPatchName()
92 <<
" has transform type "
93 << neighbPatch().transformTypeNames[neighbPatch().transform()]
106 if (rotationAngleDefined_)
108 const tensor T(rotationAxis_*rotationAxis_);
112 0, -rotationAxis_.z(), rotationAxis_.y(),
113 rotationAxis_.z(), 0, -rotationAxis_.x(),
114 -rotationAxis_.y(), rotationAxis_.x(), 0
121 +
sin(rotationAngle_)*S
128 +
sin(-rotationAngle_)*S
133 const vector transformedAreaPos =
gSum(half1Areas & revTPos);
134 const vector transformedAreaNeg =
gSum(half1Areas & revTNeg);
136 const scalar magArea0 =
mag(area0) + ROOTVSMALL;
140 const scalar errorPos =
mag(transformedAreaPos + area0);
141 const scalar errorNeg =
mag(transformedAreaNeg + area0);
143 const scalar normErrorPos = errorPos/magArea0;
144 const scalar normErrorNeg = errorNeg/magArea0;
146 if (errorPos > errorNeg && normErrorNeg < matchTolerance())
149 rotationAngle_ *= -1;
156 const scalar areaError =
min(normErrorPos, normErrorNeg);
158 if (areaError > matchTolerance())
161 <<
"Patch areas are not consistent within "
162 << 100*matchTolerance()
163 <<
" % indicating a possible error in the specified "
164 <<
"angle of rotation" <<
nl
165 <<
" owner patch : " <<
name() <<
nl
166 <<
" neighbour patch : " << neighbPatch().name()
170 <<
" area error : " << 100*areaError <<
" %"
171 <<
" match tolerance : " << matchTolerance()
177 scalar theta =
radToDeg(rotationAngle_);
179 Pout<<
"cyclicAMIPolyPatch::calcTransforms: patch:"
181 <<
" Specified rotation:"
182 <<
" swept angle: " << theta <<
" [deg]"
183 <<
" reverse transform: " << revT
191 if (half0Ctrs.size())
193 n0 = findFaceNormalMaxRadius(half0Ctrs);
195 if (half1Ctrs.size())
197 n1 = -findFaceNormalMaxRadius(half1Ctrs);
200 reduce(n0, maxMagSqrOp<point>());
201 reduce(n1, maxMagSqrOp<point>());
211 (n0 ^ rotationAxis_),
217 (-n1 ^ rotationAxis_),
226 Pout<<
"cyclicAMIPolyPatch::calcTransforms: patch:"
228 <<
" Specified rotation:"
229 <<
" n0:" << n0 <<
" n1:" << n1
230 <<
" swept angle: " << theta <<
" [deg]"
231 <<
" reverse transform: " << revT
247 Pout<<
"cyclicAMIPolyPatch::calcTransforms : patch:" <<
name()
248 <<
" Specified translation : " << separationVector_
268 <<
" Assuming cyclic AMI pairs are colocated" <<
endl;
283 <<
" forwardT = " << forwardT() <<
nl
284 <<
" reverseT = " << reverseT() <<
nl
285 <<
" separation = " << separation() <<
nl
286 <<
" collocated = " << collocated() <<
nl <<
endl;
306 neighbPatch().meshPoints()
317 transformPosition(nbrPoints);
349 AMILowWeightCorrection_,
356 Pout<<
"cyclicAMIPolyPatch : " <<
name()
357 <<
" constructed AMI with " <<
nl
358 <<
" " <<
"srcAddress:" << AMIPtr_().srcAddress().size()
360 <<
" " <<
"tgAddress :" << AMIPtr_().tgtAddress().size()
373 half0Areas[facei] = half0[facei].areaNormal(half0.
points());
380 half1Areas[facei] = half1[facei].areaNormal(half1.
points());
394 Pout<<
"calcTransforms() : patch: " <<
name() <<
nl
395 <<
" forwardT = " << forwardT() <<
nl
396 <<
" reverseT = " << reverseT() <<
nl
397 <<
" separation = " << separation() <<
nl
398 <<
" collocated = " << collocated() <<
nl <<
endl;
486 const word& patchType,
494 rotationCentre_(
Zero),
495 rotationAngleDefined_(
false),
497 separationVector_(
Zero),
501 AMIRequireMatch_(
true),
502 AMILowWeightCorrection_(-1.0),
517 const word& patchType
525 rotationCentre_(
Zero),
526 rotationAngleDefined_(
false),
528 separationVector_(
Zero),
545 AMIRequireMatch_(
true),
553 <<
"No \"neighbourPatch\" or \"coupleGroup\" provided."
557 if (nbrPatchName_ ==
name)
560 <<
"Neighbour patch name " << nbrPatchName_
561 <<
" cannot be the same as this patch " <<
name
573 rotationAngleDefined_ =
true;
574 rotationAngle_ =
degToRad(rotationAngle_);
578 Info<<
"rotationAngle: " << rotationAngle_ <<
" [rad]"
583 scalar magRot =
mag(rotationAxis_);
587 <<
"Illegal rotationAxis " << rotationAxis_ <<
endl
588 <<
"Please supply a non-zero vector."
591 rotationAxis_ /= magRot;
645 const label newStart,
646 const word& nbrPatchName
650 nbrPatchName_(nbrPatchName),
666 if (nbrPatchName_ ==
name())
669 <<
"Neighbour patch name " << nbrPatchName_
670 <<
" cannot be the same as this patch " <<
name()
717 if (nbrPatchID_ == -1)
719 nbrPatchID_ = this->
boundaryMesh().findPatchID(neighbPatchName());
721 if (nbrPatchID_ == -1)
724 <<
"Illegal neighbourPatch name " << neighbPatchName()
725 <<
nl <<
"Valid patch names are "
732 refCast<const cyclicAMIPolyPatch>
740 <<
"Patch " <<
name()
741 <<
" specifies neighbour patch " << neighbPatchName()
742 <<
nl <<
" but that in return specifies "
753 return index() < neighbPatchID();
760 return refCast<const cyclicAMIPolyPatch>(pp);
767 const word surfType(surfDict_.lookupOrDefault<
word>(
"type",
"none"));
769 if (!surfPtr_.valid() && owner() && surfType !=
"none")
771 word surfName(surfDict_.lookupOrDefault(
"name",
name()));
801 <<
"AMI interpolator only available to owner patch"
805 if (!AMIPtr_.valid())
807 resetAMI(AMIMethod_);
818 return AMI().applyLowWeightCorrection();
822 return neighbPatch().AMI().applyLowWeightCorrection();
841 else if (separated())
872 forwardT().size() == 1
886 else if (separated())
890 separation().size() == 1
892 : separation()[facei]
910 reverseT().size() == 1
924 else if (separated())
928 separation().size() == 1
930 : separation()[facei]
948 reverseT().size() == 1
1005 reverseTransformPosition(prt, facei);
1008 reverseTransformDirection(nrt, facei);
1010 label nbrFacei = -1;
1014 nbrFacei = AMI().tgtPointFace
1025 nbrFacei = neighbPatch().AMI().srcPointFace
1047 if (!nbrPatchName_.empty())
1049 os.
writeEntry(
"neighbourPatch", nbrPatchName_);
1051 coupleGroup_.
write(os);
1057 os.
writeEntry(
"rotationAxis", rotationAxis_);
1058 os.
writeEntry(
"rotationCentre", rotationCentre_);
1060 if (rotationAngleDefined_)
1069 os.
writeEntry(
"separationVector", separationVector_);
1099 if (AMILowWeightCorrection_ > 0)
1101 os.
writeEntry(
"lowWeightCorrection", AMILowWeightCorrection_);
1104 if (!surfDict_.empty())
1106 surfDict_.
writeEntry(surfDict_.dictName(), os);
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
interpolationMethod
Enumeration specifying interpolation method.
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
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,...
static const Enum< interpolationMethod > interpolationMethodNames_
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)
const Field< PointType > & points() const
Return reference to global 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))
virtual void movePoints(PstreamBuffers &, const pointField &p)
Correct patches after moving points.
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
Initialize ordering for primitivePatch. Does not.
static constexpr const zero Zero
Global zero.
cyclicAMIPolyPatch(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.
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.
Tensor< scalar > tensor
Tensor of scalars.
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.
Various functors for unary and binary operations. Can be used for parallel combine-reduce operations ...
const dictionary surfDict_
Dictionary used during projection surface construction.
Field< tensor > tensorField
Specialisation of Field<T> for tensor.
const autoPtr< searchableSurface > & surfPtr() const
Return a reference to the projection surface.
List< bool > boolList
A List of bools.
const coupleGroupIdentifier coupleGroup_
Optional patchGroup to find neighbPatch.
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 lookupOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Base class for Arbitrary Mesh Interface (AMI) methods.
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 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 initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
vector separationVector_
Translation vector.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
messageStream Info
Information stream (uses stdout - output is on the master only)
A patch is a list of labels that address the faces in the global face list.
word name(const complex &c)
Return string representation of complex.
label pointFace(const label facei, const vector &n, point &p) const
Return face index on neighbour patch which shares point p.
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
virtual void movePoints(const Field< PointType > &)
Correct patch after moving points.
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 write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual void calcTransforms()
Recalculate the transformation tensors.
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?
Macros for easy insertion into run-time selection tables.
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.
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)
Output to file stream, using an OSstream.
virtual void resetAMI(const AMIPatchToPatchInterpolation::interpolationMethod &AMIMethod=AMIPatchToPatchInterpolation::imFaceAreaWeight) const
Reset the AMI interpolator.
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.
Templated 3D Vector derived from VectorSpace adding construction from 3 components,...
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 nbr side.
dimensionedScalar acos(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label ListType::const_reference const label start
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
const word & neighbPatchName() const
Neighbour patch name.
static const word null
An empty word.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
const AMIPatchToPatchInterpolation::interpolationMethod AMIMethod_
AMI method.
const Time & time() const
Return the top-level database.
const List< Face > & localFaces() const
Return patch faces addressing into local point list.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
PrimitivePatch< face, SubList, const pointField & > primitivePatch
Addressing for a faceList slice.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
const scalar AMILowWeightCorrection_
Low weight correction threshold for AMI.
vector point
Point is a vector.
bool AMIRequireMatch_
Flag to indicate that patches should match/overlap.
void setSize(const label newSize)
Alias for resize(const label)
const word & constant() const
Return constant name.
defineTypeNameAndDebug(combustionModel, 0)
const bool AMIReverse_
Flag to indicate that slave patch should be reversed for AMI.
const word & name() const
Return the patch name.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual ~cyclicAMIPolyPatch()
Destructor.
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)