vtkPVblockMesh.C 14.8 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-2019 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 "vtkPVblockMesh.H"
#include "vtkPVblockMeshReader.h"
28 29 30

// OpenFOAM includes
#include "blockMesh.H"
31
#include "blockMeshTools.H"
32
#include "foamVersion.H"
33 34
#include "Time.H"
#include "patchZones.H"
35
#include "StringStream.H"
36 37 38 39 40 41 42

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

// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

namespace Foam
{
49
    defineTypeNameAndDebug(vtkPVblockMesh, 0);
50 51
}

52 53
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //

54
namespace Foam
55
{
56 57 58 59 60 61 62 63
    // file-scope
    //- Create a text actor
    static vtkSmartPointer<vtkTextActor> createTextActor
    (
        const std::string& s,
        const Foam::point& pt
    )
    {
64
        auto txt = vtkSmartPointer<vtkTextActor>::New();
65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82

        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;
    }
83 84 85
}


86 87
// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //

88
void Foam::vtkPVblockMesh::resetCounters()
89 90
{
    // Reset mesh part ids and sizes
91 92 93
    rangeBlocks_.reset();
    rangeEdges_.reset();
    rangeCorners_.reset();
94 95 96
}


97
void Foam::vtkPVblockMesh::updateInfoBlocks
98
(
99
    vtkDataArraySelection* select
100 101 102 103
)
{
    if (debug)
    {
104 105
        Info<< "<beg> updateInfoBlocks"
            << " [meshPtr=" << (meshPtr_ ? "set" : "null") << "]" << endl;
106 107
    }

108
    rangeBlocks_.reset(select->GetNumberOfArrays());
109 110

    const blockMesh& blkMesh = *meshPtr_;
111

112 113 114
    const int nBlocks = blkMesh.size();
    for (int blockI = 0; blockI < nBlocks; ++blockI)
    {
115
        const blockDescriptor& blockDef = blkMesh[blockI];
116

117 118
        // Display either blockI as a number or with its name
        // (looked up from blockMeshDict)
119 120
        OStringStream ostr;
        blockDescriptor::write(ostr, blockI, blkMesh.meshDict());
121 122 123 124

        // append the (optional) zone name
        if (!blockDef.zoneName().empty())
        {
125
            ostr << " - " << blockDef.zoneName();
126 127
        }

128
        // Add "blockId" or "blockId - zoneName" to GUI list
129
        select->AddArray(ostr.str().c_str());
130
        ++rangeBlocks_;
131 132 133 134
    }

    if (debug)
    {
135
        Info<< "<end> updateInfoBlocks" << endl;
136 137 138 139
    }
}


140
void Foam::vtkPVblockMesh::updateInfoEdges
141
(
142
    vtkDataArraySelection* select
143 144 145 146
)
{
    if (debug)
    {
147 148
        Info<< "<beg> updateInfoEdges"
            << " [meshPtr=" << (meshPtr_ ? "set" : "null") << "]" << endl;
149 150
    }

151
    rangeEdges_.reset(select->GetNumberOfArrays());
152 153

    const blockMesh& blkMesh = *meshPtr_;
154
    const blockEdgeList& edges = blkMesh.edges();
155 156 157 158

    forAll(edges, edgeI)
    {
        OStringStream ostr;
159 160 161 162
        blockVertex::write(ostr, edges[edgeI].start(), blkMesh.meshDict());
        ostr<< ":";
        blockVertex::write(ostr, edges[edgeI].end(), blkMesh.meshDict());
        ostr << " - " << edges[edgeI].type();
163 164

        // Add "beg:end - type" to GUI list
165
        select->AddArray(ostr.str().c_str());
166
        ++rangeEdges_;
167 168 169 170
    }

    if (debug)
    {
171
        Info<< "<end> updateInfoEdges" << endl;
172 173 174 175 176 177
    }
}


// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

178
Foam::vtkPVblockMesh::vtkPVblockMesh
179 180
(
    const char* const FileName,
181
    vtkPVblockMeshReader* reader
182 183 184
)
:
    reader_(reader),
185 186
    dbPtr_(nullptr),
    meshPtr_(nullptr),
187 188
    meshRegion_(polyMesh::defaultRegion),
    meshDir_(polyMesh::meshSubDir),
189 190 191
    rangeBlocks_("block"),
    rangeEdges_("edges"),
    rangeCorners_("corners")
192 193 194
{
    if (debug)
    {
195
        Info<< "vtkPVblockMesh - " << FileName << endl;
196 197 198 199 200 201 202 203 204 205 206 207 208 209
    }

    // avoid argList and get rootPath/caseName directly from the file
    fileName fullCasePath(fileName(FileName).path());

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

210 211 212
    // OPENFOAM API
    setEnv("FOAM_API", std::to_string(foamVersion::api), true);

213 214 215
    // The name of the executable, unless already present in the environment
    setEnv("FOAM_EXECUTABLE", "paraview", false);

216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282
    // Set the case as an environment variable - some BCs might use this
    if (fullCasePath.name().find("processor", 0) == 0)
    {
        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());
    string::size_type beg = caseName.find_last_of("/{");
    string::size_type end = caseName.find('}', beg);

    if
    (
        beg != string::npos && caseName[beg] == '{'
     && end != string::npos && end == caseName.size()-1
    )
    {
        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") << endl;
    }

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

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

    updateInfo();
}


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

283
Foam::vtkPVblockMesh::~vtkPVblockMesh()
284 285 286
{
    if (debug)
    {
287
        Info<< "~vtkPVblockMesh" << endl;
288 289 290 291 292 293 294 295
    }

    delete meshPtr_;
}


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

296
void Foam::vtkPVblockMesh::updateInfo()
297 298 299
{
    if (debug)
    {
300 301
        Info<< "<beg> updateInfo"
            << " [meshPtr=" << (meshPtr_ ? "set" : "null") << "] " << endl;
302 303 304 305 306
    }

    resetCounters();

    vtkDataArraySelection* blockSelection = reader_->GetBlockSelection();
307
    vtkDataArraySelection* edgeSelection  = reader_->GetCurvedEdgesSelection();
308

309
    const bool firstTime = (!blockSelection->GetNumberOfArrays() && !meshPtr_);
310 311 312 313

    // Preserve the enabled selections if possible
    HashSet<string> enabledParts;
    HashSet<string> enabledEdges;
314
    if (!firstTime)
315
    {
316 317
        enabledParts = getSelectedArraySet(blockSelection);
        enabledEdges = getSelectedArraySet(edgeSelection);
318 319 320 321 322 323 324 325 326 327
    }

    // Clear current mesh parts list
    blockSelection->RemoveAllArrays();
    edgeSelection->RemoveAllArrays();

    // need a blockMesh
    updateFoamMesh();

    // Update mesh parts list
328
    updateInfoBlocks(blockSelection);
329 330

    // Update curved edges list
331
    updateInfoEdges(edgeSelection);
332

333
    // Restore the enabled selections
334 335 336
    if (!firstTime)
    {
        setSelectedArrayEntries(blockSelection, enabledParts);
337
        setSelectedArrayEntries(edgeSelection,  enabledEdges);
338 339 340 341
    }

    if (debug)
    {
342
        Info<< "<end> updateInfo" << endl;
343 344 345 346
    }
}


347
void Foam::vtkPVblockMesh::updateFoamMesh()
348 349 350
{
    if (debug)
    {
351
        Info<< "<beg> updateFoamMesh" << endl;
352 353 354 355 356 357 358 359 360 361 362
    }

    // Check to see if the OpenFOAM mesh has been created
    if (!meshPtr_)
    {
        if (debug)
        {
            Info<< "Creating blockMesh at time=" << dbPtr_().timeName()
                << endl;
        }

363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379
        // Set path for the blockMeshDict
        const word dictName("blockMeshDict");
        fileName dictPath(dbPtr_().system()/dictName);

        // Check if dictionary is present in the constant directory
        if
        (
            exists
            (
                dbPtr_().path()/dbPtr_().constant()
               /polyMesh::meshSubDir/dictName
            )
        )
        {
            dictPath = dbPtr_().constant()/polyMesh::meshSubDir/dictName;
        }

380 381 382
        // Store dictionary since is used as database inside blockMesh class
        // for names of vertices and blocks
        IOdictionary* meshDictPtr = new IOdictionary
383 384 385
        (
            IOobject
            (
386
                dictPath,
387 388 389
                dbPtr_(),
                IOobject::MUST_READ_IF_MODIFIED,
                IOobject::NO_WRITE,
390
                true
391 392
            )
        );
393
        meshDictPtr->store();
394

395
        meshPtr_ = new blockMesh(*meshDictPtr, meshRegion_);
396 397 398 399 400
    }


    if (debug)
    {
401
        Info<< "<end> updateFoamMesh" << endl;
402 403 404 405
    }
}


406
void Foam::vtkPVblockMesh::Update
407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427
(
    vtkMultiBlockDataSet* output
)
{
    reader_->UpdateProgress(0.1);

    // Update the OpenFOAM mesh
    updateFoamMesh();
    reader_->UpdateProgress(0.5);

    // Convert mesh elemente
    int blockNo = 0;

    convertMeshCorners(output, blockNo);
    convertMeshBlocks(output, blockNo);
    convertMeshEdges(output, blockNo);

    reader_->UpdateProgress(0.8);
}


428
void Foam::vtkPVblockMesh::UpdateFinalize()
429 430 431 432 433
{
    reader_->UpdateProgress(1.0);
}


434
void Foam::vtkPVblockMesh::renderPatchNames
435 436 437 438 439 440
(
    vtkRenderer* renderer,
    const bool show
)
{
    // always remove old actors first
441
    forAll(patchTextActors_, actori)
442
    {
443
        renderer->RemoveViewProp(patchTextActors_[actori]);
444
    }
445
    patchTextActors_.clear();
446 447 448

    // the number of text actors
    label nActors = 0;
449 450 451

    if (show && meshPtr_)
    {
452
        const blockMesh& blkMesh = *meshPtr_;
453
        const dictionary& meshDescription = blkMesh.meshDict();
454 455
        const pointField& cornerPts = blkMesh.vertices();
        const scalar scaleFactor = blkMesh.scaleFactor();
456

457
        if (!meshDescription.found("boundary"))
458
        {
459 460
            return;
        }
461

462
        // 8 sides per block is plenty
463
        patchTextActors_.setSize(8*blkMesh.size());
464 465 466 467 468 469 470 471 472 473 474 475 476

        // Collect all variables
        dictionary varDict(meshDescription.subOrEmptyDict("namedVertices"));
        varDict.merge(meshDescription.subOrEmptyDict("namedBlocks"));

        // Read like boundary file
        const PtrList<entry> patchesInfo(meshDescription.lookup("boundary"));

        forAll(patchesInfo, patchi)
        {
            const entry& patchInfo = patchesInfo[patchi];

            if (!patchInfo.isDict())
477
            {
478 479 480 481 482
                IOWarningInFunction(meshDescription)
                    << "Entry " << patchInfo << " in boundary section is not a"
                    << " valid dictionary."
                    << endl;
                break;
483
            }
484

485
            const word& patchName = patchInfo.keyword();
486

487 488
            // Read block faces
            faceList patchFaces = blockMeshTools::read<face>
489
            (
490 491
                patchInfo.dict().lookup("faces"),
                varDict
492 493
            );

494 495 496 497 498
            forAll(patchFaces, facei)
            {
                const face& f = patchFaces[facei];

                // Into a list for later removal
499
                patchTextActors_[nActors++] = createTextActor
500 501 502 503 504
                (
                    patchName,
                    f.centre(cornerPts) * scaleFactor
                );

505
                if (nActors == patchTextActors_.size())
506 507 508 509 510
                {
                    // hit max allocated space - bail out
                    break;
                }
            }
511

512
            if (nActors == patchTextActors_.size())
513 514 515 516
            {
                // hit max allocated space - bail out
                break;
            }
517
        }
518

519
        patchTextActors_.setSize(nActors);
520 521 522
    }

    // Add text to each renderer
523
    forAll(patchTextActors_, actori)
524
    {
525
        renderer->AddViewProp(patchTextActors_[actori]);
526 527 528 529
    }
}


530 531 532 533 534 535 536 537
void Foam::vtkPVblockMesh::renderPointNumbers
(
    vtkRenderer* renderer,
    const bool show
)
{
    // always remove old actors first

538
    forAll(pointTextActors_, actori)
539
    {
540
        renderer->RemoveViewProp(pointTextActors_[actori]);
541
    }
542
    pointTextActors_.clear();
543 544 545 546 547 548 549

    if (show && meshPtr_)
    {
        const blockMesh& blkMesh = *meshPtr_;
        const pointField& cornerPts = blkMesh.vertices();
        const scalar scaleFactor = blkMesh.scaleFactor();

550
        pointTextActors_.setSize(cornerPts.size());
551 552 553 554 555 556 557 558
        forAll(cornerPts, pointi)
        {
            // Display either pointi as a number or with its name
            // (looked up from blockMeshDict)
            OStringStream os;
            blockVertex::write(os, pointi, blkMesh.meshDict());

            // Into a list for later removal
559
            pointTextActors_[pointi] = createTextActor
560 561 562 563 564 565 566 567
            (
                os.str(),
                cornerPts[pointi]*scaleFactor
            );
        }
    }

    // Add text to each renderer
568
    forAll(pointTextActors_, actori)
569
    {
570
        renderer->AddViewProp(pointTextActors_[actori]);
571 572 573
    }
}

574

575
void Foam::vtkPVblockMesh::PrintSelf(ostream& os, vtkIndent indent) const
576 577 578 579 580 581 582 583 584 585 586 587 588 589
{
#if 0
    os  << indent << "Number of nodes: "
        << (meshPtr_ ? meshPtr_->nPoints() : 0) << "\n";

    os  << indent << "Number of cells: "
        << (meshPtr_ ? meshPtr_->nCells() : 0) << "\n";

    os  << indent << "Number of available time steps: "
        << (dbPtr_.valid() ? dbPtr_().times().size() : 0) << endl;
#endif
}

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