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-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 "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.primitiveField());
50 
51  // Modify internal field
52  Field<Type>& intFld = fld.primitiveFieldRef();
53 
54  intFld.setSize(mesh.nCells());
55 
56  intFld.rmap(oldInternalField, meshMap.oldCellMap());
57  intFld.rmap(fldToAdd.primitiveField(), meshMap.addedCellMap());
58  }
59 
60 
61  // Patch fields from old mesh
62  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
63 
64  typename GeometricField<Type, fvPatchField, volMesh>::
65  Boundary& bfld = fld.boundaryFieldRef();
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, nullptr);
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(), // 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(), // 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  if (debug)
284  {
285  Pout<< "MapVolFields : Storing old time for " << fieldIter()->name()
286  << endl;
287  }
288 
289  const_cast<GeometricField<Type, fvPatchField, volMesh>*>(fieldIter())
290  ->storeOldTimes();
291  }
292 
293 
294  for
295  (
297  iterator fieldIter = fields.begin();
298  fieldIter != fields.end();
299  ++fieldIter
300  )
301  {
304  (
305  *fieldIter()
306  );
307 
308  if (fieldsToAdd.found(fld.name()))
309  {
311  *fieldsToAdd[fld.name()];
312 
313  if (debug)
314  {
315  Pout<< "MapVolFields : mapping " << fld.name()
316  << " and " << fldToAdd.name() << endl;
317  }
318 
319  MapVolField<Type>(meshMap, fld, fldToAdd);
320  }
321  else
322  {
324  << "Not mapping field " << fld.name()
325  << " since not present on mesh to add"
326  << endl;
327  }
328  }
329 }
330 
331 
332 template<class Type>
333 void Foam::fvMeshAdder::MapSurfaceField
334 (
335  const mapAddedPolyMesh& meshMap,
336 
339 )
340 {
341  const fvMesh& mesh = fld.mesh();
342  const labelList& oldPatchStarts = meshMap.oldPatchStarts();
343 
345  Boundary& bfld = fld.boundaryFieldRef();
346 
347  // Internal field
348  // ~~~~~~~~~~~~~~
349 
350  // Store old internal field
351  {
352  Field<Type> oldField(fld);
353 
354  // Modify internal field
355  Field<Type>& intFld = fld.primitiveFieldRef();
356 
357  intFld.setSize(mesh.nInternalFaces());
358 
359  intFld.rmap(oldField, meshMap.oldFaceMap());
360  intFld.rmap(fldToAdd, meshMap.addedFaceMap());
361 
362 
363  // Faces that were boundary faces but are not anymore.
364  // Use owner value (so lowest numbered cell, i.e. from 'old' not 'added'
365  // mesh)
366  forAll(bfld, patchi)
367  {
368  const fvsPatchField<Type>& pf = bfld[patchi];
369 
370  label start = oldPatchStarts[patchi];
371 
372  forAll(pf, i)
373  {
374  label newFacei = meshMap.oldFaceMap()[start + i];
375 
376  if (newFacei >= 0 && newFacei < mesh.nInternalFaces())
377  {
378  intFld[newFacei] = pf[i];
379  }
380  }
381  }
382  }
383 
384 
385  // Patch fields from old mesh
386  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
387 
388  {
389  const labelList& oldPatchMap = meshMap.oldPatchMap();
390  const labelList& oldPatchSizes = meshMap.oldPatchSizes();
391 
392  // Reorder old patches in order of new ones. Put removed patches at end.
393 
394  label unusedPatchi = 0;
395 
396  forAll(oldPatchMap, patchi)
397  {
398  label newPatchi = oldPatchMap[patchi];
399 
400  if (newPatchi != -1)
401  {
402  unusedPatchi++;
403  }
404  }
405 
406  label nUsedPatches = unusedPatchi;
407 
408  // Reorder list for patchFields
409  labelList oldToNew(oldPatchMap.size());
410 
411  forAll(oldPatchMap, patchi)
412  {
413  label newPatchi = oldPatchMap[patchi];
414 
415  if (newPatchi != -1)
416  {
417  oldToNew[patchi] = newPatchi;
418  }
419  else
420  {
421  oldToNew[patchi] = unusedPatchi++;
422  }
423  }
424 
425 
426  // Sort deleted ones last so is now in newPatch ordering
427  bfld.reorder(oldToNew);
428  // Extend to covers all patches
429  bfld.setSize(mesh.boundaryMesh().size());
430  // Delete unused patches
431  for
432  (
433  label newPatchi = nUsedPatches;
434  newPatchi < bfld.size();
435  newPatchi++
436  )
437  {
438  bfld.set(newPatchi, nullptr);
439  }
440 
441 
442  // Map old values
443  // ~~~~~~~~~~~~~~
444 
445  forAll(oldPatchMap, patchi)
446  {
447  label newPatchi = oldPatchMap[patchi];
448 
449  if (newPatchi != -1)
450  {
451  labelList newToOld
452  (
453  calcPatchMap
454  (
455  oldPatchStarts[patchi],
456  oldPatchSizes[patchi],
457  meshMap.oldFaceMap(),
458  mesh.boundaryMesh()[newPatchi],
459  -1 // unmapped value
460  )
461  );
462 
463  directFvPatchFieldMapper patchMapper(newToOld);
464 
465  // Create new patchField with same type as existing one.
466  // Note:
467  // - boundaryField already in new order so access with newPatchi
468  // - bfld[newPatchi] both used for type and old
469  // value
470  // - hope that field mapping allows aliasing since old and new
471  // are same memory!
472  bfld.set
473  (
474  newPatchi,
476  (
477  bfld[newPatchi], // old field
478  mesh.boundary()[newPatchi], // new fvPatch
479  fld(), // new internal field
480  patchMapper // mapper (new to old)
481  )
482  );
483  }
484  }
485  }
486 
487 
488 
489  // Patch fields from added mesh
490  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
491 
492  {
493  const labelList& addedPatchMap = meshMap.addedPatchMap();
494 
495  // Add addedMesh patches
496  forAll(addedPatchMap, patchi)
497  {
498  label newPatchi = addedPatchMap[patchi];
499 
500  if (newPatchi != -1)
501  {
502  const polyPatch& newPatch = mesh.boundaryMesh()[newPatchi];
503  const polyPatch& oldPatch =
504  fldToAdd.mesh().boundaryMesh()[patchi];
505 
506  if (!bfld(newPatchi))
507  {
508  // First occurrence of newPatchi. Map from existing
509  // patchField
510 
511  // From new patch faces to patch faces on added mesh.
512  labelList newToAdded
513  (
514  calcPatchMap
515  (
516  oldPatch.start(),
517  oldPatch.size(),
518  meshMap.addedFaceMap(),
519  newPatch,
520  -1 // unmapped values
521  )
522  );
523 
524  directFvPatchFieldMapper patchMapper(newToAdded);
525 
526  bfld.set
527  (
528  newPatchi,
530  (
531  fldToAdd.boundaryField()[patchi],// added field
532  mesh.boundary()[newPatchi], // new fvPatch
533  fld(), // new int. field
534  patchMapper // mapper
535  )
536  );
537  }
538  else
539  {
540  // PatchField will have correct size already. Just slot in
541  // my elements.
542 
543  labelList addedToNew(oldPatch.size(), -1);
544  forAll(addedToNew, i)
545  {
546  label addedFacei = oldPatch.start()+i;
547  label newFacei = meshMap.addedFaceMap()[addedFacei];
548  label patchFacei = newFacei-newPatch.start();
549  if (patchFacei >= 0 && patchFacei < newPatch.size())
550  {
551  addedToNew[i] = patchFacei;
552  }
553  }
554 
555  bfld[newPatchi].rmap
556  (
557  fldToAdd.boundaryField()[patchi],
558  addedToNew
559  );
560  }
561  }
562  }
563  }
564 }
565 
566 
567 template<class Type>
569 (
570  const mapAddedPolyMesh& meshMap,
571  const fvMesh& mesh,
572  const fvMesh& meshToAdd
573 )
574 {
576 
578  (
579  mesh.objectRegistry::lookupClass<fldType>()
580  );
581 
582  HashTable<const fldType*> fieldsToAdd
583  (
584  meshToAdd.objectRegistry::lookupClass<fldType>()
585  );
586 
587  // It is necessary to enforce that all old-time fields are stored
588  // before the mapping is performed. Otherwise, if the
589  // old-time-level field is mapped before the field itself, sizes
590  // will not match.
591 
592  for
593  (
594  typename HashTable<const fldType*>::
595  iterator fieldIter = fields.begin();
596  fieldIter != fields.end();
597  ++fieldIter
598  )
599  {
600  if (debug)
601  {
602  Pout<< "MapSurfaceFields : Storing old time for "
603  << fieldIter()->name() << endl;
604  }
605 
606  const_cast<fldType*>(fieldIter())->storeOldTimes();
607  }
608 
609 
610  for
611  (
612  typename HashTable<const fldType*>::
613  iterator fieldIter = fields.begin();
614  fieldIter != fields.end();
615  ++fieldIter
616  )
617  {
618  fldType& fld = const_cast<fldType&>(*fieldIter());
619 
620  if (fieldsToAdd.found(fld.name()))
621  {
622  const fldType& fldToAdd = *fieldsToAdd[fld.name()];
623 
624  if (debug)
625  {
626  Pout<< "MapSurfaceFields : mapping " << fld.name()
627  << " and " << fldToAdd.name() << endl;
628  }
629 
630  MapSurfaceField<Type>(meshMap, fld, fldToAdd);
631  }
632  else
633  {
635  << "Not mapping field " << fld.name()
636  << " since not present on mesh to add"
637  << endl;
638  }
639  }
640 }
641 
642 
643 template<class Type>
644 void Foam::fvMeshAdder::MapDimField
645 (
646  const mapAddedPolyMesh& meshMap,
647 
649  const DimensionedField<Type, volMesh>& fldToAdd
650 )
651 {
652  const fvMesh& mesh = fld.mesh();
653 
654  // Store old field
655  Field<Type> oldField(fld);
656 
657  fld.setSize(mesh.nCells());
658 
659  fld.rmap(oldField, meshMap.oldCellMap());
660  fld.rmap(fldToAdd, meshMap.addedCellMap());
661 }
662 
663 
664 template<class Type>
666 (
667  const mapAddedPolyMesh& meshMap,
668  const fvMesh& mesh,
669  const fvMesh& meshToAdd
670 )
671 {
672  typedef DimensionedField<Type, volMesh> fldType;
673 
674  // Note: use strict flag on lookupClass to avoid picking up
675  // volFields
677  (
678  mesh.objectRegistry::lookupClass<fldType>(true)
679  );
680 
681  HashTable<const fldType*> fieldsToAdd
682  (
683  meshToAdd.objectRegistry::lookupClass<fldType>(true)
684  );
685 
686  for
687  (
688  typename HashTable<const fldType*>::
689  iterator fieldIter = fields.begin();
690  fieldIter != fields.end();
691  ++fieldIter
692  )
693  {
694  fldType& fld = const_cast<fldType&>(*fieldIter());
695 
696  if (fieldsToAdd.found(fld.name()))
697  {
698  const fldType& fldToAdd = *fieldsToAdd[fld.name()];
699 
700  if (debug)
701  {
702  Pout<< "MapDimFields : mapping " << fld.name()
703  << " and " << fldToAdd.name() << endl;
704  }
705 
706  MapDimField<Type>(meshMap, fld, fldToAdd);
707  }
708  else
709  {
710  WarningIn("fvMeshAdder::MapDimFields(..)")
711  << "Not mapping field " << fld.name()
712  << " since not present on mesh to add"
713  << endl;
714  }
715  }
716 }
717 
718 
719 // ************************************************************************* //
Foam::surfaceFields.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:88
#define WarningIn(functionName)
Report a warning using Foam::Warning.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const labelList & oldPatchSizes() const
Return list of the old patch sizes.
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
const word & name() const
Return name.
Definition: IOobject.H:291
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:110
label nInternalFaces() const
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:163
label nCells() const
const labelList & oldCellMap() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
const labelList & addedPatchMap() const
From added mesh patch index to new patch index or -1 if.
const labelList & addedCellMap() const
Class containing mesh-to-mesh mapping information after a mesh addition where we add a mesh (&#39;added m...
const labelList & oldPatchStarts() const
Return list of the old patch start labels.
Generic GeometricField class.
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:104
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 found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:113
Pre-declare SubField and related Field type.
Definition: Field.H:57
word name() const
Return file name (part beyond last /)
Definition: fileName.C:180
const labelList & oldPatchMap() const
From old patch index to new patch index or -1 if patch.
const labelList & oldFaceMap() const
static tmp< fvPatchField< Type > > New(const word &, const fvPatch &, const DimensionedField< Type, volMesh > &)
Return a pointer to a new patchField created on freestore given.
List< label > labelList
A List of labels.
Definition: labelList.H:56
iterator begin()
Iterator set to the beginning of the HashTable.
Definition: HashTableI.H:419
An STL-conforming hash table.
Definition: HashTable.H:62
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
const Mesh & mesh() const
Return mesh.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
void setSize(const label)
Reset size of List.
Definition: List.C:281
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:584
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static void MapSurfaceFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all surfaceFields of Type.
static void MapDimFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all DimensionedFields of Type.
static void MapVolFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all volFields of Type.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
const labelList & addedFaceMap() const
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:545