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