pointFieldReconstructorTemplates.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-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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
33template<class Type>
36(
37 const IOobject& fieldObject,
39) const
40{
42
43 // Create the internalField
44 Field<Type> internalField(mesh_.size());
45
46 // Create the patch fields
47 PtrList<pointPatchField<Type>> patchFields(mesh_.boundary().size());
48
49
50 forAll(procMeshes_, proci)
51 {
52 const auto& procField = procFields[proci];
53
54 // Get processor-to-global addressing for use in rmap
55 const labelList& procToGlobalAddr = pointProcAddressing_[proci];
56
57 // Set the cell values in the reconstructed field
58 internalField.rmap
59 (
60 procField.primitiveField(),
61 procToGlobalAddr
62 );
63
64 // Set the boundary patch values in the reconstructed field
65 forAll(boundaryProcAddressing_[proci], patchi)
66 {
67 // Get patch index of the original patch
68 const label curBPatch = boundaryProcAddressing_[proci][patchi];
69
70 // check if the boundary patch is not a processor patch
71 if (curBPatch >= 0)
72 {
73 if (!patchFields(curBPatch))
74 {
75 patchFields.set(
76 curBPatch,
78 (
79 procField.boundaryField()[patchi],
80 mesh_.boundary()[curBPatch],
83 (
84 mesh_.boundary()[curBPatch].size()
85 )
86 )
87 );
88 }
89
90 patchFields[curBPatch].rmap
91 (
92 procField.boundaryField()[patchi],
93 patchPointAddressing_[proci][patchi]
94 );
95 }
96 }
97 }
98
99 // Construct and write the field
100 // setting the internalField and patchFields
102 (
104 (
105 fieldObject.name(),
106 mesh_.thisDb().time().timeName(),
107 mesh_.thisDb(),
110 ),
111 mesh_,
112 procFields[0].dimensions(),
113 internalField,
114 patchFields
115 );
116}
117
118
119template<class Type>
122(
123 const IOobject& fieldObject
124)
125{
127
128 // Read the field for all the processors
129 PtrList<fieldType> procFields(procMeshes_.size());
130
131 forAll(procMeshes_, proci)
132 {
133 procFields.set
134 (
135 proci,
136 new fieldType
137 (
139 (
140 fieldObject.name(),
141 procMeshes_[proci].thisDb().time().timeName(),
142 procMeshes_[proci].thisDb(),
145 ),
146 procMeshes_[proci]
147 )
148 );
149 }
150
151 return reconstructField
152 (
154 (
155 fieldObject.name(),
156 mesh_.thisDb().time().timeName(),
157 mesh_.thisDb(),
160 ),
161 procFields
162 );
163}
164
165
166template<class Type>
168(
169 const UPtrList<const IOobject>& fieldObjects
170)
171{
173
174 label nFields = 0;
175
176 for (const IOobject& io : fieldObjects)
177 {
178 if (io.isHeaderClass<fieldType>())
179 {
180 if (verbose_)
181 {
182 if (!nFields)
183 {
184 Info<< " Reconstructing "
185 << fieldType::typeName << "s\n" << nl;
186 }
187 Info<< " " << io.name() << endl;
188 }
189 ++nFields;
190
191 reconstructPointField<Type>(io)().write();
192 ++nReconstructed_;
193 }
194 }
195
196 if (verbose_ && nFields) Info<< endl;
197 return nFields;
198}
199
200
201template<class Type>
203(
204 const IOobjectList& objects,
205 const wordRes& selectedFields
206)
207{
209
210 return reconstructPointFields<Type>
211 (
212 (
213 selectedFields.empty()
214 ? objects.sorted<fieldType>()
215 : objects.sorted<fieldType>(selectedFields)
216 )
217 );
218}
219
220
221// ************************************************************************* //
Generic templated field type.
Definition: Field.H:82
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:466
Generic GeometricField class.
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:59
UPtrList< const IOobject > sorted() const
The sorted list of IOobjects.
Definition: IOobjectList.C:336
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
bool isHeaderClass() const
Check if headerClassName() equals Type::typeName.
Definition: IOobjectI.H:156
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
const T * set(const label i) const
Definition: PtrList.H:138
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
bool empty() const noexcept
True if the UList is empty (ie, size() is zero)
Definition: UListI.H:427
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
Definition: UPtrList.H:71
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
const Time & time() const noexcept
Return time registry.
tmp< GeometricField< Type, pointPatchField, pointMesh > > reconstructPointField(const IOobject &fieldObject)
Read and reconstruct point field.
label reconstructPointFields(const UPtrList< const IOobject > &fieldObjects)
Reconstruct and write specified point fields.
tmp< GeometricField< Type, pointPatchField, pointMesh > > reconstructField(const IOobject &fieldObject, const PtrList< GeometricField< Type, pointPatchField, pointMesh > > &) const
Reconstruct field.
static label size(const Mesh &mesh)
Return size. Number of points.
Definition: pointMesh.H:102
const objectRegistry & thisDb() const
Return database. For now is its polyMesh.
Definition: pointMesh.H:126
const pointBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: pointMesh.H:114
Abstract base class for point-mesh patch fields.
A class for managing temporary objects.
Definition: tmp.H:65
A List of wordRe with additional matching capabilities.
Definition: wordRes.H:54
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
runTime write()
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333