Go to the documentation of this file.
55 const UList<face>& faces,
63 ctrs[facei] = faces[facei].centre(
points);
72 const UList<face>& faces,
80 anchors[facei] =
points[faces[facei][0]];
87 Foam::label Foam::oldCyclicPolyPatch::findMaxArea
94 scalar maxAreaSqr = -GREAT;
98 scalar areaSqr =
magSqr(faces[facei].areaNormal(
points));
100 if (maxAreaSqr < areaSqr)
102 maxAreaSqr = areaSqr;
110 bool Foam::oldCyclicPolyPatch::getGeometricHalves
127 boolList regionEdge(pp.nEdges(),
false);
131 label nRegionEdges = 0;
135 const labelList& eFaces = edgeFaces[edgeI];
139 if (eFaces.size() == 2)
143 regionEdge[edgeI] =
true;
153 patchZones ppZones(pp, regionEdge);
157 Pout<<
"oldCyclicPolyPatch::getGeometricHalves : "
158 <<
"Found " << nRegionEdges <<
" edges on patch " <<
name()
159 <<
" where the cos of the angle between two connected faces"
160 <<
" was less than " << featureCos_ <<
nl
161 <<
"Patch divided by these and by single sides edges into "
162 << ppZones.nZones() <<
" parts." <<
endl;
169 for (label zoneI = 0; zoneI < ppZones.nZones(); zoneI++)
176 Pout<<
"oldCyclicPolyPatch::getGeometricHalves : Writing zone "
177 << zoneI <<
" face centres to OBJ file " << stream.
name()
184 writeOBJ(stream, pp[zoneFaces[i]].centre(pp.points()));
187 nZoneFaces[zoneI] = zoneFaces.size();
192 if (ppZones.nZones() == 2)
201 Pout<<
"oldCyclicPolyPatch::getGeometricHalves :"
202 <<
" falling back to face-normal comparison" <<
endl;
205 half0ToPatch.setSize(pp.size());
208 half1ToPatch.setSize(pp.size());
215 half0ToPatch[n0Faces++] = facei;
219 half1ToPatch[n1Faces++] = facei;
222 half0ToPatch.setSize(n0Faces);
223 half1ToPatch.setSize(n1Faces);
227 Pout<<
"oldCyclicPolyPatch::getGeometricHalves :"
228 <<
" Number of faces per zone:("
229 << n0Faces <<
' ' << n1Faces <<
')' <<
endl;
233 if (half0ToPatch.size() != half1ToPatch.size())
235 fileName casePath(boundaryMesh().
mesh().time().
path());
239 fileName nm0(casePath/
name()+
"_half0_faces.obj");
240 Pout<<
"oldCyclicPolyPatch::getGeometricHalves : Writing half0"
241 <<
" faces to OBJ file " << nm0 <<
endl;
242 writeOBJ(nm0, UIndirectList<face>(pp, half0ToPatch)(), pp.points());
244 fileName nm1(casePath/
name()+
"_half1_faces.obj");
245 Pout<<
"oldCyclicPolyPatch::getGeometricHalves : Writing half1"
246 <<
" faces to OBJ file " << nm1 <<
endl;
247 writeOBJ(nm1, UIndirectList<face>(pp, half1ToPatch)(), pp.points());
252 OFstream str0(casePath/
name()+
"_half0.obj");
253 Pout<<
"oldCyclicPolyPatch::getGeometricHalves : Writing half0"
254 <<
" face centres to OBJ file " << str0.
name() <<
endl;
258 writeOBJ(str0, pp[half0ToPatch[i]].centre(pp.points()));
261 OFstream str1(casePath/
name()+
"_half1.obj");
262 Pout<<
"oldCyclicPolyPatch::getGeometricHalves : Writing half1"
263 <<
" face centres to OBJ file " << str1.
name() <<
endl;
266 writeOBJ(str1, pp[half1ToPatch[i]].centre(pp.points()));
271 <<
"Patch " <<
name() <<
" gets decomposed in two zones of"
272 <<
"inequal size: " << half0ToPatch.size()
273 <<
" and " << half1ToPatch.size() <<
endl
274 <<
"This means that the patch is either not two separate regions"
275 <<
" or one region where the angle between the different regions"
276 <<
" is not sufficiently sharp." <<
endl
277 <<
"Please adapt the featureCos setting." <<
endl
278 <<
"Continuing with incorrect face ordering from now on!" <<
endl;
287 void Foam::oldCyclicPolyPatch::getCentresAndAnchors
301 half0Ctrs = calcFaceCentres(half0Faces, pp.points());
302 anchors0 = getAnchorPoints(half0Faces, pp.points());
303 half1Ctrs = calcFaceCentres(half1Faces, pp.points());
309 label face0 = getConsistentRotationFace(half0Ctrs);
310 label face1 = getConsistentRotationFace(half1Ctrs);
315 (half0Ctrs[face0] - rotationCentre_) ^ rotationAxis_
321 (half1Ctrs[face1] - rotationCentre_) ^ -rotationAxis_
326 Pout<<
"oldCyclicPolyPatch::getCentresAndAnchors :"
327 <<
" Specified rotation :"
328 <<
" n0:" << n0 <<
" n1:" << n1 <<
endl;
369 label max0I = findMaxArea(pp.points(), half0Faces);
370 const vector n0 = half0Faces[max0I].unitNormal(pp.points());
372 label max1I = findMaxArea(pp.points(), half1Faces);
373 const vector n1 = half1Faces[max1I].unitNormal(pp.points());
375 if (
mag(n0 & n1) < 1-matchTolerance())
379 Pout<<
"oldCyclicPolyPatch::getCentresAndAnchors :"
380 <<
" Detected rotation :"
381 <<
" n0:" << n0 <<
" n1:" << n1 <<
endl;
408 const pointField& half0Pts = half0.localPoints();
409 const point ctr0(
sum(half0Pts)/half0Pts.size());
412 const pointField& half1Pts = half1.localPoints();
413 const point ctr1(
sum(half1Pts)/half1Pts.size());
417 Pout<<
"oldCyclicPolyPatch::getCentresAndAnchors :"
418 <<
" Detected translation :"
419 <<
" n0:" << n0 <<
" n1:" << n1
420 <<
" ctr0:" << ctr0 <<
" ctr1:" << ctr1 <<
endl;
423 half0Ctrs += ctr1 - ctr0;
424 anchors0 += ctr1 - ctr0;
425 ppPoints = pp.points() + ctr1 - ctr0;
433 tols = matchTolerance()*calcFaceTol(half1Faces, pp.points(), half1Ctrs);
437 bool Foam::oldCyclicPolyPatch::matchAnchors
457 forAll(half0ToPatch, half0Facei)
460 label patchFacei = half0ToPatch[half0Facei];
462 faceMap[patchFacei] = half0Facei;
465 rotation[patchFacei] = 0;
468 bool fullMatch =
true;
470 forAll(from1To0, half1Facei)
472 label patchFacei = half1ToPatch[half1Facei];
475 label half0Facei = from1To0[half1Facei];
477 label newFacei = half0Facei + pp.size()/2;
479 faceMap[patchFacei] = newFacei;
485 const point& wantedAnchor = anchors0[half0Facei];
487 rotation[newFacei] = getRotation
490 half1Faces[half1Facei],
495 if (rotation[newFacei] == -1)
501 const face&
f = half1Faces[half1Facei];
503 <<
"Patch:" <<
name() <<
" : "
504 <<
"Cannot find point on face " <<
f
506 << UIndirectList<point>(pp.points(),
f)
507 <<
" that matches point " << wantedAnchor
508 <<
" when matching the halves of cyclic patch " <<
name()
510 <<
"Continuing with incorrect face ordering from now on!"
519 Foam::label Foam::oldCyclicPolyPatch::getConsistentRotationFace
526 magSqr((faceCentres - rotationCentre_) ^ rotationAxis_)
530 (faceCentres - rotationCentre_) & rotationAxis_
532 axisLen = axisLen -
min(axisLen);
535 magRadSqr + axisLen*axisLen
539 scalar maxMagLenSqr = -GREAT;
540 scalar maxMagRadSqr = -GREAT;
543 if (magLenSqr[i] >= maxMagLenSqr)
545 if (magRadSqr[i] > maxMagRadSqr)
548 maxMagLenSqr = magLenSqr[i];
549 maxMagRadSqr = magRadSqr[i];
556 Info<<
"getConsistentRotationFace(const pointField&)" <<
nl
557 <<
" rotFace = " << rotFace <<
nl
558 <<
" point = " << faceCentres[rotFace] <<
endl;
574 const word& patchType,
581 rotationCentre_(
Zero),
582 separationVector_(
Zero)
592 const word& patchType
598 rotationCentre_(
Zero),
599 separationVector_(
Zero)
604 <<
"Found \"neighbourPatch\" entry when reading cyclic patch "
606 <<
"Is this mesh already with split cyclics?" <<
endl
607 <<
"If so run a newer version that supports it"
608 <<
", if not comment out the \"neighbourPatch\" entry and re-run"
642 featureCos_(pp.featureCos_),
643 rotationAxis_(pp.rotationAxis_),
644 rotationCentre_(pp.rotationCentre_),
645 separationVector_(pp.separationVector_)
659 featureCos_(pp.featureCos_),
660 rotationAxis_(pp.rotationAxis_),
661 rotationCentre_(pp.rotationCentre_),
662 separationVector_(pp.separationVector_)
764 <<
"Size of cyclic " <<
name() <<
" should be a multiple of 2"
768 label halfSize = pp.size()/2;
786 half1ToPatch = half0ToPatch + halfSize;
793 pointField half0Ctrs, half1Ctrs, anchors0, ppPoints;
821 Pout<<
"oldCyclicPolyPatch::order : test if already ordered:"
822 << matchedAll <<
endl;
826 Pout<<
"oldCyclicPolyPatch::order : Writing half0"
827 <<
" faces to OBJ file " << nm0 <<
endl;
828 writeOBJ(nm0, half0Faces, ppPoints);
831 Pout<<
"oldCyclicPolyPatch::order : Writing half1"
832 <<
" faces to OBJ file " << nm1 <<
endl;
838 /
"match1_"+
name() +
"_faceCentres.obj"
840 Pout<<
"oldCyclicPolyPatch::order : "
841 <<
"Dumping currently found cyclic match as lines between"
842 <<
" corresponding face centres to file " << ccStr.
name()
856 const point& c0 = half0Ctrs[i];
857 const point&
c1 = half1Ctrs[i];
870 for (label i = 0; i < halfSize; i++)
872 half0ToPatch[i] = facei++;
873 half1ToPatch[i] = facei++;
905 Pout<<
"oldCyclicPolyPatch::order : test if pairwise ordered:"
906 << matchedAll <<
endl;
909 Pout<<
"oldCyclicPolyPatch::order : Writing half0"
910 <<
" faces to OBJ file " << nm0 <<
endl;
911 writeOBJ(nm0, half0Faces, ppPoints);
914 Pout<<
"oldCyclicPolyPatch::order : Writing half1"
915 <<
" faces to OBJ file " << nm1 <<
endl;
921 /
"match2_"+
name()+
"_faceCentres.obj"
923 Pout<<
"oldCyclicPolyPatch::order : "
924 <<
"Dumping currently found cyclic match as lines between"
925 <<
" corresponding face centres to file " << ccStr.
name()
933 if (from1To0[i] != -1)
936 const point& c0 = half0Ctrs[from1To0[i]];
937 const point&
c1 = half1Ctrs[i];
957 label matchedFacei = -1;
961 label otherFacei =
pFaces[i];
963 if (otherFacei > facei)
971 matchedFacei = otherFacei;
977 if (matchedFacei != -1)
979 half0ToPatch[baffleI] = facei;
980 half1ToPatch[baffleI] = matchedFacei;
985 if (baffleI == halfSize)
1016 Pout<<
"oldCyclicPolyPatch::order : test if baffles:"
1017 << matchedAll <<
endl;
1020 Pout<<
"oldCyclicPolyPatch::order : Writing half0"
1021 <<
" faces to OBJ file " << nm0 <<
endl;
1022 writeOBJ(nm0, half0Faces, ppPoints);
1025 Pout<<
"oldCyclicPolyPatch::order : Writing half1"
1026 <<
" faces to OBJ file " << nm1 <<
endl;
1032 /
"match3_"+
name() +
"_faceCentres.obj"
1034 Pout<<
"oldCyclicPolyPatch::order : "
1035 <<
"Dumping currently found cyclic match as lines between"
1036 <<
" corresponding face centres to file " << ccStr.
name()
1044 if (from1To0[i] != -1)
1047 const point& c0 = half0Ctrs[from1To0[i]];
1048 const point&
c1 = half1Ctrs[i];
1063 bool okSplit = getGeometricHalves(pp, half0ToPatch, half1ToPatch);
1075 getCentresAndAnchors
1100 Pout<<
"oldCyclicPolyPatch::order : automatic ordering result:"
1101 << matchedAll <<
endl;
1104 Pout<<
"oldCyclicPolyPatch::order : Writing half0"
1105 <<
" faces to OBJ file " << nm0 <<
endl;
1106 writeOBJ(nm0, half0Faces, ppPoints);
1109 Pout<<
"oldCyclicPolyPatch::order : Writing half1"
1110 <<
" faces to OBJ file " << nm1 <<
endl;
1116 /
"match4_"+
name() +
"_faceCentres.obj"
1118 Pout<<
"oldCyclicPolyPatch::order : "
1119 <<
"Dumping currently found cyclic match as lines between"
1120 <<
" corresponding face centres to file " << ccStr.
name()
1128 if (from1To0[i] != -1)
1131 const point& c0 = half0Ctrs[from1To0[i]];
1132 const point&
c1 = half1Ctrs[i];
1140 if (!matchedAll ||
debug)
1144 Pout<<
"oldCyclicPolyPatch::order : Writing half0"
1145 <<
" faces to OBJ file " << nm0 <<
endl;
1149 Pout<<
"oldCyclicPolyPatch::order : Writing half1"
1150 <<
" faces to OBJ file " << nm1 <<
endl;
1156 /
name() +
"_faceCentres.obj"
1158 Pout<<
"oldCyclicPolyPatch::order : "
1159 <<
"Dumping currently found cyclic match as lines between"
1160 <<
" corresponding face centres to file " << ccStr.
name()
1169 if (from1To0[i] != -1)
1173 const point& c0 = half0Ctrs[from1To0[i]];
1174 const point&
c1 = half1Ctrs[i];
1184 <<
"Patch:" <<
name() <<
" : "
1185 <<
"Cannot match vectors to faces on both sides of patch" <<
endl
1186 <<
" Perhaps your faces do not match?"
1187 <<
" The obj files written contain the current match." <<
endl
1188 <<
" Continuing with incorrect face ordering from now on!"
1215 if (
faceMap[facei] != facei || rotation[facei] != 0)
1228 os.
writeEntry(
"type", cyclicPolyPatch::typeName);
1239 os.
writeEntry(
"rotationAxis", rotationAxis_);
1240 os.
writeEntry(
"rotationCentre", rotationCentre_);
1245 os.
writeEntry(
"separationVector", separationVector_);
1255 <<
"Please run foamUpgradeCyclics to convert these old-style"
1256 <<
" cyclics into two separate cyclics patches."
int debug
Static debugging option.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const labelListList & pointFaces() const
Return point-face addressing.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
A class for handling words, derived from Foam::string.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
A class for handling file names.
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.
PrimitivePatch< List< face >, const pointField & > primitiveFacePatch
A PrimitivePatch with List storage for the faces, const reference for the point field.
virtual void movePoints(PstreamBuffers &, const pointField &)
Correct patches after moving points.
static constexpr const zero Zero
Global zero (0)
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
Template functions to aid in the implementation of demand driven data.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
List< bool > boolList
A List of bools.
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
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
OSstream wrapped stdout (std::cout) with parallel prefix.
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.
#define forAll(list, i)
Loop across all elements in list.
virtual ~oldCyclicPolyPatch()
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
Determine correspondence between points. See below.
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=cellModel::ref(cellModel::HEX);labelList hexLabels(8);label activeCells=0;labelList pointMap(nPoints);forAll(pointMap, i){ pointMap[i]=i;}for(label i=0;i< nPoints;i++){ if(f[i] > 0.0) { hexLabels[0]=i;hexLabels[1]=i1tab[i];hexLabels[2]=i3tab[i1tab[i]];hexLabels[3]=i3tab[i];hexLabels[4]=i8tab[i];hexLabels[5]=i1tab[i8tab[i]];hexLabels[6]=i3tab[i1tab[i8tab[i]]];hexLabels[7]=i3tab[i8tab[i]];cellShapes[activeCells]=cellShape(hex, hexLabels);edgeList edges=cellShapes[activeCells].edges();forAll(edges, ei) { if(edges[ei].mag(points)< SMALL) { label start=pointMap[edges[ei].start()];while(start !=pointMap[start]) { start=pointMap[start];} label end=pointMap[edges[ei].end()];while(end !=pointMap[end]) { end=pointMap[end];} label minLabel=min(start, end);pointMap[start]=pointMap[end]=minLabel;} } cellZoning[activeCells]=idreg[i];activeCells++;}}cellShapes.setSize(activeCells);cellZoning.setSize(activeCells);forAll(cellShapes, celli){ cellShape &cs=cellShapes[celli];forAll(cs, i) { cs[i]=pointMap[cs[i]];} cs.collapse();}label bcIDs[11]={-1, 0, 2, 4, -1, 5, -1, 6, 7, 8, 9};const label nBCs=12;const word *kivaPatchTypes[nBCs]={ &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &wallPolyPatch::typeName, &symmetryPolyPatch::typeName, &wedgePolyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &polyPatch::typeName, &symmetryPolyPatch::typeName, &oldCyclicPolyPatch::typeName};enum patchTypeNames{ PISTON, VALVE, LINER, CYLINDERHEAD, AXIS, WEDGE, INFLOW, OUTFLOW, PRESIN, PRESOUT, SYMMETRYPLANE, CYCLIC};const char *kivaPatchNames[nBCs]={ "piston", "valve", "liner", "cylinderHead", "axis", "wedge", "inflow", "outflow", "presin", "presout", "symmetryPlane", "cyclic"};List< SLList< face > > pFaces[nBCs]
surfaceVectorField faceNormals(mesh.Sf()/mesh.magSf())
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
messageStream Info
Information stream (uses stdout - output is on the master only)
word name(const complex &c)
Return string representation of complex.
virtual const fileName & name() const
Return the name of the stream.
#define SeriousErrorInFunction
Report an error message using Foam::SeriousError.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const Field< point_type > & points() const
Return reference to global points.
Macros for easy insertion into run-time selection tables.
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
const List< face_type > & localFaces() const
Return patch faces addressing into local point list.
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
errorManip< error > abort(error &err)
Vector< scalar > vector
A scalar version of the templated Vector.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Output to file stream, using an OSstream.
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
List< labelList > labelListList
A List of labelList.
'old' style cyclic polyPatch with all faces in single patch. Does ordering but cannot be used to run....
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
List< face > faceList
A List of faces.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
A List with indirect addressing.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
void write(Ostream &os) const
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field.
labelList findIndices(const ListType &input, typename ListType::const_reference val, label start=0)
Linear search to find all occurrences of given element.
A face is a list of labels corresponding to mesh vertices.
#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,...
tensor rotationTensor(const vector &n1, const vector &n2)
Rotational transformation tensor from vector n1 to n2.
vector point
Point is a vector.
void setSize(const label newSize)
Alias for resize(const label)
defineTypeNameAndDebug(combustionModel, 0)
bool matchPoints(const UList< point > &pts0, const UList< point > &pts1, const UList< scalar > &matchDistance, const bool verbose, labelList &from0To1, const point &origin=point::zero)
Determine correspondence between pointFields. Gets passed.
#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.
oldCyclicPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType, const transformType transform=UNKNOWN)
Construct from components.