FMpiTreeBuilder.hpp 25 KB
Newer Older
1
// ===================================================================================
2
3
4
5
6
7
8
9
// Logiciel initial: ScalFmm Version 0.5
// Co-auteurs : Olivier Coulaud, Bérenger Bramas.
// Propriétaires : INRIA.
// Copyright © 2011-2012, diffusé sous les termes et conditions d’une licence propriétaire.
// Initial software: ScalFmm Version 0.5
// Co-authors: Olivier Coulaud, Bérenger Bramas.
// Owners: INRIA.
// Copyright © 2011-2012, spread under the terms and conditions of a proprietary license.
10
// ===================================================================================
11
12
13
#ifndef FMPITREEBUILDER_H
#define FMPITREEBUILDER_H

14
#include "../Utils/FMpi.hpp"
15
#include "../Utils/FQuickSortMpi.hpp"
16
#include "../Utils/FBitonicSort.hpp"
berenger-bramas's avatar
berenger-bramas committed
17

18
#include "../Utils/FMemUtils.hpp"
19
#include "../Utils/FTrace.hpp"
berenger-bramas's avatar
berenger-bramas committed
20

21
22
23
24
/** 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
  * or bitonic sort.
  * After it needs to merge the leaves to finally equalize the data.
berenger-bramas's avatar
berenger-bramas committed
25
  */
26
template<class ParticleClass>
27
class FMpiTreeBuilder{
28
29
30
31
32
33
34
35
public:
    /** What sorting algorithm to use */
    enum SortingType{
        QuickSort,
        BitonicSort,
    };

private:
berenger-bramas's avatar
berenger-bramas committed
36
37
38
    /** This method has been tacken from the octree
      * it computes a tree coordinate (x or y or z) from real position
      */
39
    static int getTreeCoordinate(const FReal inRelativePosition, const FReal boxWidthAtLeafLevel) {
40
            const FReal indexFReal = inRelativePosition / boxWidthAtLeafLevel;
41
            const int index = int(FMath::dfloor(indexFReal));
42
            if( index && FMath::LookEqual(inRelativePosition, boxWidthAtLeafLevel * FReal(index) ) ){
43
44
45
46
47
48
                    return index - 1;
            }
            return index;
    }


berenger-bramas's avatar
berenger-bramas committed
49
50
    /** A particle may not have a MortonIndex Method
      * but they are sorted on their morton index
51
      * so this struct store a particle + its index
berenger-bramas's avatar
berenger-bramas committed
52
      */
53
54
55
56
57
58
59
60
61
    struct IndexedParticle{
        MortonIndex index;
        ParticleClass particle;

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

62
63
64
    //////////////////////////////////////////////////////////////////////////
    // To sort tha particles we hold
    //////////////////////////////////////////////////////////////////////////
65
66


67
68
69
70
71
72
73
74
    template <class LoaderClass>
    static void SortParticles( const FMpi::FComm& communicator, LoaderClass& loader, const SortingType type,
                        const int TreeHeight, IndexedParticle**const outputArray, FSize* const outputSize){
        FTRACE( FTrace::FFunction functionTrace(__FUNCTION__ , "Loader to Tree" , __FILE__ , __LINE__) );

        // create particles
        IndexedParticle*const realParticlesIndexed = new IndexedParticle[loader.getNumberOfParticles()];
        FMemUtils::memset(realParticlesIndexed, 0, sizeof(IndexedParticle) * loader.getNumberOfParticles());
BRAMAS Berenger's avatar
BRAMAS Berenger committed
75
        FPoint boxCorner(loader.getCenterOfBox() - (loader.getBoxWidth()/2));
76
77
78
        FTreeCoordinate host;

        const FReal boxWidthAtLeafLevel = loader.getBoxWidth() / FReal(1 << (TreeHeight - 1) );
79
        for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
80
81
82
83
84
85
86
            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 ));

            realParticlesIndexed[idxPart].index = host.getMortonIndex(TreeHeight - 1);
        }
87

88
89
        // sort particles
        if(type == QuickSort){
90
            FQuickSortMpi<IndexedParticle,MortonIndex, FSize>::QsMpi(realParticlesIndexed, loader.getNumberOfParticles(), *outputArray, *outputSize,communicator);
91
92
93
94
95
96
97
            delete [] (realParticlesIndexed);
        }
        else {
            FBitonicSort<IndexedParticle,MortonIndex, FSize>::Sort( realParticlesIndexed, loader.getNumberOfParticles(), communicator );
            *outputArray = realParticlesIndexed;
            *outputSize = loader.getNumberOfParticles();
        }
98
99
    }

100
    static void SortParticlesFromArray( const FMpi::FComm& communicator, const ParticleClass array[], const int size, const SortingType type,
BRAMAS Berenger's avatar
BRAMAS Berenger committed
101
                                       const FPoint& centerOfBox, const FReal boxWidth,
102
103
104
105
106
107
108
                        const int TreeHeight, IndexedParticle**const outputArray, FSize* const outputSize){
        FTRACE( FTrace::FFunction functionTrace(__FUNCTION__ , "Loader to Tree" , __FILE__ , __LINE__) );

        // create particles
        IndexedParticle*const realParticlesIndexed = new IndexedParticle[size];
        FMemUtils::memset(realParticlesIndexed, 0, sizeof(IndexedParticle) * size);

BRAMAS Berenger's avatar
BRAMAS Berenger committed
109
        FPoint boxCorner(centerOfBox - (boxWidth/2));
110
111
112
113
        FTreeCoordinate host;

        const FReal boxWidthAtLeafLevel = boxWidth / FReal(1 << (TreeHeight - 1) );

114
        for(int idxPart = 0 ; idxPart < size ; ++idxPart){
115
116
117
118
119
120
121
122
123
124
            realParticlesIndexed[idxPart].particle = array[idxPart];
            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 ));

            realParticlesIndexed[idxPart].index = host.getMortonIndex(TreeHeight - 1);
        }

        // sort particles
        if(type == QuickSort){
125
            FQuickSortMpi<IndexedParticle,MortonIndex, FSize>::QsMpi(realParticlesIndexed, size, *outputArray, *outputSize,communicator);
126
127
128
129
130
131
132
133
134
            delete [] (realParticlesIndexed);
        }
        else {
            FBitonicSort<IndexedParticle,MortonIndex, FSize>::Sort( realParticlesIndexed, size, communicator );
            *outputArray = realParticlesIndexed;
            *outputSize = size;
        }
    }

135

136
137
138
    //////////////////////////////////////////////////////////////////////////
    // To merge the leaves
    //////////////////////////////////////////////////////////////////////////
berenger-bramas's avatar
berenger-bramas committed
139

140
141
142
143
144
    static void MergeLeaves(const FMpi::FComm& communicator, IndexedParticle*& workingArray, FSize& workingSize,
                            char** leavesArray, FSize* const leavesSize){
        FTRACE( FTrace::FFunction functionTrace(__FUNCTION__, "Loader to Tree" , __FILE__ , __LINE__) );
        const int rank = communicator.processId();
        const int nbProcs = communicator.processCount();
145
146

        // be sure there is no splited leaves
berenger-bramas's avatar
berenger-bramas committed
147
        // to do that we exchange the first index with the left proc
148
        {
149
150
            FTRACE( FTrace::FRegion regionTrace("Remove Splited leaves", __FUNCTION__ , __FILE__ , __LINE__) );

151
            MortonIndex otherFirstIndex = -1;
152
153
            if(workingSize != 0 && rank != 0 && rank != nbProcs - 1){
                MPI_Sendrecv(&workingArray[0].index, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs,
154
155
156
157
158
159
160
                             &otherFirstIndex, 1, MPI_LONG_LONG, rank + 1, FMpi::TagExchangeIndexs,
                             MPI_COMM_WORLD, MPI_STATUS_IGNORE);
            }
            else if( rank == 0){
                MPI_Recv(&otherFirstIndex, 1, MPI_LONG_LONG, rank + 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
            }
            else if( rank == nbProcs - 1){
161
                MPI_Send( &workingArray[0].index, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD);
162
163
164
165
            }
            else {
                MPI_Recv(&otherFirstIndex, 1, MPI_LONG_LONG, rank + 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
                MPI_Send(&otherFirstIndex, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD);
166
            }
167

168
            // at this point every one know the first index of his right neighbors
169
            const bool needToRecvBeforeSend = (rank != 0 && ((workingSize && otherFirstIndex == workingArray[0].index ) || !workingSize));
170
171
172
173
            MPI_Request requestSendLeaf;

            IndexedParticle* sendBuffer = 0;
            if(rank != nbProcs - 1 && needToRecvBeforeSend == false){
174
175
                FSize idxPart = workingSize - 1 ;
                while(idxPart >= 0 && workingArray[idxPart].index == otherFirstIndex){
176
177
                    --idxPart;
                }
178
                const int particlesToSend = int(workingSize - 1 - idxPart);
179
                if(particlesToSend){
180
                    workingSize -= particlesToSend;
181
                    sendBuffer = new IndexedParticle[particlesToSend];
182
                    memcpy(sendBuffer, &workingArray[idxPart + 1], particlesToSend * sizeof(IndexedParticle));
183

184
                    MPI_Isend( sendBuffer, particlesToSend * int(sizeof(IndexedParticle)), MPI_BYTE, rank + 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &requestSendLeaf);
185
186
187
188
189
                }
                else{
                    MPI_Isend( 0, 0, MPI_BYTE, rank + 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &requestSendLeaf);
                }
            }
berenger-bramas's avatar
berenger-bramas committed
190

191
            if( rank != 0 ){
192
193
                int sendByOther = 0;

194
                MPI_Status probStatus;
195
                MPI_Probe(rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &probStatus);
196
                MPI_Get_count( &probStatus,  MPI_BYTE, &sendByOther);
197

198
                if(sendByOther){
199
                    sendByOther /= int(sizeof(IndexedParticle));
200

201
202
                    const IndexedParticle* const reallocOutputArray = workingArray;
                    const FSize reallocOutputSize = workingSize;
203

204
205
206
                    workingSize += sendByOther;
                    workingArray = new IndexedParticle[workingSize];
                    FMemUtils::memcpy(&workingArray[sendByOther], reallocOutputArray, reallocOutputSize * sizeof(IndexedParticle));
207
208
                    delete[] reallocOutputArray;

209
                    MPI_Recv(workingArray, int(sizeof(IndexedParticle)) * sendByOther, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
210
211
                }
                else{
212
                    MPI_Recv( 0, 0, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
213
214
                }
            }
215

216
            if(rank != nbProcs - 1 && needToRecvBeforeSend == true){
217
218
219
220
                MPI_Send( workingArray, int(workingSize * sizeof(IndexedParticle)), MPI_BYTE, rank + 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD);
                delete[] workingArray;
                workingArray = 0;
                workingSize  = 0;
221
            }
222
223
224
225
            else if(rank != nbProcs - 1){
                MPI_Wait( &requestSendLeaf, MPI_STATUS_IGNORE);
                delete[] sendBuffer;
                sendBuffer = 0;
226
            }
227
228
        }

229
230
231
        {
            FTRACE( FTrace::FRegion regionTrace("Remove Splited leaves", __FUNCTION__ , __FILE__ , __LINE__) );
            // We now copy the data from a sorted type into real particles array + counter
berenger-bramas's avatar
berenger-bramas committed
232

233
234
235
236
237
            (*leavesSize)  = 0;
            (*leavesArray) = 0;

            if(workingSize){
                (*leavesArray) = new char[workingSize * (sizeof(ParticleClass) + sizeof(int))];
238

239
                MortonIndex previousIndex = -1;
240
                char* writeIndex = (*leavesArray);
241
                int* writeCounter = 0;
242

243
244
245
246
                for( FSize idxPart = 0; idxPart < workingSize ; ++idxPart){
                    if( workingArray[idxPart].index != previousIndex ){
                        previousIndex = workingArray[idxPart].index;
                        ++(*leavesSize);
247

248
249
                        writeCounter = reinterpret_cast<int*>( writeIndex );
                        writeIndex += sizeof(int);
250

251
252
                        (*writeCounter) = 0;
                    }
253

254
                    memcpy(writeIndex, &workingArray[idxPart].particle, sizeof(ParticleClass));
255

256
257
258
                    writeIndex += sizeof(ParticleClass);
                    ++(*writeCounter);
                }
259
            }
260
            delete [] workingArray;
261

262
263
264
            workingArray = 0;
            workingSize = 0;
        }
265
266
    }

267
268
269
270
    //////////////////////////////////////////////////////////////////////////
    // To equalize (same number of leaves among the procs)
    //////////////////////////////////////////////////////////////////////////

berenger-bramas's avatar
berenger-bramas committed
271
    /** Put the interval into a tree */
272
    template <class OctreeClass>
273
274
275
276
277
    static void EqualizeAndFillTree(const FMpi::FComm& communicator,  OctreeClass& realTree,
                             char*& leavesArray, FSize& nbLeavesInIntervals ){
        FTRACE( FTrace::FFunction functionTrace(__FUNCTION__, "Loader to Tree" , __FILE__ , __LINE__) );
        const int rank = communicator.processId();
        const int nbProcs = communicator.processCount();
278
        const FSize currentNbLeafs = nbLeavesInIntervals;
279

280
281
        // We have to know the number of leaves each procs holds
        FSize leavesPerProcs[nbProcs];
282
        memset(leavesPerProcs, 0, sizeof(int) * nbProcs);
283
        MPI_Allgather(&nbLeavesInIntervals, 1, MPI_LONG_LONG, leavesPerProcs, 1, MPI_LONG_LONG, MPI_COMM_WORLD);
284

285
        // Count the number of leaves on each side
286
287
288
289
290
291
        FSize currentLeafsOnMyLeft  = 0;
        FSize currentLeafsOnMyRight = 0;
        for(int idxProc = 0 ; idxProc < nbProcs ; ++idxProc ){
            if(idxProc < rank) currentLeafsOnMyLeft  += leavesPerProcs[idxProc];
            if(rank < idxProc) currentLeafsOnMyRight += leavesPerProcs[idxProc];
        }
292

293
294
295
        // So we can know the number of total leafs and
        // the number of leaves each procs must have at the end
        const FSize totalNbLeaves               = (currentLeafsOnMyLeft + currentNbLeafs + currentLeafsOnMyRight);
296
297
        const FSize correctLeftLeavesNumber     = communicator.getLeft(totalNbLeaves);
        const FSize correctRightLeavesIndex     = communicator.getRight(totalNbLeaves);
298

299
        // This will be used for the all to all
300
        int leavesToSend[nbProcs];
301
        memset(leavesToSend, 0, sizeof(int) * nbProcs);
302
303
304
305
        int bytesToSend[nbProcs];
        memset(bytesToSend, 0, sizeof(int) * nbProcs);
        int bytesOffset[nbProcs];
        memset(bytesToSend, 0, sizeof(int) * nbProcs);
306

307
        // Buffer Position
308
        FSize currentIntervalPosition = 0;
309

310
311
312
313
        // I need to send left if I hold leaves that belong to procs on my left
        const bool iNeedToSendToLeft  = currentLeafsOnMyLeft < correctLeftLeavesNumber;
        // I need to send right if I hold leaves that belong to procs on my right
        const bool iNeedToSendToRight = correctRightLeavesIndex < currentLeafsOnMyLeft + currentNbLeafs;
314

315
        if(iNeedToSendToLeft){
316
            FTRACE( FTrace::FRegion regionTrace("Calcul SendToLeft", __FUNCTION__ , __FILE__ , __LINE__) );
317
            // Find the first proc that need my data
318
            int idxProc = rank - 1;
319
            while( idxProc > 0 ){
320
                const FSize thisProcRight = communicator.getOtherRight(totalNbLeaves, idxProc - 1);
321
322
                // Go to left until proc-1 has a right index lower than my current left
                if( thisProcRight < currentLeafsOnMyLeft){
323
                    break;
324
                }
325
                --idxProc;
326
327
            }

328
329
            // Count data for this proc
            int ICanGive = int(currentNbLeafs);
330
331
            leavesToSend[idxProc] = int(FMath::Min(communicator.getOtherRight(totalNbLeaves, idxProc), totalNbLeaves - currentLeafsOnMyRight)
                                        - FMath::Max( currentLeafsOnMyLeft , communicator.getOtherLeft(totalNbLeaves, idxProc)));
332
333
334
            {
                bytesOffset[idxProc] = 0;
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
335
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
336
                }
337
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
338
            }
339
            ICanGive -= leavesToSend[idxProc];
340
            ++idxProc;
341

342
343
            // count data to other proc
            while(idxProc < rank && ICanGive){
344
                leavesToSend[idxProc] = int(FMath::Min( communicator.getOtherRight(totalNbLeaves, idxProc) - communicator.getOtherLeft(totalNbLeaves, idxProc), FSize(ICanGive)));
345

346
                bytesOffset[idxProc] = int(currentIntervalPosition);
347
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
348
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
349
                }
350
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
351

352
                ICanGive -= leavesToSend[idxProc];
353
                ++idxProc;
354
355
            }
        }
356

357
358
        // Store the index of my data but we do not insert the now
        const FSize myParticlesPosition = currentIntervalPosition;
359
        {
360
361
            FTRACE( FTrace::FRegion regionTrace("Jump My particles", __FUNCTION__ , __FILE__ , __LINE__) );
            const FSize iNeedToSendLeftCount = correctLeftLeavesNumber - currentLeafsOnMyLeft;
362
            FSize endForMe = currentNbLeafs;
363
364
365
366
            if(iNeedToSendToRight){
                const FSize iNeedToSendRightCount = currentLeafsOnMyLeft + currentNbLeafs - correctRightLeavesIndex;
                endForMe -= iNeedToSendRightCount;
            }
berenger-bramas's avatar
berenger-bramas committed
367

368
            // We have to jump the correct number of leaves
369
370
            for(FSize idxLeaf = FMath::Max(iNeedToSendLeftCount,FSize(0)) ; idxLeaf < endForMe ; ++idxLeaf){
                const int nbPartInLeaf = (*(int*)&leavesArray[currentIntervalPosition]);
371
372
373
                currentIntervalPosition += (nbPartInLeaf * sizeof(ParticleClass)) + sizeof(int);
            }
        }
374

375
        // Proceed same on the right
376
        if(iNeedToSendToRight){
377
            FTRACE( FTrace::FRegion regionTrace("Calcul SendToRight", __FUNCTION__ , __FILE__ , __LINE__) );
378
            // Find the last proc on the right that need my data
379
            int idxProc = rank + 1;
380
            while( idxProc < nbProcs ){
381
382
                const FSize thisProcLeft = communicator.getOtherLeft(totalNbLeaves, idxProc);
                const FSize thisProcRight = communicator.getOtherRight(totalNbLeaves, idxProc);
383
                // Progress until the proc+1 has its left index upper to my current right
berenger-bramas's avatar
berenger-bramas committed
384
                if( thisProcLeft < currentLeafsOnMyLeft || (totalNbLeaves - currentLeafsOnMyRight) < thisProcRight){
385
                    break;
386
                }
387
                ++idxProc;
388
            }
389

390
            // Count the data
berenger-bramas's avatar
berenger-bramas committed
391
            int ICanGive = int(currentLeafsOnMyLeft + currentNbLeafs - correctRightLeavesIndex);
392
393
            leavesToSend[idxProc] = int(FMath::Min(communicator.getOtherRight(totalNbLeaves, idxProc) , (totalNbLeaves - currentLeafsOnMyRight))
                                        - FMath::Max(communicator.getOtherLeft(totalNbLeaves, idxProc), currentLeafsOnMyLeft) );
berenger-bramas's avatar
berenger-bramas committed
394

395
            {
396
                bytesOffset[idxProc] = int(currentIntervalPosition);
397
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
398
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
399
                }
400
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
401
            }
402
            ICanGive -= leavesToSend[idxProc];
berenger-bramas's avatar
berenger-bramas committed
403
            ++idxProc;
404

405
406
            // Now Count the data to other
            while(idxProc < nbProcs && ICanGive){
407
                leavesToSend[idxProc] = int(FMath::Min( communicator.getOtherRight(totalNbLeaves, idxProc) - communicator.getOtherLeft(totalNbLeaves, idxProc), FSize(ICanGive)));
408

409
                bytesOffset[idxProc] = int(currentIntervalPosition);
410
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
411
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
412
                }
413
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
414

415
416
                ICanGive -= leavesToSend[idxProc];
                ++idxProc;
417
            }
418
        }
419

420
        // Inform other about who will send/receive what
421
422
423
        int bytesToSendRecv[nbProcs * nbProcs];
        memset(bytesToSendRecv, 0, sizeof(int) * nbProcs * nbProcs);
        MPI_Allgather(bytesToSend, nbProcs, MPI_INT, bytesToSendRecv, nbProcs, MPI_INT, MPI_COMM_WORLD);
424

425
426
427
428
        int bytesToRecv[nbProcs];
        memset(bytesToRecv, 0, sizeof(int) * nbProcs);
        int bytesOffsetToRecv[nbProcs];
        memset(bytesOffsetToRecv, 0, sizeof(int) * nbProcs);
429

430
        // Prepare needed buffer
431
432
433
        FSize sumBytesToRecv = 0;
        for(int idxProc = 0 ; idxProc < nbProcs ; ++idxProc){
            if( bytesToSendRecv[idxProc * nbProcs + rank] ){
434
435
                bytesOffsetToRecv[idxProc] = int(sumBytesToRecv);
                sumBytesToRecv += FSize(bytesToSendRecv[idxProc * nbProcs + rank]);
436
437
                bytesToRecv[idxProc] = bytesToSendRecv[idxProc * nbProcs + rank];
            }
438
439
        }

440
        // Send alll to  all
441
        char* const recvbuf = new char[sumBytesToRecv];
442
        MPI_Alltoallv(leavesArray, bytesToSend, bytesOffset, MPI_BYTE,
443
444
                      recvbuf, bytesToRecv, bytesOffsetToRecv, MPI_BYTE,
                      MPI_COMM_WORLD);
445

446
        { // Insert received data
447
448
449
450
451
            FTRACE( FTrace::FRegion regionTrace("Insert Received data", __FUNCTION__ , __FILE__ , __LINE__) );
            FSize recvBufferPosition = 0;
            while( recvBufferPosition < sumBytesToRecv){
                const int nbPartInLeaf = (*reinterpret_cast<int*>(&recvbuf[recvBufferPosition]));
                ParticleClass* const particles = reinterpret_cast<ParticleClass*>(&recvbuf[recvBufferPosition] + sizeof(int));
452

453
454
455
456
                for(int idxPart = 0 ; idxPart < nbPartInLeaf ; ++idxPart){
                    realTree.insert(particles[idxPart]);
                }
                recvBufferPosition += (nbPartInLeaf * sizeof(ParticleClass)) + sizeof(int);
berenger-bramas's avatar
berenger-bramas committed
457

458
459
            }
        }
460
        delete[] recvbuf;
461

462
463
464
465
466
467
468
469
470
471
472

        { // Insert my data
            FTRACE( FTrace::FRegion regionTrace("Insert My particles", __FUNCTION__ , __FILE__ , __LINE__) );
            currentIntervalPosition = myParticlesPosition;
            const FSize iNeedToSendLeftCount = correctLeftLeavesNumber - currentLeafsOnMyLeft;
            FSize endForMe = currentNbLeafs;
            if(iNeedToSendToRight){
                const FSize iNeedToSendRightCount = currentLeafsOnMyLeft + currentNbLeafs - correctRightLeavesIndex;
                endForMe -= iNeedToSendRightCount;
            }

473
474
475
            for(FSize idxLeaf = FMath::Max(iNeedToSendLeftCount,FSize(0)) ; idxLeaf < endForMe ; ++idxLeaf){
                const int nbPartInLeaf = (*(int*)&leavesArray[currentIntervalPosition]);
                ParticleClass* const particles = reinterpret_cast<ParticleClass*>(&leavesArray[currentIntervalPosition] + sizeof(int));
476
477
478
479
480

                for(int idxPart = 0 ; idxPart < nbPartInLeaf ; ++idxPart){
                    realTree.insert(particles[idxPart]);
                }
                currentIntervalPosition += (nbPartInLeaf * sizeof(ParticleClass)) + sizeof(int);
berenger-bramas's avatar
berenger-bramas committed
481

482
483
484
            }
        }

485
486
487
488
489
490
        delete[] leavesArray;
        leavesArray = 0;
        nbLeavesInIntervals = 0;
    }

public:
berenger-bramas's avatar
berenger-bramas committed
491

492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
    //////////////////////////////////////////////////////////////////////////
    // The builder function
    //////////////////////////////////////////////////////////////////////////

    template <class LoaderClass, class OctreeClass>
    static void LoaderToTree(const FMpi::FComm& communicator, LoaderClass& loader,
                             OctreeClass& realTree, const SortingType type = QuickSort){

        IndexedParticle* particlesArray = 0;
        FSize particlesSize = 0;
        SortParticles(communicator, loader, type, realTree.getHeight(), &particlesArray, &particlesSize);

        char* leavesArray = 0;
        FSize leavesSize = 0;
        MergeLeaves(communicator, particlesArray, particlesSize, &leavesArray, &leavesSize);

        EqualizeAndFillTree(communicator, realTree, leavesArray, leavesSize);
509
    }
510

511
512
    template <class OctreeClass>
    static void ArrayToTree(const FMpi::FComm& communicator, const ParticleClass array[], const int size,
BRAMAS Berenger's avatar
BRAMAS Berenger committed
513
                            const FPoint& boxCenter, const FReal boxWidth,
514
515
516
517
518
519
520
521
522
523
524
525
526
                             OctreeClass& realTree, const SortingType type = QuickSort){

        IndexedParticle* particlesArray = 0;
        FSize particlesSize = 0;
        SortParticlesFromArray(communicator, array, size, type, boxCenter, boxWidth, realTree.getHeight(), &particlesArray, &particlesSize);

        char* leavesArray = 0;
        FSize leavesSize = 0;
        MergeLeaves(communicator, particlesArray, particlesSize, &leavesArray, &leavesSize);

        EqualizeAndFillTree(communicator, realTree, leavesArray, leavesSize);
    }

527

528
529
530
};

#endif // FMPITREEBUILDER_H