vtkPVblockMeshConvert.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-2018 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 "vtkPVblockMesh.H"
27 #include "vtkPVblockMeshReader.h"
28 
29 // OpenFOAM includes
30 #include "blockMesh.H"
31 #include "Time.H"
32 
33 #include "vtkOpenFOAMPoints.H"
34 
35 // VTK includes
36 #include "vtkCellArray.h"
37 #include "vtkDataArraySelection.h"
38 #include "vtkMultiBlockDataSet.h"
39 #include "vtkPoints.h"
40 #include "vtkPolyData.h"
41 #include "vtkUnstructuredGrid.h"
42 
43 
44 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
45 
46 void Foam::vtkPVblockMesh::convertMeshBlocks
47 (
48  vtkMultiBlockDataSet* output,
49  int& blockNo
50 )
51 {
52  vtkDataArraySelection* selection = reader_->GetBlockSelection();
53  arrayRange& range = arrayRangeBlocks_;
54  range.block(blockNo); // set output block
55  label datasetNo = 0; // restart at dataset 0
56 
57  const blockMesh& blkMesh = *meshPtr_;
58  const Foam::pointField& blockPoints = blkMesh.vertices();
59 
60  if (debug)
61  {
62  Info<< "<beg> Foam::vtkPVblockMesh::convertMeshBlocks" << endl;
63  }
64 
65  int blockI = 0;
66  const scalar scaleFactor = blkMesh.scaleFactor();
67 
68  for
69  (
70  int partId = range.start();
71  partId < range.end();
72  ++partId, ++blockI
73  )
74  {
75  if (!blockStatus_[partId])
76  {
77  continue;
78  }
79 
80  const blockDescriptor& blockDef = blkMesh[blockI];
81 
82  vtkUnstructuredGrid* vtkmesh = vtkUnstructuredGrid::New();
83 
84  // Convert OpenFOAM mesh vertices to VTK
85  vtkPoints *vtkpoints = vtkPoints::New();
86  vtkpoints->Allocate( blockDef.nPoints() );
87  const labelList& blockLabels = blockDef.blockShape();
88 
89  vtkmesh->Allocate(1);
90  vtkIdType nodeIds[8];
91 
92  forAll(blockLabels, ptI)
93  {
95  (
96  vtkpoints,
97  blockPoints[blockLabels[ptI]],
98  scaleFactor
99  );
100 
101  nodeIds[ptI] = ptI;
102  }
103 
104  vtkmesh->InsertNextCell
105  (
106  VTK_HEXAHEDRON,
107  8,
108  nodeIds
109  );
110 
111  vtkmesh->SetPoints(vtkpoints);
112  vtkpoints->Delete();
113 
114  AddToBlock
115  (
116  output, vtkmesh, range, datasetNo,
117  selection->GetArrayName(partId)
118  );
119 
120  vtkmesh->Delete();
121  datasetNo++;
122  }
123 
124 
125  // anything added?
126  if (datasetNo)
127  {
128  ++blockNo;
129  }
130 
131  if (debug)
132  {
133  Info<< "<end> Foam::vtkPVblockMesh::convertMeshBlocks" << endl;
134  }
135 }
136 
137 
138 void Foam::vtkPVblockMesh::convertMeshEdges
139 (
140  vtkMultiBlockDataSet* output,
141  int& blockNo
142 )
143 {
144  vtkDataArraySelection* selection = reader_->GetCurvedEdgesSelection();
145  arrayRange& range = arrayRangeEdges_;
146 
147  range.block(blockNo); // set output block
148  label datasetNo = 0; // restart at dataset 0
149 
150  const blockMesh& blkMesh = *meshPtr_;
151  const blockEdgeList& edges = blkMesh.edges();
152 
153  int edgeI = 0;
154  const scalar scaleFactor = blkMesh.scaleFactor();
155 
156  for
157  (
158  int partId = range.start();
159  partId < range.end();
160  ++partId, ++edgeI
161  )
162  {
163  if (!edgeStatus_[partId])
164  {
165  continue;
166  }
167 
168  // search each block
169  forAll(blkMesh, blockI)
170  {
171  const blockDescriptor& blockDef = blkMesh[blockI];
172 
173  edgeList blkEdges = blockDef.blockShape().edges();
174 
175  // List of edge point and weighting factors
176  pointField edgesPoints[12];
177  scalarList edgesWeights[12];
178  blockDef.edgesPointsWeights(edgesPoints, edgesWeights);
179 
180  // find the corresponding edge within the block
181  label foundEdgeI = -1;
182  forAll(blkEdges, blkEdgeI)
183  {
184  if (edges[edgeI].compare(blkEdges[blkEdgeI]))
185  {
186  foundEdgeI = blkEdgeI;
187  break;
188  }
189  }
190 
191  if (foundEdgeI != -1)
192  {
193  const List<point>& edgePoints = edgesPoints[foundEdgeI];
194 
195  vtkPolyData* vtkmesh = vtkPolyData::New();
196  vtkPoints* vtkpoints = vtkPoints::New();
197 
198  vtkpoints->Allocate( edgePoints.size() );
199  vtkmesh->Allocate(1);
200 
201  vtkIdType pointIds[edgePoints.size()];
202  forAll(edgePoints, ptI)
203  {
205  (
206  vtkpoints,
207  edgePoints[ptI],
208  scaleFactor
209  );
210  pointIds[ptI] = ptI;
211  }
212 
213  vtkmesh->InsertNextCell
214  (
215  VTK_POLY_LINE,
216  edgePoints.size(),
217  pointIds
218  );
219 
220  vtkmesh->SetPoints(vtkpoints);
221  vtkpoints->Delete();
222 
223  AddToBlock
224  (
225  output, vtkmesh, range, datasetNo,
226  selection->GetArrayName(partId)
227  );
228 
229  vtkmesh->Delete();
230  datasetNo++;
231 
232  break;
233  }
234  }
235  }
236 
237 
238  // anything added?
239  if (datasetNo)
240  {
241  ++blockNo;
242  }
243 
244  if (debug)
245  {
246  Info<< "<end> Foam::vtkPVblockMesh::convertMeshEdges" << endl;
247  }
248 
249 }
250 
251 
252 void Foam::vtkPVblockMesh::convertMeshCorners
253 (
254  vtkMultiBlockDataSet* output,
255  int& blockNo
256 )
257 {
258  arrayRange& range = arrayRangeCorners_;
259  range.block(blockNo); // set output block
260  label datasetNo = 0; // restart at dataset 0
261 
262  const pointField& blockPoints = meshPtr_->vertices();
263  const scalar& scaleFactor = meshPtr_->scaleFactor();
264 
265  if (debug)
266  {
267  Info<< "<beg> Foam::vtkPVblockMesh::convertMeshCorners" << endl;
268  }
269 
270  if (true) // or some flag or other condition
271  {
272  vtkPolyData* vtkmesh = vtkPolyData::New();
273  vtkPoints* vtkpoints = vtkPoints::New();
274  vtkCellArray* vtkcells = vtkCellArray::New();
275 
276  vtkpoints->Allocate( blockPoints.size() );
277  vtkcells->Allocate( blockPoints.size() );
278 
279  vtkIdType pointId = 0;
280  forAll(blockPoints, ptI)
281  {
283  (
284  vtkpoints,
285  blockPoints[ptI],
286  scaleFactor
287  );
288 
289  vtkcells->InsertNextCell(1, &pointId);
290  pointId++;
291  }
292 
293  vtkmesh->SetPoints(vtkpoints);
294  vtkpoints->Delete();
295 
296  vtkmesh->SetVerts(vtkcells);
297  vtkcells->Delete();
298 
299  AddToBlock
300  (
301  output, vtkmesh, range, datasetNo,
302  arrayRangeCorners_.name()
303  );
304  vtkmesh->Delete();
305 
306  datasetNo++;
307  }
308 
309  // anything added?
310  if (datasetNo)
311  {
312  ++blockNo;
313  }
314 
315  if (debug)
316  {
317  Info<< "<end> Foam::vtkPVblockMesh::convertMeshCorners" << endl;
318  }
319 }
320 
321 
322 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void vtkInsertNextOpenFOAMPoint(vtkPoints *points, const Foam::point &p)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
scalar scaleFactor() const
The scaling factor used to convert to meters.
const pointField & vertices() const
Reference to point field defining the blockMesh.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
List< edge > edgeList
Definition: edgeList.H:38
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
List< label > labelList
A List of labels.
Definition: labelList.H:56
PtrList< blockEdge > blockEdgeList
A PtrList of blockEdges.
Definition: blockEdgeList.H:45
messageStream Info