Skip to content

S
ScalFMM
Commit 3abae251 authored by COULAUD Olivier's avatar COULAUD Olivier
Browse files
Options

delete file

parent 7318012f
Hide whitespace changes
Inline Side-by-side
Showing with 0 additions and 357 deletions
Tests/Utils/generateDistrib.cpp deleted 100644 → 0
View file @ 7318012f
// See LICENCE file at project root
// ==== CMAKE =====
// ================
//
// make testPeriodicAlgorithm
// mpirun --oversubscribe -np 4 --tag-output Tests/Release/testPeriodicAlgorithm -h 3 -sh 1
#include <iostream>
#include <vector>
#include <array>
#include <algorithm>
#include <stdexcept>
//#include <cstdio>
#include <cstdlib>
//
// Specific part
//
#include <string>
//
// Generic part
//
#include "ScalFmmConfig.h"
#include "Files/FFmaGenericLoader.hpp"
#include "Utils/FParameters.hpp"
#include "Utils/FParameterNames.hpp"
//
// Order of the Interpolation approximation
using FReal_t = double;
/** This program show an example of use of
* the fmm basic algo
* it also check that each particles is impacted each other particles
*/
/* Mock particle structure to balance the tree over the processes. */
struct TestParticle{
FSize index; // Index of the particle in the original file.
FPoint<FReal_t> position; // Spatial position of the particle.
FReal_t physicalValue; // Physical value of the particle.
/* Returns the particle position. */
const FPoint<FReal_t>& getPosition(){
return position;
}
};
void setOneParticleInEachCell(const FReal_t &BoxWidth , const int &Height, std::vector<TestParticle> &particles,
MortonIndex & startIdx, MortonIndex & endIdx, int nbProc = 1, int myId = 0){
//
auto leavesLevel = Height- 1 ;
auto BoxCellWidth = BoxWidth/ (2 << (leavesLevel-1)) ;
auto nbLeaves = 2 << (leavesLevel-1) ;
nbLeaves = nbLeaves*nbLeaves*nbLeaves ;
std::cout << " H: " << leavesLevel << " nbCells: " << nbLeaves << " cellwidth " << BoxCellWidth << std::endl ;
auto bloc = nbLeaves / nbProc ;
startIdx = bloc*myId ;
endIdx = startIdx+bloc;
if(myId == nbProc-1){
endIdx = nbLeaves ;}
particles.resize(endIdx-startIdx) ;
FReal_t val = 0.1 *( myId %2 ==0 ? 1 : -1 );
std::cout << "procId " << myId << " MortonIdx: [ " << startIdx << " , " << endIdx <<" [" <<std::endl ;
for(MortonIndex idx = startIdx, pos=0 ; idx < endIdx ; ++idx,++pos ) {
FTreeCoordinate coor ;
coor.setPositionFromMorton(idx) ;
particles[pos].position = { (0.5+coor[0])*BoxCellWidth,(0.5+coor[1])*BoxCellWidth,(0.5+coor[2])*BoxCellWidth } ;
particles[pos].index = idx ;
particles[pos].physicalValue =val ;
val *= -1.0 ;
std::cout << idx << " " << particles[pos].index << " " << particles[pos].position << " " << particles[pos].physicalValue << std::endl;
}
}
void setOneParticleInInCellsLinesZ(const std::vector<std::array<int,2> >& Lines, const FReal_t &BoxWidth , const int &Height, std::vector<TestParticle> &particles,
MortonIndex & startIdx, MortonIndex & endIdx, int nbProc = 1, int myId = 0){
//
auto leavesLevel = Height - 1 ;
auto BoxCellWidth = BoxWidth/ (2 << (leavesLevel-1)) ;
auto nb1DLeaves = 2 << (leavesLevel-1) ;
auto nbLines = Lines.size();
std::vector<MortonIndex> cells(nb1DLeaves*nbLines) ;
// fill cells
std::cout << "Cells " << std::endl ;
for (auto v : Lines){
std::cout << " ("<< v[0]<< " "<< v[1] << ") ";
}
std::cout << std::endl;
std::size_t posI = 0 ;
for(std::size_t i = 0 ; i < nbLines ; ++i ) {
auto & index = Lines[i] ;
for(MortonIndex idx = 0; idx < nb1DLeaves ; ++idx ) {
FTreeCoordinate coor(index[0],index[1],idx) ;
cells[posI] = coor.getMortonIndex() ;
++posI ;
}
}
std::cout << "Cells " << std::endl ;
for (auto v : cells){
std::cout << " "<< v;
}
std::cout << std::endl ;
std::sort(cells.begin(), cells.end()) ;
std::cout << "Sorted Cells " << std::endl ;
for (auto v : cells){
std::cout << " "<< v;
}
std::cout << std::endl ;
// nbLeaves = nbLeaves*nbLeaves*nbLeaves ;
std::cout << " H: " << leavesLevel << " nb1DLeaves "<< nb1DLeaves
<<" nbCells: " << cells.size() << " cellwidth " << BoxCellWidth << std::endl ;
auto nbLeaves = cells.size() ;
auto bloc = nbLeaves / nbProc ;
startIdx = bloc*myId ;
endIdx = startIdx+bloc;
if(myId == nbProc-1){
endIdx = nbLeaves ;}
particles.resize(endIdx-startIdx) ;
FReal_t val ;
if(particles.size() ==1) {
val = 0.1*( myId %2 ==0 ? 1 : -1 );
}
else {
val = 0.1;
}
std::cout << "procId " << myId << " MortonIdx: [ " << startIdx << " , " << endIdx <<" [" <<std::endl ;
for(MortonIndex idx = startIdx,pos=0; idx < endIdx ; ++idx,++pos ) {
FTreeCoordinate coor(cells[idx]) ;
particles[pos].position = { (0.5+coor[0])*BoxCellWidth,(0.5+coor[1])*BoxCellWidth,(0.5+coor[2])*BoxCellWidth } ;
particles[pos].index = cells[idx];
particles[pos].physicalValue =val ;
val *= -1.0 ;
std::cout << idx << " " << coor << " "<< particles[pos].index << " "
<< particles[pos].position << " " << particles[pos].physicalValue << std::endl;
}
std::cout << " END setOneParticleInInCellsLinesZ"<<std::endl;
}
void setOneParticleInInCellsLineZ(const FReal_t &BoxWidth , const int &Height,std::vector<TestParticle> &particles,
MortonIndex & startIdx, MortonIndex & endIdx, int nbProc = 1, int myId = 0){
//
auto leavesLevel = Height - 1 ;
auto BoxCellWidth = BoxWidth/ (2 << (leavesLevel-1)) ;
auto nbLeaves = 2 << (leavesLevel-1) ;
// nbLeaves = nbLeaves*nbLeaves*nbLeaves ;
std::cout << " H: " << leavesLevel << " nbCells: " << nbLeaves << " cellwidth " << BoxCellWidth << std::endl ;
auto bloc = nbLeaves / nbProc ;
startIdx = bloc*myId ;
endIdx = startIdx+bloc;
if(myId == nbProc-1){
endIdx = nbLeaves ;}
particles.resize(endIdx-startIdx) ;
FReal_t val ;
if(particles.size() ==1) {
val = 0.1*( myId %2 ==0 ? 1 : -1 );
}
else {
val = 0.1;
}
std::cout << "procId " << myId << " MortonIdx: [ " << startIdx << " , " << endIdx <<" [" <<std::endl ;
int mid = 1 ;//nbLeaves/2;
for(MortonIndex idx = startIdx,pos=0; idx < endIdx ; ++idx,++pos ) {
FTreeCoordinate coor(mid,mid,idx) ;
particles[pos].position = { (0.5+coor[0])*BoxCellWidth,(0.5+coor[1])*BoxCellWidth,(0.5+coor[2])*BoxCellWidth } ;
particles[pos].index = coor.getMortonIndex() ;
particles[pos].physicalValue =val ;
val *= -1.0 ;
// std::cout << idx << " " << particles[pos].index << " " << particles[pos].position << " " << particles[pos].physicalValue << std::endl;
}
}
#ifdef COMPIL
void setOneParticleInInCellsLineY(OctreeClass & tree, std::vector<TestParticle> &particles,
MortonIndex & startIdx, MortonIndex & endIdx, int nbProc = 1, int myId = 0){
//
auto leavesLevel = tree.getHeight() - 1 ;
auto BoxCellWidth = tree .getBoxWidth()/ (2 << (leavesLevel-1)) ;
auto nbLeaves = 2 << (leavesLevel-1) ;
// nbLeaves = nbLeaves*nbLeaves*nbLeaves ;
std::cout << " H: " << leavesLevel << " nbCells: " << nbLeaves << " cellwidth " << BoxCellWidth << std::endl ;
auto bloc = nbLeaves / nbProc ;
startIdx = bloc*myId ;
endIdx = startIdx+bloc;
if(myId == nbProc-1){
endIdx = nbLeaves ;}
particles.resize(endIdx-startIdx) ;
FReal_t val ;
if(particles.size() ==1) {
val = 0.1*( myId %2 ==0 ? 1 : -1 );
}
else {
val = 0.1;
}
std::cout << "procId " << myId << " MortonIdx: [ " << startIdx << " , " << endIdx <<" [" <<std::endl ;
int mid = 1 ;//nbLeaves/2;
for(MortonIndex idx = startIdx,pos=0; idx < endIdx ; ++idx,++pos ) {
FTreeCoordinate coor(mid,idx,mid) ;
particles[pos].position = { (0.5+coor[0])*BoxCellWidth,(0.5+coor[1])*BoxCellWidth,(0.5+coor[2])*BoxCellWidth } ;
particles[pos].index = coor.getMortonIndex() ;
particles[pos].physicalValue =val ;
val *= -1.0 ;
// std::cout << idx << " " << particles[pos].index << " " << particles[pos].position << " " << particles[pos].physicalValue << std::endl;
}
}
void setOneParticleInInCellsLineX(OctreeClass & tree, std::vector<TestParticle> &particles,
MortonIndex & startIdx, MortonIndex & endIdx, int nbProc = 1, int myId = 0){
//
auto leavesLevel = tree.getHeight() - 1 ;
auto BoxCellWidth = tree .getBoxWidth()/ (2 << (leavesLevel-1)) ;
auto nbLeaves = 2 << (leavesLevel-1) ;
// nbLeaves = nbLeaves*nbLeaves*nbLeaves ;
std::cout << " H: " << leavesLevel << " nbCells: " << nbLeaves << " cellwidth " << BoxCellWidth << std::endl ;
auto bloc = nbLeaves / nbProc ;
startIdx = bloc*myId ;
endIdx = startIdx+bloc;
if(myId == nbProc-1){
endIdx = nbLeaves ;}
particles.resize(endIdx-startIdx) ;
FReal_t val ;
if(particles.size() ==1) {
val = 0.1*( myId %2 ==0 ? 1 : -1 );
}
else {
val = 0.1;
}
std::cout << "procId " << myId << " MortonIdx: [ " << startIdx << " , " << endIdx <<" [" <<std::endl ;
MortonIndex mid = 1 ;//nbLeaves/2;
for(MortonIndex idx = startIdx,pos=0; idx < endIdx ; ++idx,++pos ) {
FTreeCoordinate coor(idx,mid,mid) ;
particles[pos].position = { (0.5+coor[0])*BoxCellWidth,(0.5+coor[1])*BoxCellWidth,(0.5+coor[2])*BoxCellWidth } ;
particles[pos].index = coor.getMortonIndex() ;
particles[pos].physicalValue =val ;
val *= -1.0 ;
// std::cout << idx << " " << particles[pos].index << " " << particles[pos].position << " " << particles[pos].physicalValue << std::endl;
}
}
#endif
namespace Param {
const FParameterNames Charge
= {{"-charge"}, "generate physical values between -1 and 1 otherwise generate between 0 and 1"};
const FParameterNames ZM
= {{"-zeromean"}, "the average of the physical values is zero"};
}
int main(int argc, char ** argv){
FHelpDescribeAndExit(argc, argv,
"Generate particules distribution . ",
FParameterDefinitions::OutputFile,
FParameterDefinitions::OctreeHeight,
Param::Charge, Param::ZM) ;
const std::string filename = FParameters::getStr(argc,argv, FParameterDefinitions::OutputFile.options, "output.fma");
const unsigned int TreeHeight = FParameters::getValue(argc, argv, FParameterDefinitions::OctreeHeight.options, 5);
FReal_t boxWidth = 1.0 ;
FPoint<FReal_t> centerOfBox = {0.5,0.5,0.5} ;
std::vector<TestParticle> particles;
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
//
MortonIndex startIdx=0,endIdx=0;
//#ifdef SLICES
// Check a slice (Y,Z) for a given X = 0
#ifdef ss
auto nb1DLeaves = 2 << (TreeHeight-2) ;
MortonIndex mPos = 0 , mPos1 = nb1DLeaves-1 ;
std::vector<std::array<MortonIndex,2> > Lines ;
for( MortonIndex i = 0 ; i < nb1DLeaves; ++i){
// first slice
auto tmpPos = mPos;
std::array<MortonIndex,2> tmp = {mPos,i};
// if( i==0 || i== nb1DLeaves-1)
Lines.emplace_back(tmp);
// second slice
tmp[0] = ++tmpPos;
// if( i==0|| i== nb1DLeaves-1
Lines.emplace_back(tmp);
// tmp[0] = ++tmpPos ;
// // if( i==0 || i== nb1DLeaves-1)
// Lines.emplace_back(tmp);
// tmp[0] = mPos1;
// Lines.emplace_back(tmp);
}
#else
std::vector<std::array<int,2> > Lines ={{1,1}} ;
#endif
for ( auto v : Lines){
std::cout << " {" << v[0] << ", " << v[1]<< "}" ;
}
std::cout << std::endl;
setOneParticleInInCellsLinesZ(Lines, boxWidth,TreeHeight, particles,startIdx, endIdx);
//#endif
// setOneParticleInInCellsLineZ(boxWidth,TreeHeight, particles,startIdx, endIdx);
// setOneParticleInInCellsLineY(tree, particles,startIdx);
// setOneParticleInInCellsLineY(tree, particles,startIdx);
//setOneParticleInInCellsLineX(tree, particles,startIdx, endIdx);
// setOneParticleInEachCell(boxWidth,TreeHeight, particles,startIdx, endIdx);
FSize NbPoints = particles.size();
if(FParameters::existParameter(argc, argv, "-zeromean")){
// set the total charge to zero in order to have a solution
FReal_t charge = 0 ;
for ( auto v : particles){
charge += v.physicalValue ;
}
charge /= static_cast<FReal_t>(NbPoints) ;
for ( auto v : particles){
v.physicalValue -= charge ;
}
charge = 0.0 ;
for ( auto v : particles){
charge += v.physicalValue ;
}
std::cout << " GlobalCharge: "<< charge << std::endl;
}
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
// -----------------------------------------------------
std::cout << "Write "<< NbPoints <<" particles to '" << filename << "'" << std::endl;
FFmaGenericWriter<FReal_t> writer(filename);
constexpr int nbData = 4 ;
writer.writeHeader(centerOfBox, boxWidth, NbPoints, sizeof(FReal_t),nbData);
std::vector<FmaRWParticle<FReal_t, nbData ,nbData> > finalPart(NbPoints);
for (int i=0 ; i < NbPoints ; ++i){
finalPart[i].setPosition(particles[i].position) ;
finalPart[i].setPhysicalValue(particles[i].physicalValue) ;
}
writer.writeArrayOfReal(finalPart[0].getPtrFirstData(), nbData, NbPoints);
std::cout << "End of writing" <<std::endl;
return 0;
}
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