source: XIOS/dev/branch_openmp/extern/remap/src/mpi_routing.cpp @ 1289

Last change on this file since 1289 was 1289, checked in by yushan, 6 years ago

EP update part 2

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