fvFieldReconstructorReconstructFields.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-2023 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 "fvFieldReconstructor.H"
27 #include "Time.H"
28 #include "PtrList.H"
29 #include "fvPatchFields.H"
30 #include "emptyFvPatch.H"
31 #include "emptyFvPatchField.H"
32 #include "emptyFvsPatchField.H"
33 #include "processorCyclicFvPatch.H"
35 #include "stringOps.H"
36 
37 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
38 
39 template<class Type>
40 void Foam::fvFieldReconstructor::rmapFaceToFace
41 (
42  Field<Type>& toField,
43  const Field<Type>& fromField,
44  const labelUList& addressing,
45  const bool isFlux
46 )
47 {
48  forAll(addressing, i)
49  {
50  toField[mag(addressing[i]) - 1] =
51  (isFlux && addressing[i] < 0 ? -1 : +1)*fromField[i];
52  }
53 }
54 
55 
56 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
57 
58 template<class Type>
61 (
62  const IOobject& fieldIoObject,
63  const PtrList<DimensionedField<Type, volMesh>>& procFields
64 ) const
65 {
66  // Create the internalField
67  Field<Type> internalField(completeMesh_.nCells());
68 
69  forAll(procMeshes_, proci)
70  {
71  const DimensionedField<Type, volMesh>& procField = procFields[proci];
72 
73  // Set the cell values in the reconstructed field
74  internalField.rmap
75  (
76  procField.field(),
77  cellProcAddressing_[proci]
78  );
79  }
80 
81  return tmp<DimensionedField<Type, volMesh>>
82  (
83  new DimensionedField<Type, volMesh>
84  (
85  fieldIoObject,
86  completeMesh_,
87  procFields[0].dimensions(),
88  internalField
89  )
90  );
91 }
92 
93 
94 template<class Type>
97 (
98  const IOobject& fieldIoObject
99 ) const
100 {
101  PtrList<DimensionedField<Type, volMesh>>
102  procFields(procMeshes_.size());
103 
104  forAll(procMeshes_, proci)
105  {
106  procFields.set
107  (
108  proci,
109  new DimensionedField<Type, volMesh>
110  (
111  IOobject
112  (
113  fieldIoObject.name(),
114  procMeshes_[proci].time().name(),
115  procMeshes_[proci],
118  false
119  ),
120  procMeshes_[proci]
121  )
122  );
123  }
124 
125  return reconstructFvVolumeInternalField
126  (
127  IOobject
128  (
129  fieldIoObject.name(),
130  completeMesh_.time().name(),
131  completeMesh_,
134  false
135  ),
136  procFields
137  );
138 }
139 
140 
141 template<class Type>
144 (
145  const IOobject& fieldIoObject,
146  const PtrList<VolField<Type>>& procFields
147 ) const
148 {
149  // Create the internalField
150  Field<Type> internalField(completeMesh_.nCells());
151 
152  // Create the patch fields
153  PtrList<fvPatchField<Type>> patchFields(completeMesh_.boundary().size());
154 
155  forAll(procFields, proci)
156  {
157  const VolField<Type>& procField =
158  procFields[proci];
159 
160  // Set the cell values in the reconstructed field
161  internalField.rmap
162  (
163  procField.primitiveField(),
164  cellProcAddressing_[proci]
165  );
166 
167  // Set the boundary patch values in the reconstructed field
168  forAll(procField.boundaryField(), procPatchi)
169  {
170  const fvPatch& procPatch =
171  procMeshes_[proci].boundary()[procPatchi];
172 
173  const label completePatchi = completePatchID(proci, procPatchi);
174 
175  if (completePatchi == procPatchi)
176  {
177  if (!patchFields(completePatchi))
178  {
179  patchFields.set
180  (
181  completePatchi,
183  (
184  procField.boundaryField()[procPatchi],
185  completeMesh_.boundary()[completePatchi],
187  setSizeFvPatchFieldMapper
188  (
189  completeMesh_.boundary()[completePatchi].size()
190  )
191  )
192  );
193  }
194 
195  patchFields[completePatchi].map
196  (
197  procField.boundaryField()[procPatchi],
198  reverseFvPatchFieldMapper
199  (
200  faceProcAddressingBf_[proci][procPatchi] - 1
201  )
202  );
203  }
204  else if (isA<processorCyclicFvPatch>(procPatch))
205  {
206  if (!patchFields(completePatchi))
207  {
208  patchFields.set
209  (
210  completePatchi,
212  (
213  completeMesh_.boundary()[completePatchi].type(),
214  completeMesh_.boundary()[completePatchi],
216  )
217  );
218  }
219 
220  if (patchFields[completePatchi].overridesConstraint())
221  {
222  OStringStream str;
223  str << "\nThe field \"" << procFields[0].name()
224  << "\" on cyclic patch \""
225  << patchFields[completePatchi].patch().name()
226  << "\" cannot be reconstructed as it is not a cyclic "
227  << "patch field. A \"patchType cyclic;\" setting has "
228  << "been used to override the cyclic patch type.\n\n"
229  << "Cyclic patches like this with non-cyclic boundary "
230  << "conditions should be confined to a single "
231  << "processor using decomposition constraints.";
233  << stringOps::breakIntoIndentedLines(str.str()).c_str()
234  << exit(FatalError);
235  }
236 
237  patchFields[completePatchi].map
238  (
239  procField.boundaryField()[procPatchi],
240  reverseFvPatchFieldMapper
241  (
242  faceProcAddressingBf_[proci][procPatchi] - 1
243  )
244  );
245  }
246  }
247  }
248 
249  // Construct and return the field
250  return tmp<VolField<Type>>
251  (
252  new VolField<Type>
253  (
254  fieldIoObject,
255  completeMesh_,
256  procFields[0].dimensions(),
257  internalField,
258  patchFields
259  )
260  );
261 }
262 
263 
264 template<class Type>
267 (
268  const IOobject& fieldIoObject
269 ) const
270 {
271  PtrList<VolField<Type>>
272  procFields(procMeshes_.size());
273 
274  forAll(procMeshes_, proci)
275  {
276  procFields.set
277  (
278  proci,
279  new VolField<Type>
280  (
281  IOobject
282  (
283  fieldIoObject.name(),
284  procMeshes_[proci].time().name(),
285  procMeshes_[proci],
288  false
289  ),
290  procMeshes_[proci]
291  )
292  );
293  }
294 
295  return reconstructFvVolumeField
296  (
297  IOobject
298  (
299  fieldIoObject.name(),
300  completeMesh_.time().name(),
301  completeMesh_,
304  false
305  ),
306  procFields
307  );
308 }
309 
310 
311 template<class Type>
314 (
315  const IOobject& fieldIoObject,
316  const PtrList<SurfaceField<Type>>& procFields
317 ) const
318 {
319  // Create the internalField
320  Field<Type> internalField(completeMesh_.nInternalFaces());
321 
322  // Create the patch fields
323  PtrList<fvsPatchField<Type>> patchFields(completeMesh_.boundary().size());
324 
325  forAll(procMeshes_, proci)
326  {
327  const SurfaceField<Type>& procField =
328  procFields[proci];
329 
330  // Set the internal face values in the reconstructed field
331  rmapFaceToFace
332  (
333  internalField,
334  procField.primitiveField(),
335  SubList<label>
336  (
337  faceProcAddressing_[proci],
338  procMeshes_[proci].nInternalFaces()
339  ),
340  isFlux(procFields[proci])
341  );
342 
343  // Set the boundary patch values in the reconstructed field
344  forAll(procField.boundaryField(), procPatchi)
345  {
346  const fvPatch& procPatch =
347  procMeshes_[proci].boundary()[procPatchi];
348 
349  const label completePatchi = completePatchID(proci, procPatchi);
350 
351  if (completePatchi == procPatchi)
352  {
353  if (!patchFields(completePatchi))
354  {
355  patchFields.set
356  (
357  completePatchi,
359  (
360  procField.boundaryField()[procPatchi],
361  completeMesh_.boundary()[completePatchi],
363  setSizeFvPatchFieldMapper
364  (
365  completeMesh_.boundary()[completePatchi].size()
366  )
367  )
368  );
369  }
370 
371  patchFields[completePatchi].map
372  (
373  procField.boundaryField()[procPatchi],
374  reverseFvPatchFieldMapper
375  (
376  faceProcAddressingBf_[proci][procPatchi] - 1
377  )
378  );
379  }
380  else if (isA<processorCyclicFvPatch>(procPatch))
381  {
382  if (!patchFields(completePatchi))
383  {
384  patchFields.set
385  (
386  completePatchi,
388  (
389  completeMesh_.boundary()[completePatchi].type(),
390  completeMesh_.boundary()[completePatchi],
392  )
393  );
394  }
395 
396  patchFields[completePatchi].map
397  (
398  procField.boundaryField()[procPatchi],
399  reverseFvPatchFieldMapper
400  (
401  faceProcAddressingBf_[proci][procPatchi] - 1
402  )
403  );
404  }
405  else if (isA<processorFvPatch>(procPatch))
406  {
407  rmapFaceToFace
408  (
409  internalField,
410  procField.boundaryField()[procPatchi],
411  faceProcAddressingBf_[proci][procPatchi],
412  isFlux(procFields[proci])
413  );
414  }
415  }
416  }
417 
418  // Construct and return the field
419  return tmp<SurfaceField<Type>>
420  (
421  new SurfaceField<Type>
422  (
423  fieldIoObject,
424  completeMesh_,
425  procFields[0].dimensions(),
426  internalField,
427  patchFields
428  )
429  );
430 }
431 
432 
433 template<class Type>
436 (
437  const IOobject& fieldIoObject
438 ) const
439 {
440  PtrList<SurfaceField<Type>>
441  procFields(procMeshes_.size());
442 
443  forAll(procMeshes_, proci)
444  {
445  procFields.set
446  (
447  proci,
448  new SurfaceField<Type>
449  (
450  IOobject
451  (
452  fieldIoObject.name(),
453  procMeshes_[proci].time().name(),
454  procMeshes_[proci],
457  false
458  ),
459  procMeshes_[proci]
460  )
461  );
462  }
463 
464  return reconstructFvSurfaceField
465  (
466  IOobject
467  (
468  fieldIoObject.name(),
469  completeMesh_.time().name(),
470  completeMesh_,
473  false
474  ),
475  procFields
476  );
477 }
478 
479 
480 template<class Type>
482 (
483  const IOobjectList& objects,
484  const HashSet<word>& selectedFields
485 )
486 {
487  const word& fieldClassName = DimensionedField<Type, volMesh>::typeName;
488 
489  IOobjectList fields = objects.lookupClass(fieldClassName);
490 
491  if (fields.size())
492  {
493  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
494 
495  forAllConstIter(IOobjectList, fields, fieldIter)
496  {
497  if
498  (
499  selectedFields.empty()
500  || selectedFields.found(fieldIter()->name())
501  )
502  {
503  Info<< " " << fieldIter()->name() << endl;
504 
505  reconstructFvVolumeInternalField<Type>(*fieldIter())().write();
506 
507  nReconstructed_++;
508  }
509  }
510  Info<< endl;
511  }
512 }
513 
514 
515 template<class Type>
517 (
518  const IOobjectList& objects,
519  const HashSet<word>& selectedFields
520 )
521 {
522  const word& fieldClassName =
524 
525  IOobjectList fields = objects.lookupClass(fieldClassName);
526 
527  if (fields.size())
528  {
529  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
530 
531  forAllConstIter(IOobjectList, fields, fieldIter)
532  {
533  if
534  (
535  selectedFields.empty()
536  || selectedFields.found(fieldIter()->name())
537  )
538  {
539  Info<< " " << fieldIter()->name() << endl;
540 
541  reconstructFvVolumeField<Type>(*fieldIter())().write();
542 
543  nReconstructed_++;
544  }
545  }
546  Info<< endl;
547  }
548 }
549 
550 
551 template<class Type>
553 (
554  const IOobjectList& objects,
555  const HashSet<word>& selectedFields
556 )
557 {
558  const word& fieldClassName =
560 
561  IOobjectList fields = objects.lookupClass(fieldClassName);
562 
563  if (fields.size())
564  {
565  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
566 
567  forAllConstIter(IOobjectList, fields, fieldIter)
568  {
569  if
570  (
571  selectedFields.empty()
572  || selectedFields.found(fieldIter()->name())
573  )
574  {
575  Info<< " " << fieldIter()->name() << endl;
576 
577  reconstructFvSurfaceField<Type>(*fieldIter())().write();
578 
579  nReconstructed_++;
580  }
581  }
582  Info<< endl;
583  }
584 }
585 
586 
587 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
static const DimensionedField< Type, GeoMesh > & null()
Return a null DimensionedField.
static const char *const typeName
Definition: Field.H:105
tmp< SurfaceField< Type > > reconstructFvSurfaceField(const IOobject &fieldIoObject, const PtrList< SurfaceField< Type >> &) const
Reconstruct surface field.
void reconstructFvVolumeInternalFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected volume internal fields.
void reconstructFvVolumeFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected volume fields.
tmp< DimensionedField< Type, volMesh > > reconstructFvVolumeInternalField(const IOobject &fieldIoObject, const PtrList< DimensionedField< Type, volMesh >> &procFields) const
Reconstruct volume internal field.
tmp< VolField< Type > > reconstructFvVolumeField(const IOobject &fieldIoObject, const PtrList< VolField< Type >> &) const
Reconstruct volume field.
void reconstructFvSurfaceFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected surface fields.
static tmp< fvPatchField< Type > > New(const word &, const fvPatch &, const DimensionedField< Type, volMesh > &)
Return a pointer to a new patchField created on freestore given.
static tmp< fvsPatchField< Type > > New(const word &, const fvPatch &, const DimensionedField< Type, surfaceMesh > &)
Return a pointer to a new patchField created on freestore given.
A class for managing temporary objects.
Definition: tmp.H:55
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
Info<< "Calculating turbulent flame speed field St\n"<< endl;volScalarField St(IOobject("St", runTime.name(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), flameWrinkling->Xi() *Su);multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:230
string breakIntoIndentedLines(const string &str, const string::size_type nLength=80, const string::size_type nIndent=0)
Break a string up into indented lines.
Definition: stringOps.C:950
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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:251
messageStream Info
bool isFlux(const DimensionedField< Type, surfaceMesh > &df)
Check if surfaceField is a flux.
Definition: surfaceFields.H:55
dimensioned< scalar > mag(const dimensioned< Type > &)
error FatalError
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
UList< label > labelUList
Definition: UList.H:65
objects