epirecipes / sir-julia

Various implementations of the classical SIR model in Julia
MIT License
195 stars 39 forks source link

Add C++ RHS using Boost odeint #118

Open sdwfrost opened 2 months ago

sdwfrost commented 2 months ago

Write boilerplate code for using RHS written for Boost odeint, so we can directly apply models from Epiverse-Trace epidemics Rcpp package - see here for headers.

sdwfrost commented 2 months ago

Here's a simple test class to wrap using CxxWrap.jl

#include <iostream>
#include <boost/array.hpp>
#include <boost/numeric/odeint.hpp>

using namespace boost::numeric::odeint;

typedef boost::array< double , 3 > state_type;
typedef boost::array< double , 3 > param_type;

class SIR {
    param_type p;

    public:
        SIR(param_type q) : p(q) {}

        void operator() (const state_type &x , state_type &dxdt , const double t)
        {
            double N = x[0] + x[1] + x[2];
            double beta = p[0];
            double cee = p[1];
            double gamma = p[2];
            dxdt[0] = -beta * cee * x[0] * x[1]/N;
            dxdt[1] = beta * cee * x[0] * x[1]/N - gamma * x[1];
            dxdt[2] = gamma * x[1];
        }

        // Setter method for parameters
        void set_parameters(const param_type &new_p)
        {
            p = new_p;
        }

        // Method to get derivatives for given state and time
        void rhs(state_type &dxdt, const state_type &x, const param_type &p, const double t)
        {
            this->set_parameters(p); // Set the parameters
            (*this)(x, dxdt, t); // Call the functor directly
        } 
};
""";