VTK  9.2.5
vtkLagrangianParticleTracker.h
Go to the documentation of this file.
1/*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkLagrangianParticleTracker.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=========================================================================*/
91#ifndef vtkLagrangianParticleTracker_h
92#define vtkLagrangianParticleTracker_h
93
94#include "vtkBoundingBox.h" // For cached bounds
96#include "vtkFiltersFlowPathsModule.h" // For export macro
97
98#include <atomic> // for atomic
99#include <mutex> // for mutexes
100#include <queue> // for particle queue
101
102class vtkBoundingBox;
103class vtkCellArray;
104class vtkDataSet;
105class vtkDoubleArray;
106class vtkIdList;
107class vtkInformation;
113class vtkPointData;
114class vtkPoints;
115class vtkPolyData;
116class vtkPolyLine;
117struct IntegratingFunctor;
119
120class VTKFILTERSFLOWPATHS_EXPORT vtkLagrangianParticleTracker : public vtkDataObjectAlgorithm
121{
122public:
124 void PrintSelf(ostream& os, vtkIndent indent) override;
126
128 {
129 STEP_CUR_CELL_LENGTH = 1,
130 STEP_CUR_CELL_VEL_DIR = 3,
131 STEP_CUR_CELL_DIV_THEO = 5
132 } CellLengthComputation;
133
135
140 vtkGetObjectMacro(IntegrationModel, vtkLagrangianBasicIntegrationModel);
142
144
149 vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
151
153
158 vtkSetMacro(GeneratePolyVertexInteractionOutput, bool);
159 vtkGetMacro(GeneratePolyVertexInteractionOutput, bool);
161
163
176 vtkSetMacro(CellLengthComputationMode, int);
177 vtkGetMacro(CellLengthComputationMode, int);
179
181
184 vtkSetMacro(StepFactor, double);
185 vtkGetMacro(StepFactor, double);
187
189
192 vtkSetMacro(StepFactorMin, double);
193 vtkGetMacro(StepFactorMin, double);
195
197
200 vtkSetMacro(StepFactorMax, double);
201 vtkGetMacro(StepFactorMax, double);
203
205
208 vtkSetMacro(MaximumNumberOfSteps, int);
209 vtkGetMacro(MaximumNumberOfSteps, int);
211
213
217 vtkSetMacro(MaximumIntegrationTime, double);
218 vtkGetMacro(MaximumIntegrationTime, double);
220
222
228 vtkSetMacro(AdaptiveStepReintegration, bool);
229 vtkGetMacro(AdaptiveStepReintegration, bool);
230 vtkBooleanMacro(AdaptiveStepReintegration, bool);
232
234
238 vtkSetMacro(GenerateParticlePathsOutput, bool);
239 vtkGetMacro(GenerateParticlePathsOutput, bool);
240 vtkBooleanMacro(GenerateParticlePathsOutput, bool);
242
244
253
258
260
269
274
278 int FillInputPortInformation(int port, vtkInformation* info) override;
279
283 int FillOutputPortInformation(int port, vtkInformation* info) override;
284
289
294 vtkInformationVector* outputVector) override;
295
301
307
308protected:
311
312 virtual bool InitializeFlow(vtkDataObject* flow, vtkBoundingBox* bounds);
313 virtual bool InitializeParticles(const vtkBoundingBox* bounds, vtkDataSet* seeds,
314 std::queue<vtkLagrangianParticle*>& particles, vtkPointData* seedData);
315 virtual void GenerateParticles(const vtkBoundingBox* bounds, vtkDataSet* seeds,
316 vtkDataArray* initialVelocities, vtkDataArray* initialIntegrationTimes, vtkPointData* seedData,
317 int nVar, std::queue<vtkLagrangianParticle*>& particles);
319 virtual void InitializeSurface(vtkDataObject*& surfaces);
320
325 vtkPointData* seedData, vtkIdType numberOfSeeds, vtkPolyData*& particlePathsOutput);
326
331 vtkPointData* seedData, vtkDataObject* surfaces, vtkDataObject*& interractionOutput);
332
333 virtual bool FinalizeOutputs(vtkPolyData* particlePathsOutput, vtkDataObject* interactionOutput);
334
335 static void InsertPolyVertexCell(vtkPolyData* polydata);
336 static void InsertVertexCells(vtkPolyData* polydata);
337
338 virtual void GetParticleFeed(std::queue<vtkLagrangianParticle*>& particleQueue);
339
344 std::queue<vtkLagrangianParticle*>&, vtkPolyData* particlePathsOutput,
345 vtkPolyLine* particlePath, vtkDataObject* interactionOutput);
346
350 void InsertPathOutputPoint(vtkLagrangianParticle* particle, vtkPolyData* particlePathsOutput,
351 vtkIdList* particlePathPointId, bool prev = false);
352
357 unsigned int interactedSurfaceFlatIndex, vtkDataObject* interactionOutput);
358
364
368 bool ComputeNextStep(vtkInitialValueProblemSolver* integrator, double* xprev, double* xnext,
369 double t, double& delT, double& delTActual, double minStep, double maxStep, double cellLength,
370 int& integrationRes, vtkLagrangianParticle* particle);
371
376 virtual void DeleteParticle(vtkLagrangianParticle* particle);
377
380
388 bool GenerateParticlePathsOutput = true;
390 std::atomic<vtkIdType> ParticleCounter;
391 std::atomic<vtkIdType> IntegratedParticleCounter;
394
395 // internal parameters use for step computation
398
399 // Cache related parameters
403 bool FlowCacheInvalid = true;
406 bool SurfaceCacheInvalid = true;
407
408 std::mutex ProgressMutex;
409 friend struct IntegratingFunctor;
410
412
413private:
415 void operator=(const vtkLagrangianParticleTracker&) = delete;
416};
417
418#endif
Proxy object to connect input/output ports.
Fast, simple class for representing and operating on 3D bounds.
object to represent cell connectivity
Definition: vtkCellArray.h:187
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:56
Superclass for algorithms that produce only data object as output.
general representation of visualization data
Definition: vtkDataObject.h:66
abstract class to specify dataset behavior
Definition: vtkDataSet.h:63
dynamic, self-adjusting array of double
list of point or cell ids
Definition: vtkIdList.h:34
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.
vtkFunctionSet abstract implementation to be used in the vtkLagrangianParticleTracker integrator.
Filter to inject and track particles in a flow.
void SetSourceConnection(vtkAlgorithmOutput *algOutput)
Specify the source object used to generate particle initial position (seeds).
vtkLagrangianBasicIntegrationModel * IntegrationModel
virtual bool UpdateSurfaceCacheIfNeeded(vtkDataObject *&surfaces)
virtual vtkIdType GetNewParticleId()
Get an unique id for a particle This method is thread safe.
void InsertPathOutputPoint(vtkLagrangianParticle *particle, vtkPolyData *particlePathsOutput, vtkIdList *particlePathPointId, bool prev=false)
This method is thread safe.
vtkInitialValueProblemSolver * Integrator
int FillInputPortInformation(int port, vtkInformation *info) override
Declare input port type.
virtual void GetParticleFeed(std::queue< vtkLagrangianParticle * > &particleQueue)
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
Process inputs to integrate particle and generate output.
bool ComputeNextStep(vtkInitialValueProblemSolver *integrator, double *xprev, double *xnext, double t, double &delT, double &delTActual, double minStep, double maxStep, double cellLength, int &integrationRes, vtkLagrangianParticle *particle)
This method is thread safe.
virtual void InitializeSurface(vtkDataObject *&surfaces)
virtual void DeleteParticle(vtkLagrangianParticle *particle)
This method is thread safe Call the ParticleAboutToBeDeleted model method and delete the particle.
virtual void GenerateParticles(const vtkBoundingBox *bounds, vtkDataSet *seeds, vtkDataArray *initialVelocities, vtkDataArray *initialIntegrationTimes, vtkPointData *seedData, int nVar, std::queue< vtkLagrangianParticle * > &particles)
virtual bool InitializeParticles(const vtkBoundingBox *bounds, vtkDataSet *seeds, std::queue< vtkLagrangianParticle * > &particles, vtkPointData *seedData)
vtkDataObject * GetSource()
Specify the source object used to generate particle initial position (seeds).
static void InsertPolyVertexCell(vtkPolyData *polydata)
double ComputeCellLength(vtkLagrangianParticle *particle)
Computes the cell length for the next step using the method set by CellLengthComputationMode.
void SetIntegrator(vtkInitialValueProblemSolver *integrator)
Set/Get the integrator.
void SetSurfaceData(vtkDataObject *source)
Specify the source object used to compute surface interaction with Note that this method does not con...
static void InsertVertexCells(vtkPolyData *polydata)
vtkMTimeType GetMTime() override
Get the tracker modified time taking into account the integration model and the integrator.
std::atomic< vtkIdType > IntegratedParticleCounter
void SetSurfaceConnection(vtkAlgorithmOutput *algOutput)
Specify the object used to compute surface interaction with.
vtkLagrangianThreadedData * SerialThreadedData
vtkDataObject * GetSurface()
Specify the source object used to compute surface interaction with Note that this method does not con...
void SetIntegrationModel(vtkLagrangianBasicIntegrationModel *integrationModel)
Set/Get the integration model.
virtual int Integrate(vtkInitialValueProblemSolver *integrator, vtkLagrangianParticle *, std::queue< vtkLagrangianParticle * > &, vtkPolyData *particlePathsOutput, vtkPolyLine *particlePath, vtkDataObject *interactionOutput)
This method is thread safe.
virtual bool InitializeFlow(vtkDataObject *flow, vtkBoundingBox *bounds)
static vtkLagrangianParticleTracker * New()
int FillOutputPortInformation(int port, vtkInformation *info) override
Declare output port type.
int RequestDataObject(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
Create outputs objects.
virtual bool InitializePathsOutput(vtkPointData *seedData, vtkIdType numberOfSeeds, vtkPolyData *&particlePathsOutput)
This method is thread safe.
virtual bool InitializeInteractionOutput(vtkPointData *seedData, vtkDataObject *surfaces, vtkDataObject *&interractionOutput)
This method is thread safe.
void SetSourceData(vtkDataObject *source)
Specify the source object used to generate particle initial position (seeds).
void InsertInteractionOutputPoint(vtkLagrangianParticle *particle, unsigned int interactedSurfaceFlatIndex, vtkDataObject *interactionOutput)
This method is thread safe.
virtual bool FinalizeOutputs(vtkPolyData *particlePathsOutput, vtkDataObject *interactionOutput)
~vtkLagrangianParticleTracker() override
Basis class for Lagrangian particles.
Composite dataset that organizes datasets into blocks.
composite dataset to encapsulates pieces of dataset.
represent and manipulate point attribute data
Definition: vtkPointData.h:42
represent and manipulate 3D points
Definition: vtkPoints.h:40
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:91
cell represents a set of 1D lines
Definition: vtkPolyLine.h:40
struct to hold a user data
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
int vtkIdType
Definition: vtkType.h:332
vtkTypeUInt32 vtkMTimeType
Definition: vtkType.h:287