simulation.hpp

#ifndef __SIMULATIONS_UNIT_HPP__
#define __SIMULATIONS_UNIT_HPP__
 
#include "mixture/species_descriptor.hpp"
#include <Kokkos_Core.hpp>
#include <Kokkos_ScatterView.hpp>
#include <biocma_cst_config.hpp>
#include <cassert>
#include <cma_utils/alias.hpp>
#include <common/common.hpp>
#include <common/logger.hpp>
#include <cstddef>
#include <mc/domain.hpp>
#include <mc/events.hpp>
#include <mc/prng/prng.hpp>
#include <mc/unit.hpp>
#include <memory>
#include <optional>
#include <simulation/descriptors/dimensions.hpp>
#include <simulation/feed_descriptor.hpp>
#include <simulation/kernels/kernels.hpp>
#include <simulation/mass_transfer.hpp>
#include <simulation/probe.hpp>
#include <simulation/scalar_initializer.hpp>
#include <simulation/simulation_exception.hpp>
#include <simulation/simulation_getter.hpp>
#include <simulation/simulation_times.hpp>
// TODO Clean
static constexpr size_t trigger_const_particle_number = 1e6;
 
namespace CmaUtils
{
  class PreCalculatedHydroState;
} // namespace CmaUtils

namespace Simulation
{
  class ScalarSimulation;
 
  class SimulationUnit
  {
  private:
    Getter accesor;
 
  public:
    friend class Getter;
 
    // Getter get() const;
 
    SimulationUnit(
        std::unique_ptr<MC::MonteCarloUnit>&& _unit,
        ScalarInitializer&& scalar_init,
        std::optional<Feed::SimulationFeed> _feed,
        std::shared_ptr<Mixture::SpecieTable> table,
        std::optional<MassTransfer::Type::MtrTypeVariant>&& mt_model_variant,
        Sparam params);
 
    ~SimulationUnit();
 
    SimulationUnit(SimulationUnit&& other) noexcept;
    SimulationUnit(const SimulationUnit& other) = delete;
    SimulationUnit& operator=(SimulationUnit&& rhs) = delete;
    SimulationUnit& operator=(const SimulationUnit& rhs) = delete;
 
    // Setters
    // void setEndTime(double _end_time) noexcept;
    void overwriteStartTime(double _start_time) noexcept;
    void setProbes(ProbeType type, ProbeAutogeneratedBuffer&& _probes);
 
    void setLogger(std::shared_ptr<IO::Logger>);
 
    Getter getter() const noexcept;

    [[nodiscard]] std::unique_ptr<MC::MonteCarloUnit> move_mcu() noexcept;
 
    // Getters
 
    // Simulation methods
    void cycleProcess(auto& container, double d_t, auto& _functors);
    void ode_step(double d_t) const;
    void clearContribution() const noexcept;
    void update_feed(double d_t, bool update_scalar = true) noexcept;
 
    double advance(double d_t) noexcept;
 
    void updateHydro(const CmaUtils::IterationStatePtrType& newstate);
    void updateMCHydro(std::span<const double> newliquid_volume,
                       std::span<const std::size_t> neighors_flat,
                       std::span<const double> proba_flat,
                       std::span<const double> out_flows) const;
 
    void calculate_mtr();
 
    bool checkScalar() const;
 
    // Memory management
    [[deprecated]] void reset();
 
    template <typename Space, ModelType Model>
    KernelInline::CycleFunctors<Space, Model>
    init_functors(MC::ParticlesContainer<Model> container,
                  KernelDispatchOptions options);
 
  private:
    std::unique_ptr<MC::MonteCarloUnit> mc_unit; // NOLINT
 
    std::shared_ptr<Mixture::SpecieTable> m_table;
 
    void updateScalarHydro(const CmaUtils::IterationStatePtrType& newstate);
 
    void setLiquidFlow(CmaUtils::PreCalculatedHydroState* _flows_l);
    void setGasFlow(CmaUtils::PreCalculatedHydroState* _flows_g);
    void post_init_compartments();
 
    [[nodiscard]] MC::KernelConcentrationType getkernel_concentration() const;
 
    MC::ContributionView contribs_scatter;
    MapProbes probes;
    Dimensions dims;
 
    // Attributes
    Simulation::Feed::SimulationFeed m_feed;
 
    bool const_number_simulation = true;
    bool is_two_phase_flow;
    bool f_reaction;
 
    SimulatimeTimes m_times;
 
    void scatter_contribute();
 
    void setMtrModel(MassTransfer::Type::MtrTypeVariant&& variant);
 
    [[nodiscard]] MC::kernelContribution get_kernel_contribution() const;
 
    void post_init_concentration(ScalarInitializer&& scalar_init);
 
    std::shared_ptr<ScalarSimulation> liquid_scalar;
    std::shared_ptr<ScalarSimulation> gas_scalar;
    // TODO add default null value (no model)
    MassTransfer::MassTransferModel m_mt_model;
 
    std::shared_ptr<IO::Logger> logger;
 
    template <ModelType Model>
    void post_cycle(MC::ParticlesContainer<Model>& container,
                    auto& cycle_functors);
 
    template <ModelType Model>
    void pre_cycle(MC::ParticlesContainer<Model>& container,
                   double d_t,
                   auto& functors);
  };
 
  template <ModelType Model>
  void
  SimulationUnit::pre_cycle(MC::ParticlesContainer<Model>& container,
                            double d_t,
                            auto& cycle_functors)
  {
    PROFILE_SECTION("Simulation::pre_cycle")
 
    // Must precede the move kernel: it consumes move_probability, which is
    // only valid for this d_t and this flowmap.
    this->mc_unit->domain.set_move_probability(d_t);
 
    cycle_functors.update(
        d_t, container, this->mc_unit->domain.get_const_inner());
  }
 
  template <typename Space, ModelType Model>
  KernelInline::CycleFunctors<Space, Model>
  SimulationUnit::init_functors(MC::ParticlesContainer<Model> container,
                                KernelDispatchOptions options)
  {
 
    return KernelInline::CycleFunctors<Space, Model>(
        options,
        container,
        mc_unit->rng.random_pool,
        getkernel_concentration(),
        contribs_scatter,
        mc_unit->events,
        mc_unit->domain.get_const_inner(),
        probes[ProbeType ::LeavingTime],
        probes[ProbeType ::DivisionTime]);
  }
 
  void
  SimulationUnit::cycleProcess(auto& container,
                               double d_t,
                               auto& cycle_functors)
  {
    PROFILE_SECTION("cycleProcess")
    using CurrentModel =
        typename std::remove_reference<decltype(container)>::type::UsedModel;
    const size_t n_particle = container.n_particles();
    if (n_particle == 0)
    {
      return;
    }
 
    pre_cycle(container, d_t, cycle_functors);
 
    if (f_reaction)
    {
      this->contribs_scatter.reset();
      cycle_functors.launch_model(n_particle);
    }
 
    if (cycle_functors.move_kernel.need_launch())
    {
      cycle_functors.launch_move(n_particle);
    }
 
    post_cycle<CurrentModel>(container, cycle_functors);
  }
 
  template <ModelType Model>
  void
  SimulationUnit::post_cycle(MC::ParticlesContainer<Model>& container,
                             auto& cycle_functors)
  {
    PROFILE_SECTION("Simulation::post_cycle")
    Kokkos::fence();
    this->scatter_contribute();
    const auto [host_red, host_out_counter]
        = cycle_functors.get_host_reduction();
 
    container.update_and_remove_inactive(host_out_counter, host_red.dead_total);
 
    container.merge_buffer();
 
    if (host_red.waiting_allocation_particle != 0)
    {
      if (logger)
      {
        logger->alert("Simulation",
                      "Overflow of particle not implemented yet (ignore "
                      "_waiting_allocation_particle)");
      }
    }
 
    // set_kernel_contribs_to_host();
  }
 
} // namespace Simulation
 
#endif //__SIMULATIONS_UNIT_HPP__