pointFieldReconstructorFields.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2018 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
30 
31 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32 
33 template<class Type>
36 {
37  // Read the field for all the processors
39  (
40  procMeshes_.size()
41  );
42 
43  forAll(procMeshes_, proci)
44  {
45  procFields.set
46  (
47  proci,
49  (
50  IOobject
51  (
52  fieldIoObject.name(),
53  procMeshes_[proci]().time().timeName(),
54  procMeshes_[proci](),
57  ),
58  procMeshes_[proci]
59  )
60  );
61  }
62 
63 
64  // Create the internalField
65  Field<Type> internalField(mesh_.size());
66 
67  // Create the patch fields
68  PtrList<pointPatchField<Type>> patchFields(mesh_.boundary().size());
69 
70 
71  forAll(procMeshes_, proci)
72  {
74  procField = procFields[proci];
75 
76  // Get processor-to-global addressing for use in rmap
77  const labelList& procToGlobalAddr = pointProcAddressing_[proci];
78 
79  // Set the cell values in the reconstructed field
80  internalField.rmap
81  (
82  procField.primitiveField(),
83  procToGlobalAddr
84  );
85 
86  // Set the boundary patch values in the reconstructed field
87  forAll(boundaryProcAddressing_[proci], patchi)
88  {
89  // Get patch index of the original patch
90  const label curBPatch = boundaryProcAddressing_[proci][patchi];
91 
92  // check if the boundary patch is not a processor patch
93  if (curBPatch >= 0)
94  {
95  if (!patchFields(curBPatch))
96  {
97  patchFields.set(
98  curBPatch,
100  (
101  procField.boundaryField()[patchi],
102  mesh_.boundary()[curBPatch],
105  (
106  mesh_.boundary()[curBPatch].size()
107  )
108  )
109  );
110  }
111 
112  patchFields[curBPatch].rmap
113  (
114  procField.boundaryField()[patchi],
115  patchPointAddressing_[proci][patchi]
116  );
117  }
118  }
119  }
120 
121  // Construct and write the field
122  // setting the internalField and patchFields
124  (
125  IOobject
126  (
127  fieldIoObject.name(),
128  mesh_().time().timeName(),
129  mesh_(),
132  ),
133  mesh_,
134  procFields[0].dimensions(),
135  internalField,
136  patchFields
137  );
138 }
139 
140 
141 template<class Type>
143 (
144  const IOobjectList& objects,
145  const UList<word>& fieldNames
146 )
147 {
149 
150  label nFields = 0;
151  for (const word& fieldName : fieldNames)
152  {
153  const IOobject* io = objects.cfindObject<fieldType>(fieldName);
154  if (io)
155  {
156  if (!nFields++)
157  {
158  Info<< " Reconstructing "
159  << fieldType::typeName << "s\n" << nl;
160  }
161  Info<< " " << fieldName << endl;
162 
163  reconstructField<Type>(*io)().write();
164  ++nReconstructed_;
165  }
166  }
167 
168  if (nFields) Info<< endl;
169  return nFields;
170 }
171 
172 
173 template<class Type>
175 (
176  const IOobjectList& objects,
177  const wordRes& selectedFields
178 )
179 {
181 
182  const wordList fieldNames =
183  (
184  selectedFields.empty()
185  ? objects.sortedNames<fieldType>()
186  : objects.sortedNames<fieldType>(selectedFields)
187  );
188 
189  return reconstructFields<Type>(objects, fieldNames);
190 }
191 
192 
193 template<class Type>
195 (
196  const IOobjectList& objects,
197  const wordHashSet& selectedFields
198 )
199 {
201 
202  const wordList fieldNames =
203  (
204  selectedFields.empty()
205  ? objects.sortedNames<fieldType>()
206  : objects.sortedNames<fieldType>(selectedFields)
207  );
208 
209  return reconstructFields<Type>(objects, fieldNames);
210 }
211 
212 
213 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::pointMesh::boundary
const pointBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: pointMesh.H:109
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::pointMesh::size
static label size(const Mesh &mesh)
Return size. Number of points.
Definition: pointMesh.H:97
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::GeometricField::primitiveField
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
Definition: GeometricFieldI.H:53
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::IOobjectList::sortedNames
wordList sortedNames() const
The sorted names of the IOobjects.
Definition: IOobjectList.C:345
Foam::DimensionedField::null
static const DimensionedField< Type, GeoMesh > & null()
Return a NullObjectRef DimensionedField.
Definition: DimensionedFieldI.H:33
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::pointMesh::time
const Time & time() const
Return Time from polyMesh.
Definition: pointMesh.H:127
Foam::HashSet< word, Hash< word > >
Foam::pointPatchField
Abstract base class for point-mesh patch fields.
Definition: pointMVCWeight.H:60
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
pointFieldReconstructor.H
Foam::pointFieldReconstructor::reconstructFields
label reconstructFields(const IOobjectList &objects, const UList< word > &fieldNames)
Reconstruct and write specified fields.
Foam::Field
Generic templated field type.
Definition: Field.H:63
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::IOobjectList::cfindObject
const IOobject * cfindObject(const word &objName) const
Return const pointer to the object found by name.
Definition: IOobjectList.C:245
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:59
fieldNames
const wordRes fieldNames(propsDict.getOrDefault< wordRes >("fields", wordRes()))
Foam::IOobjectList
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:55
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::IOobject::name
const word & name() const noexcept
Return name.
Definition: IOobjectI.H:65
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::List< label >
Foam::UList
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:103
Foam::pointFieldReconstructor::pointPatchFieldReconstructor
Definition: pointFieldReconstructor.H:89
Foam::wordRes
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:51
Foam::vtk::write
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Definition: foamVtkOutputTemplates.C:36
Foam::GeometricField
Generic GeometricField class.
Definition: areaFieldsFwd.H:53
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::pointFieldReconstructor::reconstructField
tmp< GeometricField< Type, pointPatchField, pointMesh > > reconstructField(const IOobject &fieldIoObject)
Reconstruct field.
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::IOobject::MUST_READ
Definition: IOobject.H:185