[SciPy-User] Unexpectedly large memory usage in scipy.ode class

Per Nielsen evilper@gmail....
Mon Jun 3 08:35:35 CDT 2013


You are right, I checked the size of the working arrays and they
corresponded perfectly with my memory usage.

I checked some of the Runge-Kutta based solvers made available through the
complex_ode class, and they have a lower memory overhead compared to zvode,
about half according to my tests.

Thank you for your help :)

Per


On Fri, May 31, 2013 at 4:23 PM, Warren Weckesser <
warren.weckesser@gmail.com> wrote:

>
> On Fri, May 31, 2013 at 5:58 AM, Per Nielsen <evilper@gmail.com> wrote:
>
>> Hi all,
>>
>> I am solving large linear ODE systems using the QuTip python package (
>> https://code.google.com/p/qutip/) which uses scipy ODE solvers under the
>> hood. The system is of he form
>>
>> dydt = L*y,
>>
>> where L is a large complex sparse matrix, all pretty standard. In this
>> type of problem the matrix L is the biggest memory user, expected to be
>> much larger than the solution vector y itself.
>>
>> Below is the output of a @profile from the memory_profiler package on the
>> function setting up the ode object, no actual time-stepping is done (the
>> code can be found here:
>> https://github.com/qutip/qutip/blob/master/qutip/mesolve.py#L561).
>>
>> Line #    Mem usage    Increment   Line Contents
>> ================================================
>>    562                             @profile
>>    563                             def _mesolve_const(H, rho0, tlist,
>> c_op_list, expt_ops, args, opt,
>>    564                                                progress_bar):
>>    565                                 """!
>>    566                                 Evolve the density matrix using an
>> ODE solver, for constant hamiltonian
>>    567                                 and collapse operators.
>>    568                                 """
>>    569    61.961 MB     0.000 MB
>>    570    61.961 MB     0.000 MB       if debug:
>>     571                                     print(inspect.stack()[0][3])
>>    572
>>    573                                 #
>>    574                                 # check initial state
>>    575                                 #
>>    576    61.961 MB     0.000 MB       if isket(rho0):
>>    577                                     # if initial state is a ket
>> and no collapse operator where given,
>>    578                                     # fallback on the unitary
>> schrodinger equation solver
>>    579    61.961 MB     0.000 MB           if len(c_op_list) == 0 and
>> isoper(H):
>>     580                                         return _sesolve_const(H,
>> rho0, tlist, expt_ops, args, opt)
>>    581
>>    582                                     # Got a wave function as
>> initial state: convert to density matrix.
>>    583    61.973 MB     0.012 MB           rho0 = rho0 * rho0.dag()
>>     584
>>    585                                 #
>>    586                                 # construct liouvillian
>>    587                                 #
>>    588    61.973 MB     0.000 MB       if opt.tidy:
>>     589    61.973 MB     0.000 MB           H = H.tidyup(opt.atol)
>>    590
>>    591   327.887 MB   265.914 MB       L = liouvillian_fast(H, c_op_list)
>>    592
>>    593                                 #
>>    594                                 # setup integrator
>>    595                                 #
>>    596   343.168 MB    15.281 MB       initial_vector =
>> mat2vec(rho0.full())
>>    597   343.168 MB     0.000 MB       r = scipy.integrate.ode(cy_ode_rhs)
>>    598   343.168 MB     0.000 MB       r.set_f_params(L.data.data,
>> L.data.indices, L.data.indptr)
>>    599   343.168 MB     0.000 MB       r.set_integrator('zvode',
>> method=opt.method, order=opt.order,
>>    600   343.168 MB     0.000 MB                        atol=opt.atol,
>> rtol=opt.rtol, nsteps=opt.nsteps,
>>    601   343.168 MB     0.000 MB
>>  first_step=opt.first_step, min_step=opt.min_step,
>>    602   343.172 MB     0.004 MB
>>  max_step=opt.max_step)
>>    603   572.055 MB   228.883 MB
>> r.set_initial_value(initial_vector, tlist[0])
>>    604
>>    605                                 #
>>    606                                 # call generic ODE code
>>    607                                 #
>>    608   602.805 MB    30.750 MB       return _generic_ode_solve(r, rho0,
>> tlist, expt_ops, opt, progress_bar)
>>
>> On line 591 the L matrix generated and eats a large chunk of memory, as
>> expected. However, on line 603 setting the initial condition eats an almost
>> comparable chunk, despite the fact that the initial vector itself only
>> takes up ~ 15 MB (line 596).
>>
>> I find this strange, as I would expect that setting the initial condition
>> would at most increase the memory usage by approximately the size of the
>> initial vector.
>>
>> I have tried to reproduce the problem using a minimal script (see
>> attachment), but here the memory usage is as expected:
>>
>> Filename: test_ode2.py
>>
>> Line #    Mem usage    Increment   Line Contents
>> ================================================
>>      7                             @profile
>>      8    18.707 MB     0.000 MB   def runode():
>>      9    18.707 MB     0.000 MB       N = 5000
>>     10
>>     11                                 # M = np.random.rand(N, N)
>>     12   111.230 MB    92.523 MB       M = sparse.rand(N, N,
>> density=0.05, format='csr') \
>>     13   198.797 MB    87.566 MB           + 1j * sparse.rand(N, N,
>> density=0.05, format='csr')
>>     14   199.031 MB     0.234 MB       y0 = np.random.rand(N, 1) + 1j *
>> np.random.rand(N, 1)
>>     15
>>     16   199.031 MB     0.000 MB       t0 = 0.0
>>     17
>>     18   199.031 MB     0.000 MB       def f(t, y, M):
>>     19                                     # return np.dot(M, y)
>>     20                                     return M.dot(y)
>>     21
>>     22   199.031 MB     0.000 MB       r = ode(f)
>>      23   199.031 MB     0.000 MB       r.set_integrator('zvode',
>> atol=1e-10)
>>     24   199.035 MB     0.004 MB       r.set_f_params(M)
>>     25   199.035 MB     0.000 MB       r.set_initial_value(y0, t0)
>>
>> Does someone with more insight into the scipy.ode solver might have an
>> idea of whats going on? I looked in the file myself but didnt not see any
>> indications of large memory consumptions.
>>
>
>
> The `set_initial_value` method calls the integrator's `reset` method.  The
> `reset` method of the 'zvode' integrator allocates three "work" arrays,
> `iwork`, `rwork` and `zwork`, whose sizes depends on the size of `y0`.  To
> verify that these are the cause of the memory growth, you can access these
> arrays after calling `r.set_initial_value(y0, t0)` as
> `r._integrator.iwork`, etc.
>
> Warren
>
>
>
>> Best,
>> Per
>>
>>
>> _______________________________________________
>> SciPy-User mailing list
>> SciPy-User@scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-user
>>
>>
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User@scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mail.scipy.org/pipermail/scipy-user/attachments/20130603/ffbeb1d3/attachment-0001.html 


More information about the SciPy-User mailing list