vtkPVFoamFields.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 // VTK includes
27 #include "vtkDataArraySelection.h"
28 #include "vtkPolyData.h"
29 #include "vtkUnstructuredGrid.h"
30 
31 // OpenFOAM includes
32 #include "IOobjectList.H"
33 #include "vtkPVFoam.H"
34 #include "vtkPVFoamReader.h"
35 #include "vtkPVFoamVolFields.H"
36 #include "vtkPVFoamPointFields.H"
38 
39 
40 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
41 
42 Foam::IOobjectList Foam::vtkPVFoam::getObjects
43 (
44  const wordHashSet& selected,
45  const fvMesh& mesh,
46  const fileName& instance,
47  const fileName& local
48 )
49 {
50  // If nothing is selected then return an empty list
51  if (selected.empty())
52  {
53  return IOobjectList(0);
54  }
55 
56  // Create the list of objects at the instance
57  IOobjectList objects(mesh, instance, local);
58 
59  // Add any objects from constant that are not already present
60  IOobjectList objectsConstant(mesh, dbPtr_().constant(), local);
61  forAllIter(IOobjectList, objectsConstant, iter)
62  {
63  if (!objects.found(iter.key()))
64  {
65  objects.add
66  (
67  *objectsConstant.HashPtrTable<IOobject>::remove(iter)
68  );
69  }
70  }
71 
72  // Remove everything that is not selected
73  forAllIter(IOobjectList, objects, iter)
74  {
75  if (!selected.found(iter()->name()))
76  {
77  objects.erase(iter);
78  }
79  }
80 
81  return objects;
82 }
83 
84 
85 void Foam::vtkPVFoam::convertVolFields
86 (
87  vtkMultiBlockDataSet* output
88 )
89 {
90  const fvMesh& mesh = *meshPtr_;
91 
92  wordHashSet selectedFields = getSelected
93  (
94  reader_->GetVolFieldSelection()
95  );
96 
97  if (selectedFields.empty())
98  {
99  return;
100  }
101 
102  // Get objects (fields) for this time - only keep selected fields
103  // the region name is already in the mesh db
104  IOobjectList objects
105  (
106  getObjects
107  (
108  selectedFields,
109  mesh,
110  dbPtr_().timeName()
111  )
112  );
113 
114  if (objects.empty())
115  {
116  return;
117  }
118 
119  if (debug)
120  {
121  Info<< "<beg> Foam::vtkPVFoam::convertVolFields" << nl
122  << "converting OpenFOAM volume fields" << endl;
123  forAllConstIter(IOobjectList, objects, iter)
124  {
125  Info<< " " << iter()->name()
126  << " == " << iter()->objectPath() << nl;
127  }
128  printMemory();
129  }
130 
131 
132  PtrList<PrimitivePatchInterpolation<primitivePatch>>
133  ppInterpList(mesh.boundaryMesh().size());
134 
135  forAll(ppInterpList, i)
136  {
137  ppInterpList.set
138  (
139  i,
140  new PrimitivePatchInterpolation<primitivePatch>
141  (
142  mesh.boundaryMesh()[i]
143  )
144  );
145  }
146 
147 
148  bool interpFields = reader_->GetInterpolateVolFields();
149 
150  convertVolFields<scalar>
151  (
152  mesh, ppInterpList, objects, interpFields, output
153  );
154  convertVolFields<vector>
155  (
156  mesh, ppInterpList, objects, interpFields, output
157  );
158  convertVolFields<sphericalTensor>
159  (
160  mesh, ppInterpList, objects, interpFields, output
161  );
162  convertVolFields<symmTensor>
163  (
164  mesh, ppInterpList, objects, interpFields, output
165  );
166  convertVolFields<tensor>
167  (
168  mesh, ppInterpList, objects, interpFields, output
169  );
170 
171  if (debug)
172  {
173  Info<< "<end> Foam::vtkPVFoam::convertVolFields" << endl;
174  printMemory();
175  }
176 }
177 
178 
179 void Foam::vtkPVFoam::convertPointFields
180 (
181  vtkMultiBlockDataSet* output
182 )
183 {
184  const fvMesh& mesh = *meshPtr_;
185 
186  wordHashSet selectedFields = getSelected
187  (
188  reader_->GetPointFieldSelection()
189  );
190 
191  if (selectedFields.empty())
192  {
193  if (debug)
194  {
195  Info<< "no point fields selected" << endl;
196  }
197  return;
198  }
199 
200  // Get objects (fields) for this time - only keep selected fields
201  // the region name is already in the mesh db
202  IOobjectList objects
203  (
204  getObjects
205  (
206  selectedFields,
207  mesh,
208  dbPtr_().timeName()
209  )
210  );
211 
212  if (objects.empty())
213  {
214  return;
215  }
216 
217  if (debug)
218  {
219  Info<< "<beg> Foam::vtkPVFoam::convertPointFields" << nl
220  << "converting OpenFOAM volume fields -> point fields" << endl;
221  forAllConstIter(IOobjectList, objects, iter)
222  {
223  Info<< " " << iter()->name()
224  << " == " << iter()->objectPath() << nl;
225  }
226  printMemory();
227  }
228 
229  // Construct interpolation on the raw mesh
230  const pointMesh& pMesh = pointMesh::New(mesh);
231 
232 
233  convertPointFields<scalar>
234  (
235  mesh, pMesh, objects, output
236  );
237  convertPointFields<vector>
238  (
239  mesh, pMesh, objects, output
240  );
241  convertPointFields<sphericalTensor>
242  (
243  mesh, pMesh, objects, output
244  );
245  convertPointFields<symmTensor>
246  (
247  mesh, pMesh, objects, output
248  );
249  convertPointFields<tensor>
250  (
251  mesh, pMesh, objects, output
252  );
253 
254  if (debug)
255  {
256  Info<< "<end> Foam::vtkPVFoam::convertPointFields" << endl;
257  printMemory();
258  }
259 }
260 
261 
262 void Foam::vtkPVFoam::convertLagrangianFields
263 (
264  vtkMultiBlockDataSet* output
265 )
266 {
267  arrayRange& range = arrayRangeLagrangian_;
268  const fvMesh& mesh = *meshPtr_;
269 
270  wordHashSet selectedFields = getSelected
271  (
272  reader_->GetLagrangianFieldSelection()
273  );
274 
275  if (selectedFields.empty())
276  {
277  return;
278  }
279 
280  if (debug)
281  {
282  Info<< "<beg> Foam::vtkPVFoam::convertLagrangianFields" << endl;
283  printMemory();
284  }
285 
286  for (int partId = range.start(); partId < range.end(); ++partId)
287  {
288  const word cloudName = getPartName(partId);
289  const label datasetNo = partDataset_[partId];
290 
291  if (!partStatus_[partId] || datasetNo < 0)
292  {
293  continue;
294  }
295 
296 
297  // Get the Lagrangian fields for this time and this cloud
298  // but only keep selected fields
299  // the region name is already in the mesh db
300  IOobjectList objects
301  (
302  getObjects
303  (
304  selectedFields,
305  mesh,
306  dbPtr_().timeName(),
307  cloud::prefix/cloudName
308  )
309  );
310 
311  if (objects.empty())
312  {
313  continue;
314  }
315 
316  if (debug)
317  {
318  Info<< "converting OpenFOAM lagrangian fields" << nl;
319  forAllConstIter(IOobjectList, objects, iter)
320  {
321  Info<< " " << iter()->name()
322  << " == " << iter()->objectPath() << nl;
323  }
324  }
325 
326  convertLagrangianFields<label>
327  (
328  objects, output, datasetNo
329  );
330  convertLagrangianFields<scalar>
331  (
332  objects, output, datasetNo
333  );
334  convertLagrangianFields<vector>
335  (
336  objects, output, datasetNo
337  );
338  convertLagrangianFields<sphericalTensor>
339  (
340  objects, output, datasetNo
341  );
342  convertLagrangianFields<symmTensor>
343  (
344  objects, output, datasetNo
345  );
346  convertLagrangianFields<tensor>
347  (
348  objects, output, datasetNo
349  );
350  }
351 
352  if (debug)
353  {
354  Info<< "<end> Foam::vtkPVFoam::convertLagrangianFields" << endl;
355  printMemory();
356  }
357 }
358 
359 
360 // ************************************************************************* //
#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
static void printMemory()
Simple memory used debugging information.
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:459
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
static const pointMesh & New(const polyMesh &mesh)
Definition: MeshObject.C:44
dynamicFvMesh & mesh
word timeName
Definition: getTimeIndex.H:3
HashSet wordHashSet
A HashSet with word keys.
Definition: HashSet.H:205
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
static const char nl
Definition: Ostream.H:265
objects
static const word prefix
The prefix to local: lagrangian.
Definition: cloud.H:62
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
messageStream Info