fvMeshAdderTemplates.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 "volFields.H"
27 #include "surfaceFields.H"
28 #include "pointFields.H"
29 #include "forwardFieldMapper.H"
30 #include "reverseFieldMapper.H"
31 
32 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
33 
34 template<class Type>
35 void Foam::fvMeshAdder::MapVolField
36 (
37  const mapAddedPolyMesh& meshMap,
38 
39  VolField<Type>& fld,
40  const VolField<Type>& fldToAdd
41 )
42 {
43  const fvMesh& mesh = fld.mesh();
44 
45  // Internal field
46  // ~~~~~~~~~~~~~~
47 
48  {
49  // Store old internal field
50  Field<Type> oldInternalField(fld.primitiveField());
51 
52  // Modify internal field
53  Field<Type>& intFld = fld.primitiveFieldRef();
54 
55  intFld.setSize(mesh.nCells());
56 
57  intFld.rmap(oldInternalField, meshMap.oldCellMap());
58  intFld.rmap(fldToAdd.primitiveField(), meshMap.addedCellMap());
59  }
60 
61 
62  // Patch fields from old mesh
63  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
64 
65  typename VolField<Type>::
66  Boundary& bfld = fld.boundaryFieldRef();
67 
68  {
69  const labelList& oldPatchMap = meshMap.oldPatchMap();
70  const labelList& oldPatchStarts = meshMap.oldPatchStarts();
71  const labelList& oldPatchSizes = meshMap.oldPatchSizes();
72 
73  // Reorder old patches in order of new ones. Put removed patches at end.
74 
75  label unusedPatchi = 0;
76 
77  forAll(oldPatchMap, patchi)
78  {
79  const label newPatchi = oldPatchMap[patchi];
80 
81  if (newPatchi != -1)
82  {
83  unusedPatchi++;
84  }
85  }
86 
87  label nUsedPatches = unusedPatchi;
88 
89  // Reorder list for patchFields
90  labelList oldToNew(oldPatchMap.size());
91 
92  forAll(oldPatchMap, patchi)
93  {
94  const label newPatchi = oldPatchMap[patchi];
95 
96  if (newPatchi != -1)
97  {
98  oldToNew[patchi] = newPatchi;
99  }
100  else
101  {
102  oldToNew[patchi] = unusedPatchi++;
103  }
104  }
105 
106 
107  // Sort deleted ones last so is now in newPatch ordering
108  bfld.reorder(oldToNew);
109  // Extend to covers all patches
110  bfld.setSize(mesh.boundaryMesh().size());
111  // Delete unused patches
112  for
113  (
114  label newPatchi = nUsedPatches;
115  newPatchi < bfld.size();
116  newPatchi++
117  )
118  {
119  bfld.set(newPatchi, nullptr);
120  }
121 
122 
123  // Map old values
124  // ~~~~~~~~~~~~~~
125 
126  forAll(oldPatchMap, patchi)
127  {
128  const label newPatchi = oldPatchMap[patchi];
129 
130  if (newPatchi != -1)
131  {
132  labelList newToOld
133  (
134  calcPatchMap
135  (
136  oldPatchStarts[patchi],
137  oldPatchSizes[patchi],
138  meshMap.oldFaceMap(),
139  mesh.boundaryMesh()[newPatchi],
140  -1 // unmapped value
141  )
142  );
143 
144  // Create new patchField with same type as existing one.
145  // Note:
146  // - boundaryField already in new order so access with newPatchi
147  // - fld.boundaryField()[newPatchi] both used for type and old
148  // value
149  // - hope that field mapping allows aliasing since old and new
150  // are same memory!
151  bfld.set
152  (
153  newPatchi,
155  (
156  bfld[newPatchi], // old field
157  mesh.boundary()[newPatchi], // new fvPatch
158  fld(), // new internal field
159  forwardFieldMapper(newToOld) // mapper (new to old)
160  )
161  );
162  }
163  }
164  }
165 
166 
167 
168  // Patch fields from added mesh
169  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
170 
171  {
172  const labelList& addedPatchMap = meshMap.addedPatchMap();
173 
174  // Add addedMesh patches
175  forAll(addedPatchMap, patchi)
176  {
177  const label newPatchi = addedPatchMap[patchi];
178 
179  if (newPatchi != -1)
180  {
181  const polyPatch& newPatch = mesh.boundaryMesh()[newPatchi];
182  const polyPatch& oldPatch =
183  fldToAdd.mesh().boundaryMesh()[patchi];
184 
185  if (!bfld(newPatchi))
186  {
187  // First occurrence of newPatchi. Map from existing
188  // patchField
189 
190  // From new patch faces to patch faces on added mesh.
191  const labelList newToAdded
192  (
193  calcPatchMap
194  (
195  oldPatch.start(),
196  oldPatch.size(),
197  meshMap.addedFaceMap(),
198  newPatch,
199  -1 // unmapped values
200  )
201  );
202 
203  bfld.set
204  (
205  newPatchi,
207  (
208  fldToAdd.boundaryField()[patchi], // added field
209  mesh.boundary()[newPatchi], // new fvPatch
210  fld(), // new int. field
211  forwardFieldMapper(newToAdded) // mapper
212  )
213  );
214  }
215  else
216  {
217  // PatchField will have correct size already. Just slot in
218  // my elements.
219 
220  labelList addedToNew(oldPatch.size(), -1);
221  forAll(addedToNew, i)
222  {
223  const label addedFacei = oldPatch.start() + i;
224  const label newFacei =
225  meshMap.addedFaceMap()[addedFacei];
226  const label patchFacei = newFacei-newPatch.start();
227 
228  if (patchFacei >= 0 && patchFacei < newPatch.size())
229  {
230  addedToNew[i] = patchFacei;
231  }
232  }
233 
234  bfld[newPatchi].map
235  (
236  fldToAdd.boundaryField()[patchi],
237  reverseFieldMapper(addedToNew)
238  );
239  }
240  }
241  }
242  }
243 }
244 
245 
246 template<class Type>
248 (
249  const mapAddedPolyMesh& meshMap,
250  const fvMesh& mesh,
251  const fvMesh& meshToAdd
252 )
253 {
255  (
256  mesh.objectRegistry::lookupClass<VolField<Type>>()
257  );
258 
259  HashTable<const VolField<Type>*> fieldsToAdd
260  (
261  meshToAdd.objectRegistry::lookupClass<VolField<Type>>()
262  );
263 
264  for
265  (
266  typename HashTable<const VolField<Type>*>::
267  iterator fieldIter = fields.begin();
268  fieldIter != fields.end();
269  ++fieldIter
270  )
271  {
272  if (fvMesh::geometryFields.found(fieldIter()->name())) continue;
273 
275  const_cast<VolField<Type>&>
276  (
277  *fieldIter()
278  );
279 
280  if (fieldsToAdd.found(fld.name()))
281  {
282  const VolField<Type>& fldToAdd =
283  *fieldsToAdd[fld.name()];
284 
285  if (debug)
286  {
287  Pout<< "MapVolFields : mapping " << fld.name()
288  << " and " << fldToAdd.name() << endl;
289  }
290 
291  MapVolField<Type>(meshMap, fld, fldToAdd);
292  }
293  else
294  {
296  << "Not mapping field " << fld.name()
297  << " since not present on mesh to add"
298  << endl;
299  }
300  }
301 }
302 
303 
304 template<class Type>
305 void Foam::fvMeshAdder::MapSurfaceField
306 (
307  const mapAddedPolyMesh& meshMap,
308 
310  const SurfaceField<Type>& fldToAdd
311 )
312 {
313  const fvMesh& mesh = fld.mesh();
314  const labelList& oldPatchStarts = meshMap.oldPatchStarts();
315 
316  typename SurfaceField<Type>::
317  Boundary& bfld = fld.boundaryFieldRef();
318 
319  // Internal field
320  // ~~~~~~~~~~~~~~
321 
322  // Store old internal field
323  {
324  Field<Type> oldField(fld);
325 
326  // Modify internal field
327  Field<Type>& intFld = fld.primitiveFieldRef();
328 
329  intFld.setSize(mesh.nInternalFaces());
330 
331  intFld.rmap(oldField, meshMap.oldFaceMap());
332  intFld.rmap(fldToAdd, meshMap.addedFaceMap());
333 
334 
335  // Faces that were boundary faces but are not anymore.
336  // Use owner value (so lowest numbered cell, i.e. from 'old' not 'added'
337  // mesh)
338  forAll(bfld, patchi)
339  {
340  const fvsPatchField<Type>& pf = bfld[patchi];
341 
342  const label start = oldPatchStarts[patchi];
343 
344  forAll(pf, i)
345  {
346  const label newFacei = meshMap.oldFaceMap()[start + i];
347 
348  if (newFacei >= 0 && newFacei < mesh.nInternalFaces())
349  {
350  intFld[newFacei] = pf[i];
351  }
352  }
353  }
354  }
355 
356 
357  // Patch fields from old mesh
358  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
359 
360  {
361  const labelList& oldPatchMap = meshMap.oldPatchMap();
362  const labelList& oldPatchSizes = meshMap.oldPatchSizes();
363 
364  // Reorder old patches in order of new ones. Put removed patches at end.
365 
366  label unusedPatchi = 0;
367 
368  forAll(oldPatchMap, patchi)
369  {
370  const label newPatchi = oldPatchMap[patchi];
371 
372  if (newPatchi != -1)
373  {
374  unusedPatchi++;
375  }
376  }
377 
378  label nUsedPatches = unusedPatchi;
379 
380  // Reorder list for patchFields
381  labelList oldToNew(oldPatchMap.size());
382 
383  forAll(oldPatchMap, patchi)
384  {
385  const label newPatchi = oldPatchMap[patchi];
386 
387  if (newPatchi != -1)
388  {
389  oldToNew[patchi] = newPatchi;
390  }
391  else
392  {
393  oldToNew[patchi] = unusedPatchi++;
394  }
395  }
396 
397 
398  // Sort deleted ones last so is now in newPatch ordering
399  bfld.reorder(oldToNew);
400  // Extend to covers all patches
401  bfld.setSize(mesh.boundaryMesh().size());
402  // Delete unused patches
403  for
404  (
405  label newPatchi = nUsedPatches;
406  newPatchi < bfld.size();
407  newPatchi++
408  )
409  {
410  bfld.set(newPatchi, nullptr);
411  }
412 
413 
414  // Map old values
415  // ~~~~~~~~~~~~~~
416 
417  forAll(oldPatchMap, patchi)
418  {
419  const label newPatchi = oldPatchMap[patchi];
420 
421  if (newPatchi != -1)
422  {
423  const labelList newToOld
424  (
425  calcPatchMap
426  (
427  oldPatchStarts[patchi],
428  oldPatchSizes[patchi],
429  meshMap.oldFaceMap(),
430  mesh.boundaryMesh()[newPatchi],
431  -1 // unmapped value
432  )
433  );
434 
435  // Create new patchField with same type as existing one.
436  // Note:
437  // - boundaryField already in new order so access with newPatchi
438  // - bfld[newPatchi] both used for type and old
439  // value
440  // - hope that field mapping allows aliasing since old and new
441  // are same memory!
442  bfld.set
443  (
444  newPatchi,
446  (
447  bfld[newPatchi], // old field
448  mesh.boundary()[newPatchi], // new fvPatch
449  fld(), // new internal field
450  forwardFieldMapper(newToOld) // mapper (new to old)
451  )
452  );
453  }
454  }
455  }
456 
457 
458 
459  // Patch fields from added mesh
460  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
461 
462  {
463  const labelList& addedPatchMap = meshMap.addedPatchMap();
464 
465  // Add addedMesh patches
466  forAll(addedPatchMap, patchi)
467  {
468  const label newPatchi = addedPatchMap[patchi];
469 
470  if (newPatchi != -1)
471  {
472  const polyPatch& newPatch = mesh.boundaryMesh()[newPatchi];
473  const polyPatch& oldPatch =
474  fldToAdd.mesh().boundaryMesh()[patchi];
475 
476  if (!bfld(newPatchi))
477  {
478  // First occurrence of newPatchi. Map from existing
479  // patchField
480 
481  // From new patch faces to patch faces on added mesh.
482  const labelList newToAdded
483  (
484  calcPatchMap
485  (
486  oldPatch.start(),
487  oldPatch.size(),
488  meshMap.addedFaceMap(),
489  newPatch,
490  -1 // unmapped values
491  )
492  );
493 
494  bfld.set
495  (
496  newPatchi,
498  (
499  fldToAdd.boundaryField()[patchi],// added field
500  mesh.boundary()[newPatchi], // new fvPatch
501  fld(), // new int. field
502  forwardFieldMapper(newToAdded) // mapper
503  )
504  );
505  }
506  else
507  {
508  // PatchField will have correct size already. Just slot in
509  // my elements.
510 
511  labelList addedToNew(oldPatch.size(), -1);
512  forAll(addedToNew, i)
513  {
514  const label addedFacei = oldPatch.start() + i;
515  const label newFacei =
516  meshMap.addedFaceMap()[addedFacei];
517  const label patchFacei = newFacei-newPatch.start();
518 
519  if (patchFacei >= 0 && patchFacei < newPatch.size())
520  {
521  addedToNew[i] = patchFacei;
522  }
523  }
524 
525  bfld[newPatchi].map
526  (
527  fldToAdd.boundaryField()[patchi],
528  reverseFieldMapper(addedToNew)
529  );
530  }
531  }
532  }
533  }
534 }
535 
536 
537 template<class Type>
539 (
540  const mapAddedPolyMesh& meshMap,
541  const fvMesh& mesh,
542  const fvMesh& meshToAdd
543 )
544 {
546  (
547  mesh.objectRegistry::lookupClass<SurfaceField<Type>>()
548  );
549 
551  (
552  meshToAdd.objectRegistry::lookupClass<SurfaceField<Type>>()
553  );
554 
555  for
556  (
557  typename HashTable<const SurfaceField<Type>*>::
558  iterator fieldIter = fields.begin();
559  fieldIter != fields.end();
560  ++fieldIter
561  )
562  {
563  if (fvMesh::geometryFields.found(fieldIter()->name())) continue;
564 
565  SurfaceField<Type>& fld = const_cast<SurfaceField<Type>&>(*fieldIter());
566 
567  if (fieldsToAdd.found(fld.name()))
568  {
569  const SurfaceField<Type>& fldToAdd = *fieldsToAdd[fld.name()];
570 
571  if (debug)
572  {
573  Pout<< "MapSurfaceFields : mapping " << fld.name()
574  << " and " << fldToAdd.name() << endl;
575  }
576 
577  MapSurfaceField<Type>(meshMap, fld, fldToAdd);
578  }
579  else
580  {
582  << "Not mapping field " << fld.name()
583  << " since not present on mesh to add"
584  << endl;
585  }
586  }
587 }
588 
589 
590 template<class Type>
591 void Foam::fvMeshAdder::MapPointField
592 (
593  const pointMesh& mesh,
594  const mapAddedPolyMesh& meshMap,
595  const labelListList& oldMeshPoints,
596 
598  const PointField<Type>& fldToAdd
599 )
600 {
601  // This is a bit tricky:
602  // - mesh pointed to by fld is invalid
603  // - pointPatches pointed to be fld are invalid
604 
605  typename PointField<Type>::
606  Boundary& bfld = fld.boundaryFieldRef();
607 
608  // Internal field
609  // ~~~~~~~~~~~~~~
610 
611  // Store old internal field
612  {
613  Field<Type> oldField(fld);
614 
615  // Modify internal field
616  Field<Type>& intFld = fld.primitiveFieldRef();
617 
618  intFld.setSize(mesh.size());
619 
620  intFld.rmap(oldField, meshMap.oldPointMap());
621  intFld.rmap(fldToAdd, meshMap.addedPointMap());
622  }
623 
624 
625  // Patch fields from old mesh
626  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
627 
628  {
629  const labelList& oldPatchMap = meshMap.oldPatchMap();
630 
631  // Reorder old patches in order of new ones. Put removed patches at end.
632 
633  label unusedPatchi = 0;
634 
635  forAll(oldPatchMap, patchi)
636  {
637  const label newPatchi = oldPatchMap[patchi];
638 
639  if (newPatchi != -1)
640  {
641  unusedPatchi++;
642  }
643  }
644 
645  label nUsedPatches = unusedPatchi;
646 
647  // Reorder list for patchFields
648  labelList oldToNew(oldPatchMap.size());
649 
650  forAll(oldPatchMap, patchi)
651  {
652  const label newPatchi = oldPatchMap[patchi];
653 
654  if (newPatchi != -1)
655  {
656  oldToNew[patchi] = newPatchi;
657  }
658  else
659  {
660  oldToNew[patchi] = unusedPatchi++;
661  }
662  }
663 
664 
665  // Sort deleted ones last so is now in newPatch ordering
666  bfld.reorder(oldToNew);
667 
668  // Extend to covers all patches
669  bfld.setSize(mesh.boundary().size());
670 
671  // Delete unused patches
672  for
673  (
674  label newPatchi = nUsedPatches;
675  newPatchi < bfld.size();
676  newPatchi++
677  )
678  {
679  bfld.set(newPatchi, nullptr);
680  }
681 
682 
683  // Map old values
684  // ~~~~~~~~~~~~~~
685 
686  forAll(oldPatchMap, patchi)
687  {
688  const label newPatchi = oldPatchMap[patchi];
689 
690  if (newPatchi != -1)
691  {
692  const labelList& oldMp = oldMeshPoints[patchi];
693  const pointPatch& newPp = mesh.boundary()[newPatchi];
694  const labelList& newMeshPoints = newPp.meshPoints();
695  const labelList& oldPointMap = meshMap.oldPointMap();
696 
697  Map<label> newMeshPointMap(2*newMeshPoints.size());
698  forAll(newMeshPoints, ppi)
699  {
700  newMeshPointMap.insert(newMeshPoints[ppi], ppi);
701  }
702 
703  labelList newToOld(newPp.size(), -1);
704  forAll(oldMp, oldPointi)
705  {
706  const label newPointi = oldPointMap[oldMp[oldPointi]];
707 
709  newMeshPointMap.find(newPointi);
710 
711  if (fnd == newMeshPointMap.end())
712  {
713  // Possibly an internal point
714  }
715  else
716  {
717  // Part of new patch
718  newToOld[fnd()] = oldPointi;
719  }
720  }
721 
722  // Create new patchField with same type as existing one.
723  // Note:
724  // - boundaryField already in new order so access with newPatchi
725  // - bfld[newPatchi] both used for type and old
726  // value
727  // - hope that field mapping allows aliasing since old and new
728  // are same memory!
729  bfld.set
730  (
731  newPatchi,
733  (
734  bfld[newPatchi], // old field
735  mesh.boundary()[newPatchi], // new pointPatch
736  fld(), // new internal field
737  forwardFieldMapper(newToOld) // mapper (new to old)
738  )
739  );
740  }
741  }
742  }
743 
744 
745 
746  // Patch fields from added mesh
747  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
748 
749  {
750  const labelList& addedPatchMap = meshMap.addedPatchMap();
751 
752  // Add addedMesh patches
753  forAll(addedPatchMap, patchi)
754  {
755  const label newPatchi = addedPatchMap[patchi];
756 
757  if (newPatchi != -1)
758  {
759  const pointPatch& oldPatch = fldToAdd.mesh().boundary()[patchi];
760  const labelList& oldMp = oldPatch.meshPoints();
761  const pointPatch& newPp = mesh.boundary()[newPatchi];
762  const labelList& newMeshPoints = newPp.meshPoints();
763  const labelList& addedPointMap = meshMap.addedPointMap();
764 
765  Map<label> newMpm(2*newMeshPoints.size());
766  forAll(newMeshPoints, ppi)
767  {
768  newMpm.insert(newMeshPoints[ppi], ppi);
769  }
770 
771  if (!bfld(newPatchi))
772  {
773  // First occurrence of newPatchi. Map from existing
774  // patchField
775 
776  labelList newToAdded(newPp.size(), -1);
777  forAll(oldMp, oldPointi)
778  {
779  const label newPointi = addedPointMap[oldMp[oldPointi]];
780 
781  Map<label>::const_iterator fnd = newMpm.find(newPointi);
782  if (fnd == newMpm.end())
783  {
784  // Possibly an internal point
785  }
786  else
787  {
788  // Part of new patch
789  newToAdded[fnd()] = oldPointi;
790  }
791  }
792 
793  bfld.set
794  (
795  newPatchi,
797  (
798  fldToAdd.boundaryField()[patchi],// added field
799  mesh.boundary()[newPatchi], // new pointPatch
800  fld(), // new int. field
801  forwardFieldMapper(newToAdded) // mapper
802  )
803  );
804  }
805  else
806  {
807  // PatchField will have correct size already. Just slot in
808  // my elements.
809 
810  labelList oldToNew(oldPatch.size(), -1);
811  forAll(oldMp, oldPointi)
812  {
813  const label newPointi = addedPointMap[oldMp[oldPointi]];
814 
815  Map<label>::const_iterator fnd = newMpm.find(newPointi);
816  if (fnd != newMpm.end())
817  {
818  // Part of new patch
819  oldToNew[oldPointi] = fnd();
820  }
821  }
822 
823  bfld[newPatchi].map
824  (
825  fldToAdd.boundaryField()[patchi],
826  reverseFieldMapper(oldToNew)
827  );
828  }
829  }
830  }
831  }
832 }
833 
834 
835 template<class Type>
837 (
838  const mapAddedPolyMesh& meshMap,
839  const pointMesh& mesh,
840  const labelListList& oldMeshPoints,
841  const objectRegistry& meshToAdd
842 )
843 {
845  (
847  );
848 
850  (
851  meshToAdd.lookupClass<PointField<Type>>()
852  );
853 
854  for
855  (
856  typename HashTable<const PointField<Type>*>::
857  iterator fieldIter = fields.begin();
858  fieldIter != fields.end();
859  ++fieldIter
860  )
861  {
862  if (fvMesh::geometryFields.found(fieldIter()->name())) continue;
863 
864  PointField<Type>& fld = const_cast<PointField<Type>&>(*fieldIter());
865 
866  if (fieldsToAdd.found(fld.name()))
867  {
868  const PointField<Type>& fldToAdd = *fieldsToAdd[fld.name()];
869 
870  if (debug)
871  {
872  Pout<< "MapPointFields : mapping " << fld.name()
873  << " and " << fldToAdd.name() << endl;
874  }
875 
876  MapPointField<Type>(mesh, meshMap, oldMeshPoints, fld, fldToAdd);
877  }
878  else
879  {
881  << "Not mapping field " << fld.name()
882  << " since not present on mesh to add"
883  << endl;
884  }
885  }
886 }
887 
888 
889 template<class Type>
890 void Foam::fvMeshAdder::MapDimField
891 (
892  const mapAddedPolyMesh& meshMap,
893 
895  const DimensionedField<Type, volMesh>& fldToAdd
896 )
897 {
898  const fvMesh& mesh = fld.mesh();
899 
900  // Store old field
901  Field<Type> oldField(fld);
902 
903  fld.setSize(mesh.nCells());
904 
905  fld.rmap(oldField, meshMap.oldCellMap());
906  fld.rmap(fldToAdd, meshMap.addedCellMap());
907 }
908 
909 
910 template<class Type>
912 (
913  const mapAddedPolyMesh& meshMap,
914  const fvMesh& mesh,
915  const fvMesh& meshToAdd
916 )
917 {
918  // Note: use strict flag on lookupClass to avoid picking up volFields
920  (
921  mesh.objectRegistry::lookupClass<VolInternalField<Type>>(true)
922  );
923 
925  (
926  meshToAdd.objectRegistry::lookupClass<VolInternalField<Type>>(true)
927  );
928 
929  for
930  (
931  typename HashTable<const VolInternalField<Type>*>::
932  iterator fieldIter = fields.begin();
933  fieldIter != fields.end();
934  ++fieldIter
935  )
936  {
938  const_cast<VolInternalField<Type>&>(*fieldIter());
939 
940  if (fieldsToAdd.found(fld.name()))
941  {
942  const VolInternalField<Type>& fldToAdd = *fieldsToAdd[fld.name()];
943 
944  if (debug)
945  {
946  Pout<< "MapDimFields : mapping " << fld.name()
947  << " and " << fldToAdd.name() << endl;
948  }
949 
950  MapDimField<Type>(meshMap, fld, fldToAdd);
951  }
952  else
953  {
954  WarningIn("fvMeshAdder::MapDimFields(..)")
955  << "Not mapping field " << fld.name()
956  << " since not present on mesh to add"
957  << endl;
958  }
959  }
960 }
961 
962 
963 // ************************************************************************* //
bool found
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Mesh & mesh() const
Return mesh.
Pre-declare SubField and related Field type.
Definition: Field.H:83
void map(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 map from the given field
Definition: Field.C:289
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:386
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
An STL-conforming hash table.
Definition: HashTable.H:127
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:138
friend class const_iterator
Declare friendship with the const_iterator.
Definition: HashTable.H:197
const word & name() const
Return name.
Definition: IOobject.H:310
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void setSize(const label)
Reset size of List.
Definition: List.C:281
Map point field on topology change. This is a partial template specialisation for GeoMesh=pointMesh.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
static void MapVolFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all volFields of Type.
static void MapDimFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all DimensionedFields of Type.
static void MapSurfaceFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd)
Map all surfaceFields of Type.
static void MapPointFields(const mapAddedPolyMesh &, const pointMesh &mesh, const labelListList &oldMeshPoints, const objectRegistry &meshToAdd)
Map all surfaceFields of Type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:857
static const HashSet< word > geometryFields
Set of names of registered geometric fields.
Definition: fvMesh.H:331
static tmp< fvPatchField< Type > > New(const word &, const fvPatch &, const DimensionedField< Type, volMesh > &)
Return a pointer to a new patchField created on freestore given.
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:83
static tmp< fvsPatchField< Type > > New(const word &, const fvPatch &, const DimensionedField< Type, surfaceMesh > &)
Return a pointer to a new patchField created on freestore given.
Class containing mesh-to-mesh mapping information after a mesh addition where we add a mesh ('added m...
const labelList & addedFaceMap() const
const labelList & oldPatchMap() const
From old patch index to new patch index or -1 if patch.
const labelList & addedCellMap() const
const labelList & addedPatchMap() const
From added mesh patch index to new patch index or -1 if.
const labelList & oldPointMap() const
From old mesh point/face/cell to new mesh point/face/cell.
const labelList & oldCellMap() const
const labelList & oldPatchStarts() const
Return list of the old patch start labels.
const labelList & addedPointMap() const
From added mesh point/face/cell to new mesh point/face/cell.
const labelList & oldFaceMap() const
const labelList & oldPatchSizes() const
Return list of the old patch sizes.
Registry of regIOobjects.
HashTable< const Type * > lookupClass(const bool strict=false) const
Lookup and return all objects of the given Type.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:53
const objectRegistry & thisDb() const
Return database. For now is its polyMesh.
Definition: pointMesh.H:137
label size() const
Return number of points.
Definition: pointMesh.H:113
const pointBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: pointMesh.H:125
static autoPtr< pointPatchField< Type > > New(const word &, const pointPatch &, const DimensionedField< Type, pointMesh > &)
Return a pointer to a new patchField created on freestore given.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:404
label nCells() const
label nInternalFaces() const
label patchi
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:228
#define WarningInFunction
Report a warning using Foam::Warning.
#define WarningIn(functionName)
Report a warning using Foam::Warning.
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
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
typename VolField< Type >::Internal VolInternalField
Definition: volFieldsFwd.H:61
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
label newPointi
Definition: readKivaGrid.H:495
Foam::surfaceFields.