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-2019 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 
34 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
35 
36 template<class Type>
39 (
40  const IOobject& fieldIoObject,
41  const PtrList<DimensionedField<Type, volMesh>>& procFields
42 ) const
43 {
44  // Create the internalField
45  Field<Type> internalField(mesh_.nCells());
46 
47  forAll(procMeshes_, proci)
48  {
49  const DimensionedField<Type, volMesh>& procField = procFields[proci];
50 
51  // Set the cell values in the reconstructed field
52  internalField.rmap
53  (
54  procField.field(),
55  cellProcAddressing_[proci]
56  );
57  }
58 
60  (
62  (
63  fieldIoObject,
64  mesh_,
65  procFields[0].dimensions(),
66  internalField
67  )
68  );
69 }
70 
71 
72 template<class Type>
75 (
76  const IOobject& fieldIoObject
77 ) const
78 {
79  // Read the field for all the processors
81  (
82  procMeshes_.size()
83  );
84 
85  forAll(procMeshes_, proci)
86  {
87  procFields.set
88  (
89  proci,
91  (
92  IOobject
93  (
94  fieldIoObject.name(),
95  procMeshes_[proci].time().timeName(),
96  procMeshes_[proci],
97  IOobject::MUST_READ,
98  IOobject::NO_WRITE,
99  false
100  ),
101  procMeshes_[proci]
102  )
103  );
104  }
105 
106 
107  return reconstructFvVolumeInternalField
108  (
109  IOobject
110  (
111  fieldIoObject.name(),
112  mesh_.time().timeName(),
113  mesh_,
114  IOobject::NO_READ,
115  IOobject::NO_WRITE,
116  false
117  ),
118  procFields
119  );
120 }
121 
122 
123 template<class Type>
126 (
127  const IOobject& fieldIoObject,
129 ) const
130 {
131  // Create the internalField
132  Field<Type> internalField(mesh_.nCells());
133 
134  // Create the patch fields
135  PtrList<fvPatchField<Type>> patchFields(mesh_.boundary().size());
136 
137  forAll(procFields, proci)
138  {
140  procFields[proci];
141 
142  // Set the cell values in the reconstructed field
143  internalField.rmap
144  (
145  procField.primitiveField(),
146  cellProcAddressing_[proci]
147  );
148 
149  // Set the boundary patch values in the reconstructed field
150  forAll(boundaryProcAddressing_[proci], patchi)
151  {
152  // Get patch index of the original patch
153  const label curBPatch = boundaryProcAddressing_[proci][patchi];
154 
155  // Get addressing slice for this patch
156  const labelList::subList cp =
157  procField.mesh().boundary()[patchi].patchSlice
158  (
159  faceProcAddressing_[proci]
160  );
161 
162  // check if the boundary patch is not a processor patch
163  if (curBPatch >= 0)
164  {
165  // Regular patch. Fast looping
166 
167  if (!patchFields(curBPatch))
168  {
169  patchFields.set
170  (
171  curBPatch,
173  (
174  procField.boundaryField()[patchi],
175  mesh_.boundary()[curBPatch],
178  (
179  mesh_.boundary()[curBPatch].size()
180  )
181  )
182  );
183  }
184 
185  const label curPatchStart =
186  mesh_.boundaryMesh()[curBPatch].start();
187 
188  labelList reverseAddressing(cp.size());
189 
190  forAll(cp, facei)
191  {
192  // Check
193  if (cp[facei] <= 0)
194  {
196  << "Processor " << proci
197  << " patch "
198  << procField.mesh().boundary()[patchi].name()
199  << " face " << facei
200  << " originates from reversed face since "
201  << cp[facei]
202  << exit(FatalError);
203  }
204 
205  // Subtract one to take into account offsets for
206  // face direction.
207  reverseAddressing[facei] = cp[facei] - 1 - curPatchStart;
208  }
209 
210 
211  patchFields[curBPatch].rmap
212  (
213  procField.boundaryField()[patchi],
214  reverseAddressing
215  );
216  }
217  else
218  {
219  const Field<Type>& curProcPatch =
220  procField.boundaryField()[patchi];
221 
222  // In processor patches, there's a mix of internal faces (some
223  // of them turned) and possible cyclics. Slow loop
224  forAll(cp, facei)
225  {
226  // Subtract one to take into account offsets for
227  // face direction.
228  label curF = cp[facei] - 1;
229 
230  // Is the face on the boundary?
231  if (curF >= mesh_.nInternalFaces())
232  {
233  label curBPatch = mesh_.boundaryMesh().whichPatch(curF);
234 
235  if (!patchFields(curBPatch))
236  {
237  patchFields.set
238  (
239  curBPatch,
241  (
242  mesh_.boundary()[curBPatch].type(),
243  mesh_.boundary()[curBPatch],
245  )
246  );
247  }
248 
249  // add the face
250  label curPatchFace =
251  mesh_.boundaryMesh()
252  [curBPatch].whichFace(curF);
253 
254  patchFields[curBPatch][curPatchFace] =
255  curProcPatch[facei];
256  }
257  }
258  }
259  }
260  }
261 
262  forAll(mesh_.boundary(), patchi)
263  {
264  // add empty patches
265  if
266  (
267  isType<emptyFvPatch>(mesh_.boundary()[patchi])
268  && !patchFields(patchi)
269  )
270  {
271  patchFields.set
272  (
273  patchi,
275  (
277  mesh_.boundary()[patchi],
279  )
280  );
281  }
282  }
283 
284 
285  // Now construct and write the field
286  // setting the internalField and patchFields
288  (
290  (
291  fieldIoObject,
292  mesh_,
293  procFields[0].dimensions(),
294  internalField,
295  patchFields
296  )
297  );
298 }
299 
300 
301 template<class Type>
304 (
305  const IOobject& fieldIoObject
306 ) const
307 {
308  // Read the field for all the processors
310  (
311  procMeshes_.size()
312  );
313 
314  forAll(procMeshes_, proci)
315  {
316  procFields.set
317  (
318  proci,
320  (
321  IOobject
322  (
323  fieldIoObject.name(),
324  procMeshes_[proci].time().timeName(),
325  procMeshes_[proci],
326  IOobject::MUST_READ,
327  IOobject::NO_WRITE,
328  false
329  ),
330  procMeshes_[proci]
331  )
332  );
333  }
334 
335  return reconstructFvVolumeField
336  (
337  IOobject
338  (
339  fieldIoObject.name(),
340  mesh_.time().timeName(),
341  mesh_,
342  IOobject::NO_READ,
343  IOobject::NO_WRITE,
344  false
345  ),
346  procFields
347  );
348 }
349 
350 
351 template<class Type>
354 (
355  const IOobject& fieldIoObject,
357 ) const
358 {
359  // Create the internalField
360  Field<Type> internalField(mesh_.nInternalFaces());
361 
362  // Create the patch fields
363  PtrList<fvsPatchField<Type>> patchFields(mesh_.boundary().size());
364 
365 
366  forAll(procMeshes_, proci)
367  {
369  procFields[proci];
370 
371  // Set the face values in the reconstructed field
372 
374  {
375  // Assume all scalar surfaceFields are oriented flux fields
376  const labelList& faceMap = faceProcAddressing_[proci];
377 
378  // Correctly oriented copy of internal field
379  Field<Type> procInternalField(procField.primitiveField());
380 
381  // Addressing into original field
382  // It is necessary to create a copy of the addressing array to
383  // take care of the face direction offset trick.
384  labelList curAddr(procInternalField.size());
385 
386  forAll(procInternalField, i)
387  {
388  curAddr[i] = mag(faceMap[i]) - 1;
389  if (faceMap[i] < 0)
390  {
391  procInternalField[i] = -procInternalField[i];
392  }
393  }
394 
395  // Map
396  internalField.rmap(procInternalField, curAddr);
397  }
398  else
399  {
400  // Map
401  internalField.rmap
402  (
403  procField.primitiveField(),
404  mag(labelField(faceProcAddressing_[proci])) - 1
405  );
406  }
407 
408  // Set the boundary patch values in the reconstructed field
409  forAll(boundaryProcAddressing_[proci], patchi)
410  {
411  // Get patch index of the original patch
412  const label curBPatch = boundaryProcAddressing_[proci][patchi];
413 
414  // Get addressing slice for this patch
415  const labelList::subList cp =
416  procMeshes_[proci].boundary()[patchi].patchSlice
417  (
418  faceProcAddressing_[proci]
419  );
420 
421  // check if the boundary patch is not a processor patch
422  if (curBPatch >= 0)
423  {
424  // Regular patch. Fast looping
425 
426  if (!patchFields(curBPatch))
427  {
428  patchFields.set
429  (
430  curBPatch,
432  (
433  procField.boundaryField()[patchi],
434  mesh_.boundary()[curBPatch],
437  (
438  mesh_.boundary()[curBPatch].size()
439  )
440  )
441  );
442  }
443 
444  const label curPatchStart =
445  mesh_.boundaryMesh()[curBPatch].start();
446 
447  labelList reverseAddressing(cp.size());
448 
449  forAll(cp, facei)
450  {
451  // Subtract one to take into account offsets for
452  // face direction.
453  reverseAddressing[facei] = cp[facei] - 1 - curPatchStart;
454  }
455 
456  patchFields[curBPatch].rmap
457  (
458  procField.boundaryField()[patchi],
459  reverseAddressing
460  );
461  }
462  else
463  {
464  const Field<Type>& curProcPatch =
465  procField.boundaryField()[patchi];
466 
467  // In processor patches, there's a mix of internal faces (some
468  // of them turned) and possible cyclics. Slow loop
469  forAll(cp, facei)
470  {
471  label curF = cp[facei] - 1;
472 
473  // Is the face turned the right side round
474  if (curF >= 0)
475  {
476  // Is the face on the boundary?
477  if (curF >= mesh_.nInternalFaces())
478  {
479  label curBPatch =
480  mesh_.boundaryMesh().whichPatch(curF);
481 
482  if (!patchFields(curBPatch))
483  {
484  patchFields.set
485  (
486  curBPatch,
488  (
489  mesh_.boundary()[curBPatch].type(),
490  mesh_.boundary()[curBPatch],
492  ::null()
493  )
494  );
495  }
496 
497  // add the face
498  label curPatchFace =
499  mesh_.boundaryMesh()
500  [curBPatch].whichFace(curF);
501 
502  patchFields[curBPatch][curPatchFace] =
503  curProcPatch[facei];
504  }
505  else
506  {
507  // Internal face
508  internalField[curF] = curProcPatch[facei];
509  }
510  }
511  }
512  }
513  }
514  }
515 
516  forAll(mesh_.boundary(), patchi)
517  {
518  // add empty patches
519  if
520  (
521  isType<emptyFvPatch>(mesh_.boundary()[patchi])
522  && !patchFields(patchi)
523  )
524  {
525  patchFields.set
526  (
527  patchi,
529  (
531  mesh_.boundary()[patchi],
533  )
534  );
535  }
536  }
537 
538 
539  // Now construct and write the field
540  // setting the internalField and patchFields
542  (
544  (
545  fieldIoObject,
546  mesh_,
547  procFields[0].dimensions(),
548  internalField,
549  patchFields
550  )
551  );
552 }
553 
554 
555 template<class Type>
558 (
559  const IOobject& fieldIoObject
560 ) const
561 {
562  // Read the field for all the processors
564  (
565  procMeshes_.size()
566  );
567 
568  forAll(procMeshes_, proci)
569  {
570  procFields.set
571  (
572  proci,
574  (
575  IOobject
576  (
577  fieldIoObject.name(),
578  procMeshes_[proci].time().timeName(),
579  procMeshes_[proci],
580  IOobject::MUST_READ,
581  IOobject::NO_WRITE,
582  false
583  ),
584  procMeshes_[proci]
585  )
586  );
587  }
588 
589  return reconstructFvSurfaceField
590  (
591  IOobject
592  (
593  fieldIoObject.name(),
594  mesh_.time().timeName(),
595  mesh_,
596  IOobject::NO_READ,
597  IOobject::NO_WRITE,
598  false
599  ),
600  procFields
601  );
602 }
603 
604 
605 template<class Type>
607 (
608  const IOobjectList& objects,
609  const HashSet<word>& selectedFields
610 )
611 {
612  const word& fieldClassName = DimensionedField<Type, volMesh>::typeName;
613 
614  IOobjectList fields = objects.lookupClass(fieldClassName);
615 
616  if (fields.size())
617  {
618  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
619 
620  forAllConstIter(IOobjectList, fields, fieldIter)
621  {
622  if
623  (
624  selectedFields.empty()
625  || selectedFields.found(fieldIter()->name())
626  )
627  {
628  Info<< " " << fieldIter()->name() << endl;
629 
630  reconstructFvVolumeInternalField<Type>(*fieldIter())().write();
631 
632  nReconstructed_++;
633  }
634  }
635  Info<< endl;
636  }
637 }
638 
639 
640 template<class Type>
642 (
643  const IOobjectList& objects,
644  const HashSet<word>& selectedFields
645 )
646 {
647  const word& fieldClassName =
649 
650  IOobjectList fields = objects.lookupClass(fieldClassName);
651 
652  if (fields.size())
653  {
654  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
655 
656  forAllConstIter(IOobjectList, fields, fieldIter)
657  {
658  if
659  (
660  selectedFields.empty()
661  || selectedFields.found(fieldIter()->name())
662  )
663  {
664  Info<< " " << fieldIter()->name() << endl;
665 
666  reconstructFvVolumeField<Type>(*fieldIter())().write();
667 
668  nReconstructed_++;
669  }
670  }
671  Info<< endl;
672  }
673 }
674 
675 
676 template<class Type>
678 (
679  const IOobjectList& objects,
680  const HashSet<word>& selectedFields
681 )
682 {
683  const word& fieldClassName =
685 
686  IOobjectList fields = objects.lookupClass(fieldClassName);
687 
688  if (fields.size())
689  {
690  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
691 
692  forAllConstIter(IOobjectList, fields, fieldIter)
693  {
694  if
695  (
696  selectedFields.empty()
697  || selectedFields.found(fieldIter()->name())
698  )
699  {
700  Info<< " " << fieldIter()->name() << endl;
701 
702  reconstructFvSurfaceField<Type>(*fieldIter())().write();
703 
704  nReconstructed_++;
705  }
706  }
707  Info<< endl;
708  }
709 }
710 
711 
712 // ************************************************************************* //
A HashTable with keys but without contents.
Definition: HashSet.H:59
Field< label > labelField
Specialisation of Field<T> for label.
Definition: labelField.H:49
#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
const word & name() const
Return name.
Definition: IOobject.H:303
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
tmp< GeometricField< Type, fvPatchField, volMesh > > reconstructFvVolumeField(const IOobject &fieldIoObject, const PtrList< GeometricField< Type, fvPatchField, volMesh >> &) const
Reconstruct volume field.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void reconstructFvSurfaceFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected surface fields.
bool empty() const
Return true if the hash table is empty.
Definition: HashTableI.H:72
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:753
Traits class for primitives.
Definition: pTraits.H:50
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Generic GeometricField class.
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
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
A List obtained as a section of another List.
Definition: SubList.H:53
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
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
const Mesh & mesh() const
Return mesh.
Mapper for sizing only - does not do any actual mapping.
const Field< Type > & field() const
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
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:362
label patchi
Foam::emptyFvsPatchField.
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
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
tmp< DimensionedField< Type, volMesh > > reconstructFvVolumeInternalField(const IOobject &fieldIoObject, const PtrList< DimensionedField< Type, volMesh >> &procFields) const
Reconstruct volume internal field.
This boundary condition provides an &#39;empty&#39; condition for reduced dimensions cases, i.e. 1- and 2-D geometries. Apply this condition to patches whose normal is aligned to geometric directions that do not constitue solution directions.
void reconstructFvVolumeFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected volume fields.
A class for managing temporary objects.
Definition: PtrList.H:53
label size() const
Return the number of elements in the UList.
Definition: UListI.H:311
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
void reconstructFvVolumeInternalFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected volume internal fields.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > reconstructFvSurfaceField(const IOobject &fieldIoObject, const PtrList< GeometricField< Type, fvsPatchField, surfaceMesh >> &) const
Reconstruct surface field.