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"
36 virtual ~CommInterface(){}
37 int worldSize() const {
39 MPI_Comm_size(MPI_COMM_WORLD, &size);
41 int commSize(MPI_Comm comm, int* size) const { return MPI_Comm_size(comm,size); }
42 int commRank(MPI_Comm comm, int* rank) const { return MPI_Comm_rank(comm,rank); }
43 int commGroup(MPI_Comm comm, MPI_Group* group) const { return MPI_Comm_group(comm, group); }
44 int groupIncl(MPI_Group group, int size, int* ranks, MPI_Group* group_output) const { return MPI_Group_incl(group, size, ranks, group_output); }
45 int commCreate(MPI_Comm comm, MPI_Group group, MPI_Comm* comm_output) const { return MPI_Comm_create(comm,group,comm_output); }
46 int groupFree(MPI_Group* group) const { return MPI_Group_free(group); }
47 int commFree(MPI_Comm* comm) const { return MPI_Comm_free(comm); }
49 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); }
50 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); }
51 int sendRecv(void* sendbuf, int sendcount, MPI_Datatype sendtype,
52 int dest, int sendtag, void* recvbuf, int recvcount,
53 MPI_Datatype recvtype, int source, int recvtag, MPI_Comm comm,
54 MPI_Status* status) { return MPI_Sendrecv(sendbuf, sendcount, sendtype, dest, sendtag, recvbuf, recvcount, recvtype, source, recvtag, comm,status); }
56 int Isend(void* buffer, int count, MPI_Datatype datatype, int target,
57 int tag, MPI_Comm comm, MPI_Request *request) const { return MPI_Isend(buffer,count, datatype, target, tag, comm, request); }
58 int Irecv(void* buffer, int count, MPI_Datatype datatype, int source,
59 int tag, MPI_Comm comm, MPI_Request* request) const { return MPI_Irecv(buffer,count, datatype, source, tag, comm, request); }
61 int wait(MPI_Request *request, MPI_Status *status) const { return MPI_Wait(request, status); }
62 int test(MPI_Request *request, int *flag, MPI_Status *status) const { return MPI_Test(request, flag, status); }
63 int requestFree(MPI_Request *request) const { return MPI_Request_free(request); }
64 int waitany(int count, MPI_Request *array_of_requests, int *index, MPI_Status *status) const { return MPI_Waitany(count, array_of_requests, index, status); }
65 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); }
66 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); }
67 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); }
68 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); }
69 int testsome(int incount, MPI_Request *array_of_requests, int *outcount,
70 int *array_of_indices, MPI_Status *array_of_status) const { return MPI_Testsome(incount, array_of_requests, outcount, array_of_indices, array_of_status); }
71 int probe(int source, int tag, MPI_Comm comm, MPI_Status *status) const { return MPI_Probe(source, tag, comm, status) ; }
72 int Iprobe(int source, int tag, MPI_Comm comm, int *flag, MPI_Status *status) const { return MPI_Iprobe(source, tag, comm, flag, status) ; }
73 int cancel(MPI_Request *request) const { return MPI_Cancel(request); }
74 int testCancelled(MPI_Status *status, int *flag) const { return MPI_Test_cancelled(status, flag); }
75 int barrier(MPI_Comm comm) const { return MPI_Barrier(comm); }
76 int errorString(int errorcode, char *string, int *resultlen) const { return MPI_Error_string(errorcode, string, resultlen); }
77 int getCount(MPI_Status *status, MPI_Datatype datatype, int *count) const { return MPI_Get_count(status, datatype, count); }
79 int broadcast(void* buffer, int count, MPI_Datatype datatype, int root, MPI_Comm comm) const { return MPI_Bcast(buffer, count, datatype, root, comm); }
80 int allGather(void* sendbuf, int sendcount, MPI_Datatype sendtype,
81 void* recvbuf, int recvcount, MPI_Datatype recvtype,
82 MPI_Comm comm) const { return MPI_Allgather(sendbuf,sendcount, sendtype, recvbuf, recvcount, recvtype, comm); }
83 int allGatherV(const void *sendbuf, int sendcount, MPI_Datatype sendtype, void *recvbuf, const int recvcounts[],
84 const int displs[], MPI_Datatype recvtype, MPI_Comm comm) { return MPI_Allgatherv(sendbuf,sendcount,sendtype,recvbuf,recvcounts,displs,recvtype,comm); }
85 int allToAll(void* sendbuf, int sendcount, MPI_Datatype sendtype,
86 void* recvbuf, int recvcount, MPI_Datatype recvtype,
87 MPI_Comm comm) const { return MPI_Alltoall(sendbuf, sendcount, sendtype, recvbuf, recvcount, recvtype, comm); }
88 int allToAllV(const void* sendbuf, int* sendcounts, int* senddispls,
89 MPI_Datatype sendtype, void* recvbuf, int* recvcounts,
90 int* recvdispls, MPI_Datatype recvtype,
91 MPI_Comm comm) const { return MPI_Alltoallv(sendbuf, sendcounts, senddispls, sendtype, recvbuf, recvcounts, recvdispls, recvtype, comm); }
93 int reduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype,
94 MPI_Op op, int root, MPI_Comm comm) const { return MPI_Reduce(sendbuf, recvbuf, count, datatype, op, root, comm); }
95 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); }
99 * \a counts is expected to be an array of array length. This method returns an array of split array.
101 static std::unique_ptr<mcIdType[]> SplitArrayOfLength(const std::unique_ptr<mcIdType[]>& counts, std::size_t countsSz, int rk, int size)
103 std::unique_ptr<mcIdType[]> ret(new mcIdType[countsSz]);
104 for(std::size_t i=0;i<countsSz;++i)
107 DataArray::GetSlice(0,counts[i],1,rk,size,a,b);
114 * Helper of alltoallv and allgatherv
117 static std::unique_ptr<int []> ToIntArray(const std::unique_ptr<T []>& arr, std::size_t size)
119 std::unique_ptr<int []> ret(new int[size]);
120 std::copy(arr.get(),arr.get()+size,ret.get());
125 * Helper of alltoallv and allgatherv
127 static std::unique_ptr<int []> ComputeOffset(const std::unique_ptr<int []>& counts, std::size_t sizeOfCounts)
129 std::unique_ptr<int []> ret(new int[sizeOfCounts]);
131 for(std::size_t i = 1 ; i < sizeOfCounts ; ++i)
133 ret[i] = ret[i-1] + counts[i-1];