Possible incompatibility with joblib?

[juanZaragozaChichell]

Hello there.

I am working on a collision detection protocol and using this wonderful library for such purpose. My collision detection algorithm could work, a priori, in parallel and for that matter I employ joblib. When evaluated sequentially, it works perfectly fine, although a bit slow when then number of points to evaluate is kind of high, hence the will to parallelize.

The following piece of code is a minimum working example in sequential execution:

` import trimesh import numpy as np from joblib import Parallel, delayed import matplotlib.tri as mtri

def check_plane_side(point_to_check, normal_vector, point_in_plane): '''Evaluates in which side of a plane does a point lie'''

evaluates <v, p-a>, where v is the normal vector to a plane

# a is a point in the plane and p is the point we want to check
# since we want for multiple points, matrix implementation follows
return np.dot(normal_vector, (point_to_check - point_in_plane).T)

def fsurface(X):

parametric definition of the surface

return [X[0], X[1], X[0]*X[1]/18 + np.cos(X[0]*X[1])/9 - np.sin(X[0]*X[1])/9]

R = 0.25 # tool radius Range = [[-2,2], [-2,2]] # parametric space of the surface

triangulation of surface

n_points_u, n_points_v = 50,50 u, v = np.linspace(Range[0][0], Range[0][1], n_points_u), np.linspace(Range[1][0], Range[1][1], n_points_v) u, v = np.meshgrid(u, v) u, v = u.flatten(), v.flatten() x,y,z = fsurface([u,v]) tri = mtri.Triangulation(u, v)

surface as trimesh object

surface_trimesh = trimesh.Trimesh(vertices = np.array([x,y,z]).T, faces = tri.triangles)

ruled1 = np.array([[[ 1.34040191, 2.11004573, 0.19064014], [ 0.87044755, 1.63201049, 0.64355706]],

   [[ 0.84066793,  1.34625883,  0.15207392],
    [ 0.63746655,  1.0448646 ,  0.87483198]],

   [[ 0.16907399,  0.59061713,  0.24413452],
    [-0.05660804,  0.21906928,  0.92643087]],

   [[-0.53801992, -0.17252205,  0.28301347],
    [-0.83689162, -0.69098233,  0.82742429]],

   [[-1.29800817, -0.89476396,  0.21490749],
    [-1.70053631, -1.49566996,  0.57738727]]])

ruled2 = np.array([[[ 1.71080332, 2.04715708, 0.34155765], [ 1.13714579, 1.51165944, 0.53819905]],

   [[ 1.17063421,  1.23534972,  0.13437579],
    [ 0.92275571,  0.92828212,  0.84061529]],

   [[ 0.39273099,  0.39246898,  0.24465826],
    [ 0.15541006,  0.01755566,  0.9211398 ]],

   [[-0.41486873, -0.47868044,  0.28325126],
    [-0.73224307, -1.01866843,  0.79530105]],

   [[-1.24853497, -1.29356539,  0.18388341],
    [-1.64838651, -1.87570714,  0.57852026]]])

ruled3 = np.array([[[ 2.03525707, 1.75931073, 0.37406457], [ 1.45542057, 1.2073419 , 0.490806 ]],

   [[ 1.47199465,  1.08236089,  0.13496352],
    [ 1.17868662,  0.76144528,  0.81723331]],

   [[ 0.73171729,  0.31953616,  0.24040671],
    [ 0.48548851, -0.03991204,  0.92207882]],

   [[-0.01314475, -0.48703867,  0.28530257],
    [-0.30182164, -0.98488722,  0.85391576]],

   [[-0.7831714 , -1.30347988,  0.23025359],
    [-1.18365404, -1.90891723,  0.5874199 ]]])

ruled_surfaces = np.array([ruled1, ruled2, ruled3]) def optimal_shank_non_colliding(shank): '''Function that we wish to parallelize'''

penetration_data = {}
penetration_data['is_shank_collision'] = False
penetration_data['larges_penetration_point'] = None
penetration_data['footpoint'] = None
penetration_data['distance'] = None

p0,pN = shank[0], shank[-1]
footpoints, distances, triangles = surface_trimesh.nearest.on_surface(shank)
footpoint, distance, triangle = [_[-1] for _ in [footpoints, distances, triangles]]
## we enter the loop and only exit when save distance is considerably small
while distance - R > 1e-3:
    save_distance = distance - R
    µ = save_distance/np.linalg.norm(p0-pN)
    p = (1-µ)*pN + µ*p0
    pN = p
    footpoints, distances, triangles = surface_trimesh.nearest.on_surface([p0,pN])
    footpoint, distance, triangle = [_[-1] for _ in [footpoints, distances, triangles]]

if (np.linalg.norm(pN-p0) >= R) or (check_plane_side(footpoint, surface_trimesh.face_normals[triangle], p0) >= 0):
    penetration_data['is_shank_collision'] = True
    penetration_data['larges_penetration_point'] = pN
    penetration_data['footpoint'] = footpoint
    penetration_data['distance'] = np.sign(check_plane_side(footpoint, surface_trimesh.face_normals[triangle], p0))*distance
else:
    pass

return penetration_data

all_shanks = ruled_surfaces.reshape(-1,2,3) results = [] for shank in all_shanks: results.append(optimal_shank_non_colliding(shank)) Now, when parallelizing, the following error is returned: results = Parallel(n_jobs=-1)(delayed(optimal_shank_non_colliding)(shank) for shank in all_shanks)

Out[]

_RemoteTraceback Traceback (most recent call last) _RemoteTraceback: """ Traceback (most recent call last): File "/Library/Frameworks/Python.framework/Versions/3.11/lib/python3.11/site-packages/joblib/externals/loky/process_executor.py", line 428, in _process_worker r = call_item() ^^^^^^^^^^^ File "/Library/Frameworks/Python.framework/Versions/3.11/lib/python3.11/site-packages/joblib/externals/loky/process_executor.py", line 275, in call return self.fn(*self.args, self.kwargs) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Library/Frameworks/Python.framework/Versions/3.11/lib/python3.11/site-packages/joblib/_parallel_backends.py", line 620, in call return self.func(*args, *kwargs) ^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Library/Frameworks/Python.framework/Versions/3.11/lib/python3.11/site-packages/joblib/parallel.py", line 288, in call return [func(args, kwargs) ^^^^^^^^^^^^^^^^^^^^^^ File "/Library/Frameworks/Python.framework/Versions/3.11/lib/python3.11/site-packages/joblib/parallel.py", line 288, in return [func(*args, *kwargs) ^^^^^^^^^^^^^^^^^^^^^ File "/var/folders/vg/0lr4vt5527zffw8t0hxyk_n80000gp/T/ipykernel_334/2268496250.py", line 10, in optimal_shank_non_colliding File "/Library/Frameworks/Python.framework/Versions/3.11/lib/python3.11/site-packages/trimesh/constants.py", line 146, in timed result = method(args, **kwargs) ^^^^^^^^^^^^^^^^^^^^^^^ File "/Library/Frameworks/Python.framework/Versions/3.11/lib/python3.11/site-packages/trimesh/proximity.py", line 332, in on_surface ... 402 finally: 403 # Break a reference cycle with the exception in self._exception 404 self = None

IndexError: arrays used as indices must be of integer (or boolean) type `

Which has led me to think that there might be some incompatibility. Any info about this or a suggestion that might work?

Thank you for your nice work on this library!!

[mikedh]

Hey, it's hard to tell without a full traceback but it's probable that a c-extension (I would suspect rtree) doesn't like the stuff joblib is doing for some reason.

For what it's worth I've had pretty good luck getting a core-X speedup with ProcessPoolExecutor.map which is built in to the standard library:

import trimesh
import numpy as np
from concurrent.futures import ProcessPoolExecutor

def do_work(mesh: trimesh.Trimesh) -> float:
    # an example function doing work on a mesh
    return np.prod(mesh.extents)

if __name__ == '__main__':
    meshes = [trimesh.creation.box().permutate.transform()
              for _ in range(100)]

    # by default will use all available cores
    with ProcessPoolExecutor() as pool:
        # run the function on multiple cores
        work = list(pool.map(do_work, meshes))