]> SALOME platform Git repositories - tools/medcoupling.git/blobdiff - src/ParaMEDMEM/ParaUMesh.cxx
Salome HOME
refactor!: remove adm_local/ directory
[tools/medcoupling.git] / src / ParaMEDMEM / ParaUMesh.cxx
index 7e3060892ffd450b1d87173a3848f24cfb9c5716..b2502f9f3bbce22e9fe5e630dfc6b08eb0eaf79b 100644 (file)
@@ -1,5 +1,5 @@
 //
-// Copyright (C) 2020  CEA/DEN, EDF R&D
+// Copyright (C) 2020-2024  CEA, EDF
 //
 // This library is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
@@ -74,7 +74,13 @@ std::vector<const BigMemoryObject *> ParaUMesh::getDirectChildrenWithNull() cons
 MCAuto<DataArrayIdType> ParaUMesh::getCellIdsLyingOnNodes(const DataArrayIdType *globalNodeIds, bool fullyIn) const
 {
   if(fullyIn)
-    throw INTERP_KERNEL::Exception("ParaUMesh::getCellIdsLyingOnNodes : not implemented yet for fullyIn == True !");
+    return this->getCellIdsLyingOnNodesTrue(globalNodeIds);
+  else
+    return this->getCellIdsLyingOnNodesFalse(globalNodeIds);
+}
+
+MCAuto<DataArrayIdType> ParaUMesh::getCellIdsLyingOnNodesTrue(const DataArrayIdType *globalNodeIds) const
+{
   MPI_Comm comm(MPI_COMM_WORLD);
   CommInterface ci;
   int size;
@@ -83,17 +89,19 @@ MCAuto<DataArrayIdType> ParaUMesh::getCellIdsLyingOnNodes(const DataArrayIdType
   mcIdType nbOfNodeIdsLoc(globalNodeIds->getNumberOfTuples());
   ci.allGather(&nbOfNodeIdsLoc,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,comm);
   std::vector< MCAuto<DataArrayIdType> > tabs(size);
+  //store for each proc the local nodeids intercepted by current proc
+  int nbOfCollectiveCalls = 1;// this parameter controls the memory peak
   // loop to avoid to all procs to have all the nodes per proc
-  for(int subDiv = 0 ; subDiv < size ; ++subDiv)
+  for(int subDiv = 0 ; subDiv < nbOfCollectiveCalls ; ++subDiv)
   {
-    std::unique_ptr<mcIdType[]> nbOfElemsSp(CommInterface::SplitArrayOfLength(nbOfElems,size,subDiv,size));
+    std::unique_ptr<mcIdType[]> nbOfElemsSp(CommInterface::SplitArrayOfLength(nbOfElems,size,subDiv,nbOfCollectiveCalls));
     mcIdType nbOfNodeIdsSum(std::accumulate(nbOfElemsSp.get(),nbOfElemsSp.get()+size,0));
     std::unique_ptr<mcIdType[]> allGlobalNodeIds(new mcIdType[nbOfNodeIdsSum]);
     std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElemsSp,size) );
     std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) );
     mcIdType startGlobalNodeIds,endGlobalNodeIds;
-    DataArray::GetSlice(0,globalNodeIds->getNumberOfTuples(),1,subDiv,size,startGlobalNodeIds,endGlobalNodeIds);
-    ci.allGatherV(globalNodeIds->begin()+startGlobalNodeIds,endGlobalNodeIds-startGlobalNodeIds,MPI_ID_TYPE,allGlobalNodeIds.get(),nbOfElemsInt.get(),offsetsIn.get(),MPI_ID_TYPE,comm);
+    DataArray::GetSlice(0,globalNodeIds->getNumberOfTuples(),1,subDiv,nbOfCollectiveCalls,startGlobalNodeIds,endGlobalNodeIds);
+    ci.allGatherV(globalNodeIds->begin()+startGlobalNodeIds,FromIdType<int>(endGlobalNodeIds-startGlobalNodeIds),MPI_ID_TYPE,allGlobalNodeIds.get(),nbOfElemsInt.get(),offsetsIn.get(),MPI_ID_TYPE,comm);
     mcIdType offset(0);
     for(int curRk = 0 ; curRk < size ; ++curRk)
     {
@@ -102,14 +110,23 @@ MCAuto<DataArrayIdType> ParaUMesh::getCellIdsLyingOnNodes(const DataArrayIdType
       offset += nbOfElemsSp[curRk];
       MCAuto<DataArrayIdType> globalNodeIdsCaptured(_node_global->buildIntersection(globalNodeIdsOfCurProc));
       MCAuto<DataArrayIdType> localNodeIdsToLocate(_node_global->findIdForEach(globalNodeIdsCaptured->begin(),globalNodeIdsCaptured->end()));
-      MCAuto<DataArrayIdType> localCellCaptured(_mesh->getCellIdsLyingOnNodes(localNodeIdsToLocate->begin(),localNodeIdsToLocate->end(),false));
-      MCAuto<DataArrayIdType> localCellCapturedGlob(_cell_global->selectByTupleIdSafe(localCellCaptured->begin(),localCellCaptured->end()));
       if(tabs[curRk].isNull())
-        tabs[curRk] = localCellCapturedGlob;
+        tabs[curRk] = localNodeIdsToLocate;
       else
-        tabs[curRk]->insertAtTheEnd(localCellCapturedGlob->begin(),localCellCapturedGlob->end());
+        tabs[curRk]->insertAtTheEnd(localNodeIdsToLocate->begin(),localNodeIdsToLocate->end());
     }
   }
+
+  for(int curRk = 0 ; curRk < size ; ++curRk)
+  {
+    MCAuto<DataArrayIdType> localNodeIds(tabs[curRk]);
+    localNodeIds->sort();
+    MCAuto<DataArrayIdType> localNodeIdsUnique(localNodeIds->buildUnique());
+    MCAuto<DataArrayIdType> localCellCaptured(_mesh->getCellIdsLyingOnNodes(localNodeIdsUnique->begin(),localNodeIdsUnique->end(),true));
+    MCAuto<DataArrayIdType> localCellCapturedGlob(_cell_global->selectByTupleIdSafe(localCellCaptured->begin(),localCellCaptured->end()));
+    tabs[curRk] = localCellCapturedGlob;
+  }
+  
   for(int curRk = 0 ; curRk < size ; ++curRk)
   {
     tabs[curRk] = tabs[curRk]->buildUniqueNotSorted();
@@ -131,83 +148,83 @@ MCAuto<DataArrayIdType> ParaUMesh::getCellIdsLyingOnNodes(const DataArrayIdType
   return cellIdsFromProcs;
 }
 
-DataArrayIdType *ParaUMesh::redistributeCellField(const DataArrayIdType *globalCellIds, const DataArrayIdType *fieldValueToRed) const
+MCAuto<DataArrayIdType> ParaUMesh::getCellIdsLyingOnNodesFalse(const DataArrayIdType *globalNodeIds) const
 {
   MPI_Comm comm(MPI_COMM_WORLD);
   CommInterface ci;
-  std::unique_ptr<mcIdType[]> allGlobalCellIds,allGlobalCellIdsIndex;
-  int size(ci.allGatherArrays(comm,globalCellIds,allGlobalCellIds,allGlobalCellIdsIndex));
-  // Prepare ParaUMesh parts to be sent : compute for each proc the contribution of current rank.
-  std::vector< MCAuto<DataArrayIdType> > globalCellIdsToBeSent(size);
-  std::vector< MCAuto<DataArrayIdType> > fieldToBeSent(size);
+  int size;
+  ci.commSize(comm,&size);
+  std::unique_ptr<mcIdType[]> nbOfElems(new mcIdType[size]),nbOfElems2(new mcIdType[size]),nbOfElems3(new mcIdType[size]);
+  mcIdType nbOfNodeIdsLoc(globalNodeIds->getNumberOfTuples());
+  ci.allGather(&nbOfNodeIdsLoc,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,comm);
+  // loop to avoid to all procs to have all the nodes per proc
+  int nbOfCollectiveCalls = 1;// this parameter controls the memory peak
+  std::vector< MCAuto<DataArrayIdType> > tabs(size);
+  for(int subDiv = 0 ; subDiv < nbOfCollectiveCalls ; ++subDiv)
+  {
+    std::unique_ptr<mcIdType[]> nbOfElemsSp(CommInterface::SplitArrayOfLength(nbOfElems,size,subDiv,nbOfCollectiveCalls));
+    mcIdType nbOfNodeIdsSum(std::accumulate(nbOfElemsSp.get(),nbOfElemsSp.get()+size,0));
+    std::unique_ptr<mcIdType[]> allGlobalNodeIds(new mcIdType[nbOfNodeIdsSum]);
+    std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElemsSp,size) );
+    std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) );
+    mcIdType startGlobalNodeIds,endGlobalNodeIds;
+    DataArray::GetSlice(0,globalNodeIds->getNumberOfTuples(),1,subDiv,nbOfCollectiveCalls,startGlobalNodeIds,endGlobalNodeIds);
+    ci.allGatherV(globalNodeIds->begin()+startGlobalNodeIds,FromIdType<int>(endGlobalNodeIds-startGlobalNodeIds),MPI_ID_TYPE,allGlobalNodeIds.get(),nbOfElemsInt.get(),offsetsIn.get(),MPI_ID_TYPE,comm);
+    mcIdType offset(0);
+    for(int curRk = 0 ; curRk < size ; ++curRk)
+    {
+      MCAuto<DataArrayIdType> globalNodeIdsOfCurProc(DataArrayIdType::New());
+      globalNodeIdsOfCurProc->useArray(allGlobalNodeIds.get()+offset,false,DeallocType::CPP_DEALLOC,nbOfElemsSp[curRk],1);
+      offset += nbOfElemsSp[curRk];
+      MCAuto<DataArrayIdType> globalNodeIdsCaptured(_node_global->buildIntersection(globalNodeIdsOfCurProc));
+      MCAuto<DataArrayIdType> localNodeIdsToLocate(_node_global->findIdForEach(globalNodeIdsCaptured->begin(),globalNodeIdsCaptured->end()));
+      MCAuto<DataArrayIdType> localCellCaptured(_mesh->getCellIdsLyingOnNodes(localNodeIdsToLocate->begin(),localNodeIdsToLocate->end(),false));
+      MCAuto<DataArrayIdType> localCellCapturedGlob(_cell_global->selectByTupleIdSafe(localCellCaptured->begin(),localCellCaptured->end()));
+      if(tabs[curRk].isNull())
+        tabs[curRk] = localCellCapturedGlob;
+      else
+        tabs[curRk]->insertAtTheEnd(localCellCapturedGlob->begin(),localCellCapturedGlob->end());
+    }
+  }
   for(int curRk = 0 ; curRk < size ; ++curRk)
   {
-    mcIdType offset(allGlobalCellIdsIndex[curRk]);
-    MCAuto<DataArrayIdType> globalCellIdsOfCurProc(DataArrayIdType::New());
-    globalCellIdsOfCurProc->useArray(allGlobalCellIds.get()+offset,false,DeallocType::CPP_DEALLOC,allGlobalCellIdsIndex[curRk+1]-offset,1);
-    // the key call is here : compute for rank curRk the cells to be sent
-    MCAuto<DataArrayIdType> globalCellIdsCaptured(_cell_global->buildIntersection(globalCellIdsOfCurProc));// OK for the global cellIds
-    MCAuto<DataArrayIdType> localCellIdsCaptured(_cell_global->findIdForEach(globalCellIdsCaptured->begin(),globalCellIdsCaptured->end()));
-    globalCellIdsToBeSent[curRk] = globalCellIdsCaptured;
-    fieldToBeSent[curRk] = fieldValueToRed->selectByTupleIdSafe(localCellIdsCaptured->begin(),localCellIdsCaptured->end());
+    tabs[curRk] = tabs[curRk]->buildUniqueNotSorted();
+    nbOfElems3[curRk] = tabs[curRk]->getNumberOfTuples();
   }
-  // Receive
-  std::vector< MCAuto<DataArrayIdType> > globalCellIdsReceived;
-  ci.allToAllArrays(comm,globalCellIdsToBeSent,globalCellIdsReceived);
-  std::vector< MCAuto<DataArrayIdType> > fieldValueReceived;
-  ci.allToAllArrays(comm,fieldToBeSent,fieldValueReceived);
-  // use globalCellIdsReceived to reorganize everything
-  MCAuto<DataArrayIdType> aggregatedCellIds( DataArrayIdType::Aggregate(FromVecAutoToVecOfConst<DataArrayIdType>(globalCellIdsReceived)) );
-  MCAuto<DataArrayIdType> aggregatedCellIdsSorted(aggregatedCellIds->copySorted());
-  MCAuto<DataArrayIdType> idsIntoAggregatedIds(DataArrayIdType::FindPermutationFromFirstToSecondDuplicate(aggregatedCellIdsSorted,aggregatedCellIds));
-  MCAuto<DataArrayIdType> cellIdsOfSameNodeIds(aggregatedCellIdsSorted->indexOfSameConsecutiveValueGroups());
-  MCAuto<DataArrayIdType> n2o_cells(idsIntoAggregatedIds->selectByTupleIdSafe(cellIdsOfSameNodeIds->begin(),cellIdsOfSameNodeIds->end()-1));//new == new ordering so that global cell ids are sorted . old == coarse ordering implied by the aggregation
-  //
-  MCAuto<DataArrayIdType> fieldAggregated(DataArrayIdType::Aggregate(FromVecAutoToVecOfConst<DataArrayIdType>(fieldValueReceived)));
-  MCAuto<DataArrayIdType> ret(fieldAggregated->selectByTupleIdSafe(n2o_cells->begin(),n2o_cells->end()));
-  return ret.retn();
+  std::vector<const DataArrayIdType *> tabss(tabs.begin(),tabs.end());
+  MCAuto<DataArrayIdType> cells(DataArrayIdType::Aggregate(tabss));
+  ci.allToAll(nbOfElems3.get(),1,MPI_ID_TYPE,nbOfElems2.get(),1,MPI_ID_TYPE,comm);
+  mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems2.get(),nbOfElems2.get()+size,0));
+  MCAuto<DataArrayIdType> cellIdsFromProcs(DataArrayIdType::New());
+  cellIdsFromProcs->alloc(nbOfCellIdsSum,1);
+  {
+    std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElems3,size) ),nbOfElemsOutInt( CommInterface::ToIntArray<mcIdType>(nbOfElems2,size) );
+    std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) ), offsetsOut( CommInterface::ComputeOffset(nbOfElemsOutInt,size) );
+    ci.allToAllV(cells->begin(),nbOfElemsInt.get(),offsetsIn.get(),MPI_ID_TYPE,
+                 cellIdsFromProcs->getPointer(),nbOfElemsOutInt.get(),offsetsOut.get(),MPI_ID_TYPE,comm);
+  }
+  cellIdsFromProcs->sort();
+  return cellIdsFromProcs;
+}
+
+DataArrayIdType *ParaUMesh::redistributeCellField(const DataArrayIdType *globalCellIds, const DataArrayIdType *fieldValueToRed) const
+{
+  return this->redistributeCellFieldT<mcIdType>(globalCellIds,fieldValueToRed);
+}
+
+DataArrayDouble *ParaUMesh::redistributeCellField(const DataArrayIdType *globalCellIds, const DataArrayDouble *fieldValueToRed) const
+{
+  return this->redistributeCellFieldT<double>(globalCellIds,fieldValueToRed);
 }
 
 DataArrayIdType *ParaUMesh::redistributeNodeField(const DataArrayIdType *globalCellIds, const DataArrayIdType *fieldValueToRed) const
 {
-  MPI_Comm comm(MPI_COMM_WORLD);
-  CommInterface ci;
-  std::unique_ptr<mcIdType[]> allGlobalCellIds,allGlobalCellIdsIndex;
-  int size(ci.allGatherArrays(comm,globalCellIds,allGlobalCellIds,allGlobalCellIdsIndex));
-  // Prepare ParaUMesh parts to be sent : compute for each proc the contribution of current rank.
-  std::vector< MCAuto<DataArrayIdType> > globalNodeIdsToBeSent(size);
-  std::vector< MCAuto<DataArrayIdType> > fieldToBeSent(size);
-  for(int curRk = 0 ; curRk < size ; ++curRk)
-  {
-    mcIdType offset(allGlobalCellIdsIndex[curRk]);
-    MCAuto<DataArrayIdType> globalCellIdsOfCurProc(DataArrayIdType::New());
-    globalCellIdsOfCurProc->useArray(allGlobalCellIds.get()+offset,false,DeallocType::CPP_DEALLOC,allGlobalCellIdsIndex[curRk+1]-offset,1);
-    // the key call is here : compute for rank curRk the cells to be sent
-    MCAuto<DataArrayIdType> globalCellIdsCaptured(_cell_global->buildIntersection(globalCellIdsOfCurProc));// OK for the global cellIds
-    MCAuto<DataArrayIdType> localCellIdsCaptured(_cell_global->findIdForEach(globalCellIdsCaptured->begin(),globalCellIdsCaptured->end()));
-    MCAuto<MEDCouplingUMesh> meshPart(_mesh->buildPartOfMySelf(localCellIdsCaptured->begin(),localCellIdsCaptured->end(),true));
-    MCAuto<DataArrayIdType> o2n(meshPart->zipCoordsTraducer());// OK for the mesh
-    MCAuto<DataArrayIdType> n2o(o2n->invertArrayO2N2N2O(meshPart->getNumberOfNodes()));
-    MCAuto<DataArrayIdType> globalNodeIdsPart(_node_global->selectByTupleIdSafe(n2o->begin(),n2o->end())); // OK for the global nodeIds
-    globalNodeIdsToBeSent[curRk] = globalNodeIdsPart;
-    fieldToBeSent[curRk] = fieldValueToRed->selectByTupleIdSafe(n2o->begin(),n2o->end());
-  }
-  // Receive
-  std::vector< MCAuto<DataArrayIdType> > globalNodeIdsReceived;
-  ci.allToAllArrays(comm,globalNodeIdsToBeSent,globalNodeIdsReceived);
-  std::vector< MCAuto<DataArrayIdType> > fieldValueReceived;
-  ci.allToAllArrays(comm,fieldToBeSent,fieldValueReceived);
-  // firstly deal with nodes.
-  MCAuto<DataArrayIdType> aggregatedNodeIds( DataArrayIdType::Aggregate(FromVecAutoToVecOfConst<DataArrayIdType>(globalNodeIdsReceived)) );
-  MCAuto<DataArrayIdType> aggregatedNodeIdsSorted(aggregatedNodeIds->copySorted());
-  MCAuto<DataArrayIdType> nodeIdsIntoAggregatedIds(DataArrayIdType::FindPermutationFromFirstToSecondDuplicate(aggregatedNodeIdsSorted,aggregatedNodeIds));
-  MCAuto<DataArrayIdType> idxOfSameNodeIds(aggregatedNodeIdsSorted->indexOfSameConsecutiveValueGroups());
-  MCAuto<DataArrayIdType> n2o_nodes(nodeIdsIntoAggregatedIds->selectByTupleIdSafe(idxOfSameNodeIds->begin(),idxOfSameNodeIds->end()-1));//new == new ordering so that global node ids are sorted . old == coarse ordering implied by the aggregation
-  //
-  MCAuto<DataArrayIdType> fieldAggregated(DataArrayIdType::Aggregate(FromVecAutoToVecOfConst<DataArrayIdType>(fieldValueReceived)));
-  MCAuto<DataArrayIdType> ret(fieldAggregated->selectByTupleIdSafe(n2o_nodes->begin(),n2o_nodes->end()));
-  //
-  return ret.retn();
+  return this->redistributeNodeFieldT<mcIdType>(globalCellIds,fieldValueToRed);
+}
+
+DataArrayDouble *ParaUMesh::redistributeNodeField(const DataArrayIdType *globalCellIds, const DataArrayDouble *fieldValueToRed) const
+{
+  return this->redistributeNodeFieldT<double>(globalCellIds,fieldValueToRed);
 }
 
 /*!