1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*----------------------------------------------------------------------------*/
27 #include "fvMeshDistribute.H"
28 #include "ProcessorTopology.H"
29 #include "commSchedule.H"
30 #include "PstreamCombineReduceOps.H"
31 #include "fvMeshAdder.H"
32 #include "faceCoupleInfo.H"
33 #include "processorFvPatchField.H"
34 #include "polyTopoChange.H"
35 #include "removeCells.H"
36 #include "polyModifyFace.H"
37 #include "polyRemovePoint.H"
38 #include "mergePoints.H"
39 #include "mapDistributePolyMesh.H"
40 #include "surfaceFields.H"
42 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
44 defineTypeNameAndDebug(Foam::fvMeshDistribute, 0);
47 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
49 //Foam::List<Foam::labelPair> Foam::fvMeshDistribute::getSchedule
51 // const labelList& distribution
54 // labelList nCellsPerProc(countCells(distribution));
58 // Pout<< "getSchedule : Wanted distribution:" << nCellsPerProc << endl;
61 // // Processors I need to send data to
62 // labelListList mySendProcs(Pstream::nProcs());
65 // label nSendProcs = 0;
66 // forAll(nCellsPerProc, sendProcI)
68 // if (sendProcI != Pstream::myProcNo() && nCellsPerProc[sendProcI] > 0)
75 // mySendProcs[Pstream::myProcNo()].setSize(nSendProcs);
77 // forAll(nCellsPerProc, sendProcI)
79 // if (sendProcI != Pstream::myProcNo() && nCellsPerProc[sendProcI] > 0)
81 // mySendProcs[Pstream::myProcNo()][nSendProcs++] = sendProcI;
86 // Pstream::gatherList(mySendProcs);
87 // Pstream::scatterList(mySendProcs);
89 // // Combine into list (same on all procs) giving sending and receiving
92 // forAll(mySendProcs, procI)
94 // nComms += mySendProcs[procI].size();
97 // List<labelPair> schedule(nComms);
100 // forAll(mySendProcs, procI)
102 // const labelList& sendProcs = mySendProcs[procI];
104 // forAll(sendProcs, i)
106 // schedule[nComms++] = labelPair(procI, sendProcs[i]);
114 Foam::labelList Foam::fvMeshDistribute::select
116 const bool selectEqual,
117 const labelList& values,
125 if (selectEqual == (values[i] == value))
131 labelList indices(n);
136 if (selectEqual == (values[i] == value))
145 // Check all procs have same names and in exactly same order.
146 void Foam::fvMeshDistribute::checkEqualWordList(const wordList& lst)
148 wordList myObjects(lst);
150 // Check that all meshes have same objects
151 Pstream::listCombineGather(myObjects, checkEqualType());
152 // Below scatter only needed to balance sends and receives.
153 Pstream::listCombineScatter(myObjects);
157 Foam::wordList Foam::fvMeshDistribute::mergeWordList(const wordList& procNames)
159 List<wordList> allNames(Pstream::nProcs());
160 allNames[Pstream::myProcNo()] = procNames;
161 Pstream::gatherList(allNames);
162 Pstream::scatterList(allNames);
164 HashSet<word> mergedNames;
165 forAll(allNames, procI)
167 forAll(allNames[procI], i)
169 mergedNames.insert(allNames[procI][i]);
172 return mergedNames.toc();
176 // Print some info on mesh.
177 void Foam::fvMeshDistribute::printMeshInfo(const fvMesh& mesh)
179 Pout<< "Primitives:" << nl
180 << " points :" << mesh.nPoints() << nl
181 << " internalFaces:" << mesh.nInternalFaces() << nl
182 << " faces :" << mesh.nFaces() << nl
183 << " cells :" << mesh.nCells() << nl;
185 const fvBoundaryMesh& patches = mesh.boundary();
187 Pout<< "Patches:" << endl;
188 forAll(patches, patchI)
190 const polyPatch& pp = patches[patchI].patch();
192 Pout<< " " << patchI << " name:" << pp.name()
193 << " size:" << pp.size()
194 << " start:" << pp.start()
195 << " type:" << pp.type()
199 if (mesh.pointZones().size() > 0)
201 Pout<< "PointZones:" << endl;
202 forAll(mesh.pointZones(), zoneI)
204 const pointZone& pz = mesh.pointZones()[zoneI];
205 Pout<< " " << zoneI << " name:" << pz.name()
206 << " size:" << pz.size()
210 if (mesh.faceZones().size() > 0)
212 Pout<< "FaceZones:" << endl;
213 forAll(mesh.faceZones(), zoneI)
215 const faceZone& fz = mesh.faceZones()[zoneI];
216 Pout<< " " << zoneI << " name:" << fz.name()
217 << " size:" << fz.size()
221 if (mesh.cellZones().size() > 0)
223 Pout<< "CellZones:" << endl;
224 forAll(mesh.cellZones(), zoneI)
226 const cellZone& cz = mesh.cellZones()[zoneI];
227 Pout<< " " << zoneI << " name:" << cz.name()
228 << " size:" << cz.size()
235 void Foam::fvMeshDistribute::printCoupleInfo
237 const primitiveMesh& mesh,
238 const labelList& sourceFace,
239 const labelList& sourceProc,
240 const labelList& sourceNewProc
244 << "Current coupling info:"
247 forAll(sourceFace, bFaceI)
249 label meshFaceI = mesh.nInternalFaces() + bFaceI;
251 Pout<< " meshFace:" << meshFaceI
252 << " fc:" << mesh.faceCentres()[meshFaceI]
253 << " connects to proc:" << sourceProc[bFaceI]
254 << "/face:" << sourceFace[bFaceI]
255 << " which will move to proc:" << sourceNewProc[bFaceI]
261 // Finds (non-empty) patch that exposed internal and proc faces can be put into.
262 Foam::label Foam::fvMeshDistribute::findNonEmptyPatch() const
264 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
266 label nonEmptyPatchI = -1;
268 forAllReverse(patches, patchI)
270 const polyPatch& pp = patches[patchI];
272 if (!isA<emptyPolyPatch>(pp) && !pp.coupled())
274 nonEmptyPatchI = patchI;
279 if (nonEmptyPatchI == -1)
281 FatalErrorIn("fvMeshDistribute::findNonEmptyPatch() const")
282 << "Cannot find a patch which is neither of type empty nor"
283 << " coupled in patches " << patches.names() << endl
284 << "There has to be at least one such patch for"
285 << " distribution to work" << abort(FatalError);
290 Pout<< "findNonEmptyPatch : using patch " << nonEmptyPatchI
291 << " name:" << patches[nonEmptyPatchI].name()
292 << " type:" << patches[nonEmptyPatchI].type()
293 << " to put exposed faces into." << endl;
297 // Do additional test for processor patches intermingled with non-proc
299 label procPatchI = -1;
301 forAll(patches, patchI)
303 if (isA<processorPolyPatch>(patches[patchI]))
307 else if (procPatchI != -1)
309 FatalErrorIn("fvMeshDistribute::findNonEmptyPatch() const")
310 << "Processor patches should be at end of patch list."
312 << "Have processor patch " << procPatchI
313 << " followed by non-processor patch " << patchI
314 << " in patches " << patches.names()
315 << abort(FatalError);
319 return nonEmptyPatchI;
323 // Appends processorPolyPatch. Returns patchID.
324 Foam::label Foam::fvMeshDistribute::addProcPatch
326 const word& patchName,
330 // Clear local fields and e.g. polyMesh globalMeshData.
334 polyBoundaryMesh& polyPatches =
335 const_cast<polyBoundaryMesh&>(mesh_.boundaryMesh());
336 fvBoundaryMesh& fvPatches = const_cast<fvBoundaryMesh&>(mesh_.boundary());
338 if (polyPatches.findPatchID(patchName) != -1)
340 FatalErrorIn("fvMeshDistribute::addProcPatch(const word&, const label)")
341 << "Cannot create patch " << patchName << " since already exists."
343 << "Current patch names:" << polyPatches.names()
352 label sz = polyPatches.size();
355 polyPatches.setSize(sz+1);
359 new processorPolyPatch
365 mesh_.boundaryMesh(),
370 fvPatches.setSize(sz+1);
376 polyPatches[sz], // point to newly added polyPatch
385 // Deletes last patch
386 void Foam::fvMeshDistribute::deleteTrailingPatch()
388 // Clear local fields and e.g. polyMesh globalMeshData.
391 polyBoundaryMesh& polyPatches =
392 const_cast<polyBoundaryMesh&>(mesh_.boundaryMesh());
393 fvBoundaryMesh& fvPatches = const_cast<fvBoundaryMesh&>(mesh_.boundary());
395 if (polyPatches.size() == 0)
397 FatalErrorIn("fvMeshDistribute::deleteTrailingPatch(fvMesh&)")
398 << "No patches in mesh"
399 << abort(FatalError);
402 label sz = polyPatches.size();
404 label nFaces = polyPatches[sz-1].size();
406 // Remove last polyPatch
409 Pout<< "deleteTrailingPatch : Removing patch " << sz-1
410 << " : " << polyPatches[sz-1].name() << " size:" << nFaces << endl;
415 FatalErrorIn("fvMeshDistribute::deleteTrailingPatch()")
416 << "There are still " << nFaces << " faces in patch to be deleted "
417 << sz-1 << ' ' << polyPatches[sz-1].name()
418 << abort(FatalError);
421 // Remove actual patch
422 polyPatches.setSize(sz-1);
423 fvPatches.setSize(sz-1);
427 // Delete all processor patches. Move any processor faces into the last
428 // non-processor patch.
429 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::deleteProcPatches
431 const label destinationPatch
434 // New patchID per boundary faces to be repatched. Is -1 (no change)
436 labelList newPatchID(mesh_.nFaces() - mesh_.nInternalFaces(), -1);
438 label nProcPatches = 0;
440 forAll(mesh_.boundaryMesh(), patchI)
442 const polyPatch& pp = mesh_.boundaryMesh()[patchI];
444 if (isA<processorPolyPatch>(pp))
448 Pout<< "Moving all faces of patch " << pp.name()
449 << " into patch " << destinationPatch
453 label offset = pp.start() - mesh_.nInternalFaces();
457 newPatchID[offset+i] = destinationPatch;
464 // Note: order of boundary faces been kept the same since the
465 // destinationPatch is at the end and we have visited the patches in
466 // incremental order.
467 labelListList dummyFaceMaps;
468 autoPtr<mapPolyMesh> map = repatch(newPatchID, dummyFaceMaps);
471 // Delete (now empty) processor patches.
472 forAllReverse(mesh_.boundaryMesh(), patchI)
474 const polyPatch& pp = mesh_.boundaryMesh()[patchI];
476 if (isA<processorPolyPatch>(pp))
478 deleteTrailingPatch();
479 deleteTrailingPatchFields<volScalarField>();
480 deleteTrailingPatchFields<volVectorField>();
481 deleteTrailingPatchFields<volSphericalTensorField>();
482 deleteTrailingPatchFields<volSymmTensorField>();
483 deleteTrailingPatchFields<volTensorField>();
485 deleteTrailingPatchFields<surfaceScalarField>();
486 deleteTrailingPatchFields<surfaceVectorField>();
487 deleteTrailingPatchFields<surfaceSphericalTensorField>();
488 deleteTrailingPatchFields<surfaceSymmTensorField>();
489 deleteTrailingPatchFields<surfaceTensorField>();
498 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::repatch
500 const labelList& newPatchID, // per boundary face -1 or new patchID
501 labelListList& constructFaceMap
504 polyTopoChange meshMod(mesh_);
506 forAll(newPatchID, bFaceI)
508 if (newPatchID[bFaceI] != -1)
510 label faceI = mesh_.nInternalFaces() + bFaceI;
512 label zoneID = mesh_.faceZones().whichZone(faceI);
513 bool zoneFlip = false;
517 const faceZone& fZone = mesh_.faceZones()[zoneID];
518 zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
525 mesh_.faces()[faceI], // modified face
526 faceI, // label of face
527 mesh_.faceOwner()[faceI], // owner
530 newPatchID[bFaceI], // patch for face
531 false, // remove from zone
532 zoneID, // zone for face
533 zoneFlip // face flip in zone
540 // Do mapping of fields from one patchField to the other ourselves since
541 // is currently not supported by updateMesh.
543 // Store boundary fields (we only do this for surfaceFields)
544 PtrList<FieldField<fvsPatchField, scalar> > sFlds;
545 saveBoundaryFields<scalar, surfaceMesh>(sFlds);
546 PtrList<FieldField<fvsPatchField, vector> > vFlds;
547 saveBoundaryFields<vector, surfaceMesh>(vFlds);
548 PtrList<FieldField<fvsPatchField, sphericalTensor> > sptFlds;
549 saveBoundaryFields<sphericalTensor, surfaceMesh>(sptFlds);
550 PtrList<FieldField<fvsPatchField, symmTensor> > sytFlds;
551 saveBoundaryFields<symmTensor, surfaceMesh>(sytFlds);
552 PtrList<FieldField<fvsPatchField, tensor> > tFlds;
553 saveBoundaryFields<tensor, surfaceMesh>(tFlds);
555 // Change the mesh (no inflation). Note: parallel comms allowed.
556 autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, true);
558 // Update fields. No inflation, parallel sync.
559 mesh_.updateMesh(map);
561 // Map patch fields using stored boundary fields. Note: assumes order
562 // of fields has not changed in object registry!
563 mapBoundaryFields<scalar, surfaceMesh>(map, sFlds);
564 mapBoundaryFields<vector, surfaceMesh>(map, vFlds);
565 mapBoundaryFields<sphericalTensor, surfaceMesh>(map, sptFlds);
566 mapBoundaryFields<symmTensor, surfaceMesh>(map, sytFlds);
567 mapBoundaryFields<tensor, surfaceMesh>(map, tFlds);
570 // Move mesh (since morphing does not do this)
571 if (map().hasMotionPoints())
573 mesh_.movePoints(map().preMotionPoints());
576 // Adapt constructMaps.
580 label index = findIndex(map().reverseFaceMap(), -1);
586 "fvMeshDistribute::repatch(const labelList&, labelListList&)"
587 ) << "reverseFaceMap contains -1 at index:"
589 << "This means that the repatch operation was not just a shuffle?"
590 << abort(FatalError);
594 forAll(constructFaceMap, procI)
596 inplaceRenumber(map().reverseFaceMap(), constructFaceMap[procI]);
604 // Detect shared points. Need processor patches to be present.
605 // Background: when adding bits of mesh one can get points which
606 // share the same position but are only detectable to be topologically
607 // the same point when doing parallel analysis. This routine will
608 // merge those points.
609 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::mergeSharedPoints
611 labelListList& constructPointMap
614 // Find out which sets of points get merged and create a map from
615 // mesh point to unique point.
616 Map<label> pointToMaster
618 fvMeshAdder::findSharedPoints
625 bool merged = pointToMaster.size() > 0;
627 if (!returnReduce(merged, orOp<bool>()))
629 return autoPtr<mapPolyMesh>(NULL);
632 polyTopoChange meshMod(mesh_);
634 fvMeshAdder::mergePoints(mesh_, pointToMaster, meshMod);
636 // Change the mesh (no inflation). Note: parallel comms allowed.
637 autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, true);
639 // Update fields. No inflation, parallel sync.
640 mesh_.updateMesh(map);
642 // Move mesh (since morphing does not do this)
643 if (map().hasMotionPoints())
645 mesh_.movePoints(map().preMotionPoints());
648 // Adapt constructMaps for merged points.
649 // 1.4.1: use reversePointMap < -1 feature.
650 forAll(constructPointMap, procI)
652 labelList& constructMap = constructPointMap[procI];
654 forAll(constructMap, i)
656 label oldPointI = constructMap[i];
658 // See if merged into other point
659 Map<label>::const_iterator iter = pointToMaster.find(oldPointI);
661 if (iter != pointToMaster.end())
666 constructMap[i] = map().reversePointMap()[oldPointI];
673 // Construct the local environment of all boundary faces.
674 void Foam::fvMeshDistribute::getNeighbourData
676 const labelList& distribution,
677 labelList& sourceFace,
678 labelList& sourceProc,
679 labelList& sourceNewProc
682 sourceFace.setSize(mesh_.nFaces() - mesh_.nInternalFaces());
683 sourceProc.setSize(sourceFace.size());
684 sourceNewProc.setSize(sourceFace.size());
686 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
688 // Send meshFace labels of processor patches and destination processor.
689 forAll(patches, patchI)
691 const polyPatch& pp = patches[patchI];
693 if (isA<processorPolyPatch>(pp))
695 const processorPolyPatch& procPatch =
696 refCast<const processorPolyPatch>(pp);
698 // Labels of faces on this side
699 labelList meshFaceLabels(pp.size());
700 forAll(meshFaceLabels, i)
702 meshFaceLabels[i] = pp.start()+i;
705 // Which processor they will end up on
706 const labelList newProc
708 IndirectList<label>(distribution, pp.faceCells())
711 OPstream toNeighbour(Pstream::blocking, procPatch.neighbProcNo());
713 toNeighbour << meshFaceLabels << newProc;
717 // Receive meshFace labels and destination processors of processor faces.
718 forAll(patches, patchI)
720 const polyPatch& pp = patches[patchI];
722 label offset = pp.start() - mesh_.nInternalFaces();
724 if (isA<processorPolyPatch>(pp))
726 const processorPolyPatch& procPatch =
727 refCast<const processorPolyPatch>(pp);
730 IPstream fromNeighbour(Pstream::blocking, procPatch.neighbProcNo());
731 labelList nbrFaces(fromNeighbour);
732 labelList nbrNewProc(fromNeighbour);
734 // Check which of the two faces we store.
736 if (Pstream::myProcNo() < procPatch.neighbProcNo())
738 // Use my local face labels
741 sourceFace[offset + i] = pp.start()+i;
742 sourceProc[offset + i] = Pstream::myProcNo();
743 sourceNewProc[offset + i] = nbrNewProc[i];
748 // Use my neighbours face labels
751 sourceFace[offset + i] = nbrFaces[i];
752 sourceProc[offset + i] = procPatch.neighbProcNo();
753 sourceNewProc[offset + i] = nbrNewProc[i];
759 // Normal (physical) boundary
762 sourceFace[offset + i] = patchI;
763 sourceProc[offset + i] = -1;
764 sourceNewProc[offset + i] = -1;
771 // Subset the neighbourCell/neighbourProc fields
772 void Foam::fvMeshDistribute::subsetBoundaryData
775 const labelList& faceMap,
776 const labelList& cellMap,
778 const labelList& oldDistribution,
779 const labelList& oldFaceOwner,
780 const labelList& oldFaceNeighbour,
781 const label oldInternalFaces,
783 const labelList& sourceFace,
784 const labelList& sourceProc,
785 const labelList& sourceNewProc,
789 labelList& subNewProc
792 subFace.setSize(mesh.nFaces() - mesh.nInternalFaces());
793 subProc.setSize(mesh.nFaces() - mesh.nInternalFaces());
794 subNewProc.setSize(mesh.nFaces() - mesh.nInternalFaces());
796 forAll(subFace, newBFaceI)
798 label newFaceI = newBFaceI + mesh.nInternalFaces();
800 label oldFaceI = faceMap[newFaceI];
802 // Was oldFaceI internal face? If so which side did we get.
803 if (oldFaceI < oldInternalFaces)
805 subFace[newBFaceI] = oldFaceI;
806 subProc[newBFaceI] = Pstream::myProcNo();
808 label oldOwn = oldFaceOwner[oldFaceI];
809 label oldNei = oldFaceNeighbour[oldFaceI];
811 if (oldOwn == cellMap[mesh.faceOwner()[newFaceI]])
813 // We kept the owner side. Where does the neighbour move to?
814 subNewProc[newBFaceI] = oldDistribution[oldNei];
818 // We kept the neighbour side.
819 subNewProc[newBFaceI] = oldDistribution[oldOwn];
824 // Was boundary face. Take over boundary information
825 label oldBFaceI = oldFaceI - oldInternalFaces;
827 subFace[newBFaceI] = sourceFace[oldBFaceI];
828 subProc[newBFaceI] = sourceProc[oldBFaceI];
829 subNewProc[newBFaceI] = sourceNewProc[oldBFaceI];
835 // Find cells on mesh whose faceID/procID match the neighbour cell/proc of
836 // domainMesh. Store the matching face.
837 void Foam::fvMeshDistribute::findCouples
839 const primitiveMesh& mesh,
840 const labelList& sourceFace,
841 const labelList& sourceProc,
844 const primitiveMesh& domainMesh,
845 const labelList& domainFace,
846 const labelList& domainProc,
848 labelList& masterCoupledFaces,
849 labelList& slaveCoupledFaces
852 // Store domain neighbour as map so we can easily look for pair
853 // with same face+proc.
854 HashTable<label, labelPair, labelPairHash> map(domainFace.size());
856 forAll(domainFace, bFaceI)
858 map.insert(labelPair(domainFace[bFaceI], domainProc[bFaceI]), bFaceI);
862 // Try to match mesh data.
864 masterCoupledFaces.setSize(domainFace.size());
865 slaveCoupledFaces.setSize(domainFace.size());
868 forAll(sourceFace, bFaceI)
870 if (sourceProc[bFaceI] != -1)
872 labelPair myData(sourceFace[bFaceI], sourceProc[bFaceI]);
874 HashTable<label, labelPair, labelPairHash>::const_iterator iter =
877 if (iter != map.end())
879 label nbrBFaceI = iter();
881 masterCoupledFaces[coupledI] = mesh.nInternalFaces() + bFaceI;
882 slaveCoupledFaces[coupledI] =
883 domainMesh.nInternalFaces()
891 masterCoupledFaces.setSize(coupledI);
892 slaveCoupledFaces.setSize(coupledI);
896 Pout<< "findCouples : found " << coupledI
897 << " faces that will be stitched" << nl << endl;
902 // Map data on boundary faces to new mesh (resulting from adding two meshes)
903 Foam::labelList Foam::fvMeshDistribute::mapBoundaryData
905 const primitiveMesh& mesh, // mesh after adding
906 const mapAddedPolyMesh& map,
907 const labelList& boundaryData0, // mesh before adding
908 const label nInternalFaces1,
909 const labelList& boundaryData1 // added mesh
912 labelList newBoundaryData(mesh.nFaces() - mesh.nInternalFaces());
914 forAll(boundaryData0, oldBFaceI)
916 label newFaceI = map.oldFaceMap()[oldBFaceI + map.nOldInternalFaces()];
918 // Face still exists (is nessecary?) and still boundary face
919 if (newFaceI >= 0 && newFaceI >= mesh.nInternalFaces())
921 newBoundaryData[newFaceI - mesh.nInternalFaces()] =
922 boundaryData0[oldBFaceI];
926 forAll(boundaryData1, addedBFaceI)
928 label newFaceI = map.addedFaceMap()[addedBFaceI + nInternalFaces1];
930 if (newFaceI >= 0 && newFaceI >= mesh.nInternalFaces())
932 newBoundaryData[newFaceI - mesh.nInternalFaces()] =
933 boundaryData1[addedBFaceI];
937 return newBoundaryData;
941 // Remove cells. Add all exposed faces to patch oldInternalPatchI
942 Foam::autoPtr<Foam::mapPolyMesh> Foam::fvMeshDistribute::doRemoveCells
944 const labelList& cellsToRemove,
945 const label oldInternalPatchI
948 // Mesh change engine
949 polyTopoChange meshMod(mesh_);
951 // Cell removal topo engine. Do NOT synchronize parallel since
952 // we are doing a local cell removal.
953 removeCells cellRemover(mesh_, false);
955 // Get all exposed faces
956 labelList exposedFaces(cellRemover.getExposedFaces(cellsToRemove));
958 // Insert the topo changes
959 cellRemover.setRefinement
963 labelList(exposedFaces.size(), oldInternalPatchI), // patch for exposed
968 // Change the mesh. No inflation. Note: no parallel comms allowed.
969 autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh_, false, false);
972 mesh_.updateMesh(map);
974 // Move mesh (since morphing does not do this)
975 if (map().hasMotionPoints())
977 mesh_.movePoints(map().preMotionPoints());
984 // Delete and add processor patches. Changes mesh and returns per neighbour proc
985 // the processor patchID.
986 void Foam::fvMeshDistribute::addProcPatches
988 const labelList& neighbourNewProc, // processor that neighbour is on
989 labelList& procPatchID
992 // Now use the neighbourFace/Proc to repatch the mesh. These two lists
993 // contain for all current boundary faces the global patchID (for non-proc
994 // patch) or the processor.
996 labelList procPatchSizes(Pstream::nProcs(), 0);
998 forAll(neighbourNewProc, bFaceI)
1000 if (neighbourNewProc[bFaceI] != -1)
1002 procPatchSizes[neighbourNewProc[bFaceI]]++;
1006 // Per neighbour processor the label of the processor patch
1007 procPatchID.setSize(Pstream::nProcs());
1009 forAll(procPatchSizes, procI)
1011 if (procPatchSizes[procI] > 0)
1013 const word patchName =
1015 + name(Pstream::myProcNo())
1020 procPatchID[procI] = addProcPatch(patchName, procI);
1021 addPatchFields<volScalarField>
1023 processorFvPatchField<scalar>::typeName
1025 addPatchFields<volVectorField>
1027 processorFvPatchField<vector>::typeName
1029 addPatchFields<volSphericalTensorField>
1031 processorFvPatchField<sphericalTensor>::typeName
1033 addPatchFields<volSymmTensorField>
1035 processorFvPatchField<symmTensor>::typeName
1037 addPatchFields<volTensorField>
1039 processorFvPatchField<tensor>::typeName
1042 addPatchFields<surfaceScalarField>
1044 processorFvPatchField<scalar>::typeName
1046 addPatchFields<surfaceVectorField>
1048 processorFvPatchField<vector>::typeName
1050 addPatchFields<surfaceSphericalTensorField>
1052 processorFvPatchField<sphericalTensor>::typeName
1054 addPatchFields<surfaceSymmTensorField>
1056 processorFvPatchField<symmTensor>::typeName
1058 addPatchFields<surfaceTensorField>
1060 processorFvPatchField<tensor>::typeName
1065 procPatchID[procI] = -1;
1071 // Get boundary faces to be repatched. Is -1 or new patchID
1072 Foam::labelList Foam::fvMeshDistribute::getProcBoundaryPatch
1074 const labelList& neighbourNewProc, // new processor per boundary face
1075 const labelList& procPatchID // patchID
1078 labelList patchIDs(neighbourNewProc);
1080 forAll(neighbourNewProc, bFaceI)
1082 if (neighbourNewProc[bFaceI] != -1)
1084 label nbrProc = neighbourNewProc[bFaceI];
1086 patchIDs[bFaceI] = procPatchID[nbrProc];
1090 patchIDs[bFaceI] = -1;
1097 // Send mesh and coupling data.
1098 void Foam::fvMeshDistribute::sendMesh
1103 const wordList& pointZoneNames,
1104 const wordList& faceZoneNames,
1105 const wordList& cellZoneNames,
1107 const labelList& sourceFace,
1108 const labelList& sourceProc,
1109 const labelList& sourceNewProc
1114 Pout<< "Sending to domain " << domain << nl
1115 << " nPoints:" << mesh.nPoints() << nl
1116 << " nFaces:" << mesh.nFaces() << nl
1117 << " nCells:" << mesh.nCells() << nl
1118 << " nPatches:" << mesh.boundaryMesh().size() << nl
1122 // Assume sparse point zones. Get contents in merged-zone indices.
1123 labelListList zonePoints(pointZoneNames.size());
1125 const pointZoneMesh& pointZones = mesh.pointZones();
1127 forAll(pointZoneNames, nameI)
1129 label myZoneID = pointZones.findZoneID(pointZoneNames[nameI]);
1133 zonePoints[nameI] = pointZones[myZoneID];
1138 // Assume sparse face zones
1139 labelListList zoneFaces(faceZoneNames.size());
1140 boolListList zoneFaceFlip(faceZoneNames.size());
1142 const faceZoneMesh& faceZones = mesh.faceZones();
1144 forAll(faceZoneNames, nameI)
1146 label myZoneID = faceZones.findZoneID(faceZoneNames[nameI]);
1150 zoneFaces[nameI] = faceZones[myZoneID];
1151 zoneFaceFlip[nameI] = faceZones[myZoneID].flipMap();
1156 // Assume sparse cell zones
1157 labelListList zoneCells(cellZoneNames.size());
1159 const cellZoneMesh& cellZones = mesh.cellZones();
1161 forAll(cellZoneNames, nameI)
1163 label myZoneID = cellZones.findZoneID(cellZoneNames[nameI]);
1167 zoneCells[nameI] = cellZones[myZoneID];
1171 ////- Assume full cell zones
1172 //labelList cellZoneID;
1175 // cellZoneID.setSize(mesh.nCells());;
1178 // const cellZoneMesh& cellZones = mesh.cellZones();
1180 // forAll(cellZones, zoneI)
1182 // IndirectList<label>(cellZoneID, cellZones[zoneI]) = zoneI;
1187 OPstream toDomain(Pstream::blocking, domain);
1192 << mesh.faceNeighbour()
1193 << mesh.boundaryMesh()
1206 // Receive mesh. Opposite of sendMesh
1207 Foam::autoPtr<Foam::fvMesh> Foam::fvMeshDistribute::receiveMesh
1210 const wordList& pointZoneNames,
1211 const wordList& faceZoneNames,
1212 const wordList& cellZoneNames,
1213 const Time& runTime,
1214 labelList& domainSourceFace,
1215 labelList& domainSourceProc,
1216 labelList& domainSourceNewProc
1219 IPstream fromNbr(Pstream::blocking, domain);
1221 pointField domainPoints(fromNbr);
1222 faceList domainFaces(fromNbr);
1223 labelList domainAllOwner(fromNbr);
1224 labelList domainAllNeighbour(fromNbr);
1225 PtrList<entry> patchEntries(fromNbr);
1227 labelListList zonePoints(fromNbr);
1228 labelListList zoneFaces(fromNbr);
1229 boolListList zoneFaceFlip(fromNbr);
1230 labelListList zoneCells(fromNbr);
1235 >> domainSourceNewProc;
1238 autoPtr<fvMesh> domainMeshPtr
1244 fvMesh::defaultRegion,
1253 false // no parallel comms
1256 fvMesh& domainMesh = domainMeshPtr();
1258 List<polyPatch*> patches(patchEntries.size());
1260 forAll(patchEntries, patchI)
1262 patches[patchI] = polyPatch::New
1264 patchEntries[patchI].keyword(),
1265 patchEntries[patchI].dict(),
1267 domainMesh.boundaryMesh()
1270 // Add patches; no parallel comms
1271 domainMesh.addFvPatches(patches, false);
1275 List<pointZone*> pZonePtrs(pointZoneNames.size());
1276 forAll(pZonePtrs, i)
1278 pZonePtrs[i] = new pointZone
1283 domainMesh.pointZones()
1287 List<faceZone*> fZonePtrs(faceZoneNames.size());
1288 forAll(fZonePtrs, i)
1290 fZonePtrs[i] = new faceZone
1296 domainMesh.faceZones()
1300 List<cellZone*> cZonePtrs(cellZoneNames.size());
1301 forAll(cZonePtrs, i)
1303 cZonePtrs[i] = new cellZone
1308 domainMesh.cellZones()
1311 domainMesh.addZones(pZonePtrs, fZonePtrs, cZonePtrs);
1313 return domainMeshPtr;
1317 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
1319 // Construct from components
1320 Foam::fvMeshDistribute::fvMeshDistribute(fvMesh& mesh, const scalar mergeTol)
1327 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
1329 Foam::labelList Foam::fvMeshDistribute::countCells
1331 const labelList& distribution
1334 labelList nCells(Pstream::nProcs(), 0);
1335 forAll(distribution, cellI)
1337 label newProc = distribution[cellI];
1339 if (newProc < 0 || newProc >= Pstream::nProcs())
1341 FatalErrorIn("fvMeshDistribute::distribute(const labelList&)")
1342 << "Distribution should be in range 0.." << Pstream::nProcs()-1
1344 << "At index " << cellI << " distribution:" << newProc
1345 << abort(FatalError);
1353 Foam::autoPtr<Foam::mapDistributePolyMesh> Foam::fvMeshDistribute::distribute
1355 const labelList& distribution
1358 // Some checks on distribution
1359 if (distribution.size() != mesh_.nCells())
1361 FatalErrorIn("fvMeshDistribute::distribute(const labelList&)")
1362 << "Size of distribution:"
1363 << distribution.size() << " mesh nCells:" << mesh_.nCells()
1364 << abort(FatalError);
1368 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
1370 // Check all processors have same non-proc patches in same order.
1371 if (patches.checkParallelSync(true))
1373 FatalErrorIn("fvMeshDistribute::distribute(const labelList&)")
1374 << "This application requires all non-processor patches"
1375 << " to be present in the same order on all patches" << nl
1376 << "followed by the processor patches (which of course are unique)."
1378 << "Local patches:" << mesh_.boundaryMesh().names()
1379 << abort(FatalError);
1382 // Save some data for mapping later on
1383 const label nOldPoints(mesh_.nPoints());
1384 const label nOldFaces(mesh_.nFaces());
1385 const label nOldCells(mesh_.nCells());
1386 labelList oldPatchStarts(patches.size());
1387 labelList oldPatchNMeshPoints(patches.size());
1388 forAll(patches, patchI)
1390 oldPatchStarts[patchI] = patches[patchI].start();
1391 oldPatchNMeshPoints[patchI] = patches[patchI].nPoints();
1396 // Short circuit trivial case.
1397 if (!Pstream::parRun())
1399 // Collect all maps and return
1400 return autoPtr<mapDistributePolyMesh>
1402 new mapDistributePolyMesh
1410 oldPatchNMeshPoints,
1412 labelListList(1, identity(mesh_.nPoints())),//subPointMap
1413 labelListList(1, identity(mesh_.nFaces())), //subFaceMap
1414 labelListList(1, identity(mesh_.nCells())), //subCellMap
1415 labelListList(1, identity(patches.size())), //subPatchMap
1417 labelListList(1, identity(mesh_.nPoints())),//constructPointMap
1418 labelListList(1, identity(mesh_.nFaces())), //constructFaceMap
1419 labelListList(1, identity(mesh_.nCells())), //constructCellMap
1420 labelListList(1, identity(patches.size())) //constructPatchMap
1426 // Collect any zone names
1427 const wordList pointZoneNames(mergeWordList(mesh_.pointZones().names()));
1428 const wordList faceZoneNames(mergeWordList(mesh_.faceZones().names()));
1429 const wordList cellZoneNames(mergeWordList(mesh_.cellZones().names()));
1432 // Local environment of all boundary faces
1433 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1435 // A face is uniquely defined by
1439 // To glue the parts of meshes which can get sent from anywhere we
1440 // need to know on boundary faces what the above tuple on both sides is.
1441 // So we need to maintain:
1443 // - original processor id (= trivial)
1444 // For coupled boundaries (where the faces are 'duplicate') we take the
1445 // lowest numbered processor as the data to store.
1447 // Additionally to create the procboundaries we need to know where the owner
1448 // cell on the other side now is: newNeighbourProc.
1451 // physical boundary:
1453 // sourceNewProc = -1
1454 // sourceFace = patchID
1455 // coupled boundary:
1456 // sourceProc = proc
1457 // sourceNewProc = distribution[cell on proc]
1458 // sourceFace = face
1459 labelList sourceFace;
1460 labelList sourceProc;
1461 labelList sourceNewProc;
1462 getNeighbourData(distribution, sourceFace, sourceProc, sourceNewProc);
1465 // Remove meshPhi. Since this would otherwise dissappear anyway
1466 // during topo changes and we have to guarantee that all the fields
1469 mesh_.resetMotion();
1471 const wordList volScalars(mesh_.names(volScalarField::typeName));
1472 checkEqualWordList(volScalars);
1473 const wordList volVectors(mesh_.names(volVectorField::typeName));
1474 checkEqualWordList(volVectors);
1475 const wordList volSphereTensors
1477 mesh_.names(volSphericalTensorField::typeName)
1479 checkEqualWordList(volSphereTensors);
1480 const wordList volSymmTensors(mesh_.names(volSymmTensorField::typeName));
1481 checkEqualWordList(volSymmTensors);
1482 const wordList volTensors(mesh_.names(volTensorField::typeName));
1483 checkEqualWordList(volTensors);
1485 const wordList surfScalars(mesh_.names(surfaceScalarField::typeName));
1486 checkEqualWordList(surfScalars);
1487 const wordList surfVectors(mesh_.names(surfaceVectorField::typeName));
1488 checkEqualWordList(surfVectors);
1489 const wordList surfSphereTensors
1491 mesh_.names(surfaceSphericalTensorField::typeName)
1493 checkEqualWordList(surfSphereTensors);
1494 const wordList surfSymmTensors
1496 mesh_.names(surfaceSymmTensorField::typeName)
1498 checkEqualWordList(surfSymmTensors);
1499 const wordList surfTensors(mesh_.names(surfaceTensorField::typeName));
1500 checkEqualWordList(surfTensors);
1503 // Find patch to temporarily put exposed and processor faces into.
1504 label oldInternalPatchI = findNonEmptyPatch();
1508 // Delete processor patches, starting from the back. Move all faces into
1509 // oldInternalPatchI.
1510 labelList repatchFaceMap;
1512 autoPtr<mapPolyMesh> repatchMap = deleteProcPatches(oldInternalPatchI);
1514 // Store face map (only face ordering that changed)
1515 repatchFaceMap = repatchMap().faceMap();
1517 // Reorder all boundary face data (sourceProc, sourceFace etc.)
1522 repatchMap().reverseFaceMap(),
1523 mesh_.nFaces() - mesh_.nInternalFaces(),
1524 mesh_.nInternalFaces()
1526 - mesh_.nInternalFaces()
1529 inplaceReorder(bFaceMap, sourceFace);
1530 inplaceReorder(bFaceMap, sourceProc);
1531 inplaceReorder(bFaceMap, sourceNewProc);
1539 Pout<< nl << "MESH WITH PROC PATCHES DELETED:" << endl;
1540 printMeshInfo(mesh_);
1541 printFieldInfo<volScalarField>(mesh_);
1542 printFieldInfo<volVectorField>(mesh_);
1543 printFieldInfo<volSphericalTensorField>(mesh_);
1544 printFieldInfo<volSymmTensorField>(mesh_);
1545 printFieldInfo<volTensorField>(mesh_);
1546 printFieldInfo<surfaceScalarField>(mesh_);
1547 printFieldInfo<surfaceVectorField>(mesh_);
1548 printFieldInfo<surfaceSphericalTensorField>(mesh_);
1549 printFieldInfo<surfaceSymmTensorField>(mesh_);
1550 printFieldInfo<surfaceTensorField>(mesh_);
1556 // Maps from subsetted mesh (that is sent) back to original maps
1557 labelListList subCellMap(Pstream::nProcs());
1558 labelListList subFaceMap(Pstream::nProcs());
1559 labelListList subPointMap(Pstream::nProcs());
1560 labelListList subPatchMap(Pstream::nProcs());
1561 // Maps from subsetted mesh to reconstructed mesh
1562 labelListList constructCellMap(Pstream::nProcs());
1563 labelListList constructFaceMap(Pstream::nProcs());
1564 labelListList constructPointMap(Pstream::nProcs());
1565 labelListList constructPatchMap(Pstream::nProcs());
1570 // Find out schedule
1571 // ~~~~~~~~~~~~~~~~~
1573 labelListList nSendCells(Pstream::nProcs());
1574 nSendCells[Pstream::myProcNo()] = countCells(distribution);
1575 Pstream::gatherList(nSendCells);
1576 Pstream::scatterList(nSendCells);
1580 // What to send to neighbouring domains
1581 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1583 forAll(nSendCells[Pstream::myProcNo()], recvProc)
1587 recvProc != Pstream::myProcNo()
1588 && nSendCells[Pstream::myProcNo()][recvProc] > 0
1596 << "SUBSETTING FOR DOMAIN " << recvProc
1597 << " cells to send:"
1598 << nSendCells[Pstream::myProcNo()][recvProc]
1602 // Mesh subsetting engine
1603 fvMeshSubset subsetter(mesh_);
1605 // Subset the cells of the current domain.
1606 subsetter.setLargeCellSubset
1610 oldInternalPatchI, // oldInternalFaces patch
1611 false // no parallel sync
1614 subCellMap[recvProc] = subsetter.cellMap();
1615 subFaceMap[recvProc] = renumber
1620 subPointMap[recvProc] = subsetter.pointMap();
1621 subPatchMap[recvProc] = subsetter.patchMap();
1624 // Subset the boundary fields (owner/neighbour/processor)
1625 labelList procSourceFace;
1626 labelList procSourceProc;
1627 labelList procSourceNewProc;
1631 subsetter.subMesh(),
1632 subsetter.faceMap(), // from subMesh to mesh
1633 subsetter.cellMap(), // ,, ,,
1635 distribution, // old mesh distribution
1636 mesh_.faceOwner(), // old owner
1637 mesh_.faceNeighbour(),
1638 mesh_.nInternalFaces(),
1649 // Send to neighbour
1653 subsetter.subMesh(),
1663 sendFields<volScalarField>(recvProc, volScalars, subsetter);
1664 sendFields<volVectorField>(recvProc, volVectors, subsetter);
1665 sendFields<volSphericalTensorField>
1671 sendFields<volSymmTensorField>
1677 sendFields<volTensorField>(recvProc, volTensors, subsetter);
1679 sendFields<surfaceScalarField>(recvProc, surfScalars, subsetter);
1680 sendFields<surfaceVectorField>(recvProc, surfVectors, subsetter);
1681 sendFields<surfaceSphericalTensorField>
1687 sendFields<surfaceSymmTensorField>
1693 sendFields<surfaceTensorField>(recvProc, surfTensors, subsetter);
1699 // Subset the part that stays
1700 // ~~~~~~~~~~~~~~~~~~~~~~~~~~
1703 // Save old mesh maps before changing mesh
1704 const labelList oldFaceOwner(mesh_.faceOwner());
1705 const labelList oldFaceNeighbour(mesh_.faceNeighbour());
1706 const label oldInternalFaces = mesh_.nInternalFaces();
1709 autoPtr<mapPolyMesh> subMap
1713 select(false, distribution, Pstream::myProcNo()),
1718 // Addressing from subsetted mesh
1719 subCellMap[Pstream::myProcNo()] = subMap().cellMap();
1720 subFaceMap[Pstream::myProcNo()] = renumber
1725 subPointMap[Pstream::myProcNo()] = subMap().pointMap();
1726 subPatchMap[Pstream::myProcNo()] = identity(patches.size());
1728 // Initialize all addressing into current mesh
1729 constructCellMap[Pstream::myProcNo()] = identity(mesh_.nCells());
1730 constructFaceMap[Pstream::myProcNo()] = identity(mesh_.nFaces());
1731 constructPointMap[Pstream::myProcNo()] = identity(mesh_.nPoints());
1732 constructPatchMap[Pstream::myProcNo()] = identity(patches.size());
1734 // Subset the mesh data: neighbourCell/neighbourProc
1736 labelList domainSourceFace;
1737 labelList domainSourceProc;
1738 labelList domainSourceNewProc;
1743 subMap().faceMap(), // from new to original mesh
1746 distribution, // distribution before subsetting
1747 oldFaceOwner, // owner before subsetting
1748 oldFaceNeighbour, // neighbour ,,
1749 oldInternalFaces, // nInternalFaces ,,
1760 sourceFace.transfer(domainSourceFace);
1761 sourceProc.transfer(domainSourceProc);
1762 sourceNewProc.transfer(domainSourceNewProc);
1769 Pout<< nl << "STARTING MESH:" << endl;
1770 printMeshInfo(mesh_);
1771 printFieldInfo<volScalarField>(mesh_);
1772 printFieldInfo<volVectorField>(mesh_);
1773 printFieldInfo<volSphericalTensorField>(mesh_);
1774 printFieldInfo<volSymmTensorField>(mesh_);
1775 printFieldInfo<volTensorField>(mesh_);
1776 printFieldInfo<surfaceScalarField>(mesh_);
1777 printFieldInfo<surfaceVectorField>(mesh_);
1778 printFieldInfo<surfaceSphericalTensorField>(mesh_);
1779 printFieldInfo<surfaceSymmTensorField>(mesh_);
1780 printFieldInfo<surfaceTensorField>(mesh_);
1786 // Receive and add what was sent
1787 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1789 forAll(nSendCells, sendProc)
1791 // Did processor sendProc send anything to me?
1794 sendProc != Pstream::myProcNo()
1795 && nSendCells[sendProc][Pstream::myProcNo()] > 0
1801 << "RECEIVING FROM DOMAIN " << sendProc
1802 << " cells to receive:"
1803 << nSendCells[sendProc][Pstream::myProcNo()]
1807 // Receive from sendProc
1808 labelList domainSourceFace;
1809 labelList domainSourceProc;
1810 labelList domainSourceNewProc;
1812 // Opposite of sendMesh
1813 autoPtr<fvMesh> domainMeshPtr = receiveMesh
1820 const_cast<Time&>(mesh_.time()),
1825 fvMesh& domainMesh = domainMeshPtr();
1828 PtrList<volScalarField> vsf;
1829 receiveFields<volScalarField>
1837 PtrList<volVectorField> vvf;
1838 receiveFields<volVectorField>
1845 PtrList<volSphericalTensorField> vsptf;
1846 receiveFields<volSphericalTensorField>
1853 PtrList<volSymmTensorField> vsytf;
1854 receiveFields<volSymmTensorField>
1861 PtrList<volTensorField> vtf;
1862 receiveFields<volTensorField>
1870 PtrList<surfaceScalarField> ssf;
1871 receiveFields<surfaceScalarField>
1878 PtrList<surfaceVectorField> svf;
1879 receiveFields<surfaceVectorField>
1886 PtrList<surfaceSphericalTensorField> ssptf;
1887 receiveFields<surfaceSphericalTensorField>
1894 PtrList<surfaceSymmTensorField> ssytf;
1895 receiveFields<surfaceSymmTensorField>
1902 PtrList<surfaceTensorField> stf;
1903 receiveFields<surfaceTensorField>
1912 constructCellMap[sendProc] = identity(domainMesh.nCells());
1913 constructFaceMap[sendProc] = identity(domainMesh.nFaces());
1914 constructPointMap[sendProc] = identity(domainMesh.nPoints());
1915 constructPatchMap[sendProc] =
1916 identity(domainMesh.boundaryMesh().size());
1922 Pout<< nl << "RECEIVED MESH FROM:" << sendProc << endl;
1923 printMeshInfo(domainMesh);
1924 printFieldInfo<volScalarField>(domainMesh);
1925 printFieldInfo<volVectorField>(domainMesh);
1926 printFieldInfo<volSphericalTensorField>(domainMesh);
1927 printFieldInfo<volSymmTensorField>(domainMesh);
1928 printFieldInfo<volTensorField>(domainMesh);
1929 printFieldInfo<surfaceScalarField>(domainMesh);
1930 printFieldInfo<surfaceVectorField>(domainMesh);
1931 printFieldInfo<surfaceSphericalTensorField>(domainMesh);
1932 printFieldInfo<surfaceSymmTensorField>(domainMesh);
1933 printFieldInfo<surfaceTensorField>(domainMesh);
1937 // Now this mesh we received (from sendProc) needs to be merged
1938 // with the current mesh. On the current mesh we have for all
1939 // boundaryfaces the original face and processor. See if we can
1940 // match these up to the received domainSourceFace and
1941 // domainSourceProc.
1942 labelList masterCoupledFaces;
1943 labelList slaveCoupledFaces;
1960 // Generate additional coupling info (points, edges) from
1962 faceCoupleInfo couples
1968 mergeTol_, // merge tolerance
1969 true, // faces align
1970 true, // couples are ordered already
1975 // Add domainMesh to mesh
1976 // ~~~~~~~~~~~~~~~~~~~~~~
1978 autoPtr<mapAddedPolyMesh> map = fvMeshAdder::add
1983 false // no parallel comms
1986 // Update mesh data: sourceFace,sourceProc for added
1995 domainMesh.nInternalFaces(),
2004 domainMesh.nInternalFaces(),
2013 domainMesh.nInternalFaces(),
2017 // Update all addressing so xxProcAddressing points to correct item
2019 const labelList& oldCellMap = map().oldCellMap();
2020 const labelList& oldFaceMap = map().oldFaceMap();
2021 const labelList& oldPointMap = map().oldPointMap();
2022 const labelList& oldPatchMap = map().oldPatchMap();
2024 forAll(constructPatchMap, procI)
2026 if (procI != sendProc && constructPatchMap[procI].size() > 0)
2028 // Processor already in mesh (either myProcNo or received)
2029 inplaceRenumber(oldCellMap, constructCellMap[procI]);
2030 inplaceRenumber(oldFaceMap, constructFaceMap[procI]);
2031 inplaceRenumber(oldPointMap, constructPointMap[procI]);
2032 inplaceRenumber(oldPatchMap, constructPatchMap[procI]);
2037 inplaceRenumber(map().addedCellMap(), constructCellMap[sendProc]);
2038 inplaceRenumber(map().addedFaceMap(), constructFaceMap[sendProc]);
2039 inplaceRenumber(map().addedPointMap(), constructPointMap[sendProc]);
2040 inplaceRenumber(map().addedPatchMap(), constructPatchMap[sendProc]);
2044 Pout<< nl << "MERGED MESH FROM:" << sendProc << endl;
2045 printMeshInfo(mesh_);
2046 printFieldInfo<volScalarField>(mesh_);
2047 printFieldInfo<volVectorField>(mesh_);
2048 printFieldInfo<volSphericalTensorField>(mesh_);
2049 printFieldInfo<volSymmTensorField>(mesh_);
2050 printFieldInfo<volTensorField>(mesh_);
2051 printFieldInfo<surfaceScalarField>(mesh_);
2052 printFieldInfo<surfaceVectorField>(mesh_);
2053 printFieldInfo<surfaceSphericalTensorField>(mesh_);
2054 printFieldInfo<surfaceSymmTensorField>(mesh_);
2055 printFieldInfo<surfaceTensorField>(mesh_);
2065 Pout<< nl << "REDISTRIBUTED MESH:" << endl;
2066 printMeshInfo(mesh_);
2067 printFieldInfo<volScalarField>(mesh_);
2068 printFieldInfo<volVectorField>(mesh_);
2069 printFieldInfo<volSphericalTensorField>(mesh_);
2070 printFieldInfo<volSymmTensorField>(mesh_);
2071 printFieldInfo<volTensorField>(mesh_);
2072 printFieldInfo<surfaceScalarField>(mesh_);
2073 printFieldInfo<surfaceVectorField>(mesh_);
2074 printFieldInfo<surfaceSphericalTensorField>(mesh_);
2075 printFieldInfo<surfaceSymmTensorField>(mesh_);
2076 printFieldInfo<surfaceTensorField>(mesh_);
2082 // Add processorPatches
2083 // ~~~~~~~~~~~~~~~~~~~~
2085 // Per neighbour processor the patchID to it (or -1).
2086 labelList procPatchID;
2088 // Add processor patches.
2089 addProcPatches(sourceNewProc, procPatchID);
2091 // Put faces into correct patch. Note that we now have proper
2092 // processorPolyPatches again so repatching will take care of coupled face
2095 // Get boundary faces to be repatched. Is -1 or new patchID
2096 labelList newPatchID
2098 getProcBoundaryPatch
2105 // Change patches. Since this might change ordering of coupled faces
2106 // we also need to adapt our constructMaps.
2107 repatch(newPatchID, constructFaceMap);
2109 // See if any geometrically shared points need to be merged. Note: does
2111 mergeSharedPoints(constructPointMap);
2113 // Bit of hack: processorFvPatchField does not get reset since created
2114 // from nothing so explicitly reset.
2115 initPatchFields<volScalarField>
2117 processorFvPatchField<scalar>::typeName,
2118 pTraits<scalar>::zero
2120 initPatchFields<volVectorField>
2122 processorFvPatchField<vector>::typeName,
2123 pTraits<vector>::zero
2125 initPatchFields<volSphericalTensorField>
2127 processorFvPatchField<sphericalTensor>::typeName,
2128 pTraits<sphericalTensor>::zero
2130 initPatchFields<volSymmTensorField>
2132 processorFvPatchField<symmTensor>::typeName,
2133 pTraits<symmTensor>::zero
2135 initPatchFields<volTensorField>
2137 processorFvPatchField<tensor>::typeName,
2138 pTraits<tensor>::zero
2140 initPatchFields<surfaceScalarField>
2142 processorFvPatchField<scalar>::typeName,
2143 pTraits<scalar>::zero
2145 initPatchFields<surfaceVectorField>
2147 processorFvPatchField<vector>::typeName,
2148 pTraits<vector>::zero
2150 initPatchFields<surfaceSphericalTensorField>
2152 processorFvPatchField<sphericalTensor>::typeName,
2153 pTraits<sphericalTensor>::zero
2155 initPatchFields<surfaceSymmTensorField>
2157 processorFvPatchField<symmTensor>::typeName,
2158 pTraits<symmTensor>::zero
2160 initPatchFields<surfaceTensorField>
2162 processorFvPatchField<tensor>::typeName,
2163 pTraits<tensor>::zero
2167 mesh_.setInstance(mesh_.time().timeName());
2173 Pout<< nl << "FINAL MESH:" << endl;
2174 printMeshInfo(mesh_);
2175 printFieldInfo<volScalarField>(mesh_);
2176 printFieldInfo<volVectorField>(mesh_);
2177 printFieldInfo<volSphericalTensorField>(mesh_);
2178 printFieldInfo<volSymmTensorField>(mesh_);
2179 printFieldInfo<volTensorField>(mesh_);
2180 printFieldInfo<surfaceScalarField>(mesh_);
2181 printFieldInfo<surfaceVectorField>(mesh_);
2182 printFieldInfo<surfaceSphericalTensorField>(mesh_);
2183 printFieldInfo<surfaceSymmTensorField>(mesh_);
2184 printFieldInfo<surfaceTensorField>(mesh_);
2188 // Collect all maps and return
2189 return autoPtr<mapDistributePolyMesh>
2191 new mapDistributePolyMesh
2199 oldPatchNMeshPoints,
2210 true // reuse storage
2216 // ************************************************************************* //