VTK  9.3.0
vtkTemporalInterpolatedVelocityField.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
35#ifndef vtkTemporalInterpolatedVelocityField_h
36#define vtkTemporalInterpolatedVelocityField_h
37
38#include "vtkDeprecation.h" // For VTK_DEPRECATED_IN_9_2_0
39#include "vtkFiltersFlowPathsModule.h" // For export macro
40#include "vtkFunctionSet.h"
41#include "vtkSmartPointer.h" // For vtkSmartPointer
42
43#include <vector> // For internal structures
44
45VTK_ABI_NAMESPACE_BEGIN
49class vtkDataArray;
50class vtkDataSet;
51class vtkDoubleArray;
53class vtkGenericCell;
54class vtkLocator;
55class vtkPointData;
56
57class VTKFILTERSFLOWPATHS_EXPORT vtkTemporalInterpolatedVelocityField : public vtkFunctionSet
58{
59public:
61 void PrintSelf(ostream& os, vtkIndent indent) override;
62
68
73 {
74 INSIDE_ALL = 0,
75 OUTSIDE_ALL = 1,
76 OUTSIDE_T0 = 2,
77 OUTSIDE_T1 = 3
78 };
79
84 {
85 DIFFERENT = 0,
86 STATIC = 1,
87 LINEAR_TRANSFORMATION = 2,
88 SAME_TOPOLOGY = 3
89 };
90
92 /*
93 * Set/Get the type of variance of the mesh over time.
94 *
95 * DIFFERENT = 0,
96 * STATIC = 1,
97 * LINEAR_TRANSFORMATION = 2
98 * SAME_TOPOLOGY = 3
99 */
100 vtkSetClampMacro(MeshOverTime, int, DIFFERENT, SAME_TOPOLOGY);
101 void SetMeshOverTimeToDifferent() { this->SetMeshOverTime(DIFFERENT); }
102 void SetMeshOverTimeToStatic() { this->SetMeshOverTime(STATIC); }
103 void SetMeshOverTimeToLinearTransformation() { this->SetMeshOverTime(LINEAR_TRANSFORMATION); }
104 void SetMeshOverTimeToSameTopology() { this->SetMeshOverTime(SAME_TOPOLOGY); }
105 vtkGetMacro(MeshOverTime, int);
107
116
124
125 using Superclass::FunctionValues;
127
131 int FunctionValues(double* x, double* u) override;
132 int FunctionValuesAtT(int T, double* x, double* u);
134
140 void SelectVectors(const char* fieldName) { this->SetVectorsSelection(fieldName); }
141
143
148 void AddDataSetAtTime(int N, double T, vtkDataSet* dataset);
149 VTK_DEPRECATED_IN_9_2_0("Use AddDataSetAtTime and SetMeshOverTime instead")
150 void SetDataSetAtTime(int, int, double, vtkDataSet*, bool) {}
152
154
159 bool GetCachedCellIds(vtkIdType id[2], int ds[2]);
160 void SetCachedCellIds(vtkIdType id[2], int ds[2]);
162
168
170
174 int TestPoint(double* x);
175 int QuickTestPoint(double* x);
177
179
183 vtkGetVector3Macro(LastGoodVelocity, double);
185
187
190 vtkGetMacro(CurrentWeight, double);
192
193 bool InterpolatePoint(vtkPointData* outPD1, vtkPointData* outPD2, vtkIdType outIndex);
194
195 bool InterpolatePoint(int T, vtkPointData* outPD1, vtkIdType outIndex);
196
198 int T, double pcoords[3], double* weights, vtkGenericCell*& cell, vtkDoubleArray* cellVectors);
199
201 VTK_DEPRECATED_IN_9_2_0("Use GetMeshOverTime() instead.")
202 bool IsStatic(int) { return this->MeshOverTime == STATIC; }
203
205
207
214 vtkGetObjectMacro(FindCellStrategy, vtkFindCellStrategy);
216
217protected:
220
221 virtual void SetVectorsSelection(const char* v);
222
223 int MeshOverTime = MeshOverTimeTypes::DIFFERENT;
224
226 const std::vector<vtkDataSet*>& datasets, vtkFindCellStrategy* strategy,
227 const std::vector<vtkSmartPointer<vtkLocator>>& locators,
228 const std::vector<vtkSmartPointer<vtkAbstractCellLinks>>& links);
229
230 void CreateLocators(const std::vector<vtkDataSet*>& datasets, vtkFindCellStrategy* strategy,
231 std::vector<vtkSmartPointer<vtkLocator>>& locators);
232 void CreateLinks(const std::vector<vtkDataSet*>& datasets,
233 std::vector<vtkSmartPointer<vtkAbstractCellLinks>>& links);
235 std::vector<vtkSmartPointer<vtkLocator>>& linearCellLocators);
236
237 double Vals1[3];
238 double Vals2[3];
239 double Times[2];
240 double LastGoodVelocity[3];
241
242 static const double WEIGHT_TO_TOLERANCE;
243
244 // The weight (0.0->1.0) of the value of T between the two available
245 // time values for the current computation
247 // One minus the CurrentWeight
249 // A scaling factor used when calculating the CurrentWeight { 1.0/(T2-T1) }
251
253 std::vector<vtkSmartPointer<vtkLocator>> Locators[2];
254 std::vector<vtkSmartPointer<vtkLocator>> InitialCellLocators;
255 std::vector<vtkSmartPointer<vtkAbstractCellLinks>> Links[2];
256 std::vector<size_t> MaxCellSizes[2];
257
259
260private:
261 // Hide this since we need multiple time steps and are using a different
262 // function prototype
263 virtual void AddDataSet(vtkDataSet*) {}
264
266 void operator=(const vtkTemporalInterpolatedVelocityField&) = delete;
267};
268
269VTK_ABI_NAMESPACE_END
270#endif
abstract superclass for composite (multi-block or AMR) datasets
An abstract class for obtaining the interpolated velocity values at a point.
abstract superclass for arrays of numeric data
abstract class to specify dataset behavior
Definition vtkDataSet.h:53
dynamic, self-adjusting array of double
helper class to manage the vtkPointSet::FindCell() METHOD
Abstract interface for sets of functions.
provides thread-safe access to cells
a simple class to control print indentation
Definition vtkIndent.h:29
abstract base class for objects that accelerate spatial searches
Definition vtkLocator.h:59
represent and manipulate point attribute data
Hold a reference to a vtkObjectBase instance.
A helper class for interpolating between times during particle tracing.
int QuickTestPoint(double *x)
A utility function which evaluates the point at T1, T2 to see if it is inside the data at both times ...
int TestPoint(double *x)
A utility function which evaluates the point at T1, T2 to see if it is inside the data at both times ...
IDStates
States that define where the cell id are.
void AddDataSetAtTime(int N, double T, vtkDataSet *dataset)
In order to use this class, two sets of data must be supplied, corresponding to times T1 and T2.
MeshOverTimeTypes
Types of Variance of Mesh over time.
virtual void CopyParameters(vtkTemporalInterpolatedVelocityField *from)
Copy essential parameters between instances of this class.
bool InterpolatePoint(vtkPointData *outPD1, vtkPointData *outPD2, vtkIdType outIndex)
void SelectVectors(const char *fieldName)
If you want to work with an arbitrary vector array, then set its name here.
virtual void SetVectorsSelection(const char *v)
void CreateLinearTransformCellLocators(const std::vector< vtkSmartPointer< vtkLocator > > &locators, std::vector< vtkSmartPointer< vtkLocator > > &linearCellLocators)
bool GetVorticityData(int T, double pcoords[3], double *weights, vtkGenericCell *&cell, vtkDoubleArray *cellVectors)
std::vector< vtkSmartPointer< vtkLocator > > InitialCellLocators
bool InterpolatePoint(int T, vtkPointData *outPD1, vtkIdType outIndex)
virtual void SetFindCellStrategy(vtkFindCellStrategy *)
Set / get the strategy used to perform the FindCell() operation.
void SetCachedCellIds(vtkIdType id[2], int ds[2])
Between iterations of the Particle Tracer, Id's of the Cell are stored and then at the start of the n...
bool GetCachedCellIds(vtkIdType id[2], int ds[2])
Between iterations of the Particle Tracer, Id's of the Cell are stored and then at the start of the n...
void CreateLinks(const std::vector< vtkDataSet * > &datasets, std::vector< vtkSmartPointer< vtkAbstractCellLinks > > &links)
void ClearCache()
Set the last cell id to -1 so that the next search does not start from the previous cell.
int FunctionValues(double *x, double *u) override
Evaluate the velocity field, f, at (x, y, z, t).
void CreateLocators(const std::vector< vtkDataSet * > &datasets, vtkFindCellStrategy *strategy, std::vector< vtkSmartPointer< vtkLocator > > &locators)
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
void InitializeWithLocators(vtkCompositeInterpolatedVelocityField *ivf, const std::vector< vtkDataSet * > &datasets, vtkFindCellStrategy *strategy, const std::vector< vtkSmartPointer< vtkLocator > > &locators, const std::vector< vtkSmartPointer< vtkAbstractCellLinks > > &links)
int FunctionValuesAtT(int T, double *x, double *u)
Evaluate the velocity field, f, at (x, y, z, t).
void Initialize(vtkCompositeDataSet *t0, vtkCompositeDataSet *t1)
The Initialize() method is used to build and cache supporting structures (such as locators) which are...
static vtkTemporalInterpolatedVelocityField * New()
Construct a vtkTemporalInterpolatedVelocityField with no initial data set.
#define VTK_DEPRECATED_IN_9_2_0(reason)
int vtkIdType
Definition vtkType.h:315