Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- modelsir_cpp =function(t,x){
- S = x[1]
- I1 = x[2]
- I2 = x[3]
- N=S+I1+I2
- with(as.list(parm), {
- dS=B*I1-mu*S-beta*(S*(I1+I2)/N)
- dI1=beta*(S*(I1+I2)/N)-B*I1-lambda12*I1
- dI2=lambda12*I1
- res=c(dS,dI1,dI2)
- return(res)
- })
- }
- times = seq(0, 10, by = 1/52)
- parm=c(B=0.01,mu=0.008,beta=10,lambda12=1)
- xstart=c(S=900,I1=100,I2=0)
- out = as.data.frame(lsoda(xstart, times, modelsir, parm))
- #include <Rcpp.h>
- #include <boost/array.hpp>
- // include Boost's odeint
- #include <boost/numeric/odeint.hpp>
- // [[Rcpp::depends(BH)]]
- using namespace Rcpp;
- using namespace std;
- using namespace boost::numeric::odeint;
- typedef boost::array< double ,3 > state_type;
- void eqsir(const state_type &x, state_type &dxdt, const double t){
- Rcpp::NumericVector nvec=boost_array_to_nvec(x);
- Rcpp::NumericVector nvec2(3);
- nvec2=my_fun2(nvec,t);
- dxdt=nvec_to_boost_array(nvec2);
- }
- void write_cout_2( const state_type &x , const double t ) {
- // use Rcpp's stream
- Rcpp::Rcout << t << 't' << x[0] << 't' << x[1] << 't' << x[2] << endl;
- }
- typedef runge_kutta_dopri5< state_type > stepper_type;
- // [[Rcpp::export]]
- bool my_fun10_solver() {
- state_type x = { 900 , 100, 0 }; // initial conditions
- integrate_adaptive(make_controlled( 1E-12 , 1E-12 , stepper_type () ) ,
- eqsir , x , 1.0 , 2 , 1/52 , write_cout_2 );
- return true;
- }
- #include <Rcpp.h>
- #include <boost/array.hpp>
- // include Boost's odeint
- #include <boost/numeric/odeint.hpp>
- // [[Rcpp::depends(BH)]]
- using namespace Rcpp;
- using namespace std;
- using namespace boost::numeric::odeint;
- typedef boost::array< double ,3 > state_type;
- void rhs( const state_type &x , state_type &dxdt , const double t ) {
- dxdt[0] = 3.0/(2.0*t*t) + x[0]/(2.0*t);
- dxdt[1] = 3.0/(2.0*t*t) + x[1]/(2.0*t);
- dxdt[2] = 3.0/(2.0*t*t) + x[1]/(2.0*t);
- }
- void write_cout( const state_type &x , const double t ) {
- // use Rcpp's stream
- Rcpp::Rcout << t << 't' << x[0] << 't' << x[1] << 't' << x[2] << endl;
- }
- typedef runge_kutta_dopri5< state_type > stepper_type;
- // [[Rcpp::export]]
- bool boostExample() {
- state_type x = { 1.0 , 1.0, 1.0 }; // initial conditions
- integrate_adaptive(make_controlled( 1E-12 , 1E-12 , stepper_type () ) ,
- rhs , x , 1.0 , 10.0 , 0.1 , write_cout );
- return true;
- }
Add Comment
Please, Sign In to add comment