fvMeshDistributeTemplates.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 "polyTopoChangeMap.H"
27 #include "processorFvPatchField.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class GeoField>
33 {
34  const UPtrList<GeoField> fields(mesh.fields<GeoField>());
35 
36  forAll(fields, i)
37  {
38  const GeoField& field = fields[i];
39 
40  Pout<< "Field:" << field.name() << " internal size:" << field.size()
41  << endl;
42 
43  forAll(field.boundaryField(), patchi)
44  {
45  Pout<< " " << patchi
46  << ' ' << field.boundaryField()[patchi].patch().name()
47  << ' ' << field.boundaryField()[patchi].type()
48  << ' ' << field.boundaryField()[patchi].size()
49  << endl;
50  }
51  }
52 }
53 
54 
55 template<class Type, class Mesh>
56 void Foam::fvMeshDistribute::saveBoundaryFields
57 (
59 ) const
60 {
61  // Save whole boundary field
62 
64  (
65  mesh_.fields<SurfaceField<Type>>()
66  );
67 
68  bfields.setSize(fields.size());
69 
70  forAll(fields, i)
71  {
72  const SurfaceField<Type>& field = fields[i];
73 
74  bfields.set(i, field.boundaryField().clone().ptr());
75  }
76 }
77 
78 
79 template<class Type, class Mesh>
80 void Foam::fvMeshDistribute::mapBoundaryFields
81 (
82  const polyTopoChangeMap& map,
83  const PtrList<FieldField<fvsPatchField, Type>>& oldBfields
84 )
85 {
86  // Map boundary field
87 
88  const labelList& oldPatchStarts = map.oldPatchStarts();
89  const labelList& faceMap = map.faceMap();
90 
91  UPtrList<SurfaceField<Type>> fields
92  (
93  mesh_.fields<SurfaceField<Type>>()
94  );
95 
96  forAll(fields, i)
97  {
98  SurfaceField<Type>& field = fields[i];
99  typename SurfaceField<Type>::Boundary& bfield =
100  field.boundaryFieldRef();
101 
102  const FieldField<fvsPatchField, Type>& oldBfield = oldBfields[i];
103 
104  // Pull from old boundary field into bfield.
105 
106  forAll(bfield, patchi)
107  {
108  fvsPatchField<Type>& patchField = bfield[patchi];
109  label facei = patchField.patch().start();
110 
111  forAll(patchField, i)
112  {
113  label oldFacei = faceMap[facei++];
114 
115  // Find patch and local patch face oldFacei was in.
116  forAll(oldPatchStarts, oldPatchi)
117  {
118  label oldLocalI = oldFacei - oldPatchStarts[oldPatchi];
119 
120  if
121  (
122  oldLocalI >= 0
123  && oldLocalI < oldBfield[oldPatchi].size()
124  )
125  {
126  patchField[i] = oldBfield[oldPatchi][oldLocalI];
127  }
128  }
129  }
130  }
131  }
132 }
133 
134 
135 template<class Type>
136 void Foam::fvMeshDistribute::initMapExposedFaces
137 (
138  PtrList<Field<Type>>& ifields
139 ) const
140 {
141  const UPtrList<SurfaceField<Type>> fields
142  (
143  mesh_.fields<SurfaceField<Type>>()
144  );
145 
146  ifields.setSize(fields.size());
147 
148  forAll(fields, i)
149  {
150  ifields.set(i, Field<Type>(mesh_.nFaces()));
151 
152  const SurfaceField<Type>& field = fields[i];
153 
154  SubList<Type>(ifields[i], field.primitiveField().size()) =
155  field.primitiveField();
156 
157  forAll(field.boundaryField(), patchi)
158  {
159  const fvsPatchField<Type>& pfield = field.boundaryField()[patchi];
160 
161  SubList<Type>(ifields[i], pfield.size(), pfield.patch().start()) =
162  pfield;
163  }
164  }
165 }
166 
167 
168 template<class Type>
169 void Foam::fvMeshDistribute::mapExposedFaces
170 (
171  const polyTopoChangeMap& map,
172  const PtrList<Field<Type>>& oldFields
173 )
174 {
175  UPtrList<SurfaceField<Type>> fields
176  (
177  mesh_.fields<SurfaceField<Type>>()
178  );
179 
180  forAll(fields, i)
181  {
182  SurfaceField<Type>& field = fields[i];
183 
184  const Field<Type>& oldField = oldFields[i];
185 
186  const bool negateIfFlipped = isFlux(field);
187 
188  forAll(field.boundaryField(), patchi)
189  {
190  fvsPatchField<Type>& patchField = field.boundaryFieldRef()[patchi];
191 
192  forAll(patchField, i)
193  {
194  const label facei = patchField.patch().start()+i;
195  const label oldFacei = map.faceMap()[facei];
196 
197  if (oldFacei < map.nOldInternalFaces())
198  {
199  if (negateIfFlipped && map.flipFaceFlux().found(facei))
200  {
201  patchField[i] = flipOp()(oldField[oldFacei]);
202  }
203  else
204  {
205  patchField[i] = oldField[oldFacei];
206  }
207  }
208  else
209  {
210  patchField[i] = oldField[oldFacei];
211  }
212  }
213  }
214  }
215 }
216 
217 
218 template<class GeoField>
219 void Foam::fvMeshDistribute::correctCoupledPatchFields()
220 {
221  UPtrList<GeoField> fields(mesh_.fields<GeoField>());
222 
223  // Ensure the deltaCoeffs are available for constraint patch evaluation
224  mesh_.deltaCoeffs();
225 
226  forAll(fields, i)
227  {
228  GeoField& field = fields[i];
229 
230  typename GeoField::Boundary& bfield = field.boundaryFieldRef();
231 
232  if
233  (
236  )
237  {
238  label nReq = Pstream::nRequests();
239 
240  forAll(bfield, patchi)
241  {
242  if (bfield[patchi].coupled())
243  {
244  bfield[patchi].initEvaluate(Pstream::defaultCommsType);
245  }
246  }
247 
248  // Block for any outstanding requests
249  if
250  (
253  )
254  {
255  Pstream::waitRequests(nReq);
256  }
257 
258  forAll(bfield, patchi)
259  {
260  if (bfield[patchi].coupled())
261  {
262  bfield[patchi].evaluate(Pstream::defaultCommsType);
263  }
264  }
265  }
267  {
268  const lduSchedule& patchSchedule =
269  mesh_.globalData().patchSchedule();
270 
271  forAll(patchSchedule, patchEvali)
272  {
273  if (bfield[patchEvali].coupled())
274  {
275  if (patchSchedule[patchEvali].init)
276  {
277  bfield[patchSchedule[patchEvali].patch]
278  .initEvaluate(Pstream::commsTypes::scheduled);
279  }
280  else
281  {
282  bfield[patchSchedule[patchEvali].patch]
284  }
285  }
286  }
287  }
288  }
289 }
290 
291 
292 template<class GeoField>
293 void Foam::fvMeshDistribute::sendFields
294 (
295  const label domain,
296  const wordList& fieldNames,
297  const fvMeshSubset& subsetter,
298  Ostream& toNbr
299 )
300 {
301  // Send fields. Note order supplied so we can receive in exactly the same
302  // order.
303  // Note that field gets written as entry in dictionary so we
304  // can construct from subdictionary.
305  // (since otherwise the reading as-a-dictionary mixes up entries from
306  // consecutive fields)
307  // The dictionary constructed is:
308  // volScalarField
309  // {
310  // p {internalField ..; boundaryField ..;}
311  // k {internalField ..; boundaryField ..;}
312  // }
313  // volVectorField
314  // {
315  // U {internalField ... }
316  // }
317 
318  // volVectorField {U {internalField ..; boundaryField ..;}}
319 
320  toNbr << GeoField::typeName << token::NL << token::BEGIN_BLOCK << token::NL;
321  forAll(fieldNames, i)
322  {
323  if (debug)
324  {
325  Pout<< "Subsetting field " << fieldNames[i]
326  << " for domain:" << domain << endl;
327  }
328 
329  // Send all fieldNames. This has to be exactly the same set as is
330  // being received!
331  const GeoField& field =
332  subsetter.baseMesh().lookupObject<GeoField>(fieldNames[i]);
333 
334  tmp<GeoField> tsubfield = subsetter.interpolate(field);
335 
336  toNbr
338  << tsubfield
340  }
341  toNbr << token::END_BLOCK << token::NL;
342 }
343 
344 
345 template<class GeoField>
346 void Foam::fvMeshDistribute::receiveFields
347 (
348  const label domain,
349  const wordList& fieldNames,
350  typename GeoField::Mesh& mesh,
351  PtrList<GeoField>& fields,
352  const dictionary& fieldDicts
353 )
354 {
355  if (debug)
356  {
357  Pout<< "Receiving fields " << fieldNames
358  << " from domain:" << domain << endl;
359  }
360 
361  fields.setSize(fieldNames.size());
362 
363  forAll(fieldNames, i)
364  {
365  if (debug)
366  {
367  Pout<< "Constructing field " << fieldNames[i]
368  << " from domain:" << domain << endl;
369  }
370 
371  fields.set
372  (
373  i,
374  new GeoField
375  (
376  IOobject
377  (
378  fieldNames[i],
379  mesh.thisDb().time().name(),
380  mesh.thisDb(),
383  ),
384  mesh,
385  fieldDicts.subDict(fieldNames[i])
386  )
387  );
388  }
389 }
390 
391 
392 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Generic field type.
Definition: FieldField.H:77
tmp< FieldField< Field, Type > > clone() const
Clone.
Definition: FieldField.C:188
Generic GeometricField class.
GeometricBoundaryField< Type, PatchField, GeoMesh > Boundary
Type of the boundary field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:85
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
static label nRequests()
Get number of outstanding requests.
Definition: UPstream.C:137
static void waitRequests(const label start=0)
Wait until all requests (from start onwards) have finished.
Definition: UPstream.C:147
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static commsTypes defaultCommsType
Default commsType.
Definition: UPstream.H:272
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: UPtrList.H:66
static void printFieldInfo(const fvMesh &)
Print some field info.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
UPtrList< GeoField > fields(const bool strict=false) const
Return the list of fields of type GeoField.
@ BEGIN_BLOCK
Definition: token.H:113
@ END_BLOCK
Definition: token.H:114
label patchi
static List< word > fieldNames
Definition: globalFoam.H:46
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
List< word > wordList
A List of words.
Definition: fileName.H:54
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
bool isFlux(const DimensionedField< Type, surfaceMesh > &df)
Check if surfaceField is a flux.
Definition: surfaceFields.H:55
List< lduScheduleEntry > lduSchedule
Definition: lduSchedule.H:74
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53