1#ifndef __PARTICLES_CONTAINER_HPP__
2#define __PARTICLES_CONTAINER_HPP__
4#include "Kokkos_Core_fwd.hpp"
5#include "common/common.hpp"
7#include <Kokkos_Core.hpp>
8#include <common/has_serialize.hpp>
10#include <mc/prng/prng.hpp>
11#include <mc/traits.hpp>
18 using TeamPolicy = Kokkos::TeamPolicy<ComputeSpace>;
19 using TeamMember = TeamPolicy::member_type;
21 template <ModelType M>
struct FillGapFunctor
25 M::SelfParticle _model,
28 std::size_t _to_remove,
29 std::size_t _last_used_index)
30 : status(std::move(_status)), model(std::move(_model)), position(std::move(_position)),
31 ages(std::move(_ages)), offset(
"offset"), to_remove(_to_remove),
32 last_used_index(_last_used_index)
35 Kokkos::deep_copy(offset, 0);
38 KOKKOS_INLINE_FUNCTION
void operator()(
const int i, std::size_t& update,
const bool final)
const
42 std::size_t scan_index = update;
43 update += is_dead ? 1 : 0;
45 if (
final && is_dead && scan_index < to_remove)
48 const auto i_to_remove = i;
49 auto idx_to_move = last_used_index - Kokkos::atomic_fetch_add(&offset(), 1);
51 idx_to_move ==
static_cast<std::size_t
>(i_to_remove))
53 idx_to_move = last_used_index - Kokkos::atomic_fetch_add(&offset(), 1);
57 Kokkos::atomic_exchange(&position(i_to_remove), position(idx_to_move));
59 for (std::size_t i_properties = 0; i_properties < M::n_var; ++i_properties)
61 model(i_to_remove, i_properties) = model(idx_to_move, i_properties);
63 ages(i_to_remove, 0) = ages(idx_to_move, 0);
64 ages(i_to_remove, 1) = ages(idx_to_move, 1);
69 M::SelfParticle model;
72 Kokkos::View<std::size_t, ComputeSpace> offset;
73 std::size_t to_remove;
74 std::size_t last_used_index;
77 template <ModelType M>
struct InsertFunctor
79 InsertFunctor(std::size_t _original_size,
80 M::SelfParticle _model,
82 M::SelfParticle _buffer_model,
84 : original_size(_original_size), model(std::move(_model)), position(std::move(_position)),
85 buffer_model(std::move(_buffer_model)), buffer_position(std::move(_buffer_position))
102 KOKKOS_INLINE_FUNCTION
103 void operator()(
const TeamMember& team)
const
105 auto range = M::n_var;
106 const int i = team.league_rank();
108 position(original_size + i) = buffer_position(i);
109 Kokkos::parallel_for(Kokkos::TeamVectorRange(team, range),
110 [&](
const int& j) { model(original_size + i, j) = buffer_model(i, j); });
113 std::size_t original_size;
114 M::SelfParticle model;
116 M::SelfParticle buffer_model;
171 [[nodiscard]] KOKKOS_INLINE_FUNCTION std::size_t
n_particles()
const
182 [[nodiscard]] KOKKOS_INLINE_FUNCTION
bool
194 template <
class Archive>
void save(Archive& ar)
const
199 serialize_view(ar,
status);
200 serialize_view(ar,
model);
201 serialize_view(ar,
ages);
204 template <
class Archive>
void load(Archive& ar)
209 deserialize_view(ar,
status);
210 deserialize_view(ar,
model);
211 deserialize_view(ar,
ages);
216 [[nodiscard]] KOKKOS_INLINE_FUNCTION
double get_weight(std::size_t idx)
const;
242 template <ModelType Model>
243 KOKKOS_INLINE_FUNCTION
bool
245 std::size_t idx1)
const
247 if (Kokkos::atomic_load(&buffer_index()) < buffer_model.extent(0))
249 const auto idx2 = Kokkos::atomic_fetch_add(&buffer_index(), 1);
250 Model::division(random_pool, idx1, idx2, model, buffer_model);
251 buffer_position(idx2) = position(idx1);
255 Kokkos::printf(
"%ld %ld\r\n", Kokkos::atomic_load(&buffer_index()), buffer_model.extent(0));
261 PROFILE_SECTION(
"ParticlesContainer::merge_buffer")
262 const auto original_size = n_used_elements;
263 const auto n_add_item = buffer_index();
268 __allocate__(original_size + n_add_item);
271 auto get_policy_insert = [=]()
273 if constexpr (std::is_same_v<Kokkos::DefaultHostExecutionSpace,
274 Kokkos::DefaultExecutionSpace>)
276 return TeamPolicy(n_add_item / Kokkos::num_threads(), Kokkos::AUTO, Model::n_var);
280 return TeamPolicy(n_add_item, Kokkos::AUTO, Model::n_var);
284 Kokkos::parallel_for(
286 TeamPolicy(n_add_item, Kokkos::AUTO, Model::n_var),
287 InsertFunctor<Model>(original_size, model, position, buffer_model, buffer_position));
290 n_used_elements += n_add_item;
291 __allocate_buffer__();
296 template <ModelType Model>
299 PROFILE_SECTION(
"ParticlesContainer::__allocate__")
302 if (new_size >= n_allocated_elements)
304 const auto new_allocated_size =
305 static_cast<std::size_t
>(std::ceil(
static_cast<double>(new_size) * allocation_factor));
306 n_allocated_elements = new_allocated_size;
307 Kokkos::resize(position, n_allocated_elements);
308 Kokkos::resize(model, n_allocated_elements, Model::n_var);
309 Kokkos::resize(status, n_allocated_elements);
310 Kokkos::resize(ages, n_allocated_elements);
313 Kokkos::resize(weights, 1);
317 Kokkos::resize(weights, n_allocated_elements);
323 template <ModelType Model>
326 PROFILE_SECTION(
"ParticlesContainer::__shrink__")
327 if (new_size > 0 && (new_size > n_used_elements || force))
329 const auto new_allocated_size =
330 static_cast<std::size_t
>(std::ceil(
static_cast<double>(new_size) * allocation_factor));
331 n_allocated_elements = new_allocated_size;
332 Kokkos::resize(position, n_allocated_elements);
333 Kokkos::resize(model, n_allocated_elements, Model::n_var);
334 Kokkos::resize(status, n_allocated_elements);
335 Kokkos::resize(ages, n_allocated_elements);
338 Kokkos::resize(weights, 1);
342 Kokkos::resize(weights, n_allocated_elements);
349 PROFILE_SECTION(
"ParticlesContainer::__allocate_buffer__")
350 auto buffer_size = buffer_position.extent(0);
351 if (
static_cast<double>(buffer_size) /
static_cast<double>(n_allocated_elements) < buffer_ratio)
353 buffer_size =
static_cast<std::size_t
>(
354 std::ceil(
static_cast<double>(n_allocated_elements) * buffer_ratio));
357 Kokkos::realloc(buffer_position, buffer_size);
358 Kokkos::realloc(buffer_model, buffer_size, Model::n_var);
363 template <ModelType M>
365 : model(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_model"), 0),
366 position(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_position"), 0),
367 status(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_status"), 0),
368 weights(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_weigth"), 0),
369 ages(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_age"), 0),
370 buffer_model(
"buffer_particle_model", 0),
371 buffer_position(
"buffer_particle_position", 0),
372 buffer_index(
"buffer_index"), allocation_factor(default_allocation_factor),
373 n_allocated_elements(0), n_used_elements(n_particle)
380 template <ModelType M>
382 : model(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_model"), 0),
383 position(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_position"), 0),
384 status(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_status"), 0),
385 weights(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_weigth"), 0),
386 ages(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"particle_age"), 0),
387 buffer_model(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"buffer_particle_model"), 0),
388 buffer_position(Kokkos::view_alloc(Kokkos::WithoutInitializing,
"buffer_particle_model")),
389 buffer_index(
"buffer_index"), allocation_factor(default_allocation_factor),
390 n_allocated_elements(0), n_used_elements(0)
397 PROFILE_SECTION(
"ParticlesContainer::remove_dead")
403 if (to_remove == n_used_elements)
408 else if (to_remove > n_used_elements)
410 throw std::runtime_error(
411 "clean_dead: Error in kernel cannot remove more element than existing");
416 const auto new_used_item = n_used_elements - to_remove;
418 const auto last_used_index = n_used_elements - 1;
419 Kokkos::parallel_scan(
421 Kokkos::RangePolicy<ComputeSpace>(0, n_used_elements),
422 FillGapFunctor<M>(status, model, position, ages, to_remove, last_used_index));
425 KOKKOS_ASSERT(this->position.extent(0) == n_allocated_elements);
426 KOKKOS_ASSERT(this->model.extent(0) == n_allocated_elements);
427 KOKKOS_ASSERT(this->status.extent(0) == n_allocated_elements);
428 n_used_elements = new_used_item;
429 if (
static_cast<double>(n_used_elements) /
static_cast<double>(n_allocated_elements) <= 0.1)
432 __shrink__(n_used_elements * 2,
false);
437 template <ModelType M>
438 KOKKOS_INLINE_FUNCTION
void
443 const double weight = get_weight(idx);
444 M::contribution(idx, position(idx), weight, model, contributions);
447 template <ModelType M>
448 [[nodiscard]] KOKKOS_INLINE_FUNCTION
double
Kokkos::Random_XorShift1024_Pool< Kokkos::DefaultExecutionSpace > pool_type
Definition prng.hpp:17
Main owning object for Monte-Carlo particles.
Definition particles_container.hpp:135
~ParticlesContainer()=default
Default destructor.
void __allocate__(std::size_t new_size)
Definition particles_container.hpp:297
Model::SelfParticle model
Definition particles_container.hpp:148
KOKKOS_INLINE_FUNCTION void get_contributions(std::size_t idx, const ContributionView &contributions) const
Definition particles_container.hpp:439
KOKKOS_INLINE_FUNCTION std::size_t n_particles() const
Gets the number of particles in the container.
Definition particles_container.hpp:171
ParticlesContainer(std::size_t n_particle)
Definition particles_container.hpp:364
double allocation_factor
Definition particles_container.hpp:232
ParticleWeigths weights
Definition particles_container.hpp:151
double get_allocation_factor() const
Definition particles_container.hpp:218
Kokkos::View< uint64_t, Kokkos::SharedSpace > buffer_index
Definition particles_container.hpp:231
void __allocate_buffer__()
Definition particles_container.hpp:347
ParticlesContainer()
Definition particles_container.hpp:381
KOKKOS_INLINE_FUNCTION bool handle_division(const MC::KPRNG::pool_type &random_pool, std::size_t idx1) const
Definition particles_container.hpp:244
ParticlesContainer & operator=(ParticlesContainer &&)=default
void __shrink__(std::size_t new_size, bool force)
Definition particles_container.hpp:324
uint64_t n_used_elements
Definition particles_container.hpp:234
auto get_buffer_index() const
Definition particles_container.hpp:189
void merge_buffer()
Definition particles_container.hpp:259
Model::SelfParticle buffer_model
Definition particles_container.hpp:229
KOKKOS_INLINE_FUNCTION double get_weight(std::size_t idx) const
Definition particles_container.hpp:449
std::size_t capacity() const
Definition particles_container.hpp:223
ParticlePositions buffer_position
Definition particles_container.hpp:230
MC::ParticlePositions position
Definition particles_container.hpp:149
void save(Archive &ar) const
Definition particles_container.hpp:194
static constexpr double default_allocation_factor
Definition particles_container.hpp:236
static constexpr double buffer_ratio
Definition particles_container.hpp:137
ParticleAges ages
Definition particles_container.hpp:152
std::size_t counter
Definition particles_container.hpp:177
Model UsedModel
Alias for the model used by the container.
Definition particles_container.hpp:143
void clean_dead(std::size_t to_remove)
Definition particles_container.hpp:395
ParticlesContainer & operator=(const ParticlesContainer &)=default
ParticlesContainer(const ParticlesContainer &)=default
Default copy and move constructors and assignment operators.
void load(Archive &ar)
Definition particles_container.hpp:204
ParticlesContainer(ParticlesContainer &&)=default
std::size_t n_allocated_elements
Definition particles_container.hpp:233
MC::ParticleStatus status
Definition particles_container.hpp:150
Concept to check if a model type has uniform_weight
Definition traits.hpp:173
Namespace that contains classes and structures related to Monte Carlo (MC) simulations.
Definition alias.hpp:9
Kokkos::View< double *, ComputeSpace > ParticleWeigths
Definition alias.hpp:21
Kokkos::View< Status *, ComputeSpace > ParticleStatus
Definition alias.hpp:20
Kokkos::View< double *[2], Kokkos::LayoutLeft, ComputeSpace > ParticleAges
Definition alias.hpp:22
Kokkos::View< uint64_t *, ComputeSpace > ParticlePositions
Definition alias.hpp:19
Kokkos::Experimental::ScatterView< double **, Kokkos::LayoutRight > ContributionView
Definition alias.hpp:31