Define $\epsilon_k(x) = \epsilon(x+z_k) \chi_k(x+z_k).$ Then $$\int |\nabla \epsilon_k|^2 = \int |\nabla \epsilon|^2 \chi_k^2 - \int \epsilon^2 (\Delta \chi_k) \chi_k.$$ There is a complicated defintion of $\chi_k$ and $\epsilon$ but I don't think that it is required in proving the above identity. Here is what I have tried so far, $$|\nabla \epsilon_k|^2 = (\nabla \epsilon \chi_k + \epsilon \nabla \chi_k)^2 = |\nabla \epsilon|^2 \chi_k^2 + \epsilon^2 |\nabla \chi_k|^2 + 2\epsilon \nabla \epsilon \chi_k \nabla \chi_k.$$ I am not sure how to simplify this further. Any hints will be much appreciated.
For the sake of completeness here is the definition of $\chi_k$ for $k=1,2.$ First we define $\chi:\mathbb{R}\to \mathbb{R} $ to be a smooth function such that $$\chi = 1\text{ on } [0,1],\quad \chi = 0 \text{ on }[2,\infty], \quad \chi'\leq 0 \text{ on }\mathbb{R}.$$ Next, $$\chi_1(x) = \chi\left(\frac{|x-z_1|}{\lambda}\right),\quad \chi_2(x) = (1-\chi_1(x))^{1/2}$$ for $\lambda = |z|/4\gg 1.$
Note that $$\nabla \cdot (\epsilon^2\chi_k\nabla \chi_k)=2\epsilon\chi_k\nabla\chi_k\cdot\nabla \epsilon+\epsilon^2\nabla\chi_k\cdot\nabla\chi_k+\epsilon^2\chi_k\Delta\chi_k,$$ and so
$$|\nabla\epsilon_k|^2=|\nabla \epsilon|^2 \chi_k^2 + \epsilon^2 |\nabla \chi_k|^2 + 2\epsilon \nabla \epsilon \chi_k \nabla \chi_k=|\nabla \epsilon|^2 \chi_k^2 +\nabla \cdot (\epsilon^2\chi_k\nabla \chi_k)-\epsilon^2\chi_k\Delta\chi_k.$$ Integrating both sides and using the divergence theorem gives what you're after.