2 h-e-kB-a-1 in which e is the electronic charge and a is the latti W)and r solution of the non-interacting model in Eq.(1)has been shown [14.15]to reprodu .The the conventional Boltzmann elerc-d. with the on-site Coulomb interaction parameter as FIG =∑(-）)(du-) tric field E. (2) =4.The lin e in the syn ric limi (DMFT JIG and int ting (mode ty deviates from t 28])to treat the many-body interaction via a self istent local approximation of the self nergies.Note the self-emergy ha oth the idth W the d case.All following energis are in unit of une (w)andS(）=2ifE w)+ (w)with the Fermi-Dirac(FD)distribution fer exp(w/T川 the local retarde means and lo ions(GFs). p a homoe eous non in the Coulomb y(t) sites are quivalent. In the Coulomb gauge, this lead U2(t)2(t).The GFs are updated with this self Gu( the energy using the above Dyson's equations and the proce dure is repeated until convergence is achieved formation from the temporal gauge Wiwe gen to hig her 1 Below,we present the implementation of our DMFT scheme the Coul ts in gauge directly g ou pendicular direction and the above construction of the Dyson's equation can be carried out independen tly pe semi-infinite dissipative Hul ailed divec e S bard chains and its own reservoir)with a sel results of the model in one and three dimension E The ng We first discuss the linear response regime. Withir treated hy means of an imnurity solve the DMFT,the DC in the limit of zere For given self-energy w)r(w+E),the on site Green's functions obey the following Dyson equations the Ku G(w)-1_ w-)-2Fo() )LfED() the spectral function at a given wave-vector k p() k+订 ong a in which 2E< Thi imi to the left and lar to the one used by prange F(w +E)+F and for the (w)is the on-sit I th ain（=1whic electron-phonon interaction.Recent calculat nilar dy e FI(w)-1=w-Eiot(w)-Y2FI(w+E), (5 on is in dependent of the interaction strength i 261 both in which can b solved recursively after (a)one and (b)three-dimension The linear behav and can be 0.00 i (a). or of the impurity GFs,,are constructed using ing pa the W+-W05 mith th g)-1=w+i-2F) (6 factor=l-Re话/a g(w)Ig"(w)22irfeD(w)+2F(w)].(7) With increasing E-field,the contribution atE=U/2is2 ~ = e = kB = a = 1 in which e is the electronic charge and a is the lattice constant. In the rest of this Letter, we measure energies in units of the full TB bandwidth W = 4γ = 1 (1-d) and W = 12γ = 1 (3-d). The exact solution of the non-interacting model in Eq. (1) has been shown [14, 15] to reproduce the conventional Boltzmann transport theory despite the lack of momentum transfer scattering. The Hubbard model Hˆ = Hˆ 0 + Hˆ 1 is defined with the on-site Coulomb interaction parameter U as Hˆ 1 = U X `  d † `↑ d`↑ − 1 2  d † `↓ d`↓ − 1 2  . (2) Our calculations are in the particle-hole symmetric limit. We use the dynamical mean-field theory (DMFT [16, 28]) to treat the many-body interaction via a self￾consistent local approximation of the self-energies. Note that the self-energy has contributions from both the many-body interaction Hˆ 1 and the coupling to the reser￾voirs: Σr tot(ω) = −iΓ+Σr U (ω) and Σ< tot(ω) = 2iΓfFD(ω)+ Σ < U (ω) with the Fermi-Dirac (FD) distribution fFD(ω) ≡ [1 + exp(ω/Tb)]−1 . Once the local retarded and lesser self-energies are computed, one can access the full re￾tarded and lesser Green’s functions (GFs). Note that in a homogeneous non-equilibrium steady state, all the TB sites are equivalent. In the Coulomb gauge, this leads to G r,< ``0 (ω) = G r,< `+k,`0+k (ω + kE) and similarly for the self-energies [15, 25], as can be derived via a gauge trans￾formation from the temporal gauge. Below, we present the implementation of our DMFT scheme in the Coulomb gauge directly in the steady states. It consists in singling out one TB site – say ` = 0 – (often referred as impurity) and replacing its direct environment (i.e. semi-infinite dissipative Hub￾bard chains and its own reservoir) with a self-consistently determined non-interacting environment (often referred as Weiss “fields”). The local electronic problem is then treated by means of an impurity solver. For given self-energy [Σr,< ` (ω) ≡ Σ r,< U (ω+`E)], the on￾site Green’s functions obey the following Dyson equations G r (ω) −1 = ω − Σ r tot(ω) − γ 2F r tot(ω), (3) G <(ω) = |G r (ω)| 2 [Σ< tot(ω) + γ 2F < tot(ω)], (4) in which γ 2F r,< tot are the total hybridization functions to the left and right semi-infinite chains, F r,< tot (ω) = F r,< + (ω + E) + F r,< − (ω − E). F+(ω) is the on-site re￾tarded GF at the end of the RHS-chain (` = 1) which obeys the self-similar Dyson equation F r +(ω) −1 = ω − Σ r tot(ω) − γ 2F r +(ω + E), (5) which can be solved recursively after more than 500 iter￾ations. F−(ω) corresponds to the GF of the LHS-chain and can be obtained similarly. The non-interacting parts of the impurity GFs, G, are constructed using G r (ω) −1 = ω + iΓ − γ 2F r tot(ω) (6) G <(ω) = |Gr (ω)| 2 [2iΓfFD(ω) + γ 2F < tot(ω)]. (7) 0 0.5 1 1.5 electric-field, E/W 0 0.02 0.04 0.06 0.08 0.1 current, J U = 0 U = 1.5 1-d tight-binding chain, K/W = 0.0625, Tbath/W = 0.00125 0 0.005 0.01 E/W 0 0.005 0.01 J linear-response U/2 - peak U - peak (a) 0 0.1 0.2 0.3 0.4 0.5 0.6 electric-field, E/W 0 0.005 0.01 0.015 current, J U = 0 U = 0.33 U = 0.50 U = 0.67 U = 0.83 3-d tight-binding lattice, Γ/W = 0.0083, Tbath = 0.0 0 0.008 0.016 E/W 0 0.002 0.004 0.006 J (b) FIG. 1: (color online) Electric current (per spin) J vs. elec￾tric field E. (a) 1-d chain with damping Γ = 0.0625W and fermion bath temperature Tb = 0.00125W with the 1-d TB bandwidth W = 4γ. The linear conductance in the small field limit (magnified in the inset) is the same for non-interacting (U = 0) and interacting (U = 1.5W) models. After the con￾ductivity deviates from the linear response behavior, inelastic contributions appear at E = U/2 and E = U. (b) 3-d lat￾tice with Γ = 0.0083W and Tb = 0.00042W with the 3-d TB bandwidth W = 12γ. The main features remain similar to the 1-d case. All following energies are in unit of W, unless otherwise mentioned. The local self-energies are obtained by means of the iterative-perturbation theory (IPT) up to the second-order in the Coulomb parameter U: Σ≷ U (t) = U 2 [G ≷(t)]2G ≶(t). The GFs are updated with this self￾energy using the above Dyson’s equations and the proce￾dure is repeated until convergence is achieved. We generalize the above method to higher dimensions. With the electric-field along the principal axis direction, E = Exˆ, the lattice is translation invariant in the per￾pendicular direction and the above construction of the Dyson’s equation can be carried out independently per each perpendicular momentum vector. See Supplemen￾tary Material for a detailed discussion. Below, we present results of the model in one and three dimensions. We first discuss the linear response regime. Within the DMFT, the DC conductivity in the limit of zero temperature and zero electric field can be obtained via the Kubo formula as σDC ∝ limω→0 P k R dνρk(ν)ρk(ν + ω)[fFD(ν) − fFD(ν + ω)]/ω = P k R dν[ρk(ν)]2 δ(ν) with the spectral function at a given wave-vector k ρk(ν) = −π −1 Im[ν − k + iΓ − Σ r U (ν)]−1 . Therefore, as long as Σ r U (ν) → 0 as ν → 0, T → 0, the DC conductivity is independent of the interaction. This argument is simi￾lar to the one used by Prange and Kadanoff [29] for the electron-phonon interaction. Recent calculations did not have access to the linear response regime [21, 23, 24]. FIG. 1 confirms the validity of the linear response analysis. The initial slope of the J − E relation is in￾dependent of the interaction strength U [26] both in (a) one and (b) three-dimension. The linear behav￾ior deviates at the field Elin ≈ 0.003 in (a), orders of magnitude smaller than the renormalized bandwidth W∗ = zW ≈ 0.5 with the equilibrium renormalization factor z = [1 − Re∂Σ r U (ω)/∂ω] −1 ω=E=Tb=0. With increasing E-field, the contribution at E = U/2 is