Go to the documentation of this file.
51 Foam::vector Foam::cyclicAMIPolyPatch::findFaceNormalMaxRadius
62 label facei =
findMax(magRadSqr);
65 <<
"Patch: " <<
name() <<
nl
66 <<
" rotFace = " << facei <<
nl
68 <<
" distance = " <<
Foam::sqrt(magRadSqr[facei])
88 <<
" has transform type " << transformTypeNames[
transform()]
89 <<
", neighbour patch " << neighbPatchName()
90 <<
" has transform type "
91 << neighbPatch().transformTypeNames[neighbPatch().transform()]
104 if (rotationAngleDefined_)
106 const tensor T(rotationAxis_*rotationAxis_);
110 0, -rotationAxis_.z(), rotationAxis_.y(),
111 rotationAxis_.z(), 0, -rotationAxis_.x(),
112 -rotationAxis_.y(), rotationAxis_.x(), 0
119 +
sin(rotationAngle_)*S
126 +
sin(-rotationAngle_)*S
131 const vector transformedAreaPos =
gSum(half1Areas & revTPos);
132 const vector transformedAreaNeg =
gSum(half1Areas & revTNeg);
134 const scalar magArea0 =
mag(area0) + ROOTVSMALL;
138 const scalar errorPos =
mag(transformedAreaPos + area0);
139 const scalar errorNeg =
mag(transformedAreaNeg + area0);
141 const scalar normErrorPos = errorPos/magArea0;
142 const scalar normErrorNeg = errorNeg/magArea0;
144 if (errorPos > errorNeg && normErrorNeg < matchTolerance())
147 rotationAngle_ *= -1;
154 const scalar areaError =
min(normErrorPos, normErrorNeg);
156 if (areaError > matchTolerance())
159 <<
"Patch areas are not consistent within "
160 << 100*matchTolerance()
161 <<
" % indicating a possible error in the specified "
162 <<
"angle of rotation" <<
nl
163 <<
" owner patch : " <<
name() <<
nl
164 <<
" neighbour patch : " << neighbPatch().name()
168 <<
" area error : " << 100*areaError <<
" %"
169 <<
" match tolerance : " << matchTolerance()
175 scalar theta =
radToDeg(rotationAngle_);
177 Pout<<
"cyclicAMIPolyPatch::calcTransforms: patch:"
179 <<
" Specified rotation:"
180 <<
" swept angle: " << theta <<
" [deg]"
181 <<
" reverse transform: " << revT
189 if (half0Ctrs.size())
191 n0 = findFaceNormalMaxRadius(half0Ctrs);
193 if (half1Ctrs.size())
195 n1 = -findFaceNormalMaxRadius(half1Ctrs);
198 reduce(n0, maxMagSqrOp<point>());
199 reduce(n1, maxMagSqrOp<point>());
209 (n0 ^ rotationAxis_),
215 (-n1 ^ rotationAxis_),
224 Pout<<
"cyclicAMIPolyPatch::calcTransforms: patch:"
226 <<
" Specified rotation:"
227 <<
" n0:" << n0 <<
" n1:" << n1
228 <<
" swept angle: " << theta <<
" [deg]"
229 <<
" reverse transform: " << revT
245 Pout<<
"cyclicAMIPolyPatch::calcTransforms : patch:" <<
name()
246 <<
" Specified translation : " << separationVector_
266 <<
" Assuming cyclic AMI pairs are colocated" <<
endl;
281 <<
" forwardT = " << forwardT() <<
nl
282 <<
" reverseT = " << reverseT() <<
nl
283 <<
" separation = " << separation() <<
nl
284 <<
" collocated = " << collocated() <<
nl <<
endl;
294 const word surfType(surfDict_.getOrDefault<
word>(
"type",
"none"));
296 if (!surfPtr_ && owner() && surfType !=
"none")
298 word surfName(surfDict_.getOrDefault(
"name",
name()));
349 label patchSize0 = size();
350 label nbrPatchSize0 = nbr.size();
355 if (srcFaceIDs_.size())
357 patchSize0 = srcFaceIDs_.size();
359 if (tgtFaceIDs_.size())
361 nbrPatchSize0 = tgtFaceIDs_.size();
366 transformPosition(nbrPoints);
390 AMIPtr_->upToDate() =
false;
391 AMIPtr_->calculate(patch0, nbrPatch0, surfPtr());
395 AMIPtr_->checkSymmetricWeights(
true);
408 half0Areas[facei] = half0[facei].areaNormal(half0.
points());
415 half1Areas[facei] = half1[facei].areaNormal(half1.
points());
428 <<
"calcTransforms() : patch: " <<
name() <<
nl
429 <<
" forwardT = " << forwardT() <<
nl
430 <<
" reverseT = " << reverseT() <<
nl
431 <<
" separation = " << separation() <<
nl
432 <<
" collocated = " << collocated() <<
nl <<
endl;
441 AMIPtr_->upToDate() =
false;
481 restoreScaledGeometry();
487 AMIPtr_->upToDate() =
false;
547 AMIPtr_->upToDate() =
false;
563 const word& patchType,
565 const word& defaultAMIMethod
572 rotationCentre_(
Zero),
573 rotationAngleDefined_(
false),
575 separationVector_(
Zero),
579 createAMIFaces_(
false),
580 moveFaceCentres_(
false),
598 const word& patchType,
599 const word& defaultAMIMethod
607 rotationCentre_(
Zero),
608 rotationAngleDefined_(
false),
610 separationVector_(
Zero),
623 moveFaceCentres_(
false),
633 <<
"No \"neighbourPatch\" or \"coupleGroup\" provided."
637 if (nbrPatchName_ ==
name)
640 <<
"Neighbour patch name " << nbrPatchName_
641 <<
" cannot be the same as this patch " <<
name
653 rotationAngleDefined_ =
true;
654 rotationAngle_ =
degToRad(rotationAngle_);
658 Info<<
"rotationAngle: " << rotationAngle_ <<
" [rad]"
663 scalar magRot =
mag(rotationAxis_);
667 <<
"Illegal rotationAxis " << rotationAxis_ <<
endl
668 <<
"Please supply a non-zero vector."
671 rotationAxis_ /= magRot;
693 srcFaceIDs_.setSize(
dict.
get<label>(
"srcSize"));
694 tgtFaceIDs_.setSize(
dict.
get<label>(
"tgtSize"));
737 const label newStart,
738 const word& nbrPatchName
742 nbrPatchName_(nbrPatchName),
761 if (nbrPatchName_ ==
name())
764 <<
"Neighbour patch name " << nbrPatchName_
765 <<
" cannot be the same as this patch " <<
name()
809 if (nbrPatchID_ == -1)
811 nbrPatchID_ = this->
boundaryMesh().findPatchID(neighbPatchName());
813 if (nbrPatchID_ == -1)
816 <<
"Illegal neighbourPatch name " << neighbPatchName()
817 <<
nl <<
"Valid patch names are "
824 refCast<const cyclicAMIPolyPatch>
832 <<
"Patch " <<
name()
833 <<
" specifies neighbour patch " << neighbPatchName()
834 <<
nl <<
" but that in return specifies "
845 return index() < neighbPatchID();
852 return refCast<const cyclicAMIPolyPatch>(pp);
861 <<
"AMI interpolator only available to owner patch"
865 if (!AMIPtr_->upToDate())
878 return AMI().applyLowWeightCorrection();
882 return neighbPatch().AMI().applyLowWeightCorrection();
901 else if (separated())
932 forwardT().size() == 1
946 else if (separated())
950 separation().size() == 1
952 : separation()[facei]
970 reverseT().size() == 1
984 else if (separated())
988 separation().size() == 1
990 : separation()[facei]
1008 reverseT().size() == 1
1065 reverseTransformPosition(prt, facei);
1068 reverseTransformDirection(nrt, facei);
1070 label nbrFacei = -1;
1074 nbrFacei = AMI().tgtPointFace
1085 nbrFacei = neighbPatch().AMI().srcPointFace
1107 if (!nbrPatchName_.empty())
1109 os.
writeEntry(
"neighbourPatch", nbrPatchName_);
1111 coupleGroup_.
write(os);
1117 os.
writeEntry(
"rotationAxis", rotationAxis_);
1118 os.
writeEntry(
"rotationCentre", rotationCentre_);
1120 if (rotationAngleDefined_)
1129 os.
writeEntry(
"separationVector", separationVector_);
1144 if (!surfDict_.empty())
1146 surfDict_.writeEntry(surfDict_.dictName(), os);
1149 if (createAMIFaces_)
1151 os.
writeEntry(
"createAMIFaces", createAMIFaces_);
1152 os.
writeEntry(
"srcSize", srcFaceIDs_.size());
1153 os.
writeEntry(
"tgtSize", tgtFaceIDs_.size());
1154 os.
writeEntry(
"moveFaceCentres", moveFaceCentres_);
int debug
Static debugging option.
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
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)
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
An Ostream wrapper for parallel output to std::cout.
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 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.
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 (uses stdout - output is on the master only)
#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.
word name(const complex &c)
Return string representation of complex.
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 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?
const Field< point_type > & points() const
Return reference to global points.
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 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)
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.
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.
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.
void setSize(const label newSize)
Alias for resize(const label)
const word & constant() const
Return constant name.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
defineTypeNameAndDebug(combustionModel, 0)
const word & name() const
The patch name.
#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)