Attention une mise à jour du service Gitlab va être effectuée le mardi 18 janvier (et non lundi 17 comme annoncé précédemment) entre 18h00 et 18h30. Cette mise à jour va générer une interruption du service dont nous ne maîtrisons pas complètement la durée mais qui ne devrait pas excéder quelques minutes.

FMpiTreeBuilder.hpp 25.7 KB
Newer Older
1
2
3
#ifndef FMPITREEBUILDER_H
#define FMPITREEBUILDER_H

4
#include "../Utils/FMpi.hpp"
5
#include "../Utils/FQuickSort.hpp"
6
#include "../Utils/FMemUtils.hpp"
7
#include "../Utils/FTrace.hpp"
berenger-bramas's avatar
berenger-bramas committed
8
9
10
11
12
13

/** This class manage the loading of particles
  * for the mpi version.
  * It use a BinLoader and then sort the data
  * with a parallel quick sort.
  */
14
template<class ParticleClass>
15
class FMpiTreeBuilder{
berenger-bramas's avatar
berenger-bramas committed
16
17
18
    /** This method has been tacken from the octree
      * it computes a tree coordinate (x or y or z) from real position
      */
19
20
    static long getTreeCoordinate(const FReal inRelativePosition, const FReal boxWidthAtLeafLevel) {
            const FReal indexFReal = inRelativePosition / boxWidthAtLeafLevel;
21
22
            const long index = long(FMath::dfloor(indexFReal));
            if( index && FMath::LookEqual(inRelativePosition, boxWidthAtLeafLevel * FReal(index) ) ){
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
                    return index - 1;
            }
            return index;
    }

    /** get current rank */
    static int MpiGetRank(MPI_Comm comm = MPI_COMM_WORLD){
        int rank(0);
        MPI_Comm_rank(comm, &rank);
        return rank;
    }

    /** get current nb procs */
    static int MpiGetNbProcs(MPI_Comm comm = MPI_COMM_WORLD){
        int nb(0);
        MPI_Comm_size(comm, &nb);
        return nb;
    }

berenger-bramas's avatar
berenger-bramas committed
42
    /** receive data from a tag function */
43
44
45
46
47
48
49
50
51
    static void receiveDataFromTag(const int inSize, const int inTag, void* const inData, int* const inSource = 0, int* const inFilledSize = 0){
        MPI_Status status;
        MPI_Recv(inData, inSize, MPI_CHAR, MPI_ANY_SOURCE, inTag, MPI_COMM_WORLD, &status);
        if(inSource) *inSource = status.MPI_SOURCE;
        if(inFilledSize) MPI_Get_count(&status,MPI_CHAR,inFilledSize);
    }

    template< class T >
    static T GetLeft(const T inSize) {
52
53
        const float step = (float(inSize) / float(MpiGetNbProcs()));
        return T(FMath::Ceil(step * float(MpiGetRank())));
54
55
56
57
    }

    template< class T >
    static T GetRight(const T inSize) {
58
59
        const float step = (float(inSize) / float(MpiGetNbProcs()));
        const T res = T(FMath::Ceil(step * float(MpiGetRank()+1)));
60
61
62
63
64
65
66
67
68
69
70
71
        if(res > inSize) return inSize;
        else return res;
    }

    template< class T >
    static T GetOtherRight(const T inSize, const int other) {
        const float step = (float(inSize) / MpiGetNbProcs());
        const T res = T(FMath::Ceil(step * (other+1)));
        if(res > inSize) return inSize;
        else return res;
    }

berenger-bramas's avatar
berenger-bramas committed
72
73
74
    /** This struct is used to represent a particles group to
      * sort them easily
      */
75
76
77
78
79
80
81
82
83
84
    struct ParticlesGroup {
        int number;
        int positionInArray;
        MortonIndex index;
        ParticlesGroup(const int inNumber = 0 , const int inPositionInArray = 0, const MortonIndex inIndex = 0)
            : number(inNumber), positionInArray(inPositionInArray), index(inIndex) {
        }
    };


berenger-bramas's avatar
berenger-bramas committed
85
86
87
88
    /** A particle may not have a MortonIndex Method
      * but they are sorted on their morton index
      * so this struct store a particle and its index
      */
89
90
91
92
93
94
95
96
97
    struct IndexedParticle{
        MortonIndex index;
        ParticleClass particle;

        operator MortonIndex(){
            return this->index;
        }
    };

98
    char* intervals;
99
    FSize nbLeavesInIntervals;
100
101

private:
berenger-bramas's avatar
berenger-bramas committed
102
    // Forbid copy
103
104
105
    FMpiTreeBuilder(const FMpiTreeBuilder&){}
    FMpiTreeBuilder& operator=(const FMpiTreeBuilder&){return *this;}

106
public:
berenger-bramas's avatar
berenger-bramas committed
107
    /** Constructor */
108
109
110
111
    FMpiTreeBuilder()
        :  intervals(0), nbLeavesInIntervals(0) {
    }

berenger-bramas's avatar
berenger-bramas committed
112
    /** Destructor */
113
114
115
116
    virtual ~FMpiTreeBuilder(){
        delete [] intervals;
    }

berenger-bramas's avatar
berenger-bramas committed
117
    /** Split and sort the file */
118
119
    template <class LoaderClass>
    bool splitAndSortFile(LoaderClass& loader, const int NbLevels){
120
        FTRACE( FTrace::FFunction functionTrace(__FUNCTION__, "Loader" , __FILE__ , __LINE__) );
121
122
123
        const int rank = MpiGetRank();
        const int nbProcs = MpiGetNbProcs();

berenger-bramas's avatar
berenger-bramas committed
124
125
126
127
        // First we create the particles that belong to us (our proc)
        // we compute their morton index to be able to sort them
        //

128
        IndexedParticle* outputArray = 0;
129
        FSize outputSize = 0;
130
        {
131
            FTRACE( FTrace::FRegion regionTrace("Insert Particles", __FUNCTION__ , __FILE__ , __LINE__) );
132
133
            // create particles
            IndexedParticle*const realParticlesIndexed = new IndexedParticle[loader.getNumberOfParticles()];
134
            FMemUtils::memset(realParticlesIndexed, 0, sizeof(IndexedParticle) * loader.getNumberOfParticles());
135
136
            F3DPosition boxCorner(loader.getCenterOfBox() - (loader.getBoxWidth()/2));
            FTreeCoordinate host;
berenger-bramas's avatar
berenger-bramas committed
137

138
            const FReal boxWidthAtLeafLevel = loader.getBoxWidth() / FReal(1 << (NbLevels - 1) );
139
140
141
142
143
            for(long idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
                loader.fillParticle(realParticlesIndexed[idxPart].particle);
                host.setX( getTreeCoordinate( realParticlesIndexed[idxPart].particle.getPosition().getX() - boxCorner.getX(), boxWidthAtLeafLevel ));
                host.setY( getTreeCoordinate( realParticlesIndexed[idxPart].particle.getPosition().getY() - boxCorner.getY(), boxWidthAtLeafLevel ));
                host.setZ( getTreeCoordinate( realParticlesIndexed[idxPart].particle.getPosition().getZ() - boxCorner.getZ(), boxWidthAtLeafLevel ));
berenger-bramas's avatar
berenger-bramas committed
144

145
146
147
148
149
150
151
152
                realParticlesIndexed[idxPart].index = host.getMortonIndex(NbLevels - 1);
            }

            // sort particles
            FQuickSort::QsMpi<IndexedParticle,MortonIndex>(realParticlesIndexed, loader.getNumberOfParticles(),outputArray,outputSize);
            delete [] (realParticlesIndexed);
        }
        // be sure there is no splited leaves
berenger-bramas's avatar
berenger-bramas committed
153
        // to do that we exchange the first index with the left proc
154
        {
155
156
            FTRACE( FTrace::FRegion regionTrace("Remove Splited leaves", __FUNCTION__ , __FILE__ , __LINE__) );

157
158
159
160
            MortonIndex otherFirstIndex = -1;
            {
                FMpi::Request req[2];
                int reqiter = 0;
berenger-bramas's avatar
berenger-bramas committed
161
                // can I send my first index? == I am not rank 0 & I have data
162
                if( 0 < rank && outputSize){
163
                    MPI_Isend( &outputArray[0].index, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD, &req[reqiter++]);
164
165
                }
                if( rank != nbProcs - 1){
166
                    MPI_Irecv(&otherFirstIndex, 1, MPI_LONG_LONG, rank + 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD, &req[reqiter++]);
167
168
                }

169
                MPI_Waitall(reqiter,req,MPI_STATUSES_IGNORE);
170

berenger-bramas's avatar
berenger-bramas committed
171
172
                // I could not send because I do not have data, so I transmit the data coming
                // from my right neigbors
173
                if( 0 < rank && !outputSize){
174
                    MPI_Send( &otherFirstIndex, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD);
175
176
                }
            }
177

178
179
            MPI_Request req[2];
            int reqiter = 0;
180

181
            // at this point every one know the first index of his right neighbors
182
            const bool needToRecvBeforeSend = (rank != 0 && ((outputSize && otherFirstIndex == outputArray[0].index ) || !outputSize));
183
            if( needToRecvBeforeSend || (rank == nbProcs - 1) ){
berenger-bramas's avatar
berenger-bramas committed
184
185
186
187
188
                // Here we want to send data we do not have
                // so we first receive other data and put then into the array
                // this case happens only if I have one leaf with index MX
                // and if proc[rank - 1].last.index == MX && proc[rank + 1].first.index == MX

189
190
                int sendByOther = 0;

191
                MPI_Status probStatus;
192
                MPI_Probe(rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &probStatus);
193
                MPI_Get_count( &probStatus,  MPI_BYTE, &sendByOther);
194

195
                if(sendByOther){
196
                    sendByOther /= sizeof(IndexedParticle);
197
                    const IndexedParticle* const reallocOutputArray = outputArray;
198
                    const FSize reallocOutputSize = outputSize;
199
200
201

                    outputSize += sendByOther;
                    outputArray = new IndexedParticle[outputSize];
202
                    FMemUtils::memcpy(&outputArray[sendByOther], reallocOutputArray, reallocOutputSize * sizeof(IndexedParticle));
203
204
                    delete[] reallocOutputArray;

205
                    MPI_Recv(outputArray, sizeof(IndexedParticle) * sendByOther, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &probStatus);
206
207
                }
                else{
208
                    MPI_Irecv(0, 0, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &req[reqiter++]);
209
210
                }
            }
211

212
            if(rank != nbProcs - 1){
213

214
                FSize idxPart = outputSize - 1 ;
215
216
217
                while(idxPart >= 0 && outputArray[idxPart].index == otherFirstIndex){
                    --idxPart;
                }
218
                const int toSend = int(outputSize - 1 - idxPart);
219
                MPI_Isend( &outputArray[idxPart + 1], toSend * sizeof(IndexedParticle), MPI_BYTE, rank + 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &req[reqiter++]);
220

221
                if( rank != 0 && !needToRecvBeforeSend && (rank != nbProcs - 1)){
222
223
                    int sendByOther = 0;

224
                    MPI_Status probStatus;
225
                    MPI_Probe(rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &probStatus);
226
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sendByOther);
227

228
                    if(sendByOther){
229
230
231
                        sendByOther /= sizeof(IndexedParticle);
                        char* const tempBuffer = new char[sizeof(IndexedParticle) * sendByOther];

232
                        MPI_Irecv(tempBuffer, sizeof(IndexedParticle) * sendByOther, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &req[reqiter++]);
233

234
                        MPI_Waitall(reqiter,req, MPI_STATUSES_IGNORE);
235
236
                        reqiter = 0;

237
                        const IndexedParticle* const reallocOutputArray = outputArray;
238
                        const FSize reallocOutputSize = outputSize;
239
240
241

                        outputSize += sendByOther;
                        outputArray = new IndexedParticle[outputSize];
242
                        FMemUtils::memcpy(&outputArray[sendByOther], reallocOutputArray, reallocOutputSize * sizeof(IndexedParticle));
243
                        delete[] reallocOutputArray;
244
245
                        memcpy(outputArray, tempBuffer, sendByOther * sizeof(IndexedParticle));
                        delete[] tempBuffer;
246
                    }
247
                    else{
248
                        MPI_Irecv( 0, 0, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &req[reqiter++]);
249
250
251
                    }
                }
            }
252

253
            MPI_Waitall(reqiter,req,MPI_STATUSES_IGNORE);
254
255
        }

berenger-bramas's avatar
berenger-bramas committed
256
257
        // We now copy the data from a sorted type into real particles array + counter

258
259
260
        nbLeavesInIntervals = 0;
        if(outputSize){
            intervals = new char[outputSize * (sizeof(ParticleClass) + sizeof(int))];
261

262
263
264
            MortonIndex previousIndex = -1;
            char* writeIndex = intervals;
            int* writeCounter = 0;
265

266
            for( FSize idxPart = 0; idxPart < outputSize ; ++idxPart){
267
268
269
                if( outputArray[idxPart].index != previousIndex ){
                    previousIndex = outputArray[idxPart].index;
                    ++nbLeavesInIntervals;
270

271
272
                    writeCounter = reinterpret_cast<int*>( writeIndex );
                    writeIndex += sizeof(int);
273

274
275
                    (*writeCounter) = 0;
                }
276

277
                memcpy(writeIndex, &outputArray[idxPart].particle, sizeof(ParticleClass));
278

279
280
281
                writeIndex += sizeof(ParticleClass);
                ++(*writeCounter);
            }
282
        }
283
        delete [] outputArray;
284
285
286
287

        return true;
    }

berenger-bramas's avatar
berenger-bramas committed
288
    /** Put the interval into a tree */
289
    template <class OctreeClass>
290
    bool intervalsToTree(OctreeClass& realTree){
291
        FTRACE( FTrace::FFunction functionTrace(__FUNCTION__, "Loader" , __FILE__ , __LINE__) );
292
293
294
295
296
297
        const int rank = MpiGetRank();
        const int nbProcs = MpiGetNbProcs();

        //////////////////////////////////////////////////////////////////////////////////
        // We inform the master proc about the data we have
        //////////////////////////////////////////////////////////////////////////////////
298
        FTRACE( FTrace::FRegion preprocessTrace("Preprocess", __FUNCTION__ , __FILE__ , __LINE__) );
299

300
        FSize nbLeafs = nbLeavesInIntervals;
301
302
        FSize leftLeafs = 0;
        FSize rightLeafs = 0;
303
304
305

        // receive from left and right
        if((rank == 0)){
306
            MPI_Request requests[2];
307
308
            MPI_Isend(&nbLeafs, sizeof(FSize), MPI_BYTE , 1, FMpi::TagExchangeNbLeafs, MPI_COMM_WORLD, &requests[0]);
            MPI_Irecv(&rightLeafs, sizeof(FSize), MPI_BYTE, 1 , FMpi::TagExchangeNbLeafs , MPI_COMM_WORLD, &requests[1]);
309
            MPI_Waitall(2, requests, MPI_STATUSES_IGNORE);
310
311
        }
        else if(rank == nbProcs - 1){
312
            MPI_Request requests[2];
313
314
            MPI_Isend(&nbLeafs, sizeof(FSize), MPI_BYTE , rank - 1, FMpi::TagExchangeNbLeafs, MPI_COMM_WORLD, &requests[0]);
            MPI_Irecv(&leftLeafs, sizeof(FSize), MPI_BYTE, rank - 1 , FMpi::TagExchangeNbLeafs , MPI_COMM_WORLD, &requests[1]);
315
            MPI_Waitall(2, requests, MPI_STATUSES_IGNORE);
316
        }
317
        else { //rank != 0) && rank != nbProcs - 1
318
319
            for(int idxToReceive = 0 ; idxToReceive < 2 ; ++idxToReceive){
                int source(0);
320
                FSize temp = 0;
321
                receiveDataFromTag(sizeof(FSize), FMpi::TagExchangeNbLeafs, &temp, &source);
322
323
324
                if(source < rank){ // come from left
                    leftLeafs = temp;
                    temp += nbLeafs;
325
                    MPI_Send(&temp, sizeof(FSize), MPI_BYTE , rank + 1, FMpi::TagExchangeNbLeafs, MPI_COMM_WORLD);
326
327
328
329
                }
                else { // come from right
                    rightLeafs = temp;
                    temp += nbLeafs;
330
                    MPI_Send(&temp, sizeof(FSize), MPI_BYTE , rank - 1, FMpi::TagExchangeNbLeafs, MPI_COMM_WORLD);
331
332
333
                }
            }
        }
334
        FTRACE( preprocessTrace.end() );
335
336
337
338
        //////////////////////////////////////////////////////////////////////////////////
        // We balance the data
        //////////////////////////////////////////////////////////////////////////////////

339
340
341
        const FSize totalNbLeafs = (leftLeafs + nbLeafs + rightLeafs);
        const FSize myLeftLeaf = GetLeft(totalNbLeafs);
        const FSize myRightLeaf = GetRight(totalNbLeafs);
342
343
344
345
346
347
348
349
350
351
352

        const bool iNeedToSendToLeft = leftLeafs < myLeftLeaf;
        const bool iNeedToSendToRight = myRightLeaf < leftLeafs + nbLeafs;

        const bool iWillReceiveFromRight = leftLeafs + nbLeafs < myRightLeaf;
        const bool iWillReceiveFromLeft = leftLeafs > myLeftLeaf;

        const bool iDoNotHaveEnoughtToSendRight = myRightLeaf < leftLeafs;
        const bool iDoNotHaveEnoughtToSendLeft = leftLeafs + nbLeafs < myLeftLeaf;


353
354
        const FSize iNeedToSendLeftCount = myLeftLeaf - leftLeafs;
        const FSize iCanSendToLeft = nbLeafs;
355

356
357
        const FSize iNeedToSendRightCount = leftLeafs + nbLeafs - myRightLeaf;
        const FSize iCanSendToRight = nbLeafs;
358

359
        MPI_Request requests[2];
360
        int iterRequest = 0;
361

362
363
        FSize hasBeenSentToLeft = 0;
        FSize hasBeenSentToRight = 0;
364

365
        char* particlesToSend = 0;
366

367
368
369
        ///////////////////////////////
        // Manage data we already have
        ///////////////////////////////
370
        FTRACE( FTrace::FRegion step1Trace("Step1", __FUNCTION__ , __FILE__ , __LINE__) );
371

372
        if(nbLeafs){
373
            particlesToSend = intervals;
374

375
            FSize currentLeafPosition = 0;
376
377
378

            //Send to Left (the first leaves
            if(iNeedToSendToLeft){
379
                for(FSize idxLeaf = 0 ; idxLeaf < iNeedToSendLeftCount && idxLeaf < iCanSendToLeft ; ++idxLeaf){
380
                    currentLeafPosition += ((*(int*)&particlesToSend[currentLeafPosition]) * sizeof(ParticleClass)) + sizeof(int);
381
                }
382
                hasBeenSentToLeft = FMath::Min(iNeedToSendLeftCount, iCanSendToLeft);
383
                MPI_Isend(particlesToSend, int(currentLeafPosition), MPI_BYTE , rank - 1, FMpi::TagSandSettling, MPI_COMM_WORLD, &requests[iterRequest++]);
384
385
            }

386
            // Insert the particles I host and that belong to me
387
388
            const FSize beginForMe = (iNeedToSendToLeft ? FMath::Min(iNeedToSendLeftCount,iCanSendToLeft) : 0);
            const FSize endForMe = nbLeafs - (iNeedToSendToRight ? FMath::Min(iNeedToSendRightCount,iCanSendToRight) : 0);
389

390
            for(FSize idxLeaf = beginForMe ; idxLeaf < endForMe ; ++idxLeaf){
391

392
393
394
395
396
                const int nbPartInLeaf = (*(int*)&particlesToSend[currentLeafPosition]);
                ParticleClass* const particles = reinterpret_cast<ParticleClass*>(&particlesToSend[currentLeafPosition] + sizeof(int));

                for(int idxPart = 0 ; idxPart < nbPartInLeaf ; ++idxPart){
                    realTree.insert(particles[idxPart]);
397
                }
398

399
                currentLeafPosition += (nbPartInLeaf * sizeof(ParticleClass)) + sizeof(int);
400
401
            }

402
403
            //Send to Right (the right-est leaves
            if(iNeedToSendToRight){
404
                const FSize beginWriteIndex = currentLeafPosition;
405

406
                for(int idxLeaf = 0 ; idxLeaf < iNeedToSendRightCount && idxLeaf < iCanSendToRight ; ++idxLeaf){
407
                    currentLeafPosition += (*(int*)&particlesToSend[currentLeafPosition]* sizeof(ParticleClass)) + sizeof(int);
408
                }
409
410

                hasBeenSentToRight = FMath::Min(iNeedToSendRightCount, iCanSendToRight);
411
                MPI_Isend( &particlesToSend[beginWriteIndex], int(currentLeafPosition - beginWriteIndex), MPI_BYTE , rank + 1, FMpi::TagSandSettling,
412
                          MPI_COMM_WORLD, &requests[iterRequest++]);
413
414
415
            }
        }

416
417
418
        char* toRecvFromLeft = 0;
        char* toRecvFromRight = 0;
        int countReceive = int(iWillReceiveFromLeft) + int(iWillReceiveFromRight);
419
420
421
422
423
424
425
        int sizeOfLeftBuffer = 0;
        int sizeOfRightBuffer = 0;
        int sizeOfRightData = 0;
        int sizeOfLeftData = 0;

        int sourceToWhileRecv = MPI_ANY_SOURCE;

426
427
        // Now prepare to receive data
        while(countReceive--){
428
            MPI_Status recvStatus;
429
            MPI_Probe(sourceToWhileRecv, FMpi::TagSandSettling, MPI_COMM_WORLD, &recvStatus);
430
            // receive from left
431
432
            if(recvStatus.MPI_SOURCE == rank - 1){
                MPI_Get_count( &recvStatus,  MPI_BYTE, &sizeOfLeftBuffer);
433
434
                toRecvFromLeft = new char[sizeOfLeftBuffer];
                sizeOfLeftData = sizeOfLeftBuffer;
435
                MPI_Irecv(toRecvFromLeft, sizeOfLeftBuffer, MPI_BYTE, rank - 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, &requests[iterRequest++]);
436
                sourceToWhileRecv = rank + 1;
437
438
439
            }
            // receive from right
            else{
440
                MPI_Get_count( &recvStatus,  MPI_BYTE, &sizeOfRightBuffer);
441
442
                toRecvFromRight = new char[sizeOfRightBuffer];
                sizeOfRightData = sizeOfRightBuffer;
443
                MPI_Irecv(toRecvFromRight, sizeOfRightBuffer, MPI_BYTE, rank + 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, &requests[iterRequest++]);
444
                sourceToWhileRecv = rank - 1;
445
446
            }
        }
447
448
449
450

        ///////////////////////////////
        // Wait send receive
        ///////////////////////////////
451
        MPI_Waitall(iterRequest, requests, MPI_STATUSES_IGNORE);
452
        // We can delete the buffer use to send our particles only
453
454
        FTRACE( step1Trace.end() );
        FTRACE( FTrace::FRegion step2Trace("Step2", __FUNCTION__ , __FILE__ , __LINE__) );
455
456
457
458
459
        ///////////////////////////////
        // Process received data
        // and transfer if needed
        ///////////////////////////////
        // We have to receive from right and transfere to left
460
461
        int hasToBeReceivedFromLeft = int(leftLeafs - myLeftLeaf);
        int hasToBeReceivedFromRight = int(myRightLeaf - (leftLeafs + nbLeafs));
462
463
464
465
        int arrayIdxRight = 0;
        int arrayIdxLeft = 0;

        if(iDoNotHaveEnoughtToSendLeft){
466

467
468
            do{
                arrayIdxRight = 0;
469
                while(arrayIdxRight < sizeOfRightData && hasBeenSentToLeft < iNeedToSendLeftCount){
470
471
472
                    const int particlesInThisLeaf = *(int*)&toRecvFromRight[arrayIdxRight];
                    arrayIdxRight += sizeof(int) + sizeof(ParticleClass) * particlesInThisLeaf;
                    --hasToBeReceivedFromRight;
473
                    ++hasBeenSentToLeft;
474
                }
475

476
                MPI_Send(toRecvFromRight, arrayIdxRight, MPI_BYTE , rank - 1, FMpi::TagSandSettling, MPI_COMM_WORLD);
477
                if(hasBeenSentToLeft < iNeedToSendLeftCount){
478
                    MPI_Status probStatus;
479
                    MPI_Probe(MPI_ANY_SOURCE, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
480
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfRightData);
481
482
                    if(sizeOfRightBuffer < sizeOfRightData){
                        sizeOfRightBuffer = sizeOfRightData;
483
                        delete[] toRecvFromRight;
484
                        toRecvFromRight = new char[sizeOfRightData];
485
                    }
486

487
                    MPI_Recv(toRecvFromRight, sizeOfRightData, MPI_BYTE, rank + 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
488
                }
489
490
491
492
            } while(hasBeenSentToLeft < iNeedToSendLeftCount);
        }
        // We have to receive from left and transfere to right
        else if(iDoNotHaveEnoughtToSendRight){
493

494
495
            do{
                arrayIdxLeft = 0;
496
                while(arrayIdxLeft < sizeOfLeftData && hasBeenSentToRight < iNeedToSendRightCount){
497
498
499
                    const int particlesInThisLeaf = *(int*)&toRecvFromLeft[arrayIdxLeft];
                    arrayIdxLeft += sizeof(int) + sizeof(ParticleClass) * particlesInThisLeaf;
                    --hasToBeReceivedFromLeft;
500
                    ++hasBeenSentToRight;
501
                }
502

503
                MPI_Send(toRecvFromLeft, arrayIdxLeft, MPI_BYTE , rank + 1, FMpi::TagSandSettling, MPI_COMM_WORLD);
504
                if(hasBeenSentToRight < iNeedToSendRightCount){
505
                    MPI_Status probStatus;
506
                    MPI_Probe(MPI_ANY_SOURCE, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
507
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfLeftData);
508
509
                    if(sizeOfLeftBuffer < sizeOfLeftData){
                        sizeOfLeftBuffer = sizeOfLeftData;
510
                        delete[] toRecvFromLeft;
511
                        toRecvFromLeft = new char[sizeOfLeftData];
512
                    }
513
                    MPI_Recv(toRecvFromLeft, sizeOfLeftData, MPI_BYTE, rank - 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
514
515
516
517
518
519
                }
            } while(hasBeenSentToRight < iNeedToSendRightCount);
        }

        if(iWillReceiveFromLeft ){ // I need to wait
            do{
520
                while(arrayIdxLeft < sizeOfLeftData){
521
522
523
524
525
526
527
528
                    const int particlesInThisLeaf = *(int*)&toRecvFromLeft[arrayIdxLeft];
                    arrayIdxLeft += sizeof(int);
                    ParticleClass*const particles = reinterpret_cast<ParticleClass*>(&toRecvFromLeft[arrayIdxLeft]);
                    arrayIdxLeft += sizeof(ParticleClass) * particlesInThisLeaf;

                    for( int idxPart = 0 ; idxPart < particlesInThisLeaf ; ++idxPart){
                        realTree.insert( particles[ idxPart ] );
                    }
529

530
                    --hasToBeReceivedFromLeft;
531
                }
532
533
534
                arrayIdxLeft = 0;

                if(hasToBeReceivedFromLeft){
535
                    MPI_Status probStatus;
536
                    MPI_Probe( rank - 1, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
537
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfLeftData);
538
539
                    if(sizeOfLeftBuffer < sizeOfLeftData){
                        sizeOfLeftBuffer = sizeOfLeftData;
540
                        delete[] toRecvFromLeft;
541
                        toRecvFromLeft = new char[sizeOfLeftData];
542
                    }
543
                    MPI_Recv(toRecvFromLeft, sizeOfLeftData, MPI_BYTE, rank - 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
544
545
                }
            } while(hasToBeReceivedFromLeft);
546
        }
547

548
549
        if(iWillReceiveFromRight){
            do{
550
                while(arrayIdxRight < sizeOfRightData){
551
552
553
554
555
556
557
558
                    const int particlesInThisLeaf = *(int*)&toRecvFromRight[arrayIdxRight];
                    arrayIdxRight += sizeof(int);
                    ParticleClass*const particles = reinterpret_cast<ParticleClass*>(&toRecvFromRight[arrayIdxRight]);
                    arrayIdxRight += sizeof(ParticleClass) * particlesInThisLeaf;

                    for( int idxPart = 0 ; idxPart < particlesInThisLeaf ; ++idxPart){
                        realTree.insert( particles[ idxPart ] );
                    }
559

560
                    --hasToBeReceivedFromRight;
561
                }
562
563
564
                arrayIdxRight = 0;

                if(hasToBeReceivedFromRight){
565
                    MPI_Status probStatus;
566

567
                    MPI_Probe( rank + 1, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
568
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfRightData);
569

570
571
                    if(sizeOfRightBuffer < sizeOfRightData){
                        sizeOfRightBuffer = sizeOfRightData;
572
                        delete[] toRecvFromRight;
573
                        toRecvFromRight = new char[sizeOfRightData];
574
                    }
575
                    MPI_Recv(toRecvFromRight, sizeOfRightData, MPI_BYTE, rank + 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
576
577
                }
            } while(hasToBeReceivedFromRight);
578
579
        }

580
581
582
        // Delete reception buffers
        delete[] toRecvFromLeft;
        delete[] toRecvFromRight;
583
584
585
586
587
588

        return true;
    }
};

#endif // FMPITREEBUILDER_H