FreeFOAM The Cross-Platform CFD Toolkit
sampledSurfaces.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) 1991-2010 OpenCFD Ltd.
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 "sampledSurfaces.H"
27 #include <finiteVolume/volFields.H>
28 #include <OpenFOAM/dictionary.H>
29 #include <OpenFOAM/Time.H>
30 #include <OpenFOAM/IOmanip.H>
31 #include <OpenFOAM/ListListOps.H>
32 #include <OpenFOAM/mergePoints.H>
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  //- Used to offset faces in Pstream::combineOffset
40  template <>
41  class offsetOp<face>
42  {
43 
44  public:
45 
46  face operator()
47  (
48  const face& x,
49  const label offset
50  ) const
51  {
52  face result(x.size());
53 
54  forAll(x, xI)
55  {
56  result[xI] = x[xI] + offset;
57  }
58  return result;
59  }
60  };
61 
62 
63  defineTypeNameAndDebug(sampledSurfaces, 0);
64 }
65 
66 
67 bool Foam::sampledSurfaces::verbose_(false);
68 Foam::scalar Foam::sampledSurfaces::mergeTol_(1e-10);
69 
70 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
71 
72 Foam::label Foam::sampledSurfaces::classifyFieldTypes()
73 {
74  label nFields = 0;
75 
76  scalarFields_.clear();
77  vectorFields_.clear();
78  sphericalTensorFields_.clear();
79  symmTensorFields_.clear();
80  tensorFields_.clear();
81 
82  forAll(fieldNames_, fieldI)
83  {
84  const word& fieldName = fieldNames_[fieldI];
85  word fieldType = "";
86 
87  // check files for a particular time
88  if (loadFromFiles_)
89  {
90  IOobject io
91  (
92  fieldName,
93  mesh_.time().timeName(),
94  mesh_,
97  false
98  );
99 
100  if (io.headerOk())
101  {
102  fieldType = io.headerClassName();
103  }
104  else
105  {
106  continue;
107  }
108  }
109  else
110  {
111  // check objectRegistry
112  objectRegistry::const_iterator iter = mesh_.find(fieldName);
113 
114  if (iter != mesh_.objectRegistry::end())
115  {
116  fieldType = iter()->type();
117  }
118  else
119  {
120  continue;
121  }
122  }
123 
124 
125  if (fieldType == volScalarField::typeName)
126  {
127  scalarFields_.append(fieldName);
128  nFields++;
129  }
130  else if (fieldType == volVectorField::typeName)
131  {
132  vectorFields_.append(fieldName);
133  nFields++;
134  }
135  else if (fieldType == volSphericalTensorField::typeName)
136  {
137  sphericalTensorFields_.append(fieldName);
138  nFields++;
139  }
140  else if (fieldType == volSymmTensorField::typeName)
141  {
142  symmTensorFields_.append(fieldName);
143  nFields++;
144  }
145  else if (fieldType == volTensorField::typeName)
146  {
147  tensorFields_.append(fieldName);
148  nFields++;
149  }
150 
151  }
152 
153  return nFields;
154 }
155 
156 
157 void Foam::sampledSurfaces::writeGeometry() const
158 {
159  // Write to time directory under outputPath_
160  // skip surface without faces (eg, a failed cut-plane)
161 
162  const fileName outputDir = outputPath_/mesh_.time().timeName();
163 
164  forAll(*this, surfI)
165  {
166  const sampledSurface& s = operator[](surfI);
167 
168  if (Pstream::parRun())
169  {
170  if (Pstream::master() && mergeList_[surfI].faces.size())
171  {
172  genericFormatter_->write
173  (
174  outputDir,
175  s.name(),
176  mergeList_[surfI].points,
177  mergeList_[surfI].faces
178  );
179  }
180  }
181  else if (s.faces().size())
182  {
183  genericFormatter_->write
184  (
185  outputDir,
186  s.name(),
187  s.points(),
188  s.faces()
189  );
190  }
191  }
192 }
193 
194 
195 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
196 
197 Foam::sampledSurfaces::sampledSurfaces
198 (
199  const word& name,
200  const objectRegistry& obr,
201  const dictionary& dict,
202  const bool loadFromFiles
203 )
204 :
206  name_(name),
207  mesh_(refCast<const fvMesh>(obr)),
208  loadFromFiles_(loadFromFiles),
209  outputPath_(fileName::null),
210  fieldNames_(),
211  interpolationScheme_(word::null),
212  writeFormat_(word::null),
213  mergeList_(),
214  genericFormatter_(NULL),
215  scalarFields_(),
216  vectorFields_(),
217  sphericalTensorFields_(),
218  symmTensorFields_(),
219  tensorFields_()
220 {
221  if (Pstream::parRun())
222  {
223  outputPath_ = mesh_.time().path()/".."/name_;
224  }
225  else
226  {
227  outputPath_ = mesh_.time().path()/name_;
228  }
229 
230  read(dict);
231 }
232 
233 
234 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
235 
237 {}
238 
239 
240 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
241 
242 void Foam::sampledSurfaces::verbose(const bool verbosity)
243 {
244  verbose_ = verbosity;
245 }
246 
247 
249 {
250  // Do nothing - only valid on write
251 }
252 
253 
255 {
256  // Do nothing - only valid on write
257 }
258 
259 
261 {
262  if (size())
263  {
264  // finalize surfaces, merge points etc.
265  update();
266 
267  const label nFields = classifyFieldTypes();
268 
269  if (Pstream::master())
270  {
271  if (debug)
272  {
273  Pout<< "timeName = " << mesh_.time().timeName() << nl
274  << "scalarFields " << scalarFields_ << nl
275  << "vectorFields " << vectorFields_ << nl
276  << "sphTensorFields " << sphericalTensorFields_ << nl
277  << "symTensorFields " << symmTensorFields_ <<nl
278  << "tensorFields " << tensorFields_ <<nl;
279 
280  Pout<< "Creating directory "
281  << outputPath_/mesh_.time().timeName() << nl << endl;
282 
283  }
284 
285  mkDir(outputPath_/mesh_.time().timeName());
286  }
287 
288  // write geometry first if required, or when no fields would otherwise
289  // be written
290  if (nFields == 0 || genericFormatter_->separateFiles())
291  {
292  writeGeometry();
293  }
294 
295  sampleAndWrite(scalarFields_);
296  sampleAndWrite(vectorFields_);
297  sampleAndWrite(sphericalTensorFields_);
298  sampleAndWrite(symmTensorFields_);
299  sampleAndWrite(tensorFields_);
300  }
301 }
302 
303 
305 {
306  fieldNames_ = wordList(dict.lookup("fields"));
307 
308  const label nFields = fieldNames_.size();
309 
310  scalarFields_.reset(nFields);
311  vectorFields_.reset(nFields);
312  sphericalTensorFields_.reset(nFields);
313  symmTensorFields_.reset(nFields);
314  tensorFields_.reset(nFields);
315 
316  interpolationScheme_ = dict.lookupOrDefault<word>
317  (
318  "interpolationScheme",
319  "cell"
320  );
321  writeFormat_ = dict.lookupOrDefault<word>
322  (
323  "surfaceFormat",
324  "null"
325  );
326 
327 
328  // define the generic (geometry) writer
329  genericFormatter_ = surfaceWriter<bool>::New(writeFormat_);
330 
331 
333  (
334  dict.lookup("surfaces"),
335  sampledSurface::iNew(mesh_)
336  );
337 
338  transfer(newList);
339 
340  if (Pstream::parRun())
341  {
342  mergeList_.setSize(size());
343  }
344 
345  // ensure all surfaces and merge information are expired
346  expire();
347 
348  if (Pstream::master() && debug)
349  {
350  Pout<< "sample fields:" << fieldNames_ << nl
351  << "sample surfaces:" << nl << "(" << nl;
352 
353  forAll(*this, surfI)
354  {
355  Pout<< " " << operator[](surfI) << endl;
356  }
357  Pout<< ")" << endl;
358  }
359 }
360 
361 
363 {
364  expire();
365 
366  // pointMesh and interpolation will have been reset in mesh.update
367 }
368 
369 
371 {
372  expire();
373 }
374 
375 
377 {
378  if (state != polyMesh::UNCHANGED)
379  {
380  expire();
381  }
382 }
383 
384 
386 {
387  forAll(*this, surfI)
388  {
389  if (operator[](surfI).needsUpdate())
390  {
391  return true;
392  }
393  }
394 
395  return false;
396 }
397 
398 
400 {
401  bool justExpired = false;
402 
403  forAll(*this, surfI)
404  {
405  if (operator[](surfI).expire())
406  {
407  justExpired = true;
408  }
409 
410  // clear merge information
411  if (Pstream::parRun())
412  {
413  mergeList_[surfI].clear();
414  }
415  }
416 
417  // true if any surfaces just expired
418  return justExpired;
419 }
420 
421 
423 {
424  bool updated = false;
425 
426  if (!needsUpdate())
427  {
428  return updated;
429  }
430 
431  // serial: quick and easy, no merging required
432  if (!Pstream::parRun())
433  {
434  forAll(*this, surfI)
435  {
436  if (operator[](surfI).update())
437  {
438  updated = true;
439  }
440  }
441 
442  return updated;
443  }
444 
445  // dimension as fraction of mesh bounding box
446  scalar mergeDim = mergeTol_ * mesh_.globalData().bb().mag();
447 
448  if (Pstream::master() && debug)
449  {
450  Pout<< nl << "Merging all points within "
451  << mergeDim << " meter" << endl;
452  }
453 
454  forAll(*this, surfI)
455  {
456  sampledSurface& s = operator[](surfI);
457 
458  if (s.update())
459  {
460  updated = true;
461  }
462  else
463  {
464  continue;
465  }
466 
467 
468  // Collect points from all processors
469  List<pointField> gatheredPoints(Pstream::nProcs());
470  gatheredPoints[Pstream::myProcNo()] = s.points();
471  Pstream::gatherList(gatheredPoints);
472 
473  if (Pstream::master())
474  {
475  mergeList_[surfI].points = ListListOps::combine<pointField>
476  (
477  gatheredPoints,
479  );
480  }
481 
482  // Collect faces from all processors and renumber using sizes of
483  // gathered points
484  List<faceList> gatheredFaces(Pstream::nProcs());
485  gatheredFaces[Pstream::myProcNo()] = s.faces();
486  Pstream::gatherList(gatheredFaces);
487 
488  if (Pstream::master())
489  {
490  mergeList_[surfI].faces = static_cast<const faceList&>
491  (
492  ListListOps::combineOffset<faceList>
493  (
494  gatheredFaces,
496  (
497  gatheredPoints,
499  ),
502  )
503  );
504  }
505 
506  pointField newPoints;
507  labelList oldToNew;
508 
509  bool hasMerged = mergePoints
510  (
511  mergeList_[surfI].points,
512  mergeDim,
513  false, // verbosity
514  oldToNew,
515  newPoints
516  );
517 
518  if (hasMerged)
519  {
520  // Store point mapping
521  mergeList_[surfI].pointsMap.transfer(oldToNew);
522 
523  // Copy points
524  mergeList_[surfI].points.transfer(newPoints);
525 
526  // Relabel faces
527  faceList& faces = mergeList_[surfI].faces;
528 
529  forAll(faces, faceI)
530  {
531  inplaceRenumber(mergeList_[surfI].pointsMap, faces[faceI]);
532  }
533 
534  if (Pstream::master() && debug)
535  {
536  Pout<< "For surface " << surfI << " merged from "
537  << mergeList_[surfI].pointsMap.size() << " points down to "
538  << mergeList_[surfI].points.size() << " points" << endl;
539  }
540  }
541  }
542 
543  return updated;
544 }
545 
546 
547 // ************************ vim: set sw=4 sts=4 et: ************************ //