vtkPVFoam.C 26 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 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 27
#include "vtkPVFoam.H"
#include "vtkPVFoamReader.h"
28 29

// OpenFOAM includes
30 31
#include "areaFaMesh.H"
#include "faMesh.H"
32
#include "fvMesh.H"
33
#include "foamVersion.H"
34 35
#include "Time.H"
#include "patchZones.H"
36
#include "IOobjectList.H"
37
#include "collatedFileOperation.H"
38 39 40 41 42 43 44

// VTK includes
#include "vtkDataArraySelection.h"
#include "vtkMultiBlockDataSet.h"
#include "vtkRenderer.h"
#include "vtkTextActor.h"
#include "vtkTextProperty.h"
45
#include "vtkSmartPointer.h"
46
#include "vtkInformation.h"
47

48 49 50
// Templates (only needed here)
#include "vtkPVFoamUpdateTemplates.C"

51 52 53 54
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

namespace Foam
{
55
    defineTypeNameAndDebug(vtkPVFoam, 0);
56 57 58 59 60 61 62 63 64 65 66 67 68 69

    // file-scope
    static word updateStateName(polyMesh::readUpdateState state)
    {
        switch (state)
        {
            case polyMesh::UNCHANGED:      return "UNCHANGED";
            case polyMesh::POINTS_MOVED:   return "POINTS_MOVED";
            case polyMesh::TOPO_CHANGE:    return "TOPO_CHANGE";
            case polyMesh::TOPO_PATCH_CHANGE: return "TOPO_PATCH_CHANGE";
        };

        return "UNKNOWN";
    }
70 71
}

72 73 74

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

75
namespace Foam
76
{
77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105
    // file-scope

    //- Create a text actor
    vtkSmartPointer<vtkTextActor> createTextActor
    (
        const std::string& s,
        const Foam::point& pt
    )
    {
        vtkSmartPointer<vtkTextActor> txt =
            vtkSmartPointer<vtkTextActor>::New();

        txt->SetInput(s.c_str());

        // Set text properties
        vtkTextProperty* tprop = txt->GetTextProperty();
        tprop->SetFontFamilyToArial();
        tprop->BoldOn();
        tprop->ShadowOff();
        tprop->SetLineSpacing(1.0);
        tprop->SetFontSize(14);
        tprop->SetColor(1.0, 0.0, 1.0);
        tprop->SetJustificationToCentered();

        txt->GetPositionCoordinate()->SetCoordinateSystemToWorld();
        txt->GetPositionCoordinate()->SetValue(pt.x(), pt.y(), pt.z());

        return txt;
    }
106 107 108
}


109 110
// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //

111
void Foam::vtkPVFoam::resetCounters()
112 113
{
    // Reset array range information (ids and sizes)
114 115
    rangeVolume_.reset();
    rangePatches_.reset();
116
    rangeClouds_.reset();
117 118 119 120 121 122
    rangeCellZones_.reset();
    rangeFaceZones_.reset();
    rangePointZones_.reset();
    rangeCellSets_.reset();
    rangeFaceSets_.reset();
    rangePointSets_.reset();
123 124 125
}


126 127 128 129 130 131 132 133 134 135 136 137 138
template<class Container>
bool Foam::vtkPVFoam::addOutputBlock
(
    vtkMultiBlockDataSet* output,
    const HashTable<Container, string>& cache,
    const arrayRange& selector,
    const bool singleDataset
) const
{
    const auto blockNo = output->GetNumberOfBlocks();
    vtkSmartPointer<vtkMultiBlockDataSet> block;
    int datasetNo = 0;

139 140 141
    const List<label> partIds = selector.intersection(selectedPartIds_);

    for (const auto partId : partIds)
142
    {
143 144 145 146 147
        const auto& longName = selectedPartIds_[partId];
        const word shortName = getFoamName(longName);

        auto iter = cache.find(longName);
        if (iter.found() && iter.object().dataset)
148
        {
149
            auto dataset = iter.object().dataset;
150

151
            if (singleDataset)
152
            {
153 154
                output->SetBlock(blockNo, dataset);
                output->GetMetaData(blockNo)->Set
155 156 157 158 159 160
                (
                    vtkCompositeDataSet::NAME(),
                    shortName.c_str()
                );

                ++datasetNo;
161
                break;
162
            }
163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181
            else if (datasetNo == 0)
            {
                block = vtkSmartPointer<vtkMultiBlockDataSet>::New();
                output->SetBlock(blockNo, block);
                output->GetMetaData(blockNo)->Set
                (
                    vtkCompositeDataSet::NAME(),
                    selector.name()
                );
            }

            block->SetBlock(datasetNo, dataset);
            block->GetMetaData(datasetNo)->Set
            (
                vtkCompositeDataSet::NAME(),
                shortName.c_str()
            );

            ++datasetNo;
182 183 184 185 186 187 188
        }
    }

    return datasetNo;
}


189
int Foam::vtkPVFoam::setTime(const std::vector<double>& requestTimes)
190
{
191 192 193 194 195
    if (requestTimes.empty())
    {
        return 0;
    }

196 197
    Time& runTime = dbPtr_();

198 199
    // Get times list. Flush first to force refresh.
    fileHandler().flush();
200 201 202
    instantList Times = runTime.times();

    int nearestIndex = timeIndex_;
203
    for (const double& timeValue : requestTimes)
204
    {
205
        const int index = Time::findClosestTimeIndex(Times, timeValue);
206 207 208 209 210 211 212 213 214 215 216 217 218 219
        if (index >= 0 && index != timeIndex_)
        {
            nearestIndex = index;
            break;
        }
    }

    if (nearestIndex < 0)
    {
        nearestIndex = 0;
    }

    if (debug)
    {
220
        Info<< "<beg> setTime(";
221 222
        unsigned reqi = 0;
        for (const double& timeValue : requestTimes)
223
        {
224 225 226
            if (reqi) Info<< ", ";
            Info<< timeValue;
            ++reqi;
227 228
        }
        Info<< ") - previousIndex = " << timeIndex_
229
            << ", nearestIndex = " << nearestIndex << nl;
230 231
    }

232
    // See what has changed
233 234 235 236 237
    if (timeIndex_ != nearestIndex)
    {
        timeIndex_ = nearestIndex;
        runTime.setTime(Times[nearestIndex], nearestIndex);

238 239 240 241 242 243 244
        // When mesh changes, so do fields
        meshState_ =
        (
            volMeshPtr_
          ? volMeshPtr_->readUpdate()
          : polyMesh::TOPO_CHANGE
        );
245 246 247 248 249 250 251 252 253

        reader_->UpdateProgress(0.05);

        // this seems to be needed for catching Lagrangian fields
        updateInfo();
    }

    if (debug)
    {
254
        Info<< "<end> setTime() - selectedTime="
255 256
            << Times[nearestIndex].name() << " index=" << timeIndex_
            << "/" << Times.size()
257
            << " meshUpdateState=" << updateStateName(meshState_) << nl;
258 259 260 261 262 263
    }

    return nearestIndex;
}


264
Foam::word Foam::vtkPVFoam::getReaderPartName(const int partId) const
265
{
266
    return getFoamName(reader_->GetPartArrayName(partId));
267 268 269
}


270 271
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

272
Foam::vtkPVFoam::vtkPVFoam
273
(
274
    const char* const vtkFileName,
275
    vtkPVFoamReader* reader
276 277 278
)
:
    reader_(reader),
279
    dbPtr_(nullptr),
280 281
    volMeshPtr_(nullptr),
    areaMeshPtr_(nullptr),
282 283 284
    meshRegion_(polyMesh::defaultRegion),
    meshDir_(polyMesh::meshSubDir),
    timeIndex_(-1),
285
    decomposePoly_(false),
286
    meshState_(polyMesh::TOPO_CHANGE),
287 288
    rangeVolume_("volMesh"),
    rangeArea_("areaMesh"),
289
    rangePatches_("patch"),
290
    rangeClouds_("lagrangian"),
291 292 293 294 295 296
    rangeCellZones_("cellZone"),
    rangeFaceZones_("faceZone"),
    rangePointZones_("pointZone"),
    rangeCellSets_("cellSet"),
    rangeFaceSets_("faceSet"),
    rangePointSets_("pointSet")
297 298 299
{
    if (debug)
    {
300
        Info<< "vtkPVFoam - " << vtkFileName << nl;
301 302 303
        printMemory();
    }

304 305
    fileName FileName(vtkFileName);

306
    // avoid argList and get rootPath/caseName directly from the file
307
    fileName fullCasePath(FileName.path());
308 309 310 311 312 313 314 315 316 317

    if (!isDir(fullCasePath))
    {
        return;
    }
    if (fullCasePath == ".")
    {
        fullCasePath = cwd();
    }

318 319 320
    // OPENFOAM API
    setEnv("FOAM_API", std::to_string(foamVersion::api), true);

321 322 323
    // The name of the executable, unless already present in the environment
    setEnv("FOAM_EXECUTABLE", "paraview", false);

324
    // Set the case as an environment variable - some BCs might use this
325 326 327 328 329 330 331 332 333
    if (fullCasePath.name().find("processors", 0) == 0)
    {
        // FileName e.g. "cavity/processors256/processor1.OpenFOAM
        // Remove the processors section so it goes into processorDDD
        // checking below.
        fullCasePath = fullCasePath.path()/fileName(FileName.name()).lessExt();
    }


334 335
    if (fullCasePath.name().find("processor", 0) == 0)
    {
336 337 338
        // Give filehandler opportunity to analyse number of processors
        (void)fileHandler().filePath(fullCasePath);

339 340 341 342 343 344 345 346 347 348 349 350 351 352 353
        const fileName globalCase = fullCasePath.path();

        setEnv("FOAM_CASE", globalCase, true);
        setEnv("FOAM_CASENAME", globalCase.name(), true);
    }
    else
    {
        setEnv("FOAM_CASE", fullCasePath, true);
        setEnv("FOAM_CASENAME", fullCasePath.name(), true);
    }

    // look for 'case{region}.OpenFOAM'
    // could be stringent and insist the prefix match the directory name...
    // Note: cannot use fileName::name() due to the embedded '{}'
    string caseName(fileName(FileName).lessExt());
354 355
    const auto beg = caseName.find_last_of("/{");
    const auto end = caseName.find('}', beg);
356 357 358

    if
    (
359 360
        beg != std::string::npos && caseName[beg] == '{'
     && end != std::string::npos && end == caseName.size()-1
361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381
    )
    {
        meshRegion_ = caseName.substr(beg+1, end-beg-1);

        // some safety
        if (meshRegion_.empty())
        {
            meshRegion_ = polyMesh::defaultRegion;
        }

        if (meshRegion_ != polyMesh::defaultRegion)
        {
            meshDir_ = meshRegion_/polyMesh::meshSubDir;
        }
    }

    if (debug)
    {
        Info<< "fullCasePath=" << fullCasePath << nl
            << "FOAM_CASE=" << getEnv("FOAM_CASE") << nl
            << "FOAM_CASENAME=" << getEnv("FOAM_CASENAME") << nl
382
            << "region=" << meshRegion_ << nl;
383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403
    }

    // Create time object
    dbPtr_.reset
    (
        new Time
        (
            Time::controlDictName,
            fileName(fullCasePath.path()),
            fileName(fullCasePath.name())
        )
    );

    dbPtr_().functionObjects().off();

    updateInfo();
}


// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

404
Foam::vtkPVFoam::~vtkPVFoam()
405 406 407
{
    if (debug)
    {
408
        Info<< "~vtkPVFoam" << nl;
409 410
    }

411 412
    delete volMeshPtr_;
    delete areaMeshPtr_;
413 414 415 416 417
}


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

418
void Foam::vtkPVFoam::updateInfo()
419 420 421
{
    if (debug)
    {
422
        Info<< "<beg> updateInfo"
423 424 425
            << " [volMeshPtr=" << (volMeshPtr_ ? "set" : "nullptr")
            << "] timeIndex="
            << timeIndex_ << nl;
426 427 428 429
    }

    resetCounters();

430 431 432
    // Part selection
    {
        vtkDataArraySelection* select = reader_->GetPartSelection();
433

434 435 436 437 438 439 440
        // There are two ways to ensure we have the correct list of parts:
        // 1. remove everything and then set particular entries 'on'
        // 2. build a 'char **' list and call SetArraysWithDefault()
        //
        // Nr. 2 has the potential advantage of not touching the modification
        // time of the vtkDataArraySelection, but the qt/paraview proxy
        // layer doesn't care about that anyhow.
441

442 443 444 445 446 447 448 449 450 451 452
        HashSet<string> enabled;
        if (!select->GetNumberOfArrays() && !volMeshPtr_)
        {
            // Fake enable 'internalMesh' on the first call
            enabled = { "internalMesh" };
        }
        else
        {
            // Preserve the enabled selections
            enabled = getSelectedArraySet(select);
        }
453

454
        select->RemoveAllArrays();   // Clear existing list
455

456 457 458 459 460 461 462
        // Update mesh parts list - add Lagrangian at the bottom
        updateInfoInternalMesh(select);
        updateInfoPatches(select, enabled);
        updateInfoSets(select);
        updateInfoZones(select);
        updateInfoAreaMesh(select);
        updateInfoLagrangian(select);
463

464 465
        setSelectedArrayEntries(select, enabled);  // Adjust/restore selected
    }
466

467 468 469 470 471 472 473 474
    // Volume and area fields - includes save/restore of selected
    updateInfoContinuumFields(reader_->GetVolFieldSelection());

    // Point fields - includes save/restore of selected
    updateInfoPointFields(reader_->GetPointFieldSelection());

    // Lagrangian fields - includes save/restore of selected
    updateInfoLagrangianFields(reader_->GetLagrangianFieldSelection());
475 476 477

    if (debug)
    {
478
        Info<< "<end> updateInfo" << nl;
479 480 481 482
    }
}


483 484 485 486 487
void Foam::vtkPVFoam::Update
(
    vtkMultiBlockDataSet* output,
    vtkMultiBlockDataSet* outputLagrangian
)
488 489 490
{
    if (debug)
    {
491 492 493
        cout<< "<beg> Foam::vtkPVFoam::Update\n";
        output->Print(cout);
        if (outputLagrangian) outputLagrangian->Print(cout);
494 495
        printMemory();
    }
496
    reader_->UpdateProgress(0.1);
497

498 499 500 501
    const int caching = reader_->GetMeshCaching();
    const bool oldDecomp = decomposePoly_;
    decomposePoly_ = !reader_->GetUseVTKPolyhedron();

502
    // Set up mesh parts selection(s)
503
    // Update cached, saved, unneed values.
504
    {
505
        vtkDataArraySelection* selection = reader_->GetPartSelection();
506
        const int n = selection->GetNumberOfArrays();
507

508
        selectedPartIds_.clear();
509
        HashSet<string> nowActive;
510

511
        for (int id=0; id < n; ++id)
512
        {
513 514
            const string str(selection->GetArrayName(id));
            const bool status = selection->GetArraySetting(id);
515 516 517 518

            if (status)
            {
                selectedPartIds_.set(id, str); // id -> name
519
                nowActive.set(str);
520 521 522 523
            }

            if (debug > 1)
            {
524 525
                Info<< "  part[" << id << "] = " << status
                    << " : " << str << nl;
526 527
            }
        }
528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560

        // Dispose of unneeded components
        cachedVtp_.retain(nowActive);
        cachedVtu_.retain(nowActive);

        if
        (
            !caching
         || meshState_ == polyMesh::TOPO_CHANGE
         || meshState_ == polyMesh::TOPO_PATCH_CHANGE
        )
        {
            // Eliminate cached values that would be unreliable
            forAllIters(cachedVtp_, iter)
            {
                iter.object().clearGeom();
                iter.object().clear();
            }
            forAllIters(cachedVtu_, iter)
            {
                iter.object().clearGeom();
                iter.object().clear();
            }
        }
        else if (oldDecomp != decomposePoly_)
        {
            // poly-decompose changed - dispose of cached values
            forAllIters(cachedVtu_, iter)
            {
                iter.object().clearGeom();
                iter.object().clear();
            }
        }
561 562
    }

563 564 565
    reader_->UpdateProgress(0.15);

    // Update the OpenFOAM mesh
566
    {
567 568
        if (debug)
        {
569
            Info<< "<beg> updateFoamMesh" << nl;
570 571
            printMemory();
        }
572

573
        if (!caching)
574
        {
575 576 577 578 579
            delete volMeshPtr_;
            delete areaMeshPtr_;

            volMeshPtr_ = nullptr;
            areaMeshPtr_ = nullptr;
580
        }
581

582
        // Check to see if the OpenFOAM mesh has been created
583
        if (!volMeshPtr_)
584 585 586 587
        {
            if (debug)
            {
                Info<< "Creating OpenFOAM mesh for region " << meshRegion_
588
                    << " at time=" << dbPtr_().timeName() << nl;
589
            }
590

591
            volMeshPtr_ = new fvMesh
592 593 594 595 596 597 598 599 600 601
            (
                IOobject
                (
                    meshRegion_,
                    dbPtr_().timeName(),
                    dbPtr_(),
                    IOobject::MUST_READ
                )
            );

602
            meshState_ = polyMesh::TOPO_CHANGE; // New mesh
603 604 605 606 607
        }
        else
        {
            if (debug)
            {
608
                Info<< "Using existing OpenFOAM mesh" << nl;
609 610 611
            }
        }

612 613 614 615 616 617 618 619 620 621 622 623 624 625
        if (rangeArea_.intersects(selectedPartIds_))
        {
            if (!areaMeshPtr_)
            {
                areaMeshPtr_ = new faMesh(*volMeshPtr_);
            }
        }
        else
        {
            delete areaMeshPtr_;

            areaMeshPtr_ = nullptr;
        }

626 627
        if (debug)
        {
628
            Info<< "<end> updateFoamMesh" << nl;
629 630 631
            printMemory();
        }
    }
632 633 634

    reader_->UpdateProgress(0.4);

635 636
    convertMeshVolume();
    convertMeshPatches();
637 638 639 640
    reader_->UpdateProgress(0.6);

    if (reader_->GetIncludeZones())
    {
641 642 643
        convertMeshCellZones();
        convertMeshFaceZones();
        convertMeshPointZones();
644 645 646 647 648
        reader_->UpdateProgress(0.65);
    }

    if (reader_->GetIncludeSets())
    {
649 650 651
        convertMeshCellSets();
        convertMeshFaceSets();
        convertMeshPointSets();
652 653 654
        reader_->UpdateProgress(0.7);
    }

655 656
    convertMeshArea();

657
    convertMeshLagrangian();
658 659 660 661

    reader_->UpdateProgress(0.8);

    // Update fields
662 663
    convertVolFields();
    convertPointFields();
664
    convertAreaFields();
665

666
    convertLagrangianFields();
667 668 669 670 671 672 673 674 675 676

    // Assemble multiblock output
    addOutputBlock(output, cachedVtu_, rangeVolume_, true); // One dataset
    addOutputBlock(output, cachedVtp_, rangePatches_);
    addOutputBlock(output, cachedVtu_, rangeCellZones_);
    addOutputBlock(output, cachedVtp_, rangeFaceZones_);
    addOutputBlock(output, cachedVtp_, rangePointZones_);
    addOutputBlock(output, cachedVtu_, rangeCellSets_);
    addOutputBlock(output, cachedVtp_, rangeFaceSets_);
    addOutputBlock(output, cachedVtp_, rangePointSets_);
677
    addOutputBlock(output, cachedVtp_, rangeArea_);
678 679 680 681
    addOutputBlock
    (
        (outputLagrangian ? outputLagrangian : output),
        cachedVtp_,
682
        rangeClouds_
683
    );
684

685 686
    if (debug)
    {
687
        Info<< "done reader part" << nl << nl;
688 689 690
    }
    reader_->UpdateProgress(0.95);

691
    meshState_ = polyMesh::UNCHANGED;
692

693 694 695 696 697 698 699 700 701 702 703 704 705 706
    if (caching & 2)
    {
        // Suppress caching of Lagrangian since it normally always changes.
        cachedVtp_.filterKeys
        (
            [](const word& k){ return k.startsWith("lagrangian/"); },
            true // prune
        );
    }
    else
    {
        cachedVtp_.clear();
        cachedVtu_.clear();
    }
707 708 709
}


710
void Foam::vtkPVFoam::UpdateFinalize()
711
{
712
    if (!reader_->GetMeshCaching())
713
    {
714 715 716 717 718
        delete volMeshPtr_;
        delete areaMeshPtr_;

        volMeshPtr_ = nullptr;
        areaMeshPtr_ = nullptr;
719 720
    }

721 722 723 724
    reader_->UpdateProgress(1.0);
}


725
std::vector<double> Foam::vtkPVFoam::findTimes(const bool skipZero) const
726
{
727
    std::vector<double> times;
728 729 730

    if (dbPtr_.valid())
    {
731
        const Time& runTime = dbPtr_();
732 733
        // Get times list. Flush first to force refresh.
        fileHandler().flush();
734 735 736
        instantList timeLst = runTime.times();

        // find the first time for which this mesh appears to exist
737
        label begIndex = timeLst.size();
738
        forAll(timeLst, timei)
739 740 741
        {
            if
            (
742
                IOobject
743 744
                (
                    "points",
745
                    timeLst[timei].name(),
746 747
                    meshDir_,
                    runTime
748
                ).typeHeaderOk<pointIOField>(false, false)
749 750
            )
            {
751
                begIndex = timei;
752 753 754 755
                break;
            }
        }

756
        label nTimes = timeLst.size() - begIndex;
757 758

        // skip "constant" time whenever possible
759
        if (begIndex == 0 && nTimes > 1)
760
        {
761
            if (timeLst[begIndex].name() == runTime.constant())
762
            {
763
                ++begIndex;
764 765 766 767 768
                --nTimes;
            }
        }

        // skip "0/" time if requested and possible
769
        if (skipZero && nTimes > 1 && timeLst[begIndex].name() == "0")
770
        {
771 772
            ++begIndex;
            --nTimes;
773 774
        }

775 776
        times.reserve(nTimes);
        while (nTimes-- > 0)
777
        {
778
            times.push_back(timeLst[begIndex++].value());
779 780 781 782 783 784 785 786 787 788
        }
    }
    else
    {
        if (debug)
        {
            cout<< "no valid dbPtr:\n";
        }
    }

789
    return times;
790 791 792
}


793
void Foam::vtkPVFoam::renderPatchNames
794 795 796 797 798
(
    vtkRenderer* renderer,
    const bool show
)
{
799
    // Always remove old actors first
800

801
    for (auto& actor : patchTextActors_)
802
    {
803
        renderer->RemoveViewProp(actor);
804
    }
805
    patchTextActors_.clear();
806

807
    if (show && volMeshPtr_)
808
    {
809
        // get the display patches, strip off any prefix/suffix
810
        hashedWordList selectedPatches = getSelected
811 812
        (
            reader_->GetPartSelection(),
813
            rangePatches_
814 815 816 817 818 819 820
        );

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

821
        const polyBoundaryMesh& pbMesh = volMeshPtr_->boundaryMesh();
822 823 824 825

        // Find the total number of zones
        // Each zone will take the patch name
        // Number of zones per patch ... zero zones should be skipped
826
        labelList nZones(pbMesh.size(), Zero);
827 828

        // Per global zone number the average face centre position
829
        List<DynamicList<point>> zoneCentre(pbMesh.size());
830 831 832


        // Loop through all patches to determine zones, and centre of each zone
833
        forAll(pbMesh, patchi)
834
        {
835
            const polyPatch& pp = pbMesh[patchi];
836 837 838 839 840 841 842 843 844 845 846 847

            // Only include the patch if it is selected
            if (!selectedPatches.found(pp.name()))
            {
                continue;
            }

            const labelListList& edgeFaces = pp.edgeFaces();
            const vectorField& n = pp.faceNormals();

            boolList featEdge(pp.nEdges(), false);

848
            forAll(edgeFaces, edgei)
849
            {
850
                const labelList& eFaces = edgeFaces[edgei];
851 852 853 854 855

                if (eFaces.size() == 1)
                {
                    // Note: could also do ones with > 2 faces but this gives
                    // too many zones for baffles
856
                    featEdge[edgei] = true;
857 858 859
                }
                else if (mag(n[eFaces[0]] & n[eFaces[1]]) < 0.5)
                {
860
                    featEdge[edgei] = true;
861 862 863 864 865 866
                }
            }

            // Do topological analysis of patch, find disconnected regions
            patchZones pZones(pp, featEdge);

867
            nZones[patchi] = pZones.nZones();
868

869
            labelList zoneNFaces(pZones.nZones(), Zero);
870 871

            // Create storage for additional zone centres
872
            forAll(zoneNFaces, zonei)
873
            {
874
                zoneCentre[patchi].append(Zero);
875 876 877
            }

            // Do averaging per individual zone
878
            forAll(pp, facei)
879
            {
880 881 882
                const label zonei = pZones[facei];
                zoneCentre[patchi][zonei] += pp[facei].centre(pp.points());
                zoneNFaces[zonei]++;
883 884
            }

885
            forAll(zoneCentre[patchi], zonei)
886
            {
887
                zoneCentre[patchi][zonei] /= zoneNFaces[zonei];
888 889 890 891 892 893 894 895 896 897
            }
        }

        // Count number of zones we're actually going to display.
        // This is truncated to a max per patch

        const label MAXPATCHZONES = 20;

        label displayZoneI = 0;

898
        forAll(pbMesh, patchi)
899
        {
900
            displayZoneI += min(MAXPATCHZONES, nZones[patchi]);
901 902 903 904 905
        }

        if (debug)
        {
            Info<< "displayed zone centres = " << displayZoneI << nl
906
                << "zones per patch = " << nZones << nl;
907 908 909
        }

        // Set the size of the patch labels to max number of zones
910
        patchTextActors_.setSize(displayZoneI);
911 912 913

        if (debug)
        {
914
            Info<< "constructing patch labels" << nl;
915 916 917 918 919
        }

        // Actor index
        displayZoneI = 0;

920
        forAll(pbMesh, patchi)
921
        {
922
            const polyPatch& pp = pbMesh[patchi];
923 924

            // Only selected patches will have a non-zero number of zones
925
            const label nDisplayZones = min(MAXPATCHZONES, nZones[patchi]);
926
            label increment = 1;
927
            if (nZones[patchi] >= MAXPATCHZONES)
928
            {
929
                increment = nZones[patchi]/MAXPATCHZONES;
930 931
            }

932 933
            label globalZoneI = 0;
            for (label i = 0; i < nDisplayZones; ++i, globalZoneI += increment)
934 935 936 937
            {
                if (debug)
                {
                    Info<< "patch name = " << pp.name() << nl
938
                        << "anchor = " << zoneCentre[patchi][globalZoneI] << nl
939
                        << "globalZoneI = " << globalZoneI << nl;
940 941
                }

942
                // Into a list for later removal
943
                patchTextActors_[displayZoneI++] = createTextActor
944
                (
945 946
                    pp.name(),
                    zoneCentre[patchi][globalZoneI]
947 948 949 950 951
                );
            }
        }

        // Resize the patch names list to the actual number of patch names added
952
        patchTextActors_.setSize(displayZoneI);
953
    }
954 955

    // Add text to each renderer
956
    for (auto& actor : patchTextActors_)
957
    {
958
        renderer->AddViewProp(actor);
959
    }
960 961 962
}


963
void Foam::vtkPVFoam::PrintSelf(ostream& os, vtkIndent indent) const
964 965
{
    os  << indent << "Number of nodes: "
966
        << (volMeshPtr_ ? volMeshPtr_->nPoints() : 0) << "\n";
967 968

    os  << indent << "Number of cells: "
969
        << (volMeshPtr_ ? volMeshPtr_->nCells() : 0) << "\n";
970 971 972 973 974 975 976 977

    os  << indent << "Number of available time steps: "
        << (dbPtr_.valid() ? dbPtr_().times().size() : 0) << "\n";

    os  << indent << "mesh region: " << meshRegion_ << "\n";
}


978 979 980 981
void Foam::vtkPVFoam::printInfo() const
{
    std::cout
        << "Region: " << meshRegion_ << "\n"
982 983
        << "nPoints: " << (volMeshPtr_ ? volMeshPtr_->nPoints() : 0) << "\n"
        << "nCells: "  << (volMeshPtr_ ? volMeshPtr_->nCells() : 0) << "\n"
984 985 986 987 988 989 990 991 992 993 994 995 996 997
        << "nTimes: "
        << (dbPtr_.valid() ? dbPtr_().times().size() : 0) << "\n";

    std::vector<double> times = this->findTimes(reader_->GetSkipZeroTime());

    std::cout<<"  " << times.size() << "(";
    for (const double& val : times)
    {
        std::cout<< ' ' << val;
    }
    std::cout << " )" << std::endl;
}


998
// ************************************************************************* //