utestSphericalWithPrevious.cpp 6.64 KB
Newer Older
1
// See LICENCE file at project root
2 3


BRAMAS Berenger's avatar
BRAMAS Berenger committed
4
#include "Utils/FGlobal.hpp"
5

BRAMAS Berenger's avatar
BRAMAS Berenger committed
6 7
#include "Containers/FOctree.hpp"
#include "Containers/FVector.hpp"
8

BRAMAS Berenger's avatar
BRAMAS Berenger committed
9 10 11
#include "Kernels/Spherical/FSphericalCell.hpp"
#include "Kernels/Spherical/FSphericalKernel.hpp"
#include "Components/FSimpleLeaf.hpp"
12

BRAMAS Berenger's avatar
BRAMAS Berenger committed
13
#include "Files/FFmaGenericLoader.hpp"
14

BRAMAS Berenger's avatar
BRAMAS Berenger committed
15
#include "Files/FTreeIO.hpp"
16

BRAMAS Berenger's avatar
BRAMAS Berenger committed
17
#include "Core/FFmmAlgorithm.hpp"
18 19

#include "FUTester.hpp"
BRAMAS Berenger's avatar
BRAMAS Berenger committed
20
#include "Kernels/P2P/FP2PParticleContainerIndexed.hpp"
21

berenger-bramas's avatar
berenger-bramas committed
22
/**
23 24
 * This test compare a previous FMM result with a previous simulation result.
 */
25

26
typedef double FReal;
27
typedef FSphericalCell<FReal>           CellClass;
28
typedef FP2PParticleContainerIndexed<FReal>  ContainerClass;
29

30
typedef FSphericalKernel< FReal, CellClass, ContainerClass >          KernelClass;
31

32 33
typedef FSimpleLeaf<FReal, ContainerClass >                     LeafClass;
typedef FOctree<FReal, CellClass, ContainerClass , LeafClass >  OctreeClass;
34

35
typedef FFmmAlgorithm<OctreeClass, CellClass, ContainerClass, KernelClass, LeafClass > FmmClass;
36

37 38
/** To check if a value is correct */
bool IsSimilar(const FReal good, const FReal other){
39 40
	const FReal Epsilon = FReal(0.0001);
	return (FMath::Abs(good-other)/FMath::Abs(good)) < Epsilon;
41 42
}

berenger-bramas's avatar
berenger-bramas committed
43
/** The test class */
44
class TestSphericalWithPrevious : public FUTester<TestSphericalWithPrevious> {
45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65
	/** the test */
	void TestTree(){
		if(sizeof(FReal) == sizeof(float) ) {
			std::cerr << "No input data available for Float "<< std::endl;
			uassert(false);
		}
		//
		//  Load a Tree
		const std::string DataFile = (sizeof(FReal) == sizeof(float))?
				SCALFMMDataPath+"UTest/SphericalPrevious.data.single":
				SCALFMMDataPath+"UTest/SphericalPrevious.data.double";
		//
		// Load particles
		//

		const std::string parFile( (sizeof(FReal) == sizeof(float))?
				"UTest/DirectFloat.bfma":
				"UTest/DirectDouble.bfma");
		//
		std::string filename(SCALFMMDataPath+parFile);
		//
66
		FFmaGenericLoader<FReal> loader(filename);
67 68 69 70 71 72 73 74 75 76 77 78 79
		if(!loader.isOpen()){
			Print("Cannot open particles file.");
			uassert(false);
			return;
		}
		FSize nbParticles = loader.getNumberOfParticles() ;
		//
		const int NbLevels      = 5;
		const int SizeSubLevels = 3;
		const int DevP = 9;
		//
		// Create octree
		//
80
		FSphericalCell<FReal>::Init(DevP);
81 82 83
		//
		OctreeClass testTree(NbLevels, SizeSubLevels, loader.getBoxWidth(), loader.getCenterOfBox());
		//
84
		for(FSize idxPart = 0 ; idxPart < nbParticles ; ++idxPart){
85
            FPoint<FReal> position;
86 87 88 89 90 91 92 93 94 95 96 97 98 99
			FReal physicalValue = 0.0;
			loader.fillParticle(&position,&physicalValue);
			// put in tree
			testTree.insert(position, idxPart, physicalValue);
		}
		//
		//  Run simulation 1
		KernelClass kernels(DevP, NbLevels, loader.getBoxWidth(), loader.getCenterOfBox());
		FmmClass algo(&testTree,&kernels);
		Print("Run simulation 1 ...");

		algo.execute();

		// If needed save the result
100
        // FTreeIO<FReal>::Save<OctreeClass, CellClass, LeafClass, ContainerClass >(DataFile.c_str(), testTree);
101 102 103

		// Load previous result
		OctreeClass goodTree(NbLevels, SizeSubLevels, loader.getBoxWidth(), loader.getCenterOfBox());
104
        FTreeIO<FReal>::Load<OctreeClass, CellClass, LeafClass, ContainerClass >(DataFile.c_str(), goodTree);
105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128

		// Compare the two simulations
		Print("Check the particles...");
		{ // Check that each particle has been summed with all other
			OctreeClass::Iterator testOctreeIterator(&testTree);
			OctreeClass::Iterator goodOctreeIterator(&goodTree);

			testOctreeIterator.gotoBottomLeft();
			goodOctreeIterator.gotoBottomLeft();

			do{
				if(testOctreeIterator.getCurrentGlobalIndex() != goodOctreeIterator.getCurrentGlobalIndex()){
					uassert(false);
					break;
				}

				if(testOctreeIterator.getCurrentListSrc()->getNbParticles() != goodOctreeIterator.getCurrentListSrc()->getNbParticles()){
					uassert(false);
					break;
				}

				const ContainerClass* testLeaf = testOctreeIterator.getCurrentListSrc();
				const ContainerClass* goodLeaf = goodOctreeIterator.getCurrentListSrc();

129
				for(FSize idxPart = 0 ; idxPart < testLeaf->getNbParticles() ; ++idxPart ){
130 131 132 133
                    uassert( IsSimilar(goodLeaf->getPotentials()[idxPart], testLeaf->getPotentials()[idxPart]) );
                    uassert( IsSimilar(goodLeaf->getForcesX()[idxPart], testLeaf->getForcesX()[idxPart]) );
                    uassert( IsSimilar(goodLeaf->getForcesY()[idxPart], testLeaf->getForcesY()[idxPart]) );
                    uassert( IsSimilar(goodLeaf->getForcesZ()[idxPart], testLeaf->getForcesZ()[idxPart]) );
134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163
				}

				if(!testOctreeIterator.moveRight()){
					if(goodOctreeIterator.moveRight()){
						uassert(false);
					}
					break;
				}
				if(!goodOctreeIterator.moveRight()){
					uassert(false);
					break;
				}

			} while(true);
		}
		Print("Check the leaves...");
		{ // Ceck if there is number of NbPart summed at level 1
			OctreeClass::Iterator testOctreeIterator(&testTree);
			OctreeClass::Iterator goodOctreeIterator(&goodTree);

			testOctreeIterator.gotoBottomLeft();
			goodOctreeIterator.gotoBottomLeft();

			for(int idxLevel = NbLevels - 1 ; idxLevel > 1 ; --idxLevel ){
				do{
					if(testOctreeIterator.getCurrentGlobalIndex() != goodOctreeIterator.getCurrentGlobalIndex()){
						uassert(false);
						break;
					}

164 165 166 167 168
					for(int idxLocal = 0 ; idxLocal < CellClass::local_expansion_t::getSize() ; ++idxLocal){
						IsSimilar(testOctreeIterator.getCurrentCell()->getLocalExpansionData().get()[idxLocal].getReal(),
								goodOctreeIterator.getCurrentCell()->getLocalExpansionData().get()[idxLocal].getReal());
						IsSimilar(testOctreeIterator.getCurrentCell()->getLocalExpansionData().get()[idxLocal].getImag(),
								goodOctreeIterator.getCurrentCell()->getLocalExpansionData().get()[idxLocal].getImag());
169 170
					}

171 172 173 174 175
					for(int idxPole = 0 ; idxPole < CellClass::multipole_t::getSize() ; ++idxPole){
						IsSimilar(testOctreeIterator.getCurrentCell()->getMultipoleData().get()[idxPole].getReal(),
								goodOctreeIterator.getCurrentCell()->getMultipoleData().get()[idxPole].getReal());
						IsSimilar(testOctreeIterator.getCurrentCell()->getMultipoleData().get()[idxPole].getImag(),
								goodOctreeIterator.getCurrentCell()->getMultipoleData().get()[idxPole].getImag());
176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205
					}

					if(!testOctreeIterator.moveRight()){
						if(goodOctreeIterator.moveRight()){
							uassert(false);
						}
						break;
					}
					if(!goodOctreeIterator.moveRight()){
						uassert(false);
						break;
					}

				} while(true);

				testOctreeIterator.moveUp();
				testOctreeIterator.gotoLeft();

				goodOctreeIterator.moveUp();
				goodOctreeIterator.gotoLeft();
			}
		}
		Print("Over...");
	}


	// set test
	void SetTests(){
		AddTest(&TestSphericalWithPrevious::TestTree,"Test Simu and compare tree");
	}
206 207 208 209 210
};



// You must do this
211
TestClass(TestSphericalWithPrevious)
212 213


214