1 // Copyright (C) 2007-2024 CEA, EDF
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
22 #include "ParaIdType.hxx"
23 #include "MEDCouplingMemArray.hxx"
39 struct ParaTraits<double>
41 static MPI_Datatype MPIDataType;
45 struct ParaTraits<Int32>
47 static MPI_Datatype MPIDataType;
51 struct ParaTraits<Int64>
53 static MPI_Datatype MPIDataType;
56 /*! \anchor CommInterface-det
59 The class \a CommInterface is the gateway to the MPI library.
60 It is a wrapper around all MPI calls, thus trying to abstract the rest of the code from using the direct MPI API
61 (but this is not strictly respected overall in practice ...). It is used in all
62 the \ref parallel "DEC related classes".
64 It is typically instantiated after the MPI_Init() call in a program and is afterwards passed as a
65 parameter to the constructors of various \ref parallel "parallel objects" so that they access the
66 MPI library via this common interface.
68 As an example, the following code excerpt initializes a processor group made of the zero processor.
71 #include "CommInterface.hxx"
72 #include "ProcessorGroup.hxx"
74 int main(int argc, char** argv)
77 MPI_Init(&argc, &argv);
78 MEDCoupling::CommInterface comm_interface;
80 //setting up a processor group with proc 0
83 MEDCoupling::ProcessorGroup group(procs, comm_interface);
94 virtual ~CommInterface() { }
95 int worldSize() const {
97 MPI_Comm_size(MPI_COMM_WORLD, &size);
99 int commSize(MPI_Comm comm, int* size) const { return MPI_Comm_size(comm,size); }
100 int commRank(MPI_Comm comm, int* rank) const { return MPI_Comm_rank(comm,rank); }
101 int commGroup(MPI_Comm comm, MPI_Group* group) const { return MPI_Comm_group(comm, group); }
102 int groupIncl(MPI_Group group, int size, int* ranks, MPI_Group* group_output) const { return MPI_Group_incl(group, size, ranks, group_output); }
103 int commCreate(MPI_Comm comm, MPI_Group group, MPI_Comm* comm_output) const { return MPI_Comm_create(comm,group,comm_output); }
104 int groupFree(MPI_Group* group) const { return MPI_Group_free(group); }
105 int commFree(MPI_Comm* comm) const { return MPI_Comm_free(comm); }
107 int send(void* buffer, int count, MPI_Datatype datatype, int target, int tag, MPI_Comm comm) const { return MPI_Send(buffer,count, datatype, target, tag, comm); }
108 int recv(void* buffer, int count, MPI_Datatype datatype, int source, int tag, MPI_Comm comm, MPI_Status* status) const { return MPI_Recv(buffer,count, datatype, source, tag, comm, status); }
109 int sendRecv(void* sendbuf, int sendcount, MPI_Datatype sendtype,
110 int dest, int sendtag, void* recvbuf, int recvcount,
111 MPI_Datatype recvtype, int source, int recvtag, MPI_Comm comm,
112 MPI_Status* status) { return MPI_Sendrecv(sendbuf, sendcount, sendtype, dest, sendtag, recvbuf, recvcount, recvtype, source, recvtag, comm,status); }
114 int Isend(void* buffer, int count, MPI_Datatype datatype, int target,
115 int tag, MPI_Comm comm, MPI_Request *request) const { return MPI_Isend(buffer,count, datatype, target, tag, comm, request); }
116 int Irecv(void* buffer, int count, MPI_Datatype datatype, int source,
117 int tag, MPI_Comm comm, MPI_Request* request) const { return MPI_Irecv(buffer,count, datatype, source, tag, comm, request); }
119 int wait(MPI_Request *request, MPI_Status *status) const { return MPI_Wait(request, status); }
120 int test(MPI_Request *request, int *flag, MPI_Status *status) const { return MPI_Test(request, flag, status); }
121 int requestFree(MPI_Request *request) const { return MPI_Request_free(request); }
122 int waitany(int count, MPI_Request *array_of_requests, int *index, MPI_Status *status) const { return MPI_Waitany(count, array_of_requests, index, status); }
123 int testany(int count, MPI_Request *array_of_requests, int *index, int *flag, MPI_Status *status) const { return MPI_Testany(count, array_of_requests, index, flag, status); }
124 int waitall(int count, MPI_Request *array_of_requests, MPI_Status *array_of_status) const { return MPI_Waitall(count, array_of_requests, array_of_status); }
125 int testall(int count, MPI_Request *array_of_requests, int *flag, MPI_Status *array_of_status) const { return MPI_Testall(count, array_of_requests, flag, array_of_status); }
126 int waitsome(int incount, MPI_Request *array_of_requests,int *outcount, int *array_of_indices, MPI_Status *array_of_status) const { return MPI_Waitsome(incount, array_of_requests, outcount, array_of_indices, array_of_status); }
127 int testsome(int incount, MPI_Request *array_of_requests, int *outcount,
128 int *array_of_indices, MPI_Status *array_of_status) const { return MPI_Testsome(incount, array_of_requests, outcount, array_of_indices, array_of_status); }
129 int probe(int source, int tag, MPI_Comm comm, MPI_Status *status) const { return MPI_Probe(source, tag, comm, status) ; }
130 int Iprobe(int source, int tag, MPI_Comm comm, int *flag, MPI_Status *status) const { return MPI_Iprobe(source, tag, comm, flag, status) ; }
131 int cancel(MPI_Request *request) const { return MPI_Cancel(request); }
132 int testCancelled(MPI_Status *status, int *flag) const { return MPI_Test_cancelled(status, flag); }
133 int barrier(MPI_Comm comm) const { return MPI_Barrier(comm); }
134 int errorString(int errorcode, char *string, int *resultlen) const { return MPI_Error_string(errorcode, string, resultlen); }
135 int getCount(MPI_Status *status, MPI_Datatype datatype, int *count) const { return MPI_Get_count(status, datatype, count); }
137 int broadcast(void* buffer, int count, MPI_Datatype datatype, int root, MPI_Comm comm) const { return MPI_Bcast(buffer, count, datatype, root, comm); }
138 int gather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm) const { return MPI_Gather(const_cast<void*>(sendbuf),sendcount,sendtype,recvbuf,recvcount,recvtype,root,comm); }
139 int gatherV(const void* sendbuf, int sendcount, MPI_Datatype sendtype, void* recvbuf, const int recvcounts[], const int displs[], MPI_Datatype recvtype, int root, MPI_Comm comm) const { return MPI_Gatherv(const_cast<void*>(sendbuf),sendcount,sendtype,recvbuf,const_cast<int *>(recvcounts),const_cast<int *>(displs),recvtype,root,comm); }
140 int allGather(void* sendbuf, int sendcount, MPI_Datatype sendtype,
141 void* recvbuf, int recvcount, MPI_Datatype recvtype,
142 MPI_Comm comm) const { return MPI_Allgather(sendbuf,sendcount, sendtype, recvbuf, recvcount, recvtype, comm); }
143 int allGatherV(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int recvcounts[],
144 const int displs[], MPI_Datatype recvtype, MPI_Comm comm) const { return MPI_Allgatherv(const_cast<void*>(sendbuf),sendcount,sendtype,recvbuf,const_cast<int *>(recvcounts),const_cast<int *>(displs),recvtype,comm); }
145 int allToAll(void* sendbuf, int sendcount, MPI_Datatype sendtype,
146 void* recvbuf, int recvcount, MPI_Datatype recvtype,
147 MPI_Comm comm) const { return MPI_Alltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm); }
148 int allToAllV(const void* sendbuf, int* sendcounts, int* senddispls,
149 MPI_Datatype sendtype, void* recvbuf, int* recvcounts,
150 int* recvdispls, MPI_Datatype recvtype,
151 MPI_Comm comm) const { return MPI_Alltoallv(const_cast<void*>(sendbuf), sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm); }
153 int reduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype,
154 MPI_Op op, int root, MPI_Comm comm) const { return MPI_Reduce(sendbuf, recvbuf, count, datatype, op, root, comm); }
155 int allReduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) const { return MPI_Allreduce(sendbuf, recvbuf, count, datatype, op, comm); }
157 void gatherArrays(MPI_Comm comm, int root, const DataArrayIdType *array, std::vector< MCAuto<DataArrayIdType> >& arraysOut) const;
158 void allGatherArrays(MPI_Comm comm, const DataArrayIdType *array, std::vector< MCAuto<DataArrayIdType> >& arraysOut) const;
159 int allGatherArrays(MPI_Comm comm, const DataArrayIdType *array, std::unique_ptr<mcIdType[]>& result, std::unique_ptr<mcIdType[]>& resultIndex) const;
160 void allToAllArrays(MPI_Comm comm, const std::vector< MCAuto<DataArrayIdType> >& arrays, std::vector< MCAuto<DataArrayIdType> >& arraysOut) const;
161 void allToAllArrays(MPI_Comm comm, const std::vector< MCAuto<DataArrayDouble> >& arrays, std::vector< MCAuto<DataArrayDouble> >& arraysOut) const;
162 void allToAllArrays(MPI_Comm comm, const std::vector< MCAuto<DataArrayDouble> >& arrays, MCAuto<DataArrayDouble>& arraysOut) const;
165 int gatherArraysT(MPI_Comm comm, int root, const typename Traits<T>::ArrayType *array, std::unique_ptr<T[]>& result, std::unique_ptr<mcIdType[]>& resultIndex, int& rank) const
168 this->commSize(comm,&size);
170 this->commRank(comm,&rank);
171 std::unique_ptr<mcIdType[]> nbOfElems;
173 nbOfElems.reset(new mcIdType[size]);
174 mcIdType nbOfCellsRequested(array->getNumberOfTuples());
175 this->gather(&nbOfCellsRequested,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,root,comm);
176 std::unique_ptr<int[]> nbOfElemsInt,offsetsIn;
179 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems.get(),nbOfElems.get()+size,0));
180 result.reset(new T[nbOfCellIdsSum]);
181 nbOfElemsInt = CommInterface::ToIntArray<mcIdType>(nbOfElems,size);
182 offsetsIn = CommInterface::ComputeOffset(nbOfElemsInt,size);
184 this->gatherV(array->begin(),nbOfCellsRequested,ParaTraits<T>::MPIDataType,result.get(),nbOfElemsInt.get(),offsetsIn.get(),ParaTraits<T>::MPIDataType,root,comm);
187 resultIndex = ComputeOffsetFull<mcIdType>(nbOfElems,size);
193 void gatherArraysT2(MPI_Comm comm, int root, const typename Traits<T>::ArrayType *array, std::vector< MCAuto<typename Traits<T>::ArrayType> >& arraysOut) const
195 using DataArrayT = typename Traits<T>::ArrayType;
196 std::unique_ptr<T[]> result;
197 std::unique_ptr<mcIdType[]> resultIndex;
199 int size(this->gatherArraysT<T>(comm,root,array,result,resultIndex,rank));
200 arraysOut.resize(size);
201 for(int i = 0 ; i < size ; ++i)
203 arraysOut[i] = DataArrayT::New();
206 mcIdType nbOfEltPack(resultIndex[i+1]-resultIndex[i]);
207 arraysOut[i]->alloc(nbOfEltPack,1);
208 std::copy(result.get()+resultIndex[i],result.get()+resultIndex[i+1],arraysOut[i]->getPointer());
214 int allGatherArraysT(MPI_Comm comm, const typename Traits<T>::ArrayType *array, std::unique_ptr<T[]>& result, std::unique_ptr<mcIdType[]>& resultIndex) const
217 this->commSize(comm,&size);
218 std::unique_ptr<mcIdType[]> nbOfElems(new mcIdType[size]);
219 mcIdType nbOfCellsRequested(array->getNumberOfTuples());
220 this->allGather(&nbOfCellsRequested,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,comm);
221 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems.get(),nbOfElems.get()+size,0));
222 result.reset(new T[nbOfCellIdsSum]);
223 std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElems,size) );
224 std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) );
225 this->allGatherV(array->begin(),nbOfCellsRequested,ParaTraits<T>::MPIDataType,result.get(),nbOfElemsInt.get(),offsetsIn.get(),ParaTraits<T>::MPIDataType,comm);
226 resultIndex = ComputeOffsetFull<mcIdType>(nbOfElems,size);
231 void allGatherArraysT2(MPI_Comm comm, const typename Traits<T>::ArrayType *array, std::vector< MCAuto<typename Traits<T>::ArrayType> >& arraysOut) const
233 using DataArrayT = typename Traits<T>::ArrayType;
234 std::unique_ptr<T[]> result;
235 std::unique_ptr<mcIdType[]> resultIndex;
236 int size(this->allGatherArraysT<T>(comm,array,result,resultIndex));
237 arraysOut.resize(size);
238 for(int i = 0 ; i < size ; ++i)
240 arraysOut[i] = DataArrayT::New();
241 mcIdType nbOfEltPack(resultIndex[i+1]-resultIndex[i]);
242 arraysOut[i]->alloc(nbOfEltPack,1);
243 std::copy(result.get()+resultIndex[i],result.get()+resultIndex[i+1],arraysOut[i]->getPointer());
248 int allToAllArraysT2(MPI_Comm comm, const std::vector< MCAuto<typename Traits<T>::ArrayType> >& arrays, MCAuto<typename Traits<T>::ArrayType>& arrayOut, std::unique_ptr<mcIdType[]>& nbOfElems2, mcIdType& nbOfComponents) const
250 using DataArrayT = typename Traits<T>::ArrayType;
252 this->commSize(comm,&size);
253 if( arrays.size() != ToSizeT(size) )
254 throw INTERP_KERNEL::Exception("AllToAllArrays : internal error ! Invalid size of input array.");
256 std::vector< const DataArrayT *> arraysBis(FromVecAutoToVecOfConst<DataArrayT>(arrays));
257 std::unique_ptr<mcIdType[]> nbOfElems3(new mcIdType[size]);
258 nbOfElems2.reset(new mcIdType[size]);
259 nbOfComponents = std::numeric_limits<mcIdType>::max();
260 for(int curRk = 0 ; curRk < size ; ++curRk)
262 mcIdType curNbOfCompo( ToIdType( arrays[curRk]->getNumberOfComponents() ) );
263 if(nbOfComponents != std::numeric_limits<mcIdType>::max())
265 if( nbOfComponents != curNbOfCompo )
266 throw INTERP_KERNEL::Exception("AllToAllArrays : internal error ! Nb of components is not homogeneous !");
270 nbOfComponents = curNbOfCompo;
272 nbOfElems3[curRk] = arrays[curRk]->getNbOfElems();
274 this->allToAll(nbOfElems3.get(),1,MPI_ID_TYPE,nbOfElems2.get(),1,MPI_ID_TYPE,comm);
275 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems2.get(),nbOfElems2.get()+size,0));
276 arrayOut = DataArrayT::New();
277 arrayOut->alloc(nbOfCellIdsSum/nbOfComponents,nbOfComponents);
278 std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElems3,size) ),nbOfElemsOutInt( CommInterface::ToIntArray<mcIdType>(nbOfElems2,size) );
279 std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) ), offsetsOut( CommInterface::ComputeOffset(nbOfElemsOutInt,size) );
281 MCAuto<DataArrayT> arraysAcc(DataArrayT::Aggregate(arraysBis));
282 this->allToAllV(arraysAcc->begin(),nbOfElemsInt.get(),offsetsIn.get(),ParaTraits<T>::MPIDataType,
283 arrayOut->getPointer(),nbOfElemsOutInt.get(),offsetsOut.get(),ParaTraits<T>::MPIDataType,comm);
289 void allToAllArraysT(MPI_Comm comm, const std::vector< MCAuto<typename Traits<T>::ArrayType> >& arrays, std::vector< MCAuto<typename Traits<T>::ArrayType> >& arraysOut) const
291 using DataArrayT = typename Traits<T>::ArrayType;
292 MCAuto<DataArrayT> cellIdsFromProcs;
293 std::unique_ptr<mcIdType[]> nbOfElems2;
294 mcIdType nbOfComponents(0);
295 int size(this->allToAllArraysT2<T>(comm,arrays,cellIdsFromProcs,nbOfElems2,nbOfComponents));
296 std::unique_ptr<mcIdType[]> offsetsOutIdType( CommInterface::ComputeOffset(nbOfElems2,size) );
297 // build output arraysOut by spliting cellIdsFromProcs into parts
298 arraysOut.resize(size);
299 for(int curRk = 0 ; curRk < size ; ++curRk)
301 arraysOut[curRk] = DataArrayT::NewFromArray(cellIdsFromProcs->begin()+offsetsOutIdType[curRk],cellIdsFromProcs->begin()+offsetsOutIdType[curRk]+nbOfElems2[curRk]);
302 arraysOut[curRk]->rearrange(nbOfComponents);
308 * \a counts is expected to be an array of array length. This method returns an array of split array.
310 static std::unique_ptr<mcIdType[]> SplitArrayOfLength(const std::unique_ptr<mcIdType[]>& counts, std::size_t countsSz, int rk, int size)
312 std::unique_ptr<mcIdType[]> ret(new mcIdType[countsSz]);
313 for(std::size_t i=0;i<countsSz;++i)
316 DataArray::GetSlice(0,counts[i],1,rk,size,a,b);
323 * Helper of alltoallv and allgatherv
326 static std::unique_ptr<int []> ToIntArray(const std::unique_ptr<T []>& arr, std::size_t size)
328 std::unique_ptr<int []> ret(new int[size]);
329 std::copy(arr.get(),arr.get()+size,ret.get());
334 * Helper of alltoallv and allgatherv
337 static std::unique_ptr<T []> ComputeOffset(const std::unique_ptr<T []>& counts, std::size_t sizeOfCounts)
339 std::unique_ptr<T []> ret(new T[sizeOfCounts]);
340 ret[0] = static_cast<T>(0);
341 for(std::size_t i = 1 ; i < sizeOfCounts ; ++i)
343 ret[i] = ret[i-1] + counts[i-1];
349 * Helper of alltoallv and allgatherv
352 static std::unique_ptr<T []> ComputeOffsetFull(const std::unique_ptr<T []>& counts, std::size_t sizeOfCounts)
354 std::unique_ptr<T []> ret(new T[sizeOfCounts+1]);
355 ret[0] = static_cast<T>(0);
356 for(std::size_t i = 1 ; i < sizeOfCounts+1 ; ++i)
358 ret[i] = ret[i-1] + counts[i-1];