vtkPVFoamTemplates.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2025 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "vtkPVFoam.H"
27 #include "vtkOpenFOAMPoints.H"
28 
29 // OpenFOAM includes
30 #include "polyPatch.H"
31 #include "primitivePatch.H"
32 
33 // VTK includes
34 #include "vtkCellArray.h"
35 #include "vtkPoints.h"
36 #include "vtkPolyData.h"
37 
38 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
39 
40 template<class PatchType>
41 vtkPolyData* Foam::vtkPVFoam::patchVTKMesh(const PatchType& p)
42 {
43  vtkPolyData* vtkmesh = vtkPolyData::New();
44 
45  // Convert OpenFOAM mesh vertices to VTK
46  const Foam::pointField& points = p.localPoints();
47 
48  vtkPoints* vtkpoints = vtkPoints::New();
49  vtkpoints->Allocate(points.size());
50  forAll(points, i)
51  {
52  vtkInsertNextOpenFOAMPoint(vtkpoints, points[i]);
53  }
54 
55  vtkmesh->SetPoints(vtkpoints);
56  vtkpoints->Delete();
57 
58  // Add faces as polygons
59  const faceList& faces = p.localFaces();
60 
61  vtkCellArray* vtkcells = vtkCellArray::New();
62  vtkcells->Allocate(faces.size());
63 
64  DynamicList<vtkIdType> nodeIds(4);
65 
66  forAll(faces, facei)
67  {
68  const face& f = faces[facei];
69 
70  nodeIds.resize(f.size());
71  forAll(f, fp)
72  {
73  nodeIds[fp] = f[fp];
74  }
75 
76  vtkcells->InsertNextCell(f.size(), nodeIds.data());
77  }
78 
79  vtkmesh->SetPolys(vtkcells);
80  vtkcells->Delete();
81 
82  return vtkmesh;
83 }
84 
85 
86 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:433
const pointField & points
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
List< face > faceList
Definition: faceListFwd.H:41
labelList f(nPoints)
volScalarField & p
void vtkInsertNextOpenFOAMPoint(vtkPoints *points, const Foam::point &p)