testSphericalRotationAlgorithm.cpp 5.7 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 20 21
// ===================================================================================

#include <iostream>

#include <cstdio>
#include <cstdlib>

22 23
#include "../../Src/Utils/FTic.hpp"
#include "../../Src/Utils/FParameters.hpp"
24

25 26
#include "../../Src/Containers/FOctree.hpp"
#include "../../Src/Containers/FVector.hpp"
27

28 29 30
#include "../../Src/Core/FFmmAlgorithm.hpp"
#include "../../Src/Core/FFmmAlgorithmThread.hpp"
#include "../../Src/Core/FFmmAlgorithmTask.hpp"
31

32
#include "../../Src/Components/FSimpleLeaf.hpp"
33

34 35 36
#include "../../Src/Kernels/Spherical/FSphericalRotationKernel.hpp"
#include "../../Src/Kernels/Spherical/FSphericalCell.hpp"
#include "../../Src/Kernels/Spherical/FSphericalParticle.hpp"
37

38
#include "../../Src/Files/FFmaScanfLoader.hpp"
39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57

/** This program show an example of use of
  * the fmm basic algo
  * it also check that eachh particles is little or longer
  * related that each other
  */


// Simply create particles and try the kernels
int main(int argc, char ** argv){
    typedef FSphericalParticle             ParticleClass;
    typedef FSphericalCell                 CellClass;
    typedef FVector<ParticleClass>         ContainerClass;

    typedef FSimpleLeaf<ParticleClass, ContainerClass >                     LeafClass;
    typedef FOctree<ParticleClass, CellClass, ContainerClass , LeafClass >  OctreeClass;
    typedef FSphericalRotationKernel<ParticleClass, CellClass, ContainerClass > KernelClass;

    typedef FFmmAlgorithm<OctreeClass, ParticleClass, CellClass, ContainerClass, KernelClass, LeafClass > FmmClass;
58 59
    typedef FFmmAlgorithmThread<OctreeClass, ParticleClass, CellClass, ContainerClass, KernelClass, LeafClass > FmmClassThread;
    typedef FFmmAlgorithmTask<OctreeClass, ParticleClass, CellClass, ContainerClass, KernelClass, LeafClass > FmmClassTask;
60
    ///////////////////////What we do/////////////////////////////
61 62
    std::cout << ">> This executable has to be used to test Spherical Rotation algorithm.\n";
    std::cout << ">> You can pass -sequential or -task (thread by default).\n";
63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87
    //////////////////////////////////////////////////////////////
    const int DevP = FParameters::getValue(argc,argv,"-p", 8);
    const int NbLevels = FParameters::getValue(argc,argv,"-h", 5);
    const int SizeSubLevels = FParameters::getValue(argc,argv,"-sh", 3);
    FTic counter;

    const char* const filename = FParameters::getStr(argc,argv,"-f", "../Data/test20k.fma");
    std::cout << "Opening : " << filename << "\n";

    FFmaScanfLoader<ParticleClass> loader(filename);
    if(!loader.isOpen()){
        std::cout << "Loader Error, " << filename << " is missing\n";
        return 1;
    }

    // -----------------------------------------------------
    CellClass::Init(DevP);
    OctreeClass tree(NbLevels, SizeSubLevels,loader.getBoxWidth(),loader.getCenterOfBox());

    // -----------------------------------------------------

    std::cout << "Creating & Inserting " << loader.getNumberOfParticles() << " particles ..." << std::endl;
    std::cout << "\tHeight : " << NbLevels << " \t sub-height : " << SizeSubLevels << std::endl;
    counter.tic();

88
    loader.fillTree(tree);
89 90 91 92 93 94

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

    // -----------------------------------------------------

95
    std::cout << "Create kernel ..." << std::endl;
96 97
    counter.tic();

98
    KernelClass kernels(DevP, NbLevels, loader.getBoxWidth(), loader.getCenterOfBox());
99

100 101 102 103 104 105 106
    counter.tac();
    std::cout << "Done  " << " in " << counter.elapsed() << "s)." << std::endl;

    // -----------------------------------------------------

    std::cout << "Working on particles ..." << std::endl;

107 108
    if( FParameters::findParameter(argc,argv,"-sequential") != FParameters::NotFound){
        FmmClass algo(&tree,&kernels);
109
        counter.tic();
110 111 112 113
        algo.execute();
    }
    else if( FParameters::findParameter(argc,argv,"-task") != FParameters::NotFound){
        FmmClassTask algo(&tree,&kernels);
114
        counter.tic();
115 116 117 118
        algo.execute();
    }
    else {
        FmmClassThread algo(&tree,&kernels);
119
        counter.tic();
120 121
        algo.execute();
    }
122 123 124 125 126 127

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

    { // get sum forces&potential
        FReal potential = 0;
128
        FPoint forces;
129
        OctreeClass::Iterator octreeIterator(&tree);
130 131
        octreeIterator.gotoBottomLeft();
        do{
132
            ContainerClass::ConstBasicIterator iter(*octreeIterator.getCurrentListTargets());
133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150
            while( iter.hasNotFinished() ){
                potential += iter.data().getPotential() * iter.data().getPhysicalValue();
                forces += iter.data().getForces();

                iter.gotoNext();
            }
        } while(octreeIterator.moveRight());

        std::cout << "Foces Sum  x = " << forces.getX() << " y = " << forces.getY() << " z = " << forces.getZ() << std::endl;
        std::cout << "Potential = " << potential << std::endl;
    }

    return 0;
}