VTK
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
vtkParticleTracerBase.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtkParticleTracerBase.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 =========================================================================*/
27 #ifndef __vtkParticleTracerBase_h
28 #define __vtkParticleTracerBase_h
29 
30 #include "vtkFiltersFlowPathsModule.h" // For export macro
31 #include "vtkSmartPointer.h" // For protected ivars.
32 #include "vtkPolyDataAlgorithm.h"
33 //BTX
34 #include <vector> // STL Header
35 #include <list> // STL Header
36 //ETX
37 
39 class vtkAbstractParticleWriter;
40 class vtkCellArray;
41 class vtkCharArray;
42 class vtkCompositeDataSet;
43 class vtkDataArray;
44 class vtkDataSet;
45 class vtkDoubleArray;
46 class vtkFloatArray;
47 class vtkGenericCell;
48 class vtkInitialValueProblemSolver;
49 class vtkIntArray;
50 class vtkMultiBlockDataSet;
51 class vtkMultiProcessController;
52 class vtkPointData;
53 class vtkPoints;
54 class vtkPolyData;
56 
57 //BTX
58 namespace vtkParticleTracerBaseNamespace
59 {
60  typedef struct { double x[4]; } Position;
61  typedef struct {
62  // These are used during iteration
64  int CachedDataSetId[2];
65  vtkIdType CachedCellId[2];
67  // These are computed scalars we might display
68  int SourceID;
73  // These are useful to track for debugging etc
74  int ErrorCode;
75  float age;
76  // these are needed across time steps to compute vorticity
77  float rotation;
78  float angularVel;
79  float time;
80  float speed;
81 
82  vtkIdType PointId; //once the partice is added, PointId is valid
84 
85  typedef std::vector<ParticleInformation> ParticleVector;
86  typedef ParticleVector::iterator ParticleIterator;
87  typedef std::list<ParticleInformation> ParticleDataList;
88  typedef ParticleDataList::iterator ParticleListIterator;
89 };
90 //ETX
91 
92 class VTKFILTERSFLOWPATHS_EXPORT vtkParticleTracerBase : public vtkPolyDataAlgorithm
93 {
94 public:
95  enum Solvers
96  {
101  UNKNOWN
102  };
103 
104  vtkTypeMacro(vtkParticleTracerBase,vtkPolyDataAlgorithm)
105  void PrintSelf(ostream& os, vtkIndent indent);
106  void PrintParticleHistories();
107 
109 
111  vtkGetMacro(ComputeVorticity, bool);
112  void SetComputeVorticity(bool);
114 
116 
118  vtkGetMacro(TerminalSpeed, double);
119  void SetTerminalSpeed(double);
121 
123 
125  void SetRotationScale(double);
126  vtkGetMacro(RotationScale, double);
128 
129 
131 
133  vtkSetMacro(IgnorePipelineTime, int);
134  vtkGetMacro(IgnorePipelineTime, int);
135  vtkBooleanMacro(IgnorePipelineTime, int);
137 
139 
146  void SetForceReinjectionEveryNSteps(int);
147  vtkGetMacro(ForceReinjectionEveryNSteps,int);
149 
151 
155  void SetTerminationTime(double t);
156  vtkGetMacro(TerminationTime,double);
158 
159  void SetIntegrator(vtkInitialValueProblemSolver *);
160  vtkGetObjectMacro ( Integrator, vtkInitialValueProblemSolver );
161 
162  void SetIntegratorType(int type);
163  int GetIntegratorType();
164 
166 
170  void SetStartTime(double t);
171  vtkGetMacro(StartTime,double);
173 
174 
176 
182  vtkGetMacro(StaticSeeds,int);
184 
186 
192  vtkGetMacro(StaticMesh,int);
194 
196 
200  virtual void SetParticleWriter(vtkAbstractParticleWriter *pw);
201  vtkGetObjectMacro(ParticleWriter, vtkAbstractParticleWriter);
203 
205 
207  vtkSetStringMacro(ParticleFileName);
208  vtkGetStringMacro(ParticleFileName);
210 
212 
214  vtkSetMacro(EnableParticleWriting,int);
215  vtkGetMacro(EnableParticleWriting,int);
216  vtkBooleanMacro(EnableParticleWriting,int);
218 
220 
222  vtkSetMacro(DisableResetCache,int);
223  vtkGetMacro(DisableResetCache,int);
224  vtkBooleanMacro(DisableResetCache,int);
226 
228 
229  void AddSourceConnection(vtkAlgorithmOutput* input);
230  void RemoveAllSources();
232 
233  protected:
234  vtkSmartPointer<vtkPolyData> Output; //managed by child classes
235  vtkSmartPointer<vtkPointData> ProtoPD;
236  vtkIdType UniqueIdCounter;// global Id counter used to give particles a stamp
237  vtkParticleTracerBaseNamespace::ParticleDataList ParticleHistories;
238  vtkSmartPointer<vtkPointData> ParticlePointData; //the current particle point data consistent
239  //with particle history
240  int ReinjectionCounter;
241 
242  //Everything related to time
243  int IgnorePipelineTime; //whether to use the pipeline time for termination
244  int DisableResetCache; //whether to enable ResetCache() method
245 
246 
248  virtual ~vtkParticleTracerBase();
249 
250  //
251  // Make sure the pipeline knows what type we expect as input
252  //
253  virtual int FillInputPortInformation(int port, vtkInformation* info);
254 
255  //
256  // The usual suspects
257  //
258  virtual int ProcessRequest(vtkInformation* request,
259  vtkInformationVector** inputVector,
260  vtkInformationVector* outputVector);
261 
262  //
263  // Store any information we need in the output and fetch what we can
264  // from the input
265  //
266  virtual int RequestInformation(vtkInformation* request,
267  vtkInformationVector** inputVector,
268  vtkInformationVector* outputVector);
269 
270  //
271  // Compute input time steps given the output step
272  //
273  virtual int RequestUpdateExtent(vtkInformation* request,
274  vtkInformationVector** inputVector,
275  vtkInformationVector* outputVector);
276 
277  //
278  // what the pipeline calls for each time step
279  //
280  virtual int RequestData(vtkInformation* request,
281  vtkInformationVector** inputVector,
282  vtkInformationVector* outputVector);
283 
284  //
285  // these routines are internally called to actually generate the output
286  //
287  virtual int ProcessInput(vtkInformationVector** inputVector);
288 
289  // This is the main part of the algorithm:
290  // * move all the particles one step
291  // * Reinject particles (by adding them to this->ParticleHistories)
292  // either at the beginning or at the end of each step (modulo this->ForceReinjectionEveryNSteps)
293  // * Output a polydata representing the moved particles
294  // Note that if the starting and the ending time coincide, the polydata is still valid.
295  virtual vtkPolyData* Execute(vtkInformationVector** inputVector);
296 
297  // the RequestData will call these methods in turn
298  virtual void Initialize(){} //the first iteration
299  virtual int OutputParticles(vtkPolyData* poly)=0; //every iteration
300  virtual void Finalize(){} //the last iteration
301 
302  //
303  // Initialization of input (vector-field) geometry
304  //
305  int InitializeInterpolator();
306  int UpdateDataCache(vtkDataObject *td);
307 
308 
310 
312  void TestParticles(
315  int &count);
317 
318  void TestParticles(
319  vtkParticleTracerBaseNamespace::ParticleVector &candidates, std::vector<int> &passed);
320 
322 
326  virtual void AssignSeedsToProcessors(double time,
327  vtkDataSet *source, int sourceID, int ptId,
329  int &LocalAssignedCount);
331 
333 
335  virtual void AssignUniqueIds(
338 
340 
342  void UpdateParticleList(
345 
349 
351 
352  void IntegrateParticle(
354  double currenttime, double terminationtime,
355  vtkInitialValueProblemSolver* integrator);
357 
358  // if the particle is added to send list, then returns value is 1,
359  // if it is kept on this process after a retry return value is 0
363  {
364  return true;
365  }
366 
368 
372  bool ComputeDomainExitLocation(
373  double pos[4], double p2[4], double intersection[4],
374  vtkGenericCell *cell);
376 
377  //
378  // Scalar arrays that are generated as each particle is updated
379  //
380  void CreateProtoPD(vtkDataObject* input);
381 
382  vtkFloatArray* GetParticleAge(vtkPointData*);
383  vtkIntArray* GetParticleIds(vtkPointData*);
384  vtkCharArray* GetParticleSourceIds(vtkPointData*);
385  vtkIntArray* GetInjectedPointIds(vtkPointData*);
386  vtkIntArray* GetInjectedStepIds(vtkPointData*);
387  vtkIntArray* GetErrorCodeArr(vtkPointData*);
388  vtkFloatArray* GetParticleVorticity(vtkPointData*);
389  vtkFloatArray* GetParticleRotation(vtkPointData*);
390  vtkFloatArray* GetParticleAngularVel(vtkPointData*);
391 
392 
393  // utility function we use to test if a point is inside any of our local datasets
394  bool InsideBounds(double point[]);
395 
396 
397 
398  void CalculateVorticity( vtkGenericCell* cell, double pcoords[3],
399  vtkDoubleArray* cellVectors, double vorticity[3] );
400 
401  //------------------------------------------------------
402 
403 
404  double GetCacheDataTime(int i);
405  double GetCacheDataTime();
406 
407  virtual void ResetCache();
408  void AddParticle(vtkParticleTracerBaseNamespace::ParticleInformation &info, double* velocity);
409 
411 
414  virtual bool IsPointDataValid(vtkDataObject* input);
415  bool IsPointDataValid(vtkCompositeDataSet* input, std::vector<std::string>& arrayNames);
416  void GetPointDataArrayNames(vtkDataSet* input, std::vector<std::string>& names);
418 
419 private:
421  void SetInterpolatorPrototype(vtkAbstractInterpolatedVelocityField*) {};
422 
424 
429  bool RetryWithPush(
430  vtkParticleTracerBaseNamespace::ParticleInformation &info, double* point1,double delT, int subSteps);
432 
433  //Parameters of tracing
434  vtkInitialValueProblemSolver* Integrator;
435  double IntegrationStep;
436  double MaximumError;
437  bool ComputeVorticity;
438  double RotationScale;
439  double TerminalSpeed;
440 
441  // Important for Caching of Cells/Ids/Weights etc
442  int AllFixedGeometry;
443  int StaticMesh;
444  int StaticSeeds;
445 
446  std::vector<double> InputTimeValues;
447  double StartTime;
448  double TerminationTime;
449  double CurrentTime;
450 
451  int StartTimeStep; //InputTimeValues[StartTimeStep] <= StartTime <= InputTimeValues[StartTimeStep+1]
452  int CurrentTimeStep;
453  int TerminationTimeStep; //computed from start time
454  bool FirstIteration;
455 
456  //Innjection parameters
457  int ForceReinjectionEveryNSteps;
458  vtkTimeStamp ParticleInjectionTime;
459  bool HasCache;
460 
461  // Particle writing to disk
462  vtkAbstractParticleWriter *ParticleWriter;
463  char *ParticleFileName;
464  int EnableParticleWriting;
465 
466 
467  // The main lists which are held during operation- between time step updates
469 
470 
471  // The velocity interpolator
472  vtkSmartPointer<vtkTemporalInterpolatedVelocityField> Interpolator;
473  vtkAbstractInterpolatedVelocityField * InterpolatorPrototype;
474 
475  // Data for time step CurrentTimeStep-1 and CurrentTimeStep
476  vtkSmartPointer<vtkMultiBlockDataSet> CachedData[2];
477 
478  // Cache bounds info for each dataset we will use repeatedly
479  typedef struct {
480  double b[6];
481  } bounds;
482  std::vector<bounds> CachedBounds[2];
483 
484  // temporary variables used by Exeucte(), for convenience only
485 
486  vtkSmartPointer<vtkPoints> OutputCoordinates;
487  vtkSmartPointer<vtkFloatArray> ParticleAge;
488  vtkSmartPointer<vtkIntArray> ParticleIds;
489  vtkSmartPointer<vtkCharArray> ParticleSourceIds;
490  vtkSmartPointer<vtkIntArray> InjectedPointIds;
491  vtkSmartPointer<vtkIntArray> InjectedStepIds;
492  vtkSmartPointer<vtkIntArray> ErrorCode;
493  vtkSmartPointer<vtkFloatArray> ParticleVorticity;
494  vtkSmartPointer<vtkFloatArray> ParticleRotation;
495  vtkSmartPointer<vtkFloatArray> ParticleAngularVel;
496  vtkSmartPointer<vtkDoubleArray> CellVectors;
497  vtkSmartPointer<vtkPointData> OutputPointData;
498  vtkSmartPointer<vtkDataSet> DataReferenceT[2];
499  vtkSmartPointer<vtkCellArray> ParticleCells;
500 
501  vtkParticleTracerBase(const vtkParticleTracerBase&); // Not implemented.
502  void operator=(const vtkParticleTracerBase&); // Not implemented.
503  vtkTimeStamp ExecuteTime;
504 
505  unsigned int NumberOfParticles();
506 
509 
510  static const double Epsilon;
511 
512 };
513 
514 #endif
A helper class for interpolating between times during particle tracing.
virtual void UpdateParticleListFromOtherProcesses()
An abstract class for obtaining the interpolated velocity values at a point.
virtual bool SendParticleToAnotherProcess(vtkParticleTracerBaseNamespace::ParticleInformation &, vtkParticleTracerBaseNamespace::ParticleInformation &, vtkPointData *)
ParticleVector::iterator ParticleIterator
std::list< ParticleInformation > ParticleDataList
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
std::vector< ParticleInformation > ParticleVector
ParticleDataList::iterator ParticleListIterator
A particle tracer for vector fields.