Calculating the distance squared between all vertex-pairs of a number of 2D polygons
It would be helpful to have a simple test harness for the code.

In reviewing your code I can see some possible bottlenecks, in Cython there can be quite a large overhead in passing arrays (views) to functions, also the array operations may be convenient but when performing them on very small arrays of only a few values performance will tend to suffer considerably compared with doing the arithmetic directly — Cython can make some optimizations when it comes to calling functions into python-space and a few numpy array operations bypass python-space — but arithmetic on basic types in Cython is guaranteed to be at the speed of c as it translates directly to c code, as does array indexing.

Hence you would probably gain some performance benefits by passing floats directly instead of arrays.

cdef double calculate_2D_dist_squared(self, np.float64_t x1, np.float64_t y1, 
                                      np.float64_t x2, np.float64_t y2):
    return (x1 - x2) ** 2 + (y1 - y2) ** 2


As noted previously, indexing is very very fast in Cython, but intermediate views can result in performance losses, hence my instinct would be to call the function in this way:

dist_squared = self.calculate_2D_dist_squared(
    polygons_vertex_coords[pi_focus, vi_focus, 0], 
    polygons_vertex_coords[pi_focus, vi_focus, 1], 
    polygons_vertex_coords[pi, vi, 0],
    polygons_vertex_coords[pi, vi, 1])


Instead of the other possible way which would involve creating a variable to hold polygons_vertex_coords[pi_focus, vi_focus] (etc).