VTK  9.3.0
vtkParticleTracerBase.h
Go to the documentation of this file.
1// SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2// SPDX-License-Identifier: BSD-3-Clause
16#ifndef vtkParticleTracerBase_h
17#define vtkParticleTracerBase_h
18
19#include "vtkDeprecation.h" // For VTK_DEPRECATED_IN_9_2_0
20#include "vtkFiltersFlowPathsModule.h" // For export macro
22#include "vtkSmartPointer.h" // For vtkSmartPointer
23
24#include <list> // STL Header
25#include <mutex> // STL Header
26#include <vector> // STL Header
27
28VTK_ABI_NAMESPACE_BEGIN
30class vtkCellArray;
32class vtkDataArray;
33class vtkDataSet;
34class vtkDoubleArray;
35class vtkFloatArray;
36class vtkGenericCell;
38class vtkIntArray;
41class vtkPointData;
42class vtkPoints;
43class vtkPolyData;
46VTK_ABI_NAMESPACE_END
47
49{
50VTK_ABI_NAMESPACE_BEGIN
52{
53 double x[4];
54};
55using Position = struct Position_t;
56
58{
59 // These are used during iteration
64 // These are computed scalars we might display
66 int TimeStepAge; // amount of time steps the particle has advanced
68 int InjectedStepId; // time step the particle was injected
71 // These are useful to track for debugging etc
73 float age;
74 // these are needed across time steps to compute vorticity
75 float rotation;
77 float time;
78 float speed;
79 // once the particle is added, PointId is valid and is the tuple location
80 // in ProtoPD.
82 // if PointId is negative then in parallel this particle was just
83 // received and we need to get the tuple value from vtkPParticleTracerBase::Tail.
85};
87
88typedef std::vector<ParticleInformation> ParticleVector;
89typedef ParticleVector::iterator ParticleIterator;
90typedef std::list<ParticleInformation> ParticleDataList;
91typedef ParticleDataList::iterator ParticleListIterator;
92struct ParticleTracerFunctor;
93VTK_ABI_NAMESPACE_END
94};
95
96VTK_ABI_NAMESPACE_BEGIN
97class VTKFILTERSFLOWPATHS_EXPORT vtkParticleTracerBase : public vtkPolyDataAlgorithm
98{
99public:
100 friend struct vtkParticleTracerBaseNamespace::ParticleTracerFunctor;
102 {
107 UNKNOWN
108 };
109
111 void PrintSelf(ostream& os, vtkIndent indent) override;
113
115
120 vtkGetMacro(ComputeVorticity, bool);
123
125
128 vtkGetMacro(TerminalSpeed, double);
129 void SetTerminalSpeed(double);
131
133
137 vtkGetMacro(RotationScale, double);
138 void SetRotationScale(double);
140
142
146 vtkSetMacro(IgnorePipelineTime, vtkTypeBool);
147 vtkGetMacro(IgnorePipelineTime, vtkTypeBool);
148 vtkBooleanMacro(IgnorePipelineTime, vtkTypeBool);
150
152
161 vtkGetMacro(ForceReinjectionEveryNSteps, int);
164
166
172 void SetTerminationTime(double t);
173 vtkGetMacro(TerminationTime, double);
175
177 vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
178
179 void SetIntegratorType(int type);
181
183
187 vtkGetMacro(StartTime, double);
188 void SetStartTime(double t);
190
192
201 vtkSetMacro(StaticSeeds, vtkTypeBool);
202 vtkGetMacro(StaticSeeds, vtkTypeBool);
204
209 {
210 DIFFERENT = 0,
211 STATIC = 1,
212 LINEAR_TRANSFORMATION = 2,
213 SAME_TOPOLOGY = 3
214 };
215
217 /*
218 * Set/Get the type of variance of the mesh over time.
219 *
220 * DIFFERENT = 0,
221 * STATIC = 1,
222 * LINEAR_TRANSFORMATION = 2
223 * SAME_TOPOLOGY = 3
224 */
225 virtual void SetMeshOverTime(int meshOverTime);
226 virtual int GetMeshOverTimeMinValue() { return DIFFERENT; }
227 virtual int GetMeshOverTimeMaxValue() { return SAME_TOPOLOGY; }
228 void SetMeshOverTimeToDifferent() { this->SetMeshOverTime(DIFFERENT); }
229 void SetMeshOverTimeToStatic() { this->SetMeshOverTime(STATIC); }
230 void SetMeshOverTimeToLinearTransformation() { this->SetMeshOverTime(LINEAR_TRANSFORMATION); }
231 void SetMeshOverTimeToSameTopology() { this->SetMeshOverTime(SAME_TOPOLOGY); }
232 vtkGetMacro(MeshOverTime, int);
234
236
245 VTK_DEPRECATED_IN_9_2_0("Use SetMeshOverTime instead")
246 virtual void SetStaticMesh(vtkTypeBool staticMesh)
247 {
248 this->SetMeshOverTime(staticMesh ? STATIC : DIFFERENT);
249 }
250 VTK_DEPRECATED_IN_9_2_0("Use GetMeshOverTime instead")
251 virtual vtkTypeBool GetStaticMesh() { return this->MeshOverTime == STATIC; }
253
254 enum
255 {
257 INTERPOLATOR_WITH_CELL_LOCATOR
258 };
259
271 void SetInterpolatorType(int interpolatorType);
272
281
289
291
298 vtkGetObjectMacro(ParticleWriter, vtkAbstractParticleWriter);
300
302
306 vtkSetFilePathMacro(ParticleFileName);
307 vtkGetFilePathMacro(ParticleFileName);
309
311
315 vtkSetMacro(EnableParticleWriting, vtkTypeBool);
316 vtkGetMacro(EnableParticleWriting, vtkTypeBool);
317 vtkBooleanMacro(EnableParticleWriting, vtkTypeBool);
319
321
326 vtkSetMacro(DisableResetCache, vtkTypeBool);
327 vtkGetMacro(DisableResetCache, vtkTypeBool);
328 vtkBooleanMacro(DisableResetCache, vtkTypeBool);
330
332
338
340
344 vtkGetMacro(ForceSerialExecution, bool);
345 vtkSetMacro(ForceSerialExecution, bool);
346 vtkBooleanMacro(ForceSerialExecution, bool);
348protected:
349 vtkSmartPointer<vtkPolyData> Output; // managed by child classes
351
356 vtkIdType UniqueIdCounter; // global Id counter used to give particles a stamp
358 vtkSmartPointer<vtkPointData> ParticlePointData; // the current particle point data consistent
359 // with particle history
360 // Everything related to time
361 vtkTypeBool IgnorePipelineTime; // whether to use the pipeline time for termination
362 vtkTypeBool DisableResetCache; // whether to enable ResetCache() method
364
365 // Control execution as serial or threaded
367
370
371 //
372 // Make sure the pipeline knows what type we expect as input
373 //
374 int FillInputPortInformation(int port, vtkInformation* info) override;
375
376 //
377 // The usual suspects
378 //
380 vtkInformationVector* outputVector) override;
381
382 //
383 // Store any information we need in the output and fetch what we can
384 // from the input
385 //
387 vtkInformationVector* outputVector) override;
388
389 //
390 // Compute input time steps given the output step
391 //
393 vtkInformationVector* outputVector) override;
394
395 //
396 // what the pipeline calls for each time step
397 //
399 vtkInformationVector* outputVector) override;
400
401 //
402 // these routines are internally called to actually generate the output
403 //
404 virtual int ProcessInput(vtkInformationVector** inputVector);
405
406 // This is the main part of the algorithm:
407 // * move all the particles one step
408 // * Reinject particles (by adding them to this->ParticleHistories)
409 // either at the beginning or at the end of each step (modulo
410 // this->ForceReinjectionEveryNSteps)
411 // * Output a polydata representing the moved particles
412 // Note that if the starting and the ending time coincide, the polydata is still valid.
413 virtual vtkPolyData* Execute(vtkInformationVector** inputVector);
414
415 // the RequestData will call these methods in turn
416 virtual void Initialize() {} // the first iteration
417 virtual int OutputParticles(vtkPolyData* poly) = 0; // every iteration
418 virtual void Finalize() {} // the last iteration
419
424 virtual std::vector<vtkDataSet*> GetSeedSources(vtkInformationVector* inputVector, int timeStep);
425
426 // Initialization of input (vector-field) geometry
429
436
438 vtkParticleTracerBaseNamespace::ParticleVector& candidates, std::vector<int>& passed);
439
446 virtual void AssignSeedsToProcessors(double time, vtkDataSet* source, int sourceID, int ptId,
447 vtkParticleTracerBaseNamespace::ParticleVector& localSeedPoints, int& localAssignedCount);
448
454
460
466 virtual bool UpdateParticleListFromOtherProcesses() { return false; }
467
472 double currentTime, double targetTime, vtkInitialValueProblemSolver* integrator,
473 vtkTemporalInterpolatedVelocityField* interpolator, vtkDoubleArray* cellVectors,
474 std::atomic<vtkIdType>& particleCount, std::mutex& eraseMutex, bool sequential);
475
476 // if the particle is added to send list, then returns value is 1,
477 // if it is kept on this process after a retry return value is 0
480 {
481 return true;
482 }
483
491 double pos[4], double p2[4], double intersection[4], vtkGenericCell* cell);
492
493 //
494 // Scalar arrays that are generated as each particle is updated
495 //
497
507
508 // utility function we use to test if a point is inside any of our local datasets
509 bool InsideBounds(double point[]);
510
512 vtkGenericCell* cell, double pcoords[3], vtkDoubleArray* cellVectors, double vorticity[3]);
513
514 //------------------------------------------------------
515
516 double GetCacheDataTime(int i);
518
519 virtual void ResetCache();
521 vtkTemporalInterpolatedVelocityField* interpolator, vtkIdType particleId,
522 vtkDoubleArray* cellVectors);
523
525
530 virtual bool IsPointDataValid(vtkDataObject* input);
531 bool IsPointDataValid(vtkCompositeDataSet* input, std::vector<std::string>& arrayNames);
532 void GetPointDataArrayNames(vtkDataSet* input, std::vector<std::string>& names);
534
535 vtkGetMacro(ReinjectionCounter, int);
536 vtkGetMacro(CurrentTimeValue, double);
537
538 void ResizeArrays(vtkIdType numTuples);
539
544 virtual void InitializeExtraPointDataArrays(vtkPointData* vtkNotUsed(outputPD)) {}
545
548 {
549 }
550
552
557 virtual void AddRestartSeeds(vtkInformationVector** /*inputVector*/) {}
558
564
565private:
575 bool RetryWithPush(vtkParticleTracerBaseNamespace::ParticleInformation& info, double* point1,
576 double delT, int subSteps, vtkTemporalInterpolatedVelocityField* interpolator);
577
578 bool SetTerminationTimeNoModify(double t);
579
580 // Parameters of tracing
582 double IntegrationStep;
583 double MaximumError;
584 bool ComputeVorticity;
585 double RotationScale;
586 double TerminalSpeed;
587
588 // A counter to keep track of how many times we reinjected
589 int ReinjectionCounter;
590
591 // Important for Caching of Cells/Ids/Weights etc
592 vtkTypeBool AllFixedGeometry;
593 int MeshOverTime;
594 vtkTypeBool StaticSeeds;
595
596 std::vector<double> InputTimeValues;
597 double StartTime;
598 double TerminationTime;
599 double CurrentTimeValue;
600
601 int StartTimeStep; // InputTimeValues[StartTimeStep] <= StartTime <=
602 // InputTimeValues[StartTimeStep+1]
603 int CurrentTimeStep;
604 int TerminationTimeStep; // computed from start time
605 bool FirstIteration;
606
607 // Innjection parameters
608 int ForceReinjectionEveryNSteps;
609 vtkTimeStamp ParticleInjectionTime;
610 bool HasCache;
611
612 // Particle writing to disk
613 vtkAbstractParticleWriter* ParticleWriter;
614 char* ParticleFileName;
615 vtkTypeBool EnableParticleWriting;
616
617 // The main lists which are held during operation- between time step updates
619
620 // The velocity interpolator
622
623 // Data for time step CurrentTimeStep-1 and CurrentTimeStep
625
626 // Cache bounds info for each dataset we will use repeatedly
627 struct bounds_t
628 {
629 double b[6];
630 };
631 using bounds = struct bounds_t;
632 std::vector<bounds> CachedBounds[2];
633
634 // variables used by Execute() to produce output
635
636 vtkSmartPointer<vtkDataSet> DataReferenceT[2];
637
638 vtkSmartPointer<vtkPoints> OutputCoordinates;
639 vtkSmartPointer<vtkIdTypeArray> ParticleCellsConnectivity;
640 vtkSmartPointer<vtkIdTypeArray> ParticleCellsOffsets;
641 vtkSmartPointer<vtkCellArray> ParticleCells;
642
645 vtkSmartPointer<vtkSignedCharArray> ParticleSourceIds;
646 vtkSmartPointer<vtkIntArray> InjectedPointIds;
647 vtkSmartPointer<vtkIntArray> InjectedStepIds;
648 vtkSmartPointer<vtkIntArray> ErrorCodeArray;
649 vtkSmartPointer<vtkFloatArray> ParticleVorticity;
650 vtkSmartPointer<vtkFloatArray> ParticleRotation;
651 vtkSmartPointer<vtkFloatArray> ParticleAngularVel;
652 vtkSmartPointer<vtkPointData> OutputPointData;
653
654 // temp array
656
658 void operator=(const vtkParticleTracerBase&) = delete;
659 vtkTimeStamp ExecuteTime;
660
661 unsigned int NumberOfParticles();
662
665
666 static const double Epsilon;
667};
668
669VTK_ABI_NAMESPACE_END
670#endif
abstract class to write particle data to file
Proxy object to connect input/output ports.
object to represent cell connectivity
abstract superclass for composite (multi-block or AMR) datasets
abstract superclass for arrays of numeric data
general representation of visualization data
abstract class to specify dataset behavior
Definition vtkDataSet.h:53
dynamic, self-adjusting array of double
dynamic, self-adjusting array of float
provides thread-safe access to cells
a simple class to control print indentation
Definition vtkIndent.h:29
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:35
Composite dataset that organizes datasets into blocks.
Multiprocessing communication superclass.
A particle tracer for vector fields.
vtkTypeBool DisableResetCache
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
vtkFloatArray * GetParticleVorticity(vtkPointData *)
vtkSmartPointer< vtkPointData > ProtoPD
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
virtual void RenameGhostArray(vtkPointData *pd)
Recover the ghost array in provided point data and rename it by adding "Original_" in front of it,...
void SetComputeVorticity(bool)
Turn on/off vorticity computation at streamline points (necessary for generating proper stream-ribbon...
vtkSmartPointer< vtkPointData > ParticlePointData
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
vtkIntArray * GetInjectedPointIds(vtkPointData *)
vtkFloatArray * GetParticleAge(vtkPointData *)
void TestParticles(vtkParticleTracerBaseNamespace::ParticleVector &candidates, std::vector< int > &passed)
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)
Methods that check that the input arrays are ordered the same on all data sets.
void SetInterpolatorTypeToCellLocator()
Set the velocity field interpolator type to one that uses a cell locator to perform spatial searching...
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.
vtkIntArray * GetInjectedStepIds(vtkPointData *)
virtual bool UpdateParticleListFromOtherProcesses()
this is used during classification of seed points and also between iterations of the main loop as par...
virtual std::vector< vtkDataSet * > GetSeedSources(vtkInformationVector *inputVector, int timeStep)
Method to get the data set seed sources.
vtkIntArray * GetErrorCodeArr(vtkPointData *)
void SetParticle(vtkParticleTracerBaseNamespace::ParticleInformation &info, double *velocity, vtkTemporalInterpolatedVelocityField *interpolator, vtkIdType particleId, vtkDoubleArray *cellVectors)
void UpdateParticleList(vtkParticleTracerBaseNamespace::ParticleVector &candidates)
and sending between processors, into a list, which is used as the master list on this processor
vtkGetFilePathMacro(ParticleFileName)
Set/Get the filename to be used with the particle writer when dumping particles to disk.
virtual void SetToExtraPointDataArrays(vtkIdType, vtkParticleTracerBaseNamespace::ParticleInformation &)
vtkFloatArray * GetParticleAngularVel(vtkPointData *)
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
Upstream/Downstream requests form the generalized interface through which executives invoke a algorit...
void AddSourceConnection(vtkAlgorithmOutput *input)
Provide support for multiple seed sources.
virtual vtkPolyData * Execute(vtkInformationVector **inputVector)
void RemoveAllSources()
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)
virtual bool SendParticleToAnotherProcess(vtkParticleTracerBaseNamespace::ParticleInformation &, vtkParticleTracerBaseNamespace::ParticleInformation &, vtkPointData *)
vtkIntArray * GetParticleIds(vtkPointData *)
MeshOverTimeTypes
Types of Variance of Mesh over time.
void SetForceReinjectionEveryNSteps(int)
When animating particles, it is nice to inject new ones every Nth step to produce a continuous flow.
int UpdateDataCache(vtkDataObject *td)
void SetTerminalSpeed(double)
Specify the terminal speed value, below which integration is terminated.
virtual void SetParticleWriter(vtkAbstractParticleWriter *pw)
Set/Get the Writer associated with this Particle Tracer Ideally a parallel IO capable vtkH5PartWriter...
void SetRotationScale(double)
This can be used to scale the rate with which the streamribbons twist.
vtkParticleTracerBaseNamespace::ParticleDataList ParticleHistories
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
vtkFloatArray * GetParticleRotation(vtkPointData *)
void IntegrateParticle(vtkParticleTracerBaseNamespace::ParticleListIterator &it, double currentTime, double targetTime, vtkInitialValueProblemSolver *integrator, vtkTemporalInterpolatedVelocityField *interpolator, vtkDoubleArray *cellVectors, std::atomic< vtkIdType > &particleCount, std::mutex &eraseMutex, bool sequential)
particle between the two times supplied.
void GetPointDataArrayNames(vtkDataSet *input, std::vector< std::string > &names)
Methods that check that the input arrays are ordered the same on all data sets.
void SetInterpolatorType(int interpolatorType)
Set the type of the velocity field interpolator to determine whether INTERPOLATOR_WITH_DATASET_POINT_...
virtual void SetMeshOverTime(int meshOverTime)
int RequestUpdateExtent(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
This is called by the superclass.
void SetStartTime(double t)
Set the time value for particle tracing to begin.
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.
vtkTypeBool IgnorePipelineTime
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
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
vtkTemporalInterpolatedVelocityField * GetInterpolator()
void SetIntegratorType(int type)
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
int FillInputPortInformation(int port, vtkInformation *info) override
Fill the input port information objects for this algorithm.
vtkIdType UniqueIdCounter
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
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 *)
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()
void ResizeArrays(vtkIdType numTuples)
void SetInterpolatorTypeToDataSetPointLocator()
Set the velocity field interpolator type to one that uses a point locator to perform local spatial se...
vtkSignedCharArray * GetParticleSourceIds(vtkPointData *)
virtual void AssignUniqueIds(vtkParticleTracerBaseNamespace::ParticleVector &localSeedPoints)
give each one a unique ID.
vtkSetFilePathMacro(ParticleFileName)
Set/Get the filename to be used with the particle writer when dumping particles to disk.
virtual int OutputParticles(vtkPolyData *poly)=0
vtkSmartPointer< vtkPolyData > Output
represent and manipulate point attribute data
represent and manipulate 3D points
Definition vtkPoints.h:29
Superclass for algorithms that produce only polydata as output.
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition vtkPolyData.h:80
dynamic, self-adjusting array of signed char
Hold a reference to a vtkObjectBase instance.
A helper class for interpolating between times during particle tracing.
record modification and/or execution time
ParticleVector::iterator ParticleIterator
std::list< ParticleInformation > ParticleDataList
ParticleDataList::iterator ParticleListIterator
std::vector< ParticleInformation > ParticleVector
int vtkTypeBool
Definition vtkABI.h:64
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
#define VTK_DEPRECATED_IN_9_2_0(reason)
int vtkIdType
Definition vtkType.h:315