vtkPVFoamVolFields.H
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-2021 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 InClass
25  vtkPVFoam
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #ifndef vtkPVFoamVolFields_H
30 #define vtkPVFoamVolFields_H
31 
32 // OpenFOAM includes
33 #include "emptyFvPatchField.H"
34 #include "wallPolyPatch.H"
35 #include "faceSet.H"
36 #include "volPointInterpolation.H"
37 
38 #include "vtkPVFoamSurfaceField.H"
39 #include "vtkPVFoamPatchField.H"
40 #include "vtkOpenFOAMTupleRemap.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 template<class Type>
45 void Foam::vtkPVFoam::convertVolFields
46 (
47  const fvMesh& mesh,
48  const PtrList<PrimitivePatchInterpolation<primitivePatch>>& ppInterpList,
49  const IOobjectList& objects,
50  const bool interpFields,
51  vtkMultiBlockDataSet* output
52 )
53 {
54  const polyBoundaryMesh& patches = mesh.boundaryMesh();
55 
56  forAllConstIter(IOobjectList, objects, iter)
57  {
58  // Restrict to GeometricField<Type, ...>
59  if
60  (
61  iter()->headerClassName()
63  )
64  {
65  continue;
66  }
67 
68  // Load field
69  GeometricField<Type, fvPatchField, volMesh> tf
70  (
71  *iter(),
72  mesh
73  );
74 
75  // Interpolated field (demand driven)
76  autoPtr<GeometricField<Type, pointPatchField, pointMesh>> ptfPtr;
77  if (interpFields)
78  {
79  if (debug)
80  {
81  InfoInFunction<< "interpolating:" << tf.name() << endl;
82  }
83 
84  ptfPtr.reset
85  (
87  );
88  }
89 
90 
91  // Convert activated internalMesh regions
92  convertVolFieldBlock
93  (
94  tf,
95  ptfPtr,
96  output,
97  arrayRangeVolume_,
98  regionPolyDecomp_
99  );
100 
101  // Convert activated cellZones
102  convertVolFieldBlock
103  (
104  tf,
105  ptfPtr,
106  output,
107  arrayRangeCellZones_,
108  zonePolyDecomp_
109  );
110 
111  // Convert activated cellSets
112  convertVolFieldBlock
113  (
114  tf,
115  ptfPtr,
116  output,
117  arrayRangeCellSets_,
118  csetPolyDecomp_
119  );
120 
121 
122  // Convert patches - if activated
123  for
124  (
125  int partId = arrayRangePatches_.start();
126  partId < arrayRangePatches_.end();
127  ++partId
128  )
129  {
130  const word patchName = getPartName(partId);
131  const label datasetNo = partDataset_[partId];
132  const label patchId = patches.findPatchID(patchName);
133 
134  if (!partStatus_[partId] || datasetNo < 0 || patchId < 0)
135  {
136  continue;
137  }
138 
139  const fvPatchField<Type>& ptf = tf.boundaryField()[patchId];
140 
141  if
142  (
143  isType<emptyFvPatchField<Type>>(ptf)
144  ||
145  (
146  reader_->GetExtrapolatePatches()
147  && !polyPatch::constraintType(patches[patchId].type())
148  )
149  )
150  {
151  fvPatch p(ptf.patch().patch(), tf.mesh().boundary());
152 
153  tmp<Field<Type>> tpptf
154  (
155  fvPatchField<Type>(p, tf).patchInternalField()
156  );
157 
158  convertPatchField
159  (
160  tf.name(),
161  tpptf(),
162  output,
163  arrayRangePatches_,
164  datasetNo
165  );
166 
167  if (interpFields)
168  {
169  convertPatchPointField
170  (
171  tf.name(),
172  ppInterpList[patchId].faceToPointInterpolate(tpptf)(),
173  output,
174  arrayRangePatches_,
175  datasetNo
176  );
177  }
178  }
179  else
180  {
181  convertPatchField
182  (
183  tf.name(),
184  ptf,
185  output,
186  arrayRangePatches_,
187  datasetNo
188  );
189 
190  if (interpFields)
191  {
192  convertPatchPointField
193  (
194  tf.name(),
195  ppInterpList[patchId].faceToPointInterpolate(ptf)(),
196  output,
197  arrayRangePatches_,
198  datasetNo
199  );
200  }
201  }
202  }
203 
204  // Convert face zones - if activated
205  for
206  (
207  int partId = arrayRangeFaceZones_.start();
208  partId < arrayRangeFaceZones_.end();
209  ++partId
210  )
211  {
212  const word zoneName = getPartName(partId);
213  const label datasetNo = partDataset_[partId];
214 
215  if (!partStatus_[partId] || datasetNo < 0)
216  {
217  continue;
218  }
219 
220  const meshFaceZones& zMesh = mesh.faceZones();
221  const label zoneId = zMesh.findZoneID(zoneName);
222 
223  if (zoneId < 0)
224  {
225  continue;
226  }
227 
228  convertSurfaceField
229  (
230  tf,
231  output,
232  arrayRangeFaceZones_,
233  datasetNo,
234  mesh,
235  zMesh[zoneId]
236  );
237  }
238 
239  // Convert face sets - if activated
240  for
241  (
242  int partId = arrayRangeFaceSets_.start();
243  partId < arrayRangeFaceSets_.end();
244  ++partId
245  )
246  {
247  const word selectName = getPartName(partId);
248  const label datasetNo = partDataset_[partId];
249 
250  if (!partStatus_[partId] || datasetNo < 0)
251  {
252  continue;
253  }
254 
255  const faceSet fSet(mesh, selectName);
256 
257  convertSurfaceField
258  (
259  tf,
260  output,
261  arrayRangeFaceSets_,
262  datasetNo,
263  mesh,
264  fSet.toc()
265  );
266  }
267  }
268 }
269 
270 
271 template<class Type>
272 void Foam::vtkPVFoam::convertVolInternalFields
273 (
274  const fvMesh& mesh,
275  const IOobjectList& objects,
276  vtkMultiBlockDataSet* output
277 )
278 {
279  forAllConstIter(IOobjectList, objects, iter)
280  {
281  // Restrict to GeometricField<Type, ...>::Internal
282  if
283  (
284  iter()->headerClassName()
286  )
287  {
288  continue;
289  }
290 
291  // Load field
293  (
294  *iter(),
295  mesh
296  );
297 
298  // Convert activated internalMesh regions
299  convertVolInternalFieldBlock<Type>
300  (
301  tf,
302  output,
303  arrayRangeVolume_,
304  regionPolyDecomp_
305  );
306 
307  // Convert activated cellZones
308  convertVolInternalFieldBlock<Type>
309  (
310  tf,
311  output,
312  arrayRangeCellZones_,
313  zonePolyDecomp_
314  );
315 
316  // Convert activated cellSets
317  convertVolInternalFieldBlock<Type>
318  (
319  tf,
320  output,
321  arrayRangeCellSets_,
322  csetPolyDecomp_
323  );
324  }
325 }
326 
327 
328 template<class Type>
329 void Foam::vtkPVFoam::convertVolFieldBlock
330 (
331  const GeometricField<Type, fvPatchField, volMesh>& tf,
332  autoPtr<GeometricField<Type, pointPatchField, pointMesh>>& ptfPtr,
333  vtkMultiBlockDataSet* output,
334  const arrayRange& range,
335  const List<polyDecomp>& decompLst
336 )
337 {
338  for (int partId = range.start(); partId < range.end(); ++partId)
339  {
340  const label datasetNo = partDataset_[partId];
341 
342  if (datasetNo >= 0 && partStatus_[partId])
343  {
344  convertVolInternalField<Type>
345  (
346  tf,
347  output,
348  range,
349  datasetNo,
350  decompLst[datasetNo]
351  );
352 
353  if (ptfPtr.valid())
354  {
355  convertPointField
356  (
357  ptfPtr(),
358  tf,
359  output,
360  range,
361  datasetNo,
362  decompLst[datasetNo]
363  );
364  }
365  }
366  }
367 }
368 
369 
370 template<class Type>
371 void Foam::vtkPVFoam::convertVolInternalFieldBlock
372 (
374  vtkMultiBlockDataSet* output,
375  const arrayRange& range,
376  const List<polyDecomp>& decompLst
377 )
378 {
379  for (int partId = range.start(); partId < range.end(); ++partId)
380  {
381  const label datasetNo = partDataset_[partId];
382 
383  if (datasetNo >= 0 && partStatus_[partId])
384  {
385  convertVolInternalField<Type>
386  (
387  tf,
388  output,
389  range,
390  datasetNo,
391  decompLst[datasetNo]
392  );
393  }
394  }
395 }
396 
397 
398 template<class Type>
399 void Foam::vtkPVFoam::convertVolInternalField
400 (
402  vtkMultiBlockDataSet* output,
403  const arrayRange& range,
404  const label datasetNo,
405  const polyDecomp& decompInfo
406 )
407 {
408  const label nComp = pTraits<Type>::nComponents;
409  const labelList& superCells = decompInfo.superCells();
410 
411  vtkFloatArray* celldata = vtkFloatArray::New();
412  celldata->SetNumberOfTuples(superCells.size());
413  celldata->SetNumberOfComponents(nComp);
414  celldata->Allocate(nComp*superCells.size());
415  celldata->SetName(tf.name().c_str());
416 
417  if (debug)
418  {
420  << "converting volField::Internal: " << tf.name()
421  << " size = " << tf.size()
422  << " nComp=" << nComp
423  << " nTuples = " << superCells.size() << endl;
424  }
425 
426  float vec[nComp];
427  forAll(superCells, i)
428  {
429  const Type& t = tf[superCells[i]];
430  for (direction d=0; d<nComp; ++d)
431  {
432  vec[d] = component(t, d);
433  }
434  vtkOpenFOAMTupleRemap<Type>(vec);
435 
436  celldata->InsertTuple(i, vec);
437  }
438 
439  vtkUnstructuredGrid::SafeDownCast
440  (
441  GetDataSetFromBlock(output, range, datasetNo)
442  ) ->GetCellData()
443  ->AddArray(celldata);
444 
445  celldata->Delete();
446 }
447 
448 
449 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
450 
451 #endif
452 
453 // ************************************************************************* //
label patchId(-1)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const word & name() const
Return name.
Definition: IOobject.H:315
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
static const char *const typeName
Definition: Field.H:105
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
uint8_t direction
Definition: direction.H:45
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Traits class for primitives.
Definition: pTraits.H:50
DimensionedField< Type, GeoMesh > Internal
Type of the internal field from which this GeometricField is derived.
const tensorField & tf
scalar range
static bool constraintType(const word &pt)
Return true if the given type is a constraint type.
Definition: polyPatch.C:265
bool isType(const Type &t)
Check the typeid.
Definition: typeInfo.H:126
const Mesh & mesh() const
Return mesh.
MeshZones< faceZone, polyMesh > meshFaceZones
A MeshZones with the type faceZone.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
volScalarField & p
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
#define InfoInFunction
Report an information message using Foam::Info.