source: XIOS/dev/dev_trunk_omp/extern/remap/src/mpi_routing.cpp @ 1602

Last change on this file since 1602 was 1602, checked in by yushan, 5 years ago

branch_openmp merged with trunk r1597

File size: 18.2 KB
Line 
1#include "mpi_routing.hpp"
2#include "mpi.hpp"
3#include "node.hpp"
4#include "elt.hpp"
5#include "timerRemap.hpp"
6#include <iostream>
7using namespace ep_lib;
8
9namespace sphereRemap {
10
11const int verbose = 0;
12
13CMPIRouting::CMPIRouting(MPI_Comm comm) : communicator(comm)
14{
15        MPI_Comm_rank(comm, &mpiRank);
16        MPI_Comm_size(comm, &mpiSize);
17}
18
19/* sparse alltoallv when it is known that only a subset of ranks will communicate,
20    but message lengths are *known* to receiver */
21template <typename T>
22void alltoalls_known(const vector<vector<T> >& send, vector<vector<T> >& recv, const vector<int>& ranks, MPI_Comm communicator)
23{
24        vector<MPI_Request> request(ranks.size() * 2);
25        vector<MPI_Status>  status(ranks.size() * 2);
26
27        // communicate data
28        int nbRequest = 0;
29        for (int i = 0; i < ranks.size(); i++)
30                if (recv[i].size())
31                        MPI_Irecv(&recv[i][0], recv[i].size()*sizeof(T), MPI_CHAR, ranks[i], 0, communicator, &request[nbRequest++]);
32        for (int i = 0; i < ranks.size(); i++)
33                if (send[i].size())
34                        MPI_Isend((void *) &send[i][0], send[i].size()*sizeof(T), MPI_CHAR, ranks[i], 0, communicator, &request[nbRequest++]);
35        MPI_Waitall(nbRequest, &request[0], &status[0]);
36}
37
38/* sparse alltoallv when it is known that only a subset of ranks will communicate,
39    but message lengths are *unknown* to receiver */
40template <typename T>
41void alltoalls_unknown(const vector<vector<T> >& send, vector<vector<T> >& recv, const vector<int>& ranks, MPI_Comm communicator)
42{
43        vector<MPI_Request> request(ranks.size() * 2);
44        vector<MPI_Status>  status(ranks.size() * 2);
45
46        // communicate sizes
47        int nbRequest = 0;
48        vector<int> sendSizes(ranks.size());
49        vector<int> recvSizes(ranks.size());
50        for (int i = 0; i < ranks.size(); i++)
51                sendSizes[i] = send[i].size();
52        for (int i = 0; i < ranks.size(); i++)
53                MPI_Irecv(&recvSizes[i], 1, MPI_INT, ranks[i], 0, communicator, &request[nbRequest++]);
54        for (int i = 0; i < ranks.size(); i++)
55                MPI_Isend(&sendSizes[i], 1, MPI_INT, ranks[i], 0, communicator, &request[nbRequest++]);
56        MPI_Waitall(nbRequest, &request[0], &status[0]);
57
58        // allocate
59        for (int i = 0; i < ranks.size(); i++)
60                if (recvSizes[i]) recv[i].resize(recvSizes[i]);
61        // communicate data
62        alltoalls_known(send, recv, ranks, communicator);
63}
64
65void setElementsSendCnt(int route, vector<int>& sendCnt)
66{
67        sendCnt[route]++;
68}
69
70void setElementsSendCnt(const vector<int>& route, vector<int>& sendCnt)
71{
72        for (int i = 0; i < route.size(); i++)
73                sendCnt[route[i]]++;
74}
75
76void setTargetElementIndex(int route, vector<int>& indices, const int *rankToIndex)
77{
78        int index = rankToIndex[route];
79        indices.push_back(index);
80}
81
82void setTargetElementIndex(const vector<int>& route, vector<int>& indices, const int *rankToIndex)
83{
84        for (int i = 0; i < route.size(); i++)
85        {
86                int index = rankToIndex[route[i]];
87                indices.push_back(index);
88        }
89}
90
91template<typename T>
92void CMPIRouting::init(const vector<T>& route, CMPICascade *cascade)
93{
94        vector<int> nbElementSend(mpiSize);
95        int *toSend = new int[mpiSize];
96        int *recvCount = new int[mpiSize];
97        int *targetRankToIndex;
98
99        for (int i = 0; i < route.size(); i++)
100                setElementsSendCnt(route[i], nbElementSend);
101
102        nbTarget = 0;
103vector<int> destRanks;
104        for (int i = 0; i < mpiSize; i++)
105        {
106                if (nbElementSend[i] == 0)
107                        toSend[i] = 0;
108                else
109                {
110destRanks.push_back(i);
111                        toSend[i] = 1;
112                        nbTarget++;
113                }
114                recvCount[i] = 1;
115        }
116//int recvCntDeb = (stree) ? stree->scatter_reduce(destRanks) : -1;
117
118        MPI_Barrier(communicator);
119        CTimer::get("CMPIRouting::init(reduce_scatter)").reset();
120        CTimer::get("CMPIRouting::init(reduce_scatter)").resume();
121        MPI_Reduce_scatter(toSend, &nbSource, recvCount, MPI_INT, MPI_SUM, communicator);
122        CTimer::get("CMPIRouting::init(reduce_scatter)").suspend();
123        CTimer::get("CMPIRouting::init(reduce_scatter)").print();
124
125        MPI_Alloc_mem(nbTarget *sizeof(int), MPI_INFO_NULL, &targetRank);
126        MPI_Alloc_mem(nbSource *sizeof(int), MPI_INFO_NULL, &sourceRank);
127
128        targetRankToIndex = new int[mpiSize];
129        int index = 0;
130        for (int i = 0; i < mpiSize; i++)
131        {
132                if (toSend[i] == 1)
133                {
134                        targetRank[index] = i;
135                        targetRankToIndex[i] = index;
136                        index++;
137                }
138        }
139
140        MPI_Barrier(communicator);
141        CTimer::get("CMPIRouting::init(get_source)").reset();
142        CTimer::get("CMPIRouting::init(get_source)").resume();
143
144        MPI_Request *request = new MPI_Request[nbSource + nbTarget];
145        MPI_Status  *status = new MPI_Status[nbSource + nbTarget];
146
147        int indexRequest = 0;
148        if (verbose) cout << "CMPIRouting::init     nbSource : " << nbSource << " nbTarget : " << nbTarget << endl;
149//      assert(recvCntDeb == -1 || recvCntDeb == nbSource);
150//cout << mpiRank <<  "DEB : " << recvCntDeb << "    nbSource " << nbSource << " nbTarget : " << nbTarget << endl;
151        for (int i = 0; i < nbSource; i++)
152        {
153                MPI_Irecv(&sourceRank[i], 1, MPI_INT, -2, 0, communicator, &request[indexRequest++]);
154        }
155        MPI_Barrier(communicator);
156        for (int i = 0; i < nbTarget; i++)
157        {
158                MPI_Isend(&mpiRank, 1, MPI_INT, targetRank[i], 0, communicator, &request[indexRequest++]);
159        }
160        MPI_Waitall(indexRequest, request, status);
161        MPI_Barrier(communicator);  //needed
162        CTimer::get("CMPIRouting::init(get_source)").suspend();
163        CTimer::get("CMPIRouting::init(get_source)").print();
164
165        CTimer::get("CMPIRouting::init(get_source)").reset();
166        CTimer::get("CMPIRouting::init(get_source)").resume();
167
168        indexRequest = 0;
169        for (int i = 0; i < nbSource; i++)
170        {
171                MPI_Irecv(&sourceRank[i], 1, MPI_INT, -2, 0, communicator, &request[indexRequest]);
172                indexRequest++;
173        }
174
175        for (int i = 0; i < nbTarget; i++)
176        {
177                MPI_Isend(&mpiRank, 1, MPI_INT, targetRank[i], 0, communicator, &request[indexRequest]);
178                indexRequest++;
179        }
180        MPI_Waitall(indexRequest, request, status);
181        MPI_Barrier(communicator);
182        CTimer::get("CMPIRouting::init(get_source)").suspend();
183        CTimer::get("CMPIRouting::init(get_source)").print();
184
185        CTimer::get("CMPIRouting::init(send_element)").reset();
186        CTimer::get("CMPIRouting::init(send_element)").resume();
187
188        nbTargetElement.resize(nbTarget);
189        nbSourceElement.resize(nbSource);
190
191        for (int i = 0; i < route.size(); i++)
192                setTargetElementIndex(route[i], targetElementIndex, targetRankToIndex);
193
194        for (int i = 0; i < targetElementIndex.size(); i++)
195                nbTargetElement[targetElementIndex[i]]++;
196
197        indexRequest = 0;
198        totalSourceElement = 0;
199        totalTargetElement = 0;
200        for (int i = 0; i < nbSource; i++)
201        {
202                MPI_Irecv(&nbSourceElement[i], 1, MPI_INT, sourceRank[i], 0, communicator, &request[indexRequest]);
203                indexRequest++;
204        }
205
206        for (int i = 0; i < nbTarget; i++)
207        {
208                totalTargetElement += nbTargetElement[i];
209                MPI_Isend(&nbTargetElement[i], 1, MPI_INT, targetRank[i], 0, communicator, &request[indexRequest]);
210                indexRequest++;
211        }
212
213        MPI_Waitall(indexRequest, request, status);
214
215        CTimer::get("CMPIRouting::init(send_element)").suspend();
216        CTimer::get("CMPIRouting::init(send_element)").print();
217
218        totalSourceElement = 0;
219        for (int i = 0; i < nbSource; i++)
220                totalSourceElement += nbSourceElement[i];
221
222        sourceElementIndex.resize(totalSourceElement);
223
224        totalSourceElement = 0;
225        for (int i = 0; i < nbSource; i++)
226        {
227                for (int j = 0; j < nbSourceElement[i]; j++)
228                {
229                        sourceElementIndex[totalSourceElement] = i;
230                        totalSourceElement++;
231                }
232        }
233
234        delete[]  toSend;
235        delete[]  recvCount;
236        delete[]  request;
237        delete[]  status;
238}
239
240
241int CMPIRouting::getTotalSourceElement(void)
242{
243        return totalSourceElement;
244}
245
246template<typename T>
247void CMPIRouting::transferToTarget(T* targetElements, T* sourceElements)
248{
249        char** targetBuffer = new char*[nbTarget];
250        int* indexTargetBuffer = new int[nbTarget];
251
252        for(int i=0;i<nbTarget; i++)
253        {
254                targetBuffer[i] = new char[sizeof(T)*nbTargetElement[i]];
255                indexTargetBuffer[i]= 0;
256        }
257
258        char** sourceBuffer = new char*[nbSource];
259        int* indexSourceBuffer = new int[nbSource];
260
261        for(int i=0;i<nbSource; i++)
262        {
263                sourceBuffer[i] = new char[sizeof(T)*nbSourceElement[i]];
264                indexSourceBuffer[i]= 0;
265        }
266
267        // pack the data
268        int index;
269        for(int i=0;i<totalTargetElement;i++)
270        {
271                index=targetElementIndex[i];
272                *((T*) &(targetBuffer[index][indexTargetBuffer[index]])) = targetElements[i];
273                indexTargetBuffer[index]+=sizeof(T);
274        }
275
276
277        MPI_Request* request=new MPI_Request[nbSource+nbTarget];
278        MPI_Status*  status=new MPI_Status[nbSource+nbTarget];
279        int indexRequest=0;
280
281        MPI_Barrier(communicator);
282        CTimer::get("CMPIRouting::transferToTarget").reset();
283        CTimer::get("CMPIRouting::transferToTarget").resume();
284
285        for(int i=0; i<nbSource; i++)
286        {
287                MPI_Irecv(sourceBuffer[i],nbSourceElement[i]*sizeof(T),MPI_CHAR, sourceRank[i], 0, communicator, &request[indexRequest]);
288                indexRequest++;
289        }
290
291        for(int i=0;i<nbTarget; i++)
292        {
293                MPI_Isend(targetBuffer[i],nbTargetElement[i]*sizeof(T), MPI_CHAR, targetRank[i], 0, communicator, &request[indexRequest]);
294                indexRequest++;
295        }
296
297        MPI_Waitall(indexRequest,request,status);
298
299        CTimer::get("CMPIRouting::transferToTarget").suspend();
300        CTimer::get("CMPIRouting::transferToTarget").print();
301        MPI_Barrier(communicator);
302
303        // unpack the data
304        for(int i=0;i<totalSourceElement;i++)
305        {
306                index=sourceElementIndex[i];
307                sourceElements[i] = *((T*) &(sourceBuffer[index][indexSourceBuffer[index]]));
308                indexSourceBuffer[index]+=sizeof(T);
309        }
310
311        for(int i=0;i<nbTarget; i++) delete[] targetBuffer[i];
312        for(int i=0;i<nbSource; i++) delete[] sourceBuffer[i];
313        delete[] targetBuffer;
314        delete[] indexTargetBuffer;
315        delete[] sourceBuffer;
316        delete[] indexSourceBuffer;
317        delete[] request;
318        delete[] status;
319}
320
321
322template<typename T, typename t_pack, typename t_unpack>
323void CMPIRouting::transferToTarget(T* targetElements, T* sourceElements, t_pack pack, t_unpack unpack)
324{
325        char** targetBuffer = new char*[nbTarget];
326        int* indexTargetBuffer = new int[nbTarget];
327        int* targetMessageSize = new int[nbTarget];
328        int* sourceMessageSize = new int[nbSource];
329        int index;
330
331        // compute target buffer size
332        for (int i = 0; i < nbTarget; i++)
333                targetMessageSize[i] = 0;
334
335        for (int i = 0; i < totalTargetElement; i++)
336        {
337                index = targetElementIndex[i];
338                pack(targetElements[i], NULL, targetMessageSize[index]);
339        }
340
341        MPI_Request *request = new MPI_Request[nbSource + nbTarget];
342        MPI_Status  *status = new MPI_Status[nbSource + nbTarget];
343        int indexRequest = 0;
344
345        MPI_Barrier(communicator);
346        CTimer::get("CMPIRouting::transferToTarget(messageSize)").reset();
347        CTimer::get("CMPIRouting::transferToTarget(messageSize)").resume();
348
349        for(int i=0; i<nbSource; i++)
350        {
351                MPI_Irecv(&sourceMessageSize[i],1,MPI_INT, sourceRank[i], 0, communicator, &request[indexRequest]);
352                indexRequest++;
353        }
354
355        for(int i=0; i<nbTarget; i++)
356        {
357                MPI_Isend(&targetMessageSize[i],1, MPI_INT, targetRank[i], 0, communicator, &request[indexRequest]);
358                indexRequest++;
359        }
360
361        MPI_Waitall(indexRequest,request,status);
362
363        MPI_Barrier(communicator);
364        CTimer::get("CMPIRouting::transferToTarget(messageSize)").suspend();
365        CTimer::get("CMPIRouting::transferToTarget(messageSize)").print();
366
367        for(int i=0; i<nbTarget; i++)
368        {
369                targetBuffer[i] = new char[targetMessageSize[i]];
370                indexTargetBuffer[i] = 0;
371        }
372
373        char** sourceBuffer = new char*[nbSource];
374        int* indexSourceBuffer = new int[nbSource];
375
376        for(int i=0;i<nbSource; i++)
377        {
378                sourceBuffer[i] = new char[sourceMessageSize[i]];
379                indexSourceBuffer[i]= 0;
380        }
381
382        // pack the data
383        for(int i=0; i<totalTargetElement; i++)
384        {
385                index=targetElementIndex[i];
386                pack(targetElements[i], targetBuffer[index], indexTargetBuffer[index] );
387        }
388
389        indexRequest=0;
390
391        MPI_Barrier(communicator);
392        CTimer::get("CMPIRouting::transferToTarget(data)").reset();
393        CTimer::get("CMPIRouting::transferToTarget(data)").resume();
394        for(int i=0; i<nbSource; i++)
395        {
396                MPI_Irecv(sourceBuffer[i],sourceMessageSize[i],MPI_CHAR, sourceRank[i], 0, communicator, &request[indexRequest]);
397                indexRequest++;
398        }
399
400        for(int i=0;i<nbTarget; i++)
401        {
402                MPI_Isend(targetBuffer[i],targetMessageSize[i], MPI_CHAR, targetRank[i], 0, communicator, &request[indexRequest]);
403                indexRequest++;
404        }
405
406        MPI_Waitall(indexRequest,request,status);
407
408        MPI_Barrier(communicator);
409        CTimer::get("CMPIRouting::transferToTarget(data)").suspend();
410        CTimer::get("CMPIRouting::transferToTarget(data)").print();
411
412        // unpack the data
413        for(int i=0; i<totalSourceElement; i++)
414        {
415                index=sourceElementIndex[i];
416                unpack(sourceElements[i], sourceBuffer[index], indexSourceBuffer[index]);
417        }
418
419        for (int i=0; i<nbTarget; i++) delete[] targetBuffer[i];
420        for (int i=0; i<nbSource; i++) delete[] sourceBuffer[i];
421        delete[] targetBuffer;
422        delete[] indexTargetBuffer;
423        delete[] targetMessageSize;
424        delete[] sourceBuffer;
425        delete[] indexSourceBuffer;
426        delete[] sourceMessageSize;
427        delete[] request;
428        delete[] status;
429}
430
431template<typename T>
432void CMPIRouting::transferFromSource(T* targetElements, T* sourceElements)
433{
434        char** targetBuffer = new char*[nbTarget];
435        int* indexTargetBuffer = new int[nbTarget];
436
437        for (int i = 0; i < nbTarget; i++)
438        {
439                targetBuffer[i] = new char[sizeof(T)*nbTargetElement[i]];
440                indexTargetBuffer[i] = 0;
441        }
442
443        char** sourceBuffer = new char*[nbSource];
444        int* indexSourceBuffer = new int[nbSource];
445
446        for (int i = 0; i < nbSource; i++)
447        {
448                sourceBuffer[i] = new char[sizeof(T)*nbSourceElement[i]];
449                indexSourceBuffer[i] = 0;
450        }
451
452        // pack the data
453        int index;
454        for (int i = 0; i < totalSourceElement; i++)
455        {
456                index = sourceElementIndex[i];
457                *((T*) &(sourceBuffer[index][indexSourceBuffer[index]])) = sourceElements[i];
458                indexSourceBuffer[index] += sizeof(T);
459        }
460
461        MPI_Request* request=new MPI_Request[nbSource+nbTarget];
462        MPI_Status*  status=new MPI_Status[nbSource+nbTarget];
463        int indexRequest=0;
464
465        for(int i=0; i<nbSource; i++)
466        {
467                MPI_Isend(sourceBuffer[i],nbSourceElement[i]*sizeof(T),MPI_CHAR, sourceRank[i], 0, communicator, &request[indexRequest]);
468                indexRequest++;
469        }
470
471        for(int i=0;i<nbTarget; i++)
472        {
473                MPI_Irecv(targetBuffer[i],nbTargetElement[i]*sizeof(T), MPI_CHAR, targetRank[i], 0, communicator, &request[indexRequest]);
474                indexRequest++;
475        }
476
477        MPI_Waitall(indexRequest,request,status);
478
479        // unpack the data
480        for(int i=0;i<totalTargetElement;i++)
481        {
482                index=targetElementIndex[i];
483                targetElements[i] = *((T*) &(targetBuffer[index][indexTargetBuffer[index]]));
484                indexTargetBuffer[index]+=sizeof(T);
485        }
486
487        for(int i=0;i<nbTarget; i++) delete[] targetBuffer[i];
488        for(int i=0;i<nbSource; i++) delete[] sourceBuffer[i];
489        delete[] targetBuffer;
490        delete[] indexTargetBuffer;
491        delete[] sourceBuffer;
492        delete[] indexSourceBuffer;
493        delete[] request;
494        delete[] status;
495}
496
497
498/* number of source and target elements is known from previous call to init() */
499template<typename T, typename t_pack, typename t_unpack>
500void CMPIRouting::transferFromSource(T *targetElements, T *sourceElements, t_pack pack, t_unpack unpack)
501{
502        char **targetBuffer = new char*[nbTarget];
503        int *indexTargetBuffer = new int[nbTarget];
504        int *targetMessageSize = new int[nbTarget];
505        int *sourceMessageSize = new int[nbSource];
506        int index;
507
508        // compute target buffer size
509        for (int i = 0; i < nbSource; i++)  sourceMessageSize[i] = 0;
510
511        for (int i = 0; i < totalSourceElement; i++)
512        {
513                index = sourceElementIndex[i];
514                pack(sourceElements[i], NULL, sourceMessageSize[index]);
515        }
516
517        MPI_Request *request = new MPI_Request[nbSource + nbTarget];
518        MPI_Status  *status = new MPI_Status[nbSource + nbTarget];
519        int indexRequest = 0;
520        for (int i = 0; i < nbSource; i++)
521        {
522                MPI_Isend(&sourceMessageSize[i], 1, MPI_INT, sourceRank[i], 0, communicator, &request[indexRequest]);
523                indexRequest++;
524        }
525        for (int i = 0; i < nbTarget; i++)
526        {
527                MPI_Irecv(&targetMessageSize[i], 1, MPI_INT, targetRank[i], 0, communicator, &request[indexRequest]);
528                indexRequest++;
529        }
530        MPI_Waitall(indexRequest, request, status);
531
532        for (int i = 0; i < nbTarget; i++)
533        {
534                targetBuffer[i] = new char[targetMessageSize[i]];
535                indexTargetBuffer[i] = 0;
536        }
537
538        char **sourceBuffer = new char*[nbSource];
539        int *indexSourceBuffer = new int[nbSource];
540
541        for (int i = 0; i < nbSource; i++)
542        {
543                sourceBuffer[i] = new char[sourceMessageSize[i]];
544                indexSourceBuffer[i] = 0;
545        }
546
547
548        // pack the data
549        for (int i = 0; i < totalSourceElement; i++)
550        {
551                index = sourceElementIndex[i];
552                pack(sourceElements[i], sourceBuffer[index], indexSourceBuffer[index] );
553        }
554
555        indexRequest = 0;
556        for (int i = 0; i < nbSource; i++)
557        {
558                MPI_Isend(sourceBuffer[i], sourceMessageSize[i], MPI_CHAR, sourceRank[i], 0, communicator, &request[indexRequest]);
559                indexRequest++;
560        }
561        for (int i = 0; i < nbTarget; i++)
562        {
563                MPI_Irecv(targetBuffer[i], targetMessageSize[i], MPI_CHAR, targetRank[i], 0, communicator, &request[indexRequest]);
564                indexRequest++;
565        }
566        MPI_Waitall(indexRequest, request, status);
567
568        // unpack the data
569        for (int i = 0; i < totalTargetElement; i++)
570        {
571                index = targetElementIndex[i];
572                unpack(targetElements[i], targetBuffer[index], indexTargetBuffer[index]);
573        }
574
575        for (int i = 0; i < nbTarget; i++) delete[] targetBuffer[i];
576        for (int i = 0; i < nbSource; i++) delete[] sourceBuffer[i];
577        delete[] targetBuffer;
578        delete[] indexTargetBuffer;
579        delete[] targetMessageSize;
580        delete[] sourceBuffer;
581        delete[] indexSourceBuffer;
582        delete[] sourceMessageSize;
583        delete[] request;
584        delete[] status;
585}
586
587CMPIRouting::~CMPIRouting()
588{
589};
590
591template void CMPIRouting::init(const std::vector<int>& route, CMPICascade *cascade);
592template void CMPIRouting::init(const std::vector<vector<int> >& route, CMPICascade *cascade);
593
594template void CMPIRouting::transferToTarget(int *targetElements, int *sourceElements);
595template void CMPIRouting::transferToTarget(int *targetElements, int *sourceElements, void (*pack)(int&, char*, int&), void (* unpack)(int&, char *, int&));
596template void CMPIRouting::transferFromSource(int *targetElements, int *sourceElements);
597template void CMPIRouting::transferFromSource(int *targetElements, int *sourceElements, void (*pack)(int&, char*, int&), void (* unpack)(int&, char *, int&));
598
599template void CMPIRouting::transferToTarget(Node* targetElements, Node* sourceElements, void (*pack)(Node&, char*, int&), void (* unpack)(Node&, char *, int&));
600template void CMPIRouting::transferToTarget(Elt **targetElements, Elt **sourceElements, void (*pack)(Elt *, char*, int&), void (* unpack)(Elt *, char *, int&));
601template void CMPIRouting::transferFromSource(std::vector<int> *targetElements, std::vector<int> *sourceElements, void (*pack)(const std::vector<int>&, char*, int&), void (* unpack)(std::vector<int>&, const char *, int&));
602
603struct NES { int cnt; int risc; int rank; };
604
605template void alltoalls_unknown(const std::vector<std::vector<NES> >& send, std::vector<std::vector<NES> >& recv,
606                                const std::vector<int>& ranks, MPI_Comm communicator);
607
608template void alltoalls_known(const std::vector<std::vector<int> >& send, std::vector<std::vector<int> >& recv,
609                              const std::vector<int>& ranks, MPI_Comm communicator);
610
611}
Note: See TracBrowser for help on using the repository browser.