[SciPy-User] Number of equations and precision of integrate.odeint()

[Sachin Kumar]

Hi,

Is there any reason why integrate.odeint() should become less precise
when the number of equations increase?  I'm trying to solve these two
sets of differential equations:

  dy1/dx = x**2     (1)
  dy2/dx = x**3

and,

  dy1/dx = x**2     (2)

Since the two equations in (1) are uncoupled, y1 from integrating (1)
and (2) should be the same.  I'm implementing the two in the following
snippet:

  import numpy as np
  from scipy.integrate import odeint

  def f(y, x):
      return [x**2, x**3]

  def g(y, x):
      return x**2

  a = odeint(f, [0.0, 0.0], np.arange(0, 5, 0.0001))
  b = odeint(g, 0.0, np.arange(0, 5, 0.0001))

  print a[-1][0], b[-1][0], abs(a[-1][0] - b[-1][0])

Running the above code gives me:

  41.6641667161 41.6641667354 1.93298319573e-08

The difference seems to be very insignificant here.  But in cases where
the number of equations becomes large (for e.g. in Lyapunov exponent
calculations), it appears to cause significant differences.

What could be going on here?

-sachin