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-2022 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 {
35  (
36  mesh.objectRegistry::lookupClass<GeoField>()
37  );
38 
39  forAllConstIter(typename HashTable<const GeoField*>, flds, iter)
40  {
41  const GeoField& fld = *iter();
42 
43  Pout<< "Field:" << iter.key() << " internal size:" << fld.size()
44  << endl;
45 
46  forAll(fld.boundaryField(), patchi)
47  {
48  Pout<< " " << patchi
49  << ' ' << fld.boundaryField()[patchi].patch().name()
50  << ' ' << fld.boundaryField()[patchi].type()
51  << ' ' << fld.boundaryField()[patchi].size()
52  << endl;
53  }
54  }
55 }
56 
57 
58 template<class Type, class Mesh>
59 void Foam::fvMeshDistribute::saveBoundaryFields
60 (
62 ) const
63 {
64  // Save whole boundary field
65 
67  (
68  static_cast<const fvMesh&>(mesh_)
70  );
71 
72  bflds.setSize(flds.size());
73 
74  label i = 0;
75 
76  forAllConstIter(typename HashTable<const SurfaceField<Type>*>, flds, iter)
77  {
78  const SurfaceField<Type>& fld = *iter();
79 
80  bflds.set(i, fld.boundaryField().clone().ptr());
81 
82  i++;
83  }
84 }
85 
86 
87 template<class Type, class Mesh>
88 void Foam::fvMeshDistribute::mapBoundaryFields
89 (
90  const polyTopoChangeMap& map,
91  const PtrList<FieldField<fvsPatchField, Type>>& oldBflds
92 )
93 {
94  // Map boundary field
95 
96  const labelList& oldPatchStarts = map.oldPatchStarts();
97  const labelList& faceMap = map.faceMap();
98 
99  HashTable<SurfaceField<Type>*> flds
100  (
101  mesh_.objectRegistry::lookupClass<SurfaceField<Type>>()
102  );
103 
104  if (flds.size() != oldBflds.size())
105  {
107  << abort(FatalError);
108  }
109 
110  label fieldi = 0;
111 
112  forAllIter(typename HashTable<SurfaceField<Type>*>, flds, iter)
113  {
114  SurfaceField<Type>& fld = *iter();
115  typename SurfaceField<Type>::Boundary& bfld =
116  fld.boundaryFieldRef();
117 
118  const FieldField<fvsPatchField, Type>& oldBfld = oldBflds[fieldi++];
119 
120  // Pull from old boundary field into bfld.
121 
122  forAll(bfld, patchi)
123  {
124  fvsPatchField<Type>& patchFld = bfld[patchi];
125  label facei = patchFld.patch().start();
126 
127  forAll(patchFld, i)
128  {
129  label oldFacei = faceMap[facei++];
130 
131  // Find patch and local patch face oldFacei was in.
132  forAll(oldPatchStarts, oldPatchi)
133  {
134  label oldLocalI = oldFacei - oldPatchStarts[oldPatchi];
135 
136  if (oldLocalI >= 0 && oldLocalI < oldBfld[oldPatchi].size())
137  {
138  patchFld[i] = oldBfld[oldPatchi][oldLocalI];
139  }
140  }
141  }
142  }
143  }
144 }
145 
146 
147 template<class Type>
148 void Foam::fvMeshDistribute::initMapExposedFaces
149 (
150  PtrList<Field<Type>>& iflds
151 ) const
152 {
153  HashTable<const SurfaceField<Type>*> flds
154  (
155  static_cast<const fvMesh&>(mesh_).lookupClass<SurfaceField<Type>>()
156  );
157 
158  iflds.setSize(flds.size());
159 
160  label fieldi = 0;
161 
162  forAllConstIter(typename HashTable<const SurfaceField<Type>*>, flds, iter)
163  {
164  iflds.set(fieldi, Field<Type>(mesh_.nFaces()));
165 
166  const SurfaceField<Type>& fld = *iter();
167 
168  SubList<Type>(iflds[fieldi], fld.primitiveField().size()) =
169  fld.primitiveField();
170 
171  forAll(fld.boundaryField(), patchi)
172  {
173  const fvsPatchField<Type>& pfld = fld.boundaryField()[patchi];
174 
175  SubList<Type>(iflds[fieldi], pfld.size(), pfld.patch().start()) =
176  pfld;
177  }
178 
179  fieldi++;
180  }
181 }
182 
183 
184 template<class Type>
185 void Foam::fvMeshDistribute::mapExposedFaces
186 (
187  const polyTopoChangeMap& map,
188  const PtrList<Field<Type>>& oldFlds
189 )
190 {
191  HashTable<SurfaceField<Type>*> flds
192  (
193  mesh_.objectRegistry::lookupClass<SurfaceField<Type>>()
194  );
195 
196  label fieldi = 0;
197 
198  forAllIter(typename HashTable<SurfaceField<Type>*>, flds, iter)
199  {
200  SurfaceField<Type>& fld = *iter();
201 
202  const Field<Type>& oldFld = oldFlds[fieldi];
203 
204  const bool negateIfFlipped = isFlux(fld);
205 
206  forAll(fld.boundaryField(), patchi)
207  {
208  fvsPatchField<Type>& patchFld = fld.boundaryFieldRef()[patchi];
209 
210  forAll(patchFld, i)
211  {
212  const label facei = patchFld.patch().start()+i;
213  const label oldFacei = map.faceMap()[facei];
214 
215  if (oldFacei < map.nOldInternalFaces())
216  {
217  if (negateIfFlipped && map.flipFaceFlux().found(facei))
218  {
219  patchFld[i] = flipOp()(oldFld[oldFacei]);
220  }
221  else
222  {
223  patchFld[i] = oldFld[oldFacei];
224  }
225  }
226  else
227  {
228  patchFld[i] = oldFld[oldFacei];
229  }
230  }
231  }
232 
233  fieldi++;
234  }
235 }
236 
237 
238 template<class GeoField>
239 void Foam::fvMeshDistribute::correctCoupledPatchFields()
240 {
241  HashTable<GeoField*> flds
242  (
243  mesh_.objectRegistry::lookupClass<GeoField>()
244  );
245 
246  forAllIter(typename HashTable<GeoField*>, flds, iter)
247  {
248  GeoField& fld = *iter();
249 
250  typename GeoField::Boundary& bfld = fld.boundaryFieldRef();
251 
252  if
253  (
256  )
257  {
258  label nReq = Pstream::nRequests();
259 
260  forAll(bfld, patchi)
261  {
262  if (bfld[patchi].coupled())
263  {
264  bfld[patchi].initEvaluate(Pstream::defaultCommsType);
265  }
266  }
267 
268  // Block for any outstanding requests
269  if
270  (
273  )
274  {
275  Pstream::waitRequests(nReq);
276  }
277 
278  forAll(bfld, patchi)
279  {
280  if (bfld[patchi].coupled())
281  {
282  bfld[patchi].evaluate(Pstream::defaultCommsType);
283  }
284  }
285  }
287  {
288  const lduSchedule& patchSchedule =
289  mesh_.globalData().patchSchedule();
290 
291  forAll(patchSchedule, patchEvali)
292  {
293  if (bfld[patchEvali].coupled())
294  {
295  if (patchSchedule[patchEvali].init)
296  {
297  bfld[patchSchedule[patchEvali].patch]
298  .initEvaluate(Pstream::commsTypes::scheduled);
299  }
300  else
301  {
302  bfld[patchSchedule[patchEvali].patch]
304  }
305  }
306  }
307  }
308  }
309 }
310 
311 
312 template<class GeoField>
313 void Foam::fvMeshDistribute::sendFields
314 (
315  const label domain,
316  const wordList& fieldNames,
317  const fvMeshSubset& subsetter,
318  Ostream& toNbr
319 )
320 {
321  // Send fields. Note order supplied so we can receive in exactly the same
322  // order.
323  // Note that field gets written as entry in dictionary so we
324  // can construct from subdictionary.
325  // (since otherwise the reading as-a-dictionary mixes up entries from
326  // consecutive fields)
327  // The dictionary constructed is:
328  // volScalarField
329  // {
330  // p {internalField ..; boundaryField ..;}
331  // k {internalField ..; boundaryField ..;}
332  // }
333  // volVectorField
334  // {
335  // U {internalField ... }
336  // }
337 
338  // volVectorField {U {internalField ..; boundaryField ..;}}
339 
340  toNbr << GeoField::typeName << token::NL << token::BEGIN_BLOCK << token::NL;
341  forAll(fieldNames, i)
342  {
343  if (debug)
344  {
345  Pout<< "Subsetting field " << fieldNames[i]
346  << " for domain:" << domain << endl;
347  }
348 
349  // Send all fieldNames. This has to be exactly the same set as is
350  // being received!
351  const GeoField& fld =
352  subsetter.baseMesh().lookupObject<GeoField>(fieldNames[i]);
353 
354  tmp<GeoField> tsubfld = subsetter.interpolate(fld);
355 
356  toNbr
358  << tsubfld
360  }
361  toNbr << token::END_BLOCK << token::NL;
362 }
363 
364 
365 template<class GeoField>
366 void Foam::fvMeshDistribute::receiveFields
367 (
368  const label domain,
369  const wordList& fieldNames,
370  typename GeoField::Mesh& mesh,
371  PtrList<GeoField>& fields,
372  const dictionary& fieldDicts
373 )
374 {
375  if (debug)
376  {
377  Pout<< "Receiving fields " << fieldNames
378  << " from domain:" << domain << endl;
379  }
380 
381  fields.setSize(fieldNames.size());
382 
383  forAll(fieldNames, i)
384  {
385  if (debug)
386  {
387  Pout<< "Constructing field " << fieldNames[i]
388  << " from domain:" << domain << endl;
389  }
390 
391  fields.set
392  (
393  i,
394  new GeoField
395  (
396  IOobject
397  (
398  fieldNames[i],
399  mesh.thisDb().time().name(),
400  mesh.thisDb(),
403  ),
404  mesh,
405  fieldDicts.subDict(fieldNames[i])
406  )
407  );
408  }
409 }
410 
411 
412 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:459
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Generic field type.
Definition: FieldField.H:77
Generic GeometricField class.
GeometricBoundaryField< Type, PatchField, GeoMesh > Boundary
Type of the boundary field.
An STL-conforming hash table.
Definition: HashTable.H:127
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
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
static void printFieldInfo(const fvMesh &)
Print some field info.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:101
HashTable< const Type * > lookupClass(const bool strict=false) const
Lookup and return all objects of the given Type.
@ BEGIN_BLOCK
Definition: token.H:110
@ END_BLOCK
Definition: token.H:111
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
label patchi
static List< word > fieldNames
Definition: globalFoam.H:46
gmvFile<< "tracers "<< particles.size()<< nl;{ pointField positions(particles.size());label particlei=0;forAllConstIter(Cloud< passiveParticle >, particles, iter) { positions[particlei++]=iter().position(mesh);} for(i=0;i< pTraits< point >::nComponents;i++) { forAll(positions, particlei) { gmvFile<< component(positions[particlei], i)<< ' ';} gmvFile<< nl;}}forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
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:230
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:251
errorManip< error > abort(error &err)
Definition: errorManip.H:131
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
error FatalError