pointFieldReconstructorReconstructFields.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 
27 
28 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29 
30 template<class Type>
33 {
34  // Read the field for all the processors
36  (
37  procMeshes_.size()
38  );
39 
40  forAll(procMeshes_, proci)
41  {
42  procFields.set
43  (
44  proci,
46  (
47  IOobject
48  (
49  fieldIoObject.name(),
50  procMeshes_[proci]().time().timeName(),
51  procMeshes_[proci](),
54  ),
55  procMeshes_[proci]
56  )
57  );
58  }
59 
60 
61  // Create the internalField
62  Field<Type> internalField(mesh_.size());
63 
64  // Create the patch fields
65  PtrList<pointPatchField<Type>> patchFields(mesh_.boundary().size());
66 
67 
68  forAll(procMeshes_, proci)
69  {
71  procField = procFields[proci];
72 
73  // Get processor-to-global addressing for use in rmap
74  const labelList& procToGlobalAddr = pointProcAddressing_[proci];
75 
76  // Set the cell values in the reconstructed field
77  internalField.rmap
78  (
79  procField.primitiveField(),
80  procToGlobalAddr
81  );
82 
83  // Set the boundary patch values in the reconstructed field
84  forAll(boundaryProcAddressing_[proci], patchi)
85  {
86  // Get patch index of the original patch
87  const label curBPatch = boundaryProcAddressing_[proci][patchi];
88 
89  // check if the boundary patch is not a processor patch
90  if (curBPatch >= 0)
91  {
92  if (!patchFields(curBPatch))
93  {
94  patchFields.set(
95  curBPatch,
97  (
98  procField.boundaryField()[patchi],
99  mesh_.boundary()[curBPatch],
102  (
103  mesh_.boundary()[curBPatch].size()
104  )
105  )
106  );
107  }
108 
109  patchFields[curBPatch].rmap
110  (
111  procField.boundaryField()[patchi],
112  patchPointAddressing_[proci][patchi]
113  );
114  }
115  }
116  }
117 
118  // Construct and write the field
119  // setting the internalField and patchFields
121  (
123  (
124  IOobject
125  (
126  fieldIoObject.name(),
127  mesh_().time().timeName(),
128  mesh_(),
131  ),
132  mesh_,
133  procFields[0].dimensions(),
134  internalField,
135  patchFields
136  )
137  );
138 }
139 
140 
141 // Reconstruct and write all point fields
142 template<class Type>
144 (
145  const IOobjectList& objects,
146  const HashSet<word>& selectedFields
147 )
148 {
149  word fieldClassName
150  (
152  );
153 
154  IOobjectList fields = objects.lookupClass(fieldClassName);
155 
156  if (fields.size())
157  {
158  Info<< " Reconstructing " << fieldClassName << "s\n" << endl;
159 
160  forAllConstIter(IOobjectList, fields, fieldIter)
161  {
162  if
163  (
164  !selectedFields.size()
165  || selectedFields.found(fieldIter()->name())
166  )
167  {
168  Info<< " " << fieldIter()->name() << endl;
169 
170  reconstructField<Type>(*fieldIter())().write();
171  }
172  }
173 
174  Info<< endl;
175  }
176 }
177 
178 
179 // ************************************************************************* //
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
const word & name() const
Return name.
Definition: IOobject.H:303
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
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
Abstract base class for point-mesh patch fields.
static const DimensionedField< Type, GeoMesh > & null()
Return a null DimensionedField.
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
void reconstructFields(const IOobjectList &objects, const HashSet< word > &selectedFields)
Reconstruct and write all fields.
A class for handling words, derived from string.
Definition: word.H:59
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
label size() const
Return number of points.
Definition: pointMesh.H:91
Mapper for sizing only - does not do any actual mapping.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
IOobjectList lookupClass(const word &className) const
Return the list for all IOobjects of a given class.
Definition: IOobjectList.C:192
label patchi
const Time & time() const
Return time.
Definition: IOobject.C:318
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
messageStream Info
A class for managing temporary objects.
Definition: PtrList.H:53
const pointBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: pointMesh.H:103
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
tmp< GeometricField< Type, pointPatchField, pointMesh > > reconstructField(const IOobject &fieldIoObject)
Reconstruct field.