Commit ab595023 authored by DUFOYER Benjamin's avatar DUFOYER Benjamin
Browse files

Add Creation of GroupLetTree

parent 824bbc34
......@@ -7,7 +7,27 @@
#define _FDISTRIBUTED_GROUPTREE_BUILDER_HPP_
#include "FGroupTree.hpp"
#include "FOutOfBlockInteraction.hpp"
#include <cmath>
#include <algorithm>
#include <stdint.h>
#include <limits.h>
#if SIZE_MAX == UCHAR_MAX
#define my_MPI_SIZE_T MPI_UNSIGNED_CHAR
#elif SIZE_MAX == USHRT_MAX
#define my_MPI_SIZE_T MPI_UNSIGNED_SHORT
#elif SIZE_MAX == UINT_MAX
#define my_MPI_SIZE_T MPI_UNSIGNED
#elif SIZE_MAX == ULONG_MAX
#define my_MPI_SIZE_T MPI_UNSIGNED_LONG
#elif SIZE_MAX == ULLONG_MAX
#define my_MPI_SIZE_T MPI_UNSIGNED_LONG_LONG
#else
#error "what is happening here?"
#endif
namespace dstr_grp_tree_builder{
......@@ -155,10 +175,8 @@ void parrallel_build_block(const inria::mpi_config& conf, std::vector<node_t>* l
// We declare the new linear_tree who have the new number of cell
std::vector<node_t>* new_linear_tree = new std::vector<node_t>(array_global_nb_leaf[my_rank]);
//Posting receiv message
//For every interaction
// Array of request to know the status
inria::mpi::request tab_mpi_status[interaction_recev.size()];
for(unsigned i = 0 ; i < interaction_recev.size(); i++ ){
......@@ -185,6 +203,7 @@ void parrallel_build_block(const inria::mpi_config& conf, std::vector<node_t>* l
if(my_rank < interaction_send.at(i).process_rank){
start = (unsigned)linear_tree->size() - interaction_send.at(i).nb_leaf;
}
//sending leaf
conf.comm.isend(&linear_tree->at(start),
size_send,
......@@ -196,16 +215,824 @@ void parrallel_build_block(const inria::mpi_config& conf, std::vector<node_t>* l
// The function need to know how many cell the current process send on
// the right and on the left because the MPI request write on the same
// pointer
fill_new_linear_tree(linear_tree,new_linear_tree,nb_leaf_recev_left,nb_leaf_recev_right,nb_leaf_send_left,nb_leaf_send_right);
fill_new_linear_tree(linear_tree,
new_linear_tree,
nb_leaf_recev_left,
nb_leaf_recev_right,
nb_leaf_send_left,
nb_leaf_send_right);
// waiting for the end of MPI request
inria::mpi::request::waitall(interaction_recev.size(),tab_mpi_status);
free(array_global_nb_leaf);
// swaping linear_tree pointer
std::swap(*linear_tree,*new_linear_tree);
}
/**
* The idea of this function is to know the repartition of particle
* morton index on every proc
* Every proc send he first and his last particle morton index with a MPI
* allgather
* @author benjamin.dufoyer@inria.fr
* @param conf MPI
* @param particle : vector where particles a stock, they will be sorted
* BEFORE calling function
* @param particle_repartition : a empty vector with the size of the
* number of process
*/
template<class particle_t,
class type1_t,
class type2_t>
void know_particle_division(
const inria::mpi_config& conf,
std::vector<particle_t>& particle,
std::vector<std::pair<type1_t,type2_t>>& particle_repartition)
{
FAssert(particle_repartition.size() == (unsigned)conf.comm.size());
std::pair<type1_t,type2_t> my_idx;
my_idx.first = particle.front().morton_index;
my_idx.second = particle.back().morton_index;
// Distribute the local min max morton_index to every process
conf.comm.allgather(&my_idx,
sizeof(std::pair<type1_t,type2_t>),
MPI_CHAR,
particle_repartition.data(),
sizeof(std::pair<type1_t,type2_t>),
MPI_CHAR);
}
/**
* This function compute the morton index of every leaf needed for the P2P
* First we compute every morton index needed for every leaf
* We sort the result
* And to finish we remove multiple iteration
* @author benjamin.dufoyer@inria.fr
* @param GroupOctreeClass Actual group tree
* @param global_min_m_idx min morton index of the simulation
* @param global_max_m_idx max morton index of the simulation
* @return std::vector with morton_index of cell who is needed to compute P2P
*/
template<class GroupOctreeClass>
std::vector<size_t> get_leaf_P2P_interaction(
GroupOctreeClass& tree,
const MortonIndex& global_min_m_idx,
const MortonIndex& global_max_m_idx,
const MortonIndex& local_min_m_idx,
const MortonIndex& local_max_m_idx
){
// 26 is for every interaction
std::vector<size_t> externalInteractionsLeafLevel(tree.getTotalNbLeaf()*26,0);
// Reset interactions
// idx to know where we are in the vector
size_t idx_vector= 0;
// First leaf level
{
// We iterate on every particle group
for(int idxGroup = 0 ; idxGroup < tree.getNbParticleGroup() ; ++idxGroup){
// get the particle group
// it's ugly but, if i use template, it's not convenient
auto* containers = tree.getParticleGroup(idxGroup);
{
// Iterate on every leaf
for(int leafIdx = 0;
leafIdx < containers->getNumberOfLeavesInBlock();
++leafIdx){
// Getting the morton index of the leaf
const MortonIndex mindex = containers->getLeafMortonIndex(leafIdx);
// Define array to receive computed information
MortonIndex interactionsIndexes[26];
int interactionsPosition[26];
//Getting coordinate of the current leaf
FTreeCoordinate coord(mindex);
// Getting neigbors
int counter = coord.getNeighborsIndexes(tree.getHeight(),interactionsIndexes,interactionsPosition);
// Iterate on every neighbors
for(int idxInter = 0 ; idxInter < counter ; ++idxInter){
// Check if the current proc already have the leaf
if(interactionsIndexes[idxInter] < local_min_m_idx || interactionsIndexes[idxInter] > local_max_m_idx ){
// Check if the leaf exist
if(interactionsIndexes[idxInter] >= global_min_m_idx && interactionsIndexes[idxInter] <= global_max_m_idx ){
externalInteractionsLeafLevel[idx_vector] = interactionsIndexes[idxInter];
++idx_vector;
}
}
}
}
}
}
}
if(idx_vector != 0) {
// Sort every morton index
FQuickSort<size_t>::QsSequential(externalInteractionsLeafLevel.data(),idx_vector);
// Compute the number of different morton index
size_t nb_leaf = 1;
size_t last_leaf = externalInteractionsLeafLevel[0];
for(unsigned i = 1 ; i < idx_vector ; ++i){
if(last_leaf != externalInteractionsLeafLevel[i]){
last_leaf = externalInteractionsLeafLevel[i];
nb_leaf++;
}
}
// Alloc the returned vector
std::vector<size_t> leaf_needed(nb_leaf,0);
// Fill the returned vector
size_t current_idx = 1;
last_leaf = externalInteractionsLeafLevel[0];
leaf_needed[0] = externalInteractionsLeafLevel[0];
for(unsigned i = 1 ; i < idx_vector ; ++i){
if(last_leaf != externalInteractionsLeafLevel[i]){
last_leaf = externalInteractionsLeafLevel[i];
leaf_needed[current_idx] = externalInteractionsLeafLevel[i];
++current_idx;
}
if(current_idx == nb_leaf)
break;
}
return leaf_needed;
} else {
std::vector<size_t> leaf_needed(0,0);
return leaf_needed;
}
}
/**
* This function compute the leaf needed for the M2L operation
* We take every leaf of the tree, get her parent, get the neigbors of
* the parents and take every child of the parent's neighbors.
* After we sort the result and remove duplicate data
* @author benjamin.dufoyer@inria.fr
* @param global_min_m_idx global min morton index
* @param global_max_m_idx global max morton index
* @param local_min_m_idx local min morton index
* @param local_max_m_idx local max morton index
* @param tree local group tree
* @param dim optionnal parameter to compute in other dimension
* @return vector with leaf needed for the M2M operation
*/
template<class GroupOctreeClass>
std::vector<size_t> get_leaf_M2L_interaction(
const MortonIndex& global_min_m_idx,
const MortonIndex& global_max_m_idx,
const MortonIndex& local_min_m_idx,
const MortonIndex& local_max_m_idx,
GroupOctreeClass& tree,
int dim = 3)
{
size_t idx_vector = 0;
std::vector<size_t> external_leaf(tree.getTotalNbLeaf()*26*8,0);
// iterate on the group
for(int idxGroup = 0 ; idxGroup < tree.getNbParticleGroup() ; ++idxGroup){
auto* containers = tree.getParticleGroup(idxGroup);
size_t curr_m_idx;
// +1 is to pass the test at the first try
size_t last_m_idx = (containers->getLeafMortonIndex(0) >> dim)+1;
for(int leafIdx = 0;
leafIdx < containers->getNumberOfLeavesInBlock();
++leafIdx){
// Getting the current morton index
curr_m_idx = containers->getLeafMortonIndex(leafIdx);
// Compute the morton index of the father
curr_m_idx = curr_m_idx >> dim;
// If it's a new father
if(curr_m_idx != last_m_idx){
last_m_idx = curr_m_idx;
// Compute coordinate
MortonIndex interactionsIndexes[26];
int interactionsPosition[26];
FTreeCoordinate coord(curr_m_idx);
// Getting neigbors of the father
int counter = coord.getNeighborsIndexes(tree.getHeight(),interactionsIndexes,interactionsPosition);
for(int idxNeighbor = 0 ; idxNeighbor < counter ; ++idxNeighbor){
// Compute the first and the last child
size_t first_child =
interactionsIndexes[idxNeighbor] << dim;
size_t last_child =
((interactionsIndexes[idxNeighbor]+1) << dim)-1;
// Add child if needed
for(size_t idx_Child = first_child ; idx_Child <= last_child ; ++idx_Child){
if( idx_Child >= (size_t)global_min_m_idx
&& idx_Child <= (size_t)global_max_m_idx)
{
if(idx_Child < (size_t)local_min_m_idx || idx_Child > (size_t)local_max_m_idx){
//Stock the leaf
external_leaf[idx_vector] = idx_Child;
++idx_vector;
}
}
} // end for child
} // end for neigbors
} // end if
} // end for leaf
} // end for group
// if we have leaf in the vector
if(idx_vector != 0 ){
//sort the result
FQuickSort<size_t>::QsSequential(external_leaf.data(),idx_vector);
// compute the number of leaf
size_t nb_leaf = 1;
size_t last_leaf = external_leaf[0];
for(unsigned i = 1 ; i < idx_vector ; ++i){
if(last_leaf != external_leaf[i]){
last_leaf = external_leaf[i];
nb_leaf++;
}
}
// Alloc the returned vector
std::vector<size_t> leaf_needed(nb_leaf,0);
// Fill the returned vector
size_t current_idx = 1;
last_leaf = external_leaf[0];
leaf_needed[0] = external_leaf[0];
for(unsigned i = 1 ; i < idx_vector ; ++i){
if(last_leaf != external_leaf[i]){
last_leaf = external_leaf[i];
leaf_needed[current_idx] = external_leaf[i];
++current_idx;
}
if(current_idx == nb_leaf)
break;
}
return leaf_needed;
} else {
std::vector<size_t> leaf_needed(0,0);
return leaf_needed;
}
}
/**
* The goal of this function is to concat the two vector
* and erase duplicate data
* first we compute the number of leaf
* second we fill the vector with erase duplicate data
* @author benjamin.dufoyer@inria.fr
* @param leaf_P2P vector of P2P leaf
* @param leaf_M2L vector of M2L leaf
* @return vector with morton index of leaf needed
*/
std::vector<size_t> concat_M2L_P2P(
std::vector<size_t>& leaf_P2P,
std::vector<size_t>& leaf_M2L)
{
// Checking if one of vector is empty
if(leaf_P2P.size() == 0)
return leaf_M2L;
if(leaf_M2L.size() == 0)
return leaf_P2P;
// Compute the total number of leaf needed
unsigned idx_P2P = 0;
unsigned idx_M2L = 0;
size_t nb_leaf = 0;
while(idx_P2P != leaf_P2P.size() && idx_M2L != leaf_M2L.size()){
size_t curr_P2P = leaf_P2P[idx_P2P];
size_t curr_M2M = leaf_M2L[idx_M2L];
if(curr_P2P != curr_M2M){
if(curr_P2P < curr_M2M){
++idx_P2P;
} else {
++idx_M2L;
}
} else {
++idx_P2P;
++idx_M2L;
}
++nb_leaf;
}
if(idx_P2P == leaf_P2P.size()){
nb_leaf += (leaf_M2L.size()-1) - idx_M2L;
} else {
nb_leaf += (leaf_P2P.size()-1) - idx_P2P;
}
// Allocate the vector
std::vector<size_t> leaf_needed(nb_leaf,0);
idx_P2P = 0;
idx_M2L = 0;
size_t idx_leaf = 0;
// fill the vector
while(idx_P2P != leaf_P2P.size() && idx_M2L != leaf_M2L.size()){
size_t curr_P2P = leaf_P2P[idx_P2P];
size_t curr_M2M = leaf_M2L[idx_M2L];
if(curr_P2P != curr_M2M){
if(curr_P2P < curr_M2M){
leaf_needed[idx_leaf] = leaf_P2P[idx_P2P];
++idx_P2P;
} else {
leaf_needed[idx_leaf] = leaf_M2L[idx_M2L];
++idx_M2L;
}
++nb_leaf;
} else {
leaf_needed[idx_leaf] = leaf_M2L[idx_M2L];
++idx_P2P;
++idx_M2L;
}
idx_leaf++;
}
// Copy the rest of leaf with a memcpy
void* destination = &leaf_needed.data()[idx_leaf];
void* source;
size_t num;
if(idx_P2P == leaf_P2P.size()){
source = &leaf_M2L[idx_M2L];
num = sizeof(size_t)* ((leaf_M2L.size()-1) - idx_M2L);
} else {
source = &leaf_P2P[idx_P2P];
num = sizeof(size_t)* ((leaf_P2P.size()-1) - idx_P2P);
}
memcpy(destination,source,num);
return leaf_needed;
}
/**
* This function compute the number of leaf needed on every proc,
* she check if leaf exit, if a leaf doesn't exist this function
* remove here at the end of the function
* @author benjamin.dufoyer@inria.fr
* @param needed_leaf vector with the morton index of every leaf
* needed for my proc
* @param particle_repartition vector with the repartation of particle on
* every proc
* @param conf conf MPI
* @return vector<size_t> with a size nb_proc*nb_proc
* it's a interaction matrix
*/
std::vector<size_t> get_matrix_interaction(
std::vector<size_t>& needed_leaf,
std::vector<std::pair<unsigned long long int,unsigned long long int>>& particle_repartition,
const inria::mpi_config& conf)
{
// Getting MPI Info
const int nb_proc = conf.comm.size();
// Alloc interaction matrix
std::vector<size_t> my_matrix_interaction(nb_proc,0);
std::vector<size_t> global_matrix_interaction(nb_proc*nb_proc,0);
// Initialise idx on particle_repartition
size_t idx_part = 0;
// Interate on every leaf to know where she is
for(unsigned idx_leaf = 0; idx_leaf < needed_leaf.size(); ++idx_leaf){
size_t current_leaf = needed_leaf[idx_leaf];
// if she is on the current proc
if(current_leaf >= particle_repartition[idx_part].first
&& current_leaf <= particle_repartition[idx_part].second){
my_matrix_interaction[idx_part] += 1;
} else {
// While the current leaf is not on the good interval
while(particle_repartition[idx_part].second < current_leaf){
idx_part += 1;
}
if(particle_repartition[idx_part].first > current_leaf){
// in this case the leaf is not in interval, so she doesn't exist
needed_leaf[idx_leaf] = 0;
} else {
// In the case it's a normal case, we juste increment the
// number of leaf send at the proc idx_part
my_matrix_interaction[idx_part] += 1;
}
}
}
// now we have the number of leaf to send at every proc
// we proceed a AllGather to share this information at every proc
conf.comm.allgather(my_matrix_interaction.data(),
nb_proc,
my_MPI_SIZE_T,
global_matrix_interaction.data(),
nb_proc,
my_MPI_SIZE_T);
// removing bad leaf
needed_leaf.erase(std::remove(needed_leaf.begin(),needed_leaf.end(),0),needed_leaf.end());
return {begin(global_matrix_interaction),end(global_matrix_interaction)};
}
/**
* compute the number of block where leaf are stock
* Leafs are stock in leaf_needed
* @author benjamin.dufoyer@inria.fr
* @param tree GroupTree
* @param leaf_needed Vector where leaf are stock
* @return the number of block
*/
template<class GroupOctreeClass>
std::vector<size_t> get_nb_block(GroupOctreeClass& tree,
size_t* leaf_needed,
size_t nb_leaf)
{
std::vector<size_t> block_to_send(nb_leaf,0);
int idx_vect = 0 ;
for(unsigned i = 0 ; i < nb_leaf; ++i){
bool leaf_ok = false;
for(unsigned idxGroup = 0 ; idxGroup < (unsigned)tree.getNbParticleGroup() ; ++idxGroup){
if(leaf_ok)
break;
auto* containers = tree.getParticleGroup(idxGroup);
for(unsigned idx_leaf_block = 0; idx_leaf_block < (unsigned)containers->getNumberOfLeavesInBlock() ; ++idx_leaf_block){
if(leaf_ok)
break;
MortonIndex curr_m_idx = containers->getLeafMortonIndex(idx_leaf_block);
if(leaf_needed[i] == (unsigned)curr_m_idx){
if(idx_vect == 0){
block_to_send[idx_vect] = idxGroup;
++idx_vect;
leaf_ok = true;
} else if(block_to_send[idx_vect-1] != idxGroup){
block_to_send[idx_vect] = idxGroup;
++idx_vect;
leaf_ok = true;
}
}
}
}
}
return {block_to_send.begin(),block_to_send.begin()+idx_vect};
}
std::vector<size_t> send_get_leaf_number(
std::vector<size_t>& needed_leaf,
std::vector<size_t>& global_matrix_interaction,
size_t& nb_msg_recv,
size_t& nb_leaf_recv,
const inria::mpi_config& conf)
{
// Get MPI Info
const unsigned nb_proc = conf.comm.size();
const unsigned my_rank = conf.comm.rank();
// allocate tab of mpi status for synchronisazion
inria::mpi::request tab_mpi_status[nb_msg_recv];
// Initialiser the reception vector
std::vector<size_t> vect_recv(nb_leaf_recv,0);
// Initialize every variable
int idx_status = 0;
size_t idx_vect = 0;
int idx_proc = 0;
// Posting every recv message
for(unsigned i = my_rank; i < global_matrix_interaction.size() ; i+= nb_proc ){
if(global_matrix_interaction[i] != 0){
tab_mpi_status[idx_status] = conf.comm.irecv(
&vect_recv.data()[idx_vect],
(int)global_matrix_interaction[i],
my_MPI_SIZE_T,
idx_proc,1
);
idx_status += 1;
idx_vect += global_matrix_interaction[i];
}
idx_proc += 1;
}
// Posting every send message
idx_proc = 0;
idx_vect = 0;
for(unsigned i = my_rank*nb_proc;
i < (my_rank*nb_proc)+nb_proc;
++i)
{
if(global_matrix_interaction[i] != 0){
conf.comm.isend(
&needed_leaf.data()[idx_vect],
(int)global_matrix_interaction[i],
my_MPI_SIZE_T,
idx_proc,1
);
idx_vect += global_matrix_interaction[i];
}
idx_proc += 1;
}
if(nb_msg_recv != 0 ){
inria::mpi::request::waitall(nb_msg_recv,tab_mpi_status);
}
return{begin(vect_recv),end(vect_recv)};
}
/**
* This function send and get the number of block sended to me by other proc
* and sended from me to other proc
* In this step, i anticipe the next step, i stock the index of the block i will
* send to not reproduct the computation
*
* @author benjamin.dufoyer@inria.fr
* @param vect_recv it's the vector where the list of leaf
* needed by a proc is stock
* @param global_matrix_interaction it's the interaction matrix
* @param nb_msg_recv it's the number of message i will recev
* from other proc
* @param tree it's the local grouptree
* @param nb_block_to_receiv it's the number of block i wiil recev,
* it's here where the result is stock
* @param block_to_send it's the index of block what i will send
* in the next step
* @param conf MPI conf
*/
template<class GroupOctreeClass>
void send_get_number_of_block(
std::vector<size_t>& vect_recv,
std::vector<size_t>& global_matrix_interaction,
size_t& nb_msg_recv,
GroupOctreeClass& tree,
std::vector<std::pair<int,int>>& nb_block_to_receiv,
std::vector<std::pair<int,std::vector<size_t>>>& block_to_send,
const inria::mpi_config& conf
)
{
// Get MPI Info
const unsigned nb_proc = conf.comm.size();
const unsigned my_rank = conf.comm.rank();
int idx_status = 0;
int idx_proc = 0;
inria::mpi::request tab_mpi_status[nb_msg_recv];
// Post the number reception of the number of block