vtkPVFoamMesh.C 20.4 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-2018 OpenCFD Ltd.
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
-------------------------------------------------------------------------------
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/>.

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

26
#include "vtkPVFoam.H"
27
28
29
30
31

// OpenFOAM includes
#include "cellSet.H"
#include "faceSet.H"
#include "pointSet.H"
32
#include "faMesh.H"
33
#include "fvMeshSubset.H"
34
#include "vtkPVFoamReader.h"
35
36
37
38
39
40
41
#include "uindirectPrimitivePatch.H"

// VTK includes
#include "vtkDataArraySelection.h"
#include "vtkMultiBlockDataSet.h"
#include "vtkPolyData.h"
#include "vtkUnstructuredGrid.h"
42
#include "vtkSmartPointer.h"
43
44
45

// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

46
void Foam::vtkPVFoam::convertMeshVolume()
47
{
48
49
50
51
52
    // Convert internalMesh - actually only one part
    const List<label> partIds =
        rangeVolume_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;
53
54
55

    if (debug)
    {
56
        Info<< "<beg> " << FUNCTION_NAME << nl;
57
58
59
        printMemory();
    }

60
    for (const auto partId : partIds)
61
    {
62
63
        const auto& longName = selectedPartIds_[partId];
        foamVtuData& vtuData = cachedVtu_(longName);
64

65
66
67
68
        vtkSmartPointer<vtkUnstructuredGrid> vtkgeom;
        if (vtuData.nPoints())
        {
            if (meshState_ == polyMesh::UNCHANGED)
69
            {
70
                if (debug)
71
                {
72
                    Info<< "reuse " << longName << nl;
73
                }
74
75
                vtuData.reuse(); // reuse
                continue;
76
            }
77
            else if (meshState_ == polyMesh::POINTS_MOVED)
78
            {
79
80
81
82
83
                if (debug)
                {
                    Info<< "move points " << longName << nl;
                }
                vtkgeom = vtuData.getCopy();
84
                vtkgeom->SetPoints(vtuData.points(mesh));
85
            }
86
        }
87

88
89
90
        if (!vtkgeom)
        {
            // Nothing usable from cache - create new geometry
91
            vtkgeom = vtuData.internal(mesh, this->decomposePoly_);
92
        }
93
94

        vtuData.set(vtkgeom);
95
96
97
98
    }

    if (debug)
    {
99
        Info<< "<end> " << FUNCTION_NAME << nl;
100
101
102
103
104
        printMemory();
    }
}


105
void Foam::vtkPVFoam::convertMeshLagrangian()
106
{
107
108
109
110
    const List<label> partIds =
        rangeClouds_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;
111
112
113

    if (debug)
    {
114
        Info<< "<beg> " << FUNCTION_NAME << nl;
115
116
117
        printMemory();
    }

118
    for (const auto partId : partIds)
119
    {
120
121
122
123
124
125
126
127
        const auto& longName = selectedPartIds_[partId];
        const word cloudName = getFoamName(longName);

        // Direct conversion, no caching for Lagrangian
        cachedVtp_(longName).set
        (
            lagrangianVTKMesh(mesh, cloudName)
        );
128
129
130
131
    }

    if (debug)
    {
132
        Info<< "<end> " << FUNCTION_NAME << nl;
133
134
135
136
137
        printMemory();
    }
}


138
void Foam::vtkPVFoam::convertMeshPatches()
139
{
140
141
142
143
    const List<label> partIds =
        rangePatches_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;
144
145
146
147
    const polyBoundaryMesh& patches = mesh.boundaryMesh();

    if (debug)
    {
148
        Info<< "<beg> " << FUNCTION_NAME << nl;
149
150
151
        printMemory();
    }

152
    for (const auto partId : partIds)
153
    {
154
        const auto& longName = selectedPartIds_[partId];
155
156
        const word partName = getFoamName(longName);
        foamVtpData& vtpData = cachedVtp_(longName);
157

158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
        vtkSmartPointer<vtkPolyData> vtkgeom;
        if (vtpData.nPoints())
        {
            if (meshState_ == polyMesh::UNCHANGED)
            {
                // Without movement is easy.
                if (debug)
                {
                    Info<< "reuse " << longName << nl;
                }
                vtpData.reuse();
                continue;
            }
            else if (meshState_ == polyMesh::POINTS_MOVED)
            {
                // Point movement on single patch is OK

                const labelList& patchIds = vtpData.additionalIds();
                if (patchIds.size() == 1)
                {
                    vtkgeom = vtpData.getCopy();
179
180
181
182
                    vtkgeom->SetPoints
                    (
                        vtk::Tools::Patch::points(patches[patchIds[0]])
                    );
183
184
185
186
                    continue;
                }
            }
        }
187
188

        if (longName.startsWith("group/"))
189
190
        {
            // Patch group. Collect patch faces.
191

192
193
            vtpData.clear(); // Remove any old mappings

194
195
196
197
198
199
            const labelList& patchIds =
                patches.groupPatchIDs().lookup(partName, labelList());

            if (debug)
            {
                Info<< "Creating VTK mesh for patches [" << patchIds <<"] "
200
                    << longName << nl;
201
202
            }

203
204
205
            // Store good patch ids as additionalIds
            vtpData.additionalIds().reserve(patchIds.size());

206
            label sz = 0;
207
            for (const auto id : patchIds)
208
            {
209
210
211
212
213
214
215
                const label n = patches[id].size();

                if (n)
                {
                    vtpData.additionalIds().append(id);
                    sz += n;
                }
216
            }
217
218
219
220
221
222
            Foam::sort(vtpData.additionalIds());

            // Temporarily (mis)use cellMap for face labels
            DynamicList<label>& faceLabels = vtpData.cellMap();
            faceLabels.reserve(sz);

223
            for (const auto id : vtpData.additionalIds())
224
            {
225
                const auto& pp = patches[id];
226
227
                forAll(pp, i)
                {
228
                    faceLabels.append(pp.start()+i);
229
230
                }
            }
231

232
233
234
            if (faceLabels.size())
            {
                uindirectPrimitivePatch pp
235
                (
236
237
238
239
                    UIndirectList<face>(mesh.faces(), faceLabels),
                    mesh.points()
                );

240
                vtkgeom = vtk::Tools::Patch::mesh(pp);
241
            }
242
243

            faceLabels.clear();  // Unneeded
244
245
246
        }
        else
        {
247
248
            vtpData.clear(); // Remove any old mappings

249
250
251
252
253
            const label patchId = patches.findPatchID(partName);

            if (debug)
            {
                Info<< "Creating VTK mesh for patch [" << patchId <<"] "
254
                    << partName << nl;
255
            }
256

257
258
            if (patchId >= 0)
            {
259
260
261
                // Store good patch id as additionalIds
                vtpData.additionalIds() = {patchId};

262
                vtkgeom = vtk::Tools::Patch::mesh(patches[patchId]);
263
            }
264
265
        }

266
        if (vtkgeom)
267
        {
268
269
270
271
272
            vtpData.set(vtkgeom);
        }
        else
        {
            cachedVtp_.erase(longName);
273
274
275
276
277
        }
    }

    if (debug)
    {
278
        Info<< "<end> " << FUNCTION_NAME << nl;
279
280
281
282
283
        printMemory();
    }
}


284
void Foam::vtkPVFoam::convertMeshArea()
285
{
286
287
288
    // Convert areaMesh - actually only one part
    const List<label> partIds =
        rangeArea_.intersection(selectedPartIds_);
289

290
    if (partIds.empty())
291
292
293
    {
        return;
    }
294
295
296
297
298
299
300
301
302
303
304
    else if (!areaMeshPtr_ && volMeshPtr_)
    {
        areaMeshPtr_ = new faMesh(*volMeshPtr_);
    }

    if (!areaMeshPtr_)
    {
        return;
    }

    const faMesh& mesh = *areaMeshPtr_;
305
306
307

    if (debug)
    {
308
        Info<< "<beg> " << FUNCTION_NAME << nl;
309
310
311
        printMemory();
    }

312
    for (const auto partId : partIds)
313
    {
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
        const auto& longName = selectedPartIds_[partId];
        foamVtpData& vtpData = cachedVtp_(longName);

        vtkSmartPointer<vtkPolyData> vtkgeom;

// This needs checking
// #if 0
//         if (vtpData.nPoints())
//         {
//             if (meshState_ == polyMesh::UNCHANGED)
//             {
//                 if (debug)
//                 {
//                     Info<< "reuse " << longName << nl;
//                 }
//                 vtpData.reuse(); // reuse
//                 continue;
//             }
//             else if (meshState_ == polyMesh::POINTS_MOVED)
//             {
//                 if (debug)
//                 {
//                     Info<< "move points " << longName << nl;
//                 }
//                 vtkgeom = vtpData.getCopy();
//                 vtkgeom->SetPoints(movePoints(mesh, vtpData));
//             }
//         }
// #endif

        if (!vtkgeom)
345
        {
346
347
348
            vtpData.clear(); // Remove any old mappings

            // Nothing usable from cache - create new geometry
349
            vtkgeom = vtk::Tools::Patch::mesh(mesh.patch());
350
351
        }

352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
        vtpData.set(vtkgeom);
    }

    if (debug)
    {
        Info<< "<end> " << FUNCTION_NAME << nl;
        printMemory();
    }
}


void Foam::vtkPVFoam::convertMeshCellZones()
{
    const List<label> partIds =
        rangeCellZones_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;

    if (partIds.empty())
    {
        return;
    }

    if (debug)
    {
        Info<< "<beg> " << FUNCTION_NAME << nl;
        printMemory();
    }

    const cellZoneMesh& zMesh = mesh.cellZones();
    for (const auto partId : partIds)
    {
384
        const auto& longName = selectedPartIds_[partId];
385
        const word zoneName = getFoamName(longName);
386
387
        const label  zoneId = zMesh.findZoneID(zoneName);

388
        if (zoneId < 0)
389
390
391
392
393
394
395
        {
            continue;
        }

        if (debug)
        {
            Info<< "Creating VTK mesh for cellZone[" << zoneId << "] "
396
                << zoneName << nl;
397
398
        }

399
        foamVtuData& vtuData = cachedVtu_(longName);
400

401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
        vtkSmartPointer<vtkUnstructuredGrid> vtkgeom;
        if (vtuData.nPoints() && vtuData.pointMap().size())
        {
            if (meshState_ == polyMesh::UNCHANGED)
            {
                if (debug)
                {
                    Info<< "reuse " << longName << nl;
                }
                vtuData.reuse();
                continue;
            }
            else if (meshState_ == polyMesh::POINTS_MOVED)
            {
                if (debug)
                {
                    Info<< "move points " << longName << nl;
                }
                vtkgeom = vtuData.getCopy();
420
                vtkgeom->SetPoints(vtuData.points(mesh, vtuData.pointMap()));
421
422
423
424
425
            }
        }

        if (!vtkgeom)
        {
426
            fvMeshSubset subsetter(mesh, zMesh[zoneId]);
427

428
            vtkgeom = vtuData.subset(subsetter, this->decomposePoly_);
429
430
431
        }

        vtuData.set(vtkgeom);
432
433
434
435
    }

    if (debug)
    {
436
        Info<< "<end> " << FUNCTION_NAME << nl;
437
438
439
440
441
        printMemory();
    }
}


442
void Foam::vtkPVFoam::convertMeshCellSets()
443
{
444
445
446
447
    const List<label> partIds =
        rangeCellSets_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;
448
449
450

    if (debug)
    {
451
        Info<< "<beg> " << FUNCTION_NAME << nl;
452
453
454
        printMemory();
    }

455
    for (const auto partId : partIds)
456
    {
457
        const auto& longName = selectedPartIds_[partId];
458
        const word partName = getFoamName(longName);
459

460
461
        if (debug)
        {
462
            Info<< "Creating VTK mesh for cellSet=" << partName << nl;
463
464
        }

465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
        foamVtuData& vtuData = cachedVtu_(longName);

        vtkSmartPointer<vtkUnstructuredGrid> vtkgeom;
        if (vtuData.nPoints() && vtuData.pointMap().size())
        {
            if (meshState_ == polyMesh::UNCHANGED)
            {
                if (debug)
                {
                    Info<< "reuse " << longName << nl;
                }
                vtuData.reuse();
                continue;
            }
            else if (meshState_ == polyMesh::POINTS_MOVED)
            {
                if (debug)
                {
                    Info<< "move points " << longName << nl;
                }
                vtkgeom = vtuData.getCopy();
486
                vtkgeom->SetPoints(vtuData.points(mesh, vtuData.pointMap()));
487
488
            }
        }
489

490
491
        if (!vtkgeom)
        {
492
            fvMeshSubset subsetter(mesh, cellSet(mesh, partName));
493

494
            vtkgeom = vtuData.subset(subsetter, this->decomposePoly_);
495
496
497
        }

        vtuData.set(vtkgeom);
498
499
500
501
    }

    if (debug)
    {
502
        Info<< "<end> " << FUNCTION_NAME << nl;
503
504
505
506
507
        printMemory();
    }
}


508
void Foam::vtkPVFoam::convertMeshFaceZones()
509
{
510
511
512
513
    const List<label> partIds =
        rangeFaceZones_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;
514

515
    if (partIds.empty())
516
517
518
519
520
521
    {
        return;
    }

    if (debug)
    {
522
        Info<< "<beg> " << FUNCTION_NAME << nl;
523
524
525
526
        printMemory();
    }

    const faceZoneMesh& zMesh = mesh.faceZones();
527
    for (const auto partId : partIds)
528
    {
529
        const auto& longName = selectedPartIds_[partId];
530
        const word zoneName = getFoamName(longName);
531
532
        const label  zoneId = zMesh.findZoneID(zoneName);

533
        if (zoneId < 0)
534
535
536
537
538
539
        {
            continue;
        }
        if (debug)
        {
            Info<< "Creating VTKmesh for faceZone[" << zoneId << "] "
540
                << zoneName << nl;
541
542
        }

543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
        foamVtpData& vtpData = cachedVtp_(longName);

        vtkSmartPointer<vtkPolyData> vtkgeom;
        if (vtpData.nPoints())
        {
            if (meshState_ == polyMesh::UNCHANGED)
            {
                // Without movement is easy.
                if (debug)
                {
                    Info<<"reuse " << longName << nl;
                }
                vtpData.reuse();
                continue;
            }
            else if (meshState_ == polyMesh::POINTS_MOVED)
            {
                // Need point maps etc - not worth it at the moment
            }
        }
563

564
565
566
567
568
        if (!vtkgeom)
        {
            vtpData.clear();  // No additional ids, maps

            const primitiveFacePatch& pp = zMesh[zoneId]();
569
            vtkgeom = vtk::Tools::Patch::mesh(pp);
570
571
572
        }

        vtpData.set(vtkgeom);
573
574
575
576
    }

    if (debug)
    {
577
        Info<< "<end> " << FUNCTION_NAME << nl;
578
579
580
581
582
        printMemory();
    }
}


583
void Foam::vtkPVFoam::convertMeshFaceSets()
584
{
585
586
587
588
    const List<label> partIds =
        rangeFaceSets_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;
589
590
591

    if (debug)
    {
592
        Info<< "<beg> " << FUNCTION_NAME << nl;
593
594
595
        printMemory();
    }

596
    for (const auto partId : partIds)
597
    {
598
        const auto& longName = selectedPartIds_[partId];
599
        const word partName = getFoamName(longName);
600

601
602
        if (debug)
        {
603
            Info<< "Creating VTK mesh for faceSet=" << partName << nl;
604
605
        }

606
        foamVtpData& vtpData = cachedVtp_(longName);
607

608
609
        vtkSmartPointer<vtkPolyData> vtkgeom;
        if (vtpData.nPoints())
610
        {
611
612
613
614
615
616
617
618
619
620
621
622
623
624
            if (meshState_ == polyMesh::UNCHANGED)
            {
                // Without movement is easy.
                if (debug)
                {
                    Info<<"reuse " << longName << nl;
                }
                vtpData.reuse();
                continue;
            }
            else if (meshState_ == polyMesh::POINTS_MOVED)
            {
                // Need point maps etc - not worth it at the moment
            }
625
        }
626

627
628
629
        if (!vtkgeom)
        {
            vtpData.clear();  // No other additional ids, maps
630

631
632
633
634
635
636
637
638
639
640
641
            // Misuse cellMap for face labels - sorted order for reliability
            vtpData.cellMap() = faceSet(mesh, partName).sortedToc();

            if (vtpData.cellMap().size())
            {
                uindirectPrimitivePatch pp
                (
                    UIndirectList<face>(mesh.faces(), vtpData.cellMap()),
                    mesh.points()
                );

642
                vtkgeom = vtk::Tools::Patch::mesh(pp);
643
644
645
646
647
648
649
650
651
652
653
            }
        }

        if (vtkgeom)
        {
            vtpData.set(vtkgeom);
        }
        else
        {
            cachedVtp_.erase(longName);
        }
654
655
656
657
    }

    if (debug)
    {
658
        Info<< "<end> " << FUNCTION_NAME << nl;
659
660
661
662
663
        printMemory();
    }
}


664
void Foam::vtkPVFoam::convertMeshPointZones()
665
{
666
667
668
669
670
671
672
673
674
    const List<label> partIds =
        rangePointZones_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;

    if (partIds.empty())
    {
        return;
    }
675
676
677

    if (debug)
    {
678
        Info<< "<beg> " << FUNCTION_NAME << nl;
679
680
681
        printMemory();
    }

682
683
    const pointZoneMesh& zMesh = mesh.pointZones();
    for (const auto partId : partIds)
684
    {
685
686
687
688
689
        const auto& longName = selectedPartIds_[partId];
        const word zoneName = getFoamName(longName);
        const label zoneId = zMesh.findZoneID(zoneName);

        if (zoneId < 0)
690
        {
691
692
            continue;
        }
693

694
        foamVtpData& vtpData = cachedVtp_(longName);
695

696
697
698
699
        vtkSmartPointer<vtkPolyData> vtkgeom;
        if (vtpData.nPoints() && vtpData.pointMap().size())
        {
            if (meshState_ == polyMesh::UNCHANGED)
700
            {
701
702
703
704
705
                if (debug)
                {
                    Info<< "reusing " << longName << nl;
                }
                vtpData.reuse();
706
707
                continue;
            }
708
            else if (meshState_ == polyMesh::POINTS_MOVED)
709
            {
710
                if (debug)
711
                {
712
                    Info<< "move points " << longName << nl;
713
                }
714
                vtkgeom = vtpData.getCopy();
715
            }
716
        }
717

718
719
720
721
722
723
        if (!vtkgeom)
        {
            // First time, or topo change
            vtkgeom = vtkSmartPointer<vtkPolyData>::New();
            vtpData.pointMap() = zMesh[zoneId];
        }
724

725
726
        const pointField& points = mesh.points();
        const labelUList& pointMap = vtpData.pointMap();
727

728
        auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
729

730
731
732
733
        vtkpoints->SetNumberOfPoints(pointMap.size());
        forAll(pointMap, pointi)
        {
            vtkpoints->SetPoint(pointi, points[pointMap[pointi]].v_);
734
        }
735
736
737

        vtkgeom->SetPoints(vtkpoints);
        vtpData.set(vtkgeom);
738
739
740
741
    }

    if (debug)
    {
742
        Info<< "<end> " << FUNCTION_NAME << nl;
743
744
745
746
747
        printMemory();
    }
}


748
void Foam::vtkPVFoam::convertMeshPointSets()
749
{
750
751
752
753
    const List<label> partIds =
        rangePointSets_.intersection(selectedPartIds_);

    const fvMesh& mesh = *volMeshPtr_;
754
755
756

    if (debug)
    {
757
        Info<< "<beg> " << FUNCTION_NAME << nl;
758
759
760
        printMemory();
    }

761
    for (const auto partId : partIds)
762
    {
763
        const auto& longName = selectedPartIds_[partId];
764
        const word partName = getFoamName(longName);
765

766
767
768
769
        foamVtpData& vtpData = cachedVtp_(longName);

        vtkSmartPointer<vtkPolyData> vtkgeom;
        if (vtpData.nPoints() && vtpData.pointMap().size())
770
        {
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
            if (meshState_ == polyMesh::UNCHANGED)
            {
                if (debug)
                {
                    Info<< "reusing " << longName << nl;
                }
                vtpData.reuse();
                continue;
            }
            else if (meshState_ == polyMesh::POINTS_MOVED)
            {
                if (debug)
                {
                    Info<< "move points " << longName << nl;
                }
                vtkgeom = vtpData.getCopy();
            }
788
789
        }

790
791
792
793
794
795
        if (!vtkgeom)
        {
            // First time, or topo change
            vtkgeom = vtkSmartPointer<vtkPolyData>::New();
            vtpData.pointMap() = pointSet(mesh, partName).sortedToc();
        }
796

797
798
        const pointField& points = mesh.points();
        const labelUList& pointMap = vtpData.pointMap();
799

800
        auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
801

802
803
        vtkpoints->SetNumberOfPoints(pointMap.size());
        forAll(pointMap, pointi)
804
        {
805
            vtkpoints->SetPoint(pointi, points[pointMap[pointi]].v_);
806
807
        }

808
809
        vtkgeom->SetPoints(vtkpoints);
        vtpData.set(vtkgeom);
810
811
812
813
    }

    if (debug)
    {
814
        Info<< "<end> " << FUNCTION_NAME << nl;
815
816
817
818
819
820
        printMemory();
    }
}


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