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
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.
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.