wallBoundedStreamLine.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-2016 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 "Pstream.H"
27 #include "functionObjectList.H"
28 #include "wallBoundedStreamLine.H"
29 #include "fvMesh.H"
31 #include "ReadFields.H"
32 #include "meshSearch.H"
33 #include "sampledSet.H"
34 #include "globalIndex.H"
35 #include "mapDistribute.H"
36 #include "interpolationCellPoint.H"
37 #include "PatchTools.H"
38 #include "meshSearchMeshObject.H"
39 #include "faceSet.H"
40 #include "mapPolyMesh.H"
42 
43 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 namespace functionObjects
48 {
49  defineTypeNameAndDebug(wallBoundedStreamLine, 0);
50 
52  (
53  functionObject,
54  wallBoundedStreamLine,
55  dictionary
56  );
57 }
58 }
59 
60 
61 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
62 
64 Foam::functionObjects::wallBoundedStreamLine::wallPatch() const
65 {
66  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
67 
68  const polyBoundaryMesh& patches = mesh.boundaryMesh();
69 
70  label nFaces = 0;
71 
72  forAll(patches, patchi)
73  {
74  //if (!polyPatch::constraintType(patches[patchi].type()))
75  if (isA<wallPolyPatch>(patches[patchi]))
76  {
77  nFaces += patches[patchi].size();
78  }
79  }
80 
81  labelList addressing(nFaces);
82 
83  nFaces = 0;
84 
85  forAll(patches, patchi)
86  {
87  //if (!polyPatch::constraintType(patches[patchi].type()))
88  if (isA<wallPolyPatch>(patches[patchi]))
89  {
90  const polyPatch& pp = patches[patchi];
91 
92  forAll(pp, i)
93  {
94  addressing[nFaces++] = pp.start()+i;
95  }
96  }
97  }
98 
99  return autoPtr<indirectPrimitivePatch>
100  (
102  (
103  IndirectList<face>
104  (
105  mesh.faces(),
106  addressing
107  ),
108  mesh.points()
109  )
110  );
111 }
112 
113 
114 Foam::tetIndices Foam::functionObjects::wallBoundedStreamLine::findNearestTet
115 (
116  const PackedBoolList& isWallPatch,
117  const point& seedPt,
118  const label celli
119 ) const
120 {
121  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
122 
123  const cell& cFaces = mesh.cells()[celli];
124 
125  label minFacei = -1;
126  label minTetPtI = -1;
127  scalar minDistSqr = sqr(GREAT);
128 
129  forAll(cFaces, cFacei)
130  {
131  label facei = cFaces[cFacei];
132 
133  if (isWallPatch[facei])
134  {
135  const face& f = mesh.faces()[facei];
136  const label fp0 = mesh.tetBasePtIs()[facei];
137  const point& basePoint = mesh.points()[f[fp0]];
138 
139  label fp = f.fcIndex(fp0);
140  for (label i = 2; i < f.size(); i++)
141  {
142  const point& thisPoint = mesh.points()[f[fp]];
143  label nextFp = f.fcIndex(fp);
144  const point& nextPoint = mesh.points()[f[nextFp]];
145 
146  const triPointRef tri(basePoint, thisPoint, nextPoint);
147 
148  scalar d2 = magSqr(tri.centre() - seedPt);
149  if (d2 < minDistSqr)
150  {
151  minDistSqr = d2;
152  minFacei = facei;
153  minTetPtI = i-1;
154  }
155  fp = nextFp;
156  }
157  }
158  }
159 
160  // Put particle in tet
161  return tetIndices
162  (
163  celli,
164  minFacei,
165  minTetPtI,
166  mesh
167  );
168 }
169 
170 
171 void Foam::functionObjects::wallBoundedStreamLine::track()
172 {
173  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
174 
175 
176  // Determine the 'wall' patches
177  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
178  // These are the faces that need to be followed
179 
180  autoPtr<indirectPrimitivePatch> boundaryPatch(wallPatch());
181  PackedBoolList isWallPatch(mesh.nFaces());
182  forAll(boundaryPatch().addressing(), i)
183  {
184  isWallPatch[boundaryPatch().addressing()[i]] = 1;
185  }
186 
187 
188 
189  // Find nearest wall particle for the seedPoints
190  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
191 
192  IDLList<wallBoundedStreamLineParticle> initialParticles;
193  wallBoundedStreamLineParticleCloud particles
194  (
195  mesh,
196  cloudName_,
197  initialParticles
198  );
199 
200  {
201  // Get the seed points
202  // ~~~~~~~~~~~~~~~~~~~
203 
204  const sampledSet& seedPoints = sampledSetPtr_();
205 
206 
207  forAll(seedPoints, i)
208  {
209  const point& seedPt = seedPoints[i];
210  label celli = seedPoints.cells()[i];
211 
212  tetIndices ids(findNearestTet(isWallPatch, seedPt, celli));
213 
214  if (ids.face() != -1 && isWallPatch[ids.face()])
215  {
216  //Pout<< "Seeding particle :" << nl
217  // << " seedPt:" << seedPt << nl
218  // << " face :" << ids.face() << nl
219  // << " at :" << mesh.faceCentres()[ids.face()] << nl
220  // << " cell :" << mesh.cellCentres()[ids.cell()] << nl
221  // << endl;
222 
223  particles.addParticle
224  (
225  new wallBoundedStreamLineParticle
226  (
227  mesh,
228  ids.faceTri(mesh).centre(),
229  ids.cell(),
230  ids.face(), // tetFace
231  ids.tetPt(),
232  -1, // not on a mesh edge
233  -1, // not on a diagonal edge
234  lifeTime_ // lifetime
235  )
236  );
237  }
238  else
239  {
240  Pout<< type() << " : ignoring seed " << seedPt
241  << " since not in wall cell." << endl;
242  }
243  }
244  }
245 
246  label nSeeds = returnReduce(particles.size(), sumOp<label>());
247 
248  Info<< type() << " : seeded " << nSeeds << " particles." << endl;
249 
250 
251 
252  // Read or lookup fields
253  PtrList<volScalarField> vsFlds;
254  PtrList<interpolation<scalar>> vsInterp;
255  PtrList<volVectorField> vvFlds;
256  PtrList<interpolation<vector>> vvInterp;
257 
258  label UIndex = -1;
259 
260  label nScalar = 0;
261  label nVector = 0;
262 
263  forAll(fields_, i)
264  {
265  if (mesh.foundObject<volScalarField>(fields_[i]))
266  {
267  nScalar++;
268  }
269  else if (mesh.foundObject<volVectorField>(fields_[i]))
270  {
271  nVector++;
272  }
273  else
274  {
276  << "Cannot find field " << fields_[i] << endl
277  << "Valid scalar fields are:"
278  << mesh.names(volScalarField::typeName) << endl
279  << "Valid vector fields are:"
280  << mesh.names(volVectorField::typeName)
281  << exit(FatalError);
282  }
283  }
284 
285  vsInterp.setSize(nScalar);
286  nScalar = 0;
287  vvInterp.setSize(nVector);
288  nVector = 0;
289 
290  forAll(fields_, i)
291  {
292  if (mesh.foundObject<volScalarField>(fields_[i]))
293  {
294  const volScalarField& f = mesh.lookupObject<volScalarField>
295  (
296  fields_[i]
297  );
298  vsInterp.set
299  (
300  nScalar++,
302  (
303  interpolationScheme_,
304  f
305  )
306  );
307  }
308  else if (mesh.foundObject<volVectorField>(fields_[i]))
309  {
310  const volVectorField& f = mesh.lookupObject<volVectorField>
311  (
312  fields_[i]
313  );
314 
315  if (f.name() == UName_)
316  {
317  UIndex = nVector;
318  }
319 
320  vvInterp.set
321  (
322  nVector++,
324  (
325  interpolationScheme_,
326  f
327  )
328  );
329  }
330  }
331 
332  // Store the names
333  scalarNames_.setSize(vsInterp.size());
334  forAll(vsInterp, i)
335  {
336  scalarNames_[i] = vsInterp[i].psi().name();
337  }
338  vectorNames_.setSize(vvInterp.size());
339  forAll(vvInterp, i)
340  {
341  vectorNames_[i] = vvInterp[i].psi().name();
342  }
343 
344  // Check that we know the index of U in the interpolators.
345 
346  if (UIndex == -1)
347  {
349  << "Cannot find field to move particles with : " << UName_
350  << endl
351  << "This field has to be present in the sampled fields "
352  << fields_
353  << " and in the objectRegistry." << endl
354  << exit(FatalError);
355  }
356 
357  // Sampled data
358  // ~~~~~~~~~~~~
359 
360  // Size to maximum expected sizes.
361  allTracks_.clear();
362  allTracks_.setCapacity(nSeeds);
363  allScalars_.setSize(vsInterp.size());
364  forAll(allScalars_, i)
365  {
366  allScalars_[i].clear();
367  allScalars_[i].setCapacity(nSeeds);
368  }
369  allVectors_.setSize(vvInterp.size());
370  forAll(allVectors_, i)
371  {
372  allVectors_[i].clear();
373  allVectors_[i].setCapacity(nSeeds);
374  }
375 
376 
377  // additional particle info
378  wallBoundedStreamLineParticle::trackingData td
379  (
380  particles,
381  vsInterp,
382  vvInterp,
383  UIndex, // index of U in vvInterp
384  trackForward_, // track in +u direction?
385  trackLength_, // fixed track length
386  isWallPatch, // which faces are to follow
387 
388  allTracks_,
389  allScalars_,
390  allVectors_
391  );
392 
393 
394  // Set very large dt. Note: cannot use GREAT since 1/GREAT is SMALL
395  // which is a trigger value for the tracking...
396  const scalar trackTime = Foam::sqrt(GREAT);
397 
398  // Track
399  particles.move(td, trackTime);
400 }
401 
402 
403 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
404 
405 Foam::functionObjects::wallBoundedStreamLine::wallBoundedStreamLine
406 (
407  const word& name,
408  const Time& runTime,
409  const dictionary& dict
410 )
411 :
412  functionObject(name),
413  obr_
414  (
416  (
418  )
419  ),
420  dict_(dict)
421 {
422  if (!isA<fvMesh>(obr_))
423  {
425  << "objectRegistry is not an fvMesh" << exit(FatalError);
426  }
427 
428  read(dict_);
429 }
430 
431 
432 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
433 
435 {}
436 
437 
438 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
439 
441 {
442  //dict_ = dict;
443  dict.lookup("fields") >> fields_;
444  if (dict.found("U"))
445  {
446  dict.lookup("U") >> UName_;
447  }
448  else
449  {
450  UName_ = "U";
451  if (dict.found("U"))
452  {
454  (
455  dict
456  ) << "Using deprecated entry \"U\"."
457  << " Please use \"UName\" instead."
458  << endl;
459  dict.lookup("U") >> UName_;
460  }
461  }
462 
463  if (findIndex(fields_, UName_) == -1)
464  {
466  (
467  dict
468  ) << "Velocity field for tracking " << UName_
469  << " should be present in the list of fields " << fields_
470  << exit(FatalIOError);
471  }
472 
473 
474  dict.lookup("trackForward") >> trackForward_;
475  dict.lookup("lifeTime") >> lifeTime_;
476  if (lifeTime_ < 1)
477  {
479  << "Illegal value " << lifeTime_ << " for lifeTime"
480  << exit(FatalError);
481  }
482  trackLength_ = VGREAT;
483  if (dict.found("trackLength"))
484  {
485  dict.lookup("trackLength") >> trackLength_;
486 
487  Info<< type() << " : fixed track length specified : "
488  << trackLength_ << nl << endl;
489  }
490 
491 
492  interpolationScheme_ = dict.lookupOrDefault
493  (
494  "interpolationScheme",
496  );
497 
498  cloudName_ = dict.lookupOrDefault<word>
499  (
500  "cloudName",
501  "wallBoundedStreamLine"
502  );
503  dict.lookup("seedSampleSet") >> seedSet_;
504 
505  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
506 
507  const dictionary& coeffsDict = dict.subDict(seedSet_ + "Coeffs");
508 
509  sampledSetPtr_ = sampledSet::New
510  (
511  seedSet_,
512  mesh,
514  coeffsDict
515  );
516  coeffsDict.lookup("axis") >> sampledSetAxis_;
517 
518  scalarFormatterPtr_ = writer<scalar>::New(dict.lookup("setFormat"));
519  vectorFormatterPtr_ = writer<vector>::New(dict.lookup("setFormat"));
520 
521 
522  // Make sure that the mesh is trackable
523  if (debug)
524  {
525  // 1. positive volume decomposition tets
526  faceSet faces(mesh, "lowQualityTetFaces", mesh.nFaces()/100+1);
527  if
528  (
530  (
531  mesh,
533  true,
534  &faces
535  )
536  )
537  {
538  label nFaces = returnReduce(faces.size(), sumOp<label>());
539 
541  << "Found " << nFaces
542  <<" faces with low quality or negative volume "
543  << "decomposition tets. Writing to faceSet " << faces.name()
544  << endl;
545  }
546 
547  // 2. all edges on a cell having two faces
548  EdgeMap<label> numFacesPerEdge;
549  forAll(mesh.cells(), celli)
550  {
551  const cell& cFaces = mesh.cells()[celli];
552 
553  numFacesPerEdge.clear();
554 
555  forAll(cFaces, cFacei)
556  {
557  label facei = cFaces[cFacei];
558  const face& f = mesh.faces()[facei];
559  forAll(f, fp)
560  {
561  const edge e(f[fp], f.nextLabel(fp));
563  numFacesPerEdge.find(e);
564  if (eFnd != numFacesPerEdge.end())
565  {
566  eFnd()++;
567  }
568  else
569  {
570  numFacesPerEdge.insert(e, 1);
571  }
572  }
573  }
574 
575  forAllConstIter(EdgeMap<label>, numFacesPerEdge, iter)
576  {
577  if (iter() != 2)
578  {
580  << "problem cell:" << celli
581  << abort(FatalError);
582  }
583  }
584  }
585  }
586 
587  return true;
588 }
589 
590 
592 {
593  return true;
594 }
595 
596 
598 {
599  const Time& runTime = obr_.time();
600  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
601 
602 
603  // Do all injection and tracking
604  track();
605 
606 
607  if (Pstream::parRun())
608  {
609  // Append slave tracks to master ones
610  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
611 
612  globalIndex globalTrackIDs(allTracks_.size());
613 
614  // Construct a distribution map to pull all to the master.
615  labelListList sendMap(Pstream::nProcs());
616  labelListList recvMap(Pstream::nProcs());
617 
618  if (Pstream::master())
619  {
620  // Master: receive all. My own first, then consecutive
621  // processors.
622  label trackI = 0;
623 
624  forAll(recvMap, proci)
625  {
626  labelList& fromProc = recvMap[proci];
627  fromProc.setSize(globalTrackIDs.localSize(proci));
628  forAll(fromProc, i)
629  {
630  fromProc[i] = trackI++;
631  }
632  }
633  }
634 
635  labelList& toMaster = sendMap[0];
636  toMaster.setSize(globalTrackIDs.localSize());
637  forAll(toMaster, i)
638  {
639  toMaster[i] = i;
640  }
641 
642  const mapDistribute distMap
643  (
644  globalTrackIDs.size(),
645  sendMap.xfer(),
646  recvMap.xfer()
647  );
648 
649  // Distribute the track positions. Note: use scheduled comms
650  // to prevent buffering.
651  allTracks_.shrink();
653  (
655  distMap.schedule(),
656  distMap.constructSize(),
657  distMap.subMap(),
658  false,
659  distMap.constructMap(),
660  false,
661  allTracks_,
662  flipOp()
663  );
664 
665  // Distribute the scalars
666  forAll(allScalars_, scalarI)
667  {
668  allScalars_[scalarI].shrink();
670  (
672  distMap.schedule(),
673  distMap.constructSize(),
674  distMap.subMap(),
675  false,
676  distMap.constructMap(),
677  false,
678  allScalars_[scalarI],
679  flipOp()
680  );
681  allScalars_[scalarI].setCapacity(allScalars_[scalarI].size());
682  }
683  // Distribute the vectors
684  forAll(allVectors_, vectorI)
685  {
686  allVectors_[vectorI].shrink();
688  (
690  distMap.schedule(),
691  distMap.constructSize(),
692  distMap.subMap(),
693  false,
694  distMap.constructMap(),
695  false,
696  allVectors_[vectorI],
697  flipOp()
698  );
699  allVectors_[vectorI].setCapacity(allVectors_[vectorI].size());
700  }
701  }
702 
703 
704  label n = 0;
705  forAll(allTracks_, trackI)
706  {
707  n += allTracks_[trackI].size();
708  }
709 
710  Info<< " Tracks:" << allTracks_.size() << nl
711  << " Total samples:" << n << endl;
712 
713 
714  // Massage into form suitable for writers
715  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
716 
717  if (Pstream::master() && allTracks_.size())
718  {
719  // Make output directory
720 
721  fileName vtkPath
722  (
724  ? runTime.path()/".."/"postProcessing"/"sets"/name()
725  : runTime.path()/"postProcessing"/"sets"/name()
726  );
727  if (mesh.name() != fvMesh::defaultRegion)
728  {
729  vtkPath = vtkPath/mesh.name();
730  }
731  vtkPath = vtkPath/mesh.time().timeName();
732 
733  mkDir(vtkPath);
734 
735  // Convert track positions
736 
737  PtrList<coordSet> tracks(allTracks_.size());
738  forAll(allTracks_, trackI)
739  {
740  tracks.set
741  (
742  trackI,
743  new coordSet
744  (
745  "track" + Foam::name(trackI),
746  sampledSetAxis_ //"xyz"
747  )
748  );
749  tracks[trackI].transfer(allTracks_[trackI]);
750  }
751 
752  // Convert scalar values
753 
754  if (allScalars_.size() > 0)
755  {
756  List<List<scalarField>> scalarValues(allScalars_.size());
757 
758  forAll(allScalars_, scalarI)
759  {
760  DynamicList<scalarList>& allTrackVals =
761  allScalars_[scalarI];
762  scalarValues[scalarI].setSize(allTrackVals.size());
763 
764  forAll(allTrackVals, trackI)
765  {
766  scalarList& trackVals = allTrackVals[trackI];
767  scalarValues[scalarI][trackI].transfer(trackVals);
768  }
769  }
770 
771  fileName vtkFile
772  (
773  vtkPath
774  / scalarFormatterPtr_().getFileName
775  (
776  tracks[0],
777  scalarNames_
778  )
779  );
780 
781  Info<< " Writing data to " << vtkFile.path() << endl;
782 
783  scalarFormatterPtr_().write
784  (
785  true, // writeTracks
786  tracks,
787  scalarNames_,
788  scalarValues,
789  OFstream(vtkFile)()
790  );
791  }
792 
793  // Convert vector values
794 
795  if (allVectors_.size() > 0)
796  {
797  List<List<vectorField>> vectorValues(allVectors_.size());
798 
799  forAll(allVectors_, vectorI)
800  {
801  DynamicList<vectorList>& allTrackVals =
802  allVectors_[vectorI];
803  vectorValues[vectorI].setSize(allTrackVals.size());
804 
805  forAll(allTrackVals, trackI)
806  {
807  vectorList& trackVals = allTrackVals[trackI];
808  vectorValues[vectorI][trackI].transfer(trackVals);
809  }
810  }
811 
812  fileName vtkFile
813  (
814  vtkPath
815  / vectorFormatterPtr_().getFileName
816  (
817  tracks[0],
818  vectorNames_
819  )
820  );
821 
822  vectorFormatterPtr_().write
823  (
824  true, // writeTracks
825  tracks,
826  vectorNames_,
827  vectorValues,
828  OFstream(vtkFile)()
829  );
830  }
831  }
832 
833  return true;
834 }
835 
836 
838 (
839  const mapPolyMesh& mpm
840 )
841 {
842  const fvMesh& mesh_ = dynamic_cast<const fvMesh&>(obr_);
843 
844  if (&mpm.mesh() == &mesh_)
845  {
846  read(dict_);
847  }
848 }
849 
850 
852 (
853  const polyMesh& mesh
854 )
855 {
856  const fvMesh& mesh_ = dynamic_cast<const fvMesh&>(obr_);
857 
858  if (&mesh == &mesh_)
859  {
860  // Moving mesh affects the search tree
861  read(dict_);
862  }
863 }
864 
865 
866 // ************************************************************************* //
const Time & time() const
Return time.
Class containing functor to negate primitives. Dummy for all other types.
Definition: flipOp.H:50
virtual bool write()
Calculate and write the wall-bounded streamlines.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
A class for handling file names.
Definition: fileName.H:69
A list of face labels.
Definition: faceSet.H:48
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
const double e
Elementary charge.
Definition: doubleFloat.H:78
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
static autoPtr< sampledSet > New(const word &name, const polyMesh &mesh, const meshSearch &searchEngine, const dictionary &dict)
Return a reference to the selected sampledSet.
Definition: sampledSet.C:442
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Xfer< List< T > > xfer()
Transfer contents to the Xfer container.
Definition: ListI.H:90
static const char *const typeName
Definition: Field.H:94
static bool checkFaceTets(const polyMesh &mesh, scalar tol=minTetQuality, const bool report=false, labelHashSet *setPtr=NULL)
Check face-decomposition tet volume.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Output to file stream.
Definition: OFstream.H:81
static autoPtr< interpolation< scalar > > New(const word &interpolationType, const GeometricField< scalar, fvPatchField, volMesh > &psi)
Return a reference to the specified interpolation scheme.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
virtual bool read(const dictionary &)
Read the field average data.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:411
Field reading functions for post-processing utilities.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
Abstract base-class for Time/database function objects.
PrimitivePatch< face, IndirectList, const pointField & > indirectPrimitivePatch
Foam::indirectPrimitivePatch.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
const cellList & cells() const
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Macros for easy insertion into run-time selection tables.
static const meshSearchMeshObject & New(const polyMesh &mesh)
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:138
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
bool read(const char *, int32_t &)
Definition: int32IO.C:85
void setSize(const label)
Alter the addressed list size.
Definition: DynamicListI.H:163
Holds list of sampling positions.
Definition: coordSet.H:49
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
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
void clear()
Clear all entries from table.
Definition: HashTable.C:464
static void distribute(const Pstream::commsTypes commsType, const List< labelPair > &schedule, const label constructSize, const labelListList &subMap, const bool subHasFlip, const labelListList &constructMap, const bool constructHasFlip, List< T > &, const negateOp &negOp, const int tag=UPstream::msgType())
Distribute data. Note:schedule only used for Pstream::scheduled.
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
List< label > labelList
A List of labels.
Definition: labelList.H:56
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
Definition: tetIndices.H:81
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
static autoPtr< writer > New(const word &writeFormat)
Return a reference to the selected writer.
Definition: writer.C:35
dimensioned< scalar > magSqr(const dimensioned< Type > &)
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
static const char nl
Definition: Ostream.H:262
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:295
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label nFaces() const
void setSize(const label)
Reset size of List.
Definition: List.C:295
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:393
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:399
label patchi
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
Class containing processor-to-processor mapping information.
vector point
Point is a vector.
Definition: point.H:41
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
#define WarningInFunction
Report a warning using Foam::Warning.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
addToRunTimeSelectionTable(functionObject, blendingFactor, dictionary)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual void updateMesh(const mapPolyMesh &)
Update for changes of mesh.
triangle< point, const point & > triPointRef
Definition: triPointRef.H:44
messageStream Info
label n
label nextLabel(const label i) const
Next vertex on face.
Definition: faceI.H:117
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
virtual Ostream & write(const token &)=0
Write next token to stream.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
Registry of regIOobjects.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1004
fileName path() const
Return path.
Definition: Time.H:269
void transfer(List< T > &)
Transfer contents of the argument List into this.
Definition: DynamicListI.H:259
fileName path() const
Return directory path name (part before last /)
Definition: fileName.C:249
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
static const scalar minTetQuality
Minimum tetrahedron quality.
const word & name() const
Return reference to name.
Definition: fvMesh.H:257
const polyMesh & mesh() const
Return polyMesh.
Definition: mapPolyMesh.H:357
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451
IOerror FatalIOError