VTK  9.4.20241202
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_4_0
20#include "vtkFiltersFlowPathsModule.h" // For export macro
22#include "vtkSmartPointer.h" // For vtkSmartPointer
23#include "vtkTemporalAlgorithm.h" // Fro vtkTemporalAlgorithm
24#include "vtkWeakPointer.h" // For vtkWeakPointer
25
26#include <list> // STL Header
27#include <mutex> // STL Header
28#include <numeric> // STL Header
29#include <unordered_map> // STL Header
30#include <vector> // STL Header
31
32#ifndef __VTK_WRAP__
33#define vtkPolyDataAlgorithm vtkTemporalAlgorithm<vtkPolyDataAlgorithm>
34#endif
35
36VTK_ABI_NAMESPACE_BEGIN
39class vtkDataArray;
40class vtkDataSet;
41class vtkDoubleArray;
42class vtkFloatArray;
43class vtkGenericCell;
45class vtkIntArray;
48class vtkPointData;
49class vtkPoints;
50class vtkPolyData;
53VTK_ABI_NAMESPACE_END
54
56{
57VTK_ABI_NAMESPACE_BEGIN
59{
60 double x[4];
61};
62using Position = struct Position_t;
63
65{
66 // These are used during iteration
71 // These are computed scalars we might display
73 int TimeStepAge; // amount of time steps the particle has advanced
75 int InjectedStepId; // time step the particle was injected
77 // These are useful to track for debugging etc
79 float age;
80 // these are needed across time steps to compute vorticity
81 float rotation;
83 float time;
84 float speed;
85 // once the particle is added, PointId is valid and is the tuple location
86 // in ProtoPD.
88
89 double velocity[3];
90};
92
93typedef std::vector<ParticleInformation> ParticleVector;
94typedef ParticleVector::iterator ParticleIterator;
95typedef std::list<ParticleInformation> ParticleDataList;
96typedef ParticleDataList::iterator ParticleListIterator;
97struct ParticleTracerFunctor;
98VTK_ABI_NAMESPACE_END
99}
100
101VTK_ABI_NAMESPACE_BEGIN
102class VTKFILTERSFLOWPATHS_EXPORT vtkParticleTracerBase : public vtkPolyDataAlgorithm
103{
104public:
105 friend struct vtkParticleTracerBaseNamespace::ParticleTracerFunctor;
107 {
112 UNKNOWN
113 };
114
117
118#ifndef __VTK_WRAP__
119#undef vtkPolyDataAlgorithm
120#endif
122 void PrintSelf(ostream& os, vtkIndent indent) override;
124
125#if defined(__VTK_WRAP__) || defined(__WRAP_GCCXML)
127#endif
128
130
137
139
144 vtkGetMacro(ComputeVorticity, bool);
147
149
152 vtkGetMacro(TerminalSpeed, double);
153 void SetTerminalSpeed(double);
155
157
161 vtkGetMacro(RotationScale, double);
162 void SetRotationScale(double);
164
166
170 vtkSetMacro(IgnorePipelineTime, vtkTypeBool);
171 vtkGetMacro(IgnorePipelineTime, vtkTypeBool);
172 vtkBooleanMacro(IgnorePipelineTime, vtkTypeBool);
174
176
185 vtkGetMacro(ForceReinjectionEveryNSteps, int);
188
190
197 "Please edit the TIME_STEPS information key in vtkAlgorithm::RequestInformation() instead")
198 virtual void SetTerminationTime(double) {}
200 "Please edit the TIME_STEPS information key in vtkAlgorithm::RequestInformation() instead")
201 virtual double GetTerminationTime() { return std::numeric_limits<double>::quiet_NaN(); }
203 "Please edit the TIME_STEPS information key in vtkAlgorithm::RequestInformation() instead")
204 virtual void SetStartTime(double) {}
206 "Please edit the TIME_STEPS information key in vtkAlgorithm::RequestInformation() instead")
207 virtual double GetStartTime() { return std::numeric_limits<double>::quiet_NaN(); }
209
211
214 VTK_DEPRECATED_IN_9_4_0("Caching is now automated")
215 virtual void SetDisableResetCache(bool) {}
216 VTK_DEPRECATED_IN_9_4_0("Caching is now automated")
217 virtual bool GetDisableResetCache() { return false; }
218 VTK_DEPRECATED_IN_9_4_0("Caching is now automated")
219 virtual void DisableResetCacheOn() {}
220 VTK_DEPRECATED_IN_9_4_0("Caching is now automated")
221 virtual void DisableResetCacheOff() {}
223
225 vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
226
227 void SetIntegratorType(int type);
229
231
240 vtkSetMacro(StaticSeeds, vtkTypeBool);
241 vtkGetMacro(StaticSeeds, vtkTypeBool);
243
248 {
249 DIFFERENT = 0,
250 STATIC = 1,
251 LINEAR_TRANSFORMATION = 2,
252 SAME_TOPOLOGY = 3
253 };
254
256 /*
257 * Set/Get the type of variance of the mesh over time.
258 *
259 * DIFFERENT = 0,
260 * STATIC = 1,
261 * LINEAR_TRANSFORMATION = 2
262 * SAME_TOPOLOGY = 3
263 */
264 virtual void SetMeshOverTime(int meshOverTime);
265 virtual int GetMeshOverTimeMinValue() { return DIFFERENT; }
266 virtual int GetMeshOverTimeMaxValue() { return SAME_TOPOLOGY; }
267 void SetMeshOverTimeToDifferent() { this->SetMeshOverTime(DIFFERENT); }
268 void SetMeshOverTimeToStatic() { this->SetMeshOverTime(STATIC); }
269 void SetMeshOverTimeToLinearTransformation() { this->SetMeshOverTime(LINEAR_TRANSFORMATION); }
270 void SetMeshOverTimeToSameTopology() { this->SetMeshOverTime(SAME_TOPOLOGY); }
271 vtkGetMacro(MeshOverTime, int);
273
274 enum
275 {
277 INTERPOLATOR_WITH_CELL_LOCATOR
278 };
279
291 void SetInterpolatorType(int interpolatorType);
292
301
309
311
318 vtkGetObjectMacro(ParticleWriter, vtkAbstractParticleWriter);
320
322
326 vtkSetFilePathMacro(ParticleFileName);
327 vtkGetFilePathMacro(ParticleFileName);
329
331
335 vtkSetMacro(EnableParticleWriting, vtkTypeBool);
336 vtkGetMacro(EnableParticleWriting, vtkTypeBool);
337 vtkBooleanMacro(EnableParticleWriting, vtkTypeBool);
339
341
347
349
353 vtkGetMacro(ForceSerialExecution, bool);
354 vtkSetMacro(ForceSerialExecution, bool);
355 vtkBooleanMacro(ForceSerialExecution, bool);
357protected:
359
365 vtkTypeBool IgnorePipelineTime; // whether to use the pipeline time for termination
367
368 // Control execution as serial or threaded
370
372
375
376 //
377 // Make sure the pipeline knows what type we expect as input
378 //
379 int FillInputPortInformation(int port, vtkInformation* info) override;
380
384 int Initialize(vtkInformation* request, vtkInformationVector** inputVector,
385 vtkInformationVector* outputVector) override;
386
394 int Execute(vtkInformation* request, vtkInformationVector** inputVector,
395 vtkInformationVector* outputVector) override;
396
401
406 virtual std::vector<vtkDataSet*> GetSeedSources(vtkInformationVector* inputVector);
407
408 // Initialization of input (vector-field) geometry
410
416
423
425 vtkParticleTracerBaseNamespace::ParticleVector& candidates, std::vector<int>& passed);
426
433 virtual void AssignSeedsToProcessors(double time, vtkDataSet* source,
435
441 bool SendReceiveParticles(std::vector<vtkIdType>&);
442
448
455
460 double currentTime, double targetTime, vtkInitialValueProblemSolver* integrator,
461 vtkTemporalInterpolatedVelocityField* interpolator, vtkDoubleArray* cellVectors,
462 std::atomic<vtkIdType>& particleCount, std::mutex& eraseMutex, bool sequential);
463
471 double pos[4], double p2[4], double intersection[4], vtkGenericCell* cell);
472
473 //
474 // Scalar arrays that are generated as each particle is updated
475 //
477
487
488 // utility function we use to test if a point is inside any of our local datasets
489 bool InsideBounds(double point[]);
490
492 vtkGenericCell* cell, double pcoords[3], vtkDoubleArray* cellVectors, double vorticity[3]);
493
494 //------------------------------------------------------
495
497 vtkTemporalInterpolatedVelocityField* interpolator, vtkDoubleArray* cellVectors);
498
499 std::unordered_map<vtkIdType, int> InjectedPointIdToProcessId;
500
502
507 virtual bool IsPointDataValid(vtkDataObject* input);
508 bool IsPointDataValid(vtkCompositeDataSet* input, std::vector<std::string>& arrayNames);
509 void GetPointDataArrayNames(vtkDataSet* input, std::vector<std::string>& names);
511
512 vtkGetMacro(ReinjectionCounter, int);
513
514 void ResizeArrays(vtkIdType numTuples);
515
520 virtual void InitializeExtraPointDataArrays(vtkPointData* vtkNotUsed(outputPD)) {}
521
524 {
525 }
526
528
533 virtual void AddRestartSeeds(vtkInformationVector** /*inputVector*/) {}
534
545 double delT, int subSteps, vtkTemporalInterpolatedVelocityField* interpolator);
546
548
550
551 // Parameters of tracing
558
559 // A counter to keep track of how many times we reinjected
561
562 // Important for Caching of Cells/Ids/Weights etc
566
567 // Innjection parameters
571
572 // Particle writing to disk
576
578
579 // The velocity interpolator
581
582 // Data for time step CurrentTimeStep-1 and CurrentTimeStep
584
585 // MPI controller needed when running in parallel
587
593
601 std::unordered_map<vtkIdType, vtkParticleTracerBaseNamespace::ParticleInformation> MPIRecvList;
602
603 // Cache bounds info for each dataset we will use repeatedly
604 struct bounds_t
605 {
606 double b[6];
607 };
608 using bounds = struct bounds_t;
609 std::vector<bounds> CachedBounds[2];
610
611 // variables used by Execute() to produce output
612
614
617
628
629 // temp array
631
633 void operator=(const vtkParticleTracerBase&) = delete;
635
636 unsigned int NumberOfParticles();
637
638 friend class ParticlePathFilterInternal;
639 friend class StreaklineFilterInternal;
640
641 static const double Epsilon;
642};
643
644VTK_ABI_NAMESPACE_END
645#endif
abstract class to write particle data to file
Proxy object to connect input/output ports.
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:166
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:108
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
Composite dataset that organizes datasets into blocks.
Multiprocessing communication superclass.
Allocate and hold a VTK object.
Definition vtkNew.h:167
A particle tracer for vector fields.
vtkFloatArray * GetParticleVorticity(vtkPointData *)
vtkInitialValueProblemSolver * Integrator
vtkSmartPointer< vtkPointData > ProtoPD
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
vtkNew< vtkPoints > OutputCoordinates
vtkNew< vtkPointData > OutputPointData
void SetComputeVorticity(bool)
Turn on/off vorticity computation at streamline points (necessary for generating proper stream-ribbon...
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...
vtkNew< vtkFloatArray > ParticleRotation
vtkNew< vtkDoubleArray > CellVectors
vtkNew< vtkIdTypeArray > ParticleCellsConnectivity
vtkIntArray * GetInjectedStepIds(vtkPointData *)
virtual bool UpdateParticleListFromOtherProcesses()
this is used during classification of seed points and also between iterations of the main loop as par...
vtkIntArray * GetErrorCodeArr(vtkPointData *)
bool RetryWithPush(vtkParticleTracerBaseNamespace::ParticleInformation &info, double *point1, double delT, int subSteps, vtkTemporalInterpolatedVelocityField *interpolator)
When particles leave the domain, they must be collected and sent to the other processes for possible ...
vtkSmartPointer< vtkTemporalInterpolatedVelocityField > Interpolator
int Initialize(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
Resets internal cache for a clean start.
vtkNew< vtkIdTypeArray > InjectedPointIds
void operator=(const vtkParticleTracerBase &)=delete
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 *)
bool SendReceiveParticles(std::vector< vtkIdType > &)
this is used during classification of seed points and also between iterations of the main loop as par...
std::unordered_map< vtkIdType, int > InjectedPointIdToProcessId
void AddSourceConnection(vtkAlgorithmOutput *input)
Provide support for multiple seed sources.
vtkNew< vtkIntArray > ErrorCodeArray
void RemoveAllSources()
Provide support for multiple seed sources.
void CalculateVorticity(vtkGenericCell *cell, double pcoords[3], vtkDoubleArray *cellVectors, double vorticity[3])
void CreateProtoPD(vtkDataObject *input)
std::unordered_map< vtkIdType, vtkParticleTracerBaseNamespace::ParticleInformation > MPIRecvList
Storage of the particles we received.
virtual void SetController(vtkMultiProcessController *)
Get/Set the controller to use.
vtkIntArray * GetParticleIds(vtkPointData *)
vtkNew< vtkIntArray > InjectedStepIds
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.
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...
virtual std::vector< vtkDataSet * > GetSeedSources(vtkInformationVector *inputVector)
Method to get the data set seed sources.
void SetRotationScale(double)
This can be used to scale the rate with which the streamribbons twist.
vtkParticleTracerBaseNamespace::ParticleVector MPISendList
Storage of the particles we want to send to another rank.
vtkParticleTracerBaseNamespace::ParticleDataList ParticleHistories
ProtoPD is used just to keep track of the input array names and number of components for copy allocat...
vtkParticleTracerBase(const vtkParticleTracerBase &)=delete
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.
vtkNew< vtkFloatArray > ParticleVorticity
void SetInterpolatorType(int interpolatorType)
Set the type of the velocity field interpolator to determine whether INTERPOLATOR_WITH_DATASET_POINT_...
virtual void SetMeshOverTime(int meshOverTime)
bool SetTerminationTimeNoModify(double t)
int Finalize(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
Generates an output using the data provided after Execute was ran.
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[])
void TestParticles(vtkParticleTracerBaseNamespace::ParticleVector &candidates, vtkParticleTracerBaseNamespace::ParticleVector &passed)
inside our data.
virtual bool IsPointDataValid(vtkDataObject *input)
Methods that check that the input arrays are ordered the same on all data sets.
vtkTemporalInterpolatedVelocityField * GetInterpolator()
void SetIntegratorType(int type)
vtkNew< vtkIntArray > ParticleIds
vtkParticleTracerBaseNamespace::ParticleVector LocalSeeds
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.
unsigned int NumberOfParticles()
vtkNew< vtkSignedCharArray > ParticleSourceIds
virtual void InitializeExtraPointDataArrays(vtkPointData *vtkNotUsed(outputPD))
Methods to append values to existing point data arrays that may only be desired on specific concrete ...
void SetParticle(vtkParticleTracerBaseNamespace::ParticleInformation &info, vtkTemporalInterpolatedVelocityField *interpolator, vtkDoubleArray *cellVectors)
void EnqueueParticleToAnotherProcess(vtkParticleTracerBaseNamespace::ParticleInformation &)
vtkNew< vtkFloatArray > ParticleAngularVel
~vtkParticleTracerBase() override
static const double Epsilon
virtual void AssignSeedsToProcessors(double time, vtkDataSet *source, vtkParticleTracerBaseNamespace::ParticleVector &localSeedPoints)
all the injection/seed points according to which processor they belong to.
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.
void ResizeArrays(vtkIdType numTuples)
vtkNew< vtkFloatArray > ParticleAge
void SetInterpolatorTypeToDataSetPointLocator()
Set the velocity field interpolator type to one that uses a point locator to perform local spatial se...
vtkSmartPointer< vtkMultiProcessController > Controller
vtkSignedCharArray * GetParticleSourceIds(vtkPointData *)
vtkSetFilePathMacro(ParticleFileName)
Set/Get the filename to be used with the particle writer when dumping particles to disk.
virtual vtkMultiProcessController * GetController()
Get/Set the controller to use.
vtkAbstractParticleWriter * ParticleWriter
vtkSmartPointer< vtkDataSet > Seeds
All ranks have the same representation of the seeds.
int Execute(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
Moves the particles one time step further.
represent and manipulate point attribute data
represent and manipulate 3D points
Definition vtkPoints.h:139
Superclass for algorithms that produce only polydata as output.
concrete dataset represents vertices, lines, polygons, and triangle strips
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_4_0(reason)
#define vtkCreateWrappedTemporalAlgorithmInterface()
int vtkIdType
Definition vtkType.h:315