The following mophological method can be used to extract the contours in a 3D NumPy array. The idea comes from here.
import numpy as np
from skimage import morphology
arr = np.zeros((5,5,5))
arr[1:4, 1:4, 1:4] = 1
dila = morphology.dilation(arr!=1, morphology.cube(3))
res = np.argwhere((a==1) & (dila==True)) # indices of the contours
Note that in this case, if all elements of arr is 1, the result res will be empty.