diff --git a/Utils/src/ChebyshevInterpolationCmpAlgo.cpp b/Utils/src/ChebyshevInterpolationCmpAlgo.cpp
deleted file mode 100644
index f22105dcead1c650e1b4a7640c1acc65d29bff21..0000000000000000000000000000000000000000
--- a/Utils/src/ChebyshevInterpolationCmpAlgo.cpp
+++ /dev/null
@@ -1,33 +0,0 @@
-// ==== CMAKE =====
-// @FUSE_BLAS
-// ================
-// Keep in private GIT
-//
-//
-/** \brief Chebyshev FMM example
- *
- * \file
- * \authors B. Bramas, O. Coulaud
- *
- * This program runs the FMM Algorithm with the interpolation kernel based on
- * Chebyshev (grid points) interpolation (1/r kernel). It then compares the
- * results with a direct computation.
- */
-#include <string>
-#include "Kernels/Chebyshev/FChebCell.hpp"
-#include "Kernels/Chebyshev/FChebSymKernel.hpp"
-//
-template<typename FReal, int ORDER>
-using FInterpolationCell = FChebCell<FReal, ORDER>;
-
-template<typename FReal, typename GroupCellClass,
- typename GroupContainerClass,
- typename MatrixKernelClass, int ORDER>
-using FInterpolationKernel = FChebSymKernel<FReal,
- GroupCellClass,
- GroupContainerClass,
- MatrixKernelClass,
- ORDER> ;
-const std::string interpolationType("Chebyshev interpolation");
-
-#include "sharedMemoryInterpolatiomCmpAlgos.hpp"
diff --git a/Utils/src/UniformInterpolationCmpAlgo.cpp b/Utils/src/UniformInterpolationCmpAlgo.cpp
deleted file mode 100644
index 43cdf1bdeef8dbd1eeb75f72e3d616296f1a3cce..0000000000000000000000000000000000000000
--- a/Utils/src/UniformInterpolationCmpAlgo.cpp
+++ /dev/null
@@ -1,34 +0,0 @@
-// ==== CMAKE =====
-// @FUSE_BLAS
-// @FUSE_FFT
-// ================
-// Keep in private GIT
-//
-//
-/** \brief Uniform FMM example
- *
- * \file
- * \authors B. Bramas, O. Coulaud
- *
- * This program runs the FMM Algorithm with the interpolation kernel based on
- * uniform (grid points) interpolation (1/r kernel). It then compares the
- * results with a direct computation.
- */
-#include <string>
-#include "Kernels/Uniform/FUnifCell.hpp"
-#include "Kernels/Uniform/FUnifKernel.hpp"
-//
-template<typename FReal, int ORDER>
-using FInterpolationCell = FUnifCell<FReal, ORDER>;
-
-template<typename FReal, typename GroupCellClass,
- typename GroupContainerClass,
- typename MatrixKernelClass, int ORDER>
-using FInterpolationKernel = FUnifKernel<FReal,
- GroupCellClass,
- GroupContainerClass,
- MatrixKernelClass,
- ORDER> ;
-const std::string interpolationType("Uniform interpolation");
-
-#include "sharedMemoryInterpolatiomCmpAlgos.hpp"
diff --git a/Utils/src/sharedMemoryInterpolatiomCmpAlgos.hpp b/Utils/src/sharedMemoryInterpolatiomCmpAlgos.hpp
deleted file mode 100644
index 80d9ce8b7bb75bce851174b38d1e903c794e044e..0000000000000000000000000000000000000000
--- a/Utils/src/sharedMemoryInterpolatiomCmpAlgos.hpp
+++ /dev/null
@@ -1,387 +0,0 @@
-// -*-c++-*-
-// See LICENCE file at project root
-
-// ==== CMAKE =====
-// @FUSE_FFT
-// @FUSE_BLAS
-// ==== Git =====
-
-// ================
-
-/** \brief Uniform FMM example
- *
- * \file
- * \authors B. Bramas, O. Coulaud
- *
- * This program runs the FMM Algorithm with the interpolation kernel based on
- * either uniform (grid points) or Chebychev interpolation (1/r kernel).
- * It then compares the results with a direct computation.
- */
-
-#include <iostream>
-#include <iomanip>
-#include <memory>
-
-#include <cstdio>
-#include <cstdlib>
-#include <string>
-
-
-#include "ScalFmmConfig.h"
-#include "Utils/FGlobal.hpp"
-
-#include "Utils/FParameters.hpp"
-#include "Utils/FParameterNames.hpp"
-
-#include "Files/FFmaGenericLoader.hpp"
-// Leaves
-#include "Components/FSimpleLeaf.hpp"
-#include "Kernels/P2P/FP2PParticleContainerIndexed.hpp"
-// Octree
-#include "Containers/FOctree.hpp"
-//
-#include "Kernels/Interpolation/FInterpMatrixKernel.hpp"
-//
-#ifdef _OPENMP
-#include "Core/FFmmAlgorithmTask.hpp"
-#include "Core/FFmmAlgorithmNewTask.hpp"
-#include "Core/FFmmAlgorithmSectionTask.hpp"
-#include "Core/FFmmAlgorithmOmp4.hpp"
-#else
-#include "Core/FFmmAlgorithm.hpp"
-#endif
-//
-// Order of the Interpolation approximation
-static constexpr unsigned ORDER = 6 ;
-using FReal = double;
-// 1/r kernel
-//
-using MatrixKernelType = FInterpMatrixKernelR<FReal> ;
-
-//
-/// definition of the common tree structure
-using CellType = FInterpolationCell<FReal, ORDER>;
-//using CellUpType = typename CellType::multipole_t;
-//using CellDownType = typename CellType::local_expansion_t;
-//using CellSymbType = FSymbolicData;
-using ContainerType = FP2PParticleContainerIndexed<FReal>;
-using LeafType = FSimpleLeaf<FReal, ContainerType > ;
-using OctreeType = FOctree<FReal, CellType,ContainerType,LeafType>;
-using KernelType = FInterpolationKernel<FReal,CellType,ContainerType,MatrixKernelType,ORDER> ;
-
-
-#ifdef _OPENMP
-using TaskFmmAlgo = FFmmAlgorithmTask<OctreeType,CellType,ContainerType,KernelType,LeafType>;
-using TaskNewFmmAlgo = FFmmAlgorithmNewTask<OctreeType,CellType,ContainerType,KernelType,LeafType>;
-using SectionTaskFmmAlgo = FFmmAlgorithmSectionTask<OctreeType,CellType,ContainerType,KernelType,LeafType> ;
-using Omp4FmmAlgo = FFmmAlgorithmOmp4<OctreeType,CellType,ContainerType,KernelType,LeafType>;
-#else
-using FmmType = FFmmAlgorithm<OctreeType,CellType,ContainerType,KernelType,LeafType>;
-#endif
-
-
-
-// Simply create particles and try the kernels
-int main(int argc, char* argv[]) {
- const FParameterNames LocalOptionAlgo= { {"-algo"} , " Algorithm to run (task, newtask, sectiontask, omp4)\n"};
- const FParameterNames LocalOptionCmp = {
- {"-cmp"} , "Use to check the result with the exact solution given in the input file\n" };
- FHelpDescribeAndExit(
- argc, argv,
- "Driver for Lagrange Or Chebychev interpolation kernel (1/r kernel).",
- FParameterDefinitions::OctreeHeight,
- FParameterDefinitions::OctreeSubHeight,
- FParameterDefinitions::InputFile,
- FParameterDefinitions::OutputFile,
- FParameterDefinitions::NbThreads,
- LocalOptionAlgo,LocalOptionCmp
- );
-
-
- const std::string defaultFile(SCALFMMDataPath+"unitCubeXYZQ100.bfma" );
- const std::string filename = FParameters::getStr(argc,argv,FParameterDefinitions::InputFile.options, defaultFile.c_str());
- const int TreeHeight = FParameters::getValue(argc, argv, FParameterDefinitions::OctreeHeight.options, 5);
- const int SubTreeHeight = FParameters::getValue(argc, argv, FParameterDefinitions::OctreeSubHeight.options, 2);
-
-#ifdef _OPENMP
- const unsigned int NbThreads = FParameters::getValue(argc, argv, FParameterDefinitions::NbThreads.options, omp_get_max_threads());
- omp_set_num_threads(NbThreads);
- std::cout << "\n>> Using " << NbThreads << " threads.\n" << std::endl;
-#else
- const int NbThreads = FParameters::getValue(argc, argv, FParameterDefinitions::NbThreads.options, 1);
- std::cout << "\n>> Sequential version.\n" << std::endl;
-#endif
- //
- {
- std::string indent(" ");
- auto w = std::setw(18);
- std::cout << "Parameters" << std::endl << std::left
- << indent << w << "Octree Depth: " << TreeHeight << std::endl
- << indent << w << "SubOctree depth: " << SubTreeHeight << std::endl
- << indent << w << "Input file name: " << filename << std::endl
- << indent << w << "Thread number: " << NbThreads << std::endl
- << std::endl;
- }
- //
- // init timer
- FTic time;
-
- // open particle file
- ////////////////////////////////////////////////////////////////////
- //
- FFmaGenericLoader<FReal> loader(filename);
- if (loader.getNbRecordPerline() !=8 && FParameters::existParameter(argc, argv, "-cmp") ){
- std::cerr << "File should contain 8 data to read (x,y,z,q,p,fx,fy,fz)\n";
- std::exit(EXIT_FAILURE);
- }
- FSize nbParticles = loader.getNumberOfParticles() ;
- FmaRWParticle<FReal,8,8> * const particles = new FmaRWParticle<FReal, 8,8>[nbParticles];
- //
- //
- // open particle file
- ////////////////////////////////////////////////////////////////////
- //
- if(loader.getNbRecordPerline() == 8 ){
- loader.fillParticle(particles,nbParticles);
- }
- else {
- FmaRWParticle<FReal,4,4> * const part = new FmaRWParticle<FReal, 4,4>[nbParticles];
- loader.fillParticle(part,nbParticles);
- for(int idx = 0 ; idx < nbParticles ; ++idx){
- auto p = part[idx].getPosition() ;
- particles[idx].setPosition(p) ;
- particles[idx].setPhysicalValue(part[idx].getPhysicalValue() );
- }
- }
- time.tac();
- //
- // init oct-tree
- OctreeType tree(TreeHeight, SubTreeHeight, loader.getBoxWidth(), loader.getCenterOfBox());
- FReal energyD = 0.0 ;
- //
- // Read particles and insert them in octree
- { // -----------------------------------------------------
- std::cout << "Creating & Inserting " << loader.getNumberOfParticles()
- << " particles ..." << std::endl;
- std::cout << "\tHeight : " << TreeHeight << " \t sub-height : " << SubTreeHeight << std::endl;
- time.tic();
- //
- //
- //
- /////////////////////////////////////////////////////////////////////////////////////////////////
- // Compute direct energy
- /////////////////////////////////////////////////////////////////////////////////////////////////
-
- for(int idx = 0 ; idx < nbParticles ; ++idx){
- tree.insert(particles[idx].getPosition() , idx, particles[idx].getPhysicalValue() );
- energyD += particles[idx].getPotential()*particles[idx].getPhysicalValue() ;
- }
- time.tac();
- std::cout << "Done " << "(@Creating and Inserting Particles = "
- << time.elapsed() << " s) ." << std::endl;
- }
- ////////////////////////////////////////////////////////////////////
- //
- // Execute FMM Algorithm
- //
- ////////////////////////////////////////////////////////////////////
-
- { // -----------------------------------------------------
- std::cout << "\n" << interpolationType << " FMM (ORDER= "<< ORDER << ") ... " << std::endl;
-
- const MatrixKernelType MatrixKernel;
- time.tic();
- //
- std::unique_ptr<KernelType> kernels(new KernelType(TreeHeight, loader.getBoxWidth(),
- loader.getCenterOfBox(),&MatrixKernel));
- std::string algoStr = FParameters::getStr(argc,argv,"-algo", "task");
- //
- TaskFmmAlgo algo1(&tree, kernels.get() );
- TaskNewFmmAlgo algo2(&tree, kernels.get() );
- SectionTaskFmmAlgo algo3(&tree, kernels.get() );
- Omp4FmmAlgo algo4(&tree, kernels.get() );
-
- FAbstractAlgorithm * algo = nullptr;
- FAlgorithmTimers * timer = nullptr;
- if( "task" == algoStr) {
- algo = &algo1 ;
- timer = &algo1 ;
- } else if( "newtask" == algoStr) {
- algo = &algo2 ;
- timer = &algo2 ;
- } else if ( "sectiontask" == algoStr ) {
- algo = &algo3 ;
- timer = &algo3 ;
- } else if ( "omp4" == algoStr ) {
- algo = &algo4 ;
- timer = &algo4 ;
- } else {
- std::cout << "Unknown algorithm: " << algoStr << std::endl;
- std::exit(EXIT_FAILURE);
- }
- std::cout << "Algorithm to check: "<< algoStr <<std::endl;
- time.tic();
-
- //
- algo->execute(); // Here the call of the FMM algorithm
- //
- time.tac();
- std::cout << "Timers Far Field \n"
- << "P2M " << timer->getTime(FAlgorithmTimers::P2MTimer) << " seconds\n"
- << "M2M " << timer->getTime(FAlgorithmTimers::M2MTimer) << " seconds\n"
- << "M2L " << timer->getTime(FAlgorithmTimers::M2LTimer) << " seconds\n"
- << "L2L " << timer->getTime(FAlgorithmTimers::L2LTimer) << " seconds\n"
- << "L2P " << timer->getTime(FAlgorithmTimers::L2PTimer) << " seconds\n"
- << "P2P and L2P " << timer->getTime(FAlgorithmTimers::NearTimer) << " seconds\n"
- << std::endl;
-
-
- std::cout << "Done " << "(@Algorithm = " << time.elapsed() << " s) ." << std::endl;
- }
- // -----------------------------------------------------
- //
- // Some output
- //
- //
- { // -----------------------------------------------------
- FSize N1=0, N2= loader.getNumberOfParticles()/2, N3= loader.getNumberOfParticles() -1; ;
- FReal energy =0.0 ;
- //
- // Loop over all leaves
- //
- std::cout <<std::endl<<" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "<<std::endl;
- std::cout << std::scientific;
- std::cout.precision(10) ;
-
- tree.forEachLeaf([&](LeafType* leaf){
- const FReal*const posX = leaf->getTargets()->getPositions()[0];
- const FReal*const posY = leaf->getTargets()->getPositions()[1];
- const FReal*const posZ = leaf->getTargets()->getPositions()[2];
-
- 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 FSize nbParticlesInLeaf = leaf->getTargets()->getNbParticles();
- const FReal*const physicalValues = leaf->getTargets()->getPhysicalValues();
-
- const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
-
- for(FSize idxPart = 0 ; idxPart < nbParticlesInLeaf ; ++idxPart){
- const FSize indexPartOrig = indexes[idxPart];
- if ((indexPartOrig == N1) || (indexPartOrig == N2) || (indexPartOrig == N3) ) {
- std::cout << "Index "<< indexPartOrig <<" potential " << potentials[idxPart]
- << " Pos "<<posX[idxPart]<<" "<<posY[idxPart]<<" "<<posZ[idxPart]
- << " Forces: " << forcesX[idxPart] << " " << forcesY[idxPart] << " "<< forcesZ[idxPart] <<std::endl;
- }
- energy += potentials[idxPart]*physicalValues[idxPart] ;
- }
- });
- std::cout <<std::endl<<"Energy: "<< energy<<std::endl;
- std::cout <<std::endl<<" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "<<std::endl<<std::endl;
-
- }
- // -----------------------------------------------------
- if(FParameters::existParameter(argc, argv, FParameterDefinitions::OutputFile.options)){
- std::string name(FParameters::getStr(argc,argv,FParameterDefinitions::OutputFile.options, "output.fma"));
- FFmaGenericWriter<FReal> writer(name) ;
- //
- FSize NbPoints = loader.getNumberOfParticles();
- FReal * computedParticles = new FReal[8*NbPoints]() ;
- memset(computedParticles,0,8*NbPoints*sizeof(FReal));
- FSize j = 0 ;
- tree.forEachLeaf([&](LeafType* leaf){
- //
- // Input
- const FReal*const posX = leaf->getTargets()->getPositions()[0];
- const FReal*const posY = leaf->getTargets()->getPositions()[1];
- const FReal*const posZ = leaf->getTargets()->getPositions()[2];
- const FReal*const physicalValues = leaf->getTargets()->getPhysicalValues();
- const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
- //
- // Computed data
- 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 FSize nbParticlesInLeaf = leaf->getTargets()->getNbParticles();
- for(FSize idxPart = 0 ; idxPart < nbParticlesInLeaf ; ++idxPart){
- j = 8*indexes[idxPart];
- computedParticles[j] = posX[idxPart] ;
- computedParticles[j+1] = posY[idxPart] ;
- computedParticles[j+2] = posZ[idxPart] ;
- computedParticles[j+3] = physicalValues[idxPart] ;
- computedParticles[j+4] = potentials[idxPart] ;
- computedParticles[j+5] = forcesX[idxPart] ;
- computedParticles[j+6] = forcesY[idxPart] ;
- computedParticles[j+7] = forcesZ[idxPart] ;
- }
- });
-
- writer.writeHeader( loader.getCenterOfBox(), loader.getBoxWidth() , NbPoints, sizeof(FReal), 8) ;
- writer.writeArrayOfReal(computedParticles, 8 , NbPoints);
-
- delete[] computedParticles;
-
- //
- std::string name1( "output.fma");
- //
- FFmaGenericWriter<FReal> writer1(name1) ;
- writer1.writeDistributionOfParticlesFromOctree(&tree,NbPoints) ;
- }
- if(FParameters::existParameter(argc, argv, "-cmp") ){
- std::cout <<std::endl<<" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "<<std::endl;
- std::cout << std::scientific;
- std::cout.precision(10) ;
-
- printf("Compute Diff...");
- FMath::FAccurater<FReal> potentialDiff;
- FMath::FAccurater<FReal> fx, fy, fz, f;
- FReal energy = 0.0;
- { // Check that each particle has been summed with all other
-
- tree.forEachLeaf([&](LeafType* 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 FSize nbParticlesInLeaf = leaf->getTargets()->getNbParticles();
- const FReal*const physicalValues = leaf->getTargets()->getPhysicalValues();
-
- const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
-
- for(FSize idxPart = 0 ; idxPart < nbParticlesInLeaf ; ++idxPart){
- const FSize indexPartOrig = indexes[idxPart];
- potentialDiff.add(particles[indexPartOrig].getPotential(),potentials[idxPart]);
- fx.add(particles[indexPartOrig].getForces()[0],forcesX[idxPart]);
- fy.add(particles[indexPartOrig].getForces()[1],forcesY[idxPart]);
- fz.add(particles[indexPartOrig].getForces()[2],forcesZ[idxPart]);
- f.add(particles[indexPartOrig].getForces()[0],forcesX[idxPart]);
- f.add(particles[indexPartOrig].getForces()[1],forcesY[idxPart]);
- f.add(particles[indexPartOrig].getForces()[2],forcesZ[idxPart]);
- energy += potentials[idxPart]*physicalValues[idxPart];
- }
- });
-
- std::cout << energy << " " << energyD << std::endl;
- delete[] particles;
-
- f.setNbElements(nbParticles);
- std::cout << "FChebSymKernel Energy " << FMath::Abs(energy-energyD) << " Relative "<< FMath::Abs(energy-energyD) / FMath::Abs(energyD) <<std::endl;
- std::cout << " Potential Error " << potentialDiff << std::endl;
- std::cout << " Fx Error " << fx << std::endl;
- std::cout << " Fy Error " << fy << std::endl;
- std::cout << " Fz Error " << fz << std::endl;
- std::cout << " F Error " << f << std::endl;
- std::cout <<std::endl<<"Energy: "<< energy<<std::endl;
- std::cout <<std::endl<<" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "<<std::endl<<std::endl;
- }
-
- // -----------------------------------------------------
-
- }
-
-
- return 0;
-}