Go to the documentation of this file.
40 const scalar coupledPolyPatch::defaultMatchTol_ = 1
e-4;
50 { transformType::UNKNOWN,
"unknown" },
51 { transformType::ROTATIONAL,
"rotational" },
52 { transformType::TRANSLATIONAL,
"translational" },
53 { transformType::COINCIDENTFULLMATCH,
"coincidentFullMatch" },
54 { transformType::NOORDERING,
"noOrdering" },
62 os <<
"v " << pt.
x() <<
' ' << pt.
y() <<
' ' << pt.
z() <<
endl;
94 os <<
"l " << vertI-1 <<
' ' << vertI <<
nl;
113 const face&
f = faces[i];
117 if (foamToObj.insert(
f[fp], vertI))
127 os <<
' ' << foamToObj[
f[fp]]+1;
129 os <<
' ' << foamToObj[
f[0]]+1 <<
nl;
148 anchors[facei] =
points[faces[facei][0]];
156 const face&
f = faces[facei];
166 for (label fp2 = 0; fp2 <
f.size(); ++fp2)
168 if (
f[fp1] ==
f[fp2])
193 anchors[facei] =
points[faces[facei][0]];
214 const point& cc = faceCentres[facei];
216 const face&
f = faces[facei];
220 scalar maxLenSqr = -GREAT;
223 scalar maxCmpt = -GREAT;
228 maxLenSqr =
max(maxLenSqr,
magSqr(pt - cc));
251 scalar minDistSqr = GREAT;
257 if (distSqr < minDistSqr)
259 minDistSqr = distSqr;
264 if (anchorFp == -1 ||
Foam::sqrt(minDistSqr) > tol)
275 if (distSqr == minDistSqr && fp != anchorFp)
278 <<
"Cannot determine unique anchor point on face "
281 <<
"Both at index " << anchorFp <<
" and " << fp
282 <<
" the vertices have the same distance "
284 <<
" to the anchor " << anchor
285 <<
". Continuing but results might be wrong."
291 return (
f.size() - anchorFp) %
f.size();
309 Pout<<
"coupledPolyPatch::calcTransformTensors : " <<
name() <<
endl
310 <<
" transform:" << transformTypeNames[
transform] <<
nl
311 <<
" (half)size:" << Cf.size() <<
nl
312 <<
" absTol:" << absTol <<
nl
313 <<
" smallDist min:" <<
min(smallDist) <<
nl
314 <<
" smallDist max:" <<
max(smallDist) <<
nl
315 <<
" sum(mag(nf & nr)):" <<
sum(
mag(nf & nr)) <<
endl;
329 separation_.setSize(0);
332 collocated_.setSize(0);
357 separation_.setSize(0);
359 forwardT_.setSize(Cf.size());
360 reverseT_.setSize(Cf.size());
361 collocated_.setSize(Cf.size());
372 Pout<<
" sum(mag(forwardT_ - forwardT_[0])):"
373 <<
sum(
mag(forwardT_ - forwardT_[0]))
379 forwardT_.setSize(1);
380 reverseT_.setSize(1);
381 collocated_.setSize(1);
385 Pout<<
" difference in rotation less than"
386 <<
" local tolerance "
387 <<
error <<
". Assuming uniform rotation." <<
endl;
395 forwardT_.setSize(0);
396 reverseT_.setSize(0);
398 separation_ = Cr - Cf;
400 collocated_.setSize(separation_.size());
408 bool sameSeparation =
true;
409 bool doneWarning =
false;
411 forAll(separation_, facei)
413 scalar smallSqr =
sqr(smallDist[facei]);
415 collocated_[facei] = (
magSqr(separation_[facei]) < smallSqr);
418 if (
magSqr(separation_[facei] - separation_[0]) > smallSqr)
420 sameSeparation =
false;
422 if (!doneWarning &&
debug)
426 Pout<<
" separation " << separation_[facei]
428 <<
" differs from separation[0] " << separation_[0]
429 <<
" by more than local tolerance "
431 <<
". Assuming non-uniform separation." <<
endl;
443 Pout<<
" separation " <<
mag(separation_[0])
444 <<
" less than local tolerance " << smallDist[0]
445 <<
". Assuming zero separation." <<
endl;
448 separation_.setSize(0);
455 Pout<<
" separation " <<
mag(separation_[0])
456 <<
" more than local tolerance " << smallDist[0]
457 <<
". Assuming uniform separation." <<
endl;
460 separation_.setSize(1);
469 Pout<<
" separation_:" << separation_.size() <<
nl
470 <<
" forwardT size:" << forwardT_.size() <<
endl;
484 const word& patchType,
489 matchTolerance_(defaultMatchTol_),
500 const word& patchType
507 transformTypeNames.getOrDefault
511 transformType::UNKNOWN
524 matchTolerance_(pp.matchTolerance_),
525 transform_(pp.transform_)
538 polyPatch(pp, bm, index, newSize, newStart),
539 matchTolerance_(pp.matchTolerance_),
540 transform_(pp.transform_)
553 polyPatch(pp, bm, index, mapAddressing, newStart),
554 matchTolerance_(pp.matchTolerance_),
555 transform_(pp.transform_)
572 os.
writeEntry(
"matchTolerance", matchTolerance_);
573 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...
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)
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)
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.
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
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)