  The harmonic approximation of ionic fluctuations and the linear coupling
between phonons and electrons provide the basic approach to compute, from first
principles, the contribution that the nuclei dynamics and its interaction with
electrons have on materials' properties. These approaches are questionable, at
least, whenever quantum and anharmonic effects on the nuclei vibrational
properties are large, such as in hydrogenous systems, high-$T_c$
superconductors, and when systems are close to displacive phase transitions.
Here we propose a novel non-perturbative approach to compute the
electron-phonon interaction from first principles that includes non-linear
effects and takes into account the quantum nature of nuclei. The method is
based on the $GW^{(en)}$ approximation for the electron self-energy, given by
the effective nuclei-mediated electron-electron interaction $W^{(en)}$ and the
electron Green's function $G$. The electrons are treated at a mean-field level
and the nuclei dynamics is described with a Gaussian distribution function,
which can effectively take into account anharmonic effects at a mean-field
level, e.g. within the self-consistent harmonic approximation. Pivotal
quantities of the Gaussian $GW^{(en)}$ self-energy are the renormalized average
vertices, which are computed in supercells with a stochastic approach, using
the the self-consistent electronic potential computed for different atomic
configurations. In order to validate the method, $GW^{(en)}$ calculations are
performed on aluminum, a highly harmonic system with weak electron-phonon
coupling. The obtained results coincide with the ones obtained with standard
linear electron-phonon calculations. However, calculations performed on
palladium hydride, a very anharmonic system, show a highly non-linear
electron-phonon interaction, with $GW^{(en)}$ bringing corrections as high as
the linear-order result.
