1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34
|
require(PBSddesolve)
local(env=.PBSddeEnv, expr={
parms<-c(f1=1, f2=0.5, taui=5, tauj=10, v0=2, mu0=0.5)
myGrad<- function(t, y, parms){
with(as.list(c(parms)),{
if(t<taui)
lag1<-0
else
lag1<-pastvalue(t-taui)
if(t<tauj)
lag2<-0
else
lag2<-pastvalue(t-tauj)
y1 <- v0 - f1*(y[1]-y[2])
y2 <- f1*(lag1[1]- lag1[2])- f2*(y[2]-y[3])
y3 <- f2*(lag2[2]-lag2[3]) - mu0
return(c(y1, y2, y3))
})
}
yinit<-c(0.3,1,2)
yout<-dde(yinit,times=seq(0, 150, by=0.02), func=myGrad, parms=parms)
plot(yout$time, yout$y1, type="l", col="red", xlab="t", ylab="stocks y",ylim=c(min(yout$y1, yout$y2, yout$y3), max(yout$y1, yout$y2, yout$y3)))
lines(yout$time, yout$y2, type="l", col="blue")
lines(yout$time, yout$y3, type="l", col="green")
legend("topleft", legend=c("y1", "y2", "y3"), lwd=2, lty=1, xjust=1, yjust=1, col=c("red", "blue", "green"))
}) |
Partager