vtkPVFoam.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::vtkPVFoam
26 
27 Description
28  Provides a reader interface for OpenFOAM to VTK interaction.
29 
30 SourceFiles
31  vtkPVFoam.C
32  vtkPVFoam.H
33  vtkPVFoamFields.C
34  vtkPVFoamMesh.C
35  vtkPVFoamMeshLagrangian.C
36  vtkPVFoamTemplates.C
37  vtkPVFoamMeshSet.C
38  vtkPVFoamMeshVolume.C
39  vtkPVFoamMeshZone.C
40  vtkPVFoamSurfaceField.H
41  vtkPVFoamLagrangianFields.H
42  vtkPVFoamPatchField.H
43  vtkPVFoamPointFields.H
44  vtkPVFoamPoints.H
45  vtkPVFoamUpdateInfo.C
46  vtkPVFoamUtils.C
47  vtkPVFoamVolFields.H
48  vtkPVFoamAddToSelection.H
49 
50  // Needed by VTK:
51  vtkDataArrayTemplateImplicit.txx
52 
53 \*---------------------------------------------------------------------------*/
54 
55 #ifndef vtkPVFoam_H
56 #define vtkPVFoam_H
57 
58 // do not include legacy strstream headers
59 #ifndef VTK_EXCLUDE_STRSTREAM_HEADERS
60  #define VTK_EXCLUDE_STRSTREAM_HEADERS
61 #endif
62 
63 #include "className.H"
64 #include "fileName.H"
65 #include "stringList.H"
66 #include "wordList.H"
67 #include "primitivePatch.H"
69 #include "volPointInterpolation.H"
70 #include "surfaceFieldsFwd.H"
71 #include "instantList.H"
72 
73 // * * * * * * * * * * * * * Forward Declarations * * * * * * * * * * * * * //
74 
75 class vtkDataArraySelection;
76 class vtkDataSet;
77 class vtkPoints;
78 class vtkPVFoamReader;
79 class vtkRenderer;
80 class vtkTextActor;
81 class vtkMultiBlockDataSet;
82 class vtkPolyData;
83 class vtkUnstructuredGrid;
84 class vtkIndent;
85 
86 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
87 
88 namespace Foam
89 {
90 
91 // OpenFOAM class forward declarations
92 class argList;
93 class Time;
94 class fvMesh;
95 class IOobjectList;
96 class polyPatch;
97 class faceSet;
98 class pointSet;
99 
100 template<class Type> class IOField;
101 template<class Type> class List;
102 
103 /*---------------------------------------------------------------------------*\
104  Class vtkPVFoam Declaration
105 \*---------------------------------------------------------------------------*/
107 class vtkPVFoam
108 {
109  // Private classes
110 
111  //- Bookkeeping for GUI checklists and the multi-block organization
112  class arrayRange
113  {
114  const char *name_;
115  int block_;
116  int start_;
117  int size_;
118 
119  public:
120 
121  arrayRange(const char *name, const int blockNo=0)
122  :
123  name_(name),
124  block_(blockNo),
125  start_(0),
126  size_(0)
127  {}
128 
129  //- Return the block holding these datasets
130  int block() const
131  {
132  return block_;
133  }
134 
135  //- Assign block number, return previous value
136  int block(int blockNo)
137  {
138  int prev = block_;
139  block_ = blockNo;
140  return prev;
141  }
142 
143  //- Return block name
144  const char* name() const
145  {
146  return name_;
147  }
148 
149  //- Return array start index
150  int start() const
151  {
152  return start_;
153  }
154 
155  //- Return array end index
156  int end() const
157  {
158  return start_ + size_;
159  }
160 
161  //- Return sublist size
162  int size() const
163  {
164  return size_;
165  }
166 
167  bool empty() const
168  {
169  return !size_;
170  }
171 
172  //- Reset the size to zero and optionally assign a new start
173  void reset(const int startAt = 0)
174  {
175  start_ = startAt;
176  size_ = 0;
177  }
178 
179  //- Increment the size
180  void operator+=(const int n)
181  {
182  size_ += n;
183  }
184  };
185 
186  //- Bookkeeping for polyhedral cell decomposition
187  // hide in extra pointMap (cellSet/cellZone) for now
188  class polyDecomp
189  {
190  labelList superCells_;
191  labelList addPointCellLabels_;
192  labelList pointMap_;
193 
194  public:
195 
196  polyDecomp()
197  {}
198 
199  //- Label of original cell for decomposed cells
200  labelList& superCells()
201  {
202  return superCells_;
203  }
204 
205  //- Label of original cell for decomposed cells
206  const labelList& superCells() const
207  {
208  return superCells_;
209  }
210 
211  //- Cell-centre labels for additional points of decomposed cells
212  labelList& addPointCellLabels()
213  {
214  return addPointCellLabels_;
215  }
216 
217  //- Cell-centre labels for additional points of decomposed cells
218  const labelList& addPointCellLabels() const
219  {
220  return addPointCellLabels_;
221  }
222 
223  //- Point labels for subsetted meshes
224  labelList& pointMap()
225  {
226  return pointMap_;
227  }
228 
229  //- Point labels for subsetted meshes
230  const labelList& pointMap() const
231  {
232  return pointMap_;
233  }
234 
235 
236  //- Clear
237  void clear()
238  {
239  superCells_.clear();
240  addPointCellLabels_.clear();
241  pointMap_.clear();
242  }
243  };
244 
245 
246  // Private Data
247 
248  //- Access to the controlling vtkPVFoamReader
249  vtkPVFoamReader* reader_;
250 
251  //- OpenFOAM time control
252  autoPtr<Time> dbPtr_;
253 
254  //- Configuration dictionary
255  dictionary configDict_;
256 
257  //- OpenFOAM mesh
258  fvMesh* meshPtr_;
259 
260  //- The mesh region
261  word meshRegion_;
262 
263  //- The mesh directory for the region
264  fileName meshDir_;
265 
266  //- The time index
267  int timeIndex_;
268 
269  //- Track changes in mesh geometry
270  bool meshChanged_;
271 
272  //- Track changes in fields
273  bool fieldsChanged_;
274 
275  //- Selected geometrical parts (internalMesh, patches, ...)
276  boolList partStatus_;
277 
278  //- Datasets corresponding to selected geometrical pieces
279  // a negative number indicates that no vtkmesh exists for this piece
280  labelList partDataset_;
281 
282  //- First instance and size of various mesh parts
283  // used to index into partStatus_ and partDataset_
284  arrayRange arrayRangeVolume_;
285  arrayRange arrayRangePatches_;
286  arrayRange arrayRangeLagrangian_;
287  arrayRange arrayRangeCellZones_;
288  arrayRange arrayRangeFaceZones_;
289  arrayRange arrayRangePointZones_;
290  arrayRange arrayRangeCellSets_;
291  arrayRange arrayRangeFaceSets_;
292  arrayRange arrayRangePointSets_;
293 
294  //- Decomposed cells information (mesh regions)
295  // TODO: regions
296  List<polyDecomp> regionPolyDecomp_;
297 
298  //- Decomposed cells information (cellZone meshes)
299  List<polyDecomp> zonePolyDecomp_;
300 
301  //- Decomposed cells information (cellSet meshes)
302  List<polyDecomp> csetPolyDecomp_;
303 
304  //- List of patch names for rendering to window
305  List<vtkTextActor*> patchTextActorsPtrs_;
306 
307  // Private Member Functions
308 
309  // Convenience method use to convert the readers from VTK 5
310  // multiblock API to the current composite data infrastructure
311  static void AddToBlock
312  (
313  vtkMultiBlockDataSet* output,
314  vtkDataSet* dataset,
315  const arrayRange&,
316  const label datasetNo,
317  const std::string& datasetName
318  );
319 
320  // Convenience method use to convert the readers from VTK 5
321  // multiblock API to the current composite data infrastructure
322  static vtkDataSet* GetDataSetFromBlock
323  (
324  vtkMultiBlockDataSet* output,
325  const arrayRange&,
326  const label datasetNo
327  );
328 
329  // Convenience method use to convert the readers from VTK 5
330  // multiblock API to the current composite data infrastructure
331  static label GetNumberOfDataSets
332  (
333  vtkMultiBlockDataSet* output,
334  const arrayRange&
335  );
336 
337  //- Reset data counters
338  void resetCounters();
339 
340  // Update information helper functions
341 
342  //- Update the mesh parts selected in the GUI
343  void updateMeshPartsStatus();
344 
345  //- Internal mesh info
346  void updateInfoInternalMesh(vtkDataArraySelection*);
347 
348  //- Lagrangian info
349  void updateInfoLagrangian(vtkDataArraySelection*);
350 
351  //- Patch info
352  void updateInfoPatches
353  (
354  vtkDataArraySelection*,
355  stringList&,
356  const bool
357  );
358 
359  //- Set info
360  void updateInfoSets(vtkDataArraySelection*);
361 
362  //- Zone info
363  void updateInfoZones(vtkDataArraySelection*);
364 
365  //- Get non-empty zone names for zoneType from file
366  wordList getZoneNames(const word& zoneType) const;
367 
368  //- Get non-empty zone names from mesh info
369  template<class ZoneType>
370  wordList getZoneNames
371  (
373  ) const;
374 
375  //- Add objects of Type to paraview array selection
376  template<class Type>
377  label addToSelection
378  (
379  vtkDataArraySelection*,
380  const IOobjectList&,
381  const string& suffix=string::null
382  );
383 
384  //- Add objects of Type to paraview array selection
385  template<template<class> class patchType, class meshType>
386  void addFieldsToSelection
387  (
388  vtkDataArraySelection *select,
389  const IOobjectList& objects,
390  const string& suffix=string::null
391  );
392 
393  //- Add objects of Type to paraview array selection
394  template<template<class> class patchType, class meshType>
395  void addInternalFieldsToSelection
396  (
397  vtkDataArraySelection *select,
398  const IOobjectList& objects,
399  const string& suffix=string::null
400  );
401 
402  //- Field info
403  void updateInfoFields();
404 
405  //- Lagrangian field info
406  void updateInfoLagrangianFields();
407 
408 
409  // Update helper functions
410 
411  //- OpenFOAM mesh
412  void updateFoamMesh();
413 
414  //- Reduce memory footprint after conversion
415  void reduceMemory();
416 
417  //- Volume fields
418  void updateVolFields(vtkMultiBlockDataSet*);
419 
420  //- Point fields
421  void updatePointFields(vtkMultiBlockDataSet*);
422 
423  //- Lagrangian fields
424  void updateLagrangianFields(vtkMultiBlockDataSet*);
425 
426 
427  // Mesh conversion functions
428 
429  //- Volume mesh
430  void convertMeshVolume(vtkMultiBlockDataSet*, int& blockNo);
431 
432  //- Lagrangian mesh
433  void convertMeshLagrangian(vtkMultiBlockDataSet*, int& blockNo);
434 
435  //- Patch meshes
436  void convertMeshPatches(vtkMultiBlockDataSet*, int& blockNo);
437 
438  //- Cell zone meshes
439  void convertMeshCellZones(vtkMultiBlockDataSet*, int& blockNo);
440 
441  //- Face zone meshes
442  void convertMeshFaceZones(vtkMultiBlockDataSet*, int& blockNo);
443 
444  //- Point zone meshes
445  void convertMeshPointZones(vtkMultiBlockDataSet*, int& blockNo);
446 
447  //- Cell set meshes
448  void convertMeshCellSets(vtkMultiBlockDataSet*, int& blockNo);
449 
450  //- Face set meshes
451  void convertMeshFaceSets(vtkMultiBlockDataSet*, int& blockNo);
452 
453  //- Point set meshes
454  void convertMeshPointSets(vtkMultiBlockDataSet*, int& blockNo);
455 
456 
457  // Add mesh functions
458 
459  //- Add internal mesh/cell set meshes
460  vtkUnstructuredGrid* volumeVTKMesh(const fvMesh&, polyDecomp&);
461 
462  //- Add Lagrangian mesh
463  vtkPolyData* lagrangianVTKMesh
464  (
465  const fvMesh&,
466  const word& cloudName
467  );
468 
469  //- Add patch mesh
470  template<class PatchType>
471  vtkPolyData* patchVTKMesh(const word& name, const PatchType&);
472 
473  //- Add point zone
474  vtkPolyData* pointZoneVTKMesh
475  (
476  const fvMesh&,
477  const labelList& pointLabels
478  );
479 
480  //- Add face set mesh
481  vtkPolyData* faceSetVTKMesh
482  (
483  const fvMesh&,
484  const faceSet&
485  );
486 
487  //- Add point mesh
488  vtkPolyData* pointSetVTKMesh
489  (
490  const fvMesh&,
491  const pointSet&
492  );
493 
494  // Field conversion functions
495 
496  //- Convert fields
497  void convertFields(vtkMultiBlockDataSet*);
498 
499  //- Convert Lagrangian fields
500  void convertLagrangianFields(vtkMultiBlockDataSet*);
501 
502 
503  //- Add the fields in the selected time directory to the selection
504  // lists
505  template<class GeoField>
506  label addObjectsToSelection
507  (
508  vtkDataArraySelection*,
509  const IOobjectList&,
510  const string& suffix=string::null
511  );
512 
513 
514  // Convert OpenFOAM fields
515 
516  //- Volume fields - all types
517  template<class Type>
518  void convertVolFields
519  (
520  const fvMesh&,
522  const IOobjectList&,
523  const bool interpFields,
524  vtkMultiBlockDataSet* output
525  );
526 
527  //- Volume internal fields - all types
528  template<class Type>
529  void convertVolInternalFields
530  (
531  const fvMesh&,
532  const IOobjectList&,
533  vtkMultiBlockDataSet* output
534  );
535 
536  //- Volume field - all selected parts
537  template<class Type>
538  void convertVolFieldBlock
539  (
542  vtkMultiBlockDataSet* output,
543  const arrayRange&,
544  const List<polyDecomp>& decompLst
545  );
546 
547  //- Volume internal field - all selected parts
548  template<class Type>
549  void convertVolInternalFieldBlock
550  (
552  ::Internal&,
553  vtkMultiBlockDataSet* output,
554  const arrayRange&,
555  const List<polyDecomp>& decompLst
556  );
557 
558  //- Volume field
559  template<class Type>
560  void convertVolInternalField
561  (
563  ::Internal&,
564  vtkMultiBlockDataSet* output,
565  const arrayRange&,
566  const label datasetNo,
567  const polyDecomp&
568  );
569 
570  //- Patch field
571  template<class Type>
572  void convertPatchField
573  (
574  const word& name,
575  const Field<Type>&,
576  vtkMultiBlockDataSet* output,
577  const arrayRange&,
578  const label datasetNo
579  );
580 
581  //- Face set/zone field from volField
582  template<class Type>
583  void convertSurfaceField
584  (
586  vtkMultiBlockDataSet* output,
587  const arrayRange&,
588  const label datasetNo,
589  const fvMesh&,
590  const labelList& faceLabels
591  );
592 
593  //- Face set/zone field from surfaceField
594  template<class Type>
595  void convertSurfaceField
596  (
598  vtkMultiBlockDataSet* output,
599  const arrayRange& range,
600  const label datasetNo,
601  const fvMesh& mesh,
602  const labelList& faceLabels
603  );
604 
605  //- Surface fields - all types
606  template<class Type>
607  void convertSurfaceFields
608  (
609  const fvMesh& mesh,
610  const IOobjectList& objects,
611  vtkMultiBlockDataSet* output
612  );
613 
614  //- Lagrangian fields - all types
615  template<class Type>
616  void convertLagrangianFields
617  (
618  const IOobjectList&,
619  vtkMultiBlockDataSet* output,
620  const label datasetNo
621  );
622 
623  //- Lagrangian field
624  template<class Type>
625  void convertLagrangianField
626  (
627  const IOField<Type>&,
628  vtkMultiBlockDataSet* output,
629  const arrayRange&,
630  const label datasetNo
631  );
632 
633  //- Point fields - all types
634  template<class Type>
635  void convertPointFields
636  (
637  const fvMesh&,
638  const pointMesh&,
639  const IOobjectList&,
640  vtkMultiBlockDataSet* output
641  );
642 
643  //- Point field - all selected parts
644  template<class Type>
645  void convertPointFieldBlock
646  (
648  vtkMultiBlockDataSet* output,
649  const arrayRange&,
650  const List<polyDecomp>&
651  );
652 
653  //- Point fields
654  template<class Type>
655  void convertPointField
656  (
659  vtkMultiBlockDataSet* output,
660  const arrayRange&,
661  const label datasetNo,
662  const polyDecomp&
663  );
664 
665  //- Patch point field
666  template<class Type>
667  void convertPatchPointField
668  (
669  const word& name,
670  const Field<Type>&,
671  vtkMultiBlockDataSet* output,
672  const arrayRange&,
673  const label datasetNo
674  );
675 
676 
677  // GUI selection helper functions
678 
679  //- Only keep what is listed in hashSet
680  static void pruneObjectList
681  (
682  IOobjectList&,
683  const wordHashSet&
684  );
685 
686  //- Get the list of selected objects
687  IOobjectList getObjects
688  (
689  const wordHashSet& selected,
690  const fvMesh& mesh,
691  const fileName& instance,
692  const fileName& local = ""
693  );
694 
695  //- Retrieve the current selections
696  static wordHashSet getSelected(vtkDataArraySelection*);
697 
698  //- Retrieve a sub-list of the current selections
699  static wordHashSet getSelected
700  (
701  vtkDataArraySelection*,
702  const arrayRange&
703  );
704 
705  //- Retrieve the current selections
706  static stringList getSelectedArrayEntries(vtkDataArraySelection*);
707 
708  //- Retrieve a sub-list of the current selections
709  static stringList getSelectedArrayEntries
710  (
711  vtkDataArraySelection*,
712  const arrayRange&
713  );
714 
715  //- Set selection(s)
716  static void setSelectedArrayEntries
717  (
718  vtkDataArraySelection*,
719  const stringList&
720  );
721 
722  //- Get the first word from the mesh parts selection
723  word getPartName(const int);
724 
725 
726 public:
727 
728  //- Static data members
729 
730  ClassName("vtkPVFoam");
731 
732 
733  // Constructors
734 
735  //- Construct from components
736  vtkPVFoam
737  (
738  const char* const FileName,
739  vtkPVFoamReader* reader
740  );
741 
742  //- Disallow default bitwise copy construction
743  vtkPVFoam(const vtkPVFoam&) = delete;
744 
745 
746  //- Destructor
747  ~vtkPVFoam();
748 
749 
750  // Member Functions
751 
752  //- Update
753  void updateInfo();
754 
755  void Update
756  (
757  vtkMultiBlockDataSet* output,
758  vtkMultiBlockDataSet* lagrangianOutput
759  );
760 
761  //- Clean any storage
762  void CleanUp();
763 
764  //- Allocate and return a list of selected times
765  // returns the count via the parameter
766  double* findTimes(int& nTimeSteps);
767 
768  //- Add/remove patch names to/from the view
769  void renderPatchNames(vtkRenderer*, const bool show);
770 
771  //- Set the runTime to the first plausible request time,
772  // returns the timeIndex
773  // sets to "constant" on error
774  int setTime(int count, const double requestTimes[]);
775 
776 
777  //- The current time index
778  int timeIndex() const
779  {
780  return timeIndex_;
781  }
782 
783 
784  // Access
785 
786  //- Debug information
787  void PrintSelf(ostream&, vtkIndent) const;
788 
789  //- Simple memory used debugging information
790  static void printMemory();
791 
792 
793  // Member Operators
794 
795 
796  //- Disallow default bitwise assignment
797  void operator=(const vtkPVFoam&) = delete;
798 };
799 
800 
801 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
802 
803 } // End namespace Foam
804 
805 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
806 
807 #ifdef NoRepository
808  #include "vtkPVFoamTemplates.C"
809 #endif
810 
811 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
812 
813 #endif
814 
815 // ************************************************************************* //
A HashTable with keys but without contents.
Definition: HashSet.H:59
Interpolation class within a primitive patch. Allows interpolation from points to faces and vice vers...
int timeIndex() const
The current time index.
Definition: vtkPVFoam.H:777
tUEqn clear()
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A class for handling file names.
Definition: fileName.H:79
A list of face labels.
Definition: faceSet.H:48
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
static void printMemory()
Simple memory used debugging information.
A set of point labels.
Definition: pointSet.H:48
labelList pointLabels(nPoints, -1)
Provides a reader interface for OpenFOAM to VTK interaction.
Definition: vtkPVFoam.H:106
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
double * findTimes(int &nTimeSteps)
Allocate and return a list of selected times.
void updateInfo()
Update.
void renderPatchNames(vtkRenderer *, const bool show)
Add/remove patch names to/from the view.
A list of mesh zones.
Generic GeometricField class.
const tensorField & tf
void PrintSelf(ostream &, vtkIndent) const
Debug information.
scalar range
void Update(vtkMultiBlockDataSet *output, vtkMultiBlockDataSet *lagrangianOutput)
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
ClassName("vtkPVFoam")
Static data members.
void CleanUp()
Clean any storage.
dynamicFvMesh & mesh
Pre-declare SubField and related Field type.
Definition: Field.H:56
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
A class for handling words, derived from string.
Definition: word.H:59
Creates a single block of cells from point coordinates, numbers of cells in each direction and an exp...
Definition: block.H:63
static const string null
An empty string.
Definition: string.H:88
objects
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void operator=(const vtkPVFoam &)=delete
Disallow default bitwise assignment.
int setTime(int count, const double requestTimes[])
Set the runTime to the first plausible request time,.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
Macro definitions for declaring ClassName(), NamespaceName(), etc.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
vtkPVFoam(const char *const FileName, vtkPVFoamReader *reader)
Construct from components.
label n
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const word cloudName(propsDict.lookup("cloudName"))
~vtkPVFoam()
Destructor.
A primitive field of type <T> with automated input and output.
Definition: IOField.H:50
Namespace for OpenFOAM.