VTK  9.3.1
vtkStreamTracer.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
82 #ifndef vtkStreamTracer_h
83 #define vtkStreamTracer_h
84 
85 #include "vtkFiltersFlowPathsModule.h" // For export macro
86 #include "vtkPolyDataAlgorithm.h"
87 
88 #include "vtkDataSetAttributesFieldList.h" // Needed to identify common data arrays
89 #include "vtkInitialValueProblemSolver.h" // Needed for constants
90 
91 VTK_ABI_NAMESPACE_BEGIN
94 class vtkDataArray;
96 class vtkDoubleArray;
97 class vtkExecutive;
98 class vtkGenericCell;
99 class vtkIdList;
100 class vtkIntArray;
101 class vtkPoints;
102 
103 VTK_ABI_NAMESPACE_END
104 #include <vector> // for std::vector
105 
106 // Helper struct to convert between different length scales.
107 VTK_ABI_NAMESPACE_BEGIN
108 struct VTKFILTERSFLOWPATHS_EXPORT vtkIntervalInformation
109 {
110  double Interval;
111  int Unit;
112 
113  static double ConvertToLength(double interval, int unit, double cellLength);
114  static double ConvertToLength(vtkIntervalInformation& interval, double cellLength);
115 };
116 
128  void* clientdata, vtkPoints* points, vtkDataArray* velocity, int integrationDirection);
129 
130 class VTKFILTERSFLOWPATHS_EXPORT vtkStreamTracer : public vtkPolyDataAlgorithm
131 {
132 public:
140  static vtkStreamTracer* New();
141 
143 
147  void PrintSelf(ostream& os, vtkIndent indent) override;
149 
151 
156  vtkSetVector3Macro(StartPosition, double);
157  vtkGetVector3Macro(StartPosition, double);
159 
161 
167  void SetSourceData(vtkDataSet* source);
168  vtkDataSet* GetSource();
170 
176  void SetSourceConnection(vtkAlgorithmOutput* algOutput);
177 
178  // The previously-supported TIME_UNIT is excluded in this current
179  // enumeration definition because the underlying step size is ALWAYS in
180  // arc length unit (LENGTH_UNIT) while the 'real' time interval (virtual
181  // for steady flows) that a particle actually takes to trave in a single
182  // step is obtained by dividing the arc length by the LOCAL speed. The
183  // overall elapsed time (i.e., the life span) of the particle is the sum
184  // of those individual step-wise time intervals. The arc-length-to-time
185  // conversion only occurs for vorticity computation and for generating a
186  // point data array named 'IntegrationTime'.
187  enum Units
188  {
189  LENGTH_UNIT = 1,
190  CELL_LENGTH_UNIT = 2
191  };
192 
193  enum Solvers
194  {
200  };
201 
203  {
207  OUT_OF_LENGTH = 4,
208  OUT_OF_STEPS = 5,
209  STAGNATION = 6,
210  FIXED_REASONS_FOR_TERMINATION_COUNT
211  };
212 
214 
224  void SetIntegrator(vtkInitialValueProblemSolver*);
225  vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
226  void SetIntegratorType(int type);
227  int GetIntegratorType();
228  void SetIntegratorTypeToRungeKutta2() { this->SetIntegratorType(RUNGE_KUTTA2); }
229  void SetIntegratorTypeToRungeKutta4() { this->SetIntegratorType(RUNGE_KUTTA4); }
230  void SetIntegratorTypeToRungeKutta45() { this->SetIntegratorType(RUNGE_KUTTA45); }
232 
242  void SetInterpolatorTypeToDataSetPointLocator();
243 
250  void SetInterpolatorTypeToCellLocator();
251 
253 
256  vtkSetMacro(MaximumPropagation, double);
257  vtkGetMacro(MaximumPropagation, double);
259 
266  void SetIntegrationStepUnit(int unit);
267  int GetIntegrationStepUnit() { return this->IntegrationStepUnit; }
268 
270 
277  vtkSetMacro(InitialIntegrationStep, double);
278  vtkGetMacro(InitialIntegrationStep, double);
280 
282 
288  vtkSetMacro(MinimumIntegrationStep, double);
289  vtkGetMacro(MinimumIntegrationStep, double);
291 
293 
299  vtkSetMacro(MaximumIntegrationStep, double);
300  vtkGetMacro(MaximumIntegrationStep, double);
302 
304 
307  vtkSetMacro(MaximumError, double);
308  vtkGetMacro(MaximumError, double);
310 
312 
320  vtkSetMacro(MaximumNumberOfSteps, vtkIdType);
321  vtkGetMacro(MaximumNumberOfSteps, vtkIdType);
323 
325 
329  vtkSetMacro(TerminalSpeed, double);
330  vtkGetMacro(TerminalSpeed, double);
332 
334 
338  vtkGetMacro(SurfaceStreamlines, bool);
339  vtkSetMacro(SurfaceStreamlines, bool);
340  vtkBooleanMacro(SurfaceStreamlines, bool);
342 
343  enum
344  {
347  BOTH
348  };
349 
350  enum
351  {
353  INTERPOLATOR_WITH_CELL_LOCATOR
354  };
355 
357 
364  vtkSetClampMacro(IntegrationDirection, int, FORWARD, BOTH);
365  vtkGetMacro(IntegrationDirection, int);
366  void SetIntegrationDirectionToForward() { this->SetIntegrationDirection(FORWARD); }
367  void SetIntegrationDirectionToBackward() { this->SetIntegrationDirection(BACKWARD); }
368  void SetIntegrationDirectionToBoth() { this->SetIntegrationDirection(BOTH); }
370 
372 
377  vtkSetMacro(ComputeVorticity, bool);
378  vtkGetMacro(ComputeVorticity, bool);
380 
382 
386  vtkSetMacro(RotationScale, double);
387  vtkGetMacro(RotationScale, double);
389 
399  void SetInterpolatorPrototype(vtkAbstractInterpolatedVelocityField* ivf);
400 
410  void SetInterpolatorType(int interpType);
411 
413 
417  vtkGetMacro(ForceSerialExecution, bool);
418  vtkSetMacro(ForceSerialExecution, bool);
419  vtkBooleanMacro(ForceSerialExecution, bool);
421 
430  void AddCustomTerminationCallback(
431  CustomTerminationCallbackType callback, void* clientdata, int reasonForTermination);
432 
441  void ConvertIntervals(
442  double& step, double& minStep, double& maxStep, int direction, double cellLength);
443 
445 
449  void GenerateNormals(vtkPolyData* output, double* firstNormal, const char* vecName);
450  void CalculateVorticity(
451  vtkGenericCell* cell, double pcoords[3], vtkDoubleArray* cellVectors, double vorticity[3]);
453 
455 
465  vtkSetMacro(UseLocalSeedSource, bool);
466  vtkGetMacro(UseLocalSeedSource, bool);
467  vtkBooleanMacro(UseLocalSeedSource, bool);
469 
470 protected:
471  vtkStreamTracer();
472  ~vtkStreamTracer() override;
473 
474  // Create a default executive.
476 
477  // hide the superclass' AddInput() from the user and the compiler
479  {
480  vtkErrorMacro(<< "AddInput() must be called with a vtkDataSet not a vtkDataObject.");
481  }
482 
484  int FillInputPortInformation(int, vtkInformation*) override;
485 
486  void Integrate(vtkPointData* inputData, vtkPolyData* output, vtkDataArray* seedSource,
487  vtkIdList* seedIds, vtkIntArray* integrationDirections,
488  vtkAbstractInterpolatedVelocityField* func, int maxCellSize, int vecType,
489  const char* vecFieldName, double& propagation, vtkIdType& numSteps, double& integrationTime,
490  std::vector<CustomTerminationCallbackType>& customTerminationCallback,
491  std::vector<void*>& customTerminationClientData, std::vector<int>& customReasonForTermination);
492 
493  double SimpleIntegrate(double seed[3], double lastPoint[3], double stepSize,
495  int CheckInputs(vtkAbstractInterpolatedVelocityField*& func, int* maxCellSize);
496 
498 
499  // starting from global x-y-z position
500  double StartPosition[3];
501 
502  static const double EPSILON;
504 
505  // Used by subclasses, leave alone
507 
512 
513  int SetupOutput(vtkInformation* inInfo, vtkInformation* outInfo);
514  void InitializeSeeds(vtkDataArray*& seeds, vtkIdList*& seedIds,
515  vtkIntArray*& integrationDirections, vtkDataSet* source);
516 
519 
520  // Prototype showing the integrator type to be set by the user.
522 
523  double MaximumError;
525 
528 
529  // Compute streamlines only on surface.
531 
533 
534  // These are used to manage complex input types such as
535  // multiblock / composite datasets. Basically the filter input is
536  // converted to a composite dataset, and the point data attributes
537  // are intersected to produce a common set of output data arrays.
538  vtkCompositeDataSet* InputData; // convert input data to composite dataset
539  vtkDataSetAttributesFieldList InputPD; // intersect attributes of all datasets
540  bool
541  HasMatchingPointAttributes; // does the point data in the multiblocks have the same attributes?
542 
543  // Control execution as serial or threaded
545  bool SerialExecution; // internal use to combine information
546 
547  std::vector<CustomTerminationCallbackType> CustomTerminationCallback;
548  std::vector<void*> CustomTerminationClientData;
549  std::vector<int> CustomReasonForTermination;
550 
551  // Only relevant for this derived parallel version of vtkStreamTracer,
552  // but needs to be defined in this class to have a uniform interface
553  // between this class and the parallel override vtkPStreamTracer
555 
556  friend class PStreamTracerUtils;
557 
558 private:
559  vtkStreamTracer(const vtkStreamTracer&) = delete;
560  void operator=(const vtkStreamTracer&) = delete;
561 };
562 
563 VTK_ABI_NAMESPACE_END
564 #endif
void SetIntegrationDirectionToBoth()
Specify whether the streamline is integrated in the upstream or downstream direction, or in both directions.
static const double EPSILON
represent and manipulate point attribute data
Definition: vtkPointData.h:29
virtual vtkExecutive * CreateDefaultExecutive()
Create a default executive.
Store vtkAlgorithm input/output information.
abstract class to specify dataset behavior
Definition: vtkDataSet.h:52
virtual int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
This is called by the superclass.
helps manage arrays from multiple vtkDataSetAttributes.
An abstract class for obtaining the interpolated velocity values at a point.
void AddInput(vtkDataObject *)
bool GenerateNormalsInIntegrate
vtkCompositeDataSet * InputData
int vtkIdType
Definition: vtkType.h:315
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:79
Superclass for all pipeline executives in VTK.
Definition: vtkExecutive.h:38
double MinimumIntegrationStep
provides thread-safe access to cells
vtkIdType MaximumNumberOfSteps
Proxy object to connect input/output ports.
dynamic, self-adjusting array of double
static vtkPolyDataAlgorithm * New()
void SetIntegrationDirectionToForward()
Specify whether the streamline is integrated in the upstream or downstream direction, or in both directions.
std::vector< int > CustomReasonForTermination
dynamic, self-adjusting array of int
Definition: vtkIntArray.h:34
abstract superclass for composite (multi-block or AMR) datasets
Superclass for algorithms that produce only polydata as output.
a simple class to control print indentation
Definition: vtkIndent.h:28
int GetIntegrationStepUnit()
list of point or cell ids
Definition: vtkIdList.h:22
bool(* CustomTerminationCallbackType)(void *clientdata, vtkPoints *points, vtkDataArray *velocity, int integrationDirection)
Used to specify custom conditions which are evaluated to determine whether a streamline should be ter...
vtkInitialValueProblemSolver * Integrator
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:44
void SetIntegratorTypeToRungeKutta4()
Set/get the integrator type to be used for streamline generation.
std::vector< CustomTerminationCallbackType > CustomTerminationCallback
represent and manipulate attribute data in a dataset
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
double MaximumIntegrationStep
void SetIntegratorTypeToRungeKutta2()
Set/get the integrator type to be used for streamline generation.
void SetIntegratorTypeToRungeKutta45()
Set/get the integrator type to be used for streamline generation.
double InitialIntegrationStep
vtkAbstractInterpolatedVelocityField * InterpolatorPrototype
std::vector< void * > CustomTerminationClientData
vtkDataSetAttributesFieldList InputPD
Streamline generator.
void SetIntegrationDirectionToBackward()
Specify whether the streamline is integrated in the upstream or downstream direction, or in both directions.
int FillInputPortInformation(int port, vtkInformation *info) override
Fill the input port information objects for this algorithm.
Store zero or more vtkInformation instances.
bool HasMatchingPointAttributes
general representation of visualization data
Definition: vtkDataObject.h:54
represent and manipulate 3D points
Definition: vtkPoints.h:28
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.