[R] deSolve: Problem solving ODE including modulo-operator

[apokaly]

I have a problem integrating the 'standard map' (
http://en.wikipedia.org/wiki/Standard_map
http://en.wikipedia.org/wiki/Standard_map ) with deSolve:

By using the modulo-operator '%%' with 2*pi in the ODEs (standardmap1), the
resulting values of P and Theta, should not be greater than 2pi. Because
this was not the case, i was thinking that the function 'ode' has some
internal problems with the '%%' or integrating periodical ODEs, so I wrote a
modulo-function by myself (modulo and standardmap2). But still I get values
much higher than 2pi and I cannot find the error... Any guess?

Thanks


Full code:
----------------------------------------------------------------
library(deSolve)


iterations <- 100

Parameter <- c(k = 0.6)
State <- c(Theta = 1 ,  P = 1) 

Time <- 0:iterations/10




standardmap1 <- function(Time, State, Parameter){ 
  with(as.list(c(State, Parameter)), {
	dP  <- (P + k * sin(Theta)) %% (2 * pi)
	dTheta <- (P + Theta) %% (2 * pi)
return(list(c(dP, dTheta)))
})
}

#solve ode using standardmap1
out1 <- as.data.frame(ode(func = standardmap1, y = State,parms = Parameter,
times = Time))




# x mod y, end: maximal iterations
modulo <- function(x,y,end=1000){
	for (n in 0:end)
		if (x > (n-1)*y) 
			z = x - (n-1)*y
			if (z > 0)
				return(z)
			else break
}



standardmap2 <- function(Time, State, Parameter){
  with(as.list(c(State, Parameter)), {
	dTheta <- modulo((P + Theta),(2*pi))
	dP <- modulo(P + k *sin(Theta),(2*pi))
return(list(c(dP, dTheta)))
})
}



#solve ode using standardmap2
out2 <- as.data.frame(ode(func = standardmap2, y = State,parms = Parameter,
times = Time))

#plot the results
matplot(out1[2],out1[3], type = "p", pch = ".")
matplot(out2[2],out2[3], type = "p", pch = ".")
-- 
View this message in context: http://r.789695.n4.nabble.com/deSolve-Problem-solving-ODE-including-modulo-operator-tp3236396p3236396.html
Sent from the R help mailing list archive at Nabble.com.