MAJ terminée. Nous sommes passés en version 14.6.2 . Pour consulter les "releases notes" associées c'est ici :

https://about.gitlab.com/releases/2022/01/11/security-release-gitlab-14-6-2-released/
https://about.gitlab.com/releases/2022/01/04/gitlab-14-6-1-released/

FMpiTreeBuilder.hpp 24.9 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

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
121
122
        const int rank = MpiGetRank();
        const int nbProcs = MpiGetNbProcs();

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

127
        IndexedParticle* outputArray = 0;
128
        FSize outputSize = 0;
129
130
131
        {
            // create particles
            IndexedParticle*const realParticlesIndexed = new IndexedParticle[loader.getNumberOfParticles()];
132
            FMemUtils::memset(realParticlesIndexed, 0, sizeof(IndexedParticle) * loader.getNumberOfParticles());
133
134
            F3DPosition boxCorner(loader.getCenterOfBox() - (loader.getBoxWidth()/2));
            FTreeCoordinate host;
berenger-bramas's avatar
berenger-bramas committed
135

136
            const FReal boxWidthAtLeafLevel = loader.getBoxWidth() / FReal(1 << (NbLevels - 1) );
137
138
139
140
141
            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
142

143
144
145
146
147
148
149
150
                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
151
        // to do that we exchange the first index with the left proc
152
153
154
155
156
        {
            MortonIndex otherFirstIndex = -1;
            {
                FMpi::Request req[2];
                int reqiter = 0;
berenger-bramas's avatar
berenger-bramas committed
157
                // can I send my first index? == I am not rank 0 & I have data
158
                if( 0 < rank && outputSize){
159
                    MPI_Isend( &outputArray[0].index, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD, &req[reqiter++]);
160
161
                }
                if( rank != nbProcs - 1){
162
                    MPI_Irecv(&otherFirstIndex, 1, MPI_LONG_LONG, rank + 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD, &req[reqiter++]);
163
164
                }

165
                MPI_Waitall(reqiter,req,MPI_STATUSES_IGNORE);
166

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

174
175
            MPI_Request req[2];
            int reqiter = 0;
176

177
            // at this point every one know the first index of his right neighbors
178
            const bool needToRecvBeforeSend = (rank != 0 && ((outputSize && otherFirstIndex == outputArray[0].index ) || !outputSize));
179
            if( needToRecvBeforeSend || (rank == nbProcs - 1) ){
berenger-bramas's avatar
berenger-bramas committed
180
181
182
183
184
                // 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

185
186
                int sendByOther = 0;

187
                MPI_Status probStatus;
188
                MPI_Probe(rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &probStatus);
189
                MPI_Get_count( &probStatus,  MPI_BYTE, &sendByOther);
190

191
                if(sendByOther){
192
                    sendByOther /= sizeof(IndexedParticle);
193
                    const IndexedParticle* const reallocOutputArray = outputArray;
194
                    const FSize reallocOutputSize = outputSize;
195
196
197

                    outputSize += sendByOther;
                    outputArray = new IndexedParticle[outputSize];
198
                    FMemUtils::memcpy(&outputArray[sendByOther], reallocOutputArray, reallocOutputSize * sizeof(IndexedParticle));
199
200
                    delete[] reallocOutputArray;

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

208
            if(rank != nbProcs - 1){
209

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

217
                if( rank != 0 && !needToRecvBeforeSend && (rank != nbProcs - 1)){
218
219
                    int sendByOther = 0;

220
                    MPI_Status probStatus;
221
                    MPI_Probe(rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &probStatus);
222
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sendByOther);
223

224
                    if(sendByOther){
225
226
227
                        sendByOther /= sizeof(IndexedParticle);
                        char* const tempBuffer = new char[sizeof(IndexedParticle) * sendByOther];

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

230
                        MPI_Waitall(reqiter,req, MPI_STATUSES_IGNORE);
231
232
                        reqiter = 0;

233
                        const IndexedParticle* const reallocOutputArray = outputArray;
234
                        const FSize reallocOutputSize = outputSize;
235
236
237

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

249
            MPI_Waitall(reqiter,req,MPI_STATUSES_IGNORE);
250
251
        }

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

254
255
256
        nbLeavesInIntervals = 0;
        if(outputSize){
            intervals = new char[outputSize * (sizeof(ParticleClass) + sizeof(int))];
257

258
259
260
            MortonIndex previousIndex = -1;
            char* writeIndex = intervals;
            int* writeCounter = 0;
261

262
            for( FSize idxPart = 0; idxPart < outputSize ; ++idxPart){
263
264
265
                if( outputArray[idxPart].index != previousIndex ){
                    previousIndex = outputArray[idxPart].index;
                    ++nbLeavesInIntervals;
266

267
268
                    writeCounter = reinterpret_cast<int*>( writeIndex );
                    writeIndex += sizeof(int);
269

270
271
                    (*writeCounter) = 0;
                }
272

273
                memcpy(writeIndex, &outputArray[idxPart].particle, sizeof(ParticleClass));
274

275
276
277
                writeIndex += sizeof(ParticleClass);
                ++(*writeCounter);
            }
278
        }
279
        delete [] outputArray;
280
281
282
283

        return true;
    }

berenger-bramas's avatar
berenger-bramas committed
284
    /** Put the interval into a tree */
285
    template <class OctreeClass>
286
    bool intervalsToTree(OctreeClass& realTree){
287
288
289
290
291
292
293
        const int rank = MpiGetRank();
        const int nbProcs = MpiGetNbProcs();

        //////////////////////////////////////////////////////////////////////////////////
        // We inform the master proc about the data we have
        //////////////////////////////////////////////////////////////////////////////////

294
        FSize nbLeafs = nbLeavesInIntervals;
295
296
        FSize leftLeafs = 0;
        FSize rightLeafs = 0;
297
298
299

        // receive from left and right
        if((rank == 0)){
300
            MPI_Request requests[2];
301
302
            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]);
303
            MPI_Waitall(2, requests, MPI_STATUSES_IGNORE);
304
305
        }
        else if(rank == nbProcs - 1){
306
            MPI_Request requests[2];
307
308
            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]);
309
            MPI_Waitall(2, requests, MPI_STATUSES_IGNORE);
310
        }
311
        else { //rank != 0) && rank != nbProcs - 1
312
313
            for(int idxToReceive = 0 ; idxToReceive < 2 ; ++idxToReceive){
                int source(0);
314
                FSize temp = 0;
315
                receiveDataFromTag(sizeof(FSize), FMpi::TagExchangeNbLeafs, &temp, &source);
316
317
318
                if(source < rank){ // come from left
                    leftLeafs = temp;
                    temp += nbLeafs;
319
                    MPI_Send(&temp, sizeof(FSize), MPI_BYTE , rank + 1, FMpi::TagExchangeNbLeafs, MPI_COMM_WORLD);
320
321
322
323
                }
                else { // come from right
                    rightLeafs = temp;
                    temp += nbLeafs;
324
                    MPI_Send(&temp, sizeof(FSize), MPI_BYTE , rank - 1, FMpi::TagExchangeNbLeafs, MPI_COMM_WORLD);
325
326
327
                }
            }
        }
328

329
330
331
332
        //////////////////////////////////////////////////////////////////////////////////
        // We balance the data
        //////////////////////////////////////////////////////////////////////////////////

333
334
335
        const FSize totalNbLeafs = (leftLeafs + nbLeafs + rightLeafs);
        const FSize myLeftLeaf = GetLeft(totalNbLeafs);
        const FSize myRightLeaf = GetRight(totalNbLeafs);
336
337
338
339
340
341
342
343
344
345
346

        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;


347
348
        const FSize iNeedToSendLeftCount = myLeftLeaf - leftLeafs;
        const FSize iCanSendToLeft = nbLeafs;
349

350
351
        const FSize iNeedToSendRightCount = leftLeafs + nbLeafs - myRightLeaf;
        const FSize iCanSendToRight = nbLeafs;
352

353
        MPI_Request requests[2];
354
        int iterRequest = 0;
355

356
357
        FSize hasBeenSentToLeft = 0;
        FSize hasBeenSentToRight = 0;
358

359
        char* particlesToSend = 0;
360

361
362
363
364
        ///////////////////////////////
        // Manage data we already have
        ///////////////////////////////

365
        if(nbLeafs){
366
            particlesToSend = intervals;
367

368
            FSize currentLeafPosition = 0;
369
370
371

            //Send to Left (the first leaves
            if(iNeedToSendToLeft){
372
                for(FSize idxLeaf = 0 ; idxLeaf < iNeedToSendLeftCount && idxLeaf < iCanSendToLeft ; ++idxLeaf){
373
                    currentLeafPosition += ((*(int*)&particlesToSend[currentLeafPosition]) * sizeof(ParticleClass)) + sizeof(int);
374
                }
375
                hasBeenSentToLeft = FMath::Min(iNeedToSendLeftCount, iCanSendToLeft);
376
                MPI_Isend(particlesToSend, int(currentLeafPosition), MPI_BYTE , rank - 1, FMpi::TagSandSettling, MPI_COMM_WORLD, &requests[iterRequest++]);
377
378
            }

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

383
            for(FSize idxLeaf = beginForMe ; idxLeaf < endForMe ; ++idxLeaf){
384

385
386
387
388
389
                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]);
390
                }
391

392
                currentLeafPosition += (nbPartInLeaf * sizeof(ParticleClass)) + sizeof(int);
393
394
            }

395
396
            //Send to Right (the right-est leaves
            if(iNeedToSendToRight){
397
                const FSize beginWriteIndex = currentLeafPosition;
398

399
                for(int idxLeaf = 0 ; idxLeaf < iNeedToSendRightCount && idxLeaf < iCanSendToRight ; ++idxLeaf){
400
                    currentLeafPosition += (*(int*)&particlesToSend[currentLeafPosition]* sizeof(ParticleClass)) + sizeof(int);
401
                }
402
403

                hasBeenSentToRight = FMath::Min(iNeedToSendRightCount, iCanSendToRight);
404
                MPI_Isend( &particlesToSend[beginWriteIndex], int(currentLeafPosition - beginWriteIndex), MPI_BYTE , rank + 1, FMpi::TagSandSettling,
405
                          MPI_COMM_WORLD, &requests[iterRequest++]);
406
407
408
            }
        }

409
410
411
        char* toRecvFromLeft = 0;
        char* toRecvFromRight = 0;
        int countReceive = int(iWillReceiveFromLeft) + int(iWillReceiveFromRight);
412
413
414
415
416
417
418
        int sizeOfLeftBuffer = 0;
        int sizeOfRightBuffer = 0;
        int sizeOfRightData = 0;
        int sizeOfLeftData = 0;

        int sourceToWhileRecv = MPI_ANY_SOURCE;

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

        ///////////////////////////////
        // Wait send receive
        ///////////////////////////////
444
        MPI_Waitall(iterRequest, requests, MPI_STATUSES_IGNORE);
445
        // We can delete the buffer use to send our particles only
446

447
448
449
450
451
        ///////////////////////////////
        // Process received data
        // and transfer if needed
        ///////////////////////////////
        // We have to receive from right and transfere to left
452
453
        int hasToBeReceivedFromLeft = int(leftLeafs - myLeftLeaf);
        int hasToBeReceivedFromRight = int(myRightLeaf - (leftLeafs + nbLeafs));
454
455
456
457
        int arrayIdxRight = 0;
        int arrayIdxLeft = 0;

        if(iDoNotHaveEnoughtToSendLeft){
458

459
460
            do{
                arrayIdxRight = 0;
461
                while(arrayIdxRight < sizeOfRightData && hasBeenSentToLeft < iNeedToSendLeftCount){
462
463
464
                    const int particlesInThisLeaf = *(int*)&toRecvFromRight[arrayIdxRight];
                    arrayIdxRight += sizeof(int) + sizeof(ParticleClass) * particlesInThisLeaf;
                    --hasToBeReceivedFromRight;
465
                    ++hasBeenSentToLeft;
466
                }
467

468
                MPI_Send(toRecvFromRight, arrayIdxRight, MPI_BYTE , rank - 1, FMpi::TagSandSettling, MPI_COMM_WORLD);
469
                if(hasBeenSentToLeft < iNeedToSendLeftCount){
470
                    MPI_Status probStatus;
471
                    MPI_Probe(MPI_ANY_SOURCE, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
472
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfRightData);
473
474
                    if(sizeOfRightBuffer < sizeOfRightData){
                        sizeOfRightBuffer = sizeOfRightData;
475
                        delete[] toRecvFromRight;
476
                        toRecvFromRight = new char[sizeOfRightData];
477
                    }
478

479
                    MPI_Recv(toRecvFromRight, sizeOfRightData, MPI_BYTE, rank + 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
480
                }
481
482
483
484
            } while(hasBeenSentToLeft < iNeedToSendLeftCount);
        }
        // We have to receive from left and transfere to right
        else if(iDoNotHaveEnoughtToSendRight){
485

486
487
            do{
                arrayIdxLeft = 0;
488
                while(arrayIdxLeft < sizeOfLeftData && hasBeenSentToRight < iNeedToSendRightCount){
489
490
491
                    const int particlesInThisLeaf = *(int*)&toRecvFromLeft[arrayIdxLeft];
                    arrayIdxLeft += sizeof(int) + sizeof(ParticleClass) * particlesInThisLeaf;
                    --hasToBeReceivedFromLeft;
492
                    ++hasBeenSentToRight;
493
                }
494

495
                MPI_Send(toRecvFromLeft, arrayIdxLeft, MPI_BYTE , rank + 1, FMpi::TagSandSettling, MPI_COMM_WORLD);
496
                if(hasBeenSentToRight < iNeedToSendRightCount){
497
                    MPI_Status probStatus;
498
                    MPI_Probe(MPI_ANY_SOURCE, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
499
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfLeftData);
500
501
                    if(sizeOfLeftBuffer < sizeOfLeftData){
                        sizeOfLeftBuffer = sizeOfLeftData;
502
                        delete[] toRecvFromLeft;
503
                        toRecvFromLeft = new char[sizeOfLeftData];
504
                    }
505
                    MPI_Recv(toRecvFromLeft, sizeOfLeftData, MPI_BYTE, rank - 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
506
507
508
509
510
511
                }
            } while(hasBeenSentToRight < iNeedToSendRightCount);
        }

        if(iWillReceiveFromLeft ){ // I need to wait
            do{
512
                while(arrayIdxLeft < sizeOfLeftData){
513
514
515
516
517
518
519
520
                    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 ] );
                    }
521

522
                    --hasToBeReceivedFromLeft;
523
                }
524
525
526
                arrayIdxLeft = 0;

                if(hasToBeReceivedFromLeft){
527
                    MPI_Status probStatus;
528
                    MPI_Probe( rank - 1, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
529
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfLeftData);
530
531
                    if(sizeOfLeftBuffer < sizeOfLeftData){
                        sizeOfLeftBuffer = sizeOfLeftData;
532
                        delete[] toRecvFromLeft;
533
                        toRecvFromLeft = new char[sizeOfLeftData];
534
                    }
535
                    MPI_Recv(toRecvFromLeft, sizeOfLeftData, MPI_BYTE, rank - 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
536
537
                }
            } while(hasToBeReceivedFromLeft);
538
        }
539
540
        if(iWillReceiveFromRight){
            do{
541
                while(arrayIdxRight < sizeOfRightData){
542
543
544
545
546
547
548
549
                    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 ] );
                    }
550

551
                    --hasToBeReceivedFromRight;
552
                }
553
554
555
                arrayIdxRight = 0;

                if(hasToBeReceivedFromRight){
556
                    MPI_Status probStatus;
557

558
                    MPI_Probe( rank + 1, FMpi::TagSandSettling, MPI_COMM_WORLD, &probStatus);
559
                    MPI_Get_count( &probStatus,  MPI_BYTE, &sizeOfRightData);
560

561
562
                    if(sizeOfRightBuffer < sizeOfRightData){
                        sizeOfRightBuffer = sizeOfRightData;
563
                        delete[] toRecvFromRight;
564
                        toRecvFromRight = new char[sizeOfRightData];
565
                    }
566
                    MPI_Recv(toRecvFromRight, sizeOfRightData, MPI_BYTE, rank + 1 , FMpi::TagSandSettling , MPI_COMM_WORLD, MPI_STATUS_IGNORE);
567
568
                }
            } while(hasToBeReceivedFromRight);
569
570
        }

571
572
573
        // Delete reception buffers
        delete[] toRecvFromLeft;
        delete[] toRecvFromRight;
574
575
576
577
578
579

        return true;
    }
};

#endif // FMPITREEBUILDER_H