Commit 42c991d5 authored by PIACIBELLO Cyrille's avatar PIACIBELLO Cyrille
parents 8dad7822 e394a7e6
......@@ -57,6 +57,7 @@ OPTION( ScalFMM_USE_SSE "Set to ON to compile with SSE support"
OPTION( ScalFMM_USE_AVX "Set to ON to compile with AVX support" OFF )
OPTION( ScalFMM_USE_ASSERT "Set to ON to enable safe tests during execution" ON )
OPTION( ScalFMM_USE_MIC_NATIVE "Set to ON to compile in native mode for MIC" OFF )
OPTION( ScalFMM_BUILD_ONLY_LIB "Set to ON to compile ionly the lib (examples are not compiled) " OFF )
# Set scalfmm to default libraries
SET(SCALFMM_LIBRARIES "")
SET(ScaLFMM_CXX_FLAGS "")
......@@ -116,11 +117,11 @@ else()
IF( APPLE )
SET(SCALFMM_FLAGS_OPTI_RELEASE "-m64 -funroll-loops" CACHE STRING "Set your optimization flags for release mode.")
ELSE(APPLE)
# Not apple system - Check the compilater flags
# Not apple system - Check the compiler flags
if(CMAKE_CXX_COMPILER_ID STREQUAL "Intel")
SET(SCALFMM_FLAGS_OPTI_RELEASE "${FLAGS_64bits} -march=native -fp-model precise -fp-model source -fimf-precision=low -funroll-loops -ftree-vectorize" CACHE STRING "Set your optimization flags for release mode.")
else()
SET(SCALFMM_FLAGS_OPTI_RELEASE "${FLAGS_64bits} -ffast-math -flto -funroll-loops -ftree-vectorize" CACHE STRING "Set your optimization flags for release mode.")
SET(SCALFMM_FLAGS_OPTI_RELEASE "${FLAGS_64bits} -ffast-math -funroll-loops -ftree-vectorize" CACHE STRING "Set your optimization flags for release mode.")
endif()
ENDIF(APPLE)
#
......@@ -181,8 +182,7 @@ endif(ScalFMM_USE_FFT)
#ADD_DEFINITIONS(-Wall -Wshadow -Wpointer-arith -Wcast-qual -Wconversion -fpic )
#
# C++ 2011
#ADD_DEFINITIONS(-std=c++11)
MESSAGE( "CMAKE_SYSTEM: ${CMAKE_SYSTEM} CMAKE_CXX_COMPILER_ID ${CMAKE_CXX_COMPILER_ID} CMAKE_CXX_COMPILER ${NAME}")
MESSAGE( "CMAKE_SYSTEM: ${CMAKE_SYSTEM} CMAKE_CXX_COMPILER_ID ${CMAKE_CXX_COMPILER_ID} CMAKE_CXX_COMPILER ${NAME}")
IF(CMAKE_CXX_COMPILER_ID STREQUAL "Intel")
execute_process(COMMAND ${CMAKE_CXX_COMPILER} -dumpversion OUTPUT_VARIABLE INTEL_VERSION)
MESSAGE( "Intel: ${INTEL_VERSION}")
......@@ -301,9 +301,10 @@ add_subdirectory(Src)
# Link with scalfmm lib
set(scalfmm_lib scalfmm)
if( NOT ScalFMM_BUILD_ONLY_LIB )
# Build - Examples and drivers
add_subdirectory(Examples)
endif()
# Build - Tests
MESSAGE( STATUS "ScalFMM_BUILD_TESTS = ${ScalFMM_BUILD_TESTS}" )
......
Copyright (c) 2011-2014 Inria, All rights reserved.
===========================================================================
This file contains the main features as well as overviews of specific
bug fixes (and other actions) for each version of ScalFMM since
version 1.1
1.2 ()
-----
- New FMA format to read/write particles
- Add examples repository
- Add distributions generator
- Fix bugs in periodic model
- Fix bug and performance improvements in MPI thread algorithms
- Now we use MPI THREAD SERIALIZED rather than MPI_THREAD_MULTIPLE in MPI+OpenMP algorithm.
- Add Uniform interpolation based on Lagrange polynomials
- Improve tests -> now we use CTest and CDash (Thanks F Rue to Inria)
- Add support for Intel icl++ compiler (Intel 2015 compiler) on MacOS X and Clang on linux system
- Add SSE and AVX support
- Add SSE and AVX code for 1/r kernel
- CMake improvements
\ No newline at end of file
// ===================================================================================
// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Bérenger Bramas, Matthias Messner
// Copyright ScalFmm 2011 INRIA
// olivier.coulaud@inria.fr, berenger.bramas@inria.fr
// This software is a computer program whose purpose is to compute the FMM.
//
......@@ -18,8 +18,8 @@
#include "FAbstractSphericalKernel.hpp"
#include "../../Utils/FMemUtils.hpp"
#include "../../Utils/FBlas.hpp"
#include "Utils/FMemUtils.hpp"
#include "Utils/FBlas.hpp"
/**
* @author Berenger Bramas (berenger.bramas@inria.fr)
......@@ -123,7 +123,7 @@ public:
FF_MATRIX_ROW_DIM(Parent::harmonic.getExpSize()), FF_MATRIX_COLUMN_DIM(Parent::harmonic.getNExpSize()),
FF_MATRIX_SIZE(FF_MATRIX_ROW_DIM * FF_MATRIX_COLUMN_DIM),
temporaryMultiSource(new FComplexe[FF_MATRIX_COLUMN_DIM]),
preM2LTransitions(0){
preM2LTransitions(nullptr){
allocAndInit();
}
......
// ===================================================================================
// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Bérenger Bramas, Matthias Messner
// Copyright ScalFmm 2011 INRIA,
// olivier.coulaud@inria.fr, berenger.bramas@inria.fr
// This software is a computer program whose purpose is to compute the FMM.
//
......@@ -17,7 +17,7 @@
#define FSPHERICALROTATIONKERNEL_HPP
#include "FAbstractSphericalKernel.hpp"
#include "../../Utils/FMemUtils.hpp"
#include "Utils/FMemUtils.hpp"
/**
* @author Berenger Bramas (berenger.bramas@inria.fr)
......@@ -393,7 +393,7 @@ public:
FSphericalRotationKernel(const int inDevP, const int inTreeHeight, const FReal inBoxWidth, const FPoint& inBoxCenter)
: Parent(inDevP, inTreeHeight, inBoxWidth, inBoxCenter),
devM2lP(int(((inDevP*2)+1) * ((inDevP*2)+2) * 0.5)),
preM2LTransitions(0),
preM2LTransitions(nullptr),
rotation_Info(inDevP) {
allocAndInit();
}
......
......@@ -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,16 +16,14 @@
// ==== CMAKE =====
// @FUSE_BLAS
// ================
// ==============
#include "ScalFmmConfig.h"
#include "../Src/Utils/FGlobal.hpp"
#include "../Src/Containers/FOctree.hpp"
#include "../Src/Containers/FVector.hpp"
#include "../Src/Files/FFmaBinLoader.hpp"
#include "../Src/Files/FTreeIO.hpp"
#include "Files/FFmaGenericLoader.hpp"
#include "../Src/Core/FFmmAlgorithmThread.hpp"
#include "../Src/Core/FFmmAlgorithm.hpp"
......@@ -42,7 +40,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.
*/
......@@ -56,14 +54,24 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
///////////////////////////////////////////////////////////
template <class CellClass, class ContainerClass, class KernelClass, class MatrixKernelClass,
class LeafClass, class OctreeClass, class FmmClass, const int NbRhs>
void RunTest(const FReal epsilon) {
class LeafClass, class OctreeClass, class FmmClass, const int NVals>
void RunTest() {
// 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,119 +80,140 @@ 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 ,0.0,0.0,0.0);
}
// //
// // 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
Print("Fmm...");
KernelClass kernels(NbLevels, loader.getBoxWidth(), loader.getCenterOfBox(), epsilon);
KernelClass kernels(NbLevels, loader.getBoxWidth(), loader.getCenterOfBox());
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){
particles[idx].potential *= NVals ;
particles[idx].forces[0] *= NVals ;
particles[idx].forces[1] *= NVals ;
particles[idx].forces[2] *= NVals ;
energyD += particles[idx].potential*particles[idx].physicalValue ;
}
//
// Compare
Print("Compute Diff...");
FMath::FAccurater potentialDiff;
FMath::FAccurater fx, fy, fz;
{ // Check that each particle has been summed with all other
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();
const FVector<int>& indexes = leaf->getTargets()->getIndexes();
for(int idxPart = 0 ; idxPart < nbParticlesInLeaf ; ++idxPart){
const int indexPartOrig = indexes[idxPart];
potentialDiff.add(particles[indexPartOrig].potential,potentials[idxPart]);
fx.add(particles[indexPartOrig].forces[0],forcesX[idxPart]);
fy.add(particles[indexPartOrig].forces[1],forcesY[idxPart]);
fz.add(particles[indexPartOrig].forces[2],forcesZ[idxPart]);
}
});
}
delete[] particles;
{ // Check that each particle has been summed with all other
tree.forEachLeaf([&](LeafClass* leaf){
const FReal*const potentials = leaf->getTargets()->getPotentials();
const FReal*const physicalValues = leaf->getTargets()->getPhysicalValues();
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();
const FVector<int>& indexes = leaf->getTargets()->getIndexes();
for(int idxPart = 0 ; idxPart < nbParticlesInLeaf ; ++idxPart){
const int indexPartOrig = indexes[idxPart];
potentialDiff.add(particles[indexPartOrig].potential,potentials[idxPart]);
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];
}
});
}
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 MaximumDiffPotential = FReal(9e-3);
const FReal MaximumDiffForces = FReal(9e-2);
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;
FMath::FAccurater multiPoleDiff;
tree.forEachCell([&](CellClass* cell){
for( int idxRhs = 1 ; idxRhs < NbRhs ; ++idxRhs){
for( int idxRhs = 1 ; idxRhs < NVals ; ++idxRhs){
localDiff.add(cell->getLocal(0), cell->getLocal(idxRhs), cell->getVectorSize());
multiPoleDiff.add(cell->getMultipole(0), cell->getMultipole(idxRhs), cell->getVectorSize());
}
......@@ -222,34 +251,32 @@ class TestChebyshevDirect : public FUTester<TestChebyshevDirect> {
/** TestChebKernel */
void TestChebKernel(){
const int NbRhs = 4;
const unsigned int ORDER = 5;
const FReal epsilon = FReal(1e-5);
const int NVals = 4;
const unsigned int ORDER = 6 ;
typedef FP2PParticleContainerIndexed<> ContainerClass;
typedef FSimpleLeaf<ContainerClass> LeafClass;
typedef FInterpMatrixKernelR MatrixKernelClass;
typedef FChebCell<ORDER, 1, 1, NbRhs> CellClass;
typedef FChebCell<ORDER, 1, 1, NVals> CellClass;
typedef FOctree<CellClass,ContainerClass,LeafClass> OctreeClass;
typedef FChebKernel<CellClass,ContainerClass,MatrixKernelClass,ORDER, NbRhs> KernelClass;
typedef FChebKernel<CellClass,ContainerClass,MatrixKernelClass,ORDER, NVals> KernelClass;
typedef FFmmAlgorithm<OctreeClass,CellClass,ContainerClass,KernelClass,LeafClass> FmmClass;
// run test
RunTest<CellClass,ContainerClass,KernelClass,MatrixKernelClass,LeafClass,OctreeClass,FmmClass, NbRhs>(epsilon);
RunTest<CellClass,ContainerClass,KernelClass,MatrixKernelClass,LeafClass,OctreeClass,FmmClass, NVals>();
}
/** TestChebSymKernel */
void TestChebSymKernel(){
const int NbRhs = 4;
const unsigned int ORDER = 5;
const FReal epsilon = FReal(1e-5);
const int NVals = 4;
const unsigned int ORDER = 6;
typedef FP2PParticleContainerIndexed<> ContainerClass;
typedef FSimpleLeaf<ContainerClass> LeafClass;
typedef FInterpMatrixKernelR MatrixKernelClass;
typedef FChebCell<ORDER, 1, 1, NbRhs> CellClass;
typedef FChebCell<ORDER, 1, 1, NVals> CellClass;
typedef FOctree<CellClass,ContainerClass,LeafClass> OctreeClass;
typedef FChebSymKernel<CellClass,ContainerClass,MatrixKernelClass,ORDER, NbRhs> KernelClass;
typedef FChebSymKernel<CellClass,ContainerClass,MatrixKernelClass,ORDER, NVals> KernelClass;
typedef FFmmAlgorithm<OctreeClass,CellClass,ContainerClass,KernelClass,LeafClass> FmmClass;
// run test
RunTest<CellClass,ContainerClass,KernelClass,MatrixKernelClass,LeafClass,OctreeClass,FmmClass, NbRhs>(epsilon);
RunTest<CellClass,ContainerClass,KernelClass,MatrixKernelClass,LeafClass,OctreeClass,FmmClass, NVals>();
}
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment