1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
28 \*---------------------------------------------------------------------------*/
30 #ifndef vtkPV3FoamPointFields_H
31 #define vtkPV3FoamPointFields_H
34 #include "interpolatePointToCell.H"
36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
39 void Foam::vtkPV3Foam::convertPointFields
42 const pointMesh& pMesh,
43 const IOobjectList& objects,
44 vtkMultiBlockDataSet* output
47 const polyBoundaryMesh& patches = mesh.boundaryMesh();
49 forAllConstIter(IOobjectList, objects, iter)
51 const word& fieldName = iter()->name();
52 // restrict to this GeometricField<Type, ...>
55 iter()->headerClassName()
56 != GeometricField<Type, pointPatchField, pointMesh>::typeName
64 Info<< "Foam::vtkPV3Foam::convertPointFields : "
68 GeometricField<Type, pointPatchField, pointMesh> ptf
75 // Convert activated internalMesh regions
76 convertPointFieldBlock
84 // Convert activated cellZones
85 convertPointFieldBlock
93 // Convert activated cellSets
94 convertPointFieldBlock
104 // Convert patches - if activated
108 int partId = partInfoPatches_.start();
109 partId < partInfoPatches_.end();
113 const word patchName = getPartName(partId);
114 const label datasetNo = partDataset_[partId];
115 const label patchId = patches.findPatchID(patchName);
117 if (!partStatus_[partId] || datasetNo < 0 || patchId < 0)
122 convertPatchPointField
125 ptf.boundaryField()[patchId].patchInternalField()(),
136 void Foam::vtkPV3Foam::convertPointFieldBlock
138 const GeometricField<Type, pointPatchField, pointMesh>& ptf,
139 vtkMultiBlockDataSet* output,
140 const partInfo& selector,
141 const List<polyDecomp>& decompLst
144 for (int partId = selector.start(); partId < selector.end(); ++partId)
146 const label datasetNo = partDataset_[partId];
148 if (datasetNo >= 0 && partStatus_[partId])
153 GeometricField<Type, fvPatchField, volMesh>::null(),
165 void Foam::vtkPV3Foam::convertPointField
167 const GeometricField<Type, pointPatchField, pointMesh>& ptf,
168 const GeometricField<Type, fvPatchField, volMesh>& tf,
169 vtkMultiBlockDataSet* output,
170 const partInfo& selector,
171 const label datasetNo,
172 const polyDecomp& decomp
175 const label nComp = pTraits<Type>::nComponents;
176 const labelList& addPointCellLabels = decomp.addPointCellLabels();
177 const labelList& pointMap = decomp.pointMap();
179 // use a pointMap or address directly into mesh
183 nPoints = pointMap.size();
187 nPoints = ptf.size();
190 vtkFloatArray *pointData = vtkFloatArray::New();
191 pointData->SetNumberOfTuples(nPoints + addPointCellLabels.size());
192 pointData->SetNumberOfComponents(nComp);
193 pointData->Allocate(nComp*(nPoints + addPointCellLabels.size()));
194 pointData->SetName(ptf.name().c_str());
198 Info<< "convert convertPointField: "
200 << " size = " << nPoints
201 << " nComp=" << nComp
202 << " nTuples = " << (nPoints + addPointCellLabels.size())
212 const Type& t = ptf[pointMap[i]];
213 for (direction d=0; d<nComp; d++)
215 vec[d] = component(t, d);
217 pointData->InsertTuple(i, vec);
224 const Type& t = ptf[i];
225 for (direction d=0; d<nComp; d++)
227 vec[d] = component(t, d);
229 pointData->InsertTuple(i, vec);
233 // continue insertion from here
236 if (&tf != &GeometricField<Type, fvPatchField, volMesh>::null())
238 forAll(addPointCellLabels, apI)
240 const Type& t = tf[addPointCellLabels[apI]];
241 for (direction d=0; d<nComp; d++)
243 vec[d] = component(t, d);
245 pointData->InsertTuple(i++, vec);
250 forAll(addPointCellLabels, apI)
252 Type t = interpolatePointToCell(ptf, addPointCellLabels[apI]);
253 for (direction d=0; d<nComp; d++)
255 vec[d] = component(t, d);
257 pointData->InsertTuple(i++, vec);
261 vtkUnstructuredGrid::SafeDownCast
263 GetDataSetFromBlock(output, selector, datasetNo)
265 ->AddArray(pointData);
271 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
275 // ************************************************************************* //