Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

LinearNDInterpolator not thread safe #8856

Open
mraspaud opened this issue May 22, 2018 · 4 comments · May be fixed by #20619
Open

LinearNDInterpolator not thread safe #8856

mraspaud opened this issue May 22, 2018 · 4 comments · May be fixed by #20619
Labels
defect A clear bug or issue that prevents SciPy from being installed or used as expected scipy.interpolate
Projects
scipy.interpolate

Comments

@mraspaud
Copy link

LinearNDInterpolator is crashing with a segfault when an instance is being called concurrently.
Steps to reproduce:

  • create a LinearNDInterpolator instance
  • run the call method of the instance concurrently

I looked at the source code, but can't see anything at the python level that could cause it to not be thread safe. Can it be deeper down in the fortran routines it's calling ? Could it be related to #4763 and #4770 ?

Reproducing code example:

import threading
import numpy as np
from scipy.interpolate.interpnd import LinearNDInterpolator

r_ticks = np.arange(0, 5000, 10)
phi_ticks = np.arange(0, 5000, 10)
r_grid, phi_grid = np.meshgrid(r_ticks, phi_ticks)

def do_interp(interpolator, slice_rows, slice_cols):
    print('interpolation started')
    grid_x, grid_y = np.mgrid[slice_rows, slice_cols]
    res = interpolator((grid_x, grid_y))
    print('interpolation complete')
    return res

points = np.vstack((r_grid.ravel(), phi_grid.ravel())).T
values = (r_grid * phi_grid).ravel()
interpolator = LinearNDInterpolator(points, values)

worker_thread_1 = threading.Thread(target=do_interp, 
                                   args=(interpolator, slice(0, 2500), slice(0, 2500)))
worker_thread_2 = threading.Thread(target=do_interp, 
                                   args=(interpolator, slice(2500, 5000), slice(0, 2500)))
worker_thread_3 = threading.Thread(target=do_interp,
                                   args=(interpolator, slice(0, 2500), slice(2500, 5000)))
worker_thread_4 = threading.Thread(target=do_interp,
                                   args=(interpolator, slice(2500, 5000), slice(2500, 5000)))


worker_thread_1.start()
worker_thread_2.start()
worker_thread_3.start()
worker_thread_4.start()

worker_thread_1.join()
worker_thread_2.join()
worker_thread_3.join()
worker_thread_4.join()

Error message:

[2]    536 segmentation fault  python test_ts_interp.py

Scipy/Numpy/Python version information:

Tested on both python 2.7.5 and python 3.6.3

('1.1.0', '1.14.3', sys.version_info(major=2, minor=7, micro=5, releaselevel='final', serial=0))
(1.1.0 1.14.3 sys.version_info(major=3, minor=6, micro=3, releaselevel='final', serial=0)
@pv
Copy link
Member

pv commented May 22, 2018 via email

@mraspaud
Copy link
Author

@pv thanks for the workaround!

@ev-br ev-br added defect A clear bug or issue that prevents SciPy from being installed or used as expected scipy.interpolate labels Sep 2, 2018
@knutfrode knutfrode mentioned this issue Jan 17, 2019
Closed
@tjof2 tjof2 mentioned this issue Apr 11, 2020
Closed
@ev-br ev-br added this to N-D scattered in scipy.interpolate Apr 29, 2022
@ev-br
Copy link
Member

ev-br commented Oct 20, 2022

Having looked at this issue for a while, I still fail to see what is being cached on the first call and where. One way to make the OP crash go away is to hold the GIL in the loop over the data points, https://github.com/scipy/scipy/blob/main/scipy/interpolate/interpnd.pyx#L314
Inside that loop, the main suspect is a call to qhull._find_simplex, the rest looks benign.

So if indeed the culprit is somwhere in Qhull, I actually wonder if it is really correct to call from with nogil or without specialized locks, since the Qhull website has the note (http://qhull.org/html/qh-code.htm#reentrant):

Note: Reentrant Qhull is not thread safe. Do not invoke Qhull routines with the same qhT* pointer from multiple threads.

Concluding, ISTM locks need to be added somewhere in Qhull itself or our Cython wrappers. To properly fix this issue, it needs a look from spatial/qhull expert. @tylerjereddy if memory serves, you worked on scipy.spatial, any chance you have a pointer off the top of your head?

@Naikless
Copy link

Naikless commented Sep 7, 2023 •
edited

On first use, the interpolator computes some additional data structures, which are cached on the interpolator object. This step is not threadsafe (it should be, need to add locks around this). The crash goes away if you add interpolator((0, 0)) before the threaded calls.

This is actually very crucial information if you want to use this interpolator in any form of parallel computing. For example, if you don't "initialize" it this way when using imap from the multiprocessing library , each new iteration will suffer from the significantly increased computation time on the first call. I suppose, since this cache is saved on the object, you have to ensure that it already exists before shipping the object to any subprocesses.

This should probably be mentioned in the docs.

@Timothy-Anders0n Timothy-Anders0n mentioned this issue Dec 2, 2023
Open
@andfoy andfoy linked a pull request Apr 30, 2024 that will close this issue
Open
athob added a commit to athob/py-ananke that referenced this issue May 16, 2024
@athob
Investigating threading issue - using scipy/scipy#8856 (comment)
21a4c3a
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
defect A clear bug or issue that prevents SciPy from being installed or used as expected scipy.interpolate
Projects
scipy.interpolate
N-D scattered
Milestone
No milestone
Development

Successfully merging a pull request may close this issue.

ENH: spatial: serialize concurrent calls to QHull andfoy/scipy
4 participants
@pv @mraspaud @ev-br @Naikless