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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2015-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "volFields.H"
30#include "surfaceFields.H"
31#include "emptyFvPatchField.H"
33
34// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
35
36template<class Type>
37void Foam::fvMeshAdder::MapVolField
38(
39 const mapAddedPolyMesh& meshMap,
40
41 GeometricField<Type, fvPatchField, volMesh>& fld,
42 const GeometricField<Type, fvPatchField, volMesh>& fldToAdd,
43 const bool fullyMapped
44)
45{
46 const fvMesh& mesh = fld.mesh();
47
48 // Internal field
49 // ~~~~~~~~~~~~~~
50
51 {
52 // Store old internal field
53 const Field<Type> oldInternalField(fld.primitiveField());
54
55 // Modify internal field
56 Field<Type>& intFld = fld.primitiveFieldRef();
57
58 intFld.setSize(mesh.nCells());
59
60 intFld.rmap(oldInternalField, meshMap.oldCellMap());
61 intFld.rmap(fldToAdd.primitiveField(), meshMap.addedCellMap());
62 }
63
64
65 // Patch fields from old mesh
66 // ~~~~~~~~~~~~~~~~~~~~~~~~~~
67
68 auto& bfld = fld.boundaryFieldRef();
69
70 {
71 const labelList& oldPatchMap = meshMap.oldPatchMap();
72 const labelList& oldPatchStarts = meshMap.oldPatchStarts();
73 const labelList& oldPatchSizes = meshMap.oldPatchSizes();
74
75 // Reorder old patches in order of new ones. Put removed patches at end.
76
77 label unusedPatchi = 0;
78
79 forAll(oldPatchMap, patchi)
80 {
81 label newPatchi = oldPatchMap[patchi];
82
83 if (newPatchi != -1)
84 {
85 unusedPatchi++;
86 }
87 }
88
89 label nUsedPatches = unusedPatchi;
90
91 // Reorder list for patchFields
92 labelList oldToNew(oldPatchMap.size());
93
94 forAll(oldPatchMap, patchi)
95 {
96 const label newPatchi = oldPatchMap[patchi];
97
98 if (newPatchi != -1)
99 {
100 oldToNew[patchi] = newPatchi;
101 }
102 else
103 {
104 oldToNew[patchi] = unusedPatchi++;
105 }
106 }
107
108
109 // Sort deleted ones last so is now in newPatch ordering
110 bfld.reorder(oldToNew);
111 // Extend to covers all patches
112 bfld.setSize(mesh.boundaryMesh().size());
113 // Delete unused patches
114 for
115 (
116 label newPatchi = nUsedPatches;
117 newPatchi < bfld.size();
118 newPatchi++
119 )
120 {
121 bfld.set(newPatchi, nullptr);
122 }
123
124
125 // Map old values
126 // ~~~~~~~~~~~~~~
127
128 forAll(oldPatchMap, patchi)
129 {
130 const label newPatchi = oldPatchMap[patchi];
131
132 if (newPatchi != -1)
133 {
134 labelList newToOld
135 (
136 calcPatchMap
137 (
138 oldPatchStarts[patchi],
139 oldPatchSizes[patchi],
140 meshMap.oldFaceMap(),
141 mesh.boundaryMesh()[newPatchi],
142 -1 // unmapped value
143 )
144 );
145
146 directFvPatchFieldMapper patchMapper(newToOld);
147
148 // Override mapping (for use in e.g. fvMeshDistribute where
149 // it sorts mapping out itself)
150 if (fullyMapped)
151 {
152 patchMapper.hasUnmapped() = false;
153 }
154
155 // Create new patchField with same type as existing one.
156 // Note:
157 // - boundaryField already in new order so access with newPatchi
158 // - fld.boundaryField()[newPatchi] both used for type and old
159 // value
160 // - hope that field mapping allows aliasing since old and new
161 // are same memory!
162 bfld.set
163 (
164 newPatchi,
166 (
167 bfld[newPatchi], // old field
168 mesh.boundary()[newPatchi], // new fvPatch
169 fld(), // new internal field
170 patchMapper // mapper (new to old)
171 )
172 );
173 }
174 }
175 }
176
177
178
179 // Patch fields from added mesh
180 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
181
182 {
183 const labelList& addedPatchMap = meshMap.addedPatchMap();
184
185 // Add addedMesh patches
186 forAll(addedPatchMap, patchi)
187 {
188 const label newPatchi = addedPatchMap[patchi];
189
190 if (newPatchi != -1)
191 {
192 const polyPatch& newPatch = mesh.boundaryMesh()[newPatchi];
193 const polyPatch& oldPatch =
194 fldToAdd.mesh().boundaryMesh()[patchi];
195
196 if (!bfld(newPatchi))
197 {
198 // First occurrence of newPatchi. Map from existing
199 // patchField
200
201 // From new patch faces to patch faces on added mesh.
202 labelList newToAdded
203 (
204 calcPatchMap
205 (
206 oldPatch.start(),
207 oldPatch.size(),
208 meshMap.addedFaceMap(),
209 newPatch,
210 -1 // unmapped values
211 )
212 );
213
214 directFvPatchFieldMapper patchMapper(newToAdded);
215
216 // Override mapping (for use in e.g. fvMeshDistribute where
217 // it sorts mapping out itself)
218 if (fullyMapped)
219 {
220 patchMapper.hasUnmapped() = false;
221 }
222
223 bfld.set
224 (
225 newPatchi,
227 (
228 fldToAdd.boundaryField()[patchi], // added field
229 mesh.boundary()[newPatchi], // new fvPatch
230 fld(), // new int. field
231 patchMapper // mapper
232 )
233 );
234 }
235 else
236 {
237 // PatchField will have correct size already. Just slot in
238 // my elements.
239
240 labelList addedToNew(oldPatch.size(), -1);
241 forAll(addedToNew, i)
242 {
243 label addedFacei = oldPatch.start()+i;
244 label newFacei = meshMap.addedFaceMap()[addedFacei];
245 label patchFacei = newFacei-newPatch.start();
246 if (patchFacei >= 0 && patchFacei < newPatch.size())
247 {
248 addedToNew[i] = patchFacei;
249 }
250 }
251
252 bfld[newPatchi].rmap
253 (
254 fldToAdd.boundaryField()[patchi],
255 addedToNew
256 );
257 }
258 }
259 }
260 }
261}
262
263
264template<class Type>
266(
267 const mapAddedPolyMesh& meshMap,
268 const fvMesh& mesh,
269 const fvMesh& meshToAdd,
270 const bool fullyMapped
271)
272{
274
276 (
277 mesh.objectRegistry::lookupClass<fldType>()
278 );
279
280 HashTable<const fldType*> fieldsToAdd
281 (
282 meshToAdd.objectRegistry::lookupClass<fldType>()
283 );
284
285 // It is necessary to enforce that all old-time fields are stored
286 // before the mapping is performed. Otherwise, if the
287 // old-time-level field is mapped before the field itself, sizes
288 // will not match.
289
290 forAllIters(fields, fieldIter)
291 {
292 fldType& fld = const_cast<fldType&>(*fieldIter());
293
295 << "MapVolFields : Storing old time for " << fld.name() << endl;
296
297 fld.storeOldTimes();
298 }
299
300
301 forAllIters(fields, fieldIter)
302 {
303 fldType& fld = const_cast<fldType&>(*fieldIter());
304
305 if (fieldsToAdd.found(fld.name()))
306 {
307 const fldType& fldToAdd = *fieldsToAdd[fld.name()];
308
310 << "MapVolFields : mapping " << fld.name()
311 << " and " << fldToAdd.name() << endl;
312
313 MapVolField<Type>(meshMap, fld, fldToAdd, fullyMapped);
314 }
315 else
316 {
318 << "Not mapping field " << fld.name()
319 << " since not present on mesh to add" << endl;
320 }
321 }
322}
323
324
325template<class Type>
326void Foam::fvMeshAdder::MapSurfaceField
327(
328 const mapAddedPolyMesh& meshMap,
329
332 const bool fullyMapped
333)
334{
335 const fvMesh& mesh = fld.mesh();
336 const labelList& oldPatchStarts = meshMap.oldPatchStarts();
337
338 auto& bfld = fld.boundaryFieldRef();
339
340 // Internal field
341 // ~~~~~~~~~~~~~~
342
343 // Store old internal field
344 {
345 const Field<Type> oldField(fld);
346
347 // Modify internal field
348 Field<Type>& intFld = fld.primitiveFieldRef();
349
350 intFld.setSize(mesh.nInternalFaces());
351
352 intFld.rmap(oldField, meshMap.oldFaceMap());
353 intFld.rmap(fldToAdd, meshMap.addedFaceMap());
354
355
356 // Faces that were boundary faces but are not anymore.
357 // Use owner value (so lowest numbered cell, i.e. from 'old' not 'added'
358 // mesh)
359 forAll(bfld, patchi)
360 {
361 const fvsPatchField<Type>& pf = bfld[patchi];
362
363 label start = oldPatchStarts[patchi];
364
365 forAll(pf, i)
366 {
367 label newFacei = meshMap.oldFaceMap()[start + i];
368
369 if (newFacei >= 0 && newFacei < mesh.nInternalFaces())
370 {
371 intFld[newFacei] = pf[i];
372 }
373 }
374 }
375 }
376
377
378 // Patch fields from old mesh
379 // ~~~~~~~~~~~~~~~~~~~~~~~~~~
380
381 {
382 const labelList& oldPatchMap = meshMap.oldPatchMap();
383 const labelList& oldPatchSizes = meshMap.oldPatchSizes();
384
385 // Reorder old patches in order of new ones. Put removed patches at end.
386
387 label unusedPatchi = 0;
388
389 forAll(oldPatchMap, patchi)
390 {
391 const label newPatchi = oldPatchMap[patchi];
392
393 if (newPatchi != -1)
394 {
395 unusedPatchi++;
396 }
397 }
398
399 label nUsedPatches = unusedPatchi;
400
401 // Reorder list for patchFields
402 labelList oldToNew(oldPatchMap.size());
403
404 forAll(oldPatchMap, patchi)
405 {
406 const label newPatchi = oldPatchMap[patchi];
407
408 if (newPatchi != -1)
409 {
410 oldToNew[patchi] = newPatchi;
411 }
412 else
413 {
414 oldToNew[patchi] = unusedPatchi++;
415 }
416 }
417
418
419 // Sort deleted ones last so is now in newPatch ordering
420 bfld.reorder(oldToNew);
421 // Extend to covers all patches
422 bfld.setSize(mesh.boundaryMesh().size());
423 // Delete unused patches
424 for
425 (
426 label newPatchi = nUsedPatches;
427 newPatchi < bfld.size();
428 newPatchi++
429 )
430 {
431 bfld.set(newPatchi, nullptr);
432 }
433
434
435 // Map old values
436 // ~~~~~~~~~~~~~~
437
438 forAll(oldPatchMap, patchi)
439 {
440 const label newPatchi = oldPatchMap[patchi];
441
442 if (newPatchi != -1)
443 {
444 labelList newToOld
445 (
446 calcPatchMap
447 (
448 oldPatchStarts[patchi],
449 oldPatchSizes[patchi],
450 meshMap.oldFaceMap(),
451 mesh.boundaryMesh()[newPatchi],
452 -1 // unmapped value
453 )
454 );
455
456 directFvPatchFieldMapper patchMapper(newToOld);
457
458 // Override mapping (for use in e.g. fvMeshDistribute where
459 // it sorts mapping out itself)
460 if (fullyMapped)
461 {
462 patchMapper.hasUnmapped() = false;
463 }
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 const 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 // Override mapping (for use in e.g. fvMeshDistribute where
527 // it sorts mapping out itself)
528 if (fullyMapped)
529 {
530 patchMapper.hasUnmapped() = false;
531 }
532
533 bfld.set
534 (
535 newPatchi,
537 (
538 fldToAdd.boundaryField()[patchi],// added field
539 mesh.boundary()[newPatchi], // new fvPatch
540 fld(), // new int. field
541 patchMapper // mapper
542 )
543 );
544 }
545 else
546 {
547 // PatchField will have correct size already. Just slot in
548 // my elements.
549
550 labelList addedToNew(oldPatch.size(), -1);
551 forAll(addedToNew, i)
552 {
553 label addedFacei = oldPatch.start()+i;
554 label newFacei = meshMap.addedFaceMap()[addedFacei];
555 label patchFacei = newFacei-newPatch.start();
556 if (patchFacei >= 0 && patchFacei < newPatch.size())
557 {
558 addedToNew[i] = patchFacei;
559 }
560 }
561
562 bfld[newPatchi].rmap
563 (
564 fldToAdd.boundaryField()[patchi],
565 addedToNew
566 );
567 }
568 }
569 }
570 }
571}
572
573
574template<class Type>
576(
577 const mapAddedPolyMesh& meshMap,
578 const fvMesh& mesh,
579 const fvMesh& meshToAdd,
580 const bool fullyMapped
581)
582{
584
586 (
587 mesh.objectRegistry::lookupClass<fldType>()
588 );
589
590 HashTable<const fldType*> fieldsToAdd
591 (
592 meshToAdd.objectRegistry::lookupClass<fldType>()
593 );
594
595 // It is necessary to enforce that all old-time fields are stored
596 // before the mapping is performed. Otherwise, if the
597 // old-time-level field is mapped before the field itself, sizes
598 // will not match.
599
600 forAllIters(fields, fieldIter)
601 {
602 fldType& fld = const_cast<fldType&>(*fieldIter());
603
605 << "MapSurfaceFields : Storing old time for " << fld.name() << endl;
606
607 fld.storeOldTimes();
608 }
609
610
611 forAllIters(fields, fieldIter)
612 {
613 fldType& fld = const_cast<fldType&>(*fieldIter());
614
615 if (fieldsToAdd.found(fld.name()))
616 {
617 const fldType& fldToAdd = *fieldsToAdd[fld.name()];
618
620 << "MapSurfaceFields : mapping " << fld.name()
621 << " and " << fldToAdd.name() << endl;
622
623 MapSurfaceField<Type>(meshMap, fld, fldToAdd, fullyMapped);
624 }
625 else
626 {
628 << "Not mapping field " << fld.name()
629 << " since not present on mesh to add" << endl;
630 }
631 }
632}
633
634
635template<class Type>
636void Foam::fvMeshAdder::MapDimField
637(
638 const mapAddedPolyMesh& meshMap,
639
641 const DimensionedField<Type, volMesh>& fldToAdd
642)
643{
644 const fvMesh& mesh = fld.mesh();
645
646 // Store old field
647 const Field<Type> oldField(fld);
648
649 fld.setSize(mesh.nCells());
650
651 fld.rmap(oldField, meshMap.oldCellMap());
652 fld.rmap(fldToAdd, meshMap.addedCellMap());
653}
654
655
656template<class Type>
658(
659 const mapAddedPolyMesh& meshMap,
660 const fvMesh& mesh,
661 const fvMesh& meshToAdd
662)
663{
664 typedef DimensionedField<Type, volMesh> fldType;
665
666 // Note: use strict flag on lookupClass to avoid picking up
667 // volFields
669 (
670 mesh.objectRegistry::lookupClass<fldType>(true)
671 );
672
673 HashTable<const fldType*> fieldsToAdd
674 (
675 meshToAdd.objectRegistry::lookupClass<fldType>(true)
676 );
677
678 forAllIters(fields, fieldIter)
679 {
680 fldType& fld = const_cast<fldType&>(*fieldIter());
681
682 if (fieldsToAdd.found(fld.name()))
683 {
684 const fldType& fldToAdd = *fieldsToAdd[fld.name()];
685
687 << "MapDimFields : mapping " << fld.name()
688 << " and " << fldToAdd.name() << endl;
689
690 MapDimField<Type>(meshMap, fld, fldToAdd);
691 }
692 else
693 {
695 << "Not mapping field " << fld.name()
696 << " since not present on mesh to add" << endl;
697 }
698 }
699}
700
701
702//- Multi-mesh mapping
703
704template<class Type>
705void Foam::fvMeshAdder::MapDimField
706(
708 const labelListList& cellProcAddressing,
709 const bool fullyMapped
710)
711{
712 // Add fields to fields[0] after adding the meshes to meshes[0].
713 // Mesh[0] is the sum of all meshes. Fields are not yet mapped.
714
715 if
716 (
717 flds.size() == 0
718 || !flds.set(0)
719 || cellProcAddressing.size() != flds.size()
720 )
721 {
722 FatalErrorInFunction << "Not valid field at element 0"
723 << " in field list of size " << flds.size() << exit(FatalError);
724 }
725
726
727 // Internal field
728 // ~~~~~~~~~~~~~~
729
730 {
731 // Store old internal field
732 const Field<Type> oldInternalField(flds[0]);
733
734 // Modify internal field
735 Field<Type>& intFld = flds[0];
736
737 // Set to new mesh size
738 intFld.setSize(flds[0].mesh().nCells());
739 // Add fld0
740 intFld.rmap(oldInternalField, cellProcAddressing[0]);
741
742 for (label meshi = 1; meshi < flds.size(); meshi++)
743 {
744 if (flds.set(meshi))
745 {
746 const Field<Type>& addFld = flds[meshi];
747 intFld.rmap(addFld, cellProcAddressing[meshi]);
748 }
749 }
750 }
751}
752
753
754template<class Type>
755void Foam::fvMeshAdder::MapVolField
756(
758 const labelList& oldPatchStarts0,
759 const labelList& oldPatchSizes0,
760 const labelListList& patchProcAddressing,
761 const labelListList& cellProcAddressing,
762 const labelListList& faceProcAddressing,
763 const bool fullyMapped
764)
765{
766 // Add fields to fields[0] after adding the meshes to meshes[0].
767 // Mesh[0] is the sum of all meshes. Fields are not yet mapped.
768
769 if (flds.size() == 0 || !flds.set(0))
770 {
771 FatalErrorInFunction << "Not valid field at element 0"
772 << " in field list of size " << flds.size() << exit(FatalError);
773 }
774
775
776 // Internal field
777 // ~~~~~~~~~~~~~~
778
779 {
780 // Store old internal field
781 const Field<Type> oldInternalField(flds[0].primitiveField());
782
783 // Modify internal field
784 Field<Type>& intFld = flds[0].primitiveFieldRef();
785
786 // Set to new mesh size
787 intFld.setSize(flds[0].mesh().nCells());
788 // Add fld0
789 intFld.rmap(oldInternalField, cellProcAddressing[0]);
790
791 for (label meshi = 1; meshi < flds.size(); meshi++)
792 {
793 if (flds.set(meshi))
794 {
795 const Field<Type>& addFld = flds[meshi].primitiveFieldRef();
796 intFld.rmap(addFld, cellProcAddressing[meshi]);
797 }
798 }
799 }
800
801
802 // Patch fields from old mesh
803 // ~~~~~~~~~~~~~~~~~~~~~~~~~~
804
805 auto& bfld0 = flds[0].boundaryFieldRef();
806 //Pout<< " Adding patchfields for field " << flds[0].name() << endl;
807 forAll(bfld0, patchi)
808 {
809 //Pout<< " patch:" << patchi
810 // << " patch:" << flds[0].mesh().boundaryMesh()[patchi].name()
811 // << " size:" << flds[0].mesh().boundaryMesh()[patchi].size()
812 // << endl;
813 //Pout<< " patchField:" << bfld0[patchi].size()
814 // << endl;
815
816 labelList newToOld
817 (
818 calcPatchMap
819 (
820 oldPatchStarts0[patchi],
821 oldPatchSizes0[patchi],
822 faceProcAddressing[0],
823 bfld0[patchi].patch().patch(),
824 -1 // unmapped value
825 )
826 );
827
828 directFvPatchFieldMapper patchMapper(newToOld);
829
830 // Override mapping (for use in e.g. fvMeshDistribute where
831 // it sorts mapping out itself)
832 if (fullyMapped)
833 {
834 patchMapper.hasUnmapped() = false;
835 }
836
837 bfld0[patchi].autoMap(patchMapper);
838 }
839
840 for (label meshi = 1; meshi < flds.size(); meshi++)
841 {
842 if (flds.set(meshi))
843 {
844 const auto& bfld = flds[meshi].boundaryFieldRef();
845
846 const labelList& patchMap = patchProcAddressing[meshi];
847
848 forAll(patchMap, oldPatchi)
849 {
850 const auto& fvp = bfld[oldPatchi].patch();
851 const label newPatchi = patchMap[oldPatchi];
852
853 //Pout<< " oldPatch:" << oldPatchi
854 // << " newPatch:" << newPatchi << endl;
855
856 if (newPatchi >= 0 && newPatchi < bfld0.size())
857 {
858 const auto& fvp0 = bfld0[newPatchi].patch();
859 labelList addedToNew(bfld[oldPatchi].size(), -1);
860 forAll(addedToNew, i)
861 {
862 const label newFacei =
863 faceProcAddressing[meshi][fvp.start()+i];
864 const label patchFacei = newFacei-fvp0.start();
865 if
866 (
867 patchFacei >= 0
868 && patchFacei < fvp0.size()
869 )
870 {
871 addedToNew[i] = patchFacei;
872 }
873 }
874
875 bfld0[newPatchi].rmap(bfld[oldPatchi], addedToNew);
876 }
877 else
878 {
879 WarningInFunction << "Ignoring old patch "
880 << bfld[oldPatchi].patch().name() << " on field "
881 << flds[meshi].name() << endl; //exit(FatalError);
882 }
883 }
884 }
885 }
886}
887
888
889template<class Type>
890void Foam::fvMeshAdder::MapSurfaceField
891(
893 const labelList& oldFaceOwner0,
894 const labelList& oldPatchStarts0,
895 const labelList& oldPatchSizes0,
896 const labelListList& patchProcAddressing,
897 const labelListList& cellProcAddressing,
898 const labelListList& faceProcAddressing,
899 const bool fullyMapped
900)
901{
902 // Add fields to fields[0] after adding the meshes to meshes[0].
903 // Mesh[0] is the sum of all meshes. Fields are not yet mapped.
904
905 if (flds.size() == 0 || !flds.set(0))
906 {
907 FatalErrorInFunction << "Not valid field at element 0"
908 << " in field list of size " << flds.size() << exit(FatalError);
909 }
910
911 const fvMesh& mesh0 = flds[0].mesh();
912
913
914 // Internal field
915 // ~~~~~~~~~~~~~~
916
917 {
918 // Store old internal field
919 const Field<Type> oldInternalField(flds[0].primitiveField());
920
921 // Modify internal field
922 Field<Type>& intFld = flds[0].primitiveFieldRef();
923
924 // Set to new mesh size
925 intFld.setSize(mesh0.nInternalFaces());
926
927 // Map
928 forAll(flds, meshi)
929 {
930 if (flds.set(meshi))
931 {
932 const labelList& faceMap = faceProcAddressing[meshi];
933 const auto& fld = flds[meshi];
934
935 // Map internal field
936 if (meshi == 0)
937 {
938 intFld.rmap(oldInternalField, faceMap);
939 }
940 else
941 {
942 intFld.rmap(fld.primitiveField(), faceMap);
943 }
944
945 // Map faces that were boundary faces but are not anymore.
946 // There will be two meshes that provide the same face. Use
947 // owner one.
948 const auto& bfld = flds[meshi].boundaryField();
949
950 forAll(bfld, oldPatchi)
951 {
952 const fvsPatchField<Type>& pf = bfld[oldPatchi];
953 //Pout<< "patch:" << pf.patch().name() << endl;
954 forAll(pf, patchFacei)
955 {
956 // Get new face, mapped face owner
957 label newFacei;
958 label newOwn;
959 if (meshi == 0)
960 {
961 // Do not access mesh information since in-place
962 // modified
963 const label oldFacei =
964 oldPatchStarts0[oldPatchi]+patchFacei;
965 newFacei = faceProcAddressing[meshi][oldFacei];
966 const label oldOwn = oldFaceOwner0[oldFacei];
967 newOwn = cellProcAddressing[meshi][oldOwn];
968
969 //Pout<< "MESH0: pfi:" << patchFacei
970 // << " old face:" << oldFacei
971 // << " new face:" << newFacei
972 // << " at:" << mesh0.faceCentres()[newFacei]
973 // << " oldOwn:" << oldOwn
974 // << " newOwn:" << newOwn << endl;
975 }
976 else
977 {
978 const label oldFacei =
979 pf.patch().start()+patchFacei;
980 newFacei = faceProcAddressing[meshi][oldFacei];
981 const label oldOwn =
982 fld.mesh().faceOwner()[oldFacei];
983 newOwn = cellProcAddressing[meshi][oldOwn];
984
985 //Pout<< "MESH:" << meshi << " pfi:" << patchFacei
986 // << " old face:" << oldFacei
987 // << " new face:" << newFacei
988 // << " at:" << mesh0.faceCentres()[newFacei]
989 // << " oldOwn:" << oldOwn
990 // << " newOwn:" << newOwn << endl;
991 }
992
993 if
994 (
995 newFacei >= 0
996 && newFacei < mesh0.nInternalFaces()
997 && (newOwn == mesh0.faceOwner()[newFacei])
998 )
999 {
1000 intFld[newFacei] = pf[patchFacei];
1001 }
1002 //else
1003 //{
1004 // Pout<< " ignoring pfi:" << patchFacei
1005 // << " value:" << pf[patchFacei]
1006 // << " since newFacei:" << newFacei
1007 // << " since newOwn:" << newOwn
1008 // << " own:" << mesh0.faceOwner()[newFacei]
1009 // << endl;
1010 //}
1011 }
1012 }
1013 }
1014 }
1015 }
1016
1017
1018 // Patch fields from old mesh
1019 // ~~~~~~~~~~~~~~~~~~~~~~~~~~
1020
1021 auto& bfld0 = flds[0].boundaryFieldRef();
1022 forAll(bfld0, patchi)
1023 {
1024 labelList newToOld
1025 (
1026 calcPatchMap
1027 (
1028 oldPatchStarts0[patchi],
1029 oldPatchSizes0[patchi],
1030 faceProcAddressing[0],
1031 bfld0[patchi].patch().patch(),
1032 -1 // unmapped value
1033 )
1034 );
1035
1036 directFvPatchFieldMapper patchMapper(newToOld);
1037
1038 // Override mapping (for use in e.g. fvMeshDistribute where
1039 // it sorts mapping out itself)
1040 if (fullyMapped)
1041 {
1042 patchMapper.hasUnmapped() = false;
1043 }
1044
1045 bfld0[patchi].autoMap(patchMapper);
1046 }
1047
1048 for (label meshi = 1; meshi < flds.size(); meshi++)
1049 {
1050 if (flds.set(meshi))
1051 {
1052 const auto& bfld = flds[meshi].boundaryFieldRef();
1053
1054 const labelList& patchMap = patchProcAddressing[meshi];
1055
1056 forAll(patchMap, oldPatchi)
1057 {
1058 const auto& fvp = bfld[oldPatchi].patch();
1059 const label newPatchi = patchMap[oldPatchi];
1060 if (newPatchi >= 0 && newPatchi < bfld0.size())
1061 {
1062 const auto& fvp0 = bfld0[newPatchi].patch();
1063 labelList addedToNew(bfld[oldPatchi].size(), -1);
1064 forAll(addedToNew, i)
1065 {
1066 const label newFacei =
1067 faceProcAddressing[meshi][fvp.start()+i];
1068 const label patchFacei = newFacei-fvp0.start();
1069 if
1070 (
1071 patchFacei >= 0
1072 && patchFacei < fvp0.size()
1073 )
1074 {
1075 addedToNew[i] = patchFacei;
1076 }
1077 }
1078
1079 bfld0[newPatchi].rmap(bfld[oldPatchi], addedToNew);
1080 }
1081 else
1082 {
1083 WarningInFunction << "Ignoring old patch "
1084 << bfld[oldPatchi].patch().name() << " on field "
1085 << flds[meshi].name() << endl; //exit(FatalError);
1086 }
1087 }
1088 }
1089 }
1090}
1091
1092
1093template<class Type>
1095(
1096 const UPtrList<fvMesh>& meshes,
1097 const labelList& oldPatchStarts0,
1098 const labelList& oldPatchSizes0,
1099 const labelListList& patchProcAddressing,
1100 const labelListList& cellProcAddressing,
1101 const labelListList& faceProcAddressing,
1102 const labelListList& pointProcAddressing,
1103 const bool fullyMapped
1104)
1105{
1107
1108 if (meshes.size() == 0 || !meshes.set(0))
1109 {
1110 FatalErrorInFunction << "Not valid field at element 0"
1111 << " in mesh list of size " << meshes.size() << exit(FatalError);
1112 }
1113 const fvMesh& mesh0 = meshes[0];
1114
1116 (
1117 mesh0.objectRegistry::lookupClass<fldType>()
1118 );
1119
1120
1121 // It is necessary to enforce that all old-time fields are stored
1122 // before the mapping is performed. Otherwise, if the
1123 // old-time-level field is mapped before the field itself, sizes
1124 // will not match.
1125
1126 for (const auto& fld : fields)
1127 {
1128 DebugPout
1129 << "MapVolFields : Storing old time for " << fld->name()
1130 << endl;
1131
1132 const_cast<fldType&>(*fld).storeOldTimes();
1133 }
1134
1135
1136 for (const auto& fld : fields)
1137 {
1138 const word& name0 = fld->name();
1139
1140 DebugPout
1141 << "MapVolFields : mapping " << name0 << endl;
1142
1143 UPtrList<fldType> meshToField(meshes.size());
1144 forAll(meshes, meshi)
1145 {
1146 if (meshes.set(meshi))
1147 {
1148 auto& meshFld = meshes[meshi].
1149 objectRegistry::lookupObjectRef<fldType>(name0);
1150 meshToField.set(meshi, &meshFld);
1151 }
1152 }
1153
1154 MapVolField
1155 (
1156 meshToField,
1157 oldPatchStarts0,
1158 oldPatchSizes0,
1159 patchProcAddressing,
1160 cellProcAddressing,
1161 faceProcAddressing,
1162 fullyMapped
1163 );
1164 }
1165}
1166
1167
1168template<class Type>
1170(
1171 const UPtrList<fvMesh>& meshes,
1172 const labelListList& cellProcAddressing,
1173 const bool fullyMapped
1174)
1175{
1176 typedef DimensionedField<Type, volMesh> fldType;
1178
1179 if (meshes.size() == 0 || !meshes.set(0))
1180 {
1181 FatalErrorInFunction << "Not valid field at element 0"
1182 << " in mesh list of size " << meshes.size() << exit(FatalError);
1183 }
1184 const fvMesh& mesh0 = meshes[0];
1185
1187 (
1188 mesh0.objectRegistry::lookupClass<fldType>()
1189 );
1190
1191
1192 for (const auto& fld : fields)
1193 {
1194 if (!isA<excludeType>(*fld))
1195 {
1196 const word& name0 = fld->name();
1197
1198 DebugPout
1199 << "MapDimFields : mapping " << name0 << endl;
1200
1201 UPtrList<fldType> meshToField(meshes.size());
1202 forAll(meshes, meshi)
1203 {
1204 if (meshes.set(meshi))
1205 {
1206 auto& meshFld = meshes[meshi].
1207 objectRegistry::lookupObjectRef<fldType>(name0);
1208 meshToField.set(meshi, &meshFld);
1209 }
1210 }
1211
1212 MapDimField(meshToField, cellProcAddressing, fullyMapped);
1213 }
1214 else
1215 {
1216 DebugPout
1217 << "MapDimFields : ignoring " << fld->name() << endl;
1218 }
1219 }
1220}
1221
1222
1223template<class Type>
1225(
1226 const UPtrList<fvMesh>& meshes,
1227 const labelList& oldFaceOwner0,
1228 const labelList& oldPatchStarts0,
1229 const labelList& oldPatchSizes0,
1230 const labelListList& patchProcAddressing,
1231 const labelListList& cellProcAddressing,
1232 const labelListList& faceProcAddressing,
1233 const labelListList& pointProcAddressing,
1234 const bool fullyMapped
1235)
1236{
1238
1239 if (meshes.size() == 0 || !meshes.set(0))
1240 {
1241 FatalErrorInFunction << "Not valid field at element 0"
1242 << " in mesh list of size " << meshes.size() << exit(FatalError);
1243 }
1244 const auto& mesh0 = meshes[0];
1245
1247 (
1248 mesh0.objectRegistry::lookupClass<fldType>()
1249 );
1250
1251
1252 // It is necessary to enforce that all old-time fields are stored
1253 // before the mapping is performed. Otherwise, if the
1254 // old-time-level field is mapped before the field itself, sizes
1255 // will not match.
1256
1257 for (const auto& fld : fields)
1258 {
1259 DebugPout
1260 << "MapSurfaceFields : Storing old time for " << fld->name()
1261 << endl;
1262
1263 const_cast<fldType&>(*fld).storeOldTimes();
1264 }
1265
1266
1267 for (const auto& fld : fields)
1268 {
1269 const word& name0 = fld->name();
1270
1271 DebugPout
1272 << "MapSurfaceFields : Mapping " << fld->name() << endl;
1273
1274 UPtrList<fldType> meshToField(meshes.size());
1275 forAll(meshes, meshi)
1276 {
1277 if (meshes.set(meshi))
1278 {
1279 auto& meshFld = meshes[meshi].
1280 objectRegistry::lookupObjectRef<fldType>(name0);
1281 meshToField.set(meshi, &meshFld);
1282 }
1283 }
1284
1285 MapSurfaceField
1286 (
1287 meshToField,
1288 oldFaceOwner0,
1289 oldPatchStarts0,
1290 oldPatchSizes0,
1291 patchProcAddressing,
1292 cellProcAddressing,
1293 faceProcAddressing,
1294 fullyMapped
1295 );
1296 }
1297}
1298
1299
1300// ************************************************************************* //
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Mesh & mesh() const
Return mesh.
A templated direct mapper for the given FieldMapper type.
virtual bool hasUnmapped() const
Any unmapped values?
Generic templated field type.
Definition: Field.H:82
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:466
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A HashTable similar to std::unordered_map.
Definition: HashTable.H:123
bool found(const Key &key) const
Return true if hashed entry is found in table.
Definition: HashTableI.H:100
void setSize(const label n)
Alias for resize()
Definition: List.H:218
const T * set(const label i) const
Definition: PtrList.H:138
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A list of pointers to objects of type <T>, without allocation/deallocation management of the pointers...
Definition: UPtrList.H:71
const T * set(const label i) const
Definition: UPtrList.H:248
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
static void MapSurfaceFields(const mapAddedPolyMesh &, const fvMesh &mesh, const fvMesh &meshToAdd, const bool fullyMapped=false)
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, const bool fullyMapped=false)
Map all volFields of Type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:712
virtual label start() const
Return start label of this patch in the polyMesh face list.
Definition: fvPatch.H:179
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:79
const fvPatch & patch() const
Return patch.
Class containing mesh-to-mesh mapping information after a mesh addition where we add a mesh ('added m...
const labelList & oldPatchStarts() const
Return list of the old patch start labels.
const labelList & oldCellMap() const
const labelList & addedCellMap() const
const labelList & addedPatchMap() const
From added mesh patch index to new patch index or -1 if.
const labelList & oldPatchMap() const
From old patch index to new patch index or -1 if patch.
const labelList & oldFaceMap() const
const labelList & oldPatchSizes() const
Return list of the old patch sizes.
const labelList & addedFaceMap() const
const polyMesh & mesh() const noexcept
Return the mesh reference.
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1121
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:456
label nInternalFaces() const noexcept
Number of internal faces.
label nCells() const noexcept
Number of mesh cells.
A class for handling words, derived from Foam::string.
Definition: word.H:68
dynamicFvMesh & mesh
Foam::PtrList< Foam::fvMesh > meshes(regionNames.size())
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
#define WarningInFunction
Report a warning using Foam::Warning.
#define DebugPout
Report an information message using Foam::Pout.
DirectFieldMapper< fvPatchFieldMapper > directFvPatchFieldMapper
A fvPatchFieldMapper with direct mapping.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
List< label > labelList
A List of labels.
Definition: List.H:66
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
multivariateSurfaceInterpolationScheme< scalar >::fieldTable fields
Definition: createFields.H:97
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllIters(container, iter)
Iterate across all elements in the container object.
Definition: stdFoam.H:260
Foam::surfaceFields.