fvFieldReconstructorReconstructFields.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2012 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(),
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  ),
100  procMeshes_[procI]
101  )
102  );
103  }
104 
105 
106  return reconstructFvVolumeInternalField
107  (
108  IOobject
109  (
110  fieldIoObject.name(),
111  mesh_.time().timeName(),
112  mesh_,
113  IOobject::NO_READ,
114  IOobject::NO_WRITE
115  ),
116  procFields
117  );
118 }
119 
120 
121 template<class Type>
124 (
125  const IOobject& fieldIoObject,
127 ) const
128 {
129  // Create the internalField
130  Field<Type> internalField(mesh_.nCells());
131 
132  // Create the patch fields
133  PtrList<fvPatchField<Type> > patchFields(mesh_.boundary().size());
134 
135  forAll(procFields, procI)
136  {
138  procFields[procI];
139 
140  // Set the cell values in the reconstructed field
141  internalField.rmap
142  (
143  procField.internalField(),
144  cellProcAddressing_[procI]
145  );
146 
147  // Set the boundary patch values in the reconstructed field
148  forAll(boundaryProcAddressing_[procI], patchI)
149  {
150  // Get patch index of the original patch
151  const label curBPatch = boundaryProcAddressing_[procI][patchI];
152 
153  // Get addressing slice for this patch
154  const labelList::subList cp =
155  procField.mesh().boundary()[patchI].patchSlice
156  (
157  faceProcAddressing_[procI]
158  );
159 
160  // check if the boundary patch is not a processor patch
161  if (curBPatch >= 0)
162  {
163  // Regular patch. Fast looping
164 
165  if (!patchFields(curBPatch))
166  {
167  patchFields.set
168  (
169  curBPatch,
171  (
172  procField.boundaryField()[patchI],
173  mesh_.boundary()[curBPatch],
176  (
177  mesh_.boundary()[curBPatch].size()
178  )
179  )
180  );
181  }
182 
183  const label curPatchStart =
184  mesh_.boundaryMesh()[curBPatch].start();
185 
186  labelList reverseAddressing(cp.size());
187 
188  forAll(cp, faceI)
189  {
190  // Check
191  if (cp[faceI] <= 0)
192  {
194  (
195  "fvFieldReconstructor::reconstructFvVolumeField\n"
196  "(\n"
197  " const IOobject&,\n"
198  " const PtrList<GeometricField<Type,"
199  " fvPatchField, volMesh> >&\n"
200  ") const\n"
201  ) << "Processor " << procI
202  << " patch "
203  << procField.mesh().boundary()[patchI].name()
204  << " face " << faceI
205  << " originates from reversed face since "
206  << cp[faceI]
207  << exit(FatalError);
208  }
209 
210  // Subtract one to take into account offsets for
211  // face direction.
212  reverseAddressing[faceI] = cp[faceI] - 1 - curPatchStart;
213  }
214 
215 
216  patchFields[curBPatch].rmap
217  (
218  procField.boundaryField()[patchI],
219  reverseAddressing
220  );
221  }
222  else
223  {
224  const Field<Type>& curProcPatch =
225  procField.boundaryField()[patchI];
226 
227  // In processor patches, there's a mix of internal faces (some
228  // of them turned) and possible cyclics. Slow loop
229  forAll(cp, faceI)
230  {
231  // Subtract one to take into account offsets for
232  // face direction.
233  label curF = cp[faceI] - 1;
234 
235  // Is the face on the boundary?
236  if (curF >= mesh_.nInternalFaces())
237  {
238  label curBPatch = mesh_.boundaryMesh().whichPatch(curF);
239 
240  if (!patchFields(curBPatch))
241  {
242  patchFields.set
243  (
244  curBPatch,
246  (
247  mesh_.boundary()[curBPatch].type(),
248  mesh_.boundary()[curBPatch],
250  )
251  );
252  }
253 
254  // add the face
255  label curPatchFace =
256  mesh_.boundaryMesh()
257  [curBPatch].whichFace(curF);
258 
259  patchFields[curBPatch][curPatchFace] =
260  curProcPatch[faceI];
261  }
262  }
263  }
264  }
265  }
266 
267  forAll(mesh_.boundary(), patchI)
268  {
269  // add empty patches
270  if
271  (
272  isType<emptyFvPatch>(mesh_.boundary()[patchI])
273  && !patchFields(patchI)
274  )
275  {
276  patchFields.set
277  (
278  patchI,
280  (
282  mesh_.boundary()[patchI],
284  )
285  );
286  }
287  }
288 
289 
290  // Now construct and write the field
291  // setting the internalField and patchFields
293  (
295  (
296  fieldIoObject,
297  mesh_,
298  procFields[0].dimensions(),
300  patchFields
301  )
302  );
303 }
304 
305 
306 template<class Type>
309 (
310  const IOobject& fieldIoObject
311 ) const
312 {
313  // Read the field for all the processors
315  (
316  procMeshes_.size()
317  );
318 
319  forAll(procMeshes_, procI)
320  {
321  procFields.set
322  (
323  procI,
325  (
326  IOobject
327  (
328  fieldIoObject.name(),
329  procMeshes_[procI].time().timeName(),
330  procMeshes_[procI],
331  IOobject::MUST_READ,
332  IOobject::NO_WRITE
333  ),
334  procMeshes_[procI]
335  )
336  );
337  }
338 
339  return reconstructFvVolumeField
340  (
341  IOobject
342  (
343  fieldIoObject.name(),
344  mesh_.time().timeName(),
345  mesh_,
346  IOobject::NO_READ,
347  IOobject::NO_WRITE
348  ),
349  procFields
350  );
351 }
352 
353 
354 template<class Type>
357 (
358  const IOobject& fieldIoObject,
360 ) const
361 {
362  // Create the internalField
363  Field<Type> internalField(mesh_.nInternalFaces());
364 
365  // Create the patch fields
366  PtrList<fvsPatchField<Type> > patchFields(mesh_.boundary().size());
367 
368 
369  forAll(procMeshes_, procI)
370  {
372  procFields[procI];
373 
374  // Set the face values in the reconstructed field
375 
376  // It is necessary to create a copy of the addressing array to
377  // take care of the face direction offset trick.
378  //
379  {
380  const labelList& faceMap = faceProcAddressing_[procI];
381 
382  // Correctly oriented copy of internal field
383  Field<Type> procInternalField(procField.internalField());
384  // Addressing into original field
385  labelList curAddr(procInternalField.size());
386 
387  forAll(procInternalField, addrI)
388  {
389  curAddr[addrI] = mag(faceMap[addrI])-1;
390  if (faceMap[addrI] < 0)
391  {
392  procInternalField[addrI] = -procInternalField[addrI];
393  }
394  }
395 
396  // Map
397  internalField.rmap(procInternalField, curAddr);
398  }
399 
400  // Set the boundary patch values in the reconstructed field
401  forAll(boundaryProcAddressing_[procI], patchI)
402  {
403  // Get patch index of the original patch
404  const label curBPatch = boundaryProcAddressing_[procI][patchI];
405 
406  // Get addressing slice for this patch
407  const labelList::subList cp =
408  procMeshes_[procI].boundary()[patchI].patchSlice
409  (
410  faceProcAddressing_[procI]
411  );
412 
413  // check if the boundary patch is not a processor patch
414  if (curBPatch >= 0)
415  {
416  // Regular patch. Fast looping
417 
418  if (!patchFields(curBPatch))
419  {
420  patchFields.set
421  (
422  curBPatch,
424  (
425  procField.boundaryField()[patchI],
426  mesh_.boundary()[curBPatch],
429  (
430  mesh_.boundary()[curBPatch].size()
431  )
432  )
433  );
434  }
435 
436  const label curPatchStart =
437  mesh_.boundaryMesh()[curBPatch].start();
438 
439  labelList reverseAddressing(cp.size());
440 
441  forAll(cp, faceI)
442  {
443  // Subtract one to take into account offsets for
444  // face direction.
445  reverseAddressing[faceI] = cp[faceI] - 1 - curPatchStart;
446  }
447 
448  patchFields[curBPatch].rmap
449  (
450  procField.boundaryField()[patchI],
451  reverseAddressing
452  );
453  }
454  else
455  {
456  const Field<Type>& curProcPatch =
457  procField.boundaryField()[patchI];
458 
459  // In processor patches, there's a mix of internal faces (some
460  // of them turned) and possible cyclics. Slow loop
461  forAll(cp, faceI)
462  {
463  label curF = cp[faceI] - 1;
464 
465  // Is the face turned the right side round
466  if (curF >= 0)
467  {
468  // Is the face on the boundary?
469  if (curF >= mesh_.nInternalFaces())
470  {
471  label curBPatch =
472  mesh_.boundaryMesh().whichPatch(curF);
473 
474  if (!patchFields(curBPatch))
475  {
476  patchFields.set
477  (
478  curBPatch,
480  (
481  mesh_.boundary()[curBPatch].type(),
482  mesh_.boundary()[curBPatch],
484  ::null()
485  )
486  );
487  }
488 
489  // add the face
490  label curPatchFace =
491  mesh_.boundaryMesh()
492  [curBPatch].whichFace(curF);
493 
494  patchFields[curBPatch][curPatchFace] =
495  curProcPatch[faceI];
496  }
497  else
498  {
499  // Internal face
500  internalField[curF] = curProcPatch[faceI];
501  }
502  }
503  }
504  }
505  }
506  }
507 
508  forAll(mesh_.boundary(), patchI)
509  {
510  // add empty patches
511  if
512  (
513  isType<emptyFvPatch>(mesh_.boundary()[patchI])
514  && !patchFields(patchI)
515  )
516  {
517  patchFields.set
518  (
519  patchI,
521  (
523  mesh_.boundary()[patchI],
525  )
526  );
527  }
528  }
529 
530 
531  // Now construct and write the field
532  // setting the internalField and patchFields
534  (
536  (
537  fieldIoObject,
538  mesh_,
539  procFields[0].dimensions(),
541  patchFields
542  )
543  );
544 }
545 
546 
547 template<class Type>
550 (
551  const IOobject& fieldIoObject
552 ) const
553 {
554  // Read the field for all the processors
556  (
557  procMeshes_.size()
558  );
559 
560  forAll(procMeshes_, procI)
561  {
562  procFields.set
563  (
564  procI,
566  (
567  IOobject
568  (
569  fieldIoObject.name(),
570  procMeshes_[procI].time().timeName(),
571  procMeshes_[procI],
572  IOobject::MUST_READ,
573  IOobject::NO_WRITE
574  ),
575  procMeshes_[procI]
576  )
577  );
578  }
579 
580  return reconstructFvSurfaceField
581  (
582  IOobject
583  (
584  fieldIoObject.name(),
585  mesh_.time().timeName(),
586  mesh_,
587  IOobject::NO_READ,
588  IOobject::NO_WRITE
589  ),
590  procFields
591  );
592 }
593 
594 
595 template<class Type>
597 (
598  const IOobjectList& objects,
599  const HashSet<word>& selectedFields
600 )
601 {
602  const word& fieldClassName = DimensionedField<Type, volMesh>::typeName;
603 
604  IOobjectList fields = objects.lookupClass(fieldClassName);
605 
606  if (fields.size())
607  {
608  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
609 
610  forAllConstIter(IOobjectList, fields, fieldIter)
611  {
612  if
613  (
614  selectedFields.empty()
615  || selectedFields.found(fieldIter()->name())
616  )
617  {
618  Info<< " " << fieldIter()->name() << endl;
619 
620  reconstructFvVolumeInternalField<Type>(*fieldIter())().write();
621 
622  nReconstructed_++;
623  }
624  }
625  Info<< endl;
626  }
627 }
628 
629 
630 template<class Type>
632 (
633  const IOobjectList& objects,
634  const HashSet<word>& selectedFields
635 )
636 {
637  const word& fieldClassName =
639 
640  IOobjectList fields = objects.lookupClass(fieldClassName);
641 
642  if (fields.size())
643  {
644  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
645 
646  forAllConstIter(IOobjectList, fields, fieldIter)
647  {
648  if
649  (
650  selectedFields.empty()
651  || selectedFields.found(fieldIter()->name())
652  )
653  {
654  Info<< " " << fieldIter()->name() << endl;
655 
656  reconstructFvVolumeField<Type>(*fieldIter())().write();
657 
658  nReconstructed_++;
659  }
660  }
661  Info<< endl;
662  }
663 }
664 
665 
666 template<class Type>
668 (
669  const IOobjectList& objects,
670  const HashSet<word>& selectedFields
671 )
672 {
673  const word& fieldClassName =
675 
676  IOobjectList fields = objects.lookupClass(fieldClassName);
677 
678  if (fields.size())
679  {
680  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
681 
682  forAllConstIter(IOobjectList, fields, fieldIter)
683  {
684  if
685  (
686  selectedFields.empty()
687  || selectedFields.found(fieldIter()->name())
688  )
689  {
690  Info<< " " << fieldIter()->name() << endl;
691 
692  reconstructFvSurfaceField<Type>(*fieldIter())().write();
693 
694  nReconstructed_++;
695  }
696  }
697  Info<< endl;
698  }
699 }
700 
701 
702 // ************************************************************************* //
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.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
Info<< "Creating field dpdt\n"<< endl;volScalarField dpdt(IOobject("dpdt", runTime.timeName(), mesh), mesh, dimensionedScalar("dpdt", p.dimensions()/dimTime, 0));Info<< "Creating field kinetic energy K\n"<< endl;volScalarField K("K", 0.5 *magSqr(U));volScalarField p_rgh(IOobject("p_rgh", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);p_rgh=p-rho *gh;mesh.setFluxRequired(p_rgh.name());multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:127
const Field< Type > & field() const
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
dimensioned< scalar > mag(const dimensioned< Type > &)
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Foam::emptyFvsPatchField.
A class for handling words, derived from string.
Definition: word.H:59
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
InternalField & internalField()
Return internal field.
IOobjectList lookupClass(const word &className) const
Return the list for all IOobjects of a given class.
Definition: IOobjectList.C:199
void reconstructFvVolumeInternalFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected volume internal fields.
messageStream Info
const Mesh & mesh() const
Return mesh.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define forAll(list, i)
Definition: UList.H:421
virtual Ostream & write(const token &)=0
Write next token to stream.
tmp< GeometricField< Type, fvPatchField, volMesh > > reconstructFvVolumeField(const IOobject &fieldIoObject, const PtrList< GeometricField< Type, fvPatchField, volMesh > > &) const
Reconstruct volume field.
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
const word & name() const
Return name.
Definition: IOobject.H:260
void reconstructFvSurfaceFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected surface fields.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Generic GeometricField class.
error FatalError
A HashTable with keys but without contents.
Definition: HashSet.H:59
A List obtained as a section of another List.
Definition: SubList.H:53
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > reconstructFvSurfaceField(const IOobject &fieldIoObject, const PtrList< GeometricField< Type, fvsPatchField, surfaceMesh > > &) const
Reconstruct surface field.
bool empty() const
Return true if the hash table is empty.
tmp< DimensionedField< Type, volMesh > > reconstructFvVolumeInternalField(const IOobject &fieldIoObject, const PtrList< DimensionedField< Type, volMesh > > &procFields) const
Reconstruct volume internal field.
Mapper for sizing only - does not do any actual mapping.
void reconstructFvVolumeFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Read, reconstruct and write all/selected volume fields.
bool found(const Key &) const
Return true if hashedEntry is found in table.
A class for managing temporary objects.
Definition: PtrList.H:118
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
conserve internalField()+
const volScalarField & cp