C extension using GSL

[Nick Craig-Wood]

jesse <jberwald at gmail.com> wrote:
>  I give up. I cannot find my memory leak! I'm hoping that someone out
>  there has come across something similar. Let me lay out the basic
>  setup:
> 
>  I'm performing multiple simulations on a model. Each iteration
>  involves solving a system of differential equations. For this I use
>  the GNU Scientific Library (GSL) -- specifically the rk4imp ODE
>  solver. After the ODE is solved the array is returned to python and is
>  analyzed. As you may have guessed, my problem is that over the course
>  of the iterations the memory keeps climbing until python crashes.
> 
>  Note: *The extension does not keep running.* It returns object (a
>  list) and is done until the next iteration. AFAIK, any memory
>  allocated during execution of the extension should be released.
>  Question: Since the extension was run from within python is memory
>  allocated within an extension part of python's heap?

No, on the normal C heap

>  Would this have
>  an adverse or unpredictable affect on any memory allocated during the
>  running of the extension?

If the library passes you data and ownership of a heap block, you'll
need to free it.


> One hypothesis I have, since most other
>  possibilities have been eliminated, is that GSL's creation of it's own
>  data structures (eg., gsl_vector) is messing up Python's control of
>  the heap. Is this possible?

Only if GSL is buggy

> If so, how would one be able to fix this
>  issue?

Valgrind
 
>  It may help some nice folks out there who are good enough to look at
>  this post if I layout the flow, broken up by where stuff happens
>  (i.e., in the Python code or C code):
> 
>  1) Python: Set up the simulation and basic data structures (numpy
>  arrays with initial conditions, coupling matrices for the ODE's, dicts
>  with parameters, etc).
>  2) Python: Pass these to the C extension
>  3) C: Python objects passed in are converted to C arrays, floats,
>  etc.
>  4) C: A PyList object, L,  is created (new reference!). This will hold
>  the solution vector for the ODE
>  5) C: Initialize GSL ODE solver and stepper functions. Solve the ODE,
>  at each step use PyList_Append(L, current_state) to append the current
>  state to L.
>  6) C: After the ODE solver finishes, free GSL objects, free coupling
>  matrices, free last state vector.
>  7) C: Return L to Python with return Py_BuildValue("N", L).
>  8) Python: Convert returned list L to array A, delete L, work with A.
>  8.1) Python: Step 8) includes plotting. (I have a small suspicion that
>  matplotlib holds onto lots of data, but I've used clf() and close(fig)
>  on all plots, so I think I'm safe here. )
>  8.2) Python: save analysis results from A, save A. (At this point
>  there should be no more use of A. In fact, at point 8) in the next
>  iteration A is replaced by a new array.)
>  9) Python: Change any parameters or initial conditions and goto 1).

At every point a memory allocation is made check who owns the memory
and that it is freed through all possible error paths.

Valgrind will help you find the memory leaks.  It works well once
you've jumped the flaming hoops of fire that setting it up is!

Another thing you can try is run your process untill it leaks loads,
then make it dump core.  Examine the core dump with a hex editor and
see what it is full of!  This technique works suprisingly often.

-- 
Nick Craig-Wood <nick at craig-wood.com> -- http://www.craig-wood.com/nick