vtkPV3FoamVolFields.H
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 InClass
25  vtkPV3Foam
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #ifndef vtkPV3FoamVolFields_H
30 #define vtkPV3FoamVolFields_H
31 
32 // OpenFOAM includes
33 #include "emptyFvPatchField.H"
34 #include "wallPolyPatch.H"
35 #include "faceSet.H"
36 #include "volPointInterpolation.H"
37 
38 #include "vtkPV3FoamFaceField.H"
39 #include "vtkPV3FoamPatchField.H"
40 
41 #include "vtkOpenFOAMTupleRemap.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 template<class Type>
46 void Foam::vtkPV3Foam::convertVolFields
47 (
48  const fvMesh& mesh,
49  const PtrList<PrimitivePatchInterpolation<primitivePatch>>& ppInterpList,
50  const IOobjectList& objects,
51  const bool interpFields,
52  vtkMultiBlockDataSet* output
53 )
54 {
55  const polyBoundaryMesh& patches = mesh.boundaryMesh();
56 
57  forAllConstIter(IOobjectList, objects, iter)
58  {
59  // restrict to GeometricField<Type, ...>
60  if
61  (
62  iter()->headerClassName()
64  )
65  {
66  continue;
67  }
68 
69  // Load field
70  GeometricField<Type, fvPatchField, volMesh> tf
71  (
72  *iter(),
73  mesh
74  );
75 
76  // Interpolated field (demand driven)
77  autoPtr<GeometricField<Type, pointPatchField, pointMesh>> ptfPtr;
78  if (interpFields)
79  {
80  if (debug)
81  {
82  Info<< "convertVolFieldBlock interpolating:" << tf.name()
83  << endl;
84  }
85 
86  ptfPtr.reset
87  (
88  volPointInterpolation::New(tf.mesh()).interpolate(tf).ptr()
89  );
90  }
91 
92 
93  // Convert activated internalMesh regions
94  convertVolFieldBlock
95  (
96  tf,
97  ptfPtr,
98  output,
99  arrayRangeVolume_,
100  regionPolyDecomp_
101  );
102 
103  // Convert activated cellZones
104  convertVolFieldBlock
105  (
106  tf,
107  ptfPtr,
108  output,
109  arrayRangeCellZones_,
110  zonePolyDecomp_
111  );
112 
113  // Convert activated cellSets
114  convertVolFieldBlock
115  (
116  tf,
117  ptfPtr,
118  output,
119  arrayRangeCellSets_,
120  csetPolyDecomp_
121  );
122 
123 
124  //
125  // Convert patches - if activated
126  //
127  for
128  (
129  int partId = arrayRangePatches_.start();
130  partId < arrayRangePatches_.end();
131  ++partId
132  )
133  {
134  const word patchName = getPartName(partId);
135  const label datasetNo = partDataset_[partId];
136  const label patchId = patches.findPatchID(patchName);
137 
138  if (!partStatus_[partId] || datasetNo < 0 || patchId < 0)
139  {
140  continue;
141  }
142 
143  const fvPatchField<Type>& ptf = tf.boundaryField()[patchId];
144 
145  if
146  (
147  isType<emptyFvPatchField<Type>>(ptf)
148  ||
149  (
150  reader_->GetExtrapolatePatches()
151  && !polyPatch::constraintType(patches[patchId].type())
152  )
153  )
154  {
155  fvPatch p(ptf.patch().patch(), tf.mesh().boundary());
156 
157  tmp<Field<Type>> tpptf
158  (
159  fvPatchField<Type>(p, tf).patchInternalField()
160  );
161 
162  convertPatchField
163  (
164  tf.name(),
165  tpptf(),
166  output,
167  arrayRangePatches_,
168  datasetNo
169  );
170 
171  if (interpFields)
172  {
173  convertPatchPointField
174  (
175  tf.name(),
176  ppInterpList[patchId].faceToPointInterpolate(tpptf)(),
177  output,
178  arrayRangePatches_,
179  datasetNo
180  );
181  }
182  }
183  else
184  {
185  convertPatchField
186  (
187  tf.name(),
188  ptf,
189  output,
190  arrayRangePatches_,
191  datasetNo
192  );
193 
194  if (interpFields)
195  {
196  convertPatchPointField
197  (
198  tf.name(),
199  ppInterpList[patchId].faceToPointInterpolate(ptf)(),
200  output,
201  arrayRangePatches_,
202  datasetNo
203  );
204  }
205  }
206  }
207 
208  //
209  // Convert face zones - if activated
210  //
211  for
212  (
213  int partId = arrayRangeFaceZones_.start();
214  partId < arrayRangeFaceZones_.end();
215  ++partId
216  )
217  {
218  const word zoneName = getPartName(partId);
219  const label datasetNo = partDataset_[partId];
220 
221  if (!partStatus_[partId] || datasetNo < 0)
222  {
223  continue;
224  }
225 
226  const faceZoneMesh& zMesh = mesh.faceZones();
227  const label zoneId = zMesh.findZoneID(zoneName);
228 
229  if (zoneId < 0)
230  {
231  continue;
232  }
233 
234  convertFaceField
235  (
236  tf,
237  output,
238  arrayRangeFaceZones_,
239  datasetNo,
240  mesh,
241  zMesh[zoneId]
242  );
243 
244  // TODO: points
245  }
246 
247  //
248  // Convert face sets - if activated
249  //
250  for
251  (
252  int partId = arrayRangeFaceSets_.start();
253  partId < arrayRangeFaceSets_.end();
254  ++partId
255  )
256  {
257  const word selectName = getPartName(partId);
258  const label datasetNo = partDataset_[partId];
259 
260  if (!partStatus_[partId] || datasetNo < 0)
261  {
262  continue;
263  }
264 
265  const faceSet fSet(mesh, selectName);
266 
267  convertFaceField
268  (
269  tf,
270  output,
271  arrayRangeFaceSets_,
272  datasetNo,
273  mesh,
274  fSet.toc()
275  );
276 
277  // TODO: points
278  }
279  }
280 }
281 
282 
283 template<class Type>
284 void Foam::vtkPV3Foam::convertVolFieldBlock
285 (
286  const GeometricField<Type, fvPatchField, volMesh>& tf,
287  autoPtr<GeometricField<Type, pointPatchField, pointMesh>>& ptfPtr,
288  vtkMultiBlockDataSet* output,
289  const arrayRange& range,
290  const List<polyDecomp>& decompLst
291 )
292 {
293  for (int partId = range.start(); partId < range.end(); ++partId)
294  {
295  const label datasetNo = partDataset_[partId];
296 
297  if (datasetNo >= 0 && partStatus_[partId])
298  {
299  convertVolField
300  (
301  tf,
302  output,
303  range,
304  datasetNo,
305  decompLst[datasetNo]
306  );
307 
308  if (ptfPtr.valid())
309  {
310  convertPointField
311  (
312  ptfPtr(),
313  tf,
314  output,
315  range,
316  datasetNo,
317  decompLst[datasetNo]
318  );
319  }
320  }
321  }
322 }
323 
324 
325 template<class Type>
326 void Foam::vtkPV3Foam::convertVolField
327 (
328  const GeometricField<Type, fvPatchField, volMesh>& tf,
329  vtkMultiBlockDataSet* output,
330  const arrayRange& range,
331  const label datasetNo,
332  const polyDecomp& decompInfo
333 )
334 {
335  const label nComp = pTraits<Type>::nComponents;
336  const labelList& superCells = decompInfo.superCells();
337 
338  vtkFloatArray* celldata = vtkFloatArray::New();
339  celldata->SetNumberOfTuples(superCells.size());
340  celldata->SetNumberOfComponents(nComp);
341  celldata->Allocate(nComp*superCells.size());
342  celldata->SetName(tf.name().c_str());
343 
344  if (debug)
345  {
346  Info<< "convert volField: "
347  << tf.name()
348  << " size = " << tf.size()
349  << " nComp=" << nComp
350  << " nTuples = " << superCells.size() << endl;
351  }
352 
353  float vec[nComp];
354  forAll(superCells, i)
355  {
356  const Type& t = tf[superCells[i]];
357  for (direction d=0; d<nComp; ++d)
358  {
359  vec[d] = component(t, d);
360  }
361  vtkOpenFOAMTupleRemap<Type>(vec);
362 
363  celldata->InsertTuple(i, vec);
364  }
365 
366  vtkUnstructuredGrid::SafeDownCast
367  (
368  GetDataSetFromBlock(output, range, datasetNo)
369  ) ->GetCellData()
370  ->AddArray(celldata);
371 
372  celldata->Delete();
373 }
374 
375 
376 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
377 
378 #endif
379 
380 // ************************************************************************* //
label patchId(-1)
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
#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
uint8_t direction
Definition: direction.H:46
static const char *const typeName
Definition: Field.H:94
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
const tensorField & tf
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
static const volPointInterpolation & New(const fvMesh &mesh)
static bool constraintType(const word &pt)
Return true if the given type is a constraint type.
Definition: polyPatch.C:246
bool isType(const Type &t)
Check the typeid.
Definition: typeInfo.H:126
List< label > labelList
A List of labels.
Definition: labelList.H:56
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
messageStream Info
volScalarField & p
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)