Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 4e4825c3 authored by Olivier COULAUD's avatar Olivier COULAUD
Browse files

utestChebyshevMultiRhs and utestChebyshev have the same format

parent 70c4c178
No related branches found
No related tags found
No related merge requests found
......@@ -63,7 +63,7 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
exit(EXIT_FAILURE);
}
const std::string parFile( (sizeof(FReal) == sizeof(float))?
"Test/DirectFloatbfma":
"Test/DirectFloat.bfma":
"UTest/DirectDouble.bfma");
//
std::string filename(SCALFMMDataPath+parFile);
......@@ -75,16 +75,16 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
const int NbLevels = 4;
const int SizeSubLevels = 2;
// Create octree
FSize nbParticles = loader.getNumberOfParticles() ;
FmaR8W8Particle* const particles = new FmaR8W8Particle[nbParticles];
loader.fillParticle(particles,nbParticles);
//
// Create octree
OctreeClass tree(NbLevels, SizeSubLevels, loader.getBoxWidth(), loader.getCenterOfBox());
// Insert particle in the tree
//
OctreeClass tree(NbLevels, SizeSubLevels, loader.getBoxWidth(), loader.getCenterOfBox());
for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
tree.insert(particles[idxPart].position , idxPart, particles[idxPart].physicalValue );
}
......
......@@ -16,7 +16,7 @@
// ==== CMAKE =====
// @FUSE_BLAS
// ================
// ==============
#include "ScalFmmConfig.h"
#include "../Src/Utils/FGlobal.hpp"
......@@ -24,7 +24,7 @@
#include "../Src/Containers/FOctree.hpp"
#include "../Src/Containers/FVector.hpp"
#include "../Src/Files/FFmaBinLoader.hpp"
#include "Files/FFmaGenericLoader.hpp"
#include "../Src/Files/FTreeIO.hpp"
#include "../Src/Core/FFmmAlgorithmThread.hpp"
......@@ -42,7 +42,7 @@
#include "../Src/Kernels/P2P/FP2PParticleContainerIndexed.hpp"
/*
In this test we compare the spherical fmm results and the direct results.
In this test we compare the spherical FMM results and the direct results.
*/
......@@ -60,10 +60,20 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
void RunTest(const FReal epsilon) {
// Warning in make test the exec dir it Build/UTests
// Load particles
const char* const filename = (sizeof(FReal) == sizeof(float))?
"../../Data/utestDirect.bin.fma.single":
"../../Data/utestDirect.bin.fma.double";
FFmaBinLoader loader(filename);
//
// Load particles
//
if(sizeof(FReal) == sizeof(float) ) {
std::cerr << "No input data available for Float "<< std::endl;
exit(EXIT_FAILURE);
}
const std::string parFile( (sizeof(FReal) == sizeof(float))?
"Test/DirectFloat.bfma":
"UTest/DirectDouble.bfma");
//
std::string filename(SCALFMMDataPath+parFile);
//
FFmaGenericLoader loader(filename);
if(!loader.isOpen()){
Print("Cannot open particles file.");
uassert(false);
......@@ -72,36 +82,39 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
Print("Number of particles:");
Print(loader.getNumberOfParticles());
const int NbLevels = 4;
const int NbLevels = 4;
const int SizeSubLevels = 2;
//const FReal epsilon = FReal(1e-5);
// Create octree
struct TestParticle{
FPoint position;
FReal forces[3];
FReal physicalValue;
FReal potential;
};
TestParticle* const particles = new TestParticle[loader.getNumberOfParticles()];
// Create octree
OctreeClass tree(NbLevels, SizeSubLevels, loader.getBoxWidth(), loader.getCenterOfBox());
for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
FPoint position;
FReal physicalValue;
loader.fillParticle(&position,&physicalValue);
// put in tree
tree.insert(position, idxPart, physicalValue);
// get copy
particles[idxPart].position = position;
particles[idxPart].physicalValue = physicalValue;
particles[idxPart].potential = 0.0;
particles[idxPart].forces[0] = 0.0;
particles[idxPart].forces[1] = 0.0;
particles[idxPart].forces[2] = 0.0;
}
//
FSize nbParticles = loader.getNumberOfParticles() ;
FmaR8W8Particle* const particles = new FmaR8W8Particle[nbParticles];
loader.fillParticle(particles,nbParticles);
//
// Create octree
OctreeClass tree(NbLevels, SizeSubLevels, loader.getBoxWidth(), loader.getCenterOfBox());
// Insert particle in the tree
//
for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
tree.insert(particles[idxPart].position , idxPart, particles[idxPart].physicalValue );
}
// //
// // Create octree
// for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){
// FPoint position;
// FReal physicalValue;
// loader.fillParticle(&position,&physicalValue);
// // put in tree
// tree.insert(position, idxPart, physicalValue);
// // get copy
// particles[idxPart].position = position;
// particles[idxPart].physicalValue = physicalValue;
// particles[idxPart].potential = 0.0;
// particles[idxPart].forces[0] = 0.0;
// particles[idxPart].forces[1] = 0.0;
// particles[idxPart].forces[2] = 0.0;
// }
// Run FMM
......@@ -109,23 +122,28 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
KernelClass kernels(NbLevels, loader.getBoxWidth(), loader.getCenterOfBox(), epsilon);
FmmClass algo(&tree,&kernels);
algo.execute();
//
FReal energy= 0.0 , energyD = 0.0 ;
// Run direct computation
Print("Direct...");
for( int idxRhs = 0 ; idxRhs < NbRhs ; ++idxRhs){
for(int idxTarget = 0 ; idxTarget < loader.getNumberOfParticles() ; ++idxTarget){
for(int idxOther = idxTarget + 1 ; idxOther < loader.getNumberOfParticles() ; ++idxOther){
FP2P::MutualParticles(particles[idxTarget].position.getX(), particles[idxTarget].position.getY(),
particles[idxTarget].position.getZ(),particles[idxTarget].physicalValue,
&particles[idxTarget].forces[0],&particles[idxTarget].forces[1],
&particles[idxTarget].forces[2],&particles[idxTarget].potential,
particles[idxOther].position.getX(), particles[idxOther].position.getY(),
particles[idxOther].position.getZ(),particles[idxOther].physicalValue,
&particles[idxOther].forces[0],&particles[idxOther].forces[1],
&particles[idxOther].forces[2],&particles[idxOther].potential);
}
}
}
for(int idx = 0 ; idx < loader.getNumberOfParticles() ; ++idx){
energyD += particles[idx].potential*particles[idx].physicalValue ;
}
// for( int idxRhs = 0 ; idxRhs < NbRhs ; ++idxRhs){
// for(int idxTarget = 0 ; idxTarget < loader.getNumberOfParticles() ; ++idxTarget){
// for(int idxOther = idxTarget + 1 ; idxOther < loader.getNumberOfParticles() ; ++idxOther){
// FP2P::MutualParticles(particles[idxTarget].position.getX(), particles[idxTarget].position.getY(),
// particles[idxTarget].position.getZ(),particles[idxTarget].physicalValue,
// &particles[idxTarget].forces[0],&particles[idxTarget].forces[1],
// &particles[idxTarget].forces[2],&particles[idxTarget].potential,
// particles[idxOther].position.getX(), particles[idxOther].position.getY(),
// particles[idxOther].position.getZ(),particles[idxOther].physicalValue,
// &particles[idxOther].forces[0],&particles[idxOther].forces[1],
// &particles[idxOther].forces[2],&particles[idxOther].potential);
// }
// }
// }
// Compare
Print("Compute Diff...");
......@@ -134,6 +152,7 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
{ // Check that each particle has been summed with all other
tree.forEachLeaf([&](LeafClass* leaf){
const FReal*const physicalValues = leaf->getTargets()->getPhysicalValues();
const FReal*const potentials = leaf->getTargets()->getPotentials();
const FReal*const forcesX = leaf->getTargets()->getForcesX();
const FReal*const forcesY = leaf->getTargets()->getForcesY();
......@@ -147,6 +166,8 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
fx.add(particles[indexPartOrig].forces[0],forcesX[idxPart]);
fy.add(particles[indexPartOrig].forces[1],forcesY[idxPart]);
fz.add(particles[indexPartOrig].forces[2],forcesZ[idxPart]);
energy += potentials[idxPart]*physicalValues[idxPart];
}
});
}
......@@ -154,31 +175,53 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
delete[] particles;
// Print for information
Print("Potential diff is = ");
Print(potentialDiff.getRelativeL2Norm());
Print(potentialDiff.getRelativeInfNorm());
Print("Fx diff is = ");
Print(fx.getRelativeL2Norm());
Print(fx.getRelativeInfNorm());
Print("Fy diff is = ");
Print(fy.getRelativeL2Norm());
Print(fy.getRelativeInfNorm());
Print("Fz diff is = ");
Print(fz.getRelativeL2Norm());
Print(fz.getRelativeInfNorm());
Print("Potential diff is = ");
printf(" Pot L2Norm %e\n",potentialDiff.getL2Norm());
printf(" Pot RL2Norm %e\n",potentialDiff.getRelativeL2Norm());
printf(" Pot RMSError %e\n",potentialDiff.getRMSError());
Print("Fx diff is = ");
printf(" Fx L2Norm %e\n",fx.getL2Norm());
printf(" Fx RL2Norm %e\n",fx.getRelativeL2Norm());
printf(" Fx RMSError %e\n",fx.getRMSError());
Print("Fy diff is = ");
printf(" Fy L2Norm %e\n",fy.getL2Norm());
printf(" Fy RL2Norm %e\n",fy.getRelativeL2Norm());
printf(" Fy RMSError %e\n",fy.getRMSError());
Print("Fz diff is = ");
printf(" Fz L2Norm %e\n",fz.getL2Norm());
printf(" Fz RL2Norm %e\n",fz.getRelativeL2Norm());
printf(" Fz RMSError %e\n",fz.getRMSError());
FReal L2error = (fx.getRelativeL2Norm()*fx.getRelativeL2Norm() + fy.getRelativeL2Norm()*fy.getRelativeL2Norm() + fz.getRelativeL2Norm() *fz.getRelativeL2Norm() );
printf(" Total L2 Force Error= %e\n",FMath::Sqrt(L2error)) ;
printf(" Energy Error = %.12e\n",FMath::Abs(energy-energyD));
printf(" Energy FMM = %.12e\n",FMath::Abs(energy));
printf(" Energy DIRECT = %.12e\n",FMath::Abs(energyD));
// Assert
const FReal MaximumDiffPotential = FReal(9e-5);
const FReal MaximumDiffForces = FReal(9e-3);
uassert(potentialDiff.getRelativeL2Norm() < MaximumDiffPotential);
uassert(potentialDiff.getRelativeInfNorm() < MaximumDiffPotential);
uassert(fx.getRelativeL2Norm() < MaximumDiffForces);
uassert(fx.getRelativeInfNorm() < MaximumDiffForces);
uassert(fy.getRelativeL2Norm() < MaximumDiffForces);
uassert(fy.getRelativeInfNorm() < MaximumDiffForces);
uassert(fz.getRelativeL2Norm() < MaximumDiffForces);
uassert(fz.getRelativeInfNorm() < MaximumDiffForces);
const FReal MaximumDiffForces = FReal(9e-3);
Print("Test1 - Error Relative L2 norm Potential ");
uassert(potentialDiff.getRelativeL2Norm() < MaximumDiffPotential); //1
Print("Test2 - Error RMS L2 norm Potential ");
uassert(potentialDiff.getRMSError() < MaximumDiffPotential); //2
Print("Test3 - Error Relative L2 norm FX ");
uassert(fx.getRelativeL2Norm() < MaximumDiffForces); //3
Print("Test4 - Error RMS L2 norm FX ");
uassert(fx.getRMSError() < MaximumDiffForces); //4
Print("Test5 - Error Relative L2 norm FY ");
uassert(fy.getRelativeL2Norm() < MaximumDiffForces); //5
Print("Test6 - Error RMS L2 norm FY ");
uassert(fy.getRMSError() < MaximumDiffForces); //6
Print("Test7 - Error Relative L2 norm FZ ");
uassert(fz.getRelativeL2Norm() < MaximumDiffForces); //8
Print("Test8 - Error RMS L2 norm FZ ");
uassert(fz.getRMSError() < MaximumDiffForces); //8
Print("Test9 - Error Relative L2 norm F ");
uassert(L2error < MaximumDiffForces); //9 Total Force
Print("Test10 - Relative error Energy ");
uassert(FMath::Abs(energy-energyD) /energyD< MaximumDiffPotential); //10 Total Energy
// Compute multipole local rhs diff
FMath::FAccurater localDiff;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment