Skip to content
GitLab
Commit a200aecd authored by ESTERIE Pierre's avatar ESTERIE Pierre
Browse files

some more time loop work

parent ae27ccec
Hide whitespace changes
Inline Side-by-side
experimental/CMakeLists.txt
View file @ a200aecd
......@@ -2,7 +2,7 @@
#
# Project Declaration
#--------------------
project(scalfmm3 CXX)
project(scalfmm3 CXX C)
# check if compiling into source directories
string(COMPARE EQUAL "${CMAKE_CURRENT_SOURCE_DIR}" "${CMAKE_CURRENT_BINARY_DIR}" insource)
if(insource)
......
experimental/cmake/dependencies.cmake
View file @ a200aecd
......@@ -8,3 +8,11 @@ include(cmake/dependencies/blas.cmake)
include(cmake/dependencies/fftw.cmake)
include(cmake/dependencies/starpu.cmake)
if(SCALFMM_USE_CATALYST)
find_package(catalyst REQUIRED)
if(catalyst_FOUND)
message(STATUS "Catalyst found.")
target_link_libraries(${CMAKE_PROJECT_NAME} INTERFACE catalyst::catalyst catalyst::conduit_headers catalyst::conduit catalyst::blueprint)
endif()
endif()
experimental/examples/catalyst_related.hpp 0 → 100644
View file @ a200aecd
// --------------------------------
// See LICENCE file at project root
// File : group_tree.hpp
// --------------------------------
#ifndef CATALYST_RELATED_HPP
#define CATALYST_RELATED_HPP
#include "scalfmm/meta/const_functions.hpp"
#include "scalfmm/meta/utils.hpp"
#include "scalfmm/utils/parallel_manager.hpp"
#include <catalyst.h>
#include <conduit.hpp>
#include <conduit_blueprint_mesh.h>
#include <conduit_blueprint_exports.h>
#include <conduit_cpp_to_c.hpp>
#include <inria/tcli/tcli.hpp>
//#include <iostream>
// catalyst
namespace catalyst_adaptor
{
template<typename... Parameters>
auto initialize(inria::tcli::parser<Parameters...> const& parser) -> void
{
parallel_manager m;
m.init();
// here goes scrpits of mpi init communicator
conduit::Node node;
//node["catalyst/scripts/script" + std::to_string(0)].set_string("/home/p13ro/dev/cpp/build/scalfmm/gcc-10.2/catalyst_pipeline.py");
node["catalyst/mpi_comm"].set(0);
catalyst_initialize(conduit::c_node(&node));
}
template<typename GroupTree>
auto execute(std::size_t step, GroupTree const& tree)
{
conduit::Node exec_params;
//// add time/cycle information
auto& state = exec_params["catalyst/state"];
state["timestep"].set(step);
state["time"].set(step);
// Add channels.
// We only have 1 channel here. Let's name it 'grid'.
//auto particles_channel = exec_params["catalyst/channels/particles"];
auto& box_channel = exec_params["catalyst/channels/box"];
// Since this example is using Conduit Mesh Blueprint to define the mesh,
// we set the channel's type to "mesh".
box_channel["type"].set("mesh");
const auto n_corners = scalfmm::meta::pow(2,GroupTree::dimension);
// now create the mesh.
auto& box_mesh = box_channel["data"];
box_mesh["coordsets/coords/type"].set("explicit");
box_mesh["coordsets/coords/values/x"].set(conduit::DataType::float64(n_corners));
box_mesh["coordsets/coords/values/y"].set(conduit::DataType::float64(n_corners));
auto const& box = tree.box();
conduit::float64* box_sim_x = box_mesh["coordsets/coords/values/x"].value();
conduit::float64* box_sim_y = box_mesh["coordsets/coords/values/x"].value();
for(std::size_t n{0}; n<n_corners; ++n)
{
box_sim_x[n] = scalfmm::meta::get<0>(box.corner(n));
box_sim_y[n] = scalfmm::meta::get<1>(box.corner(n));
std::cout << box.corner(n) << '\n';
}
// Next, add topology
box_mesh["topologies/mesh/type"].set("line");
box_mesh["topologies/mesh/coordset"].set("coords");
//box_mesh["topologies/mesh/elements/shape"].set("box");
box_mesh["topologies/mesh/elements/connectivity"].set_int32_vector({0,1,2,3});
catalyst_execute(conduit::c_node(&exec_params));
}
auto finalize() -> void
{
conduit::Node node;
catalyst_finalize(conduit::c_node(&node));
}
}
#endif // CATALYST_RELATED_HPP
experimental/examples/test-time-loop.cpp
View file @ a200aecd

#include "parameters.hpp"
#include "scalfmm/algorithms/common.hpp"
#include "scalfmm/container/particle.hpp"
#include "scalfmm/container/particle.hpp"
#include "scalfmm/container/particle_container.hpp"
#include "scalfmm/container/point.hpp"
#include "scalfmm/interpolation/uniform.hpp"
#include "scalfmm/operators/fmmOperators.hpp"
#include "scalfmm/operators/laplace_operators.hpp"
#include "scalfmm/matrix_kernels/laplace.hpp"
#include "scalfmm/matrix_kernels/scalar_kernels.hpp"
//#include "scalfmm/operators/generic/l2p.hpp"
#include "scalfmm/algorithms/common.hpp"
#include "scalfmm/algorithms/sequential/sequential.hpp"
#include "scalfmm/tools/colorized.hpp"
#include "scalfmm/tools/tikz_writer.hpp"
#include "scalfmm/tree/box.hpp"
#include "scalfmm/tree/cell.hpp"
#include "scalfmm/tree/for_each.hpp"
#include "scalfmm/tree/group_tree.hpp"
#include "scalfmm/tree/leaf.hpp"
#include "scalfmm/utils/sort.hpp"
#include "scalfmm/tools/colorized.hpp"
#include "scalfmm/tools/fma_loader.hpp"
#include "scalfmm/tools/vtk_writer.hpp"
#include "scalfmm/utils/accurater.hpp"
......@@ -31,7 +26,9 @@
#include <chrono>
#include <cstdio>
#include <inria/tcli/tcli.hpp>
#include <inria/tcli/help_descriptor.hpp>
#include <iostream>
#include <limits>
#include <sstream>
#include <string>
#include <sys/types.h>
......@@ -40,6 +37,13 @@
#include <unistd.h>
#include <utility>
#include <vector>
//#define USE_CATALYST
#ifdef USE_CATALYST
#include "catalyst_related.hpp"
#endif
// example:
// ./examples/RelWithDebInfo/test-time-loop --tree-height 4 -gs 2 --order 4 --dimension 2...
namespace local_args
......@@ -90,7 +94,7 @@ namespace local_args
using type = int;
type def = 250;
};
struct dimension
struct dimension : inria::tcli::required_tag
{
inria::tcli::str_vec flags = {"--dimension", "--d"};
const char* description = "Dimension : \n 2 for dimension 2, 3 for dimension 3";
......@@ -99,13 +103,28 @@ namespace local_args
};
struct time_steps
{
inria::tcli::str_vec flags = {"--time-step", "--ts"};
inria::tcli::str_vec flags = {"--time-steps", "--ts"};
const char* description = "Number of time steps";
using type = std::size_t;
type def = 1;
};
struct delta
{
inria::tcli::str_vec flags = {"--delta"};
const char* description = "Delta to apply on forces";
using type = double;
type def = 0.1;
};
struct catalyst
{
inria::tcli::str_vec flags = {"--catalyst"};
const char* description = "Enable catalist";
using type = bool;
type def = false;
};
} // namespace local_args
template<std::size_t Dimension, typename ContainerType, typename ValueType>
void read_data(const std::string& filename, ContainerType*& container,
scalfmm::container::point<ValueType, Dimension>& Centre, ValueType& width)
......@@ -141,6 +160,58 @@ void read_data(const std::string& filename, ContainerType*& container,
}
}
template<typename GroupTree, typename ValueType>
auto update_particles(GroupTree& tree, ValueType delta) -> void
{
scalfmm::component::for_each_leaf(tree.begin(), tree.end(), [delta](auto& leaf)
{
auto& particles{leaf.particles()};
scalfmm::container::point<ValueType, GroupTree::dimension> force{};
for(std::size_t i{0}; i<particles.size(); ++i)
{
auto proxy = particles.at(i);
force = -delta*proxy.position();
proxy.position() += force;
}
});
}
template<typename GroupTree>
auto get_new_box(GroupTree const& tree) -> typename GroupTree::box_type
{
using value_type = typename GroupTree::leaf_type::particle_type::position_type::value_type;
using box_type = typename GroupTree::box_type;
auto const old_center = tree.box().center();
std::vector<value_type> max(GroupTree::dimension, 0);
std::vector<value_type> min(GroupTree::dimension, 0);
scalfmm::component::for_each_leaf(tree.begin(), tree.end(), [&min, &max, &old_center](auto const& leaf)
{
auto& particles{leaf.particles()};
for(std::size_t i{0}; i<particles.size(); ++i)
{
const auto position = particles.at(i).position();
//std::cout << particles.at(i).inputs(0) << '\n';
for(std::size_t d{0}; d<GroupTree::dimension; ++d)
{
min.at(d) = std::min(position.at(d), min.at(d));
max.at(d) = std::max(position.at(d), max.at(d));
//std::cout << particles.at(i).outputs(d) << ' ';
}
//std::cout << '\n';
}
});
value_type width{0};
for(std::size_t i{0}; i < GroupTree::dimension; ++i)
{
width = std::max(std::abs(min.at(i)), max.at(i));
}
return box_type(width + std::numeric_limits<value_type>::epsilon(), old_center);
}
template<std::size_t Dimension, typename FmmOperatorType, typename... Parameters>
auto run(inria::tcli::parser<Parameters...> const& parser) -> int
{
......@@ -148,6 +219,8 @@ auto run(inria::tcli::parser<Parameters...> const& parser) -> int
// timer
scalfmm::utils::timer time{};
static constexpr auto dimension{Dimension};
using value_type = double;
// ---------------------------------------
using near_matrix_kernel_type = typename FmmOperatorType::near_field_type::matrix_kernel_type;
......@@ -179,11 +252,14 @@ auto run(inria::tcli::parser<Parameters...> const& parser) -> int
const std::string input_file(parser.template get<local_args::input_file>());
const std::string output_file(parser.template get<local_args::output_file>());
const std::string visu_file(parser.template get<local_args::visu_file>());
const bool catalyst_enable(parser.template get<local_args::catalyst>());
const auto delta = parser.template get<local_args::delta>();
std::cout << scalfmm::colors::blue << "<params> Tree height : " << tree_height << scalfmm::colors::reset << '\n';
std::cout << scalfmm::colors::blue << "<params> Group Size : " << group_size << scalfmm::colors::reset << '\n';
std::cout << scalfmm::colors::blue << "<params> Runtime order : " << order << scalfmm::colors::reset << '\n';
std::cout << scalfmm::colors::blue << "<params> Time Steps : " << n_time_steps << scalfmm::colors::reset << '\n';
std::cout << scalfmm::colors::blue << "<params> Delta forces : " << delta << scalfmm::colors::reset << '\n';
// ---------------------------------------
// Open particle file
......@@ -232,23 +308,45 @@ auto run(inria::tcli::parser<Parameters...> const& parser) -> int
std::cout << scalfmm::colors::yellow << "Kernel and Interp created in " << time.elapsed() << "ms\n"
<< scalfmm::colors::reset;
#ifdef USE_CATALYST
if(catalyst_enable)
{
std::cout << scalfmm::colors::on_blue << "Initializing catalyst..." << scalfmm::colors::reset << '\n';
catalyst_adaptor::initialize(parser);
}
#endif
// ---------------------------------------
// start algorithm for step 0
auto operator_to_proceed = scalfmm::algorithms::all;
scalfmm::algorithms::sequential(tree, std::move(fmm_operator), operator_to_proceed);
scalfmm::algorithms::sequential(tree, fmm_operator, operator_to_proceed);
// ---------------------------------------
// time loop
for(std::size_t steps{1}; steps < n_time_steps; ++steps)
{
// ---
// compute the forces
// ---
// update particles
update_particles(tree, delta);
// ---
// visu related post treatment
if(!visu_file.empty())
{
scalfmm::tools::io::exportVTKxml(std::to_string(steps-1)+visu_file, tree, container->size());
}
#ifdef USE_CATALYST
if(catalyst_enable)
{
catalyst_adaptor::execute(steps-1, tree);
}
#endif
// ---
// compute the corners and set new box
auto const new_box = get_new_box(tree);
std::cout << scalfmm::colors::on_green << "New box width is " << new_box.width(0) << scalfmm::colors::reset << '\n';
// ---
// compute new morton indices
// ---
// move particles according to new morton indices
// ---
......@@ -256,8 +354,29 @@ auto run(inria::tcli::parser<Parameters...> const& parser) -> int
// ---
// build tree levels
// ---
// start algorithm at step : compute forces
// ---
}
// ---------------------------------------
// last step update
// ---
// update particles
// ---
update_particles(tree, delta);
// ---
// visu related post treatment
if(!visu_file.empty())
{
scalfmm::tools::io::exportVTKxml(std::to_string(n_time_steps-1) + visu_file, tree, container->size());
}
#ifdef USE_CATALYST
if(catalyst_enable)
{
}
#endif
// ---------------------------------------
// write output fma
if(!output_file.empty())
......@@ -283,10 +402,10 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
{
//
// Parameter handling
auto parser =
inria::tcli::make_parser(inria::tcli::help{}, local_args::tree_height{}, local_args::order{},
local_args::input_file{}, local_args::output_file{}, local_args::block_size{},
local_args::dimension{}, local_args::time_steps{}, local_args::visu_file{});
auto parser = inria::tcli::make_parser(inria::tcli::help{}, local_args::tree_height{}, local_args::order{},
local_args::input_file{}, local_args::output_file{},
local_args::block_size{}, local_args::dimension{}, local_args::time_steps{},
local_args::visu_file{}, local_args::delta{}, local_args::catalyst{});
parser.parse(argc, argv);
// Getting command line parameters
......@@ -297,28 +416,29 @@ auto main([[maybe_unused]] int argc, [[maybe_unused]] char* argv[]) -> int
case 2:
{
constexpr std::size_t dim = 2;
using matrix_kernel_type = scalfmm::matrix_kernels::others::one_over_r2;
using near_field_type = scalfmm::operators::near_field_operator<matrix_kernel_type>;
using near_matrix_kernel_type = scalfmm::matrix_kernels::others::grad_one_over_r2<dim>;
using near_field_type = scalfmm::operators::near_field_operator<near_matrix_kernel_type>;
//
using interpolator_type = scalfmm::interpolation::uniform_interpolator<double, dim, matrix_kernel_type>;
using far_field_type = scalfmm::operators::far_field_operator<interpolator_type>;
using far_matrix_kernel_type = scalfmm::matrix_kernels::others::one_over_r2;
using interpolator_type = scalfmm::interpolation::uniform_interpolator<double, dim, far_matrix_kernel_type>;
using far_field_type = scalfmm::operators::far_field_operator<interpolator_type, true>;
run<dim, scalfmm::operators::fmm_operators<near_field_type, far_field_type>>(parser);
break;
}
case 3:
{
constexpr std::size_t dim = 3;
using matrix_kernel_type = scalfmm::matrix_kernels::others::one_over_r2;
// using matrix_kernel_type = scalfmm::matrix_kernels::laplace::one_over_r;
using near_field_type = scalfmm::operators::near_field_operator<matrix_kernel_type>;
//
using interpolator_type = scalfmm::interpolation::uniform_interpolator<double, dim, matrix_kernel_type>;
using far_field_type = scalfmm::operators::far_field_operator<interpolator_type>;
//case 3:
//{
// constexpr std::size_t dim = 3;
// using matrix_kernel_type = scalfmm::matrix_kernels::others::one_over_r2;
// // using matrix_kernel_type = scalfmm::matrix_kernels::laplace::one_over_r;
// using near_field_type = scalfmm::operators::near_field_operator<matrix_kernel_type>;
// //
// using interpolator_type = scalfmm::interpolation::uniform_interpolator<double, dim, matrix_kernel_type>;
// using far_field_type = scalfmm::operators::far_field_operator<interpolator_type>;
run<dim, scalfmm::operators::fmm_operators<near_field_type, far_field_type>>(parser);
break;
}
// run<dim, scalfmm::operators::fmm_operators<near_field_type, far_field_type>>(parser);
// break;
//}
default:
std::cout << "check 1/r^2 Kernel for dimension 2 and 3. Value is \n"
<< " 0 for dimension 2 "
......
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