Charge-order melting in the one-dimensional Edwards model

We use infinite matrix-product-state techniques to study the time evolution of the charge-density-wave (CDW) order after a quench or a light pulse in a fundamental fermion-boson model. The motion of fermions in the model is linked to the creation of bosonic excitations, which counteracts the melting of the CDW order. For low-energy quenches corresponding to a change of the boson relaxation rate, we find behavior similar to that in an effective $t$-$V$ model. When the boson energy is quenched instead or a light pulse is applied to the system, the transient dynamics are more complex, with the CDW order first quickly decreasing to an intermediate value while the density-wave-like order of the bosons rises. In the case of pulse irradiation, the subsequent time-evolution of the CDW order depends strongly on the photon frequency. For frequencies slightly below the boson energy, we observe a temporary increase of the CDW order parameter. Our results reveal the complex physics of driven Mott insulators in low-dimensional systems with strong correlations.

Introduction.Ultrafast spectroscopy [1] and experiments on ultracold atoms in optical lattices [2] have allowed to directly observe the dynamics of quantum states in strongly-correlated materials.From a theoretical perspective, such non-equilibrium systems are intriguing, as they can give rise to exotic metastable states [3,4], and provide insight into the relaxation dynamics of quantum matter [5,6].For materials with long-ranged order, a natural question is how the order-parameter evolves after a perturbation, e.g., a quench or a light pulse [7,8].
In the context of one-dimensional (1D) systems, this has been addressed by numerical simulations for the magnetic order in spin chains [9,10] as well as for charge density waves (CDWs) of electrons [11,12].For the latter, it was found that an electron-phonon coupling modeled by a Holstein Hamiltonian not only renormalizes the electron mass via the formation of polarons, but can significantly alter the transient dynamics.This highlights that to faithfully model a system's dynamics, it is important to take into account possible couplings to environmental degrees of freedom.
The paradigmatic Edwards model describes particles whose movement requires the creation or annihilation of local bosons that parametrize the interaction with a background medium [13].This is a very generic situation in a great variety of condensed matter systems, in which the background can be seen as a deformable lattice (phonons), spins or orbitals forming, e.g., an ordered structure [14][15][16][17].Previous studies have focused on transport properties [18,19] and the ground-state phase diagram [13,20,21].At half band-filling, it was demonstrated in particular that there is a metal-insulator transition between a Tomonaga Luttinger liquid (TLL) and a CDW phase in one dimension.However, despite the rich physics captured by the Edwards model, its nonequilibrium properties are completely unexplored so far.
In this work, we utilize infinite matrix-product-state (MPS) methods to study the melting of the CDW order in the 1D Edwards model and discuss the differences compared with the previously investigated t-V and Holstein models.We consider both sudden quenches of the Hamiltonian parameters and pulse irradiation, which could be relevant for pump-probe experiments.
Model and method.The Hamiltonian of the 1D Edwards model is where f ( †) and b( †) are fermion and boson annihilation (creation) operators, respectively.Therefore Ĥ describes boson-affected quantum transport: A fermion emits or absorbs a local boson of energy ω 0 every time it hops between neighboring lattice sites.The bosons can relax via the λ-term.Obviously, large ω 0 and small λ parametrize a rather stiff background medium.A different representation of the model is obtained by the unitary transformation e i ( bi− b † i )λ/ω0 , which removes the λ-term but adds a free nearest-neighbor hopping −t f i ( f † i+1 fi + H.c.) with a renormalized transfer integral t f = 2λt b /ω 0 .
We assume a half-filled electron band, in which case the ground state of the Edwards model is either a TLL or a Mott insulator with CDW order [20].The phase diagram was derived in Ref. [21] by density-matrixrenormalization group calculations of the TLL parameter and the charge gap.For the benefit of the reader it is redrawn in Fig. 1.The CDW phase is characterized by a finite value of the electronic order parameter where N is the number of sites.Since the broken translation symmetry in the CDW phase is reflected in the boson densities as well, it makes sense to also consider the bosonic order parameter Importantly, O CDW +2O b is conserved in the limit λ = 0.
To simulate the time evolution after the system is driven out of equilibrium, we employ numerical methods based on infinite MPS.Namely, we first calculate the ground state with the variational method of Ref. [22], and then carry out the time evolution using the infinite time-evolving block decimation (iTEBD) [23] with a second-order Suzuki-Trotter decomposition and time step 0.025 t −1 b .The fermion and boson degrees of freedom are treated as separate sites in the MPS to keep the dimension of the local Hilbert spaces small.Since this increases the range of the couplings, we use the method of Ref. [24] and apply the iTEBD gates approximately using the time-dependent variational principle (TDVP) [25].
The boson Hilbert spaces need to be truncated to a finite dimension D b .Because the Mott insulator in the Edwards model typifies a CDW phase with low boson density, a relatively small value of D b is sufficient and methods to treat large local Hilbert spaces [26][27][28][29] are not necessary.We use D b = 6 for the simulations in this work, and take t b = 1 as the unit of energy.The truncation error in the two-site TDVP is kept below 5 • 10 −7 , which leads to a maxmimum MPS bond dimension of 3600.We have confirmed that the results are converged by doing additional simulations with lower bond dimensions.
Quenches.The simplest setting to study the melting of the CDW order is a quench of one of the Hamiltonian parameters.Here, we restrict ourselves to quenches from the CDW to the TLL phase [see Fig. 1], and assume that the system is initially in the ground state.and the boson densities after the boson-relaxation parameter is abruptly switched from λ = 0.1 to 0.2 at fixed boson energy ω 0 = 5.The order parameter decreases and exhibits damped oscillations around zero.Similar behavior has been previously observed for quenches in the spin-1/2 XXZ chain and the 1D Holstein model [9,12].As the boson energy ω 0 = 5 is quite large and not modified by the quench, the boson densities do not change significantly during the time evolution.The long-lived oscillations with angular frequency ≃ ω 0 already appear for isolated boson sites and can be attributed to the change of the boson eigenmodes.There is also a small staggered order O b of the bosons that, when appropriately scaled, closely follows that of the fermions.If λ is quenched to 0.4 instead [see Fig. boson densities also become more pronounced.Overall, the above observations fit into a perturbative picture, in which the dynamics of the electrons are described by an effective t-V model (see Supplemental Material [30]) and the stiff background contains almost no bosons.Because of the reduced energy scale of the effective model, the relaxation of the CDW order is slower than in a model without boson coupling and similar hopping amplitude, however.
As a second variant, we consider the quench from ω 0 = 2 to 1 at constant λ = 0.02.Here, the perturbative model should not be applicable, since ω 0 and t b are of the same order.Another difference is that the energy density relative to the ground state is much larger that in the previous case.As shown in Fig. 3(b), the CDW order O CDW again decreases after the quench.However, although the Hamiltonian parameters for the timeevolution correspond to a point in the TLL phase, this is accompanied by an increase of the boson order O b .In fact, O CDW + 2O b remains approximately constant until time t ≃ 4, and only slowly decreases afterwards.This shows that the initial dynamics can be mainly attributed to the boson-affected hopping ∝ t b , while the free hopping term ∝ t f in the alternate representation of the model becomes important for the long-time behavior.Both the fermion and boson orders appear to decay.Simulating the long-time evolution with MPS to determine the relaxation times is not feasible, however, because we can only reach times about t = 10 before the required bond dimension becomes prohibitively large.
We have attempted to push the simulation to slightly longer times by shifting a part of the time evolution to the operators that need to be evaluated [34,35].In this approach, the operators are expressed as matrixproduct operators (MPO) and evolved according to the Heisenberg picture until their bond dimension becomes too large.Unfortunately, we found that in the Edwards model the growth of the bond dimension is too fast to simulate significantly longer times.The quantity O CDW + 2O b can be evolved slightly longer, since it is conserved in the limit λ = 0. To take advantage of the small λ, we express the operator in momentum space by adapting the method of Ref. [36] to MPOs.The result of combining the MPS and MPO simulations is displayed in Fig. 3 as a dashed line.It is consistent with a continued slow decay of the density-wave order.
Pulse irradiation.A different way to melt the CDW order, which is motivated by pump-probe experiments, is to apply a light pulse to the system.The corresponding time-dependent Hamiltonian is obtained by a Peierls substitution f † j+1 fj ( b † j + bj+1 ) → e −iA(t) f † j+1 fj ( b † j + bj+1 ) in the hopping term.We consider a Gaussian pulse centered around t 0 with width σ p , amplitude A 0 and frequency ω p , i.e., the vector potential has the form In the following, σ p = 2 and A 0 = 0.2 if not stated differently.
Based on results for the optical response function [37], we expect that the pulse has the strongest effect when its frequency ω p is close to the boson energy ω 0 .For finite λ, there should also be a noticeable response for smaller ω p .If ω 0 is large, however, the physics in that case will be qualitatively similar to that in the t-V model, for which the CDW melting due to a pulse has been studied previously [11].We therefore focus on frequencies ω p ≈ ω 0 .Figure 4 shows the time evolution of the CDW order after a pulse is applied to a system with ω 0 = 5 and λ = 0.While the order parameter quickly decreases during the pulse, it stays nearly constant afterwards, except for small, persistent oscillations.These oscillations have frequency 2ω 0 and an amplitude scaling as A 2 0 for small pulse strengths, indicating that they can be attributed to two-boson excitations created by second-order transitions.That the CDW order does not melt completely is explained by the fact that for λ = 0, flipping the CDW order requires a large energy of approximately ω 0 per unit cell.Moreover, since O CDW + 2O b is conserved, the order parameter for either fermions or bosons must remain finite in the long-time limit.The effect of the pulse is thus to decrease O CDW and increase O b .It should be noted that the fermions can move by an even number of sites via effective higher-order hopping processes, so that the absence of CDW melting for λ = 0 does not apply to the CDW phase of the Edwards model at one-third filling [38].When ω p is detuned relative to ω 0 , the results are similar but the reduction of the CDW order and the strength of the oscillations become weaker.
More complex dynamics can be expected when there is a finite boson relaxation λ = 0.1 [see Fig. 5].Similarly as for λ = 0, O CDW (O b ) initially decreases (increases).
The CDW order now continues to change after the pulse, however, because of the effective free hopping with amplitude t f .Interestingly, in this transient regime O CDW increases for some pulse frequencies, and can even exceed the ground-state value before it falls off again.This timedependence of the CDW order, which is fundamentally different from the monotonic exponential decay found in Ref. [11] for pulse irradiation of the 1D Holstein model, does not require fine tuning of the pulse parameters and also occurs, e.g., for different values of the amplitude A 0 as demonstrated in Fig. 6.
A qualitative picture for the enhancement of the CDW order is provided by an effective Floquet model that is valid for large ω 0 and frequencies ω p sufficiently detuned from ω 0 [30].Like in the case without driving, the effective Hamiltonian is that of a modified t-V -model, but the nearest-neighbor repulsion is amplified for ω p < ω 0 and diminished for ω p > ω 0 .This agrees with the result in Fig. 5 that O CDW grows for low pulse frequencies.The influence of the increased repulsion can also be seen in the time evolution of O CDW + 2O b [inset of Fig. 5 (a)], which unlike O CDW does not exhibit a sharp drop-off due to the resonant excitation of bosons.
In constrast to the charge order of the fermions, both the boson density ⟨ b † i bi + b † i+1 bi+1 ⟩ and the order parameter O b remain nearly constant after the pulse.The approximate conservation of the boson density is expected because the boson energy ω 0 = 5 is large compared to the energy scale of the electronic excitations, i.e., ω 0 ≫ t f , t 2 b /ω 0 .To explain the slow evolution of O b it is helpful to look at the effective Hamiltonian [30], which shows that up to second order in the hopping bosons can not move between odd and even sites (at higher order, terms involving nearest-neighbor boson hopping that scale as λt 3 b /ω 3 0 appear).Accordingly, the time scale for the relaxation of the boson order O b should be much longer than that for the relaxation of O CDW .
The effective model furthermore suggests that it is energetically favourable for the fermion and boson orders to be aligned so that O CDW and O b have the same sign.Since O b seems to persists for long times, we may speculate that the fermion order O CDW will first settle around a reduced but finite value before it slowly approaches zero.Similar observations have been made for certain quenches in the spinless Holstein model, where the staggered displacement of phonons in the initial state can remain for long times [12].A difference in the Edwards model is that the boson order O b is rather small in the CDW ground states.Because of the boson-controlled hopping, it instead increases when the CDW order of the fermions is disturbed by a light pulse or a quench to lower boson energies ω 0 .
Conclusions.We numerically studied the chargeorder melting in terms of the 1D Edwards model, which is a minimal description for electrons interacting with a background medium.The quench parameters were chosen to drive a CDW-TLL insulator-metal transition as observed in many low-dimensional materials.In the Edwards model, this transition is of the Mott-Hubbard type rather than the Peierls type.By using the real-time simulation of infinite matrix-product states, we avoided boundary and finite-size effects in the simulations.Although the simulatable times are too short for the system to reach a true stationary state, the main effects of the fermion-boson coupling on the relaxation dynamics can be observed.For quenches that inject a large-enough energy into the system, there is a rapid growth of the boson density on a time scale t −1 b , along with a reduction of the fermion CDW order O CDW and an increase of the boson density-wave order O b .Applying a light pulse tuned to the boson energy ω 0 has a similar initial effect, but the dynamics that follow differ in the two cases.For the quench, we find a slow decrease of both O CDW and O b , while for the pulse, where the final Hamiltonian is still in the CDW regime with large ω 0 , only O CDW continues to change significantly whereas the boson order O b stays nearly constant.Assuming that thermalization occurs, both density-wave orders should eventually disappear.However, because of the different time scales for the dynamics of the fermion and boson orders, the system may first approach a quasi-stationary state with finite O CDW .A generalized Gibbs ensemble description based on the perturbative Hamiltonian as in Ref. [39] could perhaps give more insight into the properties of such a non-equilibrium state.
Remarkably, a pulse with frequency below the boson energy can lead to an increased CDW order at intermediate times.This is explained by an effective Floquet Hamiltonian that has an amplified nearest-neighbor interaction compared with the Hamiltonian from perturbation theory without driving.The enhancement of the CDW order also occurs for smaller ω 0 , however, where the effective model is no longer accurate.
For future studies, it would be interesting to investigate how the non-equilibrium charge dynamics differ in more elaborate models that take into account, e.g., a dispersion of the bosons.Understanding the effect of the fermion-boson coupling on the entanglement dynamics should also be valuable, in particular with regard to recent results on the quantum Fisher information in nonequilibrium systems [40].
The authors gratefully acknowledge the scientific support and HPC resources provided by the Er-langen National High Performance Computing Cen-

FIG. 1 .
FIG. 1. Ground-state phase diagram of the 1D Edwards model at half-filling according to Ref. [21].The red and blue arrows indicate the quenches studied in Figs.2(a) and 3 below, respectively.
Figure 2(a)   shows the time evolution of the order parameter O CDW

FIG. 2 .
FIG. 2. Time evolution of the CDW order parameter and boson densities for a quench from the CDW to the TLL phase at ω0 = 5.The initial state is the ground state for λ = 0.1, while the time evolution is according to the Hamiltonian with λ = 0.2 [panel (a)] or λ = 0.4 [panel (b)].beven ( bodd ) is a boson annihilation operator at an arbitrary even (odd) site.

O
FIG.3.CDW order parameter and boson densities for a quench from the CDW to the TLL phase at λ = 0.02.The boson energies before and after the quench are ω0 = 2 and ω0 = 1, respectively.

FIG. 4 .
FIG. 4. Time evolution of the CDW order parameter for pulses with different frequencies.Model parameters are ω0 = 5 and λ = 0.