[R-sig-dyn-mod] forcings {deSolve}
Hi Daniel, Many thanks for your detailed answer. Combining the matrices for "a" and "b" as a list was what I was missing. Thanks for your help, All the best, Tom
Daniel Kaschek <daniel.kaschek at physik.uni-freiburg.de> 04/03/2015 21:15 >>>
Hi Tom,
even though you have several forcings, you have only one function
"initforc". Consider the following example of dx/dt=-k1*x+k2*a+k3*b
where "a" and "b" are forcings. The corresponding C code would be
/** --------- Auto-generated code ------------ **/
#include <R.h>
#include <math.h>
static double parms[4];
static double forc[2];
#define k1 parms[0]
#define k2 parms[1]
#define k3 parms[2]
#define y0_0 parms[3]
#define a forc[0]
#define b forc[1]
void initmod(void (* odeparms)(int *, double *)) {
int N=4;
odeparms(&N, parms);
}
void initforc(void (* odeforcs)(int *, double *)) {
int N=2;
odeforcs(&N, forc);
}
/** Derivatives (ODE system) **/
void derivs (int *n, double *t, double *y, double *ydot, double *RPAR,
int *IPAR) {
ydot[0] = -k1*y[0]+k2*a+k3*b;
}
Then you call ode() with func = "derivs", initfunc = "initmod", initforc
= "initforc" and forcings = myforclist where myforclist may be
$a
t x
[1,] 0 1
[2,] 1 1
[3,] 2 1
[4,] 3 1
$b
t x
[1,] 0 1
[2,] 1 2
[3,] 2 3
[4,] 3 4
i.e. a list of length two. Each list entry is a matrix with two columns,
time and value. The list should be ordered according to the appearance
of forcings in the C file.
Sometimes the forcing is known as an algebraic function, e.g. an
exponential function. In this case you can define the forcing directly
in deriv():
/** Derivatives (ODE system) **/
void derivs (int *n, double *t, double *y, double *ydot, double *RPAR,
int *IPAR) {
double a = exp(-k4 * *t);
double b = exp(-k5 * *t);
ydot[0] = -k1*y[0]+k2*a+k3*b;
}
In this case add k4 and k5 to the parameters and remove a and b from the
forcings.
Hope this helps.
Cheers,
Daniel
On Mi, 2015-03-04 at 15:36 +0000, Tom Sumner wrote:
Dear list,
I'm trying to learn how to use compiled model code (specifically C) with deSolve and am struggling with passing forcing functions to the model.
I've been able to do this for a single forcing function by following the examples in the deSolve Compiled Code vignette but can't work out how to pass multiple functions - is this possible?
I think I've understand how to set up the "forcc" function in the C code but I'm not sure how to pass multiple functions in the call to ode.
Below is an example from the vignette for a single forcing function where "forcings" is a two-column matrix of times and values. If I want to pass two functions (e.g. forcings and forcings2) how do I specify these in the example below?
out <- ode(y=xstart, times, func = "derivsc",
parms = parms, dllname = "Forcing_lv",initforc = "forcc",
forcings=forcings, initfunc = "parmsc", nout = 2,
outnames = c("Sum","signal"), method = rkMethod("rk34f"))
Apologies if this is in the wrong place and thanks in advance for any help.
All the best,
Tom
_______________________________________________ R-sig-dynamic-models mailing list R-sig-dynamic-models at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models
Daniel Kaschek Institute of Physics Freiburg University Room 210 Phone: +49 761 2038531 _______________________________________________ R-sig-dynamic-models mailing list R-sig-dynamic-models at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models