You can find a fully worked-out derivation of the propagator on an elliptic curve in Appendix A of Feynman Graph Integrals and Almost Modular Forms by S. Li (2011). Basically, the propagator $P(z)$ appears upon evaluation of $\sum_{n,m=-\infty}^{\infty}(z-m-n)^{-2}$. This gives the Weierstrass elliptic function $\wp$ and Eisenstein series $E_2$ in the second equation of your posting.
You ask "What am I doing wrong?"
From what I understand, the propagator $P(z)$ is not the inverse of the Laplacian $\Delta_z$, which is the function you are computing in your posting, but instead
$$P(z)=\int_0^\infty \frac{d^2}{dz^2}e^{-t\Delta_z}\,dt.$$