streamlines.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2023 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 "streamlines.H"
29 #include "streamlinesCloud.H"
30 #include "ReadFields.H"
31 #include "meshSearch.H"
32 #include "sampledSet.H"
33 #include "globalIndex.H"
34 #include "distributionMap.H"
35 #include "interpolationCellPoint.H"
36 #include "PatchTools.H"
37 #include "writeFile.H"
38 #include "polyTopoChangeMap.H"
39 #include "polyMeshMap.H"
40 #include "polyDistributionMap.H"
41 #include "OSspecific.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 
49 template<class Type>
51 {
52  List<List<Type>> gatheredField(Pstream::nProcs());
53  gatheredField[Pstream::myProcNo()] = field;
54  Pstream::gatherList(gatheredField);
55 
56  field =
57  ListListOps::combine<List<Type>>
58  (
59  gatheredField,
61  );
62 }
63 
64 }
65 
66 
67 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
68 
69 namespace Foam
70 {
71  template<>
72  const char*
74  {"forward", "backward", "both"};
75 
76  namespace functionObjects
77  {
80 
83  }
84 }
85 
86 
87 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
88 
90 (
91  const word& name,
92  const Time& runTime,
93  const dictionary& dict
94 )
95 :
96  fvMeshFunctionObject(name, runTime, dict),
97  dict_(dict),
98  nSubCycle_(0)
99 {
100  read(dict_);
101 }
102 
103 
104 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
105 
107 {}
108 
109 
110 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
111 
113 {
114  if (dict != dict_)
115  {
116  dict_ = dict;
117  }
118 
119  Info<< type() << " " << name() << ":" << nl;
120 
121  dict.lookup("fields") >> fields_;
122 
123  UName_ = dict.lookupOrDefault("U", word("U"));
124 
125  writeAge_ = dict.lookupOrDefault<Switch>("writeAge", true);
126 
127  trackDirection_ = trackDirectionNames_[word(dict.lookup("direction"))];
128 
129  trackOutside_ = dict.lookupOrDefault<Switch>("outside", false);
130 
131  dict.lookup("lifeTime") >> lifeTime_;
132 
133  if (lifeTime_ < 1)
134  {
136  << "Illegal value " << lifeTime_ << " for lifeTime"
137  << exit(FatalError);
138  }
139 
140  bool subCycling = dict.found("nSubCycle");
141  bool fixedLength = dict.found("trackLength");
142  if (subCycling && fixedLength)
143  {
145  << "Cannot both specify automatic time stepping (through '"
146  << "nSubCycle' specification) and fixed track length (through '"
147  << "trackLength')"
148  << exit(FatalIOError);
149  }
150  if (subCycling)
151  {
152  nSubCycle_ = max(dict.lookup<scalar>("nSubCycle"), 1);
153  trackLength_ = vGreat;
154  Info<< " automatic track length specified through"
155  << " number of sub cycles : " << nSubCycle_ << nl << endl;
156  }
157  else
158  {
159  nSubCycle_ = 1;
160  dict.lookup("trackLength") >> trackLength_;
161  Info<< " fixed track length specified : "
162  << trackLength_ << nl << endl;
163  }
164 
165  interpolationScheme_ =
166  dict.lookupOrDefault
167  (
168  "interpolationScheme",
170  );
171 
172  cloudName_ = dict.lookupOrDefault<word>("cloudName", "streamlines");
173 
174  meshSearchPtr_.reset(new meshSearch(mesh_));
175 
176  sampledSetPtr_ = sampledSet::New
177  (
178  "seedSampleSet",
179  mesh_,
180  meshSearchPtr_(),
181  dict.subDict("seedSampleSet")
182  );
183 
184  formatterPtr_ = setWriter::New(dict.lookup("setFormat"), dict);
185 
186  return true;
187 }
188 
189 
191 {
192  wordList allFields(fields_);
193  allFields.append(UName_);
194 
195  return allFields;
196 }
197 
198 
200 {
201  return true;
202 }
203 
204 
206 {
207  Info<< type() << " " << name() << " write:" << nl;
208 
209  // Create list of available fields
211  forAll(fields_, fieldi)
212  {
213  if
214  (
215  false
216  #define FoundTypeField(Type, nullArg) \
217  || foundObject<VolField<Type>>(fields_[fieldi])
219  #undef FoundTypeField
220  )
221  {
222  fieldNames.append(fields_[fieldi]);
223  }
224  else
225  {
226  cannotFindObject(fields_[fieldi]);
227  }
228  }
229 
230  // Lookup fields and construct interpolators
231  #define DeclareTypeInterpolator(Type, nullArg) \
232  PtrList<interpolation<Type>> Type##Interp(fieldNames.size());
234  #undef DeclareTypeInterpolator
235  forAll(fieldNames, fieldi)
236  {
237  #define ConstructTypeInterpolator(Type, nullArg) \
238  if (mesh_.foundObject<VolField<Type>>(fieldNames[fieldi])) \
239  { \
240  Type##Interp.set \
241  ( \
242  fieldi, \
243  interpolation<Type>::New \
244  ( \
245  interpolationScheme_, \
246  mesh_.lookupObject<VolField<Type>>(fieldNames[fieldi]) \
247  ) \
248  ); \
249  }
251  #undef ConstructTypeInterpolator
252  }
253 
254  // Create a velocity interpolator if it is not already available
255  const label UIndex = findIndex(fieldNames, UName_);
256  tmpNrc<interpolation<vector>> UInterp(nullptr);
257  if (UIndex == -1)
258  {
259  UInterp =
261  (
263  (
264  interpolationScheme_,
265  mesh_.lookupObject<volVectorField>(UName_)
266  ).ptr()
267  );
268  }
269 
270  // Do tracking to create sampled data
271  DynamicField<point> allPositions;
272  DynamicField<label> allTracks;
273  DynamicField<label> allTrackParts;
274  DynamicField<scalar> allAges;
275  #define DeclareAllTypes(Type, nullArg) \
276  List<DynamicField<Type>> all##Type##s(fieldNames.size());
278  #undef DeclareAllTypes
279  {
280  // Create a cloud and initialise with points from the sampled set
281  globalIndex gi(sampledSetPtr_().size());
283  (
284  mesh_,
285  cloudName_,
287  );
288  label nLocateBoundaryHits;
289  forAll(sampledSetPtr_(), i)
290  {
291  particles.addParticle
292  (
294  (
295  mesh_,
296  sampledSetPtr_().positions()[i],
297  sampledSetPtr_().cells()[i],
298  nLocateBoundaryHits,
299  lifeTime_,
300  gi.toGlobal(i)
301  )
302  );
303  }
304 
305  // Report the number of successful seeds
306  const label nSeeds = returnReduce(particles.size(), sumOp<label>());
307  Info << " Seeded " << nSeeds << " particles" << endl;
308 
309  // Create tracking data
311  (
312  particles,
313  #define TypeInterpolatorParameter(Type, nullArg) \
314  Type##Interp,
317  UIndex != -1 ? vectorInterp[UIndex] : UInterp(),
318  trackDirection_ == trackDirection::forward,
319  trackOutside_,
320  nSubCycle_,
321  trackLength_,
322  allPositions,
323  allTracks,
324  allTrackParts,
325  allAges
326  #define AllTypesParameter(Type, nullArg) \
327  , all##Type##s
329  #undef AllTypesParameter
330  );
331 
332  // Track
333  IDLList<streamlinesParticle> initialParticles;
334  if (trackDirection_ == trackDirection::both)
335  {
336  initialParticles = particles;
337  }
338 
339  particles.move(particles, td);
340 
341  if (trackDirection_ == trackDirection::both)
342  {
343  particles.IDLList<streamlinesParticle>::operator=(initialParticles);
344  td.trackForward_ = !td.trackForward_;
345  particles.move(particles, td);
346  }
347  }
348 
349  // Gather data on the master
350  if (Pstream::parRun())
351  {
352  gatherAndFlatten(allPositions);
353  gatherAndFlatten(allTracks);
354  gatherAndFlatten(allTrackParts);
355  gatherAndFlatten(allAges);
356  forAll(fieldNames, fieldi)
357  {
358  #define GatherAndFlattenAllTypes(Type, nullArg) \
359  if (Type##Interp.set(fieldi)) \
360  { \
361  gatherAndFlatten(all##Type##s[fieldi]); \
362  }
364  #undef GatherAndFlattenAllTypes
365  }
366  }
367 
368  // Report the total number of samples
369  Info<< " Sampled " << allPositions.size() << " locations" << endl;
370 
371  // Bin-sort by track and trackPart to build an ordering
372  labelList order(allPositions.size());
373  if (Pstream::master() && allPositions.size())
374  {
375  const label nTracks = max(allTracks) + 1;
376  const label trackParti0 = min(allTrackParts);
377  const label trackParti1 = max(allTrackParts) + 1;
378 
379  labelListList trackPartCounts
380  (
381  nTracks,
382  labelList(trackParti1 - trackParti0, 0)
383  );
384  forAll(allPositions, samplei)
385  {
386  const label tracki = allTracks[samplei];
387  const label trackParti = -trackParti0 + allTrackParts[samplei];
388  trackPartCounts[tracki][trackParti] ++;
389  }
390 
391  label offset = 0;
392  labelListList trackPartOffsets
393  (
394  nTracks,
395  labelList(trackParti1 - trackParti0, 0)
396  );
397  forAll(trackPartOffsets, tracki)
398  {
399  forAll(trackPartOffsets[tracki], trackParti)
400  {
401  trackPartOffsets[tracki][trackParti] += offset;
402  offset += trackPartCounts[tracki][trackParti];
403  }
404  }
405 
406  forAll(trackPartCounts, tracki)
407  {
408  trackPartCounts[tracki] = 0;
409  }
410 
411  forAll(allPositions, samplei)
412  {
413  const label tracki = allTracks[samplei];
414  const label trackParti = -trackParti0 + allTrackParts[samplei];
415 
416  order[samplei] =
417  trackPartOffsets[tracki][trackParti]
418  + trackPartCounts[tracki][trackParti];
419 
420  trackPartCounts[tracki][trackParti] ++;
421  }
422  }
423 
424  //auto reportTrackParts = [&]()
425  //{
426  // Info<< nl;
427  // forAll(allPositions, samplei)
428  // {
429  // if
430  // (
431  // samplei == 0
432  // || allTracks[samplei] != allTracks[samplei - 1]
433  // || allTrackParts[samplei] != allTrackParts[samplei - 1]
434  // )
435  // {
436  // Info<< "track #" << allTracks[samplei]
437  // << " part #" << allTrackParts[samplei]
438  // << " from i=" << samplei << " to ";
439  // }
440  // if
441  // (
442  // samplei == allPositions.size() - 1
443  // || allTracks[samplei + 1] != allTracks[samplei]
444  // || allTrackParts[samplei + 1] != allTrackParts[samplei]
445  // )
446  // {
447  // Info<< "i=" << samplei << nl;
448  // }
449  // }
450  //};
451 
452  //reportTrackParts();
453 
454  // Reorder
455  if (Pstream::master())
456  {
457  allPositions.rmap(allPositions, order);
458  allTracks.rmap(allTracks, order);
459  allTrackParts.rmap(allTrackParts, order);
460  allAges.rmap(allAges, order);
461  forAll(fieldNames, fieldi)
462  {
463  #define RMapAllTypes(Type, nullArg) \
464  if (Type##Interp.set(fieldi)) \
465  { \
466  all##Type##s[fieldi].rmap(all##Type##s[fieldi], order); \
467  }
469  #undef RMapAllTypes
470  }
471  }
472 
473  //reportTrackParts();
474 
475  // Relabel tracks and track parts into track labels only, and join the
476  // forward and backward track parts that are connected to the seed
477  if (Pstream::master())
478  {
479  label samplei = 0, tracki = 0;
480  forAll(allPositions, samplej)
481  {
482  const label trackj = allTracks[samplej];
483  const label trackPartj = allTrackParts[samplej];
484 
485  allPositions[samplei] = allPositions[samplej];
486  allTracks[samplei] = tracki;
487  allTrackParts[samplei] = 0;
488  allAges[samplei] = allAges[samplej];
489  forAll(fieldNames, fieldi)
490  {
491  #define ShuffleUpAllTypes(Type, nullArg) \
492  if (Type##Interp.set(fieldi)) \
493  { \
494  all##Type##s[fieldi][samplei] = \
495  all##Type##s[fieldi][samplej]; \
496  }
498  #undef ShuffleUpAllTypes
499  }
500 
501  const bool joinNewParts =
502  samplej != allPositions.size() - 1
503  && trackPartj == -1
504  && allTrackParts[samplej + 1] == 0;
505 
506  if (!joinNewParts) samplei ++;
507 
508  const bool newPart =
509  samplej == allPositions.size() - 1
510  || trackj != allTracks[samplej + 1]
511  || trackPartj != allTrackParts[samplej + 1];
512 
513  if (!joinNewParts && newPart) tracki ++;
514  }
515 
516  allPositions.resize(samplei);
517  allTracks.resize(samplei);
518  allTrackParts.resize(samplei);
519  allAges.resize(samplei);
520  forAll(fieldNames, fieldi)
521  {
522  #define ResizeAllTypes(Type, nullArg) \
523  if (Type##Interp.set(fieldi)) \
524  { \
525  all##Type##s[fieldi].resize(samplei); \
526  }
528  #undef ResizeAllTypes
529  }
530  }
531 
532  //reportTrackParts();
533 
534  // Write
535  if (Pstream::master() && allPositions.size())
536  {
537  // Make output directory
538  const fileName outputPath =
539  time_.globalPath()
541  /(mesh_.name() != polyMesh::defaultRegion ? mesh_.name() : word())
542  /name()
543  /time_.name();
544  mkDir(outputPath);
545 
546  // Pass data to the formatter to write
547  const label nValueSets = fieldNames.size() + writeAge_;
548  wordList valueSetNames(nValueSets);
549  #define DeclareTypeValueSets(Type, nullArg) \
550  UPtrList<const Field<Type>> Type##ValueSets(nValueSets);
552  #undef DeclareTypeValueSets
553  if (writeAge_)
554  {
555  valueSetNames[0] = "age";
556  scalarValueSets.set(0, &allAges);
557  }
558  forAll(fieldNames, fieldi)
559  {
560  valueSetNames[fieldi + writeAge_] = fieldNames[fieldi];
561 
562  #define SetTypeValueSetPtr(Type, nullArg) \
563  if (Type##Interp.set(fieldi)) \
564  { \
565  Type##ValueSets.set \
566  ( \
567  fieldi + writeAge_, \
568  &all##Type##s[fieldi] \
569  ); \
570  }
572  #undef SetTypeValueSetPtr
573  }
574  formatterPtr_->write
575  (
576  outputPath,
577  "tracks",
578  coordSet(allTracks, word::null, allPositions),
579  valueSetNames
580  #define TypeValueSetsParameter(Type, nullArg) , Type##ValueSets
583  );
584  }
585 
586  return true;
587 }
588 
589 
591 {
592  if (&mesh == &mesh_)
593  {
594  // Moving mesh affects the search tree
595  read(dict_);
596  }
597 }
598 
599 
601 (
602  const polyTopoChangeMap& map
603 )
604 {
605  if (&map.mesh() == &mesh_)
606  {
607  read(dict_);
608  }
609 }
610 
611 
613 (
614  const polyMeshMap& map
615 )
616 {
617  if (&map.mesh() == &mesh_)
618  {
619  read(dict_);
620  }
621 }
622 
623 
625 (
626  const polyDistributionMap& map
627 )
628 {
629  if (&map.mesh() == &mesh_)
630  {
631  read(dict_);
632  }
633 }
634 
635 
636 // ************************************************************************* //
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
Field reading functions for post-processing utilities.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Dynamically sized Field.
Definition: DynamicField.H:72
void resize(const label)
Alter the addressed list size.
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:386
Generic GeometricField class.
Template class for intrusive linked lists.
Definition: ILList.H:67
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:178
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
Holds list of sampling positions.
Definition: coordSet.H:51
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
A class for handling file names.
Definition: fileName.H:82
Abstract base-class for Time/database functionObjects.
Specialisation of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
This functionObject tracks a particle cloud in the specified velocity field of an incompressible flow...
Definition: particles.H:111
Generates streamline data by sampling a set of user-specified fields along a particle track,...
Definition: streamlines.H:216
streamlines(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: streamlines.C:90
virtual wordList fields() const
Return the list of fields required.
Definition: streamlines.C:190
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
Definition: streamlines.C:601
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
Definition: streamlines.C:625
virtual ~streamlines()
Destructor.
Definition: streamlines.C:106
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
Definition: streamlines.C:613
virtual void movePoints(const polyMesh &)
Update for mesh point-motion.
Definition: streamlines.C:590
static const NamedEnum< trackDirection, 3 > trackDirectionNames_
Track direction enumeration names.
Definition: streamlines.H:230
virtual bool execute()
Do nothing.
Definition: streamlines.C:199
virtual bool write()
Calculate and write the streamlines.
Definition: streamlines.C:205
virtual bool read(const dictionary &)
Read the field average data.
Definition: streamlines.C:112
static const word outputPrefix
Directory prefix.
Definition: writeFile.H:72
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:64
label toGlobal(const label i) const
From local to global.
Definition: globalIndexI.H:82
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
static autoPtr< interpolation< Type > > New(const word &interpolationType, const VolField< Type > &psi)
Return a reference to the specified interpolation scheme.
Definition: interpolation.C:56
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search.
Definition: meshSearch.H:58
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
const polyMesh & mesh() const
Return polyMesh.
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
const polyMesh & mesh() const
Return polyMesh.
Definition: polyMeshMap.H:75
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:267
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const polyMesh & mesh() const
Return polyMesh.
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:190
static autoPtr< setWriter > New(const word &writeType, const IOstream::streamFormat writeFormat=IOstream::ASCII, const IOstream::compressionType writeCompression=IOstream::UNCOMPRESSED)
Select given write options.
Definition: setWriter.C:286
A Cloud of streamlines particles.
Particle class that samples fields as it passes through. Used in streamlines calculation.
A class for managing temporary objects without reference counting.
Definition: tmpNrc.H:53
A class for handling words, derived from string.
Definition: word.H:62
static const word null
An empty word.
Definition: word.H:77
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:346
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
static List< word > fieldNames
Definition: globalFoam.H:46
const cellShapeList & cells
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){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
bool read(const char *, int32_t &)
Definition: int32IO.C:85
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:290
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
int order(const scalar s)
messageStream Info
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
IOerror FatalIOError
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurrence of given element and return index,.
error FatalError
FOR_ALL_FIELD_TYPES(makeFieldSourceTypedef)
void offset(label &lst, const label o)
void gatherAndFlatten(DynamicField< Type > &field)
Definition: streamlines.C:50
static const char nl
Definition: Ostream.H:266
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict
#define FoundTypeField(Type, nullArg)
#define DeclareTypeValueSets(Type, nullArg)
#define SetTypeValueSetPtr(Type, nullArg)
#define DeclareAllTypes(Type, nullArg)
#define AllTypesParameter(Type, nullArg)
#define ResizeAllTypes(Type, nullArg)
#define TypeValueSetsParameter(Type, nullArg)
#define DeclareTypeInterpolator(Type, nullArg)
#define GatherAndFlattenAllTypes(Type, nullArg)
#define ConstructTypeInterpolator(Type, nullArg)
#define RMapAllTypes(Type, nullArg)
#define TypeInterpolatorParameter(Type, nullArg)
#define ShuffleUpAllTypes(Type, nullArg)