fvMeshDistributeTemplates.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 | Copyright (C) 2011-2016 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 "mapPolyMesh.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 template<class GeoField>
32 {
34  (
35  mesh.objectRegistry::lookupClass<GeoField>()
36  );
37 
38  forAllConstIter(typename HashTable<const GeoField*>, flds, iter)
39  {
40  const GeoField& fld = *iter();
41 
42  Pout<< "Field:" << iter.key() << " internalsize:" << fld.size()
43  //<< " value:" << fld
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 T, class Mesh>
59 void Foam::fvMeshDistribute::saveBoundaryFields
60 (
62 ) const
63 {
64  // Save whole boundary field
65 
67 
69  (
70  static_cast<const fvMesh&>(mesh_).objectRegistry::lookupClass<fldType>()
71  );
72 
73  bflds.setSize(flds.size());
74 
75  label i = 0;
76 
77  forAllConstIter(typename HashTable<const fldType*>, flds, iter)
78  {
79  const fldType& fld = *iter();
80 
81  bflds.set(i, fld.boundaryField().clone().ptr());
82 
83  i++;
84  }
85 }
86 
87 
88 template<class T, class Mesh>
89 void Foam::fvMeshDistribute::mapBoundaryFields
90 (
91  const mapPolyMesh& map,
92  const PtrList<FieldField<fvsPatchField, T>>& oldBflds
93 )
94 {
95  // Map boundary field
96 
97  const labelList& oldPatchStarts = map.oldPatchStarts();
98  const labelList& faceMap = map.faceMap();
99 
101 
103  (
104  mesh_.objectRegistry::lookupClass<fldType>()
105  );
106 
107  if (flds.size() != oldBflds.size())
108  {
110  << abort(FatalError);
111  }
112 
113  label fieldi = 0;
114 
115  forAllIter(typename HashTable<fldType*>, flds, iter)
116  {
117  fldType& fld = *iter();
118  typename fldType::Boundary& bfld =
119  fld.boundaryFieldRef();
120 
121  const FieldField<fvsPatchField, T>& oldBfld = oldBflds[fieldi++];
122 
123  // Pull from old boundary field into bfld.
124 
125  forAll(bfld, patchi)
126  {
127  fvsPatchField<T>& patchFld = bfld[patchi];
128  label facei = patchFld.patch().start();
129 
130  forAll(patchFld, i)
131  {
132  label oldFacei = faceMap[facei++];
133 
134  // Find patch and local patch face oldFacei was in.
135  forAll(oldPatchStarts, oldPatchi)
136  {
137  label oldLocalI = oldFacei - oldPatchStarts[oldPatchi];
138 
139  if (oldLocalI >= 0 && oldLocalI < oldBfld[oldPatchi].size())
140  {
141  patchFld[i] = oldBfld[oldPatchi][oldLocalI];
142  }
143  }
144  }
145  }
146  }
147 }
148 
149 
150 template<class T>
151 void Foam::fvMeshDistribute::saveInternalFields
152 (
153  PtrList<Field<T> >& iflds
154 ) const
155 {
157 
159  (
160  static_cast<const fvMesh&>(mesh_).objectRegistry::lookupClass<fldType>()
161  );
162 
163  iflds.setSize(flds.size());
164 
165  label i = 0;
166 
167  forAllConstIter(typename HashTable<const fldType*>, flds, iter)
168  {
169  const fldType& fld = *iter();
170 
171  iflds.set(i, fld.primitiveField().clone());
172 
173  i++;
174  }
175 }
176 
177 
178 template<class T>
179 void Foam::fvMeshDistribute::mapExposedFaces
180 (
181  const mapPolyMesh& map,
182  const PtrList<Field<T> >& oldFlds
183 )
184 {
185  // Set boundary values of exposed internal faces
186 
187  const labelList& faceMap = map.faceMap();
188 
190 
192  (
193  mesh_.objectRegistry::lookupClass<fldType>()
194  );
195 
196  if (flds.size() != oldFlds.size())
197  {
198  FatalErrorIn("fvMeshDistribute::mapExposedFaces(..)") << "problem"
199  << abort(FatalError);
200  }
201 
202 
203  label fieldI = 0;
204 
205  forAllIter(typename HashTable<fldType*>, flds, iter)
206  {
207  fldType& fld = *iter();
208  typename fldType::Boundary& bfld = fld.boundaryFieldRef();
209 
210  const Field<T>& oldInternal = oldFlds[fieldI++];
211 
212  // Pull from old internal field into bfld.
213 
214  forAll(bfld, patchi)
215  {
216  fvsPatchField<T>& patchFld = bfld[patchi];
217 
218  forAll(patchFld, i)
219  {
220  const label faceI = patchFld.patch().start()+i;
221 
222  label oldFaceI = faceMap[faceI];
223 
224  if (oldFaceI < oldInternal.size())
225  {
226  patchFld[i] = oldInternal[oldFaceI];
227 
228  if (map.flipFaceFlux().found(faceI))
229  {
230  patchFld[i] = flipOp()(patchFld[i]);
231  }
232  }
233  }
234  }
235  }
236 }
237 
238 
239 template<class GeoField, class PatchFieldType>
240 void Foam::fvMeshDistribute::initPatchFields
241 (
242  const typename GeoField::value_type& initVal
243 )
244 {
245  // Init patch fields of certain type
246 
248  (
249  mesh_.objectRegistry::lookupClass<GeoField>()
250  );
251 
252  forAllIter(typename HashTable<GeoField*>, flds, iter)
253  {
254  GeoField& fld = *iter();
255 
256  typename GeoField::Boundary& bfld = fld.boundaryFieldRef();
257 
258  forAll(bfld, patchi)
259  {
260  if (isA<PatchFieldType>(bfld[patchi]))
261  {
262  bfld[patchi] == initVal;
263  }
264  }
265  }
266 }
267 
268 
269 template<class GeoField>
270 void Foam::fvMeshDistribute::correctBoundaryConditions()
271 {
272  // correctBoundaryConditions patch fields of certain type
273 
275  (
276  mesh_.objectRegistry::lookupClass<GeoField>()
277  );
278 
279  forAllIter(typename HashTable<GeoField*>, flds, iter)
280  {
281  const GeoField& fld = *iter();
282  fld.correctBoundaryConditions();
283  }
284 }
285 
286 
287 template<class GeoField>
288 void Foam::fvMeshDistribute::sendFields
289 (
290  const label domain,
291  const wordList& fieldNames,
292  const fvMeshSubset& subsetter,
293  Ostream& toNbr
294 )
295 {
296  // Send fields. Note order supplied so we can receive in exactly the same
297  // order.
298  // Note that field gets written as entry in dictionary so we
299  // can construct from subdictionary.
300  // (since otherwise the reading as-a-dictionary mixes up entries from
301  // consecutive fields)
302  // The dictionary constructed is:
303  // volScalarField
304  // {
305  // p {internalField ..; boundaryField ..;}
306  // k {internalField ..; boundaryField ..;}
307  // }
308  // volVectorField
309  // {
310  // U {internalField ... }
311  // }
312 
313  // volVectorField {U {internalField ..; boundaryField ..;}}
314 
315  toNbr << GeoField::typeName << token::NL << token::BEGIN_BLOCK << token::NL;
316  forAll(fieldNames, i)
317  {
318  if (debug)
319  {
320  Pout<< "Subsetting field " << fieldNames[i]
321  << " for domain:" << domain << endl;
322  }
323 
324  // Send all fieldNames. This has to be exactly the same set as is
325  // being received!
326  const GeoField& fld =
327  subsetter.baseMesh().lookupObject<GeoField>(fieldNames[i]);
328 
329  tmp<GeoField> tsubfld = subsetter.interpolate(fld);
330 
331  toNbr
332  << fieldNames[i] << token::NL << token::BEGIN_BLOCK
333  << tsubfld
334  << token::NL << token::END_BLOCK << token::NL;
335  }
336  toNbr << token::END_BLOCK << token::NL;
337 }
338 
339 
340 template<class GeoField>
341 void Foam::fvMeshDistribute::receiveFields
342 (
343  const label domain,
344  const wordList& fieldNames,
345  fvMesh& mesh,
346  PtrList<GeoField>& fields,
347  const dictionary& fieldDicts
348 )
349 {
350  if (debug)
351  {
352  Pout<< "Receiving fields " << fieldNames
353  << " from domain:" << domain << endl;
354  }
355 
356  fields.setSize(fieldNames.size());
357 
358  forAll(fieldNames, i)
359  {
360  if (debug)
361  {
362  Pout<< "Constructing field " << fieldNames[i]
363  << " from domain:" << domain << endl;
364  }
365 
366  fields.set
367  (
368  i,
369  new GeoField
370  (
371  IOobject
372  (
373  fieldNames[i],
374  mesh.time().timeName(),
375  mesh,
378  ),
379  mesh,
380  fieldDicts.subDict(fieldNames[i])
381  )
382  );
383  }
384 }
385 
386 
387 // ************************************************************************* //
Class containing functor to negate primitives. Dummy for all other types.
Definition: flipOp.H:50
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const labelHashSet & flipFaceFlux() const
Map of flipped face flux faces.
Definition: mapPolyMesh.H:556
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:633
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
const labelList & oldPatchStarts() const
Return list of the old patch start labels.
Definition: mapPolyMesh.H:626
Generic GeometricField class.
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
Generic field type.
Definition: FieldField.H:51
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
dynamicFvMesh & mesh
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
bool set(const label) const
Is element set.
Definition: PtrListI.H:65
An STL-conforming hash table.
Definition: HashTable.H:61
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:109
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
static void printFieldInfo(const fvMesh &)
Print some field info.
Post-processing mesh subset tool. Given the original mesh and the list of selected cells...
Definition: fvMeshSubset.H:73
const fvPatch & patch() const
Return patch.
label patchi
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const fvMesh & baseMesh() const
Original mesh.
Definition: fvMeshSubset.H:217
const labelList & faceMap() const
Old face map.
Definition: mapPolyMesh.H:404
A class for managing temporary objects.
Definition: PtrList.H:54
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
static tmp< GeometricField< Type, fvPatchField, volMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &, const fvMesh &sMesh, const labelList &patchMap, const labelList &cellMap, const labelList &faceMap)
Map volume field.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243