VTK
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 =========================================================================*/
87 #ifndef vtkLagrangianParticleTracker_h
88 #define vtkLagrangianParticleTracker_h
89 
90 #include "vtkFiltersFlowPathsModule.h" // For export macro
91 #include "vtkDataObjectAlgorithm.h"
92 #include "vtkBoundingBox.h" // For cached bounds
93 
94 #include <queue> // for particle queue
95 
96 class vtkBoundingBox;
97 class vtkCellArray;
98 class vtkDataSet;
99 class vtkDoubleArray;
100 class vtkIdList;
101 class vtkInformation;
105 class vtkPointData;
106 class vtkPoints;
107 class vtkPolyData;
108 
109 class VTKFILTERSFLOWPATHS_EXPORT vtkLagrangianParticleTracker :
111 {
112 public:
113 
115  void PrintSelf(ostream& os, vtkIndent indent) override;
117 
118  typedef enum CellLengthComputation{
119  STEP_LAST_CELL_LENGTH = 0,
120  STEP_CUR_CELL_LENGTH = 1,
121  STEP_LAST_CELL_VEL_DIR = 2,
122  STEP_CUR_CELL_VEL_DIR = 3,
123  STEP_LAST_CELL_DIV_THEO = 4,
124  STEP_CUR_CELL_DIV_THEO = 5
125  } CellLengthComputation;
126 
128 
131  void SetIntegrationModel(vtkLagrangianBasicIntegrationModel* integrationModel);
132  vtkGetObjectMacro(IntegrationModel, vtkLagrangianBasicIntegrationModel);
134 
136 
139  void SetIntegrator(vtkInitialValueProblemSolver* integrator);
140  vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
142 
144 
149  vtkSetMacro(GeneratePolyVertexInteractionOutput, bool);
150  vtkGetMacro(GeneratePolyVertexInteractionOutput, bool);
152 
154 
178  vtkSetMacro(CellLengthComputationMode, int);
179  vtkGetMacro(CellLengthComputationMode, int);
181 
183 
186  vtkSetMacro(StepFactor, double);
187  vtkGetMacro(StepFactor, double);
189 
191 
194  vtkSetMacro(StepFactorMin, double);
195  vtkGetMacro(StepFactorMin, double);
197 
199 
202  vtkSetMacro(StepFactorMax, double);
203  vtkGetMacro(StepFactorMax, double);
205 
207 
210  vtkSetMacro(MaximumNumberOfSteps, int);
211  vtkGetMacro(MaximumNumberOfSteps, int);
213 
215 
220  vtkSetMacro(AdaptiveStepReintegration, bool);
221  vtkGetMacro(AdaptiveStepReintegration, bool);
222  vtkBooleanMacro(AdaptiveStepReintegration, bool);
224 
226 
229  vtkSetMacro(UseParticlePathsRenderingThreshold, bool);
230  vtkGetMacro(UseParticlePathsRenderingThreshold, bool);
231  vtkBooleanMacro(UseParticlePathsRenderingThreshold, bool);
233 
235 
238  vtkSetMacro(ParticlePathsRenderingPointsThreshold, int);
239  vtkGetMacro(ParticlePathsRenderingPointsThreshold, int);
241 
243 
246  vtkSetMacro(CreateOutOfDomainParticle, bool);
247  vtkGetMacro(CreateOutOfDomainParticle, bool);
248  vtkBooleanMacro(CreateOutOfDomainParticle, bool);
250 
252 
258  void SetSourceData(vtkDataObject* source);
259  vtkDataObject* GetSource();
261 
265  void SetSourceConnection(vtkAlgorithmOutput* algOutput);
266 
268 
274  void SetSurfaceData(vtkDataObject *source);
275  vtkDataObject *GetSurface();
277 
281  void SetSurfaceConnection(vtkAlgorithmOutput* algOutput);
282 
286  int FillInputPortInformation(int port, vtkInformation* info) override;
287 
292 
298  vtkInformationVector*) override;
299 
303  int RequestData(vtkInformation *request,
304  vtkInformationVector **inputVector,
305  vtkInformationVector *outputVector) override;
306 
311  vtkMTimeType GetMTime() override;
312 
316  virtual vtkIdType GetNewParticleId();
317 
318 protected:
320  ~vtkLagrangianParticleTracker() override;
321 
322  virtual bool InitializeInputs(vtkInformationVector **inputVector,
323  vtkDataObject*& flow, vtkDataObject*& seeds, vtkDataObject*& surfaces,
324  std::queue<vtkLagrangianParticle*>& particleQueue, vtkPointData* seedData);
325  virtual bool InitializeFlow(vtkDataObject* flow, vtkBoundingBox* bounds);
326  virtual bool InitializeParticles(const vtkBoundingBox* bounds, vtkDataObject* seeds,
327  std::queue<vtkLagrangianParticle*>& particles, vtkPointData* seedData);
328  virtual void GenerateParticles(const vtkBoundingBox* bounds, vtkDataSet* seeds,
329  vtkDataArray* initialVelocities, vtkDataArray* initialIntegrationTimes,
330  vtkPointData* seedData, int nVar, std::queue<vtkLagrangianParticle*>& particles);
331  virtual bool UpdateSurfaceCacheIfNeeded(vtkDataObject*& surfaces);
332  virtual void InitializeSurface(vtkDataObject*& surfaces);
333  virtual bool InitializeOutputs(vtkInformationVector *outputVector, vtkPointData* seedData,
334  vtkIdType numberOfSeeds, vtkDataObject* surfaces,
335  vtkPolyData*& particlePathsOutput, vtkDataObject*& interactionOutput);
336 
337  virtual bool InitializePathsOutput(vtkInformationVector *outputVector,
338  vtkPointData* seedData, vtkIdType numberOfSeeds,
339  vtkPolyData*& particlePathsOutput);
340 
341  virtual bool InitializeInteractionOutput(vtkInformationVector *outputVector,
342  vtkPointData* seedData, vtkDataObject* surfaces, vtkDataObject*& interractionOutput);
343 
344  virtual void InitializeParticleData(vtkFieldData* particleData, int maxTuples = 0);
345  virtual void InitializePathData(vtkFieldData* data);
346  virtual void InitializeInteractionData(vtkFieldData* data);
347 
348  virtual bool FinalizeOutputs(vtkPolyData* particlePathsOutput,
349  vtkDataObject* interractionOutput);
350 
351  static void InsertPolyVertexCell(vtkPolyData* polydata);
352  static void InsertVertexCells(vtkPolyData* polydata);
353 
354  virtual void GetParticleFeed(std::queue<vtkLagrangianParticle*>& particleQueue);
355 
356  virtual int Integrate(vtkLagrangianParticle*, std::queue<vtkLagrangianParticle*>&,
357  vtkPolyData* particlePathsOutput, vtkIdList* particlePathPointId,
358  vtkDataObject* interactionOutput);
359 
360  void InsertPathOutputPoint(vtkLagrangianParticle* particle,
361  vtkPolyData* particlePathsOutput, vtkIdList* particlePathPointId,
362  bool prev = false);
363 
364  void InsertInteractionOutputPoint(vtkLagrangianParticle* particle,
365  unsigned int interactedSurfaceFlatIndex, vtkDataObject* interactionOutput);
366 
367  void InsertSeedData(vtkLagrangianParticle* particle, vtkFieldData* data);
368  void InsertPathData(vtkLagrangianParticle* particle, vtkFieldData* data);
369  void InsertInteractionData(vtkLagrangianParticle* particle, vtkFieldData* data);
370  void InsertParticleData(vtkLagrangianParticle* particle, vtkFieldData* data, int stepEnum);
371 
372  double ComputeCellLength(vtkLagrangianParticle* particle);
373 
374  bool ComputeNextStep(
375  double* xprev, double* xnext,
376  double t, double& delT, double& delTActual,
377  double minStep, double maxStep,
378  int& integrationRes);
379 
380  virtual bool CheckParticlePathsRenderingThreshold(vtkPolyData* particlePathsOutput);
381 
384 
386  double StepFactor;
396 
397  // internal parameters use for step computation
400 
401  // Cache related parameters
407 
408 private:
410  void operator=(const vtkLagrangianParticleTracker&) = delete;
411 };
412 
413 #endif
virtual int RequestDataObject(vtkInformation *, vtkInformationVector **, vtkInformationVector *)
This is called by the superclass.
represent and manipulate point attribute data
Definition: vtkPointData.h:31
Store vtkAlgorithm input/output information.
vtkTypeUInt32 vtkMTimeType
Definition: vtkType.h:300
abstract class to specify dataset behavior
Definition: vtkDataSet.h:56
static vtkDataObjectAlgorithm * New()
vtkInitialValueProblemSolver * Integrator
int vtkIdType
Definition: vtkType.h:345
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:79
Proxy object to connect input/output ports.
dynamic, self-adjusting array of double
int FillOutputPortInformation(int port, vtkInformation *info) override
Fill the output port information objects for this algorithm.
a simple class to control print indentation
Definition: vtkIndent.h:33
vtkFunctionSet abstract implementation to be used in the vtkLagrangianParticleTracker integrator...
list of point or cell ids
Definition: vtkIdList.h:30
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:48
virtual vtkMTimeType GetMTime()
Return this object&#39;s modified time.
Basis class for Lagrangian particles.
boost::graph_traits< vtkGraph *>::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
Superclass for algorithms that produce only data object as output.
object to represent cell connectivity
Definition: vtkCellArray.h:44
virtual int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *)
vtkLagrangianBasicIntegrationModel * IntegrationModel
Store zero or more vtkInformation instances.
int FillInputPortInformation(int port, vtkInformation *info) override
Fill the input port information objects for this algorithm.
general representation of visualization data
Definition: vtkDataObject.h:58
Filter to inject and track particles in a flow.
represent and manipulate 3D points
Definition: vtkPoints.h:33
Fast Simple Class for dealing with 3D bounds.
represent and manipulate fields of data
Definition: vtkFieldData.h:53
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
Integrate a set of ordinary differential equations (initial value problem) in time.