cellLooper.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011 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  (
55  "cellLooper::New(const word&, const polyMesh&)"
56  ) << "Unknown set type "
57  << type << nl << nl
58  << "Valid cellLooper types : " << endl
59  << wordConstructorTablePtr_->sortedToc()
60  << exit(FatalError);
61  }
62 
63  return autoPtr<cellLooper>(cstrIter()(mesh));
64 }
65 
66 
67 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
68 
69 // Get faces (on cell) connected to vertI which are not using edgeI
71 (
72  const label cellI,
73  const label edgeI,
74  const label vertI
75 ) const
76 {
77  // Get faces connected to startEdge
78  label face0, face1;
79  meshTools::getEdgeFaces(mesh(), cellI, edgeI, face0, face1);
80 
81  const labelList& pFaces = mesh().pointFaces()[vertI];
82 
83  labelList vertFaces(pFaces.size());
84  label vertFaceI = 0;
85 
86  forAll(pFaces, pFaceI)
87  {
88  label faceI = pFaces[pFaceI];
89 
90  if
91  (
92  (faceI != face0)
93  && (faceI != face1)
94  && (meshTools::faceOnCell(mesh(), cellI, faceI))
95  )
96  {
97  vertFaces[vertFaceI++] = faceI;
98  }
99  }
100  vertFaces.setSize(vertFaceI);
101 
102  return vertFaces;
103 }
104 
105 
106 // Get first edge connected to vertI and on faceI
108 (
109  const label faceI,
110  const label vertI
111 ) const
112 {
113  const labelList& fEdges = mesh().faceEdges()[faceI];
114 
115  forAll(fEdges, fEdgeI)
116  {
117  label edgeI = fEdges[fEdgeI];
118 
119  const edge& e = mesh().edges()[edgeI];
120 
121  if ((e.start() == vertI) || (e.end() == vertI))
122  {
123  return edgeI;
124  }
125  }
126 
128  (
129  "getFirstVertEdge(const label, const label)"
130  ) << "Can not find edge on face " << faceI
131  << " using vertex " << vertI
132  << abort(FatalError);
133 
134  return -1;
135 }
136 
137 
138 // Get edges (on cell) connected to vertI which are not on faceI
140 (
141  const label cellI,
142  const label faceI,
143  const label vertI
144 ) const
145 {
146  const labelList& exclEdges = mesh().faceEdges()[faceI];
147 
148  const labelList& pEdges = mesh().pointEdges()[vertI];
149 
150  labelList vertEdges(pEdges.size());
151  label vertEdgeI = 0;
152 
153  forAll(pEdges, pEdgeI)
154  {
155  label edgeI = pEdges[pEdgeI];
156 
157  if
158  (
159  (findIndex(exclEdges, edgeI) == -1)
160  && meshTools::edgeOnCell(mesh(), cellI, edgeI)
161  )
162  {
163  vertEdges[vertEdgeI++] = edgeI;
164  }
165  }
166 
167  vertEdges.setSize(vertEdgeI);
168 
169  return vertEdges;
170 }
171 
172 
173 // Return edge from cellEdges that is most perpendicular
174 // to refinement direction.
176 (
177  const vector& refDir,
178  const label cellI
179 ) const
180 {
181  const labelList& cEdges = mesh().cellEdges()[cellI];
182 
183  label cutEdgeI = -1;
184  scalar maxCos = -GREAT;
185 
186  forAll(cEdges, cEdgeI)
187  {
188  label edgeI = cEdges[cEdgeI];
189 
190  scalar cosAngle = mag(refDir & meshTools::normEdgeVec(mesh(), edgeI));
191 
192  if (cosAngle > maxCos)
193  {
194  maxCos = cosAngle;
195 
196  cutEdgeI = edgeI;
197  }
198  }
199 
200  return cutEdgeI;
201 }
202 
203 
204 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
205 
206 Foam::cellLooper::cellLooper(const polyMesh& mesh)
207 :
208  edgeVertex(mesh)
209 {}
210 
211 
212 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
213 
215 {}
216 
217 
218 // ************************************************************************* //
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
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:140
vector normEdgeVec(const primitiveMesh &, const label edgeI)
Normalized edge vector.
Definition: meshTools.C:195
dimensioned< scalar > mag(const dimensioned< Type > &)
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
const labelListList & pointEdges() const
virtual ~cellLooper()
Destructor.
Definition: cellLooper.C:214
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:470
A class for handling words, derived from string.
Definition: word.H:59
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Various functions to operate on Lists.
static autoPtr< cellLooper > New(const word &type, const polyMesh &mesh)
Return a reference to the selected cellLooper.
Definition: cellLooper.C:43
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:71
dynamicFvMesh & mesh
bool faceOnCell(const primitiveMesh &, const label cellI, const label faceI)
Is face used by cell.
Definition: meshTools.C:314
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
const labelListList & faceEdges() const
const labelListList & cellEdges() const
Namespace for OpenFOAM.
bool edgeOnCell(const primitiveMesh &, const label cellI, const label edgeI)
Is edge used by cell.
Definition: meshTools.C:291
static const char nl
Definition: Ostream.H:260
const double e
Elementary charge.
Definition: doubleFloat.H:78
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define forAll(list, i)
Definition: UList.H:421
label getFirstVertEdge(const label faceI, const label vertI) const
Get first edge connected to vertI and on faceI.
Definition: cellLooper.C:108
label getMisAlignedEdge(const vector &refDir, const label cellI) const
Return edge from cellEdges that is most perpendicular.
Definition: cellLooper.C:176
const labelListList & pointFaces() const
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
error FatalError
label end() const
Return end vertex label.
Definition: edgeI.H:92
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
label start() const
Return start vertex label.
Definition: edgeI.H:81
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,&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]
Definition: readKivaGrid.H:235
Combines edge or vertex in single label. Used to specify cuts across cell circumference.
Definition: edgeVertex.H:52
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
defineTypeNameAndDebug(combustionModel, 0)