[R-sig-dyn-mod] forcings {deSolve}

[Tom Sumner]

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