vtkPVFoam.H 18.5 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
5
    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
6
     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

Class
25
    Foam::vtkPVFoam
26
27

Description
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
    The backend for the vtkPVFoamReader reader module -
    providing a paraview reader interface for OpenFOAM meshes and fields.

    Similar, and sometimes better, functionality may be provided by the
    native VTK OpenFOAM reader. OpenCFD has recently (2017) been working
    on improving the native VTK OpenFOAM reader for the benefit of everyone.

    In some areas the reader module lacks compared to the native reader
    (notably the ability to work on decomosed datasets), but provides
    additional handling of sets,zones,groups. Some features have also since
    been adapted to the native reader. Additionally, the reader module
    provides a useful platform for testing new ideas.

Note
    The reader module allows two levels of caching. The OpenFOAM fvMesh
    can be cached in memory, for faster loading of fields. Additionally,
    the translated VTK geometries are held in a local cache. The cached
    VTK geometries should incur no additional overhead since they use
    the VTK reference counting for their storage management.
47
48

SourceFiles
49
50
51
52
53
    vtkPVFoam.C
    vtkPVFoamFields.C
    vtkPVFoamMesh.C
    vtkPVFoamMeshLagrangian.C
    vtkPVFoamMeshVolume.C
54
    vtkPVFoamTemplates.C
55
    vtkPVFoamUpdateInfo.C
56
57
    vtkPVFoamFieldTemplates.C
    vtkPVFoamUpdateTemplates.C
58
59
60
61
62
63

    // Needed by VTK:
    vtkDataArrayTemplateImplicit.txx

\*---------------------------------------------------------------------------*/

64
65
#ifndef vtkPVFoam_H
#define vtkPVFoam_H
66
67
68
69
70
71
72
73

#include "className.H"
#include "fileName.H"
#include "stringList.H"
#include "wordList.H"
#include "primitivePatch.H"
#include "PrimitivePatchInterpolation.H"
#include "volPointInterpolation.H"
74
#include "foamPvCore.H"
75
#include "foamVtkMeshMaps.H"
76

77
#include "vtkPoints.h"
78
#include "vtkPolyData.h"
79
#include "vtkSmartPointer.h"
80
81
#include "vtkUnstructuredGrid.h"

82
83
// * * * * * * * * * * * * * Forward Declarations  * * * * * * * * * * * * * //

84
class vtkCellArray;
85
86
class vtkDataArraySelection;
class vtkDataSet;
87
class vtkFloatArray;
88
89
class vtkIndent;
class vtkMultiBlockDataSet;
90
class vtkPVFoamReader;
91
92
class vtkRenderer;
class vtkTextActor;
93

94
95
96
97
98
99
100
101
102
103
104
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

// OpenFOAM class forward declarations
class argList;
class Time;
class fvMesh;
class IOobjectList;
class polyPatch;
105
class fvMeshSubset;
106
107

template<class Type> class IOField;
108
template<class Type> class Field;
109
110
111
template<class Type> class List;

/*---------------------------------------------------------------------------*\
112
                        Class vtkPVFoam Declaration
113
114
\*---------------------------------------------------------------------------*/

115
class vtkPVFoam
116
117
:
    private foamPvCore
118
{
119
120
    // Convenience typedefs
    typedef PrimitivePatchInterpolation<primitivePatch> patchInterpolator;
121
122


123
    // Private classes
124

125
126
127
128
129
130
        //- Bookkeeping for internal caching.
        //  Retain an original copy of the geometry as well as a shallow copy
        //  with the output fields.
        //  The original copy is reused for different timestep
        template<class DataType>
        class foamVtkCaching
131
        {
132
        public:
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
            typedef DataType dataType;

            //- The geometry, without any cell/point data
            vtkSmartPointer<dataType> vtkgeom;

            //- The shallow-copy of geometry, plus additional data
            vtkSmartPointer<dataType> dataset;

            //- Number of points associated with geometry
            inline uint64_t nPoints() const
            {
                return vtkgeom ? vtkgeom->GetNumberOfPoints() : 0;
            }

            //- Clear geometry and dataset
            void clearGeom()
            {
                vtkgeom = nullptr;
                dataset = nullptr;
            }

            //- Return a shallow copy of vtkgeom for manipulation
            vtkSmartPointer<dataType> getCopy() const
            {
                auto copy = vtkSmartPointer<dataType>::New();
158
159
160
161
                if (vtkgeom)
                {
                    copy->ShallowCopy(vtkgeom);
                }
162
163
164
165
166
167
168
                return copy;
            }

            //- Make a shallow copy of vtkgeom into dataset
            void reuse()
            {
                dataset = vtkSmartPointer<dataType>::New();
169
170
171
172
                if (vtkgeom)
                {
                    dataset->ShallowCopy(vtkgeom);
                }
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
            }

            //- Set the geometry and make a shallow copy to dataset
            void set(vtkSmartPointer<dataType> geom)
            {
                vtkgeom = geom;
                reuse();
            }

            //- Report basic information to output
            void PrintSelf(std::ostream& os) const
            {
                os << "geom" << nl;
                if (vtkgeom)
                {
                    vtkgeom->PrintSelf(std::cout, vtkIndent(2));
                }
                else
                {
                    os << "nullptr";
                }
                os << nl;

                os << "copy" << nl;
                if (dataset)
                {
                    dataset->PrintSelf(std::cout, vtkIndent(2));
                }
                else
                {
                    os << "nullptr";
                }
                os << nl;
            }
207
208
        };

209

210
211
212
213
214
215
216
217
        //- Bookkeeping for vtkPolyData
        class foamVtpData
        :
            public foamVtkCaching<vtkPolyData>,
            public foamVtkMeshMaps
        {};


218
        //- Bookkeeping for vtkUnstructuredGrid
219
220
221
222
223
        class foamVtuData
        :
            public foamVtkCaching<vtkUnstructuredGrid>,
            public foamVtkMeshMaps
        {};
224

225
226
227

    // Private Data

228
229
        //- Access to the controlling vtkPVFoamReader
        vtkPVFoamReader* reader_;
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245

        //- OpenFOAM time control
        autoPtr<Time> dbPtr_;

        //- OpenFOAM mesh
        fvMesh* meshPtr_;

        //- The mesh region
        word meshRegion_;

        //- The mesh directory for the region
        fileName meshDir_;

        //- The time index
        int timeIndex_;

246
247
248
        //- Previous/current decomposition request
        bool decomposePoly_;

249
        //- Track changes in mesh geometry
250
        enum polyMesh::readUpdateState meshState_;
251

252
253
        //- The index of selected parts mapped to their names
        Map<string> selectedPartIds_;
254

255
256
257
258
259
        //- Any information for 2D (VTP) geometries
        HashTable<foamVtpData, string> cachedVtp_;

        //- Cell maps and other information for 3D (VTU) geometries
        HashTable<foamVtuData, string> cachedVtu_;
260
261

        //- First instance and size of various mesh parts
262
        //  used to index into selectedPartIds and thus indirectly into
263
        //  cachedVtp, cachedVtu
264
265
266
267
268
269
270
271
272
        arrayRange rangeVolume_;
        arrayRange rangePatches_;
        arrayRange rangeLagrangian_;
        arrayRange rangeCellZones_;
        arrayRange rangeFaceZones_;
        arrayRange rangePointZones_;
        arrayRange rangeCellSets_;
        arrayRange rangeFaceSets_;
        arrayRange rangePointSets_;
273
274

        //- List of patch names for rendering to window
275
        List<vtkSmartPointer<vtkTextActor>> patchTextActors_;
276

277

278
279
    // Private Member Functions

280
281
282
283
284
285
286
287
288
289
        template<class Container>
        bool addOutputBlock
        (
            vtkMultiBlockDataSet* output,
            const HashTable<Container, string>& cache,
            const arrayRange& selector,
            const bool singleDataset = false
        ) const;


290
291
292
293
294
295
        //- Reset data counters
        void resetCounters();

      // Update information helper functions

        //- Internal mesh info
296
        void updateInfoInternalMesh(vtkDataArraySelection* select);
297
298

        //- Lagrangian info
299
        void updateInfoLagrangian(vtkDataArraySelection* select);
300

301
302
303
304
305
306
        //- Patch info, modifies enabledEntries
        void updateInfoPatches
        (
            vtkDataArraySelection* select,
            HashSet<string>& enabledEntries
        );
307
308

        //- Set info
309
        void updateInfoSets(vtkDataArraySelection* select);
310
311

        //- Zone info
312
313
        void updateInfoZones(vtkDataArraySelection* select);

314
315
316
317
318
319
320

        //- Get non-empty zone names for zoneType from file
        wordList getZoneNames(const word& zoneType) const;

        //- Get names of on non-empty zones from the mesh info
        template<class ZoneType>
        static wordList getZoneNames
321
        (
322
            const ZoneMesh<ZoneType, polyMesh>& zmesh
323
324
        );

325
326
327
        //- Field info
        template<template<class> class patchType, class meshType>
        void updateInfoFields
328
        (
329
            vtkDataArraySelection* select
330
331
        );

332
333
        //- Lagrangian field info
        void updateInfoLagrangianFields(vtkDataArraySelection* select);
334
335


336
      // Mesh conversion functions
337

338
        //- Convert InternalMesh
339
        void convertMeshVolume();
340

341
        //- Convert Lagrangian points
342
        void convertMeshLagrangian();
343

344
345
346
        //- Convert mesh patches.
        //  The additionalIds (cached data) contain the patch Ids.
        //  There will be several for groups, but only one for regular patches.
347
348
        void convertMeshPatches();

349
        //- Convert cell zones
350
        void convertMeshCellZones();
351

352
        //- Convert cell sets
353
        void convertMeshCellSets();
354

355
356
357
        //- Convert face zones
        void convertMeshFaceZones();

358
        //- Convert face sets
359
        //  The cellMap (cached data) contains the face-labels.
360
        void convertMeshFaceSets();
361

362
363
364
365
        //- Convert point zones
        //  The pointMap (cached data) contains the point-labels.
        void convertMeshPointZones();

366
        //- Convert point sets
367
        //  The pointMap (cached data) contains the point-labels.
368
        void convertMeshPointSets();
369
370


371
      // Add mesh functions
372

373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
        //- Generate vtk points for the current mesh points/decomposition
        static vtkSmartPointer<vtkPoints> movePoints
        (
            const fvMesh& mesh,
            const foamVtuData& vtuData
        );

        //- Generate vtk points for the current mesh points/decomposition,
        //  using the provided pointMap
        static vtkSmartPointer<vtkPoints> movePoints
        (
            const fvMesh& mesh,
            const foamVtuData& vtuData,
            const labelUList& pointMap
        );

        //- Volume mesh as vtkUnstructuredGrid
        static vtkSmartPointer<vtkUnstructuredGrid> volumeVTKMesh
        (
            const fvMesh& mesh,
            foamVtuData& vtuData,
            const bool decompPoly
        );

        //- Subsetted mesh as vtkUnstructuredGrid
        static vtkSmartPointer<vtkUnstructuredGrid> volumeVTKSubsetMesh
        (
            const fvMeshSubset& subsetter,
            foamVtuData& vtuData,
            const bool decompPoly
        );

        //- Volume mesh as vtkUnstructuredGrid
406
        vtkSmartPointer<vtkUnstructuredGrid> volumeVTKMesh
407
408
        (
            const fvMesh& mesh,
409
            foamVtuData& vtuData
410
411
412
413
414
415
416
417
        ) const;

        //- Subsetted mesh as vtkUnstructuredGrid
        vtkSmartPointer<vtkUnstructuredGrid> volumeVTKSubsetMesh
        (
            const fvMeshSubset& subsetter,
            foamVtuData& vtuData
        ) const;
418

419
        //- Lagrangian positions as vtkPolyData
420
        vtkSmartPointer<vtkPolyData> lagrangianVTKMesh
421
422
423
        (
            const polyMesh& mesh,
            const word& cloudName
424
425
426
427
428
429
430
431
432
433
434
435
436
437
        ) const;

        //- Patch points
        template<class PatchType>
        static vtkSmartPointer<vtkPoints> movePatchPoints
        (
            const PatchType& p
        );

        //- Patch faces as vtk-cells
        template<class PatchType>
        static vtkSmartPointer<vtkCellArray> patchFacesVTKCells
        (
            const PatchType& p
438
        );
439

440
441
        //- Patches (mesh or primitive) as vtkPolyData
        template<class PatchType>
442
        static vtkSmartPointer<vtkPolyData> patchVTKMesh
443
444
445
        (
            const PatchType& p
        );
446
447


448
      // Field conversion functions
449

450
451
        //- Copy list to pre-allocated vtk array.
        //  \return number of input items copied
452
        template<class Type>
453
454
455
456
457
458
459
460
461
462
        static label transcribeFloatData
        (
            vtkFloatArray* array,
            const UList<Type>& input,
            const label start = 0
        );

        //- Create named field initialized to zero
        template<class Type>
        static vtkSmartPointer<vtkFloatArray> zeroVTKField
463
464
        (
            const word& name,
465
            const label size
466
        );
467

468
469
470
471
472
473
474
475
        //- Convert float data to VTK field
        template<class Type>
        vtkSmartPointer<vtkFloatArray> convertFieldToVTK
        (
            const word& name,
            const UList<Type>& fld
        ) const;

476
477
        //- Face set/zone field
        template<class Type>
478
        vtkSmartPointer<vtkFloatArray> convertFaceFieldToVTK
479
480
481
        (
            const GeometricField<Type, fvPatchField, volMesh>& fld,
            const labelUList& faceLabels
482
        ) const;
483
484
485

        //- Volume field
        template<class Type>
486
        vtkSmartPointer<vtkFloatArray> convertVolFieldToVTK
487
488
        (
            const GeometricField<Type, fvPatchField, volMesh>& fld,
489
            const foamVtuData& vtuData
490
        ) const;
491
492


493
        //- Convert volume fields
494
        void convertVolFields();
495

496
        //- Convert point fields
497
        void convertPointFields();
498

499
        //- Convert Lagrangian fields
500
        void convertLagrangianFields();
501
502


503
      // Convert OpenFOAM fields
504

505
506
507
508
509
        //- Volume field - all types
        template<class Type>
        void convertVolField
        (
            const PtrList<patchInterpolator>& patchInterpList,
510
            const GeometricField<Type, fvPatchField, volMesh>& fld
511
        );
512

513
514
515
516
517
518
        //- Volume fields - all types
        template<class Type>
        void convertVolFields
        (
            const fvMesh& mesh,
            const PtrList<patchInterpolator>& patchInterpList,
519
            const IOobjectList& objects
520
        );
521

522
523
524
525
526
527
        //- Volume internal fields (DimensionedField)- all types
        template<class Type>
        void convertDimFields
        (
            const fvMesh& mesh,
            const PtrList<patchInterpolator>& patchInterpList,
528
            const IOobjectList& objects
529
        );
530

531
532
533
534
535
536
        //- Volume field - all selected parts
        template<class Type>
        void convertVolFieldBlock
        (
            const GeometricField<Type, fvPatchField, volMesh>& fld,
            autoPtr<GeometricField<Type, pointPatchField, pointMesh>>& ptfPtr,
537
            const arrayRange& range
538
        );
539

540
541
542
543
544
        //- Lagrangian fields - all types
        template<class Type>
        void convertLagrangianFields
        (
            const IOobjectList& objects,
545
            vtkPolyData* vtkmesh
546
        );
547

548
549
550
551
552
        //- Point fields - all types
        template<class Type>
        void convertPointFields
        (
            const pointMesh& pMesh,
553
            const IOobjectList& objectst
554
        );
555

556
557
558
559
560
        //- Point field - all selected parts
        template<class Type>
        void convertPointFieldBlock
        (
            const GeometricField<Type, pointPatchField, pointMesh>& pfld,
561
            const arrayRange& range
562
        );
563

564
565
        //- Point field
        template<class Type>
566
        vtkSmartPointer<vtkFloatArray> convertPointField
567
568
569
        (
            const GeometricField<Type, pointPatchField, pointMesh>& pfld,
            const GeometricField<Type, fvPatchField, volMesh>& vfld,
570
            const foamVtuData& vtuData
571
        );
572
573


574
      // GUI selection helper functions
575

576
        //- Get the first word from the reader 'parts' selection
577
        word getReaderPartName(const int partId) const;
578
579


580
    // Constructors
581
582

        //- Disallow default bitwise copy construct
583
        vtkPVFoam(const vtkPVFoam&) = delete;
584
585

        //- Disallow default bitwise assignment
586
        void operator=(const vtkPVFoam&) = delete;
587
588
589
590
591
592


public:

    //- Static data members

593
        ClassName("vtkPVFoam");
594
595
596
597
598


    // Constructors

        //- Construct from components
599
        vtkPVFoam
600
601
        (
            const char* const FileName,
602
            vtkPVFoamReader* reader
603
604
605
606
        );


    //- Destructor
607
    ~vtkPVFoam();
608
609
610
611
612
613
614
615
616
617


    // Member Functions

        //- Update
        void updateInfo();

        void Update
        (
            vtkMultiBlockDataSet* output,
618
            vtkMultiBlockDataSet* outputLagrangian
619
620
        );

621
622
623
624
625
        //- Final part of Update(), after any last minute rendering.
        void UpdateFinalize();

        //- Add/remove patch names to/from the view
        void renderPatchNames(vtkRenderer* renderer, const bool show);
626

627
628
629
        //- Return a list of selected times.
        //  Use STL container since these values are used by the plugin
        std::vector<double> findTimes(const bool skipZero = false) const;
630

631
        //- Set the runTime to the first plausible request time,
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
        //  returns the timeIndex
        //  sets to "constant" on error
        int setTime(int count, const double requestTimes[]);


        //- The current time index
        int timeIndex() const
        {
           return timeIndex_;
        }


     // Access

        //- Debug information
        void PrintSelf(ostream&, vtkIndent) const;

649
        void printInfo() const;
650
651
652
653
654
655
656
657
658
659
};


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#ifdef NoRepository
660
    #include "vtkPVFoamTemplates.C"
661
662
663
664
665
666
667
#endif

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

// ************************************************************************* //