  In solid state physics, the electron-phonon interaction (EPI) is central to
many phenomena. The theory of the renormalization of electronic properties due
to EPIs became well established with the theory of Allen-Heine-Cardona, usually
applied to second order in perturbation theory (P2). However, this is only
valid in the weak coupling regime, while strong EPIs have been reported in many
materials. Although non-perturbative (NP) methods have started to arise in the
last years, they are usually not well justified, and it is not clear to what
degree they reproduce the exact theory. To address this issue, we present a
stochastic approach for the evaluation of the non-perturbative interacting
Green's function in the adiabatic limit, and show it is equivalent to the
Feynman expansion to all orders in the perturbation. Also, by defining a
self-energy, we can reduce the effect of broadening needed in numerical
calculations, improving convergence in the supercell size. In addition, we
clarify whether it is better to average the Green's function or self-energy.
Then we apply the method to a graphene tight-binding model, and obtain several
interesting results: (i) The Debye-Waller term, which is normally neglected,
does affect the change of the Fermi velocity. (ii) The P2 and NP self-energies
differ even at room temperature for some k-points, raising the question of how
well P2 works in other materials. (iii) Close to the Dirac point, positive and
negative energy peaks merge. (iv) In the strong coupling regime, a peak appears
at energy E=0, which is consistent with previous works on disorder and
localization in graphene. (v) The spectral function becomes more asymmetric at
stronger coupling and higher temperatures. Finally, in the Appendix we show
that the method has better convergent properties when the coupling is strong
relative to when it is weak, and discuss other technical aspects.
