[R-sig-dyn-mod] running C code with multiple - plain text
Hi Andras, nice puzzle ;-) The "error" in the code has nothing to do with the error message. In the argument "outname" you just forgot putting c() around your list of names. Omitting this, every name was interpreted as another argument in the chain of possible arguments to lsoda(). Best regards, Daniel
On Sun, 2017-03-12 at 18:51 +0000, Andras Farkas wrote:
Dear All,
I have two questions about an error message I received when runnin
the following deSolve routine in C compiled code.?
#C Code?
/* file modelFME3.c */?
#include <R.h>?
static double parms[9];?
static double forc[5];?
#define K parms[0]?
#define V parms[1]?
#define K12 parms[2]?
#define K21 parms[3]?
#define K13 parms[4]?
#define K31 parms[5]?
#define VP parms[6]?
#define DOSE parms[7]?
#define TAU parms[8]?
#define importr2 forc[0]?
#define importr3 forc[1]?
#define importr4 forc[2]?
#define importr5 forc[3]?
#define importr6 forc[4]?
/* initializer??*/?
void initmod(void (* odeparms)(int *, double *))?
{?
int N=9;?
odeparms(&N, parms);?
}?
void forcc(void (* odeforcs)(int *, double *))?
{?
int N=5;?
odeforcs(&N, forc);?
}?
/* Derivatives and 6 output variable */?
void derivs (int *neq, double *t, double *y, double *ydot,?
double *yout, int *ip)?
{?
double In;?
if ( fmod(*t,TAU) < 0.25 ){?
In =??DOSE/0.25;?
}?
else {?
In??=??0;?
}?
if (ip[0] <1) error("nout should be at least 1");?
ydot[0] = -(K12+K+K13*importr2+K13*importr3)*y[0]+(K21*y[1])+?
(K31*y[2]*importr2)+(K31*y[3]*importr3)+(K31*y[4]*importr4)?
+(K31*y[5]*importr5)+(K31*y[6]*importr6);?
ydot[1] = K12*y[0]-K21*y[1];?
ydot[2] = In*importr2 + K13*y[0]*importr2 - K31*y[2]*importr2;?
ydot[3] = K13*y[0]*importr3-K31*y[3]*importr3;?
ydot[4] = K13*y[0]*importr4-K31*y[4]*importr4;?
ydot[5] = K13*y[0]*importr5-K31*y[5]*importr5;?
ydot[6] = In*importr6 + K13*y[0]*importr6-K31*y[6]*importr6;?
yout[0] = y[0]/V;?
yout[1] = y[2]/VP;?
yout[2] = y[3]/VP;?
yout[3] = y[4]/VP;?
yout[4] = y[5]/VP;?
yout[5] = y[6]/VP;?
}?
/* END file modelFME3.c */?
#R code
DOSE <-45?
TAU <-6?
pars <-c(K=0.01,V=12,K12=5,K21=5,K13=5,?
K31=5,VP=2.5,DOSE=DOSE,TAU=TAU)?
t????<- seq(0, 31,by=0.01)?
rtimes<-c(0.00 , 6.00, 12.00, 18.00, 20.25, 20.25,?
24.00, 24.15, 24.50, 25.50, 27.58, 27.60, 29.50 ,29.67)?
r2<-c(1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)?
r3<-c(0 ,1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)?
r4<-c(0 ,0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)?
r5<-c(0 ,0, 0,1,1,1, 0, 0, 0, 0, 0, 0, 0, 0)?
r6<-c(0 ,0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1)?
times??<- t?
signal <- as.data.frame(list(times = times,?
importr2 = rep(0, length(times))))?
border1 <-ifelse(rtimes[sum(r2)]==max(rtimes),?
rtimes[sum(r2)],rtimes[sum(r2)+1])?
signal$importr2[signal$times <=border1] <- 1?
ftime??<- seq(0, 31, 0.01)?
sigimp <- approxfun(signal$times, signal$importr2, rule = 2)?
Sigimp <- approx(signal$times, signal$importr2,?
xout=ftime ,rule = 2)$y?
forcingsr2 <- cbind(ftime, Sigimp)?
signal <- as.data.frame(list(times = times,?
importr3 = rep(0, length(times))))?
border2 <-ifelse(rtimes[sum(r2,r3)]==max(rtimes),?
rtimes[sum(r2,r3)],rtimes[sum(r2,r3)+1])?
signal$importr3[signal$times >border1&signal$times <=border2] <- 1?
ftime??<- seq(0, 31, 0.01)?
sigimp <- approxfun(signal$times, signal$importr3, rule = 2)?
Sigimp <- approx(signal$times, signal$importr3, xout=ftime ,rule =
2)$y?
forcingsr3 <- cbind(ftime, Sigimp)?
signal <- as.data.frame(list(times = times,?
importr4 = rep(0, length(times))))?
border3 <-ifelse(rtimes[sum(r2,r3,r4)]==max(rtimes),?
rtimes[sum(r2,r3,r4)],?
rtimes[sum(r2,r3,r4)+1])?
signal$importr4[signal$times >border2&signal$times <=border3] <- 1?
ftime??<- seq(0, 31, 0.01)?
sigimp <- approxfun(signal$times, signal$importr4, rule = 2)?
Sigimp <- approx(signal$times, signal$importr4,?
xout=ftime ,rule = 2)$y?
forcingsr4 <- cbind(ftime, Sigimp)?
signal <- as.data.frame(list(times = times,?
importr5 = rep(0, length(times))))?
border4 <-ifelse(rtimes[sum(r2,r3,r4,r5)]==max(rtimes),?
rtimes[sum(r2,r3,r4,r5)],?
rtimes[sum(r2,r3,r4,r5)+1])?
signal$importr5[signal$times >border3&signal$times <=border4] <- 1?
ftime??<- seq(0, 31, 0.01)?
sigimp <- approxfun(signal$times, signal$importr5, rule = 2)?
Sigimp <- approx(signal$times, signal$importr5,?
xout=ftime ,rule = 2)$y?
forcingsr5 <- cbind(ftime, Sigimp)?
signal <- as.data.frame(list(times = times,?
importr6 = rep(0, length(times))))?
border5 <-ifelse(rtimes[sum(r2,r3,r4,r5,r6)]==max(rtimes),?
rtimes[sum(r2,r3,r4,r5,r6)],?
rtimes[sum(r2,r3,r4,r5,r6)+1])?
signal$importr6[signal$times >border4&signal$times <=border5] <- 1?
ftime??<- seq(0, 31, 0.01)?
sigimp <- approxfun(signal$times, signal$importr6, rule = 2)?
Sigimp <- approx(signal$times, signal$importr6,?
xout=ftime ,rule = 2)$y?
forcingsr6 <- cbind(ftime, Sigimp)?
forcings<-list(forcingsr2,forcingsr3,forcingsr4,?
forcingsr5,forcingsr6)?
setwd("C:...\\fme")?
system("R CMD SHLIB modelFME3.c")?
dyn.load(paste("modelFME3", .Platform$dynlib.ext, sep = ""))?
#dyn.unload(paste("modelFME3", .Platform$dynlib.ext, sep = ""))?
model <- function(pars, times=seq(0,31, by = 0.01)) {?
state <- c(a = 0,b=0,c=0,d=0,e=0,f=0,g=0)?
return(ode(y = state,?
times = times,?
func = "derivs",?
parms = pars,?
dllname = "modelFME3",?
initforc = "forcc",?
forcings=forcings,?
method="lsoda",?
rtol = 1e-8,?
atol =1e-8,?
hmax=0.01,?
initfunc = "initmod",?
nout=6,?
outnames="cp.a","p2.c","p3.d","p4.e","p5.f","p6.g"))?
}?
data <-model(pars)
once we run the code an error message will be returned as follows:
Error in lsodar(y, times, func, parms, rtol, atol, jacfunc, jactype,?
rootfunc,??: ,jactype' must be one of 'fullint', 'fullusr', 'bandusr'
or 'bandint....
1. if familiar with the error, could you please explain why the error
?message is prompted using the code above?
2. what specific suggestions would you have to correct the code
above?
to allow the estimation of the variables of interest??
Appreciate your input,?
Andras Farkas
_______________________________________________ R-sig-dynamic-models mailing list R-sig-dynamic-models at r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models