fvMeshAdderTemplates.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-2014 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 "volFields.H"
27 #include "surfaceFields.H"
28 #include "emptyFvPatchField.H"
30 
31 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
32 
33 template<class Type>
34 void Foam::fvMeshAdder::MapVolField
35 (
36  const mapAddedPolyMesh& meshMap,
37 
38  GeometricField<Type, fvPatchField, volMesh>& fld,
39  const GeometricField<Type, fvPatchField, volMesh>& fldToAdd
40 )
41 {
42  const fvMesh& mesh = fld.mesh();
43 
44  // Internal field
45  // ~~~~~~~~~~~~~~
46 
47  {
48  // Store old internal field
49  Field<Type> oldInternalField(fld.internalField());
50 
51  // Modify internal field
52  Field<Type>& intFld = fld.internalField();
53 
54  intFld.setSize(mesh.nCells());
55 
56  intFld.rmap(oldInternalField, meshMap.oldCellMap());
57  intFld.rmap(fldToAdd.internalField(), meshMap.addedCellMap());
58  }
59 
60 
61  // Patch fields from old mesh
62  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
63 
64  typename GeometricField<Type, fvPatchField, volMesh>::
65  GeometricBoundaryField& bfld = fld.boundaryField();
66 
67  {
68  const labelList& oldPatchMap = meshMap.oldPatchMap();
69  const labelList& oldPatchStarts = meshMap.oldPatchStarts();
70  const labelList& oldPatchSizes = meshMap.oldPatchSizes();
71 
72  // Reorder old patches in order of new ones. Put removed patches at end.
73 
74  label unusedPatchI = 0;
75 
76  forAll(oldPatchMap, patchI)
77  {
78  label newPatchI = oldPatchMap[patchI];
79 
80  if (newPatchI != -1)
81  {
82  unusedPatchI++;
83  }
84  }
85 
86  label nUsedPatches = unusedPatchI;
87 
88  // Reorder list for patchFields
89  labelList oldToNew(oldPatchMap.size());
90 
91  forAll(oldPatchMap, patchI)
92  {
93  label newPatchI = oldPatchMap[patchI];
94 
95  if (newPatchI != -1)
96  {
97  oldToNew[patchI] = newPatchI;
98  }
99  else
100  {
101  oldToNew[patchI] = unusedPatchI++;
102  }
103  }
104 
105 
106  // Sort deleted ones last so is now in newPatch ordering
107  bfld.reorder(oldToNew);
108  // Extend to covers all patches
109  bfld.setSize(mesh.boundaryMesh().size());
110  // Delete unused patches
111  for
112  (
113  label newPatchI = nUsedPatches;
114  newPatchI < bfld.size();
115  newPatchI++
116  )
117  {
118  bfld.set(newPatchI, NULL);
119  }
120 
121 
122  // Map old values
123  // ~~~~~~~~~~~~~~
124 
125  forAll(oldPatchMap, patchI)
126  {
127  label newPatchI = oldPatchMap[patchI];
128 
129  if (newPatchI != -1)
130  {
131  labelList newToOld
132  (
133  calcPatchMap
134  (
135  oldPatchStarts[patchI],
136  oldPatchSizes[patchI],
137  meshMap.oldFaceMap(),
138  mesh.boundaryMesh()[newPatchI],
139  -1 // unmapped value
140  )
141  );
142 
143  directFvPatchFieldMapper patchMapper(newToOld);
144 
145 
146  // Create new patchField with same type as existing one.
147  // Note:
148  // - boundaryField already in new order so access with newPatchI
149  // - fld.boundaryField()[newPatchI] both used for type and old
150  // value
151  // - hope that field mapping allows aliasing since old and new
152  // are same memory!
153  bfld.set
154  (
155  newPatchI,
157  (
158  bfld[newPatchI], // old field
159  mesh.boundary()[newPatchI], // new fvPatch
160  fld.dimensionedInternalField(), // new internal field
161  patchMapper // mapper (new to old)
162  )
163  );
164  }
165  }
166  }
167 
168 
169 
170  // Patch fields from added mesh
171  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
172 
173  {
174  const labelList& addedPatchMap = meshMap.addedPatchMap();
175 
176  // Add addedMesh patches
177  forAll(addedPatchMap, patchI)
178  {
179  label newPatchI = addedPatchMap[patchI];
180 
181  if (newPatchI != -1)
182  {
183  const polyPatch& newPatch = mesh.boundaryMesh()[newPatchI];
184  const polyPatch& oldPatch =
185  fldToAdd.mesh().boundaryMesh()[patchI];
186 
187  if (!bfld(newPatchI))
188  {
189  // First occurrence of newPatchI. Map from existing
190  // patchField
191 
192  // From new patch faces to patch faces on added mesh.
193  labelList newToAdded
194  (
195  calcPatchMap
196  (
197  oldPatch.start(),
198  oldPatch.size(),
199  meshMap.addedFaceMap(),
200  newPatch,
201  -1 // unmapped values
202  )
203  );
204 
205  directFvPatchFieldMapper patchMapper(newToAdded);
206 
207  bfld.set
208  (
209  newPatchI,
211  (
212  fldToAdd.boundaryField()[patchI], // added field
213  mesh.boundary()[newPatchI], // new fvPatch
214  fld.dimensionedInternalField(), // new int. field
215  patchMapper // mapper
216  )
217  );
218  }
219  else
220  {
221  // PatchField will have correct size already. Just slot in
222  // my elements.
223 
224  labelList addedToNew(oldPatch.size(), -1);
225  forAll(addedToNew, i)
226  {
227  label addedFaceI = oldPatch.start()+i;
228  label newFaceI = meshMap.addedFaceMap()[addedFaceI];
229  label patchFaceI = newFaceI-newPatch.start();
230  if (patchFaceI >= 0 && patchFaceI < newPatch.size())
231  {
232  addedToNew[i] = patchFaceI;
233  }
234  }
235 
236  bfld[newPatchI].rmap
237  (
238  fldToAdd.boundaryField()[patchI],
239  addedToNew
240  );
241  }
242  }
243  }
244  }
245 }
246 
247 
248 template<class Type>
250 (
251  const mapAddedPolyMesh& meshMap,
252  const fvMesh& mesh,
253  const fvMesh& meshToAdd
254 )
255 {
257  (
258  mesh.objectRegistry::lookupClass
260  ()
261  );
262 
264  (
265  meshToAdd.objectRegistry::lookupClass
267  ()
268  );
269 
270  // It is necessary to enforce that all old-time fields are stored
271  // before the mapping is performed. Otherwise, if the
272  // old-time-level field is mapped before the field itself, sizes
273  // will not match.
274 
275  for
276  (
278  iterator fieldIter = fields.begin();
279  fieldIter != fields.end();
280  ++fieldIter
281  )
282  {
283  const_cast<GeometricField<Type, fvPatchField, volMesh>*>(fieldIter())
284  ->storeOldTimes();
285  }
286 
287 
288  for
289  (
291  iterator fieldIter = fields.begin();
292  fieldIter != fields.end();
293  ++fieldIter
294  )
295  {
298  (
299  *fieldIter()
300  );
301 
302  if (fieldsToAdd.found(fld.name()))
303  {
305  *fieldsToAdd[fld.name()];
306 
307  MapVolField<Type>(meshMap, fld, fldToAdd);
308  }
309  else
310  {
311  WarningIn("fvMeshAdder::MapVolFields(..)")
312  << "Not mapping field " << fld.name()
313  << " since not present on mesh to add"
314  << endl;
315  }
316  }
317 }
318 
319 
320 template<class Type>
321 void Foam::fvMeshAdder::MapSurfaceField
322 (
323  const mapAddedPolyMesh& meshMap,
324 
327 )
328 {
329  const fvMesh& mesh = fld.mesh();
330  const labelList& oldPatchStarts = meshMap.oldPatchStarts();
331 
333  GeometricBoundaryField& bfld = fld.boundaryField();
334 
335  // Internal field
336  // ~~~~~~~~~~~~~~
337 
338  // Store old internal field
339  {
340  Field<Type> oldField(fld);
341 
342  // Modify internal field
343  Field<Type>& intFld = fld.internalField();
344 
345  intFld.setSize(mesh.nInternalFaces());
346 
347  intFld.rmap(oldField, meshMap.oldFaceMap());
348  intFld.rmap(fldToAdd, meshMap.addedFaceMap());
349 
350 
351  // Faces that were boundary faces but are not anymore.
352  // Use owner value (so lowest numbered cell, i.e. from 'old' not 'added'
353  // mesh)
354  forAll(bfld, patchI)
355  {
356  const fvsPatchField<Type>& pf = bfld[patchI];
357 
358  label start = oldPatchStarts[patchI];
359 
360  forAll(pf, i)
361  {
362  label newFaceI = meshMap.oldFaceMap()[start + i];
363 
364  if (newFaceI >= 0 && newFaceI < mesh.nInternalFaces())
365  {
366  intFld[newFaceI] = pf[i];
367  }
368  }
369  }
370  }
371 
372 
373  // Patch fields from old mesh
374  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
375 
376  {
377  const labelList& oldPatchMap = meshMap.oldPatchMap();
378  const labelList& oldPatchSizes = meshMap.oldPatchSizes();
379 
380  // Reorder old patches in order of new ones. Put removed patches at end.
381 
382  label unusedPatchI = 0;
383 
384  forAll(oldPatchMap, patchI)
385  {
386  label newPatchI = oldPatchMap[patchI];
387 
388  if (newPatchI != -1)
389  {
390  unusedPatchI++;
391  }
392  }
393 
394  label nUsedPatches = unusedPatchI;
395 
396  // Reorder list for patchFields
397  labelList oldToNew(oldPatchMap.size());
398 
399  forAll(oldPatchMap, patchI)
400  {
401  label newPatchI = oldPatchMap[patchI];
402 
403  if (newPatchI != -1)
404  {
405  oldToNew[patchI] = newPatchI;
406  }
407  else
408  {
409  oldToNew[patchI] = unusedPatchI++;
410  }
411  }
412 
413 
414  // Sort deleted ones last so is now in newPatch ordering
415  bfld.reorder(oldToNew);
416  // Extend to covers all patches
417  bfld.setSize(mesh.boundaryMesh().size());
418  // Delete unused patches
419  for
420  (
421  label newPatchI = nUsedPatches;
422  newPatchI < bfld.size();
423  newPatchI++
424  )
425  {
426  bfld.set(newPatchI, NULL);
427  }
428 
429 
430  // Map old values
431  // ~~~~~~~~~~~~~~
432 
433  forAll(oldPatchMap, patchI)
434  {
435  label newPatchI = oldPatchMap[patchI];
436 
437  if (newPatchI != -1)
438  {
439  labelList newToOld
440  (
441  calcPatchMap
442  (
443  oldPatchStarts[patchI],
444  oldPatchSizes[patchI],
445  meshMap.oldFaceMap(),
446  mesh.boundaryMesh()[newPatchI],
447  -1 // unmapped value
448  )
449  );
450 
451  directFvPatchFieldMapper patchMapper(newToOld);
452 
453  // Create new patchField with same type as existing one.
454  // Note:
455  // - boundaryField already in new order so access with newPatchI
456  // - bfld[newPatchI] both used for type and old
457  // value
458  // - hope that field mapping allows aliasing since old and new
459  // are same memory!
460  bfld.set
461  (
462  newPatchI,
464  (
465  bfld[newPatchI], // old field
466  mesh.boundary()[newPatchI], // new fvPatch
467  fld.dimensionedInternalField(), // new internal field
468  patchMapper // mapper (new to old)
469  )
470  );
471  }
472  }
473  }
474 
475 
476 
477  // Patch fields from added mesh
478  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
479 
480  {
481  const labelList& addedPatchMap = meshMap.addedPatchMap();
482 
483  // Add addedMesh patches
484  forAll(addedPatchMap, patchI)
485  {
486  label newPatchI = addedPatchMap[patchI];
487 
488  if (newPatchI != -1)
489  {
490  const polyPatch& newPatch = mesh.boundaryMesh()[newPatchI];
491  const polyPatch& oldPatch =
492  fldToAdd.mesh().boundaryMesh()[patchI];
493 
494  if (!bfld(newPatchI))
495  {
496  // First occurrence of newPatchI. Map from existing
497  // patchField
498 
499  // From new patch faces to patch faces on added mesh.
500  labelList newToAdded
501  (
502  calcPatchMap
503  (
504  oldPatch.start(),
505  oldPatch.size(),
506  meshMap.addedFaceMap(),
507  newPatch,
508  -1 // unmapped values
509  )
510  );
511 
512  directFvPatchFieldMapper patchMapper(newToAdded);
513 
514  bfld.set
515  (
516  newPatchI,
518  (
519  fldToAdd.boundaryField()[patchI],// added field
520  mesh.boundary()[newPatchI], // new fvPatch
521  fld.dimensionedInternalField(), // new int. field
522  patchMapper // mapper
523  )
524  );
525  }
526  else
527  {
528  // PatchField will have correct size already. Just slot in
529  // my elements.
530 
531  labelList addedToNew(oldPatch.size(), -1);
532  forAll(addedToNew, i)
533  {
534  label addedFaceI = oldPatch.start()+i;
535  label newFaceI = meshMap.addedFaceMap()[addedFaceI];
536  label patchFaceI = newFaceI-newPatch.start();
537  if (patchFaceI >= 0 && patchFaceI < newPatch.size())
538  {
539  addedToNew[i] = patchFaceI;
540  }
541  }
542 
543  bfld[newPatchI].rmap
544  (
545  fldToAdd.boundaryField()[patchI],
546  addedToNew
547  );
548  }
549  }
550  }
551  }
552 }
553 
554 
555 template<class Type>
557 (
558  const mapAddedPolyMesh& meshMap,
559  const fvMesh& mesh,
560  const fvMesh& meshToAdd
561 )
562 {
564 
566  (
567  mesh.objectRegistry::lookupClass<fldType>()
568  );
569 
570  HashTable<const fldType*> fieldsToAdd
571  (
572  meshToAdd.objectRegistry::lookupClass<fldType>()
573  );
574 
575  // It is necessary to enforce that all old-time fields are stored
576  // before the mapping is performed. Otherwise, if the
577  // old-time-level field is mapped before the field itself, sizes
578  // will not match.
579 
580  for
581  (
582  typename HashTable<const fldType*>::
583  iterator fieldIter = fields.begin();
584  fieldIter != fields.end();
585  ++fieldIter
586  )
587  {
588  const_cast<fldType*>(fieldIter())
589  ->storeOldTimes();
590  }
591 
592 
593  for
594  (
595  typename HashTable<const fldType*>::
596  iterator fieldIter = fields.begin();
597  fieldIter != fields.end();
598  ++fieldIter
599  )
600  {
601  fldType& fld = const_cast<fldType&>(*fieldIter());
602 
603  if (fieldsToAdd.found(fld.name()))
604  {
605  const fldType& fldToAdd = *fieldsToAdd[fld.name()];
606 
607  MapSurfaceField<Type>(meshMap, fld, fldToAdd);
608  }
609  else
610  {
611  WarningIn("fvMeshAdder::MapSurfaceFields(..)")
612  << "Not mapping field " << fld.name()
613  << " since not present on mesh to add"
614  << endl;
615  }
616  }
617 }
618 
619 
620 // ************************************************************************* //
static void MapSurfaceFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all surfaceFields of Type.
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:142
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::surfaceFields.
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
bool set(const label) const
Is element set.
Definition: PtrListI.H:107
Info<< "Creating field dpdt\n"<< endl;volScalarField dpdt(IOobject("dpdt", runTime.timeName(), mesh), mesh, dimensionedScalar("dpdt", p.dimensions()/dimTime, 0));Info<< "Creating field kinetic energy K\n"<< endl;volScalarField K("K", 0.5 *magSqr(U));volScalarField p_rgh(IOobject("p_rgh", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);p_rgh=p-rho *gh;mesh.setFluxRequired(p_rgh.name());multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:127
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
const labelList & oldPatchSizes() const
Return list of the old patch sizes.
An STL-conforming hash table.
Definition: HashTable.H:61
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
InternalField & internalField()
Return internal field.
DimensionedInternalField & dimensionedInternalField()
Return dimensioned internal field.
const Mesh & mesh() const
Return mesh.
const labelList & addedFaceMap() const
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
static void MapVolFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all volFields of Type.
static tmp< fvPatchField< Type > > New(const word &, const fvPatch &, const DimensionedField< Type, volMesh > &)
Return a pointer to a new patchField created on freestore given.
const labelList & oldFaceMap() const
void setSize(const label)
Reset size of List.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define WarningIn(functionName)
Report a warning using Foam::Warning.
const labelList & oldPatchMap() const
From old patch index to new patch index or -1 if patch.
#define forAll(list, i)
Definition: UList.H:421
const labelList & oldPatchStarts() const
Return list of the old patch start labels.
Pre-declare SubField and related Field type.
Definition: Field.H:57
const word & name() const
Return name.
Definition: IOobject.H:260
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 ))
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
label nInternalFaces() const
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:505
Generic GeometricField class.
List< label > labelList
A List of labels.
Definition: labelList.H:56
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
void reorder(const labelUList &)
Reorders elements. Ordering does not have to be done in.
Definition: PtrList.C:208
iterator begin()
Iterator set to the beginning of the HashTable.
Definition: HashTableI.H:419
const labelList & addedPatchMap() const
From added mesh patch index to new patch index or -1 if.
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:109
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:552
Class containing mesh-to-mesh mapping information after a mesh addition where we add a mesh (&#39;added m...