testTsmAlgorithm.cpp 6.76 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

#include <iostream>

19
20
#include <cstdio>
#include <cstdlib>
21

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

25
#include "../../Src/Components/FTypedLeaf.hpp"
26

27
28
#include "../../Src/Containers/FOctree.hpp"
#include "../../Src/Containers/FVector.hpp"
29

30
#include "../../Src/Core/FFmmAlgorithmTsm.hpp"
31

32
33
#include "../../Src/Kernels/Spherical/FSphericalKernel.hpp"
#include "../../Src/Kernels/Spherical/FSphericalCell.hpp"
34

BRAMAS Berenger's avatar
BRAMAS Berenger committed
35
36
37
#include "../../Src/Kernels/Rotation/FRotationKernel.hpp"
#include "../../Src/Kernels/Rotation/FRotationCell.hpp"

38
#include "../../Src/Files/FFmaTsmLoader.hpp"
39

40
41
#include "../../Src/Kernels/P2P/FP2PParticleContainer.hpp"

42
43
/** This program show an example of use of
  * the fmm basic algo
44
  * it also check that eachh particles is little or longer
45
46
47
48
  * related that each other
  */


49
// Simply create particles and try the kernels
BRAMAS Berenger's avatar
BRAMAS Berenger committed
50
51
52
53
54
template <class CellClass, class ContainerClass, class LeafClass, class OctreeClass,
          class KernelClass, class FmmClass, typename... Args>
int testFunction(int argc, char ** argv, Args... kernelPreArgs){
    FTic counter;
    // Retrieve parameters
55
    const int NbLevels = FParameters::getValue(argc,argv,"-h", 5);
56
    const int SizeSubLevels = FParameters::getValue(argc,argv,"-sh", 3);
BRAMAS Berenger's avatar
BRAMAS Berenger committed
57
    // Get working file
58
59
    const char* const filename = FParameters::getStr(argc,argv,"-f", "../Data/test20k.tsm.fma");
    std::cout << "Opening : " << filename << "\n";
BRAMAS Berenger's avatar
BRAMAS Berenger committed
60
    // Create particles loader
61
    FFmaTsmLoader loader(filename);
62
    if(!loader.isOpen()){
63
64
65
66
67
        std::cout << "Loader Error, " << filename << " is missing\n";
        return 1;
    }

    // -----------------------------------------------------
BRAMAS Berenger's avatar
BRAMAS Berenger committed
68
    // Build the tree
berenger-bramas's avatar
berenger-bramas committed
69
    OctreeClass tree(NbLevels, SizeSubLevels,loader.getBoxWidth(),loader.getCenterOfBox());
70
71
72

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

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

77
78
79
    for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
        FPoint particlePosition;
        FReal physicalValue = 0.0;
80
81
82
        FParticleType particleType;
        loader.fillParticle(&particlePosition,&physicalValue,&particleType);
        tree.insert(particlePosition, particleType, physicalValue );
83
    }
84
85

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

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

90
    std::cout << "Create kernel ..." << std::endl;
91
92
    counter.tic();

BRAMAS Berenger's avatar
BRAMAS Berenger committed
93
    KernelClass kernels( kernelPreArgs... , NbLevels, loader.getBoxWidth(), loader.getCenterOfBox());
berenger-bramas's avatar
berenger-bramas committed
94

95
96
97
98
99
100
101
    counter.tac();
    std::cout << "Done  " << " in " << counter.elapsed() << "s)." << std::endl;

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

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

berenger-bramas's avatar
berenger-bramas committed
102
    FmmClass algo(&tree,&kernels);
103

104
105
    counter.tic();
    algo.execute();
106
    counter.tac();
107

108
    std::cout << "Done  " << "(@Algorithm = " << counter.elapsed() << "s)." << std::endl;
109
110
111

    { // get sum forces&potential
        FReal potential = 0;
112
113
114
115
116
117
118
119
120
121
122
123
124
125
        FReal fx = 0.0, fy = 0.0, fz = 0.0;

        tree.forEachLeaf([&](LeafClass* leaf){
            const FReal*const potentials = leaf->getTargets()->getPotentials();
            const FReal*const forcesX = leaf->getTargets()->getForcesX();
            const FReal*const forcesY = leaf->getTargets()->getForcesY();
            const FReal*const forcesZ = leaf->getTargets()->getForcesZ();
            const int nbParticlesInLeaf = leaf->getTargets()->getNbParticles();

            for(int idxPart = 0 ; idxPart < nbParticlesInLeaf ; ++idxPart){
                potential += potentials[idxPart];
                fx += forcesX[idxPart];
                fy += forcesY[idxPart];
                fz += forcesZ[idxPart];
126
            }
127
        });
128

129
        std::cout << "Foces Sum  x = " << fx << " y = " << fy << " z = " << fz << std::endl;
130
131
132
133
134
135
136
137
138
        std::cout << "Potential = " << potential << std::endl;
    }


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

    return 0;
}

BRAMAS Berenger's avatar
BRAMAS Berenger committed
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
// This is the real main!
int main(int argc, char ** argv){
    std::cout << "[PARAM] Use Parameters -spherical -rotation -chebyshev\n";

    if( FParameters::existParameter(argc,argv,"-spherical") ){
        std::cout << "[INFO] -spherical is used\n";
        // Create template
        typedef FTypedSphericalCell            CellClass;
        typedef FP2PParticleContainer         ContainerClass;

        typedef FTypedLeaf< ContainerClass >                      LeafClass;
        typedef FOctree< CellClass, ContainerClass , LeafClass >  OctreeClass;
        typedef FSphericalKernel< CellClass, ContainerClass >          KernelClass;

        typedef FFmmAlgorithmTsm<OctreeClass, CellClass, ContainerClass, KernelClass, LeafClass > FmmClass;
154

BRAMAS Berenger's avatar
BRAMAS Berenger committed
155
156
        const int DevP = FParameters::getValue(argc,argv,"-p", 8);
        CellClass::Init(DevP);
157

BRAMAS Berenger's avatar
BRAMAS Berenger committed
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
        // Call Main function
        testFunction< CellClass, ContainerClass, LeafClass, OctreeClass, KernelClass, FmmClass>(argc, argv, DevP);
    }

    if( FParameters::existParameter(argc,argv,"-rotation") ){
        std::cout << "[INFO] -rotation is used\n";
        // Create template
        static const int P = 9;
        typedef FTypedRotationCell<P>            CellClass;
        typedef FP2PParticleContainer         ContainerClass;

        typedef FTypedLeaf< ContainerClass >                      LeafClass;
        typedef FOctree< CellClass, ContainerClass , LeafClass >  OctreeClass;
        typedef FRotationKernel< CellClass, ContainerClass, P >          KernelClass;

        typedef FFmmAlgorithmTsm<OctreeClass, CellClass, ContainerClass, KernelClass, LeafClass > FmmClass;

        // Call Main function
        testFunction< CellClass, ContainerClass, LeafClass, OctreeClass, KernelClass, FmmClass>(argc, argv);
    }

    return 0;
}