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.4 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/FBitonicSort.hpp"
berenger-bramas's avatar
berenger-bramas committed
7

8
#include "../Utils/FMemUtils.hpp"
9
#include "../Utils/FTrace.hpp"
berenger-bramas's avatar
berenger-bramas committed
10

11
12
13
14
/** 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
15
  */
16
template<class ParticleClass>
17
class FMpiTreeBuilder{
18
19
20
21
22
23
24
25
public:
    /** What sorting algorithm to use */
    enum SortingType{
        QuickSort,
        BitonicSort,
    };

private:
berenger-bramas's avatar
berenger-bramas committed
26
27
28
    /** This method has been tacken from the octree
      * it computes a tree coordinate (x or y or z) from real position
      */
29
    static int getTreeCoordinate(const FReal inRelativePosition, const FReal boxWidthAtLeafLevel) {
30
            const FReal indexFReal = inRelativePosition / boxWidthAtLeafLevel;
31
            const int index = int(FMath::dfloor(indexFReal));
32
            if( index && FMath::LookEqual(inRelativePosition, boxWidthAtLeafLevel * FReal(index) ) ){
33
34
35
36
37
38
                    return index - 1;
            }
            return index;
    }


berenger-bramas's avatar
berenger-bramas committed
39
40
    /** A particle may not have a MortonIndex Method
      * but they are sorted on their morton index
41
      * so this struct store a particle + its index
berenger-bramas's avatar
berenger-bramas committed
42
      */
43
44
45
46
47
48
49
50
51
    struct IndexedParticle{
        MortonIndex index;
        ParticleClass particle;

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

52
53
54
    //////////////////////////////////////////////////////////////////////////
    // To sort tha particles we hold
    //////////////////////////////////////////////////////////////////////////
55
56


57
58
59
60
61
62
63
64
65
66
67
68
    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());
        F3DPosition boxCorner(loader.getCenterOfBox() - (loader.getBoxWidth()/2));
        FTreeCoordinate host;

        const FReal boxWidthAtLeafLevel = loader.getBoxWidth() / FReal(1 << (TreeHeight - 1) );
69
        for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
70
71
72
73
74
75
76
            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);
        }
77

78
79
80
81
82
83
84
85
86
87
        // sort particles
        if(type == QuickSort){
            FQuickSort<IndexedParticle,MortonIndex, FSize>::QsMpi(realParticlesIndexed, loader.getNumberOfParticles(), *outputArray, *outputSize,communicator);
            delete [] (realParticlesIndexed);
        }
        else {
            FBitonicSort<IndexedParticle,MortonIndex, FSize>::Sort( realParticlesIndexed, loader.getNumberOfParticles(), communicator );
            *outputArray = realParticlesIndexed;
            *outputSize = loader.getNumberOfParticles();
        }
88
89
    }

90
91
92
93
94
95
96
97
98
99
100
101
102
103
    static void SortParticlesFromArray( const FMpi::FComm& communicator, const ParticleClass array[], const int size, const SortingType type,
                                       const F3DPosition& centerOfBox, const FReal boxWidth,
                        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);

        F3DPosition boxCorner(centerOfBox - (boxWidth/2));
        FTreeCoordinate host;

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

104
        for(int idxPart = 0 ; idxPart < size ; ++idxPart){
105
106
107
108
109
110
111
112
113
114
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){
            FQuickSort<IndexedParticle,MortonIndex, FSize>::QsMpi(realParticlesIndexed, size, *outputArray, *outputSize,communicator);
            delete [] (realParticlesIndexed);
        }
        else {
            FBitonicSort<IndexedParticle,MortonIndex, FSize>::Sort( realParticlesIndexed, size, communicator );
            *outputArray = realParticlesIndexed;
            *outputSize = size;
        }
    }

125

126
127
128
    //////////////////////////////////////////////////////////////////////////
    // To merge the leaves
    //////////////////////////////////////////////////////////////////////////
berenger-bramas's avatar
berenger-bramas committed
129

130
131
132
133
134
    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();
135
136

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

141
            MortonIndex otherFirstIndex = -1;
142
143
            if(workingSize != 0 && rank != 0 && rank != nbProcs - 1){
                MPI_Sendrecv(&workingArray[0].index, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs,
144
145
146
147
148
149
150
                             &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){
151
                MPI_Send( &workingArray[0].index, 1, MPI_LONG_LONG, rank - 1, FMpi::TagExchangeIndexs, MPI_COMM_WORLD);
152
153
154
155
            }
            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);
156
            }
157

158
            // at this point every one know the first index of his right neighbors
159
            const bool needToRecvBeforeSend = (rank != 0 && ((workingSize && otherFirstIndex == workingArray[0].index ) || !workingSize));
160
161
162
163
            MPI_Request requestSendLeaf;

            IndexedParticle* sendBuffer = 0;
            if(rank != nbProcs - 1 && needToRecvBeforeSend == false){
164
165
                FSize idxPart = workingSize - 1 ;
                while(idxPart >= 0 && workingArray[idxPart].index == otherFirstIndex){
166
167
                    --idxPart;
                }
168
                const int particlesToSend = int(workingSize - 1 - idxPart);
169
                if(particlesToSend){
170
                    workingSize -= particlesToSend;
171
                    sendBuffer = new IndexedParticle[particlesToSend];
172
                    memcpy(sendBuffer, &workingArray[idxPart + 1], particlesToSend * sizeof(IndexedParticle));
173

174
                    MPI_Isend( sendBuffer, particlesToSend * int(sizeof(IndexedParticle)), MPI_BYTE, rank + 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &requestSendLeaf);
175
176
177
178
179
                }
                else{
                    MPI_Isend( 0, 0, MPI_BYTE, rank + 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &requestSendLeaf);
                }
            }
berenger-bramas's avatar
berenger-bramas committed
180

181
            if( rank != 0 ){
182
183
                int sendByOther = 0;

184
                MPI_Status probStatus;
185
                MPI_Probe(rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, &probStatus);
186
                MPI_Get_count( &probStatus,  MPI_BYTE, &sendByOther);
187

188
                if(sendByOther){
189
                    sendByOther /= int(sizeof(IndexedParticle));
190

191
192
                    const IndexedParticle* const reallocOutputArray = workingArray;
                    const FSize reallocOutputSize = workingSize;
193

194
195
196
                    workingSize += sendByOther;
                    workingArray = new IndexedParticle[workingSize];
                    FMemUtils::memcpy(&workingArray[sendByOther], reallocOutputArray, reallocOutputSize * sizeof(IndexedParticle));
197
198
                    delete[] reallocOutputArray;

199
                    MPI_Recv(workingArray, int(sizeof(IndexedParticle)) * sendByOther, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
200
201
                }
                else{
202
                    MPI_Recv( 0, 0, MPI_BYTE, rank - 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
203
204
                }
            }
205

206
            if(rank != nbProcs - 1 && needToRecvBeforeSend == true){
207
208
209
210
                MPI_Send( workingArray, int(workingSize * sizeof(IndexedParticle)), MPI_BYTE, rank + 1, FMpi::TagSplittedLeaf, MPI_COMM_WORLD);
                delete[] workingArray;
                workingArray = 0;
                workingSize  = 0;
211
            }
212
213
214
215
            else if(rank != nbProcs - 1){
                MPI_Wait( &requestSendLeaf, MPI_STATUS_IGNORE);
                delete[] sendBuffer;
                sendBuffer = 0;
216
            }
217
218
        }

219
220
221
        {
            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
222

223
224
225
226
227
            (*leavesSize)  = 0;
            (*leavesArray) = 0;

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

229
                MortonIndex previousIndex = -1;
230
                char* writeIndex = (*leavesArray);
231
                int* writeCounter = 0;
232

233
234
235
236
                for( FSize idxPart = 0; idxPart < workingSize ; ++idxPart){
                    if( workingArray[idxPart].index != previousIndex ){
                        previousIndex = workingArray[idxPart].index;
                        ++(*leavesSize);
237

238
239
                        writeCounter = reinterpret_cast<int*>( writeIndex );
                        writeIndex += sizeof(int);
240

241
242
                        (*writeCounter) = 0;
                    }
243

244
                    memcpy(writeIndex, &workingArray[idxPart].particle, sizeof(ParticleClass));
245

246
247
248
                    writeIndex += sizeof(ParticleClass);
                    ++(*writeCounter);
                }
249
            }
250
            delete [] workingArray;
251

252
253
254
            workingArray = 0;
            workingSize = 0;
        }
255
256
    }

257
258
259
260
    //////////////////////////////////////////////////////////////////////////
    // To equalize (same number of leaves among the procs)
    //////////////////////////////////////////////////////////////////////////

berenger-bramas's avatar
berenger-bramas committed
261
    /** Put the interval into a tree */
262
    template <class OctreeClass>
263
264
265
266
267
    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();
268
        const FSize currentNbLeafs = nbLeavesInIntervals;
269

270
271
        // We have to know the number of leaves each procs holds
        FSize leavesPerProcs[nbProcs];
272
        memset(leavesPerProcs, 0, sizeof(int) * nbProcs);
273
        MPI_Allgather(&nbLeavesInIntervals, 1, MPI_LONG_LONG, leavesPerProcs, 1, MPI_LONG_LONG, MPI_COMM_WORLD);
274

275
        // Count the number of leaves on each side
276
277
278
279
280
281
        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];
        }
282

283
284
285
        // 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);
286
287
        const FSize correctLeftLeavesNumber     = communicator.getLeft(totalNbLeaves);
        const FSize correctRightLeavesIndex     = communicator.getRight(totalNbLeaves);
288

289
        // This will be used for the all to all
290
        int leavesToSend[nbProcs];
291
        memset(leavesToSend, 0, sizeof(int) * nbProcs);
292
293
294
295
        int bytesToSend[nbProcs];
        memset(bytesToSend, 0, sizeof(int) * nbProcs);
        int bytesOffset[nbProcs];
        memset(bytesToSend, 0, sizeof(int) * nbProcs);
296

297
        // Buffer Position
298
        FSize currentIntervalPosition = 0;
299

300
301
302
303
        // 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;
304

305
        if(iNeedToSendToLeft){
306
            FTRACE( FTrace::FRegion regionTrace("Calcul SendToLeft", __FUNCTION__ , __FILE__ , __LINE__) );
307
            // Find the first proc that need my data
308
            int idxProc = rank - 1;
309
            while( idxProc > 0 ){
310
                const FSize thisProcRight = communicator.getOtherRight(totalNbLeaves, idxProc - 1);
311
312
                // Go to left until proc-1 has a right index lower than my current left
                if( thisProcRight < currentLeafsOnMyLeft){
313
                    break;
314
                }
315
                --idxProc;
316
317
            }

318
319
            // Count data for this proc
            int ICanGive = int(currentNbLeafs);
320
321
            leavesToSend[idxProc] = int(FMath::Min(communicator.getOtherRight(totalNbLeaves, idxProc), totalNbLeaves - currentLeafsOnMyRight)
                                        - FMath::Max( currentLeafsOnMyLeft , communicator.getOtherLeft(totalNbLeaves, idxProc)));
322
323
324
            {
                bytesOffset[idxProc] = 0;
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
325
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
326
                }
327
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
328
            }
329
            ICanGive -= leavesToSend[idxProc];
330
            ++idxProc;
331

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

336
                bytesOffset[idxProc] = int(currentIntervalPosition);
337
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
338
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
339
                }
340
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
341

342
                ICanGive -= leavesToSend[idxProc];
343
                ++idxProc;
344
345
            }
        }
346

347
348
        // Store the index of my data but we do not insert the now
        const FSize myParticlesPosition = currentIntervalPosition;
349
        {
350
351
            FTRACE( FTrace::FRegion regionTrace("Jump My particles", __FUNCTION__ , __FILE__ , __LINE__) );
            const FSize iNeedToSendLeftCount = correctLeftLeavesNumber - currentLeafsOnMyLeft;
352
            FSize endForMe = currentNbLeafs;
353
354
355
356
            if(iNeedToSendToRight){
                const FSize iNeedToSendRightCount = currentLeafsOnMyLeft + currentNbLeafs - correctRightLeavesIndex;
                endForMe -= iNeedToSendRightCount;
            }
berenger-bramas's avatar
berenger-bramas committed
357

358
            // We have to jump the correct number of leaves
359
360
            for(FSize idxLeaf = FMath::Max(iNeedToSendLeftCount,FSize(0)) ; idxLeaf < endForMe ; ++idxLeaf){
                const int nbPartInLeaf = (*(int*)&leavesArray[currentIntervalPosition]);
361
362
363
                currentIntervalPosition += (nbPartInLeaf * sizeof(ParticleClass)) + sizeof(int);
            }
        }
364

365
        // Proceed same on the right
366
        if(iNeedToSendToRight){
367
            FTRACE( FTrace::FRegion regionTrace("Calcul SendToRight", __FUNCTION__ , __FILE__ , __LINE__) );
368
            // Find the last proc on the right that need my data
369
            int idxProc = rank + 1;
370
            while( idxProc < nbProcs ){
371
372
                const FSize thisProcLeft = communicator.getOtherLeft(totalNbLeaves, idxProc);
                const FSize thisProcRight = communicator.getOtherRight(totalNbLeaves, idxProc);
373
                // Progress until the proc+1 has its left index upper to my current right
berenger-bramas's avatar
berenger-bramas committed
374
                if( thisProcLeft < currentLeafsOnMyLeft || (totalNbLeaves - currentLeafsOnMyRight) < thisProcRight){
375
                    break;
376
                }
377
                ++idxProc;
378
            }
379

380
            // Count the data
berenger-bramas's avatar
berenger-bramas committed
381
            int ICanGive = int(currentLeafsOnMyLeft + currentNbLeafs - correctRightLeavesIndex);
382
383
            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
384

385
            {
386
                bytesOffset[idxProc] = int(currentIntervalPosition);
387
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
388
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
389
                }
390
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
391
            }
392
            ICanGive -= leavesToSend[idxProc];
berenger-bramas's avatar
berenger-bramas committed
393
            ++idxProc;
394

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

399
                bytesOffset[idxProc] = int(currentIntervalPosition);
400
                for(FSize idxLeaf = 0 ; idxLeaf < leavesToSend[idxProc] ; ++idxLeaf){
401
                    currentIntervalPosition += ((*(int*)&leavesArray[currentIntervalPosition]) * sizeof(ParticleClass)) + sizeof(int);
402
                }
403
                bytesToSend[idxProc] = int(currentIntervalPosition - bytesOffset[idxProc]);
404

405
406
                ICanGive -= leavesToSend[idxProc];
                ++idxProc;
407
            }
408
        }
409

410
        // Inform other about who will send/receive what
411
412
413
        int bytesToSendRecv[nbProcs * nbProcs];
        memset(bytesToSendRecv, 0, sizeof(int) * nbProcs * nbProcs);
        MPI_Allgather(bytesToSend, nbProcs, MPI_INT, bytesToSendRecv, nbProcs, MPI_INT, MPI_COMM_WORLD);
414

415
416
417
418
        int bytesToRecv[nbProcs];
        memset(bytesToRecv, 0, sizeof(int) * nbProcs);
        int bytesOffsetToRecv[nbProcs];
        memset(bytesOffsetToRecv, 0, sizeof(int) * nbProcs);
419

420
        // Prepare needed buffer
421
422
423
        FSize sumBytesToRecv = 0;
        for(int idxProc = 0 ; idxProc < nbProcs ; ++idxProc){
            if( bytesToSendRecv[idxProc * nbProcs + rank] ){
424
425
                bytesOffsetToRecv[idxProc] = int(sumBytesToRecv);
                sumBytesToRecv += FSize(bytesToSendRecv[idxProc * nbProcs + rank]);
426
427
                bytesToRecv[idxProc] = bytesToSendRecv[idxProc * nbProcs + rank];
            }
428
429
        }

430
        // Send alll to  all
431
        char* const recvbuf = new char[sumBytesToRecv];
432
        MPI_Alltoallv(leavesArray, bytesToSend, bytesOffset, MPI_BYTE,
433
434
                      recvbuf, bytesToRecv, bytesOffsetToRecv, MPI_BYTE,
                      MPI_COMM_WORLD);
435

436
        { // Insert received data
437
438
439
440
441
            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));
442

443
444
445
446
                for(int idxPart = 0 ; idxPart < nbPartInLeaf ; ++idxPart){
                    realTree.insert(particles[idxPart]);
                }
                recvBufferPosition += (nbPartInLeaf * sizeof(ParticleClass)) + sizeof(int);
berenger-bramas's avatar
berenger-bramas committed
447

448
449
            }
        }
450
        delete[] recvbuf;
451

452
453
454
455
456
457
458
459
460
461
462

        { // 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;
            }

463
464
465
            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));
466
467
468
469
470

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

472
473
474
            }
        }

475
476
477
478
479
480
        delete[] leavesArray;
        leavesArray = 0;
        nbLeavesInIntervals = 0;
    }

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

482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
    //////////////////////////////////////////////////////////////////////////
    // 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);
499
    }
500

501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
    template <class OctreeClass>
    static void ArrayToTree(const FMpi::FComm& communicator, const ParticleClass array[], const int size,
                            const F3DPosition& boxCenter, const FReal boxWidth,
                             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);
    }

517

518
519
520
};

#endif // FMPITREEBUILDER_H