1 // Copyright (C) 2007-2020 CEA/DEN, EDF R&D
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;
60 virtual ~CommInterface() { }
61 int worldSize() const {
63 MPI_Comm_size(MPI_COMM_WORLD, &size);
65 int commSize(MPI_Comm comm, int* size) const { return MPI_Comm_size(comm,size); }
66 int commRank(MPI_Comm comm, int* rank) const { return MPI_Comm_rank(comm,rank); }
67 int commGroup(MPI_Comm comm, MPI_Group* group) const { return MPI_Comm_group(comm, group); }
68 int groupIncl(MPI_Group group, int size, int* ranks, MPI_Group* group_output) const { return MPI_Group_incl(group, size, ranks, group_output); }
69 int commCreate(MPI_Comm comm, MPI_Group group, MPI_Comm* comm_output) const { return MPI_Comm_create(comm,group,comm_output); }
70 int groupFree(MPI_Group* group) const { return MPI_Group_free(group); }
71 int commFree(MPI_Comm* comm) const { return MPI_Comm_free(comm); }
73 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); }
74 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); }
75 int sendRecv(void* sendbuf, int sendcount, MPI_Datatype sendtype,
76 int dest, int sendtag, void* recvbuf, int recvcount,
77 MPI_Datatype recvtype, int source, int recvtag, MPI_Comm comm,
78 MPI_Status* status) { return MPI_Sendrecv(sendbuf, sendcount, sendtype, dest, sendtag, recvbuf, recvcount, recvtype, source, recvtag, comm,status); }
80 int Isend(void* buffer, int count, MPI_Datatype datatype, int target,
81 int tag, MPI_Comm comm, MPI_Request *request) const { return MPI_Isend(buffer,count, datatype, target, tag, comm, request); }
82 int Irecv(void* buffer, int count, MPI_Datatype datatype, int source,
83 int tag, MPI_Comm comm, MPI_Request* request) const { return MPI_Irecv(buffer,count, datatype, source, tag, comm, request); }
85 int wait(MPI_Request *request, MPI_Status *status) const { return MPI_Wait(request, status); }
86 int test(MPI_Request *request, int *flag, MPI_Status *status) const { return MPI_Test(request, flag, status); }
87 int requestFree(MPI_Request *request) const { return MPI_Request_free(request); }
88 int waitany(int count, MPI_Request *array_of_requests, int *index, MPI_Status *status) const { return MPI_Waitany(count, array_of_requests, index, status); }
89 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); }
90 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); }
91 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); }
92 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); }
93 int testsome(int incount, MPI_Request *array_of_requests, int *outcount,
94 int *array_of_indices, MPI_Status *array_of_status) const { return MPI_Testsome(incount, array_of_requests, outcount, array_of_indices, array_of_status); }
95 int probe(int source, int tag, MPI_Comm comm, MPI_Status *status) const { return MPI_Probe(source, tag, comm, status) ; }
96 int Iprobe(int source, int tag, MPI_Comm comm, int *flag, MPI_Status *status) const { return MPI_Iprobe(source, tag, comm, flag, status) ; }
97 int cancel(MPI_Request *request) const { return MPI_Cancel(request); }
98 int testCancelled(MPI_Status *status, int *flag) const { return MPI_Test_cancelled(status, flag); }
99 int barrier(MPI_Comm comm) const { return MPI_Barrier(comm); }
100 int errorString(int errorcode, char *string, int *resultlen) const { return MPI_Error_string(errorcode, string, resultlen); }
101 int getCount(MPI_Status *status, MPI_Datatype datatype, int *count) const { return MPI_Get_count(status, datatype, count); }
103 int broadcast(void* buffer, int count, MPI_Datatype datatype, int root, MPI_Comm comm) const { return MPI_Bcast(buffer, count, datatype, root, comm); }
104 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(sendbuf,sendcount,sendtype,recvbuf,recvcount,recvtype,root,comm); }
105 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(sendbuf,sendcount,sendtype,recvbuf,recvcounts,displs,recvtype,root,comm); }
106 int allGather(void* sendbuf, int sendcount, MPI_Datatype sendtype,
107 void* recvbuf, int recvcount, MPI_Datatype recvtype,
108 MPI_Comm comm) const { return MPI_Allgather(sendbuf,sendcount, sendtype, recvbuf, recvcount, recvtype, comm); }
109 int allGatherV(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int recvcounts[],
110 const int displs[], MPI_Datatype recvtype, MPI_Comm comm) const { return MPI_Allgatherv(sendbuf,sendcount,sendtype,recvbuf,recvcounts,displs,recvtype,comm); }
111 int allToAll(void* sendbuf, int sendcount, MPI_Datatype sendtype,
112 void* recvbuf, int recvcount, MPI_Datatype recvtype,
113 MPI_Comm comm) const { return MPI_Alltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm); }
114 int allToAllV(const void* sendbuf, int* sendcounts, int* senddispls,
115 MPI_Datatype sendtype, void* recvbuf, int* recvcounts,
116 int* recvdispls, MPI_Datatype recvtype,
117 MPI_Comm comm) const { return MPI_Alltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm); }
119 int reduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype,
120 MPI_Op op, int root, MPI_Comm comm) const { return MPI_Reduce(sendbuf, recvbuf, count, datatype, op, root, comm); }
121 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); }
123 void gatherArrays(MPI_Comm comm, int root, const DataArrayIdType *array, std::vector< MCAuto<DataArrayIdType> >& arraysOut) const;
124 void allGatherArrays(MPI_Comm comm, const DataArrayIdType *array, std::vector< MCAuto<DataArrayIdType> >& arraysOut) const;
125 int allGatherArrays(MPI_Comm comm, const DataArrayIdType *array, std::unique_ptr<mcIdType[]>& result, std::unique_ptr<mcIdType[]>& resultIndex) const;
126 void allToAllArrays(MPI_Comm comm, const std::vector< MCAuto<DataArrayIdType> >& arrays, std::vector< MCAuto<DataArrayIdType> >& arraysOut) const;
127 void allToAllArrays(MPI_Comm comm, const std::vector< MCAuto<DataArrayDouble> >& arrays, std::vector< MCAuto<DataArrayDouble> >& arraysOut) const;
128 void allToAllArrays(MPI_Comm comm, const std::vector< MCAuto<DataArrayDouble> >& arrays, MCAuto<DataArrayDouble>& arraysOut) const;
131 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
133 using DataArrayT = typename Traits<T>::ArrayType;
135 this->commSize(comm,&size);
137 this->commRank(comm,&rank);
138 std::unique_ptr<mcIdType[]> nbOfElems;
140 nbOfElems.reset(new mcIdType[size]);
141 mcIdType nbOfCellsRequested(array->getNumberOfTuples());
142 this->gather(&nbOfCellsRequested,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,root,comm);
143 std::unique_ptr<int[]> nbOfElemsInt,offsetsIn;
146 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems.get(),nbOfElems.get()+size,0));
147 result.reset(new T[nbOfCellIdsSum]);
148 nbOfElemsInt = CommInterface::ToIntArray<mcIdType>(nbOfElems,size);
149 offsetsIn = CommInterface::ComputeOffset(nbOfElemsInt,size);
151 this->gatherV(array->begin(),nbOfCellsRequested,ParaTraits<T>::MPIDataType,result.get(),nbOfElemsInt.get(),offsetsIn.get(),ParaTraits<T>::MPIDataType,root,comm);
154 resultIndex = ComputeOffsetFull<mcIdType>(nbOfElems,size);
160 void gatherArraysT2(MPI_Comm comm, int root, const typename Traits<T>::ArrayType *array, std::vector< MCAuto<typename Traits<T>::ArrayType> >& arraysOut) const
162 using DataArrayT = typename Traits<T>::ArrayType;
163 std::unique_ptr<T[]> result;
164 std::unique_ptr<mcIdType[]> resultIndex;
166 int size(this->gatherArraysT<T>(comm,root,array,result,resultIndex,rank));
167 arraysOut.resize(size);
168 for(int i = 0 ; i < size ; ++i)
170 arraysOut[i] = DataArrayT::New();
173 mcIdType nbOfEltPack(resultIndex[i+1]-resultIndex[i]);
174 arraysOut[i]->alloc(nbOfEltPack,1);
175 std::copy(result.get()+resultIndex[i],result.get()+resultIndex[i+1],arraysOut[i]->getPointer());
181 int allGatherArraysT(MPI_Comm comm, const typename Traits<T>::ArrayType *array, std::unique_ptr<T[]>& result, std::unique_ptr<mcIdType[]>& resultIndex) const
183 using DataArrayT = typename Traits<T>::ArrayType;
185 this->commSize(comm,&size);
186 std::unique_ptr<mcIdType[]> nbOfElems(new mcIdType[size]);
187 mcIdType nbOfCellsRequested(array->getNumberOfTuples());
188 this->allGather(&nbOfCellsRequested,1,MPI_ID_TYPE,nbOfElems.get(),1,MPI_ID_TYPE,comm);
189 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems.get(),nbOfElems.get()+size,0));
190 result.reset(new T[nbOfCellIdsSum]);
191 std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElems,size) );
192 std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) );
193 this->allGatherV(array->begin(),nbOfCellsRequested,ParaTraits<T>::MPIDataType,result.get(),nbOfElemsInt.get(),offsetsIn.get(),ParaTraits<T>::MPIDataType,comm);
194 resultIndex = ComputeOffsetFull<mcIdType>(nbOfElems,size);
199 void allGatherArraysT2(MPI_Comm comm, const typename Traits<T>::ArrayType *array, std::vector< MCAuto<typename Traits<T>::ArrayType> >& arraysOut) const
201 using DataArrayT = typename Traits<T>::ArrayType;
202 std::unique_ptr<T[]> result;
203 std::unique_ptr<mcIdType[]> resultIndex;
204 int size(this->allGatherArraysT<T>(comm,array,result,resultIndex));
205 arraysOut.resize(size);
206 for(int i = 0 ; i < size ; ++i)
208 arraysOut[i] = DataArrayT::New();
209 mcIdType nbOfEltPack(resultIndex[i+1]-resultIndex[i]);
210 arraysOut[i]->alloc(nbOfEltPack,1);
211 std::copy(result.get()+resultIndex[i],result.get()+resultIndex[i+1],arraysOut[i]->getPointer());
216 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
218 using DataArrayT = typename Traits<T>::ArrayType;
220 this->commSize(comm,&size);
221 if( arrays.size() != ToSizeT(size) )
222 throw INTERP_KERNEL::Exception("AllToAllArrays : internal error ! Invalid size of input array.");
224 std::vector< const DataArrayT *> arraysBis(FromVecAutoToVecOfConst<DataArrayT>(arrays));
225 std::unique_ptr<mcIdType[]> nbOfElems3(new mcIdType[size]);
226 nbOfElems2.reset(new mcIdType[size]);
227 nbOfComponents = std::numeric_limits<mcIdType>::max();
228 for(int curRk = 0 ; curRk < size ; ++curRk)
230 mcIdType curNbOfCompo( ToIdType( arrays[curRk]->getNumberOfComponents() ) );
231 if(nbOfComponents != std::numeric_limits<mcIdType>::max())
233 if( nbOfComponents != curNbOfCompo )
234 throw INTERP_KERNEL::Exception("AllToAllArrays : internal error ! Nb of components is not homogeneous !");
238 nbOfComponents = curNbOfCompo;
240 nbOfElems3[curRk] = arrays[curRk]->getNbOfElems();
242 this->allToAll(nbOfElems3.get(),1,MPI_ID_TYPE,nbOfElems2.get(),1,MPI_ID_TYPE,comm);
243 mcIdType nbOfCellIdsSum(std::accumulate(nbOfElems2.get(),nbOfElems2.get()+size,0));
244 arrayOut = DataArrayT::New();
245 arrayOut->alloc(nbOfCellIdsSum/nbOfComponents,nbOfComponents);
246 std::unique_ptr<int[]> nbOfElemsInt( CommInterface::ToIntArray<mcIdType>(nbOfElems3,size) ),nbOfElemsOutInt( CommInterface::ToIntArray<mcIdType>(nbOfElems2,size) );
247 std::unique_ptr<int[]> offsetsIn( CommInterface::ComputeOffset(nbOfElemsInt,size) ), offsetsOut( CommInterface::ComputeOffset(nbOfElemsOutInt,size) );
249 MCAuto<DataArrayT> arraysAcc(DataArrayT::Aggregate(arraysBis));
250 this->allToAllV(arraysAcc->begin(),nbOfElemsInt.get(),offsetsIn.get(),ParaTraits<T>::MPIDataType,
251 arrayOut->getPointer(),nbOfElemsOutInt.get(),offsetsOut.get(),ParaTraits<T>::MPIDataType,comm);
257 void allToAllArraysT(MPI_Comm comm, const std::vector< MCAuto<typename Traits<T>::ArrayType> >& arrays, std::vector< MCAuto<typename Traits<T>::ArrayType> >& arraysOut) const
259 using DataArrayT = typename Traits<T>::ArrayType;
260 MCAuto<DataArrayT> cellIdsFromProcs;
261 std::unique_ptr<mcIdType[]> nbOfElems2;
262 mcIdType nbOfComponents(0);
263 int size(this->allToAllArraysT2<T>(comm,arrays,cellIdsFromProcs,nbOfElems2,nbOfComponents));
264 std::unique_ptr<mcIdType[]> offsetsOutIdType( CommInterface::ComputeOffset(nbOfElems2,size) );
265 // build output arraysOut by spliting cellIdsFromProcs into parts
266 arraysOut.resize(size);
267 for(int curRk = 0 ; curRk < size ; ++curRk)
269 arraysOut[curRk] = DataArrayT::NewFromArray(cellIdsFromProcs->begin()+offsetsOutIdType[curRk],cellIdsFromProcs->begin()+offsetsOutIdType[curRk]+nbOfElems2[curRk]);
270 arraysOut[curRk]->rearrange(nbOfComponents);
276 * \a counts is expected to be an array of array length. This method returns an array of split array.
278 static std::unique_ptr<mcIdType[]> SplitArrayOfLength(const std::unique_ptr<mcIdType[]>& counts, std::size_t countsSz, int rk, int size)
280 std::unique_ptr<mcIdType[]> ret(new mcIdType[countsSz]);
281 for(std::size_t i=0;i<countsSz;++i)
284 DataArray::GetSlice(0,counts[i],1,rk,size,a,b);
291 * Helper of alltoallv and allgatherv
294 static std::unique_ptr<int []> ToIntArray(const std::unique_ptr<T []>& arr, std::size_t size)
296 std::unique_ptr<int []> ret(new int[size]);
297 std::copy(arr.get(),arr.get()+size,ret.get());
302 * Helper of alltoallv and allgatherv
305 static std::unique_ptr<T []> ComputeOffset(const std::unique_ptr<T []>& counts, std::size_t sizeOfCounts)
307 std::unique_ptr<T []> ret(new T[sizeOfCounts]);
308 ret[0] = static_cast<T>(0);
309 for(std::size_t i = 1 ; i < sizeOfCounts ; ++i)
311 ret[i] = ret[i-1] + counts[i-1];
317 * Helper of alltoallv and allgatherv
320 static std::unique_ptr<T []> ComputeOffsetFull(const std::unique_ptr<T []>& counts, std::size_t sizeOfCounts)
322 std::unique_ptr<T []> ret(new T[sizeOfCounts+1]);
323 ret[0] = static_cast<T>(0);
324 for(std::size_t i = 1 ; i < sizeOfCounts+1 ; ++i)
326 ret[i] = ret[i-1] + counts[i-1];