1

I am using package deSolve to run some models that include an external forcing. To gain speed, I produced compiled code following the vignette of the package (see 6.2 in https://cran.r-project.org/web/packages/deSolve/vignettes/compiledCode.pdf). My problem is that now I want to introduce two external forcings in the compiled code. Does anyone have a working example/knows how to do it?

#include <R.h>

static double parms[6];
static double forc[1]; 

/* A trick to keep up with the parameters and forcings */
#define b parms[0]
#define c parms[1]
#define d parms[2]
#define e parms[3]
#define f parms[4]
#define g parms[5]
#define import forc[0] 

/* initializers: */
void odec(void (* odeparms)(int *, double *))
{
int N=6;
odeparms(&N, parms);
}

void forcc(void (* odeforcs)(int *, double *))
{
int N=1;
odeforcs(&N, forc);
}

/* derivative function */
void derivsc(int *neq, double *t, double *y, double *ydot,double *yout, int*ip)
{
if (ip[0] <2) error("nout should be at least 2");
ydot[0] = import - b*y[0]*y[1] + g*y[2];
ydot[1] =          c*y[0]*y[1] - d*y[2]*y[1];
ydot[2] =          e*y[1]*y[2] - f*y[2];
yout[0] = y[0] + y[1] + y[2];
yout[1] = import;
}

Thanks

V. J. Ontiveros
  • 21
  • 2

1 Answers1

1

Oh, it was so easy to solve. I just defined a second forcing and actualize the counters. I modified the previous example here:

static double parms[6];
static double forc[2]; 

/* A trick to keep up with the parameters and forcings */
#define b parms[0]
#define c parms[1]
#define d parms[2]
#define e parms[3]
#define f parms[4]
#define g parms[5]
#define import forc[0] 
#define import2 forc[1] 

/* initializers: */
void odec(void (* odeparms)(int *, double *))
{
int N=6;
odeparms(&N, parms);
}

void forcc(void (* odeforcs)(int *, double *))
{
int N=2;
odeforcs(&N, forc);
}

/* derivative function */
void derivsc(int *neq, double *t, double *y, double *ydot,double *yout, int*ip)
{
if (ip[0] <2) error("nout should be at least 2");
ydot[0] = import2 - b*y[0]*y[1] + g*y[2];
ydot[1] =          c*y[0]*y[1] - d*y[2]*y[1];
ydot[2] =          e*y[1]*y[2] - f*y[2];
yout[0] = y[0] + y[1] + y[2];
yout[1] = import;
}
V. J. Ontiveros
  • 21
  • 2
  • I'm trying to do a similar thing where I have a different forcing parameter for each patch in a model. Would you mind sharing your R code you used for this to provide multiple forcing parameters to the "ode" function? Anytime I try to provide anything with more than one time column and one forcing column, I get an error, so I'm not sure what the format should be – m.evans Dec 05 '22 at 12:00