testFmmAlgorithmProc.cpp 18.4 KB
Newer Older
1
2
// /!\ Please, you must read the license at the bottom of this page

3
#include "../Src/Utils/FMpi.hpp"
4
#include "../Src/Utils/FTic.hpp"
5

6
#include "../Src/Containers/FOctree.hpp"
7
#include "../Src/Containers/FVector.hpp"
8
#include "../Src/Utils/FParameters.hpp"
9
#include "../Src/Utils/FGlobal.hpp"
10

11
#include "../Src/Components/FSimpleLeaf.hpp"
12

13
#include "../Src/Utils/F3DPosition.hpp"
14
#include "../Src/Utils/FAbstractSendable.hpp"
15

16
17
18
19
#include "../Src/Components/FFmaParticle.hpp"
#include "../Src/Components/FTestParticle.hpp"
#include "../Src/Components/FTestCell.hpp"
#include "../Src/Components/FTestKernels.hpp"
20
#include "../Src/Extensions/FExtendPhysicalValue.hpp"
21

22
#include "../Src/Core/FFmmAlgorithmThreadProc.hpp"
berenger-bramas's avatar
berenger-bramas committed
23
24
#include "../Src/Core/FFmmAlgorithmThread.hpp"

25
#include "../Src/Files/FFmaBinLoader.hpp"
26
27
#include "../Src/Files/FMpiFmaLoader.hpp"
#include "../Src/Files/FMpiTreeBuilder.hpp"
28

29
#include "../Src/Components/FBasicKernels.hpp"
30

berenger-bramas's avatar
berenger-bramas committed
31
32
33
34
35
36
#include <iostream>

#include <stdio.h>
#include <stdlib.h>


37
// Compile by : g++ testFmmAlgorithmProc.cpp ../Src/Utils/FDebug.cpp ../Src/Utils/FTrace.cpp -lgomp -fopenmp -O2 -o testFmmAlgorithmProc.exe
38
39

/** This program show an example of use of
berenger-bramas's avatar
berenger-bramas committed
40
  * the fmm threaded + mpi algo
41
42
43
  * it also check that each particles is impacted each other particles
  */

44
45
46
/////////////////////////////////////////////////////////////////////////////
// Test function
/////////////////////////////////////////////////////////////////////////////
47

48
49
50
// Check if tree is built correctly
template<class OctreeClass>
void ValidateTree(OctreeClass& realTree,
51
                        OctreeClass& treeValide, const FMpi::FComm& comm){
52
    FSize totalNbLeafs = 0;
53
    {
54

55
56
57
58
59
        typename OctreeClass::Iterator octreeIterator(&treeValide);
        octreeIterator.gotoBottomLeft();
        do {
            ++totalNbLeafs;
        } while(octreeIterator.moveRight());
60
61
    }

62
63
    const FSize myLeftLeaf = comm.getLeft(totalNbLeafs);
    const FSize myRightLeaf = comm.getRight(totalNbLeafs);
64

65
    //printf("%d should go from %d to %d leaf (on %d total leafs)\n", comm.processId(), myLeftLeaf, myRightLeaf, totalNbLeafs);
66
67
68

    typename OctreeClass::Iterator octreeIteratorValide(&treeValide);
    octreeIteratorValide.gotoBottomLeft();
69
    for(FSize idxLeaf = 0 ; idxLeaf < myLeftLeaf ; ++idxLeaf){
70
        if(!octreeIteratorValide.moveRight()){
71
            printf("Error cannot access to the left leaf %lld in the valide tree\n", myLeftLeaf);
72
        }
73
    }
74
75
76
77

    typename OctreeClass::Iterator octreeIterator(&realTree);
    octreeIterator.gotoBottomLeft();

78
    for(FSize idxLeaf = myLeftLeaf ; idxLeaf < myRightLeaf ; ++idxLeaf){
79
80
81
82
83
84
85
86
87
88
89
90
91
        if(octreeIteratorValide.getCurrentGlobalIndex() != octreeIterator.getCurrentGlobalIndex()){
            printf("Error index are different, valide %lld invalid %lld\n",octreeIteratorValide.getCurrentGlobalIndex(),
                   octreeIterator.getCurrentGlobalIndex());
            break;
        }
        if(octreeIteratorValide.getCurrentListSrc()->getSize() != octreeIterator.getCurrentListSrc()->getSize()){
            printf("Error leafs do not have the same number of particles, valide %d, invalide %d\n",
                   octreeIteratorValide.getCurrentListSrc()->getSize(), octreeIterator.getCurrentListSrc()->getSize() );
        }

        //printf("index %lld with %d particles\n", octreeIteratorValide.getCurrentGlobalIndex(), octreeIteratorValide.getCurrentListSrc()->getSize());

        if(!octreeIteratorValide.moveRight() && idxLeaf != myRightLeaf - 1){
92
            printf("Error cannot valide tree end to early, idxLeaf %lld myRightLeaf %lld\n", idxLeaf, myRightLeaf);
93
94
95
96
            break;
        }

        if(!octreeIterator.moveRight() && idxLeaf != myRightLeaf - 1){
97
            printf("Error cannot test tree end to early, idxLeaf %lld myRightLeaf %lld\n", idxLeaf, myRightLeaf);
98
99
            break;
        }
100
    }
101
102

}
103
104


105

berenger-bramas's avatar
berenger-bramas committed
106
/** This function tests the octree to be sure that the fmm algorithm
107
108
  * has worked completly.
  */
109
template<class OctreeClass, class ContainerClass, class FmmClassProc>
berenger-bramas's avatar
berenger-bramas committed
110
void ValidateFMMAlgoProc(OctreeClass* const badTree,
111
112
                         OctreeClass* const valideTree,
                         FmmClassProc* const fmm){
113
114
    const int OctreeHeight = badTree->getHeight();
    {
berenger-bramas's avatar
berenger-bramas committed
115
        typename OctreeClass::Iterator octreeIterator(badTree);
116
117
        octreeIterator.gotoBottomLeft();

berenger-bramas's avatar
berenger-bramas committed
118
        typename OctreeClass::Iterator octreeIteratorValide(valideTree);
119
120
        octreeIteratorValide.gotoBottomLeft();

121
        for(int level = OctreeHeight - 1 ; level > 0 && fmm->hasWorkAtLevel(level) ; --level){
122

123
            while(octreeIteratorValide.getCurrentGlobalIndex() != octreeIterator.getCurrentGlobalIndex()) {
124
125
126
                octreeIteratorValide.moveRight();
            }

127
128
129
130
131
            while(octreeIteratorValide.getCurrentGlobalIndex() != fmm->getWorkingInterval(level).min){
                octreeIteratorValide.moveRight();
                octreeIterator.moveRight();
            }

132
            FSize countCheck = 0;
133
            do{
134
135
136
137
                if(octreeIterator.getCurrentGlobalIndex() != octreeIteratorValide.getCurrentGlobalIndex()){
                    std::cout << "Error index are not equal!" << std::endl;
                }
                else{
138
                    if(octreeIterator.getCurrentCell()->getDataUp() != octreeIteratorValide.getCurrentCell()->getDataUp()){
139
                        std::cout << "M2M error at level " << level << " up bad " << octreeIterator.getCurrentCell()->getDataUp()
140
                                << " good " << octreeIteratorValide.getCurrentCell()->getDataUp() << " index " << octreeIterator.getCurrentGlobalIndex() << std::endl;
141
                    }
142
143
                    if(octreeIterator.getCurrentCell()->getDataDown() != octreeIteratorValide.getCurrentCell()->getDataDown()){
                        std::cout << "L2L error at level " << level << " down bad " << octreeIterator.getCurrentCell()->getDataDown()
144
                                << " good " << octreeIteratorValide.getCurrentCell()->getDataDown() << " index " << octreeIterator.getCurrentGlobalIndex() << std::endl;
145
146
                    }
                }
147
148
                ++countCheck;
            } while(octreeIteratorValide.moveRight() && octreeIterator.moveRight());
149

150
            // Check that each particle has been summed with all other
151
152
153
154
155
156
157

            octreeIterator.moveUp();
            octreeIterator.gotoLeft();

            octreeIteratorValide.moveUp();
            octreeIteratorValide.gotoLeft();
        }
158
    }
159

160
    {
161
162
        FSize NbPart = 0;
        FSize NbLeafs = 0;
163
        { // Check that each particle has been summed with all other
164
            typename OctreeClass::Iterator octreeIterator(valideTree);
165
166
167
168
169
170
171
172
            octreeIterator.gotoBottomLeft();
            do{
                NbPart += octreeIterator.getCurrentListSrc()->getSize();
                ++NbLeafs;
            } while(octreeIterator.moveRight());
        }
        {
            // Check that each particle has been summed with all other
berenger-bramas's avatar
berenger-bramas committed
173
            typename OctreeClass::Iterator octreeIterator(badTree);
174
175
            octreeIterator.gotoBottomLeft();

176
            do {
berenger-bramas's avatar
berenger-bramas committed
177
                typename ContainerClass::BasicIterator iter(*octreeIterator.getCurrentListTargets());
178
                const bool isUsingTsm = (octreeIterator.getCurrentListTargets() != octreeIterator.getCurrentListSrc());
179
                for(FSize idxPart = 0 ; idxPart < octreeIterator.getCurrentListTargets()->getSize() ; ++idxPart){
180
181
                    // If a particles has been impacted by less than NbPart - 1 (the current particle)
                    // there is a problem
berenger-bramas's avatar
berenger-bramas committed
182
183
                    if( (!isUsingTsm && iter.data().getDataDown() != NbPart - 1) ||
                        (isUsingTsm && iter.data().getDataDown() != NbPart) ){
184
185
                        std::cout << "Problem L2P + P2P, value on particle is : " << iter.data().getDataDown() <<
                                     " at pos " << idxPart << " index is " << octreeIterator.getCurrentGlobalIndex() << "\n";
186
                    }
berenger-bramas's avatar
berenger-bramas committed
187
                    iter.gotoNext();
188
                }
189
            } while( octreeIterator.moveRight());
190
191
        }
    }
192
193
194
195
196
197
198
199
    {
        {
            // Check that each particle has been summed with all other
            typename OctreeClass::Iterator octreeIterator(badTree);
            octreeIterator.gotoBottomLeft();

            do {
                if(octreeIterator.getCurrentListSrc()->getSize() != octreeIterator.getCurrentCell()->getDataUp()){
200
                    printf("P2M problem nb part %d data up %lld \n",
201
202
                           octreeIterator.getCurrentListSrc()->getSize(), octreeIterator.getCurrentCell()->getDataUp());
                }
203
            } while( octreeIterator.moveRight() );
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
        }
    }

    {
        // Check that each particle has been summed with all other
        typename OctreeClass::Iterator octreeIterator(badTree);
        octreeIterator.gotoBottomLeft();

        typename OctreeClass::Iterator valideOctreeIterator(valideTree);
        valideOctreeIterator.gotoBottomLeft();
        while(valideOctreeIterator.getCurrentGlobalIndex() != octreeIterator.getCurrentGlobalIndex()){
            valideOctreeIterator.moveRight();
        }

        do {
            if(valideOctreeIterator.getCurrentGlobalIndex() != octreeIterator.getCurrentGlobalIndex()){
                printf("Do not have the same index valide %lld invalide %lld \n",
                       valideOctreeIterator.getCurrentGlobalIndex(), octreeIterator.getCurrentGlobalIndex());
                break;
            }

            if(octreeIterator.getCurrentListTargets()->getSize() != valideOctreeIterator.getCurrentListTargets()->getSize()){
                printf("Do not have the same number of particle at leaf id %lld, valide %d invalide %d \n",
                       octreeIterator.getCurrentGlobalIndex(), valideOctreeIterator.getCurrentListTargets()->getSize(), octreeIterator.getCurrentListTargets()->getSize());
            }
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
            else {
                typename ContainerClass::BasicIterator iter(*octreeIterator.getCurrentListTargets());
                typename ContainerClass::BasicIterator iterValide(*valideOctreeIterator.getCurrentListTargets());

                for(int idxPart = 0 ; idxPart < octreeIterator.getCurrentListTargets()->getSize() ; ++idxPart){
                    // If a particles has been impacted by less than NbPart - 1 (the current particle)
                    // there is a problem
                    if( iter.data().getDataDown() != iterValide.data().getDataDown()){
                        std::cout << "Problem on leaf " << octreeIterator.getCurrentGlobalIndex() <<
                                     " part " << idxPart << " valide data down " << iterValide.data().getDataDown() <<
                                     " invalide " << iter.data().getDataDown() << "\n";
                        std::cout << "Data down for leaf is: valide " << valideOctreeIterator.getCurrentCell()->getDataDown()
                                  << " invalide " << octreeIterator.getCurrentCell()->getDataDown()
                                  << " size is: valide " <<  valideOctreeIterator.getCurrentListTargets()->getSize()
                                  << " invalide " << octreeIterator.getCurrentListTargets()->getSize() << std::endl;
                    }
                    iter.gotoNext();
                    iterValide.gotoNext();
                }
            }

250
251
        }while( octreeIterator.moveRight() && valideOctreeIterator.moveRight());
    }
252
253
254

}

255

256
257
258
/** To print an octree
  * used to debug and understand how the values were passed
  */
259
260
261
262
263
264
265
266
267
268
template<class OctreeClass>
void print(OctreeClass* const valideTree){
    typename OctreeClass::Iterator octreeIterator(valideTree);
    for(int idxLevel = valideTree->getHeight() - 1 ; idxLevel > 1 ; --idxLevel ){
        do{
            std::cout << "[" << octreeIterator.getCurrentGlobalIndex() << "] up:" << octreeIterator.getCurrentCell()->getDataUp() << " down:" << octreeIterator.getCurrentCell()->getDataDown() << "\t";
        } while(octreeIterator.moveRight());
        std::cout << "\n";
        octreeIterator.gotoLeft();
        octreeIterator.moveDown();
269
    }
270
271
272
273
274
275
276
277
278
279
}

/////////////////////////////////////////////////////////////////////
// Types
/////////////////////////////////////////////////////////////////////

/** Fmb class has to extend {FExtendForces,FExtendPotential,FExtendPhysicalValue}
  * Because we use fma loader it needs {FExtendPhysicalValue}
  */
class TestParticle : public FTestParticle, public FExtendPhysicalValue {
280
281
};

282

283
284
class TestCell : public FTestCell , public FAbstractSendable {
public:
285
    static const int SerializedSizeUp = sizeof(long long int);
286
    void serializeUp(void* const buffer) const {
287
        *(long long int*)buffer = this->dataUp;
288
289
    }
    void deserializeUp(const void* const buffer){
290
        this->dataUp = *(long long int*)buffer;
291
292
    }

293
    static const int SerializedSizeDown = sizeof(long long int);
294
    void serializeDown(void* const buffer) const {
295
        *(long long int*)buffer = this->dataDown;
296
297
    }
    void deserializeDown(const void* const buffer){
298
        this->dataDown = *(long long int*)buffer;
299
300
301
    }
};

berenger-bramas's avatar
berenger-bramas committed
302
303
304
305
/////////////////////////////////////////////////////////////////////
// Define the classes to use
/////////////////////////////////////////////////////////////////////

306
typedef TestParticle               ParticleClass;
307
typedef TestCell                   CellClass;
308
309
310
311
312
typedef FVector<ParticleClass>     ContainerClass;

typedef FSimpleLeaf<ParticleClass, ContainerClass >                     LeafClass;
typedef FOctree<ParticleClass, CellClass, ContainerClass , LeafClass >  OctreeClass;
typedef FTestKernels<ParticleClass, CellClass, ContainerClass >         KernelClass;
berenger-bramas's avatar
berenger-bramas committed
313

314
315
typedef FFmmAlgorithmThread<OctreeClass, ParticleClass, CellClass, ContainerClass, KernelClass, LeafClass >     FmmClass;
typedef FFmmAlgorithmThreadProc<OctreeClass, ParticleClass, CellClass, ContainerClass, KernelClass, LeafClass >     FmmClassProc;
berenger-bramas's avatar
berenger-bramas committed
316

317
318
319
/////////////////////////////////////////////////////////////////////
// Main
/////////////////////////////////////////////////////////////////////
320
321
322

// Simply create particles and try the kernels
int main(int argc, char ** argv){
323
324
325
326
    ///////////////////////What we do/////////////////////////////
    std::cout << ">> This executable has to be used to test the FMM algorithm.\n";
    //////////////////////////////////////////////////////////////

327
328
    FMpi app( argc, argv);

329
330
    const int NbLevels = FParameters::getValue(argc,argv,"-h", 9);
    const int SizeSubLevels = FParameters::getValue(argc,argv,"-sh", 3);
331
    char defaultFilename[] = "testLoaderFMA.fma"; //../../Data/ "testLoaderFMA.fma" "testFMAlgorithm.fma" Sphere.fma
332
    char* filename;
333
334
    FTic counter;

335
336
337
338
339
340
341
342
343
    if(argc == 1){
        std::cout << "You have to give a .fma file in argument.\n";
        std::cout << "The program will try a default file : " << defaultFilename << "\n";
        filename = defaultFilename;
    }
    else{
        filename = argv[1];
        std::cout << "Opening : " << filename << "\n";
    }
344

345
    FMpiFmaLoader<ParticleClass> loader(filename,app.global());
346
    if(!loader.isOpen()){
347
348
349
        std::cout << "Loader Error, " << filename << " is missing\n";
        return 1;
    }
350

351
    // The real tree to work on
352
353
    OctreeClass realTree(NbLevels, SizeSubLevels,loader.getBoxWidth(),loader.getCenterOfBox());

354
    if( app.global().processCount() != 1){
355
        //////////////////////////////////////////////////////////////////////////////////
356
        // Build tree from mpi loader
357
        //////////////////////////////////////////////////////////////////////////////////
358
        std::cout << "Build Tree ..." << std::endl;
359
        counter.tic();
360

361
        FMpiTreeBuilder<ParticleClass>::LoaderToTree(app.global(), loader, realTree);
362

363
364
        counter.tac();
        std::cout << "Done  " << "(" << counter.elapsed() << "s)." << std::endl;
365

366
        //////////////////////////////////////////////////////////////////////////////////
367
368
369
    }    
    else{
        ParticleClass partToInsert;
370
371
        for(FSize idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
            loader.fillParticle(partToInsert);
372
373
            realTree.insert(partToInsert);
        }
374
375
    }

376
    //////////////////////////////////////////////////////////////////////////////////
377
    // Create real tree
378
379
    //////////////////////////////////////////////////////////////////////////////////

berenger-bramas's avatar
berenger-bramas committed
380
    OctreeClass treeValide(NbLevels, SizeSubLevels,loader.getBoxWidth(),loader.getCenterOfBox());
381
382
383
    {
        FFmaBinLoader<ParticleClass> loaderSeq(filename);
        ParticleClass partToInsert;
384
        for(FSize idxPart = 0 ; idxPart < loaderSeq.getNumberOfParticles() ; ++idxPart){
385
386
387
388
            loaderSeq.fillParticle(partToInsert);
            treeValide.insert(partToInsert);
        }
    }
389

390
    //////////////////////////////////////////////////////////////////////////////////
391
    // Check particles in tree
392
    //////////////////////////////////////////////////////////////////////////////////
393
394
    std::cout << "Validate tree ..." << std::endl;
    counter.tic();
395

396
    ValidateTree(realTree, treeValide, app.global());
397

398
399
    counter.tac();
    std::cout << "Done  " << "(" << counter.elapsed() << "s)." << std::endl;
400
401
402

    //////////////////////////////////////////////////////////////////////////////////

403
    std::cout << "Working parallel particles ..." << std::endl;
404
405
    counter.tic();

berenger-bramas's avatar
berenger-bramas committed
406
    KernelClass kernels;
berenger-bramas's avatar
berenger-bramas committed
407

408
    FmmClassProc algo(app.global(),&realTree,&kernels);
409
410
    algo.execute();

411
412
413
414
415
416
417
418
    counter.tac();
    std::cout << "Done  " << "(@Algorithm Particles = " << counter.elapsed() << "s)." << std::endl;

    //////////////////////////////////////////////////////////////////////////////////

    std::cout << "Working sequential particles ..." << std::endl;
    counter.tic();

419
420
    FmmClass algoValide(&treeValide,&kernels);
    algoValide.execute();
berenger-bramas's avatar
berenger-bramas committed
421

422
    counter.tac();
423
    std::cout << "Done  " << "(@Algorithm Particles = " << counter.elapsed() << "s)." << std::endl;
424
425
426
427

    //////////////////////////////////////////////////////////////////////////////////
    //////////////////////////////////////////////////////////////////////////////////

428
429
430
    std::cout << "Checking data ..." << std::endl;
    counter.tic();

431
    ValidateFMMAlgoProc<OctreeClass,ContainerClass, FmmClassProc>(&realTree,&treeValide,&algo);
432

433
434
435
    counter.tac();
    std::cout << "Done  " << "(" << counter.elapsed() << "s)." << std::endl;

436
437
438
439
440
441
442
443
    //////////////////////////////////////////////////////////////////////////////////
    //////////////////////////////////////////////////////////////////////////////////

    return 0;
}


// [--LICENSE--]