testFmmAlgorithmProcPeriodic.cpp 12.4 KB
Newer Older
1
// ===================================================================================
2 3 4 5 6 7 8 9 10 11 12 13 14
// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Bérenger Bramas, Matthias Messner
// olivier.coulaud@inria.fr, berenger.bramas@inria.fr
// This software is a computer program whose purpose is to compute the FMM.
//
// This software is governed by the CeCILL-C and LGPL licenses and
// abiding by the rules of distribution of free software.  
// 
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public and CeCILL-C Licenses for more details.
// "http://www.cecill.info". 
// "http://www.gnu.org/licenses".
15
// ===================================================================================
16

17 18 19
// ==== CMAKE =====
// @FUSE_MPI
// ================
20 21 22

#include <iostream>

23 24
#include <cstdio>
#include <cstdlib>
25

26 27
#include "../../Src/Utils/FParameters.hpp"
#include "../../Src/Utils/FTic.hpp"
28

29 30
#include "../../Src/Containers/FOctree.hpp"
#include "../../Src/Containers/FVector.hpp"
31

32
#include "../../Src/Components/FSimpleLeaf.hpp"
33
#include "../../Src/Files/FRandomLoader.hpp"
34

35
#include "../../Src/Utils/FPoint.hpp"
36

37
#include "../../Src/Components/FTestParticleContainer.hpp"
38 39
#include "../../Src/Components/FTestCell.hpp"
#include "../../Src/Components/FTestKernels.hpp"
40

41
#include "../../Src/Core/FFmmAlgorithmThreadProcPeriodic.hpp"
42
#include "../../Src/Core/FFmmAlgorithmPeriodic.hpp"
43

44
#include "../../Src/Files/FMpiTreeBuilder.hpp"
45

46
#include "../../Src/BalanceTree/FLeafBalance.hpp"
47 48 49 50 51 52 53 54 55

/** This program show an example of use of
  * the fmm basic algo
  * it also check that each particles is impacted each other particles
  */


// Simply create particles and try the kernels
int main(int argc, char ** argv){
berenger-bramas's avatar
berenger-bramas committed
56
    typedef FTestCell                   CellClass;
57
    typedef FTestParticleContainer      ContainerClass;
58

59 60 61
    typedef FSimpleLeaf< ContainerClass >                     LeafClass;
    typedef FOctree< CellClass, ContainerClass , LeafClass >  OctreeClass;
    typedef FTestKernels< CellClass, ContainerClass >         KernelClass;
62

63 64
    typedef FFmmAlgorithmThreadProcPeriodic<OctreeClass, CellClass, ContainerClass, KernelClass, LeafClass >     FmmClass;
    typedef FFmmAlgorithmPeriodic<OctreeClass, CellClass, ContainerClass, KernelClass, LeafClass >     FmmClassSeq;
65 66 67 68 69 70
    ///////////////////////What we do/////////////////////////////
    std::cout << ">> This executable has to be used to test the FMM algorithm.\n";
    //////////////////////////////////////////////////////////////

    const int NbLevels          = FParameters::getValue(argc,argv,"-h", 7);
    const int SizeSubLevels     = FParameters::getValue(argc,argv,"-sh", 3);
71
    const long NbParticles      = FParameters::getValue(argc,argv,"-nb", 5);
72
    const int PeriodicDeep      = FParameters::getValue(argc,argv,"-per", 2);
BRAMAS Berenger's avatar
BRAMAS Berenger committed
73

74 75 76 77 78 79 80 81

    FMpi app(argc, argv);

    FTic counter;

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

82
    std::cout << "Creating & Inserting " << NbParticles << " particles per boxes ..." << std::endl;
83 84 85
    std::cout << "\tHeight : " << NbLevels << " \t sub-height : " << SizeSubLevels << std::endl;
    counter.tic();

86
    FRandomLoader loader(NbParticles,FReal(1.0),FPoint(0,0,0), app.global().processId());
87 88
    OctreeClass tree(NbLevels, SizeSubLevels, loader.getBoxWidth(), loader.getCenterOfBox());

89
    {
90 91 92 93 94 95 96 97
        struct TestParticle{
            FPoint position;
            const FPoint& getPosition(){
                return position;
            }
        };

        TestParticle*const particles = new TestParticle[NbParticles];
98 99

        for(int idx = 0 ; idx < NbParticles ; ++idx){
100
            loader.fillParticle(&particles[idx].position);
101 102
        }

103
        FVector<TestParticle> finalParticles;
104
        FLeafBalance balancer;
105
        // FMpiTreeBuilder<TestParticle>::ArrayToTree(app.global(), particles, NbParticles, loader.getCenterOfBox(),
106 107 108 109 110 111
        // 					   loader.getBoxWidth(), tree.getHeight(), &finalParticles, &balancer);
        FMpiTreeBuilder< TestParticle >::DistributeArrayToContainer(app.global(),particles,
                                                                    NbParticles,
                                                                    loader.getCenterOfBox(),
                                                                    loader.getBoxWidth(),tree.getHeight(),
                                                                    &finalParticles, &balancer);
112 113 114 115

        for(int idx = 0 ; idx < finalParticles.getSize(); ++idx){
            tree.insert(finalParticles[idx].position);
        }
116 117 118 119 120 121 122 123 124 125 126 127 128
        delete[] particles;
    }

    counter.tac();
    std::cout << "Done  " << "(@Creating and Inserting Particles = " << counter.elapsed() << "s)." << std::endl;

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

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

    KernelClass kernels;
BRAMAS Berenger's avatar
BRAMAS Berenger committed
129
    FmmClass algo( app.global(), &tree, PeriodicDeep);
130
    algo.setKernel(&kernels);
131 132 133 134 135 136 137 138 139
    algo.execute();

    counter.tac();
    std::cout << "Done  " << "(@Algorithm = " << counter.elapsed() << "s)." << std::endl;

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

    {
BRAMAS Berenger's avatar
BRAMAS Berenger committed
140
        long long totalRepeatedBox = algo.theoricalRepetition();
141
        std::cout << "totalRepeatedBox in each dim is = " << totalRepeatedBox << "\n";
BRAMAS Berenger's avatar
BRAMAS Berenger committed
142
        totalRepeatedBox = (totalRepeatedBox*totalRepeatedBox*totalRepeatedBox);
143
        const long long NbParticlesEntireSystem = (NbParticles * app.global().processCount()) * totalRepeatedBox;
144
        std::cout << "The total number of particles is "  << NbParticlesEntireSystem << "\n";
145 146 147
        FTreeCoordinate min, max;
        algo.repetitionsIntervals(&min, &max);
        std::cout << "Min is " << min << " Max is " << max << std::endl;
148

149
        OctreeClass::Iterator octreeIterator(&tree);
150 151
        octreeIterator.gotoBottomLeft();
        do{
152 153
            ContainerClass* container = (octreeIterator.getCurrentListTargets());
            const long long int*const dataDown = container->getDataDown();
154

155 156
            for(int idxPart = 0 ; idxPart < container->getNbParticles() ; ++idxPart){
                if( NbParticlesEntireSystem - 1 != dataDown[idxPart]){
157
                    std::cout << "P2P probleme, should be " << NbParticlesEntireSystem - 1 <<
158
                                 " iter.data().getDataDown() "<< dataDown[idxPart] << std::endl;
159 160 161 162 163 164 165 166
                }
            }
        } while(octreeIterator.moveRight());
    }

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

167 168 169
    {
        OctreeClass treeSeq(NbLevels, SizeSubLevels, FReal(1.0), FPoint(0,0,0));
        for(int idx = 0 ; idx < app.global().processCount() ; ++idx ){
170 171 172 173 174 175
            FPoint position;
            FRandomLoader loaderSeq(NbParticles,FReal(1.0),FPoint(0,0,0), idx);
            for(int idxPart = 0 ; idxPart < loaderSeq.getNumberOfParticles() ; ++idxPart){
                loaderSeq.fillParticle(&position);
                treeSeq.insert(position);
            }
176 177
        }

BRAMAS Berenger's avatar
BRAMAS Berenger committed
178
        FmmClassSeq algoSeq( &treeSeq, PeriodicDeep);
179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200
        algoSeq.setKernel(&kernels);
        algoSeq.execute();

        { // Ceck if there is number of NbPart summed at level 1
            typename OctreeClass::Iterator octreeIterator(&tree);
            octreeIterator.gotoBottomLeft();

            typename OctreeClass::Iterator octreeIteratorSeq(&treeSeq);
            octreeIteratorSeq.gotoBottomLeft();

            for(int idxLevel = tree.getHeight() - 1 ; idxLevel >= 1 ; --idxLevel ){
                std::cout << "Process level " << idxLevel << std::endl;

                while( octreeIterator.getCurrentGlobalIndex() != octreeIteratorSeq.getCurrentGlobalIndex() ){
                    octreeIteratorSeq.moveRight();
                }

                do{
                    if( octreeIterator.getCurrentGlobalIndex() != octreeIteratorSeq.getCurrentGlobalIndex()){
                        std::cout << "Index problem !!!!!" << std::endl;
                    }

201
                    if( algo.getWorkingInterval(idxLevel).leftIndex <= octreeIteratorSeq.getCurrentGlobalIndex()){
202 203 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 229 230 231 232 233 234
                        if( octreeIterator.getCurrentCell()->getDataUp() != octreeIteratorSeq.getCurrentCell()->getDataUp() ){
                            std::cout << "Up problem at " << octreeIterator.getCurrentGlobalIndex() <<
                                         " Good is " << octreeIteratorSeq.getCurrentCell()->getDataUp() <<
                                         " Bad is " << octreeIterator.getCurrentCell()->getDataUp() << std::endl;
                        }
                        if( octreeIterator.getCurrentCell()->getDataDown() != octreeIteratorSeq.getCurrentCell()->getDataDown() ){
                            std::cout << "Down problem at " << octreeIterator.getCurrentGlobalIndex() <<
                                         " Good is " << octreeIteratorSeq.getCurrentCell()->getDataDown() <<
                                         " Bad is " << octreeIterator.getCurrentCell()->getDataDown() << std::endl;
                        }
                    }
                } while(octreeIterator.moveRight() && octreeIteratorSeq.moveRight());

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

                octreeIteratorSeq.moveUp();
                octreeIteratorSeq.gotoLeft();
            }
        }
        { // Check that each particle has been summed with all other
            typename OctreeClass::Iterator octreeIterator(&tree);
            octreeIterator.gotoBottomLeft();

            typename OctreeClass::Iterator octreeIteratorSeq(&treeSeq);
            octreeIteratorSeq.gotoBottomLeft();


            while( octreeIterator.getCurrentGlobalIndex() != octreeIteratorSeq.getCurrentGlobalIndex() ){
                octreeIteratorSeq.moveRight();
            }

            do{
235 236
                ContainerClass* container = (octreeIterator.getCurrentListTargets());
                const long long int*const dataDown = container->getDataDown();
237

238 239 240
                ContainerClass* containerValide = (octreeIteratorSeq.getCurrentListTargets());
                const long long int*const dataDownValide = containerValide->getDataDown();

241 242 243 244 245 246 247 248 249 250
                if( octreeIterator.getCurrentGlobalIndex() != octreeIteratorSeq.getCurrentGlobalIndex()){
                    std::cout << "Index problem !!!!!" << std::endl;
                }

                if(container->getNbParticles() != containerValide->getNbParticles()){
                    std::cout << "Not the same number of particles on the leaf " << octreeIterator.getCurrentGlobalIndex() << "\n";
                    std::cout << "\t Correct is " << containerValide->getNbParticles() << "\n";
                    std::cout << "\t Not Correct is " << container->getNbParticles() << "\n";
                }

251
                for(int idxPart = 0 ; idxPart < container->getNbParticles() ; ++idxPart){
252 253
                    // If a particles has been impacted by less than NbPart - 1 (the current particle)
                    // there is a problem
254 255 256 257 258 259 260 261
                    if( dataDown[idxPart] != dataDownValide[idxPart]){
                        std::cout << "Problem on leaf " << octreeIterator.getCurrentGlobalIndex() <<
                                     " part " << idxPart << " valide data down " << dataDownValide[idxPart] <<
                                     " invalide " << dataDown[idxPart] << "\n";
                        std::cout << "Data down for leaf is: valide " << octreeIteratorSeq.getCurrentCell()->getDataDown()
                                  << " invalide " << octreeIterator.getCurrentCell()->getDataDown()
                                  << " size is: valide " <<  octreeIteratorSeq.getCurrentListTargets()->getNbParticles()
                                  << " invalide " << octreeIterator.getCurrentListTargets()->getNbParticles() << std::endl;
262 263
                    }
                }
264
            } while(octreeIterator.moveRight() && octreeIteratorSeq.moveRight());
265 266 267 268
        }
    }
    std::cout << "Test is over...\n";

269 270 271 272
    return 0;
}


273