I have been attempting to model a moderately simple system where there are different types of particles. Particles of a certain type can change type depending on certain rules which can be represented by ODEs. I do not model their position, only the number of particles of each type at some particular locations.

The main thing is that the total number of particles must be conserved.

I was initially using scipy.integrate.odeint, but ran into the issue I describe here: http://scicomp.stackexchange.com/questions/12938/scipy-odeint-sum-of-conservative-ode-equations-does-not-remain-zero-as-it-is-be

Would PyDSTool be able to help by allowing me to perform symplectic integration with it?