VTK  9.0.2
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 =========================================================================*/
88 #ifndef vtkLagrangianParticleTracker_h
89 #define vtkLagrangianParticleTracker_h
90 
91 #include "vtkBoundingBox.h" // For cached bounds
92 #include "vtkDataObjectAlgorithm.h"
93 #include "vtkFiltersFlowPathsModule.h" // For export macro
94 
95 #include <atomic> // for atomic
96 #include <mutex> // for mutexes
97 #include <queue> // for particle queue
98 
99 class vtkBoundingBox;
100 class vtkCellArray;
101 class vtkDataSet;
102 class vtkDoubleArray;
103 class vtkIdList;
104 class vtkInformation;
110 class vtkPointData;
111 class vtkPoints;
112 class vtkPolyData;
113 class vtkPolyLine;
114 struct IntegratingFunctor;
115 
116 class VTKFILTERSFLOWPATHS_EXPORT vtkLagrangianParticleTracker : public vtkDataObjectAlgorithm
117 {
118 public:
120  void PrintSelf(ostream& os, vtkIndent indent) override;
122 
124  {
125  STEP_LAST_CELL_LENGTH = 0,
126  STEP_CUR_CELL_LENGTH = 1,
127  STEP_LAST_CELL_VEL_DIR = 2,
128  STEP_CUR_CELL_VEL_DIR = 3,
129  STEP_LAST_CELL_DIV_THEO = 4,
130  STEP_CUR_CELL_DIV_THEO = 5
131  } CellLengthComputation;
132 
134 
139  vtkGetObjectMacro(IntegrationModel, vtkLagrangianBasicIntegrationModel);
141 
143 
148  vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
150 
152 
157  vtkSetMacro(GeneratePolyVertexInteractionOutput, bool);
158  vtkGetMacro(GeneratePolyVertexInteractionOutput, bool);
160 
162 
185  vtkSetMacro(CellLengthComputationMode, int);
186  vtkGetMacro(CellLengthComputationMode, int);
188 
190 
193  vtkSetMacro(StepFactor, double);
194  vtkGetMacro(StepFactor, double);
196 
198 
201  vtkSetMacro(StepFactorMin, double);
202  vtkGetMacro(StepFactorMin, double);
204 
206 
209  vtkSetMacro(StepFactorMax, double);
210  vtkGetMacro(StepFactorMax, double);
212 
214 
217  vtkSetMacro(MaximumNumberOfSteps, int);
218  vtkGetMacro(MaximumNumberOfSteps, int);
220 
222 
226  vtkSetMacro(MaximumIntegrationTime, double);
227  vtkGetMacro(MaximumIntegrationTime, double);
229 
231 
237  vtkSetMacro(AdaptiveStepReintegration, bool);
238  vtkGetMacro(AdaptiveStepReintegration, bool);
239  vtkBooleanMacro(AdaptiveStepReintegration, bool);
241 
243 
247  vtkSetMacro(GenerateParticlePathsOutput, bool);
248  vtkGetMacro(GenerateParticlePathsOutput, bool);
249  vtkBooleanMacro(GenerateParticlePathsOutput, bool);
251 
253 
262 
267 
269 
278 
283 
288 
293 
298 
302  int RequestData(vtkInformation* request, vtkInformationVector** inputVector,
303  vtkInformationVector* outputVector) override;
304 
309  vtkMTimeType GetMTime() override;
310 
316 
317 protected:
320 
321  virtual bool InitializeFlow(vtkDataObject* flow, vtkBoundingBox* bounds);
322  virtual bool InitializeParticles(const vtkBoundingBox* bounds, vtkDataSet* seeds,
323  std::queue<vtkLagrangianParticle*>& particles, vtkPointData* seedData);
324  virtual void GenerateParticles(const vtkBoundingBox* bounds, vtkDataSet* seeds,
325  vtkDataArray* initialVelocities, vtkDataArray* initialIntegrationTimes, vtkPointData* seedData,
326  int nVar, std::queue<vtkLagrangianParticle*>& particles);
327  virtual bool UpdateSurfaceCacheIfNeeded(vtkDataObject*& surfaces);
328  virtual void InitializeSurface(vtkDataObject*& surfaces);
329 
333  virtual bool InitializePathsOutput(
334  vtkPointData* seedData, vtkIdType numberOfSeeds, vtkPolyData*& particlePathsOutput);
335 
340  vtkPointData* seedData, vtkDataObject* surfaces, vtkDataObject*& interractionOutput);
341 
342  virtual bool FinalizeOutputs(vtkPolyData* particlePathsOutput, vtkDataObject* interactionOutput);
343 
344  static void InsertPolyVertexCell(vtkPolyData* polydata);
345  static void InsertVertexCells(vtkPolyData* polydata);
346 
347  virtual void GetParticleFeed(std::queue<vtkLagrangianParticle*>& particleQueue);
348 
353  std::queue<vtkLagrangianParticle*>&, vtkPolyData* particlePathsOutput,
354  vtkPolyLine* particlePath, vtkDataObject* interactionOutput);
355 
359  void InsertPathOutputPoint(vtkLagrangianParticle* particle, vtkPolyData* particlePathsOutput,
360  vtkIdList* particlePathPointId, bool prev = false);
361 
366  unsigned int interactedSurfaceFlatIndex, vtkDataObject* interactionOutput);
367 
369 
373  bool ComputeNextStep(vtkInitialValueProblemSolver* integrator, double* xprev, double* xnext,
374  double t, double& delT, double& delTActual, double minStep, double maxStep, double cellLength,
375  int& integrationRes, vtkLagrangianParticle* particle);
376 
379 
381  double StepFactor;
387  bool GenerateParticlePathsOutput = true;
389  std::atomic<vtkIdType> ParticleCounter;
390  std::atomic<vtkIdType> IntegratedParticleCounter;
393 
394  // internal parameters use for step computation
397 
398  // Cache related parameters
404 
405  std::mutex ProgressMutex;
406  friend struct IntegratingFunctor;
407 
408 private:
410  void operator=(const vtkLagrangianParticleTracker&) = delete;
411 };
412 
413 #endif
Proxy object to connect input/output ports.
Fast, simple class for dealing with 3D bounds.
object to represent cell connectivity
Definition: vtkCellArray.h:180
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:50
Superclass for algorithms that produce only data object as output.
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
list of point or cell ids
Definition: vtkIdList.h:31
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.
vtkFunctionSet abstract implementation to be used in the vtkLagrangianParticleTracker integrator.
Filter to inject and track particles in a flow.
vtkDataObject * GetSurface()
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 GenerateParticles(const vtkBoundingBox *bounds, vtkDataSet *seeds, vtkDataArray *initialVelocities, vtkDataArray *initialIntegrationTimes, vtkPointData *seedData, int nVar, std::queue< vtkLagrangianParticle * > &particles)
static vtkLagrangianParticleTracker * New()
virtual bool InitializeParticles(const vtkBoundingBox *bounds, vtkDataSet *seeds, std::queue< vtkLagrangianParticle * > &particles, vtkPointData *seedData)
static void InsertPolyVertexCell(vtkPolyData *polydata)
double ComputeCellLength(vtkLagrangianParticle *particle)
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.
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)
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.
vtkDataObject * GetSource()
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:32
represent and manipulate 3D points
Definition: vtkPoints.h:34
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:85
cell represents a set of 1D lines
Definition: vtkPolyLine.h:37
@ info
Definition: vtkX3D.h:382
@ port
Definition: vtkX3D.h:453
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
int vtkIdType
Definition: vtkType.h:338
vtkTypeUInt32 vtkMTimeType
Definition: vtkType.h:293