cellLooper.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "cellLooper.H"
27 #include "polyMesh.H"
28 #include "ListOps.H"
29 #include "meshTools.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  defineTypeNameAndDebug(cellLooper, 0);
36  defineRunTimeSelectionTable(cellLooper, word);
37 }
38 
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
43 (
44  const word& type,
45  const polyMesh& mesh
46 )
47 {
48  wordConstructorTable::iterator cstrIter =
49  wordConstructorTablePtr_->find(type);
50 
51  if (cstrIter == wordConstructorTablePtr_->end())
52  {
54  << "Unknown set type "
55  << type << nl << nl
56  << "Valid cellLooper types : " << endl
57  << wordConstructorTablePtr_->sortedToc()
58  << exit(FatalError);
59  }
60 
61  return autoPtr<cellLooper>(cstrIter()(mesh));
62 }
63 
64 
65 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
66 
68 (
69  const label celli,
70  const label edgeI,
71  const label vertI
72 ) const
73 {
74  // Get faces connected to startEdge
75  label face0, face1;
76  meshTools::getEdgeFaces(mesh(), celli, edgeI, face0, face1);
77 
78  const labelList& pFaces = mesh().pointFaces()[vertI];
79 
80  labelList vertFaces(pFaces.size());
81  label vertFacei = 0;
82 
83  forAll(pFaces, pFacei)
84  {
85  label facei = pFaces[pFacei];
86 
87  if
88  (
89  (facei != face0)
90  && (facei != face1)
91  && (meshTools::faceOnCell(mesh(), celli, facei))
92  )
93  {
94  vertFaces[vertFacei++] = facei;
95  }
96  }
97  vertFaces.setSize(vertFacei);
98 
99  return vertFaces;
100 }
101 
102 
104 (
105  const label facei,
106  const label vertI
107 ) const
108 {
109  const labelList& fEdges = mesh().faceEdges()[facei];
110 
111  forAll(fEdges, fEdgeI)
112  {
113  label edgeI = fEdges[fEdgeI];
114 
115  const edge& e = mesh().edges()[edgeI];
116 
117  if ((e.start() == vertI) || (e.end() == vertI))
118  {
119  return edgeI;
120  }
121  }
122 
124  << "Can not find edge on face " << facei
125  << " using vertex " << vertI
126  << abort(FatalError);
127 
128  return -1;
129 }
130 
131 
133 (
134  const label celli,
135  const label facei,
136  const label vertI
137 ) const
138 {
139  const labelList& exclEdges = mesh().faceEdges()[facei];
140 
141  const labelList& pEdges = mesh().pointEdges()[vertI];
142 
143  labelList vertEdges(pEdges.size());
144  label vertEdgeI = 0;
145 
146  forAll(pEdges, pEdgeI)
147  {
148  label edgeI = pEdges[pEdgeI];
149 
150  if
151  (
152  (findIndex(exclEdges, edgeI) == -1)
153  && meshTools::edgeOnCell(mesh(), celli, edgeI)
154  )
155  {
156  vertEdges[vertEdgeI++] = edgeI;
157  }
158  }
159 
160  vertEdges.setSize(vertEdgeI);
161 
162  return vertEdges;
163 }
164 
165 
167 (
168  const vector& refDir,
169  const label celli
170 ) const
171 {
172  const labelList& cEdges = mesh().cellEdges()[celli];
173 
174  label cutEdgeI = -1;
175  scalar maxCos = -great;
176 
177  forAll(cEdges, cEdgeI)
178  {
179  label edgeI = cEdges[cEdgeI];
180 
181  scalar cosAngle = mag(refDir & meshTools::normEdgeVec(mesh(), edgeI));
182 
183  if (cosAngle > maxCos)
184  {
185  maxCos = cosAngle;
186 
187  cutEdgeI = edgeI;
188  }
189  }
190 
191  return cutEdgeI;
192 }
193 
194 
195 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
196 
198 :
199  edgeVertex(mesh)
200 {}
201 
202 
203 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
204 
206 {}
207 
208 
209 // ************************************************************************* //
const labelListList & cellEdges() const
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
vector normEdgeVec(const primitiveMesh &, const label edgeI)
Normalised edge vector.
Definition: meshTools.C:189
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
const labelListList & faceEdges() const
label getFirstVertEdge(const label facei, const label vertI) const
Get first edge connected to vertI and on facei.
Definition: cellLooper.C:104
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
bool edgeOnCell(const primitiveMesh &, const label celli, const label edgeI)
Is edge used by cell.
Definition: meshTools.C:285
const labelListList & pointEdges() const
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
fvMesh & mesh
labelList getVertFacesNonEdge(const label celli, const label edgeI, const label vertI) const
Get faces (on cell) connected to vertI which are not using edgeI.
Definition: cellLooper.C:68
bool faceOnCell(const primitiveMesh &, const label celli, const label facei)
Is face used by cell.
Definition: meshTools.C:308
Various functions to operate on Lists.
label getMisAlignedEdge(const vector &refDir, const label celli) const
Return edge from cellEdges that is most perpendicular.
Definition: cellLooper.C:167
cellLooper(const polyMesh &mesh)
Construct from components.
Definition: cellLooper.C:197
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
A class for handling words, derived from string.
Definition: word.H:59
Combines edge or vertex in single label. Used to specify cuts across cell circumference.
Definition: edgeVertex.H:52
labelList getVertEdgesNonFace(const label celli, const label facei, const label vertI) const
Get edges (on cell) connected to vertI which are not on facei.
Definition: cellLooper.C:133
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Info<< "Finished reading KIVA file"<< endl;cellShapeList cellShapes(nPoints);labelList cellZoning(nPoints, -1);const cellModel &hex=*(cellModeller::lookup("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, &mergedCyclicPolyPatch::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]
Definition: readKivaGrid.H:235
static const char nl
Definition: Ostream.H:260
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
defineTypeNameAndDebug(combustionModel, 0)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
virtual ~cellLooper()
Destructor.
Definition: cellLooper.C:205
label end() const
Return end vertex label.
Definition: edgeI.H:92
const labelListList & pointFaces() const
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
static autoPtr< cellLooper > New(const word &type, const polyMesh &mesh)
Return a reference to the selected cellLooper.
Definition: cellLooper.C:43
void getEdgeFaces(const primitiveMesh &, const label celli, const label edgeI, label &face0, label &face1)
Get faces on cell using edgeI. Throws error if no two found.
Definition: meshTools.C:458
label start() const
Return start vertex label.
Definition: edgeI.H:81
Namespace for OpenFOAM.