testCompareIOTree.cpp 7.64 KB
Newer Older
1
// ===================================================================================
2
// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Berenger Bramas, Matthias Messner
3 4 5 6 7 8 9 10 11 12 13 14
// 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
#include "../../Src/Files/FTreeIO.hpp"
29

30 31
#include "../../Src/Kernels/Spherical/FSphericalCell.hpp"
#include "../../Src/Components/FSimpleLeaf.hpp"
32

33 34
#include "../../Src/Components/FSimpleLeaf.hpp"
#include "../../Src/Kernels/P2P/FP2PParticleContainer.hpp"
35

36 37
#include "../../Src/Utils/FParameterNames.hpp"

38 39
// Simply create particles and try the kernels
int main(int argc, char ** argv){
40 41 42 43 44 45 46 47

    FHelpDescribeAndExit(argc, argv,
                         "Load octrees that have been saved and compare everything from leaves to cells.\n"
                         "Using the Spherical Harmonics old kernel.",
                         FParameterDefinitions::SHDevelopment, FParameterDefinitions::InputFileOne,
                         FParameterDefinitions::InputFileTwow
                         );

48
    typedef double FReal;
49
    typedef FSphericalCell                 CellClass;
50
    typedef FP2PParticleContainer<FReal>         ContainerClass;
51

52 53
    typedef FSimpleLeaf<FReal, ContainerClass >                     LeafClass;
    typedef FOctree<FReal, CellClass, ContainerClass , LeafClass >  OctreeClass;
54 55 56
    ///////////////////////What we do/////////////////////////////
    std::cout << ">> This executable has to be used to compare two trees.\n";
    //////////////////////////////////////////////////////////////
57
    const int DevP = FParameters::getValue(argc,argv,FParameterDefinitions::SHDevelopment.options, 8);
58 59 60

    // -----------------------------------------------------
    CellClass::Init(DevP, true);
61 62
    OctreeClass tree1(5, 3, 0, FPoint<FReal>());
    OctreeClass tree2(5, 3, 0, FPoint<FReal>());
63 64

    // -----------------------------------------------------
65 66
    const char* const filename1 = FParameters::getStr(argc,argv,FParameterDefinitions::InputFileOne.options, "tree.data");
    const char* const filename2 = FParameters::getStr(argc,argv,FParameterDefinitions::InputFileOne.options, "dtree.data");
67 68
    std::cout << "Compare tree " << filename1 << " and " << filename2 << std::endl;

69 70
    FTreeIO::Load<OctreeClass, CellClass, LeafClass, ContainerClass >(filename1, tree1);
    FTreeIO::Load<OctreeClass, CellClass, LeafClass, ContainerClass >(filename2, tree2);
71 72 73 74

    // -----------------------------------------------------
    std::cout << "Check Result\n";
    { // Check that each particle has been summed with all other
75
        OctreeClass::Iterator octreeIterator1(&tree1);
76 77
        octreeIterator1.gotoBottomLeft();

78
        OctreeClass::Iterator octreeIterator2(&tree2);
79 80 81 82 83 84 85 86 87 88
        octreeIterator2.gotoBottomLeft();

        int nbLeaves = 0;

        do{
            if( octreeIterator1.getCurrentGlobalIndex() != octreeIterator2.getCurrentGlobalIndex()){
                std::cout << "Index is different\n";
                break;
            }

89
            if( octreeIterator1.getCurrentListSrc()->getNbParticles() != octreeIterator2.getCurrentListSrc()->getNbParticles()){
90
                std::cout << "Number of particles different on leaf " << octreeIterator1.getCurrentGlobalIndex() <<
91 92
                             " tree1 " << octreeIterator1.getCurrentListSrc()->getNbParticles() <<
                             " tree2 " << octreeIterator2.getCurrentListSrc()->getNbParticles() << std::endl;
93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114
            }

            nbLeaves += 1;

            if( octreeIterator1.moveRight() ){
                if( !octreeIterator2.moveRight() ){
                    std::cout << "Not the same number of leaf, tree2 end before tree1\n";
                    break;
                }
            }
            else {
                if( octreeIterator2.moveRight() ){
                    std::cout << "Not the same number of leaf, tree1 end before tree2\n";
                }
                break;
            }

        } while(true);

        std::cout << "There are " << nbLeaves << " leaves ...\n";
    }
    { // Ceck if there is number of NbPart summed at level 1
115
        OctreeClass::Iterator octreeIterator1(&tree1);
116 117
        octreeIterator1.gotoBottomLeft();

118
        OctreeClass::Iterator octreeIterator2(&tree2);
119 120 121 122
        octreeIterator2.gotoBottomLeft();

        for(int idxLevel = tree1.getHeight() - 1 ; idxLevel > 1 ; --idxLevel ){
            int nbCells = 0;
123 124 125 126 127
            do{
                if( octreeIterator1.getCurrentGlobalIndex() != octreeIterator2.getCurrentGlobalIndex()){
                    std::cout << "Index is different\n";
                    break;
                }
128

129 130
                const CellClass*const cell1 = octreeIterator1.getCurrentCell();
                const CellClass*const cell2 = octreeIterator2.getCurrentCell();
131

132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149
                FReal cumul = 0;
                for(int idx = 0; idx < FSphericalCell::GetPoleSize(); ++idx){
                    cumul += FMath::Abs( cell1->getMultipole()[idx].getImag() - cell2->getMultipole()[idx].getImag() );
                    cumul += FMath::Abs( cell1->getMultipole()[idx].getReal() - cell2->getMultipole()[idx].getReal() );
                }
                if( cumul > 0.00001 || FMath::IsNan(cumul)){
                    std::cout << "Pole Data are different. Cumul " << cumul << " at level " << idxLevel
                              << " index is " << octreeIterator1.getCurrentGlobalIndex() << std::endl;
                }
                cumul = 0;
                for(int idx = 0; idx < FSphericalCell::GetLocalSize(); ++idx){
                    cumul += FMath::Abs( cell1->getLocal()[idx].getImag() - cell2->getLocal()[idx].getImag() );
                    cumul += FMath::Abs( cell1->getLocal()[idx].getReal() - cell2->getLocal()[idx].getReal() );
                }
                if( cumul > 0.00001 || FMath::IsNan(cumul)){
                    std::cout << "Local Data are different. Cumul " << cumul << " at level " << idxLevel
                              << " index is " << octreeIterator1.getCurrentGlobalIndex() << std::endl;
                }
150

151 152 153 154
                nbCells += 1;
                if( octreeIterator1.moveRight() ){
                    if( !octreeIterator2.moveRight() ){
                        std::cout << "Not the same number of leaf, tree2 end before tree1\n";
155 156
                        break;
                    }
157 158 159 160 161 162 163 164 165
                }
                else {
                    if( octreeIterator2.moveRight() ){
                        std::cout << "Not the same number of leaf, tree1 end before tree2\n";
                    }
                    break;
                }
            } while(true);

166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183
            octreeIterator1.moveUp();
            octreeIterator1.gotoLeft();

            octreeIterator2.moveUp();
            octreeIterator2.gotoLeft();

            std::cout << "There are " << nbCells << " cells at level " << idxLevel << " ...\n";
        }
    }

    std::cout << "Done\n";

    return 0;
}