In the main script, more precisely when computing each component of the force, there's an abs(rel_velo).
dragforce[i] = 0.5*coeff*math.pi * particle.data[CUBA.RADIUS]**2 * density * abs(rel_velo)*rel_velo[i]
rel_velo being a 3 dimensions array, abs(rel_velo) is a 3D array too, so instead of having a scalar component of the dragforce we have a 3D component, which makes no sense.
We assume that the it must be the norm numpy.linalg.norm(rel_velo) of the relative velocity, but we need confirmation.
We will use it as a temporary fix.