reconstructLagrangianFields.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 "IOField.H"
27 #include "CompactIOField.H"
28 #include "Time.H"
29 
30 // * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
31 
32 template<class Type>
34 (
35  const word& cloudName,
36  const polyMesh& mesh,
37  const PtrList<fvMesh>& meshes,
38  const word& fieldName
39 )
40 {
41  // Construct empty field on mesh
42  tmp<IOField<Type>> tfield
43  (
44  new IOField<Type>
45  (
46  IOobject
47  (
48  fieldName,
49  mesh.time().timeName(),
50  cloud::prefix/cloudName,
51  mesh,
52  IOobject::NO_READ,
53  IOobject::NO_WRITE
54  ),
55  Field<Type>(0)
56  )
57  );
58  Field<Type>& field = tfield.ref();
59 
60  forAll(meshes, i)
61  {
62  // Check object on local mesh
63  IOobject localIOobject
64  (
65  fieldName,
66  meshes[i].time().timeName(),
67  cloud::prefix/cloudName,
68  meshes[i],
69  IOobject::MUST_READ,
70  IOobject::NO_WRITE
71  );
72 
73  if (localIOobject.typeHeaderOk<IOField<Type>>(true))
74  {
75  IOField<Type> fieldi(localIOobject);
76 
77  label offset = field.size();
78  field.setSize(offset + fieldi.size());
79 
80  forAll(fieldi, j)
81  {
82  field[offset + j] = fieldi[j];
83  }
84  }
85  }
86 
87  return tfield;
88 }
89 
90 
91 template<class Type>
94 (
95  const word& cloudName,
96  const polyMesh& mesh,
97  const PtrList<fvMesh>& meshes,
98  const word& fieldName
99 )
100 {
101  // Construct empty field on mesh
102  tmp<CompactIOField<Field<Type>, Type >> tfield
103  (
104  new CompactIOField<Field<Type>, Type>
105  (
106  IOobject
107  (
108  fieldName,
109  mesh.time().timeName(),
110  cloud::prefix/cloudName,
111  mesh,
112  IOobject::NO_READ,
113  IOobject::NO_WRITE
114  ),
115  Field<Field<Type>>(0)
116  )
117  );
118  Field<Field<Type>>& field = tfield.ref();
119 
120  forAll(meshes, i)
121  {
122  // Check object on local mesh
123  IOobject localIOobject
124  (
125  fieldName,
126  meshes[i].time().timeName(),
127  cloud::prefix/cloudName,
128  meshes[i],
129  IOobject::MUST_READ,
130  IOobject::NO_WRITE
131  );
132 
133  if
134  (
135  localIOobject.typeHeaderOk<CompactIOField<Field<Type>, Type>>
136  (
137  false
138  )
139  || localIOobject.typeHeaderOk<IOField<Field<Type>>>(false)
140  )
141  {
142  CompactIOField<Field<Type>, Type> fieldi(localIOobject);
143 
144  label offset = field.size();
145  field.setSize(offset + fieldi.size());
146 
147  forAll(fieldi, j)
148  {
149  field[offset + j] = fieldi[j];
150  }
151  }
152  }
153 
154  return tfield;
155 }
156 
157 
158 
159 template<class Type>
161 (
162  const word& cloudName,
163  const polyMesh& mesh,
164  const PtrList<fvMesh>& meshes,
165  const IOobjectList& objects,
166  const HashSet<word>& selectedFields
167 )
168 {
169  const word fieldClassName(IOField<Type>::typeName);
170 
171  IOobjectList fields = objects.lookupClass(fieldClassName);
172 
173  if (fields.size())
174  {
175  Info<< " Reconstructing lagrangian "
176  << fieldClassName << "s\n" << endl;
177 
178  forAllConstIter(IOobjectList, fields, fieldIter)
179  {
180  if
181  (
182  selectedFields.empty()
183  || selectedFields.found(fieldIter()->name())
184  )
185  {
186  Info<< " " << fieldIter()->name() << endl;
187  reconstructLagrangianField<Type>
188  (
189  cloudName,
190  mesh,
191  meshes,
192  fieldIter()->name()
193  )().write();
194  }
195  }
196 
197  Info<< endl;
198  }
199 }
200 
201 
202 template<class Type>
204 (
205  const word& cloudName,
206  const polyMesh& mesh,
207  const PtrList<fvMesh>& meshes,
208  const IOobjectList& objects,
209  const HashSet<word>& selectedFields
210 )
211 {
212  {
213  const word fieldClassName(CompactIOField<Field<Type>, Type>::typeName);
214 
215  IOobjectList fields = objects.lookupClass(fieldClassName);
216 
217  if (fields.size())
218  {
219  Info<< " Reconstructing lagrangian "
220  << fieldClassName << "s\n" << endl;
221 
222  forAllConstIter(IOobjectList, fields, fieldIter)
223  {
224  if
225  (
226  selectedFields.empty()
227  || selectedFields.found(fieldIter()->name())
228  )
229  {
230  Info<< " " << fieldIter()->name() << endl;
231  reconstructLagrangianFieldField<Type>
232  (
233  cloudName,
234  mesh,
235  meshes,
236  fieldIter()->name()
237  )().write();
238  }
239  }
240 
241  Info<< endl;
242  }
243  }
244 
245  {
246  const word fieldClassName(IOField<Field<Type>>::typeName);
247 
248  IOobjectList fields = objects.lookupClass(fieldClassName);
249 
250  if (fields.size())
251  {
252  Info<< " Reconstructing lagrangian "
253  << fieldClassName << "s\n" << endl;
254 
255  forAllConstIter(IOobjectList, fields, fieldIter)
256  {
257  if
258  (
259  selectedFields.empty()
260  || selectedFields.found(fieldIter()->name())
261  )
262  {
263  Info<< " " << fieldIter()->name() << endl;
264  reconstructLagrangianFieldField<Type>
265  (
266  cloudName,
267  mesh,
268  meshes,
269  fieldIter()->name()
270  )().write();
271  }
272  }
273 
274  Info<< endl;
275  }
276  }
277 }
278 
279 
280 // ************************************************************************* //
A HashTable with keys but without contents.
Definition: HashSet.H:59
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual Ostream & write(const char)=0
Write character.
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
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
bool empty() const
Return true if the hash table is empty.
Definition: HashTableI.H:72
void reconstructLagrangianFields(const word &cloudName, const polyMesh &mesh, const PtrList< fvMesh > &meshes, const IOobjectList &objects, const HashSet< word > &selectedFields)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Info<< "Calculating turbulent flame speed field St\"<< endl;volScalarField St(IOobject("St", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), flameWrinkling->Xi() *Su);multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:228
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
tmp< CompactIOField< Field< Type >, Type > > reconstructLagrangianFieldField(const word &cloudName, const polyMesh &mesh, const PtrList< fvMesh > &meshes, const word &fieldName)
dynamicFvMesh & mesh
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:113
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
tmp< IOField< Type > > reconstructLagrangianField(const word &cloudName, const polyMesh &mesh, const PtrList< fvMesh > &meshes, const word &fieldName)
word timeName
Definition: getTimeIndex.H:3
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
IOobjectList lookupClass(const word &className) const
Return the list for all IOobjects of a given class.
Definition: IOobjectList.C:192
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
A Field of objects of type <T> with automated input and output using a compact storage. Behaves like IOField except when binary output in case it writes a CompactListList.
messageStream Info
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const word cloudName(propsDict.lookup("cloudName"))
A class for managing temporary objects.
Definition: PtrList.H:53
void reconstructLagrangianFieldFields(const word &cloudName, const polyMesh &mesh, const PtrList< fvMesh > &meshes, const IOobjectList &objects, const HashSet< word > &selectedFields)
A primitive field of type <T> with automated input and output.
Definition: IOField.H:50