VTK  9.2.6
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 "vtkDeprecation.h" // For VTK_DEPRECATED_IN_9_2_0
32#include "vtkFiltersFlowPathsModule.h" // For export macro
34#include "vtkSmartPointer.h" // For vtkSmartPointer
35
36#include <list> // STL Header
37#include <mutex> // STL Header
38#include <vector> // STL Header
39
41class vtkCellArray;
43class vtkDataArray;
44class vtkDataSet;
45class vtkDoubleArray;
46class vtkFloatArray;
47class vtkGenericCell;
49class vtkIntArray;
52class vtkPointData;
53class vtkPoints;
54class vtkPolyData;
57
59{
61{
62 double x[4];
63};
64using Position = struct Position_t;
65
67{
68 // These are used during iteration
73 // These are computed scalars we might display
75 int TimeStepAge; // amount of time steps the particle has advanced
77 int InjectedStepId; // time step the particle was injected
80 // These are useful to track for debugging etc
82 float age;
83 // these are needed across time steps to compute vorticity
84 float rotation;
86 float time;
87 float speed;
88 // once the particle is added, PointId is valid and is the tuple location
89 // in ProtoPD.
91 // if PointId is negative then in parallel this particle was just
92 // received and we need to get the tuple value from vtkPParticleTracerBase::Tail.
94};
96
97typedef std::vector<ParticleInformation> ParticleVector;
98typedef ParticleVector::iterator ParticleIterator;
99typedef std::list<ParticleInformation> ParticleDataList;
100typedef ParticleDataList::iterator ParticleListIterator;
101struct ParticleTracerFunctor;
102};
103
104class VTKFILTERSFLOWPATHS_EXPORT vtkParticleTracerBase : public vtkPolyDataAlgorithm
105{
106public:
107 friend struct vtkParticleTracerBaseNamespace::ParticleTracerFunctor;
109 {
114 UNKNOWN
115 };
116
118 void PrintSelf(ostream& os, vtkIndent indent) override;
120
122
127 vtkGetMacro(ComputeVorticity, bool);
130
132
135 vtkGetMacro(TerminalSpeed, double);
136 void SetTerminalSpeed(double);
138
140
144 vtkGetMacro(RotationScale, double);
145 void SetRotationScale(double);
147
149
153 vtkSetMacro(IgnorePipelineTime, vtkTypeBool);
154 vtkGetMacro(IgnorePipelineTime, vtkTypeBool);
155 vtkBooleanMacro(IgnorePipelineTime, vtkTypeBool);
157
159
168 vtkGetMacro(ForceReinjectionEveryNSteps, vtkTypeBool);
171
173
179 void SetTerminationTime(double t);
180 vtkGetMacro(TerminationTime, double);
182
184 vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
185
186 void SetIntegratorType(int type);
188
190
194 vtkGetMacro(StartTime, double);
195 void SetStartTime(double t);
197
199
208 vtkSetMacro(StaticSeeds, vtkTypeBool);
209 vtkGetMacro(StaticSeeds, vtkTypeBool);
211
216 {
217 DIFFERENT = 0,
218 STATIC = 1,
219 LINEAR_TRANSFORMATION = 2,
220 SAME_TOPOLOGY = 3
221 };
222
224 /*
225 * Set/Get the type of variance of the mesh over time.
226 *
227 * DIFFERENT = 0,
228 * STATIC = 1,
229 * LINEAR_TRANSFORMATION = 2
230 * SAME_TOPOLOGY = 3
231 */
232 virtual void SetMeshOverTime(int meshOverTime);
233 virtual int GetMeshOverTimeMinValue() { return DIFFERENT; }
234 virtual int GetMeshOverTimeMaxValue() { return SAME_TOPOLOGY; }
235 void SetMeshOverTimeToDifferent() { this->SetMeshOverTime(DIFFERENT); }
236 void SetMeshOverTimeToStatic() { this->SetMeshOverTime(STATIC); }
237 void SetMeshOverTimeToLinearTransformation() { this->SetMeshOverTime(LINEAR_TRANSFORMATION); }
238 void SetMeshOverTimeToSameTopology() { this->SetMeshOverTime(SAME_TOPOLOGY); }
239 vtkGetMacro(MeshOverTime, int);
241
243
252 VTK_DEPRECATED_IN_9_2_0("Use SetMeshOverTime instead")
253 virtual void SetStaticMesh(vtkTypeBool staticMesh)
254 {
255 this->SetMeshOverTime(staticMesh ? STATIC : DIFFERENT);
256 }
257 VTK_DEPRECATED_IN_9_2_0("Use GetMeshOverTime instead")
258 virtual vtkTypeBool GetStaticMesh() { return this->MeshOverTime == STATIC; }
260
261 enum
262 {
264 INTERPOLATOR_WITH_CELL_LOCATOR
265 };
266
278 void SetInterpolatorType(int interpolatorType);
279
288
296
298
305 vtkGetObjectMacro(ParticleWriter, vtkAbstractParticleWriter);
307
309
313 vtkSetFilePathMacro(ParticleFileName);
314 vtkGetFilePathMacro(ParticleFileName);
316
318
322 vtkSetMacro(EnableParticleWriting, vtkTypeBool);
323 vtkGetMacro(EnableParticleWriting, vtkTypeBool);
324 vtkBooleanMacro(EnableParticleWriting, vtkTypeBool);
326
328
333 vtkSetMacro(DisableResetCache, vtkTypeBool);
334 vtkGetMacro(DisableResetCache, vtkTypeBool);
335 vtkBooleanMacro(DisableResetCache, vtkTypeBool);
337
339
345
347
351 vtkGetMacro(ForceSerialExecution, bool);
352 vtkSetMacro(ForceSerialExecution, bool);
353 vtkBooleanMacro(ForceSerialExecution, bool);
355protected:
356 vtkSmartPointer<vtkPolyData> Output; // managed by child classes
358
363 vtkIdType UniqueIdCounter; // global Id counter used to give particles a stamp
365 vtkSmartPointer<vtkPointData> ParticlePointData; // the current particle point data consistent
366 // with particle history
367 // Everything related to time
368 vtkTypeBool IgnorePipelineTime; // whether to use the pipeline time for termination
369 vtkTypeBool DisableResetCache; // whether to enable ResetCache() method
371
372 // Control execution as serial or threaded
374
377
378 //
379 // Make sure the pipeline knows what type we expect as input
380 //
381 int FillInputPortInformation(int port, vtkInformation* info) override;
382
383 //
384 // The usual suspects
385 //
387 vtkInformationVector* outputVector) override;
388
389 //
390 // Store any information we need in the output and fetch what we can
391 // from the input
392 //
394 vtkInformationVector* outputVector) override;
395
396 //
397 // Compute input time steps given the output step
398 //
400 vtkInformationVector* outputVector) override;
401
402 //
403 // what the pipeline calls for each time step
404 //
406 vtkInformationVector* outputVector) override;
407
408 //
409 // these routines are internally called to actually generate the output
410 //
411 virtual int ProcessInput(vtkInformationVector** inputVector);
412
413 // This is the main part of the algorithm:
414 // * move all the particles one step
415 // * Reinject particles (by adding them to this->ParticleHistories)
416 // either at the beginning or at the end of each step (modulo
417 // this->ForceReinjectionEveryNSteps)
418 // * Output a polydata representing the moved particles
419 // Note that if the starting and the ending time coincide, the polydata is still valid.
420 virtual vtkPolyData* Execute(vtkInformationVector** inputVector);
421
422 // the RequestData will call these methods in turn
423 virtual void Initialize() {} // the first iteration
424 virtual int OutputParticles(vtkPolyData* poly) = 0; // every iteration
425 virtual void Finalize() {} // the last iteration
426
431 virtual std::vector<vtkDataSet*> GetSeedSources(vtkInformationVector* inputVector, int timeStep);
432
433 // Initialization of input (vector-field) geometry
436
443
445 vtkParticleTracerBaseNamespace::ParticleVector& candidates, std::vector<int>& passed);
446
453 virtual void AssignSeedsToProcessors(double time, vtkDataSet* source, int sourceID, int ptId,
454 vtkParticleTracerBaseNamespace::ParticleVector& localSeedPoints, int& localAssignedCount);
455
461
467
473 virtual bool UpdateParticleListFromOtherProcesses() { return false; }
474
479 double currentTime, double targetTime, vtkInitialValueProblemSolver* integrator,
480 vtkTemporalInterpolatedVelocityField* interpolator, vtkDoubleArray* cellVectors,
481 std::atomic<vtkIdType>& particleCount, std::mutex& eraseMutex, bool sequential);
482
483 // if the particle is added to send list, then returns value is 1,
484 // if it is kept on this process after a retry return value is 0
487 {
488 return true;
489 }
490
498 double pos[4], double p2[4], double intersection[4], vtkGenericCell* cell);
499
500 //
501 // Scalar arrays that are generated as each particle is updated
502 //
504
514
515 // utility function we use to test if a point is inside any of our local datasets
516 bool InsideBounds(double point[]);
517
519 vtkGenericCell* cell, double pcoords[3], vtkDoubleArray* cellVectors, double vorticity[3]);
520
521 //------------------------------------------------------
522
523 double GetCacheDataTime(int i);
525
526 virtual void ResetCache();
528 vtkTemporalInterpolatedVelocityField* interpolator, vtkIdType particleId,
529 vtkDoubleArray* cellVectors);
530
532
537 virtual bool IsPointDataValid(vtkDataObject* input);
538 bool IsPointDataValid(vtkCompositeDataSet* input, std::vector<std::string>& arrayNames);
539 void GetPointDataArrayNames(vtkDataSet* input, std::vector<std::string>& names);
541
542 vtkGetMacro(ReinjectionCounter, int);
543 vtkGetMacro(CurrentTimeValue, double);
544
545 void ResizeArrays(vtkIdType numTuples);
546
551 virtual void InitializeExtraPointDataArrays(vtkPointData* vtkNotUsed(outputPD)) {}
552
555 {
556 }
557
559
564 virtual void AddRestartSeeds(vtkInformationVector** /*inputVector*/) {}
565
566private:
576 bool RetryWithPush(vtkParticleTracerBaseNamespace::ParticleInformation& info, double* point1,
577 double delT, int subSteps, vtkTemporalInterpolatedVelocityField* interpolator);
578
579 bool SetTerminationTimeNoModify(double t);
580
581 // Parameters of tracing
583 double IntegrationStep;
584 double MaximumError;
585 bool ComputeVorticity;
586 double RotationScale;
587 double TerminalSpeed;
588
589 // A counter to keep track of how many times we reinjected
590 int ReinjectionCounter;
591
592 // Important for Caching of Cells/Ids/Weights etc
593 vtkTypeBool AllFixedGeometry;
594 int MeshOverTime;
595 vtkTypeBool StaticSeeds;
596
597 std::vector<double> InputTimeValues;
598 double StartTime;
599 double TerminationTime;
600 double CurrentTimeValue;
601
602 int StartTimeStep; // InputTimeValues[StartTimeStep] <= StartTime <=
603 // InputTimeValues[StartTimeStep+1]
604 int CurrentTimeStep;
605 int TerminationTimeStep; // computed from start time
606 bool FirstIteration;
607
608 // Innjection parameters
609 vtkTypeBool ForceReinjectionEveryNSteps;
610 vtkTimeStamp ParticleInjectionTime;
611 bool HasCache;
612
613 // Particle writing to disk
614 vtkAbstractParticleWriter* ParticleWriter;
615 char* ParticleFileName;
616 vtkTypeBool EnableParticleWriting;
617
618 // The main lists which are held during operation- between time step updates
620
621 // The velocity interpolator
623
624 // Data for time step CurrentTimeStep-1 and CurrentTimeStep
626
627 // Cache bounds info for each dataset we will use repeatedly
628 struct bounds_t
629 {
630 double b[6];
631 };
632 using bounds = struct bounds_t;
633 std::vector<bounds> CachedBounds[2];
634
635 // variables used by Execute() to produce output
636
637 vtkSmartPointer<vtkDataSet> DataReferenceT[2];
638
639 vtkSmartPointer<vtkPoints> OutputCoordinates;
640 vtkSmartPointer<vtkIdTypeArray> ParticleCellsConnectivity;
641 vtkSmartPointer<vtkIdTypeArray> ParticleCellsOffsets;
642 vtkSmartPointer<vtkCellArray> ParticleCells;
643
646 vtkSmartPointer<vtkSignedCharArray> ParticleSourceIds;
647 vtkSmartPointer<vtkIntArray> InjectedPointIds;
648 vtkSmartPointer<vtkIntArray> InjectedStepIds;
649 vtkSmartPointer<vtkIntArray> ErrorCodeArray;
650 vtkSmartPointer<vtkFloatArray> ParticleVorticity;
651 vtkSmartPointer<vtkFloatArray> ParticleRotation;
652 vtkSmartPointer<vtkFloatArray> ParticleAngularVel;
653 vtkSmartPointer<vtkPointData> OutputPointData;
654
655 // temp array
657
659 void operator=(const vtkParticleTracerBase&) = delete;
660 vtkTimeStamp ExecuteTime;
661
662 unsigned int NumberOfParticles();
663
666
667 static const double Epsilon;
668};
669
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:63
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:40
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:46
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...
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.
void SetForceReinjectionEveryNSteps(vtkTypeBool)
When animating particles, it is nice to inject new ones every Nth step to produce a continuous flow.
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.
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:40
Superclass for algorithms that produce only polydata as output.
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition vtkPolyData.h:91
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:69
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:332