writeFunsTemplates.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 #include "writeFuns.H"
27 #include "interpolatePointToCell.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 // Store List in dest
32 template<class Type>
34 (
35  const List<Type>& source,
36  DynamicList<floatScalar>& dest
37 )
38 {
39  forAll(source, i)
40  {
41  insert(source[i], dest);
42  }
43 }
44 
45 
47 //template<class Type>
48 //void Foam::writeFuns::insert
49 //(
50 // const labelList& map,
51 // const List<Type>& source,
52 // DynamicList<floatScalar>& dest
53 //)
54 //{
55 // forAll(map, i)
56 // {
57 // insert(source[map[i]], dest);
58 // }
59 //}
60 
61 
62 template<class Type>
64 (
65  std::ostream& os,
66  const bool binary,
67  const GeometricField<Type, fvPatchField, volMesh>& vvf,
68  const vtkMesh& vMesh
69 )
70 {
71  const fvMesh& mesh = vMesh.mesh();
72 
73  const labelList& superCells = vMesh.topo().superCells();
74 
75  label nValues = mesh.nCells() + superCells.size();
76 
77  os << vvf.name() << ' ' << pTraits<Type>::nComponents << ' '
78  << nValues << " float" << std::endl;
79 
80  DynamicList<floatScalar> fField(pTraits<Type>::nComponents*nValues);
81 
82  insert(vvf.primitiveField(), fField);
83 
84  forAll(superCells, superCelli)
85  {
86  label origCelli = superCells[superCelli];
87 
88  insert(vvf[origCelli], fField);
89  }
90  write(os, binary, fField);
91 }
92 
93 
94 template<class Type>
96 (
97  std::ostream& os,
98  const bool binary,
99  const GeometricField<Type, pointPatchField, pointMesh>& pvf,
100  const vtkMesh& vMesh
101 )
102 {
103  const fvMesh& mesh = vMesh.mesh();
104  const vtkTopo& topo = vMesh.topo();
105 
106  const labelList& addPointCellLabels = topo.addPointCellLabels();
107  const label nTotPoints = mesh.nPoints() + addPointCellLabels.size();
108 
109  os << pvf.name() << ' ' << pTraits<Type>::nComponents << ' '
110  << nTotPoints << " float" << std::endl;
111 
112  DynamicList<floatScalar> fField(pTraits<Type>::nComponents*nTotPoints);
113 
114  insert(pvf, fField);
115 
116  forAll(addPointCellLabels, api)
117  {
118  label origCelli = addPointCellLabels[api];
119 
120  insert(interpolatePointToCell(pvf, origCelli), fField);
121  }
122  write(os, binary, fField);
123 }
124 
125 
126 template<class Type>
128 (
129  std::ostream& os,
130  const bool binary,
131  const GeometricField<Type, fvPatchField, volMesh>& vvf,
132  const GeometricField<Type, pointPatchField, pointMesh>& pvf,
133  const vtkMesh& vMesh
134 )
135 {
136  const fvMesh& mesh = vMesh.mesh();
137  const vtkTopo& topo = vMesh.topo();
138 
139  const labelList& addPointCellLabels = topo.addPointCellLabels();
140  const label nTotPoints = mesh.nPoints() + addPointCellLabels.size();
141 
142  os << vvf.name() << ' ' << pTraits<Type>::nComponents << ' '
143  << nTotPoints << " float" << std::endl;
144 
145  DynamicList<floatScalar> fField(pTraits<Type>::nComponents*nTotPoints);
146 
147  insert(pvf, fField);
148 
149  forAll(addPointCellLabels, api)
150  {
151  label origCelli = addPointCellLabels[api];
152 
153  insert(vvf[origCelli], fField);
154  }
155  write(os, binary, fField);
156 }
157 
158 
159 template<class Type, template<class> class PatchField, class GeoMesh>
161 (
162  std::ostream& os,
163  const bool binary,
164  const PtrList<GeometricField<Type, PatchField, GeoMesh>>& flds,
165  const vtkMesh& vMesh
166 )
167 {
168  forAll(flds, i)
169  {
170  write(os, binary, flds[i], vMesh);
171  }
172 }
173 
174 
175 template<class Type>
177 (
178  std::ostream& os,
179  const bool binary,
180  const volPointInterpolation& pInterp,
181  const PtrList<GeometricField<Type, fvPatchField, volMesh>>& flds,
182  const vtkMesh& vMesh
183 )
184 {
185  forAll(flds, i)
186  {
187  write(os, binary, flds[i], pInterp.interpolate(flds[i])(), vMesh);
188  }
189 }
190 
191 
192 // ************************************************************************* //
#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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
static void insert(const point &, DynamicList< floatScalar > &dest)
Append point to given DynamicList.
List< label > labelList
A List of labels.
Definition: labelList.H:56
static void write(std::ostream &, const bool, DynamicList< floatScalar > &)
Write floats ascii or binary.
Type interpolatePointToCell(const GeometricField< Type, pointPatchField, pointMesh > &ptf, const label celli)
Interpolates (averages) the vertex values to the cell center.