Commit 4ec8fdd1 authored by berenger-bramas's avatar berenger-bramas

Create a arranger to work in parallel.

Now the tree is post processed and the mpi process exchange the particles that change
not only in leaf but also in proc host.

git-svn-id: svn+ssh://scm.gforge.inria.fr/svn/scalfmm/scalfmm/trunk@264 2616d619-271b-44dc-8df4-d4a8f33a7222
parent d1e9111a
#ifndef FOCTREEARRANGERPROC_HPP
#define FOCTREEARRANGERPROC_HPP
#include "../Utils/FGlobal.hpp"
#include "../Containers/FVector.hpp"
#include "../Utils/FAssertable.hpp"
#include "../Utils/FMpi.hpp"
template <class OctreeClass, class ContainerClass, class ParticleClass>
class FOctreeArrangerProc : FAssertable {
struct Interval{
MortonIndex min;
MortonIndex max;
};
static void mpiassert(const int test, const unsigned line, const char* const message = 0){
if(test != MPI_SUCCESS){
printf("[ERROR] Test failled at line %d, result is %d", line, test);
if(message) printf(", message: %s",message);
printf("\n");
fflush(stdout);
MPI_Abort(MPI_COMM_WORLD, int(line) );
}
}
int getInterval(const MortonIndex mindex, const int size, const Interval intervals[]) const{
for(int idxProc = 0 ; idxProc < size ; ++idxProc){
if( intervals[idxProc].min <= mindex && mindex < intervals[idxProc].max){
return idxProc;
}
}
return size - 1;
}
OctreeClass* const tree;
public:
FOctreeArrangerProc(OctreeClass* const inTree) : tree(inTree) {
fassert(tree, "Tree cannot be null", __LINE__ , __FILE__ );
}
/** return false if tree is empty after processing */
bool rearrange(const FMpi::FComm& comm){
Interval*const intervals = new Interval[comm.processCount()];
memset(intervals, 0, sizeof(Interval) * comm.processCount());
{ // We need to give an interval to each process, this interval
// will be based on the current interval
Interval myLastInterval;
typename OctreeClass::Iterator octreeIterator(tree);
octreeIterator.gotoBottomLeft();
myLastInterval.min = octreeIterator.getCurrentGlobalIndex();
octreeIterator.gotoRight();
myLastInterval.max = octreeIterator.getCurrentGlobalIndex();
// We get the min/max indexes from each procs
mpiassert( MPI_Allgather( &myLastInterval, sizeof(Interval), MPI_BYTE, intervals, sizeof(Interval), MPI_BYTE, comm.getComm()), __LINE__ );
intervals[0].min = 0;
for(int idxProc = 1 ; idxProc < comm.processCount() ; ++idxProc){
intervals[idxProc].min = ((intervals[idxProc].min - intervals[idxProc-1].max)/2) + intervals[idxProc-1].max;
intervals[idxProc-1].max = intervals[idxProc].min;
}
intervals[comm.processCount() - 1].max = ((1 << (3*(tree->getHeight()-1))) - 1);
}
FVector<ParticleClass>*const toMove = new FVector<ParticleClass>[comm.processCount()];
{ // iterate on the leafs and found particle to remove or to send
typename OctreeClass::Iterator octreeIterator(tree);
octreeIterator.gotoBottomLeft();
do{
const MortonIndex currentIndex = octreeIterator.getCurrentGlobalIndex();
typename ContainerClass::BasicIterator iter(*octreeIterator.getCurrentListTargets());
while( iter.hasNotFinished() ){
const MortonIndex particuleIndex = tree->getMortonFromPosition(iter.data().getPosition());
// is this particle need to be changed from its leaf
if(particuleIndex != currentIndex){
// find the right interval
const int procConcerned = getInterval( particuleIndex, comm.processCount(), intervals);
toMove[procConcerned].push(iter.data());
iter.remove();
}
else {
iter.gotoNext();
}
}
} while(octreeIterator.moveRight());
}
// To send and recv
ParticleClass* toReceive = 0;
MPI_Request*const requests = new MPI_Request[comm.processCount()*2];
memset(requests, 0, sizeof(MPI_Request) * comm.processCount() * 2);
long long int*const indexToReceive = new long long int[comm.processCount() + 1];
memset(indexToReceive, 0, sizeof(long long int) * comm.processCount() + 1);
int iterRequests = 0;
int limitRecvSend = 0;
int hasToRecvFrom = 0;
{ // gather what to send to who + isend data
int*const counter = new int[comm.processCount()];
memset(counter, 0, sizeof(int) * comm.processCount());
for(int idxProc = 0 ; idxProc < comm.processCount() ; ++idxProc){
counter[idxProc] = toMove[idxProc].getSize();
}
// say who send to who
int*const allcounter = new int[comm.processCount()*comm.processCount()];
mpiassert( MPI_Allgather( counter, comm.processCount(), MPI_INT, allcounter, comm.processCount(), MPI_INT, comm.getComm()), __LINE__ );
// prepare buffer to receive
long long int sumToRecv = 0;
indexToReceive[0] = 0;
for(int idxProc = 0 ; idxProc < comm.processCount() ; ++idxProc){
if( idxProc != comm.processId()){
sumToRecv += allcounter[idxProc * comm.processCount() + comm.processId()];
}
indexToReceive[idxProc + 1] = sumToRecv;
}
toReceive = new ParticleClass[sumToRecv];
// send
for(int idxProc = 0 ; idxProc < comm.processCount() ; ++idxProc){
if(idxProc != comm.processId() && allcounter[idxProc * comm.processCount() + comm.processId()]){
mpiassert( MPI_Irecv(&toReceive[indexToReceive[idxProc]], allcounter[idxProc * comm.processCount() + comm.processId()] * sizeof(ParticleClass), MPI_BYTE,
idxProc, 0, comm.getComm(), &requests[iterRequests++]), __LINE__ );
hasToRecvFrom += 1;
}
}
limitRecvSend = iterRequests;
// recv
for(int idxProc = 0 ; idxProc < comm.processCount() ; ++idxProc){
if(idxProc != comm.processId() && toMove[idxProc].getSize()){
mpiassert( MPI_Isend(toMove[idxProc].data(), toMove[idxProc].getSize() * sizeof(ParticleClass), MPI_BYTE,
idxProc, 0, comm.getComm(), &requests[iterRequests++]), __LINE__ );
}
}
delete[] allcounter;
delete[] counter;
}
{ // insert particles that moved
for(int idxPart = 0 ; idxPart < toMove[comm.processId()].getSize() ; ++idxPart){
tree->insert(toMove[comm.processId()][idxPart]);
}
}
{ // wait any recv or send
// if it is a recv then insert particles
MPI_Status status;
while( hasToRecvFrom ){
int done = 0;
mpiassert( MPI_Waitany( iterRequests, requests, &done, &status ), __LINE__ );
if( done < limitRecvSend ){
const int source = status.MPI_SOURCE;
for(long long int idxPart = indexToReceive[source] ; idxPart < indexToReceive[source+1] ; ++idxPart){
tree->insert(toReceive[idxPart]);
}
hasToRecvFrom -= 1;
}
}
}
int counterLeavesAlive = 0;
{ // Remove empty leaves
typename OctreeClass::Iterator octreeIterator(tree);
octreeIterator.gotoBottomLeft();
bool workOnNext = true;
do{
// Empty leaf
if( octreeIterator.getCurrentListTargets()->getSize() == 0 ){
const MortonIndex currentIndex = octreeIterator.getCurrentGlobalIndex();
workOnNext = octreeIterator.moveRight();
tree->removeLeaf( currentIndex );
}
// Not empty, just continue
else {
// todo remove
if( octreeIterator.getCurrentGlobalIndex() < intervals[comm.processId()].min
|| intervals[comm.processId()].max <= octreeIterator.getCurrentGlobalIndex()){
std::cout << comm.processId() << " Must delete this leaf at " << octreeIterator.getCurrentGlobalIndex()
<< " size " << octreeIterator.getCurrentListTargets()->getSize() <<std::endl;
}
workOnNext = octreeIterator.moveRight();
counterLeavesAlive += 1;
}
} while( workOnNext );
}
// wait all send
mpiassert( MPI_Waitall( iterRequests, requests, MPI_STATUSES_IGNORE), __LINE__ );
delete[] intervals;
delete[] toMove;
delete[] requests;
delete[] toReceive;
delete[] indexToReceive;
// return false if tree is empty
return counterLeavesAlive != 0;
}
};
#endif // FOCTREEARRANGERPROC_HPP
......@@ -87,6 +87,41 @@ private:
}
}
static void SortParticlesFromArray( const FMpi::FComm& communicator, const ParticleClass array[], const int size, const SortingType type,
const F3DPosition& centerOfBox, const FReal boxWidth,
const int TreeHeight, IndexedParticle**const outputArray, FSize* const outputSize){
FTRACE( FTrace::FFunction functionTrace(__FUNCTION__ , "Loader to Tree" , __FILE__ , __LINE__) );
// create particles
IndexedParticle*const realParticlesIndexed = new IndexedParticle[size];
FMemUtils::memset(realParticlesIndexed, 0, sizeof(IndexedParticle) * size);
F3DPosition boxCorner(centerOfBox - (boxWidth/2));
FTreeCoordinate host;
const FReal boxWidthAtLeafLevel = boxWidth / FReal(1 << (TreeHeight - 1) );
for(long idxPart = 0 ; idxPart < size ; ++idxPart){
realParticlesIndexed[idxPart].particle = array[idxPart];
host.setX( getTreeCoordinate( realParticlesIndexed[idxPart].particle.getPosition().getX() - boxCorner.getX(), boxWidthAtLeafLevel ));
host.setY( getTreeCoordinate( realParticlesIndexed[idxPart].particle.getPosition().getY() - boxCorner.getY(), boxWidthAtLeafLevel ));
host.setZ( getTreeCoordinate( realParticlesIndexed[idxPart].particle.getPosition().getZ() - boxCorner.getZ(), boxWidthAtLeafLevel ));
realParticlesIndexed[idxPart].index = host.getMortonIndex(TreeHeight - 1);
}
// sort particles
if(type == QuickSort){
FQuickSort<IndexedParticle,MortonIndex, FSize>::QsMpi(realParticlesIndexed, size, *outputArray, *outputSize,communicator);
delete [] (realParticlesIndexed);
}
else {
FBitonicSort<IndexedParticle,MortonIndex, FSize>::Sort( realParticlesIndexed, size, communicator );
*outputArray = realParticlesIndexed;
*outputSize = size;
}
}
//////////////////////////////////////////////////////////////////////////
// To merge the leaves
......@@ -465,6 +500,22 @@ public:
EqualizeAndFillTree(communicator, realTree, leavesArray, leavesSize);
}
template <class OctreeClass>
static void ArrayToTree(const FMpi::FComm& communicator, const ParticleClass array[], const int size,
const F3DPosition& boxCenter, const FReal boxWidth,
OctreeClass& realTree, const SortingType type = QuickSort){
IndexedParticle* particlesArray = 0;
FSize particlesSize = 0;
SortParticlesFromArray(communicator, array, size, type, boxCenter, boxWidth, realTree.getHeight(), &particlesArray, &particlesSize);
char* leavesArray = 0;
FSize leavesSize = 0;
MergeLeaves(communicator, particlesArray, particlesSize, &leavesArray, &leavesSize);
EqualizeAndFillTree(communicator, realTree, leavesArray, leavesSize);
}
};
......
......@@ -229,6 +229,10 @@ public:
return MPI_LONG_LONG;
}
static MPI_Datatype GetType(long int&){
return MPI_LONG;
}
static MPI_Datatype GetType(double&){
return MPI_DOUBLE;
}
......
// /!\ Please, you must read the license at the bottom of this page
#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include "../Src/Utils/FParameters.hpp"
#include "../Src/Utils/FTic.hpp"
#include "../Src/Containers/FOctree.hpp"
#include "../Src/Containers/FVector.hpp"
#include "../Src/Components/FSimpleLeaf.hpp"
#include "../Src/Utils/F3DPosition.hpp"
#include "../Src/Components/FTestParticle.hpp"
#include "../Src/Components/FTestCell.hpp"
#include "../Src/Arranger/FOctreeArrangerProc.hpp"
#include "../Src/Files/FMpiTreeBuilder.hpp"
// Simply create particles and try the kernels
int main(int argc, char ** argv){
typedef FTestParticle ParticleClass;
typedef FTestCell CellClass;
typedef FVector<ParticleClass> ContainerClass;
typedef FSimpleLeaf<ParticleClass, ContainerClass > LeafClass;
typedef FOctree<ParticleClass, CellClass, ContainerClass , LeafClass > OctreeClass;
///////////////////////What we do/////////////////////////////
std::cout << ">> This executable has to be used to test the FMM algorithm.\n";
//////////////////////////////////////////////////////////////
FMpi app(argc, argv);
const int NbLevels = FParameters::getValue(argc,argv,"-h", 7);
const int SizeSubLevels = FParameters::getValue(argc,argv,"-sh", 3);
const long NbPart = FParameters::getValue(argc,argv,"-nb", 20000);
const FReal FRandMax = FReal(RAND_MAX);
FTic counter;
srand ( 1 ); // volontary set seed to constant
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
const FReal BoxWidth = 1.0;
const FReal BoxCenter = 0.5;
OctreeClass tree(NbLevels, SizeSubLevels, BoxWidth, F3DPosition(BoxCenter,BoxCenter,BoxCenter));
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
std::cout << "Creating & Inserting " << NbPart << " particles ..." << std::endl;
counter.tic();
{
FTestParticle particles[NbPart];
for(int idxPart = 0 ; idxPart < NbPart ; ++idxPart){
particles[idxPart].setPosition(
(BoxWidth*FReal(rand())/FRandMax) + (BoxCenter-(BoxWidth/2)),
(BoxWidth*FReal(rand())/FRandMax) + (BoxCenter-(BoxWidth/2)),
(BoxWidth*FReal(rand())/FRandMax) + (BoxCenter-(BoxWidth/2)));
}
FMpiTreeBuilder<ParticleClass>::ArrayToTree(app.global(), particles, NbPart, F3DPosition(BoxCenter,BoxCenter,BoxCenter), BoxWidth, tree);
}
counter.tac();
std::cout << "Done " << "(@Creating and Inserting Particles = " << counter.elapsed() << "s)." << std::endl;
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
std::cout << "Working on particles ..." << std::endl;
counter.tic();
{ // Check that each particle has been summed with all other
typename OctreeClass::Iterator octreeIterator(&tree);
octreeIterator.gotoBottomLeft();
do{
typename ContainerClass::BasicIterator iter(*octreeIterator.getCurrentListTargets());
while( iter.hasNotFinished() ){
if(FReal(rand())/FRandMax > 0.5){
iter.data().setPosition(
(BoxWidth*FReal(rand())/FRandMax) + (BoxCenter-(BoxWidth/2)),
(BoxWidth*FReal(rand())/FRandMax) + (BoxCenter-(BoxWidth/2)),
(BoxWidth*FReal(rand())/FRandMax) + (BoxCenter-(BoxWidth/2)));
}
iter.gotoNext();
}
} while(octreeIterator.moveRight());
}
counter.tac();
std::cout << "Done " << "(@Moving = " << counter.elapsed() << "s)." << std::endl;
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
std::cout << "Arrange ..." << std::endl;
counter.tic();
FOctreeArrangerProc<OctreeClass, ContainerClass, ParticleClass> arrange(&tree);
arrange.rearrange(app.global());
counter.tac();
std::cout << "Done " << "(@Arrange = " << counter.elapsed() << "s)." << std::endl;
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
std::cout << "Test ..." << std::endl;
counter.tic();
{ // Check that each particle has been put into the right leaf
long counterPart = 0;
typename OctreeClass::Iterator octreeIterator(&tree);
octreeIterator.gotoBottomLeft();
MortonIndex interval[2] = {0,0};
interval[0] = octreeIterator.getCurrentGlobalIndex();
do{
const MortonIndex leafIndex = octreeIterator.getCurrentGlobalIndex();
typename ContainerClass::BasicIterator iter(*octreeIterator.getCurrentListTargets());
while( iter.hasNotFinished() ){
const MortonIndex particleIndex = tree.getMortonFromPosition( iter.data().getPosition() );
if( leafIndex != particleIndex){
std::cout << "Index problem, should be " << leafIndex <<
" particleIndex "<< particleIndex << std::endl;
}
iter.gotoNext();
}
counterPart += octreeIterator.getCurrentListTargets()->getSize();
if(octreeIterator.getCurrentListTargets()->getSize() == 0){
std::cout << "Problem, leaf is empty at index " << leafIndex << std::endl;
}
} while(octreeIterator.moveRight());
interval[1] = octreeIterator.getCurrentGlobalIndex();
counterPart = app.global().reduceSum(counterPart);
if(app.global().processId() == 0 && counterPart != NbPart * app.global().processCount() ){
std::cout <<"Wrong particles number, should be " << (NbPart * app.global().processCount()) << " but is " << counterPart << std::endl;
}
MortonIndex*const allintervals = new MortonIndex[ 2 * app.global().processCount() ];
MPI_Allgather( interval, sizeof(MortonIndex) * 2, MPI_BYTE, allintervals, sizeof(MortonIndex) * 2, MPI_BYTE, MPI_COMM_WORLD);
if(app.global().processId() == 0){
for(int idxProc = 1 ; idxProc < app.global().processCount() ; ++idxProc){
if( allintervals[idxProc*2-1] > allintervals[idxProc*2] ){
std::cout << "Interval problem for [" << idxProc-1 << "].max = " << allintervals[idxProc*2-1] <<
" [" << idxProc << "].min = "<< allintervals[idxProc*2] << std::endl;
}
}
}
delete[] allintervals;
}
{ // Check that each particle has been summed with all other
typename OctreeClass::Iterator octreeIterator(&tree);
typename OctreeClass::Iterator avoidGotoLeftIterator(octreeIterator);
const int heightMinusOne = NbLevels - 1;
for(int idxLevel = 1 ; idxLevel < heightMinusOne ; ++idxLevel ){
// for each cells
do{
int countChild = 0;
CellClass** const child = octreeIterator.getCurrentChild();
for(int idxChild = 0 ; idxChild < 8 ; ++idxChild ){
if( child[idxChild] ){
countChild += 1;
}
}
if(countChild == 0){
std::cout << "Problem at level " << idxLevel << " cell has no child " << octreeIterator.getCurrentGlobalIndex() << std::endl;
}
} while(octreeIterator.moveRight());
avoidGotoLeftIterator.moveDown();
octreeIterator = avoidGotoLeftIterator;
}
}
counter.tac();
std::cout << "Done " << "(@Test = " << counter.elapsed() << "s)." << std::endl;
//////////////////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////////////////
return 0;
}
// [--LICENSE--]
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