PatchToolsNormals.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-2014 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 "PatchTools.H"
27 #include "polyMesh.H"
28 #include "indirectPrimitivePatch.H"
29 #include "globalMeshData.H"
30 
31 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32 
33 template
34 <
35  class Face,
36  template<class> class FaceList,
37  class PointField,
38  class PointType
39 >
40 
43 (
44  const polyMesh& mesh,
46 )
47 {
48  const globalMeshData& globalData = mesh.globalData();
49  const indirectPrimitivePatch& coupledPatch = globalData.coupledPatch();
50  const Map<label>& coupledPatchMP = coupledPatch.meshPointMap();
51  const mapDistribute& map = globalData.globalPointSlavesMap();
52  const globalIndexAndTransform& transforms =
53  globalData.globalTransforms();
54 
55 
56 
57 
58  // Combine normals. Note: do on all master points. Cannot just use
59  // patch points since the master point does not have to be on the
60  // patch!
61 
62  pointField coupledPointNormals(map.constructSize(), vector::zero);
63 
64  {
65  // Collect local pointFaces (sized on patch points only)
66  List<List<point> > pointFaceNormals(map.constructSize());
67  forAll(p.meshPoints(), patchPointI)
68  {
69  label meshPointI = p.meshPoints()[patchPointI];
70  Map<label>::const_iterator fnd = coupledPatchMP.find(meshPointI);
71  if (fnd != coupledPatchMP.end())
72  {
73  label coupledPointI = fnd();
74 
75  List<point>& pNormals = pointFaceNormals[coupledPointI];
76  const labelList& pFaces = p.pointFaces()[patchPointI];
77  pNormals.setSize(pFaces.size());
78  forAll(pFaces, i)
79  {
80  pNormals[i] = p.faceNormals()[pFaces[i]];
81  }
82  }
83  }
84 
85 
86  // Pull remote data into local slots
87  map.distribute
88  (
89  transforms,
90  pointFaceNormals,
92  );
93 
94 
95  // Combine all face normals (-local, -remote,untransformed,
96  // -remote,transformed)
97 
98  const labelListList& slaves = globalData.globalPointSlaves();
99  const labelListList& transformedSlaves =
100  globalData.globalPointTransformedSlaves();
101 
102  forAll(slaves, coupledPointI)
103  {
104  const labelList& slaveSlots = slaves[coupledPointI];
105  const labelList& transformedSlaveSlots =
106  transformedSlaves[coupledPointI];
107 
108  point& n = coupledPointNormals[coupledPointI];
109 
110  // Local entries
111  const List<point>& local = pointFaceNormals[coupledPointI];
112 
113  label nFaces =
114  local.size()
115  + slaveSlots.size()
116  + transformedSlaveSlots.size();
117 
118  n = sum(local);
119 
120  // Add any remote face normals
121  forAll(slaveSlots, i)
122  {
123  n += sum(pointFaceNormals[slaveSlots[i]]);
124  }
125  forAll(transformedSlaveSlots, i)
126  {
127  n += sum(pointFaceNormals[transformedSlaveSlots[i]]);
128  }
129 
130  if (nFaces >= 1)
131  {
132  n /= mag(n)+VSMALL;
133  }
134 
135  // Put back into slave slots
136  forAll(slaveSlots, i)
137  {
138  coupledPointNormals[slaveSlots[i]] = n;
139  }
140  forAll(transformedSlaveSlots, i)
141  {
142  coupledPointNormals[transformedSlaveSlots[i]] = n;
143  }
144  }
145 
146 
147  // Send back
149  (
150  transforms,
151  coupledPointNormals.size(),
152  coupledPointNormals,
154  );
155  }
156 
157 
158  // 1. Start off with local normals (note:without calculating pointNormals
159  // to avoid them being stored)
160 
161  tmp<pointField> textrudeN(new pointField(p.nPoints(), vector::zero));
162  pointField& extrudeN = textrudeN();
163  {
164  const faceList& localFaces = p.localFaces();
165  const vectorField& faceNormals = p.faceNormals();
166 
167  forAll(localFaces, faceI)
168  {
169  const face& f = localFaces[faceI];
170  const vector& n = faceNormals[faceI];
171  forAll(f, fp)
172  {
173  extrudeN[f[fp]] += n;
174  }
175  }
176  extrudeN /= mag(extrudeN)+VSMALL;
177  }
178 
179 
180  // 2. Override patch normals on coupled points
181  forAll(p.meshPoints(), patchPointI)
182  {
183  label meshPointI = p.meshPoints()[patchPointI];
184  Map<label>::const_iterator fnd = coupledPatchMP.find(meshPointI);
185  if (fnd != coupledPatchMP.end())
186  {
187  label coupledPointI = fnd();
188  extrudeN[patchPointI] = coupledPointNormals[coupledPointI];
189  }
190  }
191 
192  return textrudeN;
193 }
194 
195 
196 template
197 <
198  class Face,
199  template<class> class FaceList,
200  class PointField,
201  class PointType
202 >
203 
206 (
207  const polyMesh& mesh,
209  const labelList& patchEdges,
210  const labelList& coupledEdges
211 )
212 {
213  // 1. Start off with local normals
214 
215  tmp<pointField> tedgeNormals(new pointField(p.nEdges(), vector::zero));
216  pointField& edgeNormals = tedgeNormals();
217  {
218  const labelListList& edgeFaces = p.edgeFaces();
219  const vectorField& faceNormals = p.faceNormals();
220 
221  forAll(edgeFaces, edgeI)
222  {
223  const labelList& eFaces = edgeFaces[edgeI];
224  forAll(eFaces, i)
225  {
226  edgeNormals[edgeI] += faceNormals[eFaces[i]];
227  }
228  }
229  edgeNormals /= mag(edgeNormals)+VSMALL;
230  }
231 
232 
233 
234  const globalMeshData& globalData = mesh.globalData();
235  const mapDistribute& map = globalData.globalEdgeSlavesMap();
236 
237 
238  // Convert patch-edge data into cpp-edge data
239  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
240 
241  //- Construct with all data in consistent orientation
242  pointField cppEdgeData(map.constructSize(), vector::zero);
243 
244  forAll(patchEdges, i)
245  {
246  label patchEdgeI = patchEdges[i];
247  label coupledEdgeI = coupledEdges[i];
248  cppEdgeData[coupledEdgeI] = edgeNormals[patchEdgeI];
249  }
250 
251 
252  // Synchronise
253  // ~~~~~~~~~~~
254 
255  globalData.syncData
256  (
257  cppEdgeData,
258  globalData.globalEdgeSlaves(),
259  globalData.globalEdgeTransformedSlaves(),
260  map,
261  globalData.globalTransforms(),
262  plusEqOp<point>(), // add since normalised later on
264  );
265  cppEdgeData /= mag(cppEdgeData)+VSMALL;
266 
267 
268  // Back from cpp-edge to patch-edge data
269  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
270 
271  forAll(patchEdges, i)
272  {
273  label patchEdgeI = patchEdges[i];
274  label coupledEdgeI = coupledEdges[i];
275  edgeNormals[patchEdgeI] = cppEdgeData[coupledEdgeI];
276  }
277 
278  return tedgeNormals;
279 }
280 
281 
282 // ************************************************************************* //
Foam::pointField
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Foam::globalMeshData::globalPointTransformedSlaves
const labelListList & globalPointTransformedSlaves() const
Definition: globalMeshData.C:2180
Foam::globalMeshData::globalPointSlaves
const labelListList & globalPointSlaves() const
Definition: globalMeshData.C:2170
p
p
Definition: pEqn.H:62
Foam::globalMeshData::globalEdgeSlavesMap
const mapDistribute & globalEdgeSlavesMap() const
Definition: globalMeshData.C:2245
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:406
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:118
globalMeshData.H
Foam::Map< label >
Foam::globalMeshData::globalPointSlavesMap
const mapDistribute & globalPointSlavesMap() const
Definition: globalMeshData.C:2191
Foam::mag
dimensioned< scalar > mag(const dimensioned< Type > &)
polyMesh.H
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
pFaces
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
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::label
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
Foam::Field
Pre-declare SubField and related Field type.
Definition: Field.H:57
Foam::globalMeshData::globalTransforms
const globalIndexAndTransform & globalTransforms() const
Global transforms numbering.
Definition: globalMeshData.C:2160
Foam::mapDistribute
Class containing processor-to-processor mapping information.
Definition: mapDistribute.H:152
Foam::globalMeshData::globalEdgeSlaves
const labelListList & globalEdgeSlaves() const
Definition: globalMeshData.C:2214
Foam::mapDistribute::distribute
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
Definition: mapDistributeTemplates.C:155
PatchTools.H
Foam::globalMeshData::syncData
static void syncData(List< Type > &pointData, const labelListList &slaves, const labelListList &transformedSlaves, const mapDistribute &slavesMap, const globalIndexAndTransform &, const CombineOp &cop, const TransformOp &top)
Helper: synchronise data with transforms.
Definition: globalMeshDataTemplates.C:34
Foam::globalMeshData
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
Definition: globalMeshData.H:106
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:18
Foam::plusEqOp
Definition: ops.H:71
Foam::mapDistribute::transform
Default transformation behaviour.
Definition: mapDistribute.H:196
Foam::mapDistributeBase::constructSize
label constructSize() const
Constructed data size.
Definition: mapDistributeBase.H:244
Foam::List::setSize
void setSize(const label)
Reset size of List.
Foam::globalMeshData::globalEdgeTransformedSlaves
const labelListList & globalEdgeTransformedSlaves() const
Definition: globalMeshData.C:2224
f
labelList f(nPoints)
Foam::PatchTools::pointNormals
static tmp< pointField > pointNormals(const polyMesh &, const PrimitivePatch< Face, FaceList, PointField, PointType > &)
Return parallel consistent point normals for patches using mesh points.
Foam::Vector< scalar >
Foam::List
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
Foam::globalMeshData::coupledPatch
const indirectPrimitivePatch & coupledPatch() const
Return patch of all coupled faces.
Definition: globalMeshData.C:2046
Foam::sum
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
Definition: DimensionedFieldFunctions.C:333
Foam::mapDistribute::reverseDistribute
void reverseDistribute(const label constructSize, List< T > &, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Reverse distribute data using default commsType.
Definition: mapDistributeTemplates.C:187
Foam::face
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
Foam::List::size
void size(const label)
Override size to be inconsistent with allocated storage.
Foam::PrimitivePatch::meshPointMap
const Map< label > & meshPointMap() const
Mesh point map. Given the global point index find its.
Definition: PrimitivePatchTemplate.C:412
Foam::PatchTools::edgeNormals
static tmp< pointField > edgeNormals(const polyMesh &, const PrimitivePatch< Face, FaceList, PointField, PointType > &, const labelList &patchEdges, const labelList &coupledEdges)
Return parallel consistent edge normals for patches using mesh points.
Foam::globalIndexAndTransform
Determination and storage of the possible independent transforms introduced by coupledPolyPatches,...
Definition: globalIndexAndTransform.H:60
Foam::PrimitivePatch
A list of faces which address into the list of points.
Definition: PrimitivePatchTemplate.H:88