pointFieldDecomposerTemplates.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-2023 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 "pointFieldDecomposer.H"
27 #include "fvMesh.H"
29 
30 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
31 
32 template<class Type>
34 Foam::pointFieldDecomposer::decomposeField
35 (
36  const IOobject& fieldIoObject
37 ) const
38 {
39  // Read the field
40  const PointField<Type> field
41  (
42  IOobject
43  (
44  fieldIoObject.name(),
45  completeMesh_.db().time().name(),
46  completeMesh_.db(),
49  false
50  ),
51  completeMesh_
52  );
53 
54  // Construct the processor fields
55  PtrList<PointField<Type>> procFields(procMeshes_.size());
56  forAll(procMeshes_, proci)
57  {
58  const pointMesh& procMesh = pointMesh::New(procMeshes_[proci]);
59 
60  // Create and map the internal field values
61  Field<Type> internalField
62  (
63  field.primitiveField(),
64  pointProcAddressing_[proci]
65  );
66 
67  // Create a list of pointers for the patchFields
68  PtrList<pointPatchField<Type>> patchFields
69  (
70  procMesh.boundary().size()
71  );
72 
73  // Create and map the patch field values
74  forAll(procMesh.boundary(), patchi)
75  {
76  if (patchi < completeMesh_.boundary().size())
77  {
78  patchFields.set
79  (
80  patchi,
82  (
83  field.boundaryField()[patchi],
84  procMesh.boundary()[patchi],
86  patchFieldDecomposers_[proci][patchi]
87  )
88  );
89  }
90  else
91  {
92  patchFields.set
93  (
94  patchi,
95  new processorPointPatchField<Type>
96  (
97  procMesh.boundary()[patchi],
99  )
100  );
101  }
102  }
103 
104  // Create the field for the processor
105  procFields.set
106  (
107  proci,
108  new PointField<Type>
109  (
110  IOobject
111  (
112  field.name(),
113  procMesh().time().name(),
114  procMesh(),
117  false
118  ),
119  procMesh,
120  field.dimensions(),
121  internalField,
122  patchFields
123  )
124  );
125  }
126 
127  return procFields;
128 }
129 
130 
131 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
132 
133 template<class Type>
135 (
136  const IOobjectList& objects
137 )
138 {
139  const word& fieldClassName = PointField<Type>::typeName;
140 
141  IOobjectList fields = objects.lookupClass(fieldClassName);
142 
143  if (fields.size())
144  {
145  Info<< nl << " Decomposing " << fieldClassName << "s" << nl << endl;
146 
147  forAllConstIter(IOobjectList, fields, fieldIter)
148  {
149  Info<< " " << fieldIter()->name() << endl;
150 
151  PtrList<PointField<Type>> procFields =
152  decomposeField<Type>(*fieldIter());
153 
154  forAll(procFields, proci)
155  {
156  procFields[proci].write();
157  }
158  }
159  }
160 }
161 
162 
163 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
static pointMesh & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
static const DimensionedField< Type, GeoMesh > & null()
Return a null DimensionedField.
static const char *const typeName
Definition: Field.H:106
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
const word & name() const
Return const reference to name.
const Time & time() const
Return time.
void decomposeFields(const IOobjectList &objects)
Read, decompose and write all fields.
const pointBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: pointMesh.H:125
static autoPtr< pointPatchField< Type > > New(const word &, const pointPatch &, const DimensionedField< Type, pointMesh > &)
Return a pointer to a new patchField created on freestore given.
label patchi
Info<< "Calculating turbulent flame speed field St\n"<< endl;volScalarField St(IOobject("St", runTime.name(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), flameWrinkling->Xi() *Su);multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:228
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
messageStream Info
static const char nl
Definition: Ostream.H:266
objects