Go to the documentation of this file.
39 const scalar coupledPolyPatch::defaultMatchTol_ = 1
e-4;
49 { transformType::UNKNOWN,
"unknown" },
50 { transformType::ROTATIONAL,
"rotational" },
51 { transformType::TRANSLATIONAL,
"translational" },
52 { transformType::COINCIDENTFULLMATCH,
"coincidentFullMatch" },
53 { transformType::NOORDERING,
"noOrdering" },
61 os <<
"v " << pt.
x() <<
' ' << pt.
y() <<
' ' << pt.
z() <<
endl;
93 os <<
"l " << vertI-1 <<
' ' << vertI <<
nl;
112 const face&
f = faces[i];
116 if (foamToObj.insert(
f[fp], vertI))
126 os <<
' ' << foamToObj[
f[fp]]+1;
128 os <<
' ' << foamToObj[
f[0]]+1 <<
nl;
147 anchors[facei] =
points[faces[facei][0]];
155 const face&
f = faces[facei];
165 for (
label fp2 = 0; fp2 <
f.size(); ++fp2)
167 if (
f[fp1] ==
f[fp2])
192 anchors[facei] =
points[faces[facei][0]];
213 const point& cc = faceCentres[facei];
215 const face&
f = faces[facei];
219 scalar maxLenSqr = -GREAT;
222 scalar maxCmpt = -GREAT;
227 maxLenSqr =
max(maxLenSqr,
magSqr(pt - cc));
250 scalar minDistSqr = GREAT;
256 if (distSqr < minDistSqr)
258 minDistSqr = distSqr;
263 if (anchorFp == -1 ||
Foam::sqrt(minDistSqr) > tol)
274 if (distSqr == minDistSqr && fp != anchorFp)
277 <<
"Cannot determine unique anchor point on face "
280 <<
"Both at index " << anchorFp <<
" and " << fp
281 <<
" the vertices have the same distance "
283 <<
" to the anchor " << anchor
284 <<
". Continuing but results might be wrong."
290 return (
f.size() - anchorFp) %
f.size();
308 Pout<<
"coupledPolyPatch::calcTransformTensors : " <<
name() <<
endl
309 <<
" transform:" << transformTypeNames[
transform] <<
nl
310 <<
" (half)size:" << Cf.size() <<
nl
311 <<
" absTol:" << absTol <<
nl
312 <<
" smallDist min:" <<
min(smallDist) <<
nl
313 <<
" smallDist max:" <<
max(smallDist) <<
nl
314 <<
" sum(mag(nf & nr)):" <<
sum(
mag(nf & nr)) <<
endl;
328 separation_.setSize(0);
331 collocated_.setSize(0);
356 separation_.setSize(0);
358 forwardT_.setSize(Cf.size());
359 reverseT_.setSize(Cf.size());
360 collocated_.setSize(Cf.size());
371 Pout<<
" sum(mag(forwardT_ - forwardT_[0])):"
372 <<
sum(
mag(forwardT_ - forwardT_[0]))
378 forwardT_.setSize(1);
379 reverseT_.setSize(1);
380 collocated_.setSize(1);
384 Pout<<
" difference in rotation less than"
385 <<
" local tolerance "
386 <<
error <<
". Assuming uniform rotation." <<
endl;
394 forwardT_.setSize(0);
395 reverseT_.setSize(0);
397 separation_ = Cr - Cf;
399 collocated_.setSize(separation_.size());
407 bool sameSeparation =
true;
408 bool doneWarning =
false;
410 forAll(separation_, facei)
412 scalar smallSqr =
sqr(smallDist[facei]);
414 collocated_[facei] = (
magSqr(separation_[facei]) < smallSqr);
417 if (
magSqr(separation_[facei] - separation_[0]) > smallSqr)
419 sameSeparation =
false;
421 if (!doneWarning &&
debug)
425 Pout<<
" separation " << separation_[facei]
427 <<
" differs from separation[0] " << separation_[0]
428 <<
" by more than local tolerance "
430 <<
". Assuming non-uniform separation." <<
endl;
442 Pout<<
" separation " <<
mag(separation_[0])
443 <<
" less than local tolerance " << smallDist[0]
444 <<
". Assuming zero separation." <<
endl;
447 separation_.setSize(0);
454 Pout<<
" separation " <<
mag(separation_[0])
455 <<
" more than local tolerance " << smallDist[0]
456 <<
". Assuming uniform separation." <<
endl;
459 separation_.setSize(1);
468 Pout<<
" separation_:" << separation_.size() <<
nl
469 <<
" forwardT size:" << forwardT_.size() <<
endl;
483 const word& patchType,
488 matchTolerance_(defaultMatchTol_),
499 const word& patchType
506 transformTypeNames.lookupOrDefault
510 transformType::UNKNOWN
523 matchTolerance_(pp.matchTolerance_),
524 transform_(pp.transform_)
537 polyPatch(pp, bm, index, newSize, newStart),
538 matchTolerance_(pp.matchTolerance_),
539 transform_(pp.transform_)
552 polyPatch(pp, bm, index, mapAddressing, newStart),
553 matchTolerance_(pp.matchTolerance_),
554 transform_(pp.transform_)
571 os.
writeEntry(
"matchTolerance", matchTolerance_);
572 os.
writeEntry(
"transform", transformTypeNames[transform_]);
static pointField getAnchorPoints(const UList< face > &, const pointField &, const transformType)
Get a unique anchor point for all faces.
int debug
Static debugging option.
const Cmpt & x() const
Access to the vector x component.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
A class for handling words, derived from Foam::string.
A class for handling file names.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
static const Enum< transformType > transformTypeNames
List< bool > boolList
A List of bools.
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...
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.
const Cmpt & z() const
Access to the vector z component.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
#define forAll(list, i)
Loop across all elements in list.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
void cmptMax(FieldField< Field, typename FieldField< Field, Type >::cmptType > &cf, const FieldField< Field, Type > &f)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
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.
static label getRotation(const pointField &points, const face &f, const point &anchor, const scalar tol)
Get the number of vertices face f needs to be rotated such that.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
static scalarField calcFaceTol(const UList< face > &faces, const pointField &points, const pointField &faceCentres)
Calculate typical tolerance per face. Is currently max distance.
Output to file stream, using an OSstream.
const Cmpt & y() const
Access to the vector y component.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual ~coupledPolyPatch()
Destructor.
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label ListType::const_reference const label start
const dimensionedScalar e
Elementary charge.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
void calcTransformTensors(const vectorField &Cf, const vectorField &Cr, const vectorField &nf, const vectorField &nr, const scalarField &smallDist, const scalar absTol, const transformType=UNKNOWN) const
Calculate the transformation tensors.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
void size(const label n) noexcept
Override size to be inconsistent with allocated storage.
A List with indirect addressing.
A face is a list of labels corresponding to mesh vertices.
Various functions to operate on Lists.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
const volScalarField & p0
tensor rotationTensor(const vector &n1, const vector &n2)
Rotational transformation tensor from vector n1 to n2.
static void writeOBJ(Ostream &os, const point &pt)
Write point in OBJ format.
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)
#define WarningInFunction
Report a warning using Foam::Warning.
Class to handle errors and exceptions in a simple, consistent stream-based manner.
labelList pointLabels(nPoints, -1)