History.hpp

// Copyright 2021 Arnaud Becheler    <abechele@umich.edu>
 
/*********************************************************************** * This program is free software; you can
 *redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free
 *Software Foundation; either version 2 of the License, or    * (at your option) any later version. *
 *                                                                      *
 ***************************************************************************/
 
#ifndef __HISTORY_H_INCLUDED__
#define __HISTORY_H_INCLUDED__
 
#include "HistoryBase.hpp"
#include "demographic_policy.hpp"
 
namespace quetzal
{
namespace demography
{
template <typename Space, typename DispersalPolicy, typename MemoryPolicy>
class History : public HistoryBase<Space, DispersalPolicy, MemoryPolicy>
{
};
 
template <typename Space, typename Memory>
class History<Space, demographic_policy::individual_based, Memory>
    : public HistoryBase<Space, demographic_policy::individual_based, Memory>
{
    // Using the HistoryBase class constructor
    using HistoryBase<Space, demographic_policy::individual_based, Memory>::HistoryBase;
 
  public:
    template <typename T, typename U, typename V> void simulate_forward(T sim_growth, U kernel, V &gen)
    {
        for (unsigned int t = 0; t < this->nb_generations(); ++t)
        {
            unsigned int landscape_individuals_count = 0;
            for (auto x : this->distribution_area(t))
            {
                auto N_tilde = sim_growth(gen, x, t);
                landscape_individuals_count += N_tilde;
                if (N_tilde >= 1)
                {
                    for (unsigned int ind = 1; ind <= N_tilde; ++ind)
                    {
                        auto y = kernel(gen, x);
                        this->m_flows->add_to_flow_from_to(x, y, t, 1);
                        this->pop_size(y, t + 1) += 1;
                    }
                }
            }
            if (landscape_individuals_count == 0)
            {
                throw std::domain_error("Landscape populations went extinct before sampling");
            }
        }
    }
};
template <typename Space, typename Memory>
class History<Space, demographic_policy::mass_based, Memory>
    : public HistoryBase<Space, demographic_policy::mass_based, Memory>
{
    // Using the HistoryBase constructor
    using HistoryBase<Space, demographic_policy::mass_based, Memory>::HistoryBase;
 
  public:
    template <typename T, typename U, typename V> void simulate_forward(T sim_growth, U kernel, V &gen)
    {
        for (unsigned int t = 0; t < this->nb_generations(); ++t)
        {
            unsigned int landscape_individuals_count = 0;
            for (auto x : this->distribution_area(t))
            {
                auto N_tilde = sim_growth(gen, x, t);
                for (auto const &y : kernel.arrival_space(x))
                {
                    auto m = kernel(x, y);
                    assert(m >= 0.0 && m <= 1.0);
                    // Use of std::ceil to allow oversea migration
                    double nb_migrants = std::ceil(m * static_cast<double>(N_tilde));
                    if (nb_migrants >= 0)
                    {
                        landscape_individuals_count += nb_migrants;
                        this->m_flows->set_flow_from_to(x, y, t, nb_migrants);
                        this->pop_size(y, t + 1) += nb_migrants;
                    }
                }
            }
            if (landscape_individuals_count == 0)
            {
                throw std::domain_error("Landscape populations went extinct before sampling");
            }
        }
    }
};
} // namespace demography
} // namespace quetzal
 
#endif