foamVtkFvMeshAdaptorGeom.C 8.58 KB
Newer Older
1 2 3 4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
OpenFOAM bot's avatar
OpenFOAM bot committed
5
    \\  /    A nd           | www.openfoam.com
6
     \\/     M anipulation  |
OpenFOAM bot's avatar
OpenFOAM bot committed
7
-------------------------------------------------------------------------------
8
    Copyright (C) 2017-2020 OpenCFD Ltd.
9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
-------------------------------------------------------------------------------
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/>.

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

#include "foamVtkFvMeshAdaptor.H"
29
#include "cellCellStencilObject.H"
30 31

// VTK includes
32 33 34 35
#include "vtkMultiBlockDataSet.h"
#include "vtkPolyData.h"
#include "vtkUnstructuredGrid.h"
#include "vtkSmartPointer.h"
36 37 38 39 40

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

void Foam::vtk::fvMeshAdaptor::convertGeometryInternal()
{
41 42
    // INTERNAL
    if (!usingInternal())
43 44 45 46 47
    {
        cachedVtu_.clear();
        return;
    }

48
    const auto& longName = internalName();
49 50 51 52 53 54 55 56

    foamVtuData& vtuData = cachedVtu_(longName);

    vtkSmartPointer<vtkUnstructuredGrid> vtkgeom;
    if (vtuData.nPoints())
    {
        if (meshState_ == polyMesh::UNCHANGED)
        {
57 58 59
            DebugInfo
                << "reuse " << longName << nl;

60
            vtuData.reuse();  // No movement - simply reuse
61 62 63 64
            return;
        }
        else if (meshState_ == polyMesh::POINTS_MOVED)
        {
65 66 67
            DebugInfo
                << "move points " << longName << nl;

68 69 70 71 72 73 74
            vtkgeom = vtuData.getCopy();
            vtkgeom->SetPoints(vtuData.points(mesh_));
        }
    }

    if (!vtkgeom)
    {
75 76
        DebugInfo
            << "Nothing usable from cache - create new geometry" << nl;
77 78 79 80 81 82 83 84 85

        // Nothing usable from cache - create new geometry
        vtkgeom = vtuData.internal(mesh_, decomposePoly_);
    }

    vtuData.set(vtkgeom);
}


86
void Foam::vtk::fvMeshAdaptor::convertGeometryBoundary()
87
{
88
    // BOUNDARY
89

90 91
    const polyBoundaryMesh& patches = mesh_.boundaryMesh();

92
    for (const label patchId : patchIds_)
93 94 95 96 97 98 99 100 101 102 103
    {
        const polyPatch& pp = patches[patchId];
        const word& longName = pp.name();

        foamVtpData& vtpData = cachedVtp_(longName);

        vtkSmartPointer<vtkPolyData> vtkgeom;
        if (vtpData.nPoints())
        {
            if (meshState_ == polyMesh::UNCHANGED)
            {
104 105 106
                DebugInfo
                    << "reuse " << longName << nl;

107
                vtpData.reuse();  // No movement - simply reuse
108 109 110 111
                continue;
            }
            else if (meshState_ == polyMesh::POINTS_MOVED)
            {
112
                // Point movement on single patch
113 114 115
                DebugInfo
                    << "move points " << longName << nl;

116 117
                vtkgeom = vtpData.getCopy();
                vtkgeom->SetPoints(vtk::Tools::Patch::points(pp));
118 119 120 121 122
            }
        }

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

123 124 125
        DebugInfo
            << "Creating VTK mesh for patch [" << patchId <<"] "
            << longName << endl;
126

127 128
        // Unused information
        vtpData.additionalIds().clear();
129 130 131 132

        // This is somewhat inconsistent, since we currently only have
        // normal (non-grouped) patches but this may change in the future.

133
        vtkgeom = vtk::Tools::Patch::mesh(pp);
134 135 136 137 138 139 140 141 142 143 144 145 146 147

        if (vtkgeom)
        {
            vtpData.set(vtkgeom);
        }
        else
        {
            // Catch any problems
            cachedVtp_.erase(longName);
        }
    }
}


148
void Foam::vtk::fvMeshAdaptor::applyGhostingInternal(const labelUList& types)
149
{
150
    // INTERNAL
151

152
    if (types.empty() || !usingInternal())
153 154 155 156
    {
        return;
    }

157
    const auto& longName = internalName();
158 159

    auto iter = cachedVtu_.find(longName);
160
    if (!iter.found() || !iter.val().dataset)
161 162
    {
        // Should not happen, but for safety require a vtk geometry
163
        Pout<<"Cache miss for VTU " << longName << endl;
164 165
        return;
    }
166
    foamVtuData& vtuData = iter.val();
167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182
    auto dataset = vtuData.dataset;

    const labelUList& cellMap = vtuData.cellMap();

    auto vtkgcell = dataset->GetCellGhostArray();
    if (!vtkgcell)
    {
        vtkgcell = dataset->AllocateCellGhostArray();
    }

    // auto vtkgpoint = dataset->GetPointGhostArray();
    // if (!vtkgpoint)
    // {
    //     vtkgpoint = dataset->AllocatePointGhostArray();
    // }

183
    UList<uint8_t> gcell = vtk::Tools::asUList(vtkgcell, cellMap.size());
184 185 186 187 188 189 190

    vtkgcell->FillValue(0); // Initialize to zero

    forAll(cellMap, i)
    {
        const label cellType = types[cellMap[i]];

191
        if (cellType == cellCellStencil::HOLE)
192
        {
193 194
            // Need duplicate (not just HIDDENCELL) for it to be properly
            // recognized
195 196 197
            gcell[i] |=
            (
                vtkDataSetAttributes::DUPLICATECELL
198
              | vtkDataSetAttributes::HIDDENCELL
199 200
            );
        }
201 202 203 204 205
        #if 0
        // No special treatment for INTERPOLATED.
        // This causes duplicate/overlapping values, but avoids holes
        // in the results
        else if (cellType == cellCellStencil::INTERPOLATED)
206 207 208 209 210 211
        {
            gcell[i] |=
            (
                vtkDataSetAttributes::DUPLICATECELL
            );
        }
212
        #endif
213 214 215 216 217 218
    }

    dataset->GetCellData()->AddArray(vtkgcell);
}


219
void Foam::vtk::fvMeshAdaptor::applyGhostingBoundary(const labelUList& types)
220
{
221
    // BOUNDARY
222

223
    if (types.empty() || !usingBoundary())
224 225 226 227 228 229
    {
        return;
    }

    const polyBoundaryMesh& patches = mesh_.boundaryMesh();

230
    for (const label patchId : patchIds_)
231 232 233 234 235
    {
        const polyPatch& pp = patches[patchId];
        const word& longName = pp.name();

        auto iter = cachedVtp_.find(longName);
236
        if (!iter.found() || !iter.val().dataset)
237 238 239 240 241 242
        {
            // Should not happen, but for safety require a vtk geometry
            Pout<<"Cache miss for VTP patch " << longName << endl;
            continue;
        }

243
        foamVtpData& vtpData = iter.val();
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 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301
        auto dataset = vtpData.dataset;

        auto vtkgcell = dataset->GetCellGhostArray();
        if (!vtkgcell)
        {
            vtkgcell = dataset->AllocateCellGhostArray();
        }


        // Determine face ghosting based on interior cells
        const labelUList& bCells = pp.faceCells();

        const label len = bCells.size();

        UList<uint8_t> gcell = vtk::Tools::asUList(vtkgcell, len);

        vtkgcell->FillValue(0); // Initialize to zero

        for (label i=0; i < len; ++i)
        {
            const label celli = bCells[i];

            const label cellType = types[celli];

            if (cellType == cellCellStencil::HOLE)
            {
                // Need duplicate (not just HIDDENCELL) for it to be properly
                // recognized
                gcell[i] |=
                (
                    vtkDataSetAttributes::DUPLICATECELL
                  | vtkDataSetAttributes::HIDDENCELL
                );
            }
            #if 0
            // No special treatment for INTERPOLATED.
            // This causes duplicate/overlapping values, but avoids holes
            // in the results
            else if (cellType == cellCellStencil::INTERPOLATED)
            {
                gcell[i] |=
                (
                    vtkDataSetAttributes::DUPLICATECELL
                );
            }
            #endif
        }

        dataset->GetCellData()->AddArray(vtkgcell);
    }
}


// These need to be rebuild here, since the component mesh may just have point
// motion without topology changes.
void Foam::vtk::fvMeshAdaptor::applyGhosting()
{
    const auto* stencilPtr =
302
        mesh_.cfindObject<cellCellStencilObject>
303 304 305 306 307 308 309 310 311
        (
            cellCellStencilObject::typeName
        );

    if (stencilPtr)
    {
        const labelUList& types = stencilPtr->cellTypes();

        applyGhostingInternal(types);
312
        applyGhostingBoundary(types);
313 314 315 316
    }
}


317
// ************************************************************************* //