2 // Copyright (C) 2020 CEA/DEN, EDF R&D
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 // Author : Anthony Geay (EDF R&D)
22 #include "ParaUMesh.hxx"
23 #include "ProcessorGroup.hxx"
24 #include "MPIProcessorGroup.hxx"
25 #include "Topology.hxx"
26 #include "BlockTopology.hxx"
27 #include "CommInterface.hxx"
28 #include "MEDCouplingMemArray.hxx"
39 using namespace MEDCoupling;
41 ParaUMesh::ParaUMesh(MEDCouplingUMesh *mesh, DataArrayIdType *globalCellIds, DataArrayIdType *globalNodeIds)
44 _cell_global.takeRef(globalCellIds);
45 _node_global.takeRef(globalNodeIds);
47 _cell_global.checkNotNull();
48 _node_global.checkNotNull();
49 _mesh->checkConsistencyLight();
50 if(_mesh->getNumberOfNodes() != _node_global->getNumberOfTuples())
51 throw INTERP_KERNEL::Exception("ParaUMesh constructor : mismatch between # nodes and len of global # nodes.");
52 if(_mesh->getNumberOfCells() != _cell_global->getNumberOfTuples())
53 throw INTERP_KERNEL::Exception("ParaUMesh constructor : mismatch between # cells and len of global # cells.");
57 * This method computes the cells part of distributed mesh lying on \a globalNodeIds nodes.
58 * The input \a globalNodeIds are not supposed to reside on the current process.
60 MCAuto<DataArrayIdType> ParaUMesh::getCellIdsLyingOnNodes(const DataArrayIdType *globalNodeIds, bool fullyIn) const
63 throw INTERP_KERNEL::Exception("ParaUMesh::getCellIdsLyingOnNodes : not implemented yet for fullyIn == True !");
64 MPI_Comm comm(MPI_COMM_WORLD);
67 ci.commSize(comm,&size);
68 std::unique_ptr<mcIdType[]> nbOfElems(new mcIdType[size]),nbOfElems2(new mcIdType[size]),nbOfElems3(new mcIdType[size]);
69 mcIdType nbOfNodeIdsLoc(globalNodeIds->getNumberOfTuples());
70 ci.allGather(&nbOfNodeIdsLoc,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,comm);
71 std::vector< MCAuto<DataArrayIdType> > tabs(size);
72 // loop to avoid to all procs to have all the nodes per proc
73 for(int subDiv = 0 ; subDiv < size ; ++subDiv)
75 std::unique_ptr<mcIdType[]> nbOfElemsSp(CommInterface::SplitArrayOfLength(nbOfElems,size,subDiv,size));
76 mcIdType nbOfNodeIdsSum(std::accumulate(nbOfElemsSp.get(),nbOfElemsSp.get()+size,0));
77 std::unique_ptr<mcIdType[]> allGlobalNodeIds(new mcIdType[nbOfNodeIdsSum]);
78 std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElemsSp,size) );
79 std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) );
80 mcIdType startGlobalNodeIds,endGlobalNodeIds;
81 DataArray::GetSlice(0,globalNodeIds->getNumberOfTuples(),1,subDiv,size,startGlobalNodeIds,endGlobalNodeIds);
82 ci.allGatherV(globalNodeIds->begin()+startGlobalNodeIds,endGlobalNodeIds-startGlobalNodeIds,MPI_ID_TYPE,allGlobalNodeIds.get(),nbOfElemsInt.get(),offsetsIn.get(),MPI_ID_TYPE,comm);
84 for(int curRk = 0 ; curRk < size ; ++curRk)
86 MCAuto<DataArrayIdType> globalNodeIdsOfCurProc(DataArrayIdType::New());
87 globalNodeIdsOfCurProc->useArray(allGlobalNodeIds.get()+offset,false,DeallocType::CPP_DEALLOC,nbOfElemsSp[curRk],1);
88 offset += nbOfElemsSp[curRk];
89 MCAuto<DataArrayIdType> globalNodeIdsCaptured(_node_global->buildIntersection(globalNodeIdsOfCurProc));
90 MCAuto<DataArrayIdType> localNodeIdsToLocate(_node_global->findIdForEach(globalNodeIdsCaptured->begin(),globalNodeIdsCaptured->end()));
91 MCAuto<DataArrayIdType> localCellCaptured(_mesh->getCellIdsLyingOnNodes(localNodeIdsToLocate->begin(),localNodeIdsToLocate->end(),false));
92 MCAuto<DataArrayIdType> localCellCapturedGlob(_cell_global->selectByTupleIdSafe(localCellCaptured->begin(),localCellCaptured->end()));
93 if(tabs[curRk].isNull())
94 tabs[curRk] = localCellCapturedGlob;
96 tabs[curRk]->insertAtTheEnd(localCellCapturedGlob->begin(),localCellCapturedGlob->end());
99 for(int curRk = 0 ; curRk < size ; ++curRk)
101 tabs[curRk] = tabs[curRk]->buildUniqueNotSorted();
102 nbOfElems3[curRk] = tabs[curRk]->getNumberOfTuples();
104 std::vector<const DataArrayIdType *> tabss(tabs.begin(),tabs.end());
105 MCAuto<DataArrayIdType> cells(DataArrayIdType::Aggregate(tabss));
106 ci.allToAll(nbOfElems3.get(),1,MPI_ID_TYPE,nbOfElems2.get(),1,MPI_ID_TYPE,comm);
107 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems2.get(),nbOfElems2.get()+size,0));
108 MCAuto<DataArrayIdType> cellIdsFromProcs(DataArrayIdType::New());
109 cellIdsFromProcs->alloc(nbOfCellIdsSum,1);
111 std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElems3,size) ),nbOfElemsOutInt( CommInterface::ToIntArray<mcIdType>(nbOfElems2,size) );
112 std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) ), offsetsOut( CommInterface::ComputeOffset(nbOfElemsOutInt,size) );
113 ci.allToAllV(cells->begin(),nbOfElemsInt.get(),offsetsIn.get(),MPI_ID_TYPE,
114 cellIdsFromProcs->getPointer(),nbOfElemsOutInt.get(),offsetsOut.get(),MPI_ID_TYPE,comm);
116 cellIdsFromProcs->sort();
117 return cellIdsFromProcs;
122 MCAuto<ParaUMesh> ParaUMesh::redistributeCells(const DataArrayIdType *globalCellIds) const
124 MPI_Comm comm(MPI_COMM_WORLD);
127 ci.commSize(comm,&size);
128 std::unique_ptr<mcIdType[]> nbOfElems(new mcIdType[size]),nbOfElems2(new mcIdType[size]);
129 mcIdType nbOfNodeIdsLoc(globalCellIds->getNumberOfTuples());
130 ci.allGather(&nbOfNodeIdsLoc,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,comm);
131 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems.get(),nbOfElems.get()+size,0));
132 std::unique_ptr<mcIdType[]> allGlobalCellIds(new mcIdType[nbOfCellIdsSum]);
133 std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElems,size) );
134 std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) );
135 ci.allGatherV(globalCellIds->begin(),nbOfNodeIdsLoc,MPI_ID_TYPE,allGlobalCellIds.get(),nbOfElemsInt.get(),offsetsIn.get(),MPI_ID_TYPE,comm);
137 for(int curRk = 0 ; curRk < size ; ++curRk)
139 MCAuto<DataArrayIdType> globalCellIdsOfCurProc(DataArrayIdType::New());
140 globalCellIdsOfCurProc->useArray(allGlobalCellIds.get()+offset,false,DeallocType::CPP_DEALLOC,nbOfElems[curRk],1);
141 offset += nbOfElems[curRk];
142 // prepare all2all session
143 MCAuto<DataArrayIdType> globalCellIdsCaptured(_cell_global->buildIntersection(globalCellIdsOfCurProc));// OK for the global cellIds
144 MCAuto<DataArrayIdType> localCellIdsCaptured(_node_global->findIdForEach(globalCellIdsCaptured->begin(),globalCellIdsCaptured->end()));
145 MCAuto<MEDCouplingUMesh> meshPart(_mesh->buildPartOfMySelf(localCellIdsCaptured->begin(),localCellIdsCaptured->end(),true));
146 MCAuto<DataArrayIdType> o2n(meshPart->zipCoordsTraducer());// OK for the mesh
147 MCAuto<DataArrayIdType> n2o(o2n->invertArrayO2N2N2O(meshPart->getNumberOfNodes()));
148 MCAuto<DataArrayIdType> globalNodeIdsPart(_node_global->selectByTupleIdSafe(n2o->begin(),n2o->end())); // OK for the global nodeIds