VTK  9.0.1
vtkParticleTracerBase.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtkParticleTracerBase.h
5 
6  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7  All rights reserved.
8  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9 
10  This software is distributed WITHOUT ANY WARRANTY; without even
11  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12  PURPOSE. See the above copyright notice for more information.
13 
14 =========================================================================*/
28 #ifndef vtkParticleTracerBase_h
29 #define vtkParticleTracerBase_h
30 
31 #include "vtkFiltersFlowPathsModule.h" // For export macro
32 #include "vtkPolyDataAlgorithm.h"
33 #include "vtkSmartPointer.h" // For protected ivars.
34 
35 #include <list> // STL Header
36 #include <vector> // STL Header
37 
40 class vtkCellArray;
41 class vtkCharArray;
43 class vtkDataArray;
44 class vtkDataSet;
45 class vtkDoubleArray;
46 class vtkFloatArray;
47 class vtkGenericCell;
49 class vtkIntArray;
52 class vtkPointData;
53 class vtkPoints;
54 class vtkPolyData;
56 
58 {
59 typedef struct
60 {
61  double x[4];
62 } Position;
63 typedef struct
64 {
65  // These are used during iteration
67  int CachedDataSetId[2];
68  vtkIdType CachedCellId[2];
70  // These are computed scalars we might display
71  int SourceID;
72  int TimeStepAge; // amount of time steps the particle has advanced
74  int InjectedStepId; // time step the particle was injected
77  // These are useful to track for debugging etc
78  int ErrorCode;
79  float age;
80  // these are needed across time steps to compute vorticity
81  float rotation;
82  float angularVel;
83  float time;
84  float speed;
85  // once the partice is added, PointId is valid and is the tuple location
86  // in ProtoPD.
88  // if PointId is negative then in parallel this particle was just
89  // received and we need to get the tuple value from vtkPParticleTracerBase::Tail.
92 
93 typedef std::vector<ParticleInformation> ParticleVector;
94 typedef ParticleVector::iterator ParticleIterator;
95 typedef std::list<ParticleInformation> ParticleDataList;
96 typedef ParticleDataList::iterator ParticleListIterator;
97 };
98 
99 class VTKFILTERSFLOWPATHS_EXPORT vtkParticleTracerBase : public vtkPolyDataAlgorithm
100 {
101 public:
102  enum Solvers
103  {
108  UNKNOWN
109  };
110 
112  void PrintSelf(ostream& os, vtkIndent indent) override;
114 
116 
121  vtkGetMacro(ComputeVorticity, bool);
124 
126 
129  vtkGetMacro(TerminalSpeed, double);
130  void SetTerminalSpeed(double);
132 
134 
138  vtkGetMacro(RotationScale, double);
139  void SetRotationScale(double);
141 
143 
147  vtkSetMacro(IgnorePipelineTime, vtkTypeBool);
148  vtkGetMacro(IgnorePipelineTime, vtkTypeBool);
149  vtkBooleanMacro(IgnorePipelineTime, vtkTypeBool);
151 
153 
162  vtkGetMacro(ForceReinjectionEveryNSteps, int);
165 
167 
173  void SetTerminationTime(double t);
174  vtkGetMacro(TerminationTime, double);
176 
178  vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
179 
182 
184 
188  vtkGetMacro(StartTime, double);
189  void SetStartTime(double t);
191 
193 
202  vtkSetMacro(StaticSeeds, int);
203  vtkGetMacro(StaticSeeds, int);
205 
207 
216  vtkSetMacro(StaticMesh, int);
217  vtkGetMacro(StaticMesh, int);
219 
221 
228  vtkGetObjectMacro(ParticleWriter, vtkAbstractParticleWriter);
230 
232 
236  vtkSetStringMacro(ParticleFileName);
237  vtkGetStringMacro(ParticleFileName);
239 
241 
245  vtkSetMacro(EnableParticleWriting, vtkTypeBool);
246  vtkGetMacro(EnableParticleWriting, vtkTypeBool);
247  vtkBooleanMacro(EnableParticleWriting, vtkTypeBool);
249 
251 
256  vtkSetMacro(DisableResetCache, vtkTypeBool);
257  vtkGetMacro(DisableResetCache, vtkTypeBool);
258  vtkBooleanMacro(DisableResetCache, vtkTypeBool);
260 
262 
268 
269 protected:
270  vtkSmartPointer<vtkPolyData> Output; // managed by child classes
272 
277  vtkIdType UniqueIdCounter; // global Id counter used to give particles a stamp
279  vtkSmartPointer<vtkPointData> ParticlePointData; // the current particle point data consistent
280  // with particle history
281  // Everything related to time
282  vtkTypeBool IgnorePipelineTime; // whether to use the pipeline time for termination
283  vtkTypeBool DisableResetCache; // whether to enable ResetCache() method
285 
288 
289  //
290  // Make sure the pipeline knows what type we expect as input
291  //
293 
294  //
295  // The usual suspects
296  //
298  vtkInformationVector* outputVector) override;
299 
300  //
301  // Store any information we need in the output and fetch what we can
302  // from the input
303  //
305  vtkInformationVector* outputVector) override;
306 
307  //
308  // Compute input time steps given the output step
309  //
311  vtkInformationVector* outputVector) override;
312 
313  //
314  // what the pipeline calls for each time step
315  //
316  int RequestData(vtkInformation* request, vtkInformationVector** inputVector,
317  vtkInformationVector* outputVector) override;
318 
319  //
320  // these routines are internally called to actually generate the output
321  //
322  virtual int ProcessInput(vtkInformationVector** inputVector);
323 
324  // This is the main part of the algorithm:
325  // * move all the particles one step
326  // * Reinject particles (by adding them to this->ParticleHistories)
327  // either at the beginning or at the end of each step (modulo
328  // this->ForceReinjectionEveryNSteps)
329  // * Output a polydata representing the moved particles
330  // Note that if the starting and the ending time coincide, the polydata is still valid.
331  virtual vtkPolyData* Execute(vtkInformationVector** inputVector);
332 
333  // the RequestData will call these methods in turn
334  virtual void Initialize() {} // the first iteration
335  virtual int OutputParticles(vtkPolyData* poly) = 0; // every iteration
336  virtual void Finalize() {} // the last iteration
337 
342  virtual std::vector<vtkDataSet*> GetSeedSources(vtkInformationVector* inputVector, int timeStep);
343 
344  //
345  // Initialization of input (vector-field) geometry
346  //
349 
356 
358  vtkParticleTracerBaseNamespace::ParticleVector& candidates, std::vector<int>& passed);
359 
366  virtual void AssignSeedsToProcessors(double time, vtkDataSet* source, int sourceID, int ptId,
367  vtkParticleTracerBaseNamespace::ParticleVector& localSeedPoints, int& localAssignedCount);
368 
374 
380 
386  virtual bool UpdateParticleListFromOtherProcesses() { return false; }
387 
392  double currenttime, double terminationtime, vtkInitialValueProblemSolver* integrator);
393 
394  // if the particle is added to send list, then returns value is 1,
395  // if it is kept on this process after a retry return value is 0
398  {
399  return true;
400  }
401 
409  double pos[4], double p2[4], double intersection[4], vtkGenericCell* cell);
410 
411  //
412  // Scalar arrays that are generated as each particle is updated
413  //
415 
425 
426  // utility function we use to test if a point is inside any of our local datasets
427  bool InsideBounds(double point[]);
428 
430  vtkGenericCell* cell, double pcoords[3], vtkDoubleArray* cellVectors, double vorticity[3]);
431 
432  //------------------------------------------------------
433 
434  double GetCacheDataTime(int i);
436 
437  virtual void ResetCache();
439 
441 
446  virtual bool IsPointDataValid(vtkDataObject* input);
447  bool IsPointDataValid(vtkCompositeDataSet* input, std::vector<std::string>& arrayNames);
448  void GetPointDataArrayNames(vtkDataSet* input, std::vector<std::string>& names);
450 
451  vtkGetMacro(ReinjectionCounter, int);
452  vtkGetMacro(CurrentTimeValue, double);
453 
458  virtual void InitializeExtraPointDataArrays(vtkPointData* vtkNotUsed(outputPD)) {}
459 
461 
463 
468  virtual void AddRestartSeeds(vtkInformationVector** /*inputVector*/) {}
469 
470 private:
474  void SetInterpolatorPrototype(vtkAbstractInterpolatedVelocityField*) {}
475 
485  bool RetryWithPush(vtkParticleTracerBaseNamespace::ParticleInformation& info, double* point1,
486  double delT, int subSteps);
487 
488  bool SetTerminationTimeNoModify(double t);
489 
490  // Parameters of tracing
491  vtkInitialValueProblemSolver* Integrator;
492  double IntegrationStep;
493  double MaximumError;
494  bool ComputeVorticity;
495  double RotationScale;
496  double TerminalSpeed;
497 
498  // A counter to keep track of how many times we reinjected
499  int ReinjectionCounter;
500 
501  // Important for Caching of Cells/Ids/Weights etc
502  int AllFixedGeometry;
503  int StaticMesh;
504  int StaticSeeds;
505 
506  std::vector<double> InputTimeValues;
507  double StartTime;
508  double TerminationTime;
509  double CurrentTimeValue;
510 
511  int StartTimeStep; // InputTimeValues[StartTimeStep] <= StartTime <=
512  // InputTimeValues[StartTimeStep+1]
513  int CurrentTimeStep;
514  int TerminationTimeStep; // computed from start time
515  bool FirstIteration;
516 
517  // Innjection parameters
518  int ForceReinjectionEveryNSteps;
519  vtkTimeStamp ParticleInjectionTime;
520  bool HasCache;
521 
522  // Particle writing to disk
523  vtkAbstractParticleWriter* ParticleWriter;
524  char* ParticleFileName;
525  vtkTypeBool EnableParticleWriting;
526 
527  // The main lists which are held during operation- between time step updates
529 
530  // The velocity interpolator
532  vtkAbstractInterpolatedVelocityField* InterpolatorPrototype;
533 
534  // Data for time step CurrentTimeStep-1 and CurrentTimeStep
536 
537  // Cache bounds info for each dataset we will use repeatedly
538  typedef struct
539  {
540  double b[6];
541  } bounds;
542  std::vector<bounds> CachedBounds[2];
543 
544  // temporary variables used by Exeucte(), for convenience only
545 
546  vtkSmartPointer<vtkPoints> OutputCoordinates;
547  vtkSmartPointer<vtkFloatArray> ParticleAge;
548  vtkSmartPointer<vtkIntArray> ParticleIds;
549  vtkSmartPointer<vtkCharArray> ParticleSourceIds;
550  vtkSmartPointer<vtkIntArray> InjectedPointIds;
551  vtkSmartPointer<vtkIntArray> InjectedStepIds;
552  vtkSmartPointer<vtkIntArray> ErrorCodeArray;
553  vtkSmartPointer<vtkFloatArray> ParticleVorticity;
554  vtkSmartPointer<vtkFloatArray> ParticleRotation;
555  vtkSmartPointer<vtkFloatArray> ParticleAngularVel;
557  vtkSmartPointer<vtkPointData> OutputPointData;
558  vtkSmartPointer<vtkDataSet> DataReferenceT[2];
559  vtkSmartPointer<vtkCellArray> ParticleCells;
560 
562  void operator=(const vtkParticleTracerBase&) = delete;
563  vtkTimeStamp ExecuteTime;
564 
565  unsigned int NumberOfParticles();
566 
569 
570  static const double Epsilon;
571 };
572 
573 #endif
An abstract class for obtaining the interpolated velocity values at a point.
abstract class to write particle data to file
Proxy object to connect input/output ports.
object to represent cell connectivity
Definition: vtkCellArray.h:180
dynamic, self-adjusting array of char
Definition: vtkCharArray.h:36
abstract superclass for composite (multi-block or AMR) datasets
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:50
general representation of visualization data
Definition: vtkDataObject.h:60
abstract class to specify dataset behavior
Definition: vtkDataSet.h:57
dynamic, self-adjusting array of double
dynamic, self-adjusting array of float
Definition: vtkFloatArray.h:36
provides thread-safe access to cells
a simple class to control print indentation
Definition: vtkIndent.h:34
Store zero or more vtkInformation instances.
Store vtkAlgorithm input/output information.
Integrate a set of ordinary differential equations (initial value problem) in time.
dynamic, self-adjusting array of int
Definition: vtkIntArray.h:40
Composite dataset that organizes datasets into blocks.
Multiprocessing communication superclass.
A particle tracer for vector fields.
vtkSmartPointer< vtkPointData > ProtoPD
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
vtkCharArray * GetParticleSourceIds(vtkPointData *)
void IntegrateParticle(vtkParticleTracerBaseNamespace::ParticleListIterator &it, double currenttime, double terminationtime, vtkInitialValueProblemSolver *integrator)
particle between the two times supplied.
void SetComputeVorticity(bool)
vtkSmartPointer< vtkPointData > ParticlePointData
void TestParticles(vtkParticleTracerBaseNamespace::ParticleVector &candidates, std::vector< int > &passed)
vtkIntArray * GetErrorCodeArr(vtkPointData *)
virtual void AddRestartSeeds(vtkInformationVector **)
For restarts of particle paths, we add in the ability to add in particles from a previous computation...
bool IsPointDataValid(vtkCompositeDataSet *input, std::vector< std::string > &arrayNames)
void SetTerminationTime(double t)
Setting TerminationTime to a positive value will cause particles to terminate when the time is reache...
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
This is called by the superclass.
virtual bool UpdateParticleListFromOtherProcesses()
this is used during classification of seed points and also between iterations of the main loop as par...
vtkFloatArray * GetParticleVorticity(vtkPointData *)
vtkFloatArray * GetParticleRotation(vtkPointData *)
void UpdateParticleList(vtkParticleTracerBaseNamespace::ParticleVector &candidates)
and sending between processors, into a list, which is used as the master list on this processor
virtual void AppendToExtraPointDataArrays(vtkParticleTracerBaseNamespace::ParticleInformation &)
double GetCacheDataTime(int i)
void TestParticles(vtkParticleTracerBaseNamespace::ParticleVector &candidates, vtkParticleTracerBaseNamespace::ParticleVector &passed, int &count)
inside our data.
vtkTypeBool ProcessRequest(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
see vtkAlgorithm for details
void AddSourceConnection(vtkAlgorithmOutput *input)
Provide support for multiple seed sources.
void CalculateVorticity(vtkGenericCell *cell, double pcoords[3], vtkDoubleArray *cellVectors, double vorticity[3])
void CreateProtoPD(vtkDataObject *input)
virtual int ProcessInput(vtkInformationVector **inputVector)
vtkTemporalInterpolatedVelocityField * GetInterpolator()
vtkIntArray * GetInjectedStepIds(vtkPointData *)
virtual bool SendParticleToAnotherProcess(vtkParticleTracerBaseNamespace::ParticleInformation &, vtkParticleTracerBaseNamespace::ParticleInformation &, vtkPointData *)
void SetForceReinjectionEveryNSteps(int)
int UpdateDataCache(vtkDataObject *td)
void SetTerminalSpeed(double)
virtual void SetParticleWriter(vtkAbstractParticleWriter *pw)
Set/Get the Writer associated with this Particle Tracer Ideally a parallel IO capable vtkH5PartWriter...
void SetRotationScale(double)
void AddParticle(vtkParticleTracerBaseNamespace::ParticleInformation &info, double *velocity)
vtkParticleTracerBaseNamespace::ParticleDataList ParticleHistories
void GetPointDataArrayNames(vtkDataSet *input, std::vector< std::string > &names)
virtual vtkPolyData * Execute(vtkInformationVector **inputVector)
int RequestUpdateExtent(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
This is called by the superclass.
void SetStartTime(double t)
virtual void AssignSeedsToProcessors(double time, vtkDataSet *source, int sourceID, int ptId, vtkParticleTracerBaseNamespace::ParticleVector &localSeedPoints, int &localAssignedCount)
all the injection/seed points according to which processor they belong to.
virtual std::vector< vtkDataSet * > GetSeedSources(vtkInformationVector *inputVector, int timeStep)
Method to get the data set seed sources.
bool InsideBounds(double point[])
virtual bool IsPointDataValid(vtkDataObject *input)
Methods that check that the input arrays are ordered the same on all data sets.
int RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
vtkFloatArray * GetParticleAge(vtkPointData *)
vtkIntArray * GetInjectedPointIds(vtkPointData *)
void SetIntegratorType(int type)
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
vtkFloatArray * GetParticleAngularVel(vtkPointData *)
int FillInputPortInformation(int port, vtkInformation *info) override
Fill the input port information objects for this algorithm.
virtual void InitializeExtraPointDataArrays(vtkPointData *vtkNotUsed(outputPD))
Methods to append values to existing point data arrays that may only be desired on specific concrete ...
~vtkParticleTracerBase() override
void SetIntegrator(vtkInitialValueProblemSolver *)
vtkIntArray * GetParticleIds(vtkPointData *)
bool ComputeDomainExitLocation(double pos[4], double p2[4], double intersection[4], vtkGenericCell *cell)
This is an old routine kept for possible future use.
virtual void ResetCache()
virtual void AssignUniqueIds(vtkParticleTracerBaseNamespace::ParticleVector &localSeedPoints)
give each one a unique ID.
virtual int OutputParticles(vtkPolyData *poly)=0
vtkSmartPointer< vtkPolyData > Output
represent and manipulate point attribute data
Definition: vtkPointData.h:32
represent and manipulate 3D points
Definition: vtkPoints.h:34
Superclass for algorithms that produce only polydata as output.
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:85
A helper class for interpolating between times during particle tracing.
record modification and/or execution time
Definition: vtkTimeStamp.h:33
ParticleVector::iterator ParticleIterator
std::list< ParticleInformation > ParticleDataList
ParticleDataList::iterator ParticleListIterator
std::vector< ParticleInformation > ParticleVector
@ point
Definition: vtkX3D.h:242
@ info
Definition: vtkX3D.h:382
@ port
Definition: vtkX3D.h:453
@ time
Definition: vtkX3D.h:503
@ type
Definition: vtkX3D.h:522
int vtkTypeBool
Definition: vtkABI.h:69
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
int vtkIdType
Definition: vtkType.h:338