diff --git a/Src/Files/FFmaBinLoaderResult.hpp b/Src/Files/FFmaBinLoaderResult.hpp new file mode 100644 index 0000000000000000000000000000000000000000..241395b301038f011a2da2ac22231ee77bc4a832 --- /dev/null +++ b/Src/Files/FFmaBinLoaderResult.hpp @@ -0,0 +1,73 @@ +// =================================================================================== +// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Bérenger Bramas, Matthias Messner +// 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". +// =================================================================================== + +#ifndef FFMABINLOADERRESULT_HPP +#define FFMABINLOADERRESULT_HPP + +#include "./FFmaBinLoader.hpp" + +/** + * @author Cyrille Piacibello (cyrille.piacibello@inria.fr) + * @class FFmaBinLoaderResult + * + * Please read the license + * + * @brief The aim of this class is to load particles and forces from a + * file, to compare with results from any Fmm Kernel. + */ + +class FFmaBinLoaderResult : public FFmaBinLoader { + +protected: + + typedef FFmaBinLoader Parent; + + size_t removeWarning; + +public: + + FFmaBinLoaderResult(const char * filename): + Parent::FFmaBinLoader(filename) + { + } + + void fillParticle(FPoint*const inParticlePosition, FReal*const physicalValue, FReal* forceX, FReal* forceY, FReal* forceZ, FReal * const potential){ + FReal x,y,z,data,fx,fy,fz,pot; + + //Same datas as particle files + removeWarning += fread(&x, sizeof(FReal), 1, file); + removeWarning += fread(&y, sizeof(FReal), 1, file); + removeWarning += fread(&z, sizeof(FReal), 1, file); + removeWarning += fread(&data, sizeof(FReal), 1, file); + + //results datas + removeWarning += fread(&fx, sizeof(FReal), 1, file); + removeWarning += fread(&fy, sizeof(FReal), 1, file); + removeWarning += fread(&fz, sizeof(FReal), 1, file); + removeWarning += fread(&pot, sizeof(FReal), 1, file); + + //return + inParticlePosition->setPosition(x,y,z); + (*physicalValue) = data; + *forceX = fx; + *forceY = fy; + *forceZ = fz; + *potential = pot; + } +}; + + +#endif //FFMABINLOADERRESULT_HPP diff --git a/Tests/Utils/testDirectStored.cpp b/Tests/Utils/testDirectStored.cpp new file mode 100644 index 0000000000000000000000000000000000000000..06e22a91a60a66385de545372e8724119db25b4c --- /dev/null +++ b/Tests/Utils/testDirectStored.cpp @@ -0,0 +1,242 @@ +// =================================================================================== +// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Bérenger Bramas, Matthias Messner +// 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". +// =================================================================================== + + +/** + *@author Cyrille Piacibello + * + * This exec will be used to store the result of direct computation in + * Binary Files. Then, a kernel (Chebyshev) is called on the same file, and + * results are compared. + * + * Format of result file : + * Each data is a FReal : posX,posY,posZ,physicalValue,forceX,forceY,forceZ,potential + */ + + +#include <iostream> + +#include <cstdio> +#include <cstdlib> +#include "ScalFmmConfig.h" + +#include "../../Src/Utils/FTic.hpp" +#include "../../Src/Utils/FParameters.hpp" + +#include "../../Src/Containers/FOctree.hpp" + +#include "../../Src/Components/FSimpleLeaf.hpp" +#include "../../Src/Kernels/P2P/FP2P.hpp" +#include "../../Src/Kernels/P2P/FP2PParticleContainerIndexed.hpp" + +#include "../../Src/Files/FFmaBinLoader.hpp" +#include "../../Src/Files/FFmaBinLoaderResult.hpp" + +#include "../../Src/Utils/FTic.hpp" + +#include "../../Src/Kernels/Chebyshev/FChebCell.hpp" +#include "../../Src/Kernels/Interpolation/FInterpMatrixKernel.hpp" +#include "../../Src/Kernels/Chebyshev/FChebKernel.hpp" +#include "../../Src/Kernels/Chebyshev/FChebSymKernel.hpp" + +#include "../../Src/Core/FFmmAlgorithm.hpp" + +int main(int argc, char* argv[]) +{ + // get info from commandline + const char* const defaultFilename = (sizeof(FReal) == sizeof(float))? + "../Data/test20k.bin.fma.single": + "../Data/test20k.bin.fma.double"; + const char* const filename = FParameters::getStr(argc,argv,"-f", defaultFilename); + const char* const fileresult = FParameters::getStr(argc,argv,"-fr","../Data/test20k.bin.fma.double.result");; + + //For Fmm Computation + const unsigned int TreeHeight = FParameters::getValue(argc, argv, "-h", 5); + const unsigned int SubTreeHeight = FParameters::getValue(argc, argv, "-sh", 2); + + FTic time; + struct TestParticle{ + FPoint position; + FReal forces[3]; + FReal physicalValue; + FReal potential; + }; + + printf("Input is %s, \n Results will be stored in %s \n",filename,fileresult); + + { + //First Part + //Direct Computation and Storage of result + + //open particle file + FFmaBinLoader loader(filename); + if(!loader.isOpen()){ + std::cout << "Loader Error, " << filename << " is missing\n"; + return 1; + } + + + time.tic(); + TestParticle* const particles = new TestParticle[loader.getNumberOfParticles()]; + for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){ + FPoint position; + FReal physicalValue = 0.0; + loader.fillParticle(&position,&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; + } + + time.tac(); + printf("Elapsed Time for Loading %lld particles from File: \t %f\n",loader.getNumberOfParticles(),time.elapsed()); + + //Direct Computation + time.tic(); + 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); + } + } + time.tac(); + printf("Elapsed Time for Direct Computation: \t %f\n",time.elapsed()); + + //Write results in output file + FILE * fd = fopen(fileresult,"w"); + + //Size of elements to be written + + + int realSize = sizeof(FReal); + FSize nbPart = loader.getNumberOfParticles(); + FReal boxWidth = loader.getBoxWidth() / 2; + FReal centerX = loader.getCenterOfBox().getX(); + FReal centerY = loader.getCenterOfBox().getY(); + FReal centerZ = loader.getCenterOfBox().getZ(); + + fwrite(&realSize, sizeof(int), 1, fd); + fwrite(&nbPart, sizeof(FSize), 1, fd); + fwrite(&boxWidth, sizeof(FReal), 1, fd); + fwrite(¢erX, sizeof(FReal), 1, fd); + fwrite(¢erY, sizeof(FReal), 1, fd); + fwrite(¢erZ, sizeof(FReal), 1, fd); + + + for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart) + { + fwrite( particles[idxPart].position.getDataValue(),realSize,3,fd); + fwrite(&particles[idxPart].physicalValue, realSize,1,fd); + fwrite( particles[idxPart].forces, realSize,3,fd); + fwrite(&particles[idxPart].potential, realSize,1,fd); + } + + delete [] particles; + fclose(fd); + } + + //Second Part + { + //Fmm Computation and comparison to Stored results + + // begin Chebyshev kernel + + // accuracy + const unsigned int ORDER = 7; + + // typedefs + typedef FP2PParticleContainerIndexed ContainerClass; + typedef FSimpleLeaf< ContainerClass > LeafClass; + + typedef FInterpMatrixKernelR MatrixKernelClass; + typedef FChebCell<ORDER> CellClass; + typedef FOctree<CellClass,ContainerClass,LeafClass> OctreeClass; + + typedef FChebSymKernel<CellClass,ContainerClass,MatrixKernelClass,ORDER> KernelClass; + typedef FFmmAlgorithm<OctreeClass,CellClass,ContainerClass,KernelClass,LeafClass> FmmClass; + + + time.tic(); + //Open result file to fill octree for Chebyshev kernel AND compare results + FFmaBinLoaderResult loader2(fileresult); + + // init oct-tree + OctreeClass tree(TreeHeight, SubTreeHeight, loader2.getBoxWidth(), loader2.getCenterOfBox()); + + TestParticle* const particles2 = new TestParticle[loader2.getNumberOfParticles()]; + for(int idxPart = 0 ; idxPart < loader2.getNumberOfParticles() ; ++idxPart) + { + loader2.fillParticle(&particles2[idxPart].position, + &particles2[idxPart].physicalValue, + &particles2[idxPart].forces[0], + &particles2[idxPart].forces[1], + &particles2[idxPart].forces[2], + &particles2[idxPart].potential); + tree.insert(FPoint(particles2[idxPart].position),idxPart,particles2[idxPart].physicalValue); + } + time.tac(); + printf("Elapsed Time for Reading File: \t %f\n",time.elapsed()); + + time.tic(); + KernelClass kernels(TreeHeight, loader2.getBoxWidth(), loader2.getCenterOfBox()); + FmmClass algorithm(&tree, &kernels); + algorithm.execute(); + time.tac(); + printf("Elapsed Time for Fmm Computation: \t %f\n",time.elapsed()); + + + FMath::FAccurater potentialDiff; + FMath::FAccurater fx, fy, fz; + + //Compare the kernel to the stored elements + { + 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(particles2[indexPartOrig].potential,potentials[idxPart]); + fx.add(particles2[indexPartOrig].forces[0],forcesX[idxPart]); + fy.add(particles2[indexPartOrig].forces[1],forcesY[idxPart]); + fz.add(particles2[indexPartOrig].forces[2],forcesZ[idxPart]); + } + }); + // Print for information + std::cout << "Potential " << potentialDiff << std::endl; + std::cout << "Fx " << fx << std::endl; + std::cout << "Fy " << fy << std::endl; + std::cout << "Fz " << fz << std::endl; + + } + + } + return 0; +} diff --git a/Tests/Utils/testDirectStoredCreate.cpp b/Tests/Utils/testDirectStoredCreate.cpp new file mode 100644 index 0000000000000000000000000000000000000000..af92dea00435390b118f656d96d9d0781e1d54bf --- /dev/null +++ b/Tests/Utils/testDirectStoredCreate.cpp @@ -0,0 +1,150 @@ +// =================================================================================== +// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Bérenger Bramas, Matthias Messner +// 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". +// =================================================================================== + + +/** + *@author Cyrille Piacibello + * + * @brief This exec can be used to create a result file from direct + * computation of a data file. The result file will contain forces, + * positions, physical values, and potential for each particles. + * + * This file can be loaded by FFmaBinLoaderResult. + * + * Format of result file : + * Each data is a FReal : posX,posY,posZ,physicalValue,forceX,forceY,forceZ,potential + */ + + +#include <iostream> + +#include <cstdio> +#include <cstdlib> +#include "ScalFmmConfig.h" + +#include "../../Src/Utils/FTic.hpp" +#include "../../Src/Utils/FParameters.hpp" + +#include "../../Src/Components/FSimpleLeaf.hpp" +#include "../../Src/Kernels/P2P/FP2P.hpp" +#include "../../Src/Kernels/P2P/FP2PParticleContainerIndexed.hpp" + +#include "../../Src/Files/FFmaBinLoader.hpp" + +#include "../../Src/Utils/FTic.hpp" + +int main(int argc, char* argv[]) +{ + // get info from commandline + const char* const defaultFilename = (sizeof(FReal) == sizeof(float))? + "../Data/test20k.bin.fma.single": + "../Data/test20k.bin.fma.double"; + const char* const filename = FParameters::getStr(argc,argv,"-f", defaultFilename); + const char* const fileresult = FParameters::getStr(argc,argv,"-fr","../Data/test20k.bin.fma.double.result");; + + printf("Input is %s \n, Results will be stored in %s \n",filename,fileresult); + + + //Direct Computation and Storage of result + //open particle file + FFmaBinLoader loader(filename); + if(!loader.isOpen()){ + std::cout << "Loader Error, " << filename << " is missing\n"; + return 1; + } + + struct TestParticle{ + FPoint position; + FReal forces[3]; + FReal physicalValue; + FReal potential; + }; + + FTic time; + time.tic(); + + TestParticle* const particles = new TestParticle[loader.getNumberOfParticles()]; + for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){ + FPoint position; + FReal physicalValue = 0.0; + loader.fillParticle(&position,&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; + } + + time.tac(); + printf("Elapsed Time for Loading File: \t %f\n",time.elapsed()); + + //Direct Computation + time.tic(); + 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); + } + } + time.tac(); + printf("Elapsed Time for Direct Computation: \t %f\n",time.elapsed()); + + time.tic(); + //Write results in output file + FILE * fd = fopen(fileresult,"w"); + + //Size of elements to be written + size_t realSize = sizeof(FReal); + + fwrite(&realSize, sizeof(int), 1, fd); + FSize nbPart = loader.getNumberOfParticles(); + FReal boxWidth = loader.getBoxWidth(); + FReal centerX = loader.getCenterOfBox().getX(); + FReal centerY = loader.getCenterOfBox().getY(); + FReal centerZ = loader.getCenterOfBox().getZ(); + + fwrite(&nbPart, sizeof(FSize), 1, fd); + fwrite(&boxWidth, sizeof(FReal), 1, fd); + fwrite(¢erX, sizeof(FReal), 1, fd); + fwrite(¢erY, sizeof(FReal), 1, fd); + fwrite(¢erZ, sizeof(FReal), 1, fd); + + + for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart) + { + fwrite( particles[idxPart].position.getDataValue(),realSize,3,fd); + fwrite( particles[idxPart].forces,realSize,3,fd); + fwrite(&particles[idxPart].physicalValue,realSize,1,fd); + fwrite(&particles[idxPart].potential,realSize,1,fd); + } + + time.tac(); + printf("Elapsed Time for Writing results into file: \t %f\n",time.elapsed()); + + delete [] particles; + fclose(fd); + + + return 0; +} diff --git a/Tests/Utils/testLoaderFMABinResult.cpp b/Tests/Utils/testLoaderFMABinResult.cpp new file mode 100644 index 0000000000000000000000000000000000000000..d3f8e35a1dfedcca30964e04f18478839c7a5735 --- /dev/null +++ b/Tests/Utils/testLoaderFMABinResult.cpp @@ -0,0 +1,104 @@ +// =================================================================================== +// Copyright ScalFmm 2011 INRIA, Olivier Coulaud, Bérenger Bramas, Matthias Messner +// 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". +// =================================================================================== + +#include <iostream> + +#include <cstdio> +#include <cstdlib> +#include <time.h> + +#include "../../Src/Utils/FParameters.hpp" +#include "../../Src/Utils/FTic.hpp" + +#include "../../Src/Containers/FOctree.hpp" +#include "../../Src/Containers/FVector.hpp" + +#include "../../Src/Utils/FAssert.hpp" +#include "../../Src/Utils/FPoint.hpp" + +#include "../../Src/Components/FBasicCell.hpp" + +#include "../../Src/Components/FSimpleLeaf.hpp" + +#include "../../Src/Files/FFmaBinLoaderResult.hpp" + +#include "../../Src/Components/FBasicParticleContainer.hpp" + + +/** + * In this file we show an example of FBasicLoader use +* Inserting 2000000 particles ... +* Done (5.77996). +* Deleting particles ... +* Done (0.171918). + */ + +int main(int argc, char ** argv ){ + typedef FBasicParticleContainer<1> ContainerClass; + typedef FSimpleLeaf< ContainerClass > LeafClass; + typedef FOctree< FBasicCell, ContainerClass , LeafClass > OctreeClass; + ///////////////////////What we do///////////////////////////// + std::cout << ">> This executable is useless to execute.\n"; + std::cout << ">> It is only interesting to wath the code to understand\n"; + std::cout << ">> how to use the FMA loader\n"; + ////////////////////////////////////////////////////////////// + + // Use testLoaderFMABinCreate.exe to create this file + FTic counter; + const char* const filename = FParameters::getStr(argc,argv,"-f", "../Data/test20k.bin.fma.double.res"); + std::cout << "Opening : " << filename << "\n"; + + // open basic result particles loader + FFmaBinLoaderResult loader(filename); + if(!loader.isOpen()){ + std::cout << "Loader Error, " << filename << "is missing\n"; + return 1; + } + + { + // octree + OctreeClass tree(FParameters::getValue(argc,argv,"-h", 5), FParameters::getValue(argc,argv,"-sh", 3), + loader.getBoxWidth(), loader.getCenterOfBox()); + + // ----------------------------------------------------- + std::cout << "Inserting " << loader.getNumberOfParticles() << " particles ..." << std::endl; + counter.tic(); + + FPoint particlePosition; + FReal physicalValue = 0.0; + FReal forces[3]; + FReal potential; + for(int idxPart = 0 ; idxPart < loader.getNumberOfParticles() ; ++idxPart){ + loader.fillParticle(&particlePosition,&physicalValue,forces,&forces[1],&forces[2],&potential); + tree.insert(particlePosition,physicalValue); + } + + counter.tac(); + std::cout << "Done " << "(" << counter.elapsed() << ")." << std::endl; + + // ----------------------------------------------------- + std::cout << "Deleting particles ..." << std::endl; + counter.tic(); + } + counter.tac(); + std::cout << "Done " << "(" << counter.elapsed() << ")." << std::endl; + // ----------------------------------------------------- + + return 0; +} + + +