streamLine.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 "Pstream.H"
27 #include "functionObjectList.H"
28 #include "streamLine.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 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43 defineTypeNameAndDebug(streamLine, 0);
44 }
45 
46 
47 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
48 
50 Foam::streamLine::wallPatch() const
51 {
52  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
53 
54  const polyBoundaryMesh& patches = mesh.boundaryMesh();
55 
56  label nFaces = 0;
57 
58  forAll(patches, patchI)
59  {
60  //if (!polyPatch::constraintType(patches[patchI].type()))
61  if (isA<wallPolyPatch>(patches[patchI]))
62  {
63  nFaces += patches[patchI].size();
64  }
65  }
66 
67  labelList addressing(nFaces);
68 
69  nFaces = 0;
70 
71  forAll(patches, patchI)
72  {
73  //if (!polyPatch::constraintType(patches[patchI].type()))
74  if (isA<wallPolyPatch>(patches[patchI]))
75  {
76  const polyPatch& pp = patches[patchI];
77 
78  forAll(pp, i)
79  {
80  addressing[nFaces++] = pp.start()+i;
81  }
82  }
83  }
84 
85  return autoPtr<indirectPrimitivePatch>
86  (
88  (
89  IndirectList<face>
90  (
91  mesh.faces(),
92  addressing
93  ),
94  mesh.points()
95  )
96  );
97 }
98 
99 
100 void Foam::streamLine::track()
101 {
102  const Time& runTime = obr_.time();
103  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
104 
105  IDLList<streamLineParticle> initialParticles;
106  streamLineParticleCloud particles
107  (
108  mesh,
109  cloudName_,
110  initialParticles
111  );
112 
113  const sampledSet& seedPoints = sampledSetPtr_();
114 
115  forAll(seedPoints, i)
116  {
117  particles.addParticle
118  (
119  new streamLineParticle
120  (
121  mesh,
122  seedPoints[i],
123  seedPoints.cells()[i],
124  lifeTime_ // lifetime
125  )
126  );
127  }
128 
129  label nSeeds = returnReduce(particles.size(), sumOp<label>());
130 
131  Info << " seeded " << nSeeds << " particles" << endl;
132 
133  // Read or lookup fields
134  PtrList<volScalarField> vsFlds;
135  PtrList<interpolation<scalar> > vsInterp;
136  PtrList<volVectorField> vvFlds;
137  PtrList<interpolation<vector> > vvInterp;
138 
139  label UIndex = -1;
140 
141  if (loadFromFiles_)
142  {
143  IOobjectList allObjects(mesh, runTime.timeName());
144 
145  IOobjectList objects(2*fields_.size());
146  forAll(fields_, i)
147  {
148  objects.add(*allObjects[fields_[i]]);
149  }
150 
151  ReadFields(mesh, objects, vsFlds);
152  vsInterp.setSize(vsFlds.size());
153  forAll(vsFlds, i)
154  {
155  vsInterp.set
156  (
157  i,
159  (
160  interpolationScheme_,
161  vsFlds[i]
162  )
163  );
164  }
165  ReadFields(mesh, objects, vvFlds);
166  vvInterp.setSize(vvFlds.size());
167  forAll(vvFlds, i)
168  {
169  vvInterp.set
170  (
171  i,
173  (
174  interpolationScheme_,
175  vvFlds[i]
176  )
177  );
178  }
179  }
180  else
181  {
182  label nScalar = 0;
183  label nVector = 0;
184 
185  forAll(fields_, i)
186  {
187  if (mesh.foundObject<volScalarField>(fields_[i]))
188  {
189  nScalar++;
190  }
191  else if (mesh.foundObject<volVectorField>(fields_[i]))
192  {
193  nVector++;
194  }
195  else
196  {
197  FatalErrorIn("streamLine::track()")
198  << "Cannot find field " << fields_[i] << nl
199  << "Valid scalar fields are:"
200  << mesh.names(volScalarField::typeName) << nl
201  << "Valid vector fields are:"
202  << mesh.names(volVectorField::typeName)
203  << exit(FatalError);
204  }
205  }
206  vsInterp.setSize(nScalar);
207  nScalar = 0;
208  vvInterp.setSize(nVector);
209  nVector = 0;
210 
211  forAll(fields_, i)
212  {
213  if (mesh.foundObject<volScalarField>(fields_[i]))
214  {
215  const volScalarField& f = mesh.lookupObject<volScalarField>
216  (
217  fields_[i]
218  );
219  vsInterp.set
220  (
221  nScalar++,
223  (
224  interpolationScheme_,
225  f
226  )
227  );
228  }
229  else if (mesh.foundObject<volVectorField>(fields_[i]))
230  {
231  const volVectorField& f = mesh.lookupObject<volVectorField>
232  (
233  fields_[i]
234  );
235 
236  if (f.name() == UName_)
237  {
238  UIndex = nVector;
239  }
240 
241  vvInterp.set
242  (
243  nVector++,
245  (
246  interpolationScheme_,
247  f
248  )
249  );
250  }
251  }
252  }
253 
254  // Store the names
255  scalarNames_.setSize(vsInterp.size());
256  forAll(vsInterp, i)
257  {
258  scalarNames_[i] = vsInterp[i].psi().name();
259  }
260  vectorNames_.setSize(vvInterp.size());
261  forAll(vvInterp, i)
262  {
263  vectorNames_[i] = vvInterp[i].psi().name();
264  }
265 
266  // Check that we know the index of U in the interpolators.
267 
268  if (UIndex == -1)
269  {
270  FatalErrorIn("streamLine::track()")
271  << "Cannot find field to move particles with : " << UName_ << nl
272  << "This field has to be present in the sampled fields " << fields_
273  << " and in the objectRegistry."
274  << exit(FatalError);
275  }
276 
277  // Sampled data
278  // ~~~~~~~~~~~~
279 
280  // Size to maximum expected sizes.
281  allTracks_.clear();
282  allTracks_.setCapacity(nSeeds);
283  allScalars_.setSize(vsInterp.size());
284  forAll(allScalars_, i)
285  {
286  allScalars_[i].clear();
287  allScalars_[i].setCapacity(nSeeds);
288  }
289  allVectors_.setSize(vvInterp.size());
290  forAll(allVectors_, i)
291  {
292  allVectors_[i].clear();
293  allVectors_[i].setCapacity(nSeeds);
294  }
295 
296 
297  // additional particle info
298  streamLineParticle::trackingData td
299  (
300  particles,
301  vsInterp,
302  vvInterp,
303  UIndex, // index of U in vvInterp
304  trackForward_, // track in +u direction?
305  nSubCycle_, // automatic track control:step through cells in steps?
306  trackLength_, // fixed track length
307 
308  allTracks_,
309  allScalars_,
310  allVectors_
311  );
312 
313 
314  // Set very large dt. Note: cannot use GREAT since 1/GREAT is SMALL
315  // which is a trigger value for the tracking...
316  const scalar trackTime = Foam::sqrt(GREAT);
317 
318  // Track
319  particles.move(td, trackTime);
320 }
321 
322 
323 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
324 
325 Foam::streamLine::streamLine
326 (
327  const word& name,
328  const objectRegistry& obr,
329  const dictionary& dict,
330  const bool loadFromFiles
331 )
332 :
333  dict_(dict),
334  name_(name),
335  obr_(obr),
336  loadFromFiles_(loadFromFiles),
337  active_(true),
338  nSubCycle_(0)
339 {
340  // Only active if a fvMesh is available
341  if (isA<fvMesh>(obr_))
342  {
343  read(dict_);
344  }
345  else
346  {
347  active_ = false;
348  WarningIn
349  (
350  "streamLine::streamLine\n"
351  "(\n"
352  "const word&,\n"
353  "const objectRegistry&,\n"
354  "const dictionary&,\n"
355  "const bool\n"
356  ")"
357  ) << "No fvMesh available, deactivating."
358  << nl << endl;
359  }
360 }
361 
362 
363 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
364 
366 {}
367 
368 
369 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
370 
372 {
373  if (active_)
374  {
375  Info<< type() << " " << name_ << ":" << nl;
376 
377  //dict_ = dict;
378  dict.lookup("fields") >> fields_;
379  if (dict.found("UName"))
380  {
381  dict.lookup("UName") >> UName_;
382  }
383  else
384  {
385  UName_ = "U";
386  if (dict.found("U"))
387  {
388  IOWarningIn("streamLine::read(const dictionary&)", dict)
389  << "Using deprecated entry \"U\"."
390  << " Please use \"UName\" instead."
391  << endl;
392  dict.lookup("U") >> UName_;
393  }
394  }
395 
396  if (findIndex(fields_, UName_) == -1)
397  {
398  FatalIOErrorIn("streamLine::read(const dictionary&)", dict)
399  << "Velocity field for tracking " << UName_
400  << " should be present in the list of fields " << fields_
401  << exit(FatalIOError);
402  }
403 
404 
405  dict.lookup("trackForward") >> trackForward_;
406  dict.lookup("lifeTime") >> lifeTime_;
407  if (lifeTime_ < 1)
408  {
409  FatalErrorIn(":streamLine::read(const dictionary&)")
410  << "Illegal value " << lifeTime_ << " for lifeTime"
411  << exit(FatalError);
412  }
413 
414 
415  bool subCycling = dict.found("nSubCycle");
416  bool fixedLength = dict.found("trackLength");
417 
418  if (subCycling && fixedLength)
419  {
420  FatalIOErrorIn("streamLine::read(const dictionary&)", dict)
421  << "Cannot both specify automatic time stepping (through '"
422  << "nSubCycle' specification) and fixed track length (through '"
423  << "trackLength')"
424  << exit(FatalIOError);
425  }
426 
427 
428  nSubCycle_ = 1;
429  if (dict.readIfPresent("nSubCycle", nSubCycle_))
430  {
431  trackLength_ = VGREAT;
432  if (nSubCycle_ < 1)
433  {
434  nSubCycle_ = 1;
435  }
436  Info<< " automatic track length specified through"
437  << " number of sub cycles : " << nSubCycle_ << nl << endl;
438  }
439  else
440  {
441  dict.lookup("trackLength") >> trackLength_;
442 
443  Info<< " fixed track length specified : "
444  << trackLength_ << nl << endl;
445  }
446 
447 
448  interpolationScheme_ = dict.lookupOrDefault
449  (
450  "interpolationScheme",
452  );
453 
454  //Info<< " using interpolation " << interpolationScheme_
455  // << endl;
456 
457  cloudName_ = dict.lookupOrDefault<word>("cloudName", "streamLine");
458  dict.lookup("seedSampleSet") >> seedSet_;
459 
460  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
461 
462  meshSearchPtr_.reset(new meshSearch(mesh));
463 
464  const dictionary& coeffsDict = dict.subDict(seedSet_ + "Coeffs");
465  sampledSetPtr_ = sampledSet::New
466  (
467  seedSet_,
468  mesh,
469  meshSearchPtr_(),
470  coeffsDict
471  );
472  coeffsDict.lookup("axis") >> sampledSetAxis_;
473 
474  scalarFormatterPtr_ = writer<scalar>::New(dict.lookup("setFormat"));
475  vectorFormatterPtr_ = writer<vector>::New(dict.lookup("setFormat"));
476  }
477 }
478 
479 
481 {
482 // const Time& runTime = obr_.time();
483 // Pout<< "**streamLine::execute : time:" << runTime.timeName() << endl;
484 //
485 // bool isOutputTime = false;
486 //
487 // const functionObjectList& fobs = runTime.functionObjects();
488 //
489 // forAll(fobs, i)
490 // {
491 // if (isA<streamLineFunctionObject>(fobs[i]))
492 // {
493 // const streamLineFunctionObject& fo =
494 // dynamic_cast<const streamLineFunctionObject&>(fobs[i]);
495 //
496 // if (fo.name() == name_)
497 // {
498 // Pout<< "found me:" << i << endl;
499 // if (fo.outputControl().output())
500 // {
501 // isOutputTime = true;
502 // break;
503 // }
504 // }
505 // }
506 // }
507 //
508 //
509 // if (active_ && isOutputTime)
510 // {
511 // track();
512 // }
513 }
514 
515 
517 {}
518 
519 
521 {}
522 
523 
525 {
526  if (active_)
527  {
528  Info<< type() << " " << name_ << " output:" << nl;
529 
530  const Time& runTime = obr_.time();
531  const fvMesh& mesh = dynamic_cast<const fvMesh&>(obr_);
532 
533 
534  // Do all injection and tracking
535  track();
536 
537 
538  if (Pstream::parRun())
539  {
540  // Append slave tracks to master ones
541  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
542 
543  globalIndex globalTrackIDs(allTracks_.size());
544 
545  // Construct a distribution map to pull all to the master.
546  labelListList sendMap(Pstream::nProcs());
547  labelListList recvMap(Pstream::nProcs());
548 
549  if (Pstream::master())
550  {
551  // Master: receive all. My own first, then consecutive
552  // processors.
553  label trackI = 0;
554 
555  forAll(recvMap, procI)
556  {
557  labelList& fromProc = recvMap[procI];
558  fromProc.setSize(globalTrackIDs.localSize(procI));
559  forAll(fromProc, i)
560  {
561  fromProc[i] = trackI++;
562  }
563  }
564  }
565 
566  labelList& toMaster = sendMap[0];
567  toMaster.setSize(globalTrackIDs.localSize());
568  forAll(toMaster, i)
569  {
570  toMaster[i] = i;
571  }
572 
573  const mapDistribute distMap
574  (
575  globalTrackIDs.size(),
576  sendMap.xfer(),
577  recvMap.xfer()
578  );
579 
580 
581  // Distribute the track positions. Note: use scheduled comms
582  // to prevent buffering.
584  (
586  distMap.schedule(),
587  distMap.constructSize(),
588  distMap.subMap(),
589  distMap.constructMap(),
590  allTracks_
591  );
592 
593  // Distribute the scalars
594  forAll(allScalars_, scalarI)
595  {
597  (
599  distMap.schedule(),
600  distMap.constructSize(),
601  distMap.subMap(),
602  distMap.constructMap(),
603  allScalars_[scalarI]
604  );
605  }
606  // Distribute the vectors
607  forAll(allVectors_, vectorI)
608  {
610  (
612  distMap.schedule(),
613  distMap.constructSize(),
614  distMap.subMap(),
615  distMap.constructMap(),
616  allVectors_[vectorI]
617  );
618  }
619  }
620 
621 
622  label n = 0;
623  forAll(allTracks_, trackI)
624  {
625  n += allTracks_[trackI].size();
626  }
627 
628  Info<< " Tracks:" << allTracks_.size() << nl
629  << " Total samples:" << n
630  << endl;
631 
632 
633  // Massage into form suitable for writers
634  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
635 
636  if (Pstream::master() && allTracks_.size())
637  {
638  // Make output directory
639 
640  fileName vtkPath
641  (
643  ? runTime.path()/".."/"postProcessing"/"sets"/name()
644  : runTime.path()/"postProcessing"/"sets"/name()
645  );
646  if (mesh.name() != fvMesh::defaultRegion)
647  {
648  vtkPath = vtkPath/mesh.name();
649  }
650  vtkPath = vtkPath/mesh.time().timeName();
651 
652  mkDir(vtkPath);
653 
654  // Convert track positions
655 
656  PtrList<coordSet> tracks(allTracks_.size());
657  forAll(allTracks_, trackI)
658  {
659  tracks.set
660  (
661  trackI,
662  new coordSet
663  (
664  "track" + Foam::name(trackI),
665  sampledSetAxis_ //"xyz"
666  )
667  );
668  tracks[trackI].transfer(allTracks_[trackI]);
669  }
670 
671  // Convert scalar values
672 
673  if (allScalars_.size() > 0)
674  {
675  List<List<scalarField> > scalarValues(allScalars_.size());
676 
677  forAll(allScalars_, scalarI)
678  {
679  DynamicList<scalarList>& allTrackVals =
680  allScalars_[scalarI];
681  scalarValues[scalarI].setSize(allTrackVals.size());
682 
683  forAll(allTrackVals, trackI)
684  {
685  scalarList& trackVals = allTrackVals[trackI];
686  scalarValues[scalarI][trackI].transfer(trackVals);
687  }
688  }
689 
690  fileName vtkFile
691  (
692  vtkPath
693  / scalarFormatterPtr_().getFileName
694  (
695  tracks[0],
696  scalarNames_
697  )
698  );
699 
700  Info<< " Writing data to " << vtkFile.path() << endl;
701 
702  scalarFormatterPtr_().write
703  (
704  true, // writeTracks
705  tracks,
706  scalarNames_,
707  scalarValues,
708  OFstream(vtkFile)()
709  );
710  }
711 
712  // Convert vector values
713 
714  if (allVectors_.size() > 0)
715  {
716  List<List<vectorField> > vectorValues(allVectors_.size());
717 
718  forAll(allVectors_, vectorI)
719  {
720  DynamicList<vectorList>& allTrackVals =
721  allVectors_[vectorI];
722  vectorValues[vectorI].setSize(allTrackVals.size());
723 
724  forAll(allTrackVals, trackI)
725  {
726  vectorList& trackVals = allTrackVals[trackI];
727  vectorValues[vectorI][trackI].transfer(trackVals);
728  }
729  }
730 
731  fileName vtkFile
732  (
733  vtkPath
734  / vectorFormatterPtr_().getFileName
735  (
736  tracks[0],
737  vectorNames_
738  )
739  );
740 
741  //Info<< " Writing vector data to " << vtkFile << endl;
742 
743  vectorFormatterPtr_().write
744  (
745  true, // writeTracks
746  tracks,
747  vectorNames_,
748  vectorValues,
749  OFstream(vtkFile)()
750  );
751  }
752  }
753  }
754 }
755 
756 
758 {
759  read(dict_);
760 }
761 
762 
764 {
765  // Moving mesh affects the search tree
766  read(dict_);
767 }
768 
769 
770 //void Foam::streamLine::readUpdate(const polyMesh::readUpdateState state)
771 //{
772 // if (state != UNCHANGED)
773 // {
774 // read(dict_);
775 // }
776 //}
777 
778 
779 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
Output to file stream.
Definition: OFstream.H:81
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:380
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
static autoPtr< writer > New(const word &writeFormat)
Return a reference to the selected writer.
Definition: writer.C:35
void transfer(List< T > &)
Transfer contents of the argument List into this.
Definition: DynamicListI.H:277
virtual void end()
Execute the averaging at the final time-loop, currently does nothing.
Definition: streamLine.C:516
static const char *const typeName
Definition: Field.H:94
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const word & name() const
Return reference to name.
Definition: fvMesh.H:257
virtual void read(const dictionary &)
Read the field average data.
Definition: streamLine.C:371
virtual void updateMesh(const mapPolyMesh &)
Update for changes of mesh.
Definition: streamLine.C:757
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
fileName path() const
Return path.
Definition: Time.H:281
static void distribute(const Pstream::commsTypes commsType, const List< labelPair > &schedule, const label constructSize, const labelListList &subMap, const labelListList &constructMap, List< T > &, const int tag=UPstream::msgType())
Distribute data. Note:schedule only used for Pstream::scheduled.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
Holds list of sampling positions.
Definition: coordSet.H:49
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
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:448
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Xfer< List< T > > xfer()
Transfer contents to the Xfer container.
Definition: ListI.H:90
Class containing processor-to-processor mapping information.
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
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
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
PrimitivePatch< face, IndirectList, const pointField & > indirectPrimitivePatch
Foam::indirectPrimitivePatch.
messageStream Info
fileName path() const
Return directory path name (part before last /)
Definition: fileName.C:287
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:638
Namespace for OpenFOAM.
#define IOWarningIn(functionName, ios)
Report an IO warning using Foam::Warning.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
label n
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:386
static const char nl
Definition: Ostream.H:260
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
IOerror FatalIOError
#define WarningIn(functionName)
Report a warning using Foam::Warning.
wordList ReadFields(const Mesh &mesh, const IOobjectList &objects, PtrList< GeoField > &fields, const bool syncPar=true)
Helper routine to read fields.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
#define forAll(list, i)
Definition: UList.H:421
virtual Ostream & write(const token &)=0
Write next token to stream.
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
Definition: meshSearch.H:57
void setSize(const label)
Alter the addressed list size.
Definition: DynamicListI.H:167
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
Helper routine to read fields.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual ~streamLine()
Destructor.
Definition: streamLine.C:365
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
virtual void execute()
Execute the averaging.
Definition: streamLine.C:480
virtual void timeSet()
Called when time was set at the end of the Time::operator++.
Definition: streamLine.C:520
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Registry of regIOobjects.
List< label > labelList
A List of labels.
Definition: labelList.H:56
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
Definition: streamLine.C:763
A class for handling file names.
Definition: fileName.H:69
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
virtual void write()
Calculate the field average data and write.
Definition: streamLine.C:524
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:398
bool read(const char *, int32_t &)
Definition: int32IO.C:87
const Time & time() const
Return time.
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:420
static autoPtr< interpolation< scalar > > New(const word &interpolationType, const GeometricField< scalar, fvPatchField, volMesh > &psi)
Return a reference to the specified interpolation scheme.
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
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)