Polaron Mass and Electron-Phonon Correlations in the Holstein Model

The Holstein Molecular Crystal Model is investigated by a strong coupling perturbative method which, unlike the standard Lang-Firsov approach, accounts for retardation effects due to the spreading of the polaron size. The effective mass is calculated to the second perturbative order in any lattice dimensionality for a broad range of (anti)adiabatic regimes and electron-phonon couplings. The crossover from a large to a small polaron state is found in all dimensionalities for adiabatic and intermediate adiabatic regimes. The phonon dispersion largely smooths such crossover which is signalled by polaron mass enhancement and on site localization of the correlation function. The notion of self-trapping together with the conditions for the existence of light polarons, mainly in two- and three-dimensions, are discussed. By the imaginary time path integral formalism I show how non local electron-phonon correlations, due to dispersive phonons, renormalize downwards the {\it e-ph} coupling justifying the possibility for light and essentially small 2D Holstein polarons.


Introduction
The interest for phonons and lattice distortions in High Temperature Superconductors (HTSc) is today more than alive [1,2]. While the microscopic origin of the pairing mechanism in cuprate HTSc has not yet been unravelled, evidence has been provided [3,4,5,6,7,8] that electron-lattice interactions and local lattice fluctuations are correlated with the onset of the superconducting transition [9,10,11]. Recent investigations [12] would suggest a similar role for the lattice also in the newly discovered layered pnictide-oxyde quaternary superconducting compounds [13,14,15]. In fact, the discovery of HTSc in copper oxides [16] was motivated by the knowledge that the Jahn-Teller effect [17,18,19] is strong in Cu 2+ . In nonlinear molecules a lattice distortion lifts the degeneracy of the electronic states and lowers the overall ground state energy: the energy gain is the Jahn-Teller stabilization energy which can be measured optically. As vibrational and electronic energies are of the same order, the nuclear motion cannot be separated from the electrons and the combined electron-lattice object becomes a mobile Jahn-Teller polaron [20] which can be described in terms of the nonlinear Holstein Hamiltonian [21]. Over the last twenty years, the focus on the HTSc has largely contributed to trigger the study of the polaron properties [22,23,24,25]. Earlier and more generally, polarons had become a significant branch in condensed matter physics [26,27,28,29] as Landau introduced the concept of an electron which can be trapped by digging its own hole in an ionic crystal [30]. A strong coupling of the electron to its own lattice deformation implies violation of the Migdal theorem [31] and polaron collapse of the electron band [32] with the appearance of time retarded interactions in the system. A great advance in the field came after Feynman used the path integrals [33] to calculate energy and effective mass of the Fröhlich polaron [34] by a variational method [35]. In the statistical path integral, the quantum mechanical electron motion is represented by a real space path x(τ ) being function of an imaginary time τ whose range has an upper bound given by the inverse temperature. As a key feature of the Feynman treatment, the phonon degrees of freedom can be integrated out exactly (being the action quadratic and linear in the oscillator coordinates) but leave a substantial trace as a time retarded potential naturally emerges in the exact integral action. The electron, at any imaginary time, interacts with its position at a past time. This self interaction mirrors the fact that the lattice needs some time to readjust to the deformation induced by the electron motion and follow it. For decades, the Fröhlich polaron problem has been extensively treated by path integral techniques [36,37] which, remarkably, can be applied for any value of the coupling constant [38]. A review work on the Fröhlich polaron is in Ref. [39].
With regard to small polarons, path integral investigations started with the groundbreaking numerical work by De Raedt and Lagendijk [40,41] who derived fundamental properties for the Holstein polaron. While a sizable electron-phonon coupling is a requisite for polaron formation [42,43,44], also the dimensionality and degree of adiabaticity of the physical system could essentially determine the stability and behavior of the polaron states [45,46,47,48,49,50,51,52,53,54,55,56,57,58,59]. When the characteristic phonon energy is smaller than the electronic energy the system is set in the adiabatic regime. In materials such as the HTSc, colossal magnetoresistance oxides [60,61], organic molecular crystals [62], DNA molecules [63,64] and finite quantum structures [65,66], intermediate adiabatic polarons are relevant as carriers are strongly coupled to high energy optical phonons. For the HTSc systems, a path integral description had been proposed for polaron scattering by anharmonic potentials, due to lattice structure instabilities, as a possible mechanism for the non metallic behavior of the c-axis resistivity [67,68].
It has been questioned whether small (bi)polarons could indeed account for high T c due to their large effective mass [69] but, later on, it was recognized that dispersive phonons renormalize the effective coupling in the Holstein model yielding much lighter masses [70]. Also long range Fröhlich e-ph interactions provide a pairing mechanism and bipolaron mass consistent with high T c [71,72,73,74] while accounting, among others, for angleresolved photoemission spectroscopy data in cuprates HTSc [75,76].
In view of the relevance of the polaron mass issue, I review in this article some work [77] on the Holstein polaron with dispersive phonons which examines the notion of self-trapping and estimates some polaron properties in the parameter space versus lattice dimensionality. The latter is introduced in the formalism by modelling the phonon spectrum through a force constant approach which weighs the first neighbors intermolecular shell. This is done in the spirit of the work by Holstein [78] on the small polaron motion for a one dimensional narrow band system. Calculating the polaron hopping between localized states due to multiphonon scattering processes in perturbation theory, Holstein first pointed out that the phonon dispersion was in fact a vital ingredient of the theory whereas a dispersionless spectrum would yield a meaningless divergent hopping probability. Time dependent perturbative theory, the hopping integral being the perturbation, applies as long as the interaction time is shorter than the polaron lifetime on a state. Computing the time integral for the hopping probability, I had earlier shown [50] that such interaction time is shorter in higher dimensionality thus making the time dependent perturbative method more appropriate in the latter conditions. For a given dimensionality, the interaction time is also reduced by increasing the strength of the intermolecular forces hence, the phonon dispersion. Consistently also the crossover temperature between band-like motion (at low T ) and hopping motion (at high T ) [79] shifts upwards, along the T axis, in higher dimensional systems with larger coordination numbers. Hopping polaron motion is the relevant transport mechanism in several systems including DNA chains [80,81].
It is worth emphasizing that our previous investigations on the Holstein polaron and the present paper assume dimensional effects as driven by the lattice while the electron transfer integral is a scalar quantity. In a different picture, first proposed by Emin [82] for quasi 1D solids in the adiabatic limit and recently generalized by variational exact diagonalization techniques [83], the electron transfer integral is instead a vector whose components are switched on in order to treat anisotropic polaron properties. The system dimensionality is then increased by tuning the electronic subsystem parameters, the 3D case corresponding to an isotropic hopping integral, whereas the phonon spectrum is dispersionless. While the method to introduce the system dimensionality in Ref. [83] differs from ours, some trends regarding the polaron mass in 1D and 3D are apparently at variance with ours as it will be pointed out in the following discussion.
Being central to polaron studies, the effective mass behavior ultimately reflects the abovementioned property of the polaron problem: when the electron moves through the lattice, it induces and drags a lattice deformation which however does not follow instantaneously, the retardation causing a spread in the size of the quasiparticle. The electron-phonon correlation function provides a measure of the deformation around the instantaneous position of the electron and may be used to quantify the polaron size. Several refined theoretical tools, including quantum Monte Carlo calculations [40], density matrix renormalization group [84], variational methods [85,86] and exact diagonalization techniques [87,88] have been applied to the Holstein Hamiltonian to compute such polaron characteristics covering various regimes of electronphonon coupling and adiabatic ratio. While a qualitative and often quantitative [84,85] agreement among several numerical methods has emerged especially with regard to ground state polaron properties, we feel that analytical investigations may provide a useful insight also in the specific regime of intermediate adiabaticity and e-ph coupling for which the leading orders of perturbation theory traditionally fail in the weak coupling [89] or become less accurate in the strong coupling [90,91].
Extremely interesting is that range of e-ph couplings which are sufficiently strong to allow for polaron formation but not too strong to prevent the polaron from moving through the lattice. When the relevant energy scales for electrons and phonons are not well separated, such range of couplings produces a composite particle which is not trapped by its lattice deformation. It is this intermediate regime which forms the focus of the present work.
The retardation effect and its consequences are here treated by means of a variational analytical method based on the Modified Lang Firsov (MLF) transformation [92] which considerably improves the standard Lang Firsov (LF) [93] approach on which strong coupling perturbation theory (SCPT) is based. Section 2 provides the generalities of the MLF method applied to the dispersive Holstein Hamiltonian. In Section 3, the polaron mass is calculated in the (anti)adiabatic regimes covering a broad range of cases in parameter space. The spreading of the adiabatic polaron size over a few lattice sites is shown in Section 4 through computation of the electron-phonon correlations. In Section 5, I give further physical motivations for the existence of light Holstein polarons: using a path-integral method, I show how the electron-phonon coupling is renormalized downwards in momentum space by nonlocal correlations arising from the dispersive nature of the phonon spectrum. Some final remarks are in Section 6.

Modified Lang-Firsov Method
The Holstein diatomic molecular model was originally cast [21] in the form of a discrete nonlinear Schrödinger equation for electrons whose probability amplitude at a molecular site depends on the interatomic vibration coordinates [94]. The nonlinearities are tuned by the electron-phonon coupling g, whose strength drives the crossover between a large and a small polaron for a given value of the adiabatic parameter. The polaron radius is measured with respect to the lattice constant [95,96]. When the size of the lattice distortion is of the same order of (or less than) the lattice constant the polaron has a small radius and the discreteness of the lattice must be taken into account.
In second quantization the dimension dependent Holstein Hamiltonian with dispersive harmonic optical phonons reads: t is the hopping integral and the first sum is over z nearest neighbors. c † i and c i are the real space electron creation and annihilation operators at the i-site, n i (= c † i c i ) is the number operator, b † i and b i are the phonons creation and annihilation operators. b † q is the Fourier transform of b † i and ω q is the frequency of the phonon with vector momentum q. Unlike the Su-Schrieffer-Heeger Hamiltonian [97,98] a paradigmatic model in polymer physics, in Eq. (1) the electron hopping does not depend on the relative displacement between adjacent molecular sites hence, the phonon created by b † i is locally coupled to the electronic density.
The LF transformation uses a phonon basis of fixed displacements (at the electron residing site) which diagonalizes the Hamiltonian in Eq. (1) in absence of hopping. The hopping term is then treated as a perturbation [70,90]. However the standard LF approach does not account for the retardation between the electron and the lattice deformations which induces a spread in the size of the polaron. Precisely this effect becomes important for the intermediate e-ph coupling values which may be appropriate for some HTSc.
The idea underlying the MLF transformation [99] is that to consider the displacements of the oscillators at different sites around an electron in order to describe the retardation effect. For the present case of dispersive phonon the MLF transformation, applied to the Hamiltonian in Eq. (1), reads: where R i are the lattice vectors and λ q are the variational parameters which represent the shifts of the equilibrium positions of the oscillators (quantized ion vibrations) with momentum q. The conventional Lang-Firsov transformation is recovered by setting: λ q = g/ω q . Explicitly, the MLF transformed Holstein Hamiltonian in Eq. (2) is: where the polaron self-energy ǫ p is: and the polaronic hopping is: The coordination number z is twice the system dimensionality.
Looking at Eq. (3), it is clear that the unperturbed Hamiltonian H 0 can be taken as while the remaining part of the Hamiltonian (H − H 0 ) in the MLF basis is considered as the perturbation part.
The energy eigenstates of H 0 are given by where, i is the electron site and n q1 , n q2 , n q3 are the phonon occupation numbers in the phonon momentum states q 1 , q 2 , q 3 , respectively. The lowest energy eigenstate of the unperturbed Hamiltonian has no phonon excitations, i.e. n q = 0 for all q. The ground state has an energy E 0 0 = −ǫ p and is N -fold degenerate, where N is the number of sites in the system. The perturbation lifts the degeneracy and to first order in t the ground state energy of the 3D-polaron with momentum k is given by The second order correction to the ground-state energy of the polaron with momentum k is given by (2) where, in principle, intermediate states having all possible phonon numbers contribute to Eq. (9).
By minimizing the zone center ground state energy, the variational parameters λ q are obtained as: and, by Eq. (10), the one phonon matrix element between the ground state |k = 0, {n q = 0} and the first excited state 1 q , k ′ | vanishes. Then, the one phonon excitation process yields no contribution to Eq. (9). By Eq. (10) one gets the λ q 's for the 1D, 2D and 3D systems and, via Eq. (5), the narrowing of the polaron band [100].

Polaron Mass
Given the formal background, the polaron mass is calculated both for the Lang-Firsov and for the Modified Lang-Firsov method to the second order in SCPT. Generally the second order order correction i) makes a relevant contribution to the ground state energy, ii) does not affect the bandwidth, iii) introduces the mass dependence on the adiabatic parameter which would be absent in the first order. Altogether the polaron landscape introduced by the second order SCPT is much more articulated than the simple picture suggested by the first order of SCPT which nonetheless retains its validity towards the antiadiabatic limit. In first order SCPT, ground state energy, bandwidth and effective mass appear as equivalent, interchangeable properties to describe the polaron state while band narrowing and abrupt mass enhancement are coincident signatures of small polaron formation in parameter space. However, there is no compelling physical reason for such coincidence to occur as the bandwidth is in fact a zone edge quantity whereas the effective mass is a zone center quantity. Thus, by decoupling onset of the band narrowing and self trapping of the polaron mass, the second order of SPCT provides the framework for a richer polaron structure in momentum space.
As emphasized in the Introduction, dispersive phonons have a relevant role in the Holstein model: a dispersionless spectrum would in fact predict larger polaron bandwidths in lower dimensionality and yield a divergent site jump probability for the small polaron in time dependent perturbation theory [78,101]. Here, a lattice model is assumed in which first neighbors molecular sites interact through a force constants pair potential. Then, the optical phonon spectrum is given in 1D, 2D (square lattice) and 3D (simple cubic lattice) respectively, by where α is the intra-molecular force constant and γ is the inter-molecular first neighbors force constant. M is the reduced molecular mass. Thus, the intra-and intermolecular energies are ω 0 = 2α/M and ω 1 = γ/M respectively. Some care should be taken in setting the phonon energies as the second order perturbative term grows faster than the first order by increasing ω 1 [70]. Hence too large dispersions may cause a breakdown of the SCPT. In terms of ω 0 , the dimensionless parameter zt/ω 0 defines the adiabatic (zt/ω 0 > 1) and the antiadiabatic (zt/ω 0 < 1) regime. Some other choices are found in the literature with t/ω 0 > 1 or t/ω 0 > 1/4 defining the adiabatic regime. As shown below, such discrepancies may lead to significantly different interpretations of the polaron behavior in parameter space mainly for higher dimensions.
Hereafter I take ω 0 = 100meV and tune t to select a set of zt/ω 0 which sample both antiadiabatic and adiabatic regime without reaching the adiabatic limit. The dynamics of the lattice is in fact central to our investigation. A moderate to strong range of e-ph couplings is assumed (g/ω 0 1) so that the general conditions for polaron formation are fulfilled throughout the range of adiabatic ratios [22,91].  Fig. 1(a). While at very strong couplings the MLF plots converge towards the LF predictions, a remarkably different behavior between the LF and the MLF mass shows up for moderate g. The LF method overestimates the polaron mass for g ∈ [∼ 1 − 2] and mostly, it does not capture the rapid mass increase found instead in the MLF description. Note that, around the crossover, the MLF polaron mass is of order ten times the bare band mass choosing ω 1 = 60meV . Large intermolecular energies enhance the phonon spectrum thus reducing the effective masses in both methods. In the MLF method, large ω 1 tend also to smooth the mass behavior in the crossover region. Fig. 1(b) presents the case of an adiabatic regime: the discrepancies between LF and MLF plots are relevant for a broad range of e-ph couplings. The latter is more extended for larger intermolecular energies as also seen in Fig. 1(a). Mass renormalization is poor in the MLF curves up to the crossover which is clearly signalled by a sudden although continuous mass enhancement whose abruptness is significantly smoothed for the largest values of intermolecular energies. Fig. 1(c) shows a fully antiadiabatic case in which the LF and MLF plots practically overlap throughout the whole range of couplings. This confirms that the LF method is appropriate in the antiadiabatic limit which is essentially free from retardation effects. Unlike the previous cases the MLF plots are always smooth versus g indicating that no self-trapping event occurs in the antiadiabatic regime. Here the polaron does not trap as it is always small for the whole range of couplings.
It should be reminded that the concept of self-trapping traditionally indicates an abrupt transition between an infinite size state at weak e-ph couplings and a finite (small) size polaron at strong e-ph couplings. In one dimension, according to the traditional adiabatic polaron theory [102,103,104], the polaron solution is always the ground state of the system and no self-trapping occurs. Instead, in higher dimensionality a minimum coupling strength is required to form finite size polarons, hence self-trapped polarons can exist at couplings larger than that minimum.
However, these conclusions have been critically reexamined in the recent polaron literature and the same notion of polaron size has been questioned due to the complexity of the polaron quasiparticle itself [105,106]. Certainly, as a shrinking of the polaron size yields a weight increase, the polaron mass appears as the most reliable indicator of the self-trapping transition. The latter, as the results in Fig. 1 suggest, is rather a crossover essentially dependent on the degree of adiabaticity of the system and crucially shaped by the internal structure of the phonon cloud which I have modelled by tuning the intermolecular forces. In this view, self trapping events can be found also in the parameter space of 1D systems. This amounts to say that also finite size polarons can self-trap if a sudden, although continuous, change in their effective mass occurs for some values of the e-ph couplings in some portions of the intermediate/adiabatic regime. The continuity of such event follows from the fact that the ground state energy is analytic function of g for optical phonon dispersions [107] hence the possibility of phase transitions is ruled out in the Holstein model. Neither at finite temperature can phase transitions occur as the free energy is smooth for any phonon dispersion [108,109].
As fluctuations in the lattice displacements around the electron site are included in the MLF variational wavefunction, the calculated polaron mass should not display discontinuities by varying the Hamiltonian parameters through the crossover [110]. Mathematically the crossover points are selected through the simultaneous occurence of a maximum in the first logarithmic derivative and a zero in the second logarithmic derivative of the MLF polaron mass with respect to the coupling parameter. Such inflection points, corresponding to the points of most rapid increase for m * , are marked in Fig. 2 on the plots of the mass ratios computed for a broad range of (anti)adiabatic parameters both in one, two and three dimensions.
In 1D, see Fig. 2(a), the crossover occurs for g values between ∼ 1.8 − 2.3 and the corresponding self trapped polaron masses are of order ∼ 5 − 50 times the bare band mass thus suggesting that relatively light small polarons can exist in 1D molecular solids with high phonon spectrum. The onset of the self-trapping line is set at the intermediate value 2t = ω 0 and the self-trapped mass values grow versus g by increasing the degree of adiabaticity. There is no self-trapping in the fully antiadiabatic regime as the electron and the dragged phonon cloud form a compact unit also at moderate e-ph couplings. Then, the mass increase is always smooth in the antiadiabatic regime. Some significant results are found in 2D as shown in Fig. 2(b): i) for a given g and adiabaticity ratio, the 2D mass is lighter than the 1D mass and the 2D LF limit is attained at a value which is roughly one order of magnitude smaller than in 1D; ii) the crossover region is shifted upwards along the g axis with the self trapping events taking place in the range, g ∼ 2.2 − 2.6. The masses are of order ∼ 5−10 times the bare band mass; iii) the line of self-trapping events changes considerably with respect to the 1D plot. The marked curve is parabolic in 2D with an extended descending branch starting at the intermediate value 4t/ω 0 = 1; iv) in the deep adiabatic regime, the lattice dimensionality smoothens the mass increase versus g.
This effect is even more evident in 3D, see Fig. 2(c), as there are no signs of abrupt mass increase even for the largest values of the adiabatic parameter. At the crossover, 3D masses are of order ∼ 5 − 10 times the bare band mass with the self trapping points lying in the range, g ∼ 2.5 − 2.9. At very large couplings the mass ratio becomes independent of t and converges towards the LF value. In this region (and for the choice ω 1 = 60meV ) the 3D Lang-Firsov mass is one order of magnitude smaller than the 2D mass. As the coordination number grows versus dimensionality, large intermolecular forces are more effective in hardening the 3D phonon spectrum thus leading to lighter 3D polaron masses than 2D ones.
The self trapping transition appears in Fig. 2 as smoother in higher dimensions thus contradicting the trend found by previous investigations, markedly by exact diagonalization techniques [57,83] and variational calculations [111]. Two reasons may be invoked to resolve the discrepancy, the first being more physical and the second more technical. 1) From Eq. (11), it can be easily seen that: ω 2 d [q = (0, 0, 0)] − ω 2 d [q = (π, π, π)] = 2dω 2 1 . Then, for a given ω 1 , the phonon band is more disper- sive in higher dimensionality d. Hence increasing the system dimension corresponds to attribute larger weight to the intermolecular forces which ultimately smoothen the crossover as made evident in Fig. 1. 2) The works in Refs. [57,83,111] take t/ω 0 as adiabatic ratio in any dimension whereas our calculation spans the same range of zt/ω 0 values in any dimension. Thus, for instance, the fully adiabatic zt/ω 0 = 3 plots in Fig. 2 would correspond to an antiadiabatic case in 3D (t/ω 0 = 1/2) and a slightly adiabatic case in 1D (t/ω 0 = 3/2) following the definition in Refs. [57,83,111]. Consistently, the 3D antiadiabatic polaron mass is smoother than the 1D adiabatic polaron mass. Viceversa, assuming the same ratio i.e. t/ω 0 = 1/2, one should compare the rightmost plot in Fig. 2(c) with the fourth from left in Fig. 2(a): the latter is smoother than the former as the 1D case is now antiadiabatic while the 3D case is fully adiabatic.
Nonetheless, the findings displayed in Fig. 2 agree with density matrix renormalization group [84], variational Hilbert space [57] methods and weak coupling perturbation theory [54] in predicting a lighter mass in higher dimensionality.

Electron-Phonon Correlations
Within the MLF formalism one may also compute the electron-phonon correlation functions in the polaron ground state. This offers a measure of the polaron size as electron and phonons displacements can be taken at different neighbors sites. The on site χ 0 , the first neighbor site χ 1 and the second neighbor site χ 2 correlation functions are plotted in Fig. 3 against the e-ph coupling both in 1D and 2D. For g ∼ 1, there is some residual quasiparticle weight (not appreciable on this scale) associated with χ 3 and χ 4 which however tend to vanish by increasing the coupling. An adiabatic regime is assumed to point out the self-trapping event which occurs when χ 0 ∼ 1 and χ n ∼ 0 for n ≥ 1. This happens for g ∼ 2 in 1D and g ∼ 2.5 in 2D for the lowest case of ω 1 here considered. By studying the correlation functions for two values of ω 1 it is seen how the intermolecular forces smooth the crossover and extend the polaron size over a larger range of g. The case ω 1 = 60meV , which allows a comparison with Fig. 2, is intermediate between the two plots presented in Fig. 3. In such case the self trapping appears at g ∼ 2.3 in 1D in fair agreement with the corresponding inflection point in the 1D effective mass (see fifth curve from right in Fig. 2(a)). In 2D the polaron size shrinks at slightly larger g than in 1D and the on site localization occurs at g somewhat larger than the mass inflection point.
Altogether the diamond marked loci displayed in Fig. 2 indicate that strong mass renormalization is accompanied by on site (or two sites) polaron localization thus identifying the self trapping line with the formation of small polarons. Before the crossover takes place, for instance at 1 ≤ g ≤ 2 in Fig. 3(a), the 1D polaron spreads over a few (two to four) lattice sites. This is the also the range of couplings in which the polaron band shows the largest deviations [87] from the cosine-like band of standard LF perturbation theory [90] due to the importance of longer (than nearest neighbors) hopping processes. The notion of intermediate polaron seems to us as the most appropriate for such mobile and relatively light objects.

Path Integral Method
The results obtained so far can be put on sound physical bases by applying the space-time path integral method to the dispersive Holstein Hamiltonian. The method permits to incorporate the effect of the electronphonon correlations in a momentum dependent effective e-ph coupling.
The phonons operators in Eq. (1) can be generally written in terms of the isotropic displacement field u n as: Then, the e-ph term in Eq. (1) transforms as follows: The sum over n spans all n th neighbors of the R i lattice site in any dimensionality. It is clear that the phonon dispersion introduces e-ph real space correlations which would be absent in a dispersionless model. Note that Eq. (13) contains the same physics as the e-ph term in Eq. (1). Fourier transforming the atomic displacement field and taking the lattice constant |a| = 1, from Eq. (13), I get for a linear chain and a square lattice respectively: cos(n∆q x ) + cos(n∆q y ) cos m∆q x + n∆q y + cos m∆q x − n∆q y While, in principle, the sum over n should cover all the N sites in the lattice, the cutoff n * permits to monitor the behavior of the coupling term as a function of the range of the e-ph correlations. In particular, in 1D the integer n numbers the neighbors shells up to n * while in the square lattice, the n * = 1 term includes the second neighbors shell, the sum up to n * = 2 includes the fifth neighbors shells, n * = 3 covers the nineth shell and so on. Switching off the interatomic forces, ω(q) = ω 0 , one would recover from Eq. (14) a local e-ph coupling model with S d ≡ 1. As no approximation has been done at this stage Eq. (14) is general.

A. Semiclassical Holstein Model
I apply to the Holstein Hamiltonian space-time mapping techniques [112,113,114] to write the general path integral for an electron particle in a bath of dispersive phonons. The method has been used to treat also e-ph polymer models [115,116] in which the electron hopping causes a shift in the atomic displacements and, as a consequence, the vertex function depends both on the electronic and the phononic wave vector. Such shift is however absent in the Holstein model as Eq. (1) makes clear [117,118].
I introduce x(τ ) and y(τ ′ ) as the electron coordinates at the i and j lattice sites respectively, and the electronic Hamiltonian (first term in Eq. (1)) transforms into τ and τ ′ are continuous variables ∈ [0, β] in the Matsubara Green's functions formalism with β being the inverse temperature hence the electron hops are not constrained to first neighbors sites. After setting τ ′ = 0, y(0) ≡ 0, the electron operators are thermally averaged over the ground state of the Hamiltonian. As a result, the average energy per lattice site, h e (τ ) ≡ < H e (τ ) > /N , associated to electron hopping reads (in d dimensions): ν n are the fermionic Matsubara frequencies and ǫ k is the electron dispersion relation.
Consider now the e-ph term. The spatial e-ph correlations contained in Eq. (14) are mapped onto the time axis introducing the τ dependence in the displacement field: u q → u q (τ ). Assuming periodic atomic particle paths: u q (τ + β) = u q (τ ), the displacement is expanded in N F Fourier components: 2 (ℜu n ) q cos(ω n τ )− (ℑu n ) q sin(ω n τ ) (17) with ω n = 2nπ/β. I take a semiclassical version of the Holstein Hamiltonian [119], assuming that the phonon coordinates in Eq. (14) as classical variables interacting with quantum mechanical fermion operators. Such approximation may affect the thermodynamics of the system as the quantum lattice fluctuations in fact play a role mainly for intermediate values of the e-ph coupling [85,120].
Averaging Eq. (14) on the electronic ground state, the e-ph energy per lattice site is defined as Then, on the base of Eq. (17) and Eq. (18), I identify the perturbing source current [121] for the Holstein model with the τ dependent averaged e-ph Hamiltonian term: Note that the Holstein source current does not depend on the electron path coordinates. The time dependence is incorporated only in the atomic displacements. This property will allow us to disentangle phonon and electron degrees of freedom in the path integral and in the total partition function.
After these premises, one can proceed to write the general path integral for an Holstein electron in a bath of dispersive phonons. Assuming a mixed representation, the electron paths are taken in real space while the phonon paths are in momentum space. Thus, the electron path integral reads: where m is the electron mass. The perturbing current in Eq. (20) is integrated over τ using Eq. (17) and Eq. (19). The result is: where g d (q) is thus a time averaged e-ph potential. The total partition function can be derived from Eq. (20) by imposing the closure condition both on the phonons (Eq. (17)) and on the electron paths, x(β) = x(0). Using Eq. (21), I obtain: Du q and Dx are the measures of integration which normalize the kinetic terms in the phonon field and electron actions respectively [121].
The phonon degrees of freedom in Eq. (22) can be integrated out analytically yielding: being λ M = π 2 β/M . Eq. (23) is the final analytical result from which the thermodynamics of the model can be computed [118]. The exponential function in P (q) embodies the effect of the non local correlations due to the dispersive nature of the phonon spectrum. Phonon and electron contributions to the partition function are decoupled although the effective potential g d (q) carries a dependence on the electron density profile in momentum space through the function ρ q .

B. Electron-Phonon Coupling
The time (temperature) averaged e-ph potential in Eq. (21) is computed in the case of a linear chain and of a square lattice. As an example, I take the electron density profile as ρ q = ρ o cos(q) in 1D and ρ q = ρ o cos(q x ) cos(q y ) for the 2D system. Since the momentum integration runs over q i ∈ [0, π/2], ρ o represents in both cases the total electron density. More structured density functions containing longer range oscillations are not expected to change the trend of the results hereafter presented. Low energy phonon spectra parameters are here assumed, ω 0 = 20meV and ω 1 = 10meV . Setting g = 3, I take a strong Holstein coupling although the general trend of the results holds for any value of g.
In Fig. 3(a), g 1D (q)/ρ o is plotted for three choices of the cutoff n * to emphasize how the potential depends on the e-ph interaction range. The constant value of the potential obtained for ω 1 = 0 is also reported on. For short correlations (n * = 4) there is a range of wave vectors in which the effective coupling becomes even larger than in the dispersionless case: this is due to the fact that the approximation related to a short cutoff is not sufficiently accurate. In fact, extending the range of the e-ph correlations, the effective potential is progressively and substantially reduced for all momenta with respect to the dispersionless Holstein model. Numerical convergence for g 1D (q) is found at the value n * = 24 which corresponds to 48 lattice sites along the chain. Then, the 1D coupling renormalization is q-dependent and generally larger for q ∼ π/2 where the phonon dispersion generates stronger e-ph correlations. Exact diagonalization techniques [87] had found that multiphonon states with longer range hopping processes substantially soften the polaron band narrowing with respect to the predictions of standard SCPT with the largest deviations occuring at the polaron momentum ∼ π/2. In fact strong, nonlocal e-ph correlations should favour hopping on further than neighbors sites occupied by phonons belonging to the lattice deformation.
The projections of the two dimensional e-ph potential along the y component of the wave vector is shown in Fig. 3(b) for three values of the cutoff. For n * = 1 the correlation range is extended to the second neighbor shell thus including 8 lattice sites. For n * = 2 and n * = 3 the normalization is over 24 and 48 lattice sites respectively. Then, the three values of n * in 2D span as many lattice sites as the three values of n * in 1D respectively. I have made this choice to normalize consistently the potential for the linear chain and the square lattice. There is a strong renormalization for the 2D effective potential with respect to the dispersionless case for any value of the wave vector. Remarkably, the reduction is already evident for n * = 1. By extending the correlation range this tendency becomes more pronounced for q y close to the center and to the edge of the reduced Brillouin zone. An analogous behavior is found by projecting the 2D potential along the q x axis. The 2D potential stabilizes by including the 9th neighbor shell (n * = 3) in the correlation range.
Then, an increased range for the e-ph correlations leads to a very strong reduction of the effective coupling, an effect which is qualitatively analogous to that one would get by hardening the phonon spectrum. This is the physical reason underlying the fact that the Holstein polaron mass can be light.

Conclusions
According to the traditional notion of small polaron, the strong electron coupling to the lattice deformation implies a polaron collapse of the electron bandwidth (zt) and remarkable mass renormalization. Thus, bipolaronic theories have been thought of being inconsistent with superconductivity at high T c as the latter is inversely proportional to the bipolaron effective mass. Motivated by these issues, I have investigated the conditions for the existence of light polarons in the Holstein model and examined the concept of self-trapping versus dimensionality for a broad range of (anti)adiabatic regimes. The self trapping events correspond to the points of most rapide increase for the polaron effective mass versus the e-ph coupling. A modified version of the Lang-Firsov transformation accounts for the spreading of the polaron size due to retardation effects which are relevant mainly for moderate e-ph couplings. Assuming large optical phonon frequencies (ω 0 ), it is found that polaron masses display an abrupt although continuous increase driven by the eph coupling in any dimension, included the 1D case. Such self trapping events occur in the adiabatic and intermediate (zt/ω 0 ∼ 1) regimes provided that the e-ph coupling is very strong (g/ω 0 ≥ 2). The second order of strong coupling perturbation theory permits to decouple two fundamental properties of the polaron landscape, namely band narrowing and effective mass enhancement. I have considerd a range of e-ph couplings (g/ω 0 ≥ 1) in which the conditions for polaron formation are fulfilled. For a given value of the adiabatic ratio, the onset of the band narrowing anticipates, along the g axis, the crossover to the heavy polaron state which pins the quasiparticle essentially on one lattice site. In such buffer zone, between band narrowing onset and self trapping, the charge carri-ers can be appropriately defined as intermediate polarons.
In all dimensions, there is room for non trapped intermediate polarons in the region of the moderate couplings. Such polarons spreads over a few lattice sites and their real space extension grows with dimensionality. Thus, analytical methods can suitably describe some polaron properties also in that interesting intermediate window of e-ph couplings (i.e. 1 ≤ g/ω 0 ≤ 2 for the 1D system) in which adiabatic polarons are not large but not yet on site localized. In the square lattice, at the onset of the self trapping state, the polaron over bare electron mass ratios are ∼ 5 − 10 and even smaller for lower couplings.
The physical origins for the possibility of light polarons have been analysed treating the dispersive Holstein Hamiltonian model by the imaginary time path integral method. The perturbing source current is the time averaged e-ph Hamiltonian term. The partition function of the model has been derived integrating out the phonon degrees of freedom in semiclassical approximation. Focusing on the role of the intermolecular forces, I have quantitatively determined the effects of the dispersive spectrum on the time-averaged e-ph coupling which incorporates the nonlocality effects. It is found that, in the square lattice, the e-ph coupling is strongly renormalized downwards and this phenomenon is already pronounced by the (minimal) inclusion of the second neighbors intermolecular shell. This may explain why 2D Holstein polarons can be light while the real space range of the e-ph correlations is relatively short.