meshRefinementTemplates.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 "meshRefinement.H"
27 #include "fvMesh.H"
28 #include "globalIndex.H"
29 #include "syncTools.H"
30 
31 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32 
33 // Add a T entry
34 template<class T> void Foam::meshRefinement::updateList
35 (
36  const labelList& newToOld,
37  const T& nullValue,
38  List<T>& elems
39 )
40 {
41  List<T> newElems(newToOld.size(), nullValue);
42 
43  forAll(newElems, i)
44  {
45  label oldI = newToOld[i];
46 
47  if (oldI >= 0)
48  {
49  newElems[i] = elems[oldI];
50  }
51  }
52 
53  elems.transfer(newElems);
54 }
55 
56 
57 template<class T>
59 (
60  const PackedBoolList& isMasterElem,
61  const UList<T>& values
62 )
63 {
64  if (values.size() != isMasterElem.size())
65  {
67  (
68  "meshRefinement::gAverage\n"
69  "(\n"
70  " const PackedBoolList& isMasterElem,\n"
71  " const UList<T>& values\n"
72  ")\n"
73  ) << "Number of elements in list " << values.size()
74  << " does not correspond to number of elements in isMasterElem "
75  << isMasterElem.size()
76  << exit(FatalError);
77  }
78 
79  T sum = pTraits<T>::zero;
80  label n = 0;
81 
82  forAll(values, i)
83  {
84  if (isMasterElem[i])
85  {
86  sum += values[i];
87  n++;
88  }
89  }
90 
91  reduce(sum, sumOp<T>());
92  reduce(n, sumOp<label>());
93 
94  if (n > 0)
95  {
96  return sum/n;
97  }
98  else
99  {
100  return pTraits<T>::max;
101  }
102 }
103 
104 
105 // Compare two lists over all boundary faces
106 template<class T>
108 (
109  const scalar tol,
110  const string& msg,
111  const UList<T>& faceData,
112  const UList<T>& syncedFaceData
113 ) const
114 {
115  label nBFaces = mesh_.nFaces() - mesh_.nInternalFaces();
116 
117  if (faceData.size() != nBFaces || syncedFaceData.size() != nBFaces)
118  {
120  (
121  "meshRefinement::testSyncBoundaryFaceList"
122  "(const scalar, const string&, const List<T>&, const List<T>&)"
123  ) << "Boundary faces:" << nBFaces
124  << " faceData:" << faceData.size()
125  << " syncedFaceData:" << syncedFaceData.size()
126  << abort(FatalError);
127  }
128 
129  const polyBoundaryMesh& patches = mesh_.boundaryMesh();
130 
131  forAll(patches, patchI)
132  {
133  const polyPatch& pp = patches[patchI];
134 
135  label bFaceI = pp.start() - mesh_.nInternalFaces();
136 
137  forAll(pp, i)
138  {
139  const T& data = faceData[bFaceI];
140  const T& syncData = syncedFaceData[bFaceI];
141 
142  if (mag(data - syncData) > tol)
143  {
144  label faceI = pp.start()+i;
145 
146  FatalErrorIn("testSyncFaces")
147  << msg
148  << "patchFace:" << i
149  << " face:" << faceI
150  << " fc:" << mesh_.faceCentres()[faceI]
151  << " patch:" << pp.name()
152  << " faceData:" << data
153  << " syncedFaceData:" << syncData
154  << " diff:" << mag(data - syncData)
155  << abort(FatalError);
156  }
157 
158  bFaceI++;
159  }
160  }
161 }
162 
163 
164 // Print list sorted by coordinates. Used for comparing non-parallel v.s.
165 // parallel operation
166 template<class T>
168 (
169  const UList<point>& points,
170  const UList<T>& data
171 )
172 {
173  globalIndex globalPoints(points.size());
174 
176  globalPoints.gather
177  (
178  Pstream::worldComm,
179  identity(Pstream::nProcs()),
180  points,
181  allPoints,
182  UPstream::msgType(),
183  Pstream::blocking
184  );
185 
186  List<T> allData;
187  globalPoints.gather
188  (
189  Pstream::worldComm,
190  identity(Pstream::nProcs()),
191  data,
192  allData,
193  UPstream::msgType(),
194  Pstream::blocking
195  );
196 
197 
198  scalarField magAllPoints(mag(allPoints-point(-0.317, 0.117, 0.501)));
199 
200  labelList visitOrder;
201  sortedOrder(magAllPoints, visitOrder);
202  forAll(visitOrder, i)
203  {
204  label allPointI = visitOrder[i];
205  Info<< allPoints[allPointI] << " : " << allData[allPointI]
206  << endl;
207  }
208 }
209 
210 
211 //template<class T, class Mesh>
212 template<class GeoField>
213 void Foam::meshRefinement::addPatchFields
214 (
215  fvMesh& mesh,
216  const word& patchFieldType
217 )
218 {
220  (
221  mesh.objectRegistry::lookupClass<GeoField>()
222  );
223 
224  forAllIter(typename HashTable<GeoField*>, flds, iter)
225  {
226  GeoField& fld = *iter();
227  typename GeoField::GeometricBoundaryField& bfld = fld.boundaryField();
228 
229  label sz = bfld.size();
230  bfld.setSize(sz+1);
231  bfld.set
232  (
233  sz,
235  (
236  patchFieldType,
237  mesh.boundary()[sz],
238  fld.dimensionedInternalField()
239  )
240  );
241  }
242 }
243 
244 
245 // Reorder patch field
246 template<class GeoField>
247 void Foam::meshRefinement::reorderPatchFields
248 (
249  fvMesh& mesh,
250  const labelList& oldToNew
251 )
252 {
254  (
255  mesh.objectRegistry::lookupClass<GeoField>()
256  );
257 
258  forAllIter(typename HashTable<GeoField*>, flds, iter)
259  {
260  GeoField& fld = *iter();
261  typename GeoField::GeometricBoundaryField& bfld = fld.boundaryField();
262 
263  bfld.reorder(oldToNew);
264  }
265 }
266 
267 
268 template<class Enum>
270 (
271  const Enum& namedEnum,
272  const wordList& words
273 )
274 {
275  int flags = 0;
276 
277  forAll(words, i)
278  {
279  int index = namedEnum[words[i]];
280  int val = 1<<index;
281  flags |= val;
282  }
283  return flags;
284 }
285 
286 
287 template<class Type>
289 (
290  const polyMesh& mesh,
291  const PackedBoolList& isMasterEdge,
292  const labelList& meshPoints,
293  const edgeList& edges,
294  const scalarField& edgeWeights,
295  const Field<Type>& pointData,
296  Field<Type>& sum
297 )
298 {
299  if
300  (
301  edges.size() != isMasterEdge.size()
302  || edges.size() != edgeWeights.size()
303  || meshPoints.size() != pointData.size()
304  )
305  {
306  FatalErrorIn("medialAxisMeshMover::weightedSum(..)")
307  << "Inconsistent sizes for edge or point data:"
308  << " isMasterEdge:" << isMasterEdge.size()
309  << " edgeWeights:" << edgeWeights.size()
310  << " edges:" << edges.size()
311  << " pointData:" << pointData.size()
312  << " meshPoints:" << meshPoints.size()
313  << abort(FatalError);
314  }
315 
316  sum.setSize(meshPoints.size());
317  sum = pTraits<Type>::zero;
318 
319  forAll(edges, edgeI)
320  {
321  if (isMasterEdge[edgeI])
322  {
323  const edge& e = edges[edgeI];
324 
325  scalar eWeight = edgeWeights[edgeI];
326 
327  label v0 = e[0];
328  label v1 = e[1];
329 
330  sum[v0] += eWeight*pointData[v1];
331  sum[v1] += eWeight*pointData[v0];
332  }
333  }
334 
336  (
337  mesh,
338  meshPoints,
339  sum,
340  plusEqOp<Type>(),
341  pTraits<Type>::zero // null value
342  );
343 }
344 
345 
346 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
static void updateList(const labelList &newToOld, const T &nullValue, List< T > &elems)
Helper: reorder list according to map.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
vector point
Point is a vector.
Definition: point.H:41
const word & name() const
Return name.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:390
dimensioned< scalar > mag(const dimensioned< Type > &)
An STL-conforming hash table.
Definition: HashTable.H:61
static T gAverage(const PackedBoolList &isMasterElem, const UList< T > &values)
Helper: calculate average.
void sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
#define forAllIter(Container, container, iter)
Definition: UList.H:440
A bit-packed bool list.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
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
label size() const
Number of entries.
Definition: PackedListI.H:722
messageStream Info
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
const volScalarField & T
Definition: createFields.H:25
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
patches[0]
label n
static void collectAndPrint(const UList< point > &points, const UList< T > &data)
Print list according to (collected and) sorted coordinate.
const double e
Elementary charge.
Definition: doubleFloat.H:78
Calculates points shared by more than two processor patches or cyclic patches.
Definition: globalPoints.H:100
void setSize(const label)
Reset size of List.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Foam::polyBoundaryMesh.
static void weightedSum(const polyMesh &mesh, const PackedBoolList &isMasterEdge, const labelList &meshPoints, const edgeList &edges, const scalarField &edgeWeights, const Field< Type > &data, Field< Type > &sum)
Helper: weighted sum (over all subset of mesh points) by.
void testSyncBoundaryFaceList(const scalar mergeDistance, const string &, const UList< T > &, const UList< T > &) const
Compare two lists over all boundary faces.
#define forAll(list, i)
Definition: UList.H:421
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
static int readFlags(const Enum &namedEnum, const wordList &)
Helper: convert wordList into bit pattern using provided.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
static void syncPointList(const polyMesh &, List< T > &, const CombineOp &cop, const T &nullValue, const TransformOp &top)
Synchronize values on all mesh points.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
tmp< pointField > allPoints(const Triangulation &t)
Extract all points in vertex-index order.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
Traits class for primitives.
Definition: pTraits.H:50
error FatalError
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:552