Commit 8bff5bfd authored by COULAUD Olivier's avatar COULAUD Olivier
Browse files

Check distribute routines

parent 165dc472

#include <array>
#include "parameters.hpp"
#include "scalfmm/container/particle.hpp"
......@@ -14,8 +14,12 @@
#include "scalfmm/utils/parallel_manager.hpp"
auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
{
parallel_manager para;
para.init();
#ifndef LET
constexpr int dimension = 3;
constexpr int nb_inputs_near = 1;
constexpr int nb_outputs_near = 1;
......@@ -24,7 +28,8 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
using value_type = double;
using particle_type = scalfmm::container::particle<value_type, dimension, value_type, nb_inputs_near, value_type,
nb_outputs_near, std::size_t>;
using container_type = scalfmm::container::particle_container<particle_type>;
using read_particle_type = scalfmm::container::particle<value_type, dimension, value_type, nb_inputs_near, value_type,
0, std::size_t>; using container_type = scalfmm::container::particle_container<particle_type>;
using position_type = typename particle_type::position_type;
using cell_type =
scalfmm::component::cell<value_type, dimension, nb_inputs_far, nb_outputs_far, std::complex<value_type>>;
......@@ -60,8 +65,7 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
// const auto order{parser.get<args::order>()};
//
//
parallel_manager para;
para.init();
///////////////////////////////////////////////////////////////////////////////////////////////////////
/// Read the data in parallel
......@@ -79,6 +83,7 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
box_type box(width, centre);
//
container_type container(local_number_of_particles);
std::vector <particle_type> particles_set(local_number_of_particles);
for(std::size_t idx = 0; idx < local_number_of_particles; ++idx)
{
loader.fillParticle(values_to_read, nb_val_to_red_per_part);
......@@ -88,7 +93,8 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
{
e = values_to_read[ii++];
}
container.push_particle(idx, p);
particles_set[idx] = p;
//container.push_particle(idx, p);
}
///
///////////////////////////////////////////////////////////////////////////////////////////////////////
......@@ -107,7 +113,11 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
std::vector<MortonIndex> mortonCellDistribution;
group_tree_type* letGroupTree = nullptr;
//
scalfmm::tree::let::buildLetTree(para, loader, container, box, tree_height, group_size, letGroupTree,
scalfmm::tree::let::buildLetTree(para, loader, particles_set, box, tree_height, group_size, letGroupTree,
mortonCellDistribution, nb_block);
///
......@@ -124,7 +134,78 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
/// writer.write( .... todo)
///
///////////////////////////////////////////////////////////////////////////////////////////////////////
#else
// MPI communicator
// auto world = inria::mpi::communicator::world();
inria::mpi_config conf(para.get_communicator());
const auto rank = conf.comm.rank();
const int N = 30;
// Input is an empty vector, no allocation is done
std::vector<int> input{};
int block =
rank == para.get_num_processes() - 1 ? N - int(N / 3) * rank : int(N / 3);
input.resize(block);
std::cout << "block: " << block << std::endl;
std::iota(begin(input), end(input), block * rank);
std::cout << " proc " << conf.comm.rank() << " in (" << input.size()
<< "): ";
for (auto& v : input) {
std::cout << v << " ";
}
std::cout << std::endl;
std::array<int, 3> proportions{100, 200, 300};
int total = std::accumulate(proportions.begin(), proportions.end(), 0);
int N_local = N;
// int(N * proportions[rank] / total + 1);
// All processes have a large enough buffer
std::vector<int> output;
output.resize(N_local);
int none_value = 0;
try {
std::array<int, 3> proportions{100, 200, 300};
// Distribution object, keeps all even elements on 0, all odd elements on 1
struct dist {
std::vector<int> mapping{15, 20, 30};
int operator()(const int& i) {
for (int k = 0; k < mapping.size(); ++k) {
if (i <= mapping[k]) {
return k;
}
}
}
};
inria::proportional_distribution<> distrib(conf, input, proportions);
// inria::distribute(conf, input, output, distrib);
inria::distribute(conf, input, output, dist{});
// inria::distribute(conf, begin(input), end(input), begin(output),
// end(output), distrib);
} catch (std::out_of_range& e) {
std::cerr << e.what() << '\n';
}
std::cout <<" output.size() " << output.size() << std::endl;
for (int i = 1; i < output.size(); ++i) {
N_local = i;
if (output[i] == none_value) {
N_local = i;
break;
}
}
output.resize(N_local);
std::cout <<" proc "<< rank<< " out (" << output.size() << "): ";
for (auto& v : output) {
std::cout << v << " ";
}
std::cout << std::endl;
#endif
para.end();
return 0;
}
#ifndef SCALFMM_TREE_LET_HPP
#define SCALFMM_TREE_LET_HPP
#include <algorithm>
#include <array>
#include <iostream>
#include <string>
#include <tuple>
#include <inria/algorithm/distributed/sort.hpp>
#include <inria/algorithm/distributed/unique.hpp>
#include <inria/linear_tree/balance_tree.hpp>
//#include <scalfmm/tree/group_linear_tree.hpp>
#include <scalfmm/tree/utils.hpp>
#include <scalfmm/utils/parallel_manager.hpp>
#ifdef SCALFMM_USE_MPI
#include "inria/algorithm/distributed/distribute.hpp"
#include <inria/algorithm/distributed/mpi.hpp>
#include <mpi.h>
#endif
namespace scalfmm::tree
{
namespace distrib
{
// template<typename index_type>
// index_type send_get_max_morton_idx(parallel_manager& manager, const index_type& morton_idx)
// {
// // Setting parametter
// auto comm = manager.get_communicator();
// const int nb_proc = comm.size();
// const int my_rank = comm.rank();
// // compute the buffer size
// index_type buff_recev{-1};
// if(nb_proc != 1)
// {
// inria::mpi::request tab_mpi_status;
// // if i'm not the last proc
// const int receiver = (my_rank + 1 == nb_proc) ? MPI_PROC_NULL : my_rank + 1;
// const int sender = (my_rank == 0) ? MPI_PROC_NULL : my_rank - 1;
// if(my_rank != nb_proc - 1)
// {
// // Sending my max
// comm.isend(&morton_idx, 1, inria::mpi::get_datatype<index_type>(), my_rank + 1, 1);
// }
// // if i'm not the first proc
// if(my_rank != 0)
// {
// // Receive the max of the left proc
// tab_mpi_status = comm.irecv(&buff_recev, 1, inria::mpi::get_datatype<index_type>(),
// my_rank - 1, 1);
// }
// // Waiting for the result of the request, if my rank is 0
// // I don't need to wait
// if(my_rank != 0)
// inria::mpi::request::waitall(1, &tab_mpi_status);
// }
// return buff_recev;
// }
///
/// \brief balanced_leaves
///
/// 1) we distribute the particles according to their morton index. The the
/// leaves
/// are split on the processor by their morton index
/// 2) balanced the morton index by some criteria to define
///
/// input[inout] mortonArray the morton index located on the processor
template<typename MortonArray_type>
void balanced_leaves(parallel_manager& manager, MortonArray_type& mortonArray)
{
auto minIndex{mortonArray[0][0]};
auto maxIndex{mortonArray[mortonArray.size() - 1][0]};
std::cout << "Mordon Index in [" << minIndex << ", " << maxIndex << "]\n";
auto maxIndexNextProc = send_get_max_morton_idx(manager, maxIndex);
std::cout << " Index on next prox " << maxIndexNextProc << std::endl;
// exchange max
// send max to the next proc
// send & receive number of values then the values
// add/remove the values on the morton index
}
} // namespace distrib
namespace let
{
template<typename Vector_type>
void print(const std::string&& title, Vector_type& v)
{
std::cout << title;
for(int i: v)
std::cout << i << ' ';
std::cout << '\n';
};
//
// @param[in] mpi_comm the MPI communicator
// @param[inout] myParticleslocal array of particles on my node. On output the
// array is sorted
// @param[in] total number of particles in the simulation
// @param[in] box size of the simulation box
// @param[in] TreeHeight Height of the tree
// @param[inout] localGroupTree the LET of the octree
// @param[out] m_idx_distribution Distribution of the leaves on the
// processors
// @param[out] nb_blocks
template<typename Loader_type, typename Index_type, typename Container_type, typename Tree_type,
typename Box_type>
void testbuildLetTree(parallel_manager& manager, Loader_type& loader, Container_type& particle_container,
const Box_type& box, const int tree_height, const int groupSize,
Tree_type*& localGroupTree, std::vector<Index_type>& m_idx_distribution, int& nb_blocks)
{
inria::mpi_config conf(manager.get_communicator());
int level = tree_height - 1;
auto nb_part = particle_container.size();
int rank = manager.get_process_id();
auto minMorton = scalfmm::index::get_morton_index(particle_container[0].position(), box, level);
auto maxMorton = scalfmm::index::get_morton_index(particle_container[nb_part - 1].position(), box, level);
std::cout << " proc :" << rank << " particles :" << nb_part << " minMorton " << minMorton << " maxMorton "
<< maxMorton << std::endl;
inria::distribute(conf, particle_container, inria::uniform_distribution{conf, particle_container});
nb_part = particle_container.size();
minMorton = scalfmm::index::get_morton_index(particle_container[0].position(), box, level);
maxMorton = scalfmm::index::get_morton_index(particle_container[nb_part - 1].position(), box, level);
std::cout << " proc :" << rank << " Linear tree :" << nb_part << " minMorton " << minMorton << " maxMorton "
<< maxMorton << std::endl;
/// The particles are distributed so we can construct the octree as in sequential case
///
/// ToDO
// // create GroupLinearTree
// scalfmm::tree::linear::FGroupLinearTree<typename decltype(linear_tree)::value_type,
// Index_type>
// group_linear_tree{manager.get_communicator()};
// group_linear_tree.create_local_group_linear_tree(&linear_tree, groupSize);
// std::cout << "%%%%%%%%%% " << group_linear_tree << std::endl;
// m_idx_distribution = group_linear_tree.get_index_particle_distribution_implicit();
/// Now we construct the Let by adding the cells used in the different operators (P2P, M2M, M2L, L2L)
///
/// TO DO
///
/// The Local essential tree (Let) is build
}
//
// @param[in] mpi_comm the MPI communicator
// @param[inout] myParticleslocal array of particles on my node. On output the
......@@ -30,10 +157,15 @@ namespace scalfmm::tree
std::vector<Index_type>& m_idx_distribution, int& nb_blocks)
{
//
using indexing_structure = std::array<std::size_t, 3>;
//
using indexing_structure = std::array<std::int64_t, 3>;
constexpr int dimension = Box_type::dimension;
const std::size_t localNumberOfParticles = particle_container.number_elements();
//
inria::mpi_config conf(manager.get_communicator());
int rank = manager.get_process_id();
int leaf_level = tree_height - 1;
auto nb_part = particle_container.size();
const std::size_t localNumberOfParticles = particle_container.size();
//
// use an array (mortonIndex, position, processor) to construct the distribution on the processes
//
......@@ -48,40 +180,53 @@ namespace scalfmm::tree
tuple_of_indexes[part][2] = manager.get_process_id();
tuple_of_indexes[part][1] = part;
tuple_of_indexes[part][0] =
scalfmm::index::get_morton_index(particle_container.position(part), box, level);
scalfmm::index::get_morton_index(particle_container[part].position(), box, level);
}
std::vector<Index_type> tmpleafMortonIdx(tuple_of_indexes.size());
#pragma omp parallel for shared(tmpleafMortonIdx, tuple_of_indexes, localNumberOfParticles, dimension)
for(std::size_t part = 0; part < localNumberOfParticles; ++part)
{
tmpleafMortonIdx[part] = tuple_of_indexes[part][0] /*>> dimension*/;
}
// Now i have all of my particles in a vector, they all have a morton index
// now we will sort them
print("rank(" + std::to_string(rank) + ") before sort: ", tmpleafMortonIdx);
/// Now i have all of my particles in a vector, they all have a morton index
/// now we will sort them
inria::sort(manager.get_communicator(), tuple_of_indexes,
[](const auto& p1, const auto& p2) { return p1[0] < p2[0]; });
//
// set the morton index at the leaf level (we just have to compute the parent morton index
//#pragma omp parallel for shared(tuple_of_indexes,localNumberOfParticles, dimension)
for(std::size_t part = 0; part < localNumberOfParticles; ++part)
// set the morton index at the leaf level (we just have to compute
// the parent morton index
std::vector<Index_type> leafMortonIdx(tuple_of_indexes.size());
#pragma omp parallel for shared(leafMortonIdx, tuple_of_indexes, localNumberOfParticles, dimension)
for(std::size_t part = 0; part < tuple_of_indexes.size(); ++part)
{
tuple_of_indexes[part][0] = tuple_of_indexes[part][0] >> dimension;
leafMortonIdx[part] = tuple_of_indexes[part][0] /* >> dimension*/;
}
// if(manager.io_master())
// {
// for(std::size_t part = 0; part < localNumberOfParticles; ++part)
// {
// std::cout << part << " " << tuple_of_indexes[part][0] << " " <<
// tuple_of_indexes[part][1]
// << std::endl;
// }
// }
//
/// A morton index should be own by only one process
///
print("rank(" + std::to_string(rank) + ") leafMortonIdx: ", leafMortonIdx);
#ifdef TT
auto last = std::unique(leafMortonIdx.begin(), leafMortonIdx.end());
leafMortonIdx.erase(last, leafMortonIdx.end());
// if(leafMortonIdx[leafMortonIdx.size])
print("rank(" + std::to_string(rank) + ") leafMortonIdx U: ", leafMortonIdx);
// Whe check
std::cout << "Local size before balanced step " << tuple_of_indexes.size() << std::endl;
// scalfmm::tree::distrib::balanced_leaves(manager, tuple_of_indexes);
// But we have lost the link between particles and the distribution
// create the linear tree
// a linear tree is a tree, with only the leaf
level = tree_height - 1;
auto linear_tree =
// scalfmm::tree::linear::create_balanced_linear_tree_at_level(manager.get_communicator(),
// tuple_of_indexes);
inria::linear_tree::create_balanced_linear_tree_at_leaf(manager.get_communicator(), level, box,
tuple_of_indexes);
// auto linear_tree =
// scalfmm::tree::linear::create_balanced_linear_tree_at_level(manager.get_communicator(),
// tuple_of_indexes);
// inria::linear_tree::create_balanced_linear_tree_at_leaf(manager.get_communicator(),
// level, box,
// tuple_of_indexes);
// Now we have a balanced tree (leaf per processor) according to the morton index
std::cout << "Linear tree \n";
// std::cout << "Linear tree \n";
// if(manager.io_master())
// {
// for(auto& e: linear_tree)
......@@ -90,11 +235,12 @@ namespace scalfmm::tree
// std::endl;
// }
// }
auto start = linear_tree.front()[0], end = linear_tree.back()[0];
std::cout << "Local size " << linear_tree.size() << " " << start << " " << end << " "
<< tuple_of_indexes.front()[0] << " " << tuple_of_indexes.back()[0] << std::endl;
// auto start = linear_tree.front()[0], end = linear_tree.back()[0];
// std::cout << "Local size " << linear_tree.size() << " " << start << " " << end << " "
// << tuple_of_indexes.front()[0] << " " << tuple_of_indexes.back()[0] <<
// std::endl;
// Now we have to redistribute the particles
// // Now we have to redistribute the particles
//
// int nbPart = 0;
// for(std::size_t part = 0; part < particle_container.number_elements(); ++part)
......@@ -126,6 +272,8 @@ namespace scalfmm::tree
/// TO DO
///
/// The Local essential tree (Let) is build
//
#endif
}
} // namespace let
......
Supports Markdown
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