Commit 790b6440 authored by COULAUD Olivier's avatar COULAUD Olivier

Bug with BPC is fixed

parent f6d6f6e8
...@@ -32,7 +32,6 @@ ...@@ -32,7 +32,6 @@
#include "Utils/FParameters.hpp" #include "Utils/FParameters.hpp"
#include "Utils/FParameterNames.hpp" #include "Utils/FParameterNames.hpp"
// //
// Order of the Interpolation approximation // Order of the Interpolation approximation
static constexpr unsigned ORDER = 6 ; static constexpr unsigned ORDER = 6 ;
...@@ -115,6 +114,7 @@ int main(int argc, char* argv[]) ...@@ -115,6 +114,7 @@ int main(int argc, char* argv[])
// Initialize values for MPI // Initialize values for MPI
FMpi app(argc,argv); FMpi app(argc,argv);
const bool masterIO = ( app.global().processId() == 0 );
// //
// Initialize timer // Initialize timer
FTic time; FTic time;
...@@ -137,7 +137,7 @@ int main(int argc, char* argv[]) ...@@ -137,7 +137,7 @@ int main(int argc, char* argv[])
FSize localParticlesNumber = 0 ; FSize localParticlesNumber = 0 ;
// ----------------------------------------------------- // -----------------------------------------------------
if(app.global().processId() == 0){ if(masterIO){
std::cout << "Loading & Inserting " << loader.getNumberOfParticles() std::cout << "Loading & Inserting " << loader.getNumberOfParticles()
<< " particles ..." << std::endl << " particles ..." << std::endl
<<" Box: "<< std::endl <<" Box: "<< std::endl
...@@ -214,7 +214,7 @@ int main(int argc, char* argv[]) ...@@ -214,7 +214,7 @@ int main(int argc, char* argv[])
MPI_Reduce(&timeUsed,&minTime,1,MPI_DOUBLE,MPI_MIN,0,app.global().getComm()); MPI_Reduce(&timeUsed,&minTime,1,MPI_DOUBLE,MPI_MIN,0,app.global().getComm());
MPI_Reduce(&timeUsed,&maxTime,1,MPI_DOUBLE,MPI_MAX,0,app.global().getComm()); MPI_Reduce(&timeUsed,&maxTime,1,MPI_DOUBLE,MPI_MAX,0,app.global().getComm());
if(app.global().processId() == 0){ if(masterIO){
std::cout << "readinsert-time-min:" << minTime std::cout << "readinsert-time-min:" << minTime
<< " readinsert-time-max:" << maxTime << " readinsert-time-max:" << maxTime
<< std::endl; << std::endl;
...@@ -225,8 +225,9 @@ int main(int argc, char* argv[]) ...@@ -225,8 +225,9 @@ int main(int argc, char* argv[])
FAbstractAlgorithm * algorithm = nullptr; FAbstractAlgorithm * algorithm = nullptr;
FAlgorithmTimers * timer = nullptr; FAlgorithmTimers * timer = nullptr;
{ // ----------------------------------------------------- { // -----------------------------------------------------
std::cout << "\n"<<interpolationType<<" FMM Proc (ORDER="<< ORDER << ") ... " << std::endl; if(masterIO) {
std::cout << "\n"<<interpolationType<<" FMM Proc (ORDER="<< ORDER << ") ... " << std::endl;
}
time.tic(); time.tic();
// Kernels to use (pointer because of the limited size of the stack) // Kernels to use (pointer because of the limited size of the stack)
...@@ -239,8 +240,10 @@ int main(int argc, char* argv[]) ...@@ -239,8 +240,10 @@ int main(int argc, char* argv[])
// //
// periodic FMM algorithm // periodic FMM algorithm
FmmClassProcPER algoPer(app.global(),&tree, aboveTree); FmmClassProcPER algoPer(app.global(),&tree, aboveTree);
KernelClass kernelsPer(algoPer.extendedTreeHeight(), algoPer.extendedBoxWidth(), KernelClass kernelsPer(algoPer.extendedTreeHeight(), algoPer.extendedBoxWidth(),
algoPer.extendedBoxCenter(),&MatrixKernel); algoPer.extendedBoxCenter(),&MatrixKernel);
algoPer.setKernel(&kernelsPer); algoPer.setKernel(&kernelsPer);
/////////////////////////////////////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////////////////////////////////////
if(! periodicCondition) {// Non periodic case if(! periodicCondition) {// Non periodic case
...@@ -252,12 +255,12 @@ int main(int argc, char* argv[]) ...@@ -252,12 +255,12 @@ int main(int argc, char* argv[])
timer = &algoPer ; timer = &algoPer ;
} }
// //
// FMM exectution FFmmFarField FFmmNearField // FMM exectution FFmmFarField FFmmNearField FFmmP2M|FFmmM2M|FFmmM2L|FFmmL2L
algorithm->execute(); algorithm->execute();
// //
algorithm->getMortonLeafDistribution(mortonLeafDistribution); algorithm->getMortonLeafDistribution(mortonLeafDistribution);
// if(app.global().processId() == 0) if(masterIO)
{ {
std::cout << app.global().processId() <<" Morton distribution "<< std::endl ; std::cout << app.global().processId() <<" Morton distribution "<< std::endl ;
for(auto v : mortonLeafDistribution) for(auto v : mortonLeafDistribution)
...@@ -268,10 +271,10 @@ int main(int argc, char* argv[]) ...@@ -268,10 +271,10 @@ int main(int argc, char* argv[])
app.global().barrier(); app.global().barrier();
time.tac(); time.tac();
timeUsed = time.elapsed(); timeUsed = time.elapsed();
std::cout << "Done " << "(@Algorithm = " << time.elapsed() << " s)." << std::endl;
MPI_Reduce(&timeUsed,&minTime,1,MPI_DOUBLE,MPI_MIN,0,app.global().getComm()); MPI_Reduce(&timeUsed,&minTime,1,MPI_DOUBLE,MPI_MIN,0,app.global().getComm());
MPI_Reduce(&timeUsed,&maxTime,1,MPI_DOUBLE,MPI_MAX,0,app.global().getComm()); MPI_Reduce(&timeUsed,&maxTime,1,MPI_DOUBLE,MPI_MAX,0,app.global().getComm());
if(app.global().processId() == 0){ if(masterIO){
std::cout << "Done " << "(@Algorithm = " << time.elapsed() << " s)." << std::endl;
std::cout << "exec-time-min:" << minTime std::cout << "exec-time-min:" << minTime
<< " exec-time-max:" << maxTime << " exec-time-max:" << maxTime
<< std::endl; << std::endl;
...@@ -331,7 +334,7 @@ int main(int argc, char* argv[]) ...@@ -331,7 +334,7 @@ int main(int argc, char* argv[])
}); });
FReal gloEnergy = app.global().reduceSum(energy); FReal gloEnergy = app.global().reduceSum(energy);
FReal TotalPhysicalValue = app.global().reduceSum(locTotalPhysicalValue); FReal TotalPhysicalValue = app.global().reduceSum(locTotalPhysicalValue);
if(0 == app.global().processId()){ if(masterIO){
std::cout <<std::endl<<"aboveRoot: " << aboveTree << " Energy: "<< gloEnergy <<" TotalPhysicalValue: " << TotalPhysicalValue<< std::endl; std::cout <<std::endl<<"aboveRoot: " << aboveTree << " Energy: "<< gloEnergy <<" TotalPhysicalValue: " << TotalPhysicalValue<< std::endl;
std::cout <<std::endl <<" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "<<std::endl<<std::endl; std::cout <<std::endl <<" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "<<std::endl<<std::endl;
} }
......
...@@ -198,6 +198,7 @@ int main(int argc, char* argv[]) ...@@ -198,6 +198,7 @@ int main(int argc, char* argv[])
timer = &algoPer ; timer = &algoPer ;
} }
// //
algorithm->execute(); // Here the call of the FMM algorithm algorithm->execute(); // Here the call of the FMM algorithm
// //
...@@ -250,10 +251,6 @@ int main(int argc, char* argv[]) ...@@ -250,10 +251,6 @@ int main(int argc, char* argv[])
<< " Pos "<<posX[idxPart]<<" "<<posY[idxPart]<<" "<<posZ[idxPart] << " Pos "<<posX[idxPart]<<" "<<posY[idxPart]<<" "<<posZ[idxPart]
<< " Forces: " << forcesX[idxPart] << " " << forcesY[idxPart] << " "<< forcesZ[idxPart] <<std::endl; << " Forces: " << forcesX[idxPart] << " " << forcesY[idxPart] << " "<< forcesZ[idxPart] <<std::endl;
} }
// if(numLeaf==288){
// std::cout << count << " idxPart " << idxPart << potentials[idxPart]<<" "<<
// physicalValues[idxPart] << " energy " << energy << std::endl;
// }
energy += potentials[idxPart]*physicalValues[idxPart] ; energy += potentials[idxPart]*physicalValues[idxPart] ;
TotalPhysicalValue += physicalValues[idxPart] ; TotalPhysicalValue += physicalValues[idxPart] ;
// ++count ; // ++count ;
......
...@@ -19,7 +19,7 @@ version 1.1 ...@@ -19,7 +19,7 @@ version 1.1
- Improvement in CmakeLists - FUSE is working again - - Improvement in CmakeLists - FUSE is working again -
- Compile with Intel and CLang (with a lot of warnings) - Compile with Intel and CLang (with a lot of warnings)
- Now morse_cmake is a git submodule - Now morse_cmake is a git submodule
- - Fix several bugs in pseudo-periodic algorithms (sequential and parallel)
1.5 1.5
......
...@@ -35,8 +35,17 @@ To obtain ScalFMM (develop branch) and its git submodules do ...@@ -35,8 +35,17 @@ To obtain ScalFMM (develop branch) and its git submodules do
git clone --recursive git@gitlab.inria.fr:solverstack/ScalFMM.git -b develop git clone --recursive git@gitlab.inria.fr:solverstack/ScalFMM.git -b develop
``` ```
``` bash or
```bash
git clone git@gitlab.inria.fr:solverstack/ScalFMM.git
cd ScalFMM
git submodule init
git submodule update
```
# Move to the build folder # Move to the build folder
``` bash
cd scalfmm/Build cd scalfmm/Build
# Use cmake, with relevant options # Use cmake, with relevant options
cmake .. # -DSCALFMM_USE_MPI=ON cmake .. # -DSCALFMM_USE_MPI=ON
......
...@@ -198,12 +198,14 @@ protected: ...@@ -198,12 +198,14 @@ protected:
if(operationsToProceed & FFmmP2M) bottomPass(); if(operationsToProceed & FFmmP2M) bottomPass();
if(operationsToProceed & FFmmM2M) upwardPass(); if(operationsToProceed & FFmmM2M) {
upwardPass();
// before downward pass we have to perform the periodicity
this->processPeriodicLevels();
}
if(operationsToProceed & FFmmM2L){ if(operationsToProceed & FFmmM2L){
transferPass(); transferPass();
// before downward pass we have to perform the periodicity
processPeriodicLevels();
} }
if(operationsToProceed & FFmmL2L) downardPass(); if(operationsToProceed & FFmmL2L) downardPass();
...@@ -561,10 +563,10 @@ protected: ...@@ -561,10 +563,10 @@ protected:
// upperCells[offsetRealTree-1] is root cell // upperCells[offsetRealTree-1] is root cell
CellClass*const upperCells = new CellClass[offsetRealTree]; CellClass*const upperCells = new CellClass[offsetRealTree];
for(int i = 0; i < offsetRealTree; ++i) { for(int i = 0; i < offsetRealTree; ++i) {
upperCells[i].setLevel(i+1); // TODO: check that the level is correct (see other todos under) upperCells[i].setLevel(i+1);
} }
{ { // Build Expansion at level 0
typename OctreeClass::Iterator octreeIterator(tree); typename OctreeClass::Iterator octreeIterator(tree);
octreeIterator.gotoLeft(); octreeIterator.gotoLeft();
...@@ -584,7 +586,7 @@ protected: ...@@ -584,7 +586,7 @@ protected:
std::array<const symbolic_data_t*, 8> level_1_symbolics; std::array<const symbolic_data_t*, 8> level_1_symbolics;
std::transform(children, children+8, level_1_symbolics.begin(), std::transform(children, children+8, level_1_symbolics.begin(),
[](CellClass* c) {return c;}); [](CellClass* c) {return c;});
// TODO check that parent symbolic has the right level //
kernels->M2M( kernels->M2M(
real_tree_root_multipole, real_tree_root_multipole,
real_tree_root_symbolic, real_tree_root_symbolic,
...@@ -598,6 +600,7 @@ protected: ...@@ -598,6 +600,7 @@ protected:
// Copy virtual cell at given level into virtual children // Copy virtual cell at given level into virtual children
FMemUtils::setall(children,&upperCells[idxLevel],8); FMemUtils::setall(children,&upperCells[idxLevel],8);
multipole_t* const virtual_parent_multipole multipole_t* const virtual_parent_multipole
= &((&upperCells[idxLevel-1])->getMultipoleData()); = &((&upperCells[idxLevel-1])->getMultipoleData());
const symbolic_data_t* const virtual_parent_symbolic const symbolic_data_t* const virtual_parent_symbolic
...@@ -618,7 +621,9 @@ protected: ...@@ -618,7 +621,9 @@ protected:
); );
} }
} }
//
// M2L PASS
//
CellClass*const downerCells = new CellClass[offsetRealTree]; CellClass*const downerCells = new CellClass[offsetRealTree];
for(int i = 0; i < offsetRealTree; ++i) { for(int i = 0; i < offsetRealTree; ++i) {
downerCells[i].setLevel(i+1); // TODO: check that the level is correct (see other todos under) downerCells[i].setLevel(i+1); // TODO: check that the level is correct (see other todos under)
...@@ -627,8 +632,8 @@ protected: ...@@ -627,8 +632,8 @@ protected:
{ {
const int idxUpperLevel = 2; const int idxUpperLevel = 2;
const CellClass* neighbors[342]; const CellClass* neighbors[342]{};
int neighborPositions[342]; int neighborPositions[342]{};
int counter = 0; int counter = 0;
// The cells are all the same, duplicate the central cell and // The cells are all the same, duplicate the central cell and
...@@ -637,7 +642,7 @@ protected: ...@@ -637,7 +642,7 @@ protected:
for(int idxY = -3 ; idxY <= 2 ; ++idxY){ for(int idxY = -3 ; idxY <= 2 ; ++idxY){
for(int idxZ = -3 ; idxZ <= 2 ; ++idxZ){ for(int idxZ = -3 ; idxZ <= 2 ; ++idxZ){
if( FMath::Abs(idxX) > 1 || FMath::Abs(idxY) > 1 || FMath::Abs(idxZ) > 1){ if( FMath::Abs(idxX) > 1 || FMath::Abs(idxY) > 1 || FMath::Abs(idxZ) > 1){
neighbors[counter] = &upperCells[idxUpperLevel-1]; neighbors[counter] = &upperCells[idxUpperLevel-1];
neighborPositions[counter] = neighIndex(idxX,idxY,idxZ); neighborPositions[counter] = neighIndex(idxX,idxY,idxZ);
++counter; ++counter;
} }
...@@ -649,13 +654,13 @@ protected: ...@@ -649,13 +654,13 @@ protected:
= &(downerCells[idxUpperLevel-1].getLocalExpansionData()); = &(downerCells[idxUpperLevel-1].getLocalExpansionData());
const symbolic_data_t* const target_symbolic const symbolic_data_t* const target_symbolic
= &downerCells[idxUpperLevel-1]; = &downerCells[idxUpperLevel-1];
std::array<const multipole_t*, 342> neighbor_multipoles;
std::array<const symbolic_data_t*, 342> neighbor_symbolics; std::array<const multipole_t*, 342> neighbor_multipoles{};
std::array<const symbolic_data_t*, 342> neighbor_symbolics{};
std::transform(neighbors, neighbors+counter, neighbor_multipoles.begin(), std::transform(neighbors, neighbors+counter, neighbor_multipoles.begin(),
[](const CellClass* c) {return &(c->getMultipoleData());}); [](const CellClass* c) {return &(c->getMultipoleData());});
std::transform(neighbors, neighbors+counter, neighbor_symbolics.begin(), std::transform(neighbors, neighbors+counter, neighbor_symbolics.begin(),
[](const CellClass* c) {return c;}); [](const CellClass* c) {return c;});
// compute M2L // compute M2L
kernels->M2L( kernels->M2L(
target_local_exp, target_local_exp,
...@@ -665,8 +670,11 @@ protected: ...@@ -665,8 +670,11 @@ protected:
neighborPositions, neighborPositions,
counter); counter);
} }
// note: the triple loop bounds are not the same than for the previous piece of code // note: the triple loop bounds are not the same than for the previous piece of code
// which handles the topmost virtual cell // which handles the topmost virtual cell
for(int idxUpperLevel = 3 ; idxUpperLevel <= offsetRealTree ; ++idxUpperLevel){ for(int idxUpperLevel = 3 ; idxUpperLevel <= offsetRealTree ; ++idxUpperLevel){
const CellClass* neighbors[342]; const CellClass* neighbors[342];
int neighborPositions[342]; int neighborPositions[342];
...@@ -675,7 +683,7 @@ protected: ...@@ -675,7 +683,7 @@ protected:
for(int idxY = -2 ; idxY <= 3 ; ++idxY){ for(int idxY = -2 ; idxY <= 3 ; ++idxY){
for(int idxZ = -2 ; idxZ <= 3 ; ++idxZ){ for(int idxZ = -2 ; idxZ <= 3 ; ++idxZ){
if( FMath::Abs(idxX) > 1 || FMath::Abs(idxY) > 1 || FMath::Abs(idxZ) > 1){ if( FMath::Abs(idxX) > 1 || FMath::Abs(idxY) > 1 || FMath::Abs(idxZ) > 1){
neighbors[counter] = &upperCells[idxUpperLevel-1]; neighbors[counter] = &upperCells[idxUpperLevel-1];
neighborPositions[counter] = neighIndex(idxX,idxY,idxZ); neighborPositions[counter] = neighIndex(idxX,idxY,idxZ);
++counter; ++counter;
} }
...@@ -687,6 +695,7 @@ protected: ...@@ -687,6 +695,7 @@ protected:
= &(downerCells[idxUpperLevel-1].getLocalExpansionData()); = &(downerCells[idxUpperLevel-1].getLocalExpansionData());
const symbolic_data_t* const target_symbolic const symbolic_data_t* const target_symbolic
= &downerCells[idxUpperLevel-1]; = &downerCells[idxUpperLevel-1];
std::array<const multipole_t*, 342> neighbor_multipoles; std::array<const multipole_t*, 342> neighbor_multipoles;
std::array<const symbolic_data_t*, 342> neighbor_symbolics; std::array<const symbolic_data_t*, 342> neighbor_symbolics;
std::transform(neighbors, neighbors+counter, neighbor_multipoles.begin(), std::transform(neighbors, neighbors+counter, neighbor_multipoles.begin(),
...@@ -706,10 +715,11 @@ protected: ...@@ -706,10 +715,11 @@ protected:
neighborPositions, neighborPositions,
counter); counter);
} }
// Run the L2L for all but the lowest virtual levels // Run the L2L for all but the lowest virtual levels
// 2 is the highest level
// offsetRealTree-1 is the root level
{ {
std::array<local_expansion_t*, 8> virtual_child_local_exps = {}; std::array<local_expansion_t*, 8> virtual_child_local_exps = {};
std::array<symbolic_data_t*, 8> virtual_child_symbolics = {}; std::array<symbolic_data_t*, 8> virtual_child_symbolics = {};
for(int idxLevel = 2 ; idxLevel < offsetRealTree-1 ; ++idxLevel){ for(int idxLevel = 2 ; idxLevel < offsetRealTree-1 ; ++idxLevel){
...@@ -732,16 +742,17 @@ protected: ...@@ -732,16 +742,17 @@ protected:
} }
} }
// Run the L2L for the lowest virtual level // Run the L2L for the lowest virtual level
{ if(offsetRealTree >2){
std::array<local_expansion_t*, 8> virtual_child_local_exps = {}; std::array<local_expansion_t*, 8> virtual_child_local_exps {};
std::array<symbolic_data_t*, 8> virtual_child_symbolics = {}; std::array<symbolic_data_t*, 8> virtual_child_symbolics {};
const int idxLevel = offsetRealTree-1; const int idxLevel = offsetRealTree-1;
local_expansion_t* const virtual_parent_local_exp local_expansion_t* const virtual_parent_local_exp
= &(downerCells[idxLevel-1].getLocalExpansionData()); = &(downerCells[idxLevel-1].getLocalExpansionData());
const symbolic_data_t* const virtual_parent_symbolic const symbolic_data_t* const virtual_parent_symbolic
= &downerCells[idxLevel-1]; = &downerCells[idxLevel-1];
virtual_child_local_exps[7] = &(downerCells[idxLevel].getLocalExpansionData()); virtual_child_local_exps[7] = &(downerCells[idxLevel].getLocalExpansionData());
virtual_child_symbolics[7] = &downerCells[idxLevel]; virtual_child_symbolics[7] = &downerCells[idxLevel];
kernels->L2L( kernels->L2L(
virtual_parent_local_exp, virtual_parent_local_exp,
virtual_parent_symbolic, virtual_parent_symbolic,
...@@ -750,7 +761,7 @@ protected: ...@@ -750,7 +761,7 @@ protected:
); );
} }
// Run L2L from the lowest vitual level to the highest real tree level // Run L2L from the lowest virtual level to the highest real tree level
{ {
typename OctreeClass::Iterator octreeIterator(tree); typename OctreeClass::Iterator octreeIterator(tree);
octreeIterator.gotoLeft(); octreeIterator.gotoLeft();
...@@ -777,7 +788,6 @@ protected: ...@@ -777,7 +788,6 @@ protected:
child_symbolics.data() child_symbolics.data()
); );
} }
delete[] upperCells; delete[] upperCells;
delete[] downerCells; delete[] downerCells;
} }
......
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