Blog

Low-energy effective Hamiltonian

We take in the present study the hole-doped Mott insulator on a square lattice, which is irradiated by an intense circularly-polarised light (CPL). This can be minimally modelled by the repulsive Hubbard model, with a Hamiltonian

where \({\hat{c}}_{j\sigma }\) annihilates an electron on the site j at R_j with spin σ = ↑,↓, and \({\hat{n}}_{i}={\sum }_{\sigma }{\hat{n}}_{i\sigma }={\sum }_{\sigma }{\hat{c}}_{i\sigma }^{{{{\dagger}}} }{\hat{c}}_{i\sigma }\) is the density operator. For the hopping amplitude t_ij, here we consider the second-neighbour (\({t}_{0}^{\prime}\)) as well as the nearest-neighbour hopping (t₀), as necessitated for realising a CPL-induced topological superconductivity [as we shall see below Eq. (17)]. The onsite repulsion is denoted as U (>0).

The circularly-polarised laser field is introduced via the Peierls phase, which involves R_ij ≡ R_i − R_j, and the vector potential,

for the right-circularly-polarised case; replace E with E* for the left circulation. We set ℏ = e = 1 hereafter.

Since the laser electric field is time-periodic, \(\hat{H}(t+2\pi /\omega )=\hat{H}(t)\), we can employ the Floquet formalism²¹, where the eigenvalue, called quasienergy, of the discrete time translation plays a role of energy. Because the quasienergy spectrum is invariant under the time-periodic unitary transformation, \({e}^{-i\hat{{{\Lambda }}}(t)}\) with \(\hat{{{\Lambda }}}(t+2\pi /\omega )=\hat{{{\Lambda }}}(t)\), we can introduce an effective static Hamiltonian \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\) as

where \(\hat{{{\Lambda }}}(t)\) is determined such that \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\) becomes time-independent. With such a transformation, we can analyse the time-dependent original problem using various approaches that are designed for static Hamiltonians.

While it is difficult to determine \(\hat{{{\Lambda }}}(t)\) and \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\) in an exact manner, we can obtain their perturbative expansions in various situations. The high-frequency expansion^{8, 22, 23} is a seminal example of such perturbative methods, where the effective Hamiltonian \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\) is obtained as an asymptotic series in 1/ω, and widely used for describing the Floquet topological phase transition. However, since the onsite interaction U in Mott insulators is typically much greater than the photon energy ω, the high-frequency expansion is expected to be invalidated for the present case.

Still, in the present system in a strongly-correlated regime t₀ ≪ U, we can employ the Floquet extension of the strong-coupling expansion (Schrieffer-Wolff transformation)^{11, 12, 19, 24, 25}, when the photon energy ω exceeds the band width ~ t₀ of the doped holes (and chosen to be off-resonant with the Mott gap ~ U), i.e., t₀ ≪ ω < U. As we describe the detailed derivation in Methods, we can obtain \(\hat{{{\Lambda }}}(t)\) and \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\) in perturbative series in the hopping amplitudes \({t}_{0},{t}_{0}^{\prime}\). While the t₀ ≪ U assumption yields a low-energy effective model known in the undriven case as a t–J model²⁶ which describes the dynamics of holes in the background of localised spins, an essential difference in the Floquet states under CPL is that extra terms emerge. Namely, the Hamiltonian of the irradiated case reads

where \({\hat{P}}_{G}={\prod }_{i}(1-{\hat{n}}_{i\uparrow }{\hat{n}}_{i\downarrow })\) projects out doubly-occupied configurations. Here, we have retained all the processes up to the second order (in the hopping amplitudes \({t}_{0},{t}_{0}^{\prime}\)), and additionally taken account of the fourth-order processes for the last line of the above equation. Thus the photon-modified exchange interaction, J_ij, and the photon-induced correlated hopping, Γ_i,j; k, are of second order, while the photon-generated chiral spin coupling, \({J}_{ijk}^{\chi }\), is of fourth order.

Let us look at the effective model in the laser field term by term. First, the hopping amplitude of the holes is renormalised from t_ij into \({\tilde{t}}_{ij}={t}_{ij}{{{{{{{{\mathcal{J}}}}}}}}}_{0}({A}_{ij})\) due to the time-averaged Peierls phase^{27, 28}, where A_ij = E∣R_ij∣/ω and \({{{{{{{{\mathcal{J}}}}}}}}}_{m}\) is the m-th Bessel function. Then we come to the interactions. The spin–spin interaction has a coupling strength dramatically affected by the intense electric field, since it is mediated by kinetic motion of electrons. Namely, the static kinetic-exchange interaction J is modulated as²⁴

which is a sum over m-photon processes, and appears in the second term on the first line of Eq. (5).

While these modulations of t and J also occur for the case of linearly-polarised lasers (i.e., for time-reversal symmetric modulations), a notable feature of the CPL is the emergence of the time-reversal breaking many-body interactions specific to the strong correlation. In the one-body Floquet physics, the two-step hopping (i.e., a perturbative process with hopping twice) in CPL becomes imaginary (and thus breaks the time-reversal symmetry) already on the noninteracting level if we take a honeycomb lattice^{2, 6, 8}, while such an imaginary hopping does not arise for the (noninteracting) square lattice due to a cancellation of contributions from different paths as depicted in Fig. 2a, which is why a honeycomb lattice fits with the FTI.

If we move on to correlated systems, on the other hand, the hopping process involves interactions for spinful electrons on the path, which lifts such cancellations, as exemplified in Fig. 2b. This implies that such a time-reversal breaking term is present in (interacting) square lattice as well. The two-step correlated hopping Γ appears in the second-order perturbation in the presence of holes [the second line of Eq. (5)], which involves three sites and sometimes referred to as the “three-site term”²⁶. While the three-site term also appears in the undriven case with the amplitude Γ_i,j; k = J_ij/4, this term involving the dynamics of holes (i.e., the change of holes’ positions) is usually neglected because its contribution is small in the low-doping regime [see the Gutzwiller factor in the next subsection, ~ g²δ in Eqs. (10), (11)]. In sharp contrast, when the system is driven by a CPL, Γ_i,j; k acquires an important three-site imaginary part as

with θ_ij defined as \({{{{{{{{\boldsymbol{R}}}}}}}}}_{ij}=| {{{{{{{{\boldsymbol{R}}}}}}}}}_{ij}| (\cos {\theta }_{ij},\sin {\theta }_{ij})\). The term triggers a topological superconductivity even when it is small, as we shall see.

In addition to the two-step correlated hopping Γ, another time-reversal breaking term arises if we go over to higher-order perturbations. In the previous studies of the half-filled case^{11, 12}, an emergent three-spin interaction (i.e., the scalar spin chirality term J^χ) is shown to appear in the fourth-order perturbation. This appears on the last line in Eq. (5) with a coefficient \({J}_{ijk}^{\chi }\). While this term is basically much smaller than the second-order terms, it may become comparable with the contribution from \({{{{{{{\rm{Im}}}}}}}}\,{{\Gamma }}\), as the spin chirality term does not accompany the dynamics of holes and has a larger Gutzwiller factor [~g³ in Eqs. (10), (11). See also Fig. 3].

Mean-field decomposition under Gutzwiller ansatz

Let us investigate the fate of the d-wave superconductivity in the presence of the time-reversal breaking terms. To this end, we here adopt the Gutzwiller ansatz^{26, 29} for the ground-state wavefunction and perform a mean-field analysis, as elaborated in Methods. After the approximate evaluation of the Gutzwiller projection, the problem of determining the ground state is turned into an energy minimisation of Eq. (5) without the Gutzwiller factor \({\hat{P}}_{G}\) if we renormalise the parameters. The mean-field solution is then formulated as a diagonalisation of the Bogoliubov-de Gennes Hamiltonian in the momentum space,

Here we have

where the expressions are independent of the site index i due to the spatial translational symmetry, and we have denoted the doping level as

with N being the number of lattice sites. Here τ runs over τ = mx + ny with \(m,n\in {\mathbb{Z}}\), and the label i + τ represents the site at R_i+τ = R_i + (m, n) in units of the lattice constant. The dependence on δ and a factor g = 2/(1 + δ) appears as a result of the Gutzwiller projection, which suppresses the contribution from charge dynamics as represented by \(\tilde{t}\) and Γ in the small δ regime. For the detailed derivation, see Methods. Throughout the present study, we take δ = 0.2. The mean-field Hamiltonian is self-consistently determined for the bond order parameter χ and the pairing amplitude Δ as

In the absence of the external field, the system undergoes a phase transition from normal to the \({d}_{{x}^{2}-{y}^{2}}\)-wave superconductivity when the temperature is sufficiently low²⁶. The corresponding order parameter is

for which the gap function F(k) becomes

as derived from the first term in Eq. (11) with τ = ± x, ± y. Here we have dropped a small correction due to δ. Without a loss of generality, we assume that Δ is real. The \({d}_{{x}^{2}-{y}^{2}}\) gap function has nodal lines along k_y = ± k_x, across which the gap function changes sign (see Fig. 4a–c below).

Now let us look into the possibility for the laser field converting this ground state into topological superconductivity from the time-reversal breaking terms [the second line in Eq. (11)]. The leading term should be those for \({\tau }^{\prime}=\pm x,\pm y\) (with \({{{\Delta }}}_{{\tau }^{\prime}}=\pm {{\Delta }}\)), which is nonzero even for the original \({d}_{{x}^{2}-{y}^{2}}\)-wave solution with no driving field and results in an imaginary gap function ∝ iΔ. In particular, gathering the terms with τ = ± (x + y), ± (x − y), we obtain the leading modulation to the gap function as

where we have defined

Thus we are indeed led to an emergence of a topological (chiral) \({d}_{{x}^{2}-{y}^{2}}+i{d}_{xy}\) superconductivity with \(F({{{{{{{\boldsymbol{k}}}}}}}})\simeq {F}^{{x}^{2}-{y}^{2}}+i{F}^{xy}\) with a full gap. For the full expression of the modulated gap function, see Eq. (65) in “Methods”. The key interactions, γ and J_χ, take nonzero values when the original Hubbard Hamiltonian Eq. (1) has the second-neighbour hopping \({t}_{0}^{\prime}\): Since the two-step correlated hopping Γ_i,j; k is composed of two hoppings, k → i and j → k, we can see that the imaginary coefficient above has \(\gamma \propto {t}_{0}{t}_{0}^{\prime}\). The chiral spin-coupling \({J}_{\chi }\propto {({t}_{0}{t}_{0}^{\prime})}^{2}\) also necessitates the next-nearest-neighbour hopping, which can be deduced from Eq. (41) in “Methods”.

Having revealed the essential couplings for the chiral superconductivity, let us now explore how we can optimise the field amplitude and frequency for realising larger topological gaps. For this, we set here \({t}_{0}^{\prime}=-0.2{t}_{0}\) and U = 12t₀ having cuprates with t₀ ≃ 0.4 eV in mind as a typical example.

We plot the essential γ and J_χ against the driving amplitude E and frequency ω in Fig. 3a, b. If we look at the overall picture, we can see that the dynamical time-reversal breaking is strongly enhanced along some characteristic frequencies in a resonant fashion at ω/U = 1/integer, which we can capture from the expressions for the coupling constants in the small-amplitude regime as follows.

The two-step correlated hopping γ is given, in the leading order in the amplitude, as

This expression, as a function of E, takes the maximal value around E ≃ 2.45ω/a, while diverges for ω → 0 or ω → U/2 for each value of E, as seen from the energy denominator. The enhancement occurs even in the low-frequency regime, which should be advantageous for experimental feasibility, since the required field amplitude (E ≃ 2.4ω/a) can be small.

If we turn to the spin-chirality term, J_χ has a complicated form involving Bessel functions [see Eq. (41) below]. If we Taylor-expand it in E, we have

which takes large values around ω = U due to the energy denominator in the above expression, while γ is small in this regime. The expression reveals that a dynamical time-reversal breaking also occurs as a fourth-order nonlinear effect with respect to the field strength E.

We can see in Fig. 3c–e how J_χ, as well as the exchange interaction J and the renormalised nearest-neighbour hopping \(\tilde{t}\), vary with the CPL intensity for several representative frequencies for which the dynamical time-reversal breaking becomes prominent. The hopping amplitude \(\tilde{t}\) vanishes around the peak of γ (E ≃ 2.45ω/a), as seen in Fig. 3c, d. This involves a zero of the Bessel function (at E ~ 2.40ω/a), and known as the dynamical localisation²⁷. In Fig. 3d, e, we can see a strong enhancement of J. Panel d has ω = 0.45U (slightly red-detuned from U/2), while panel e has ω = 0.88U (slightly red-detuned from U). If we go back to Eq. (6), the former has to do with the energy denominator for m = 2, while the latter for m = 1. While the driving frequency is set to be red-detuned from the resonance (which we can call “U − mω resonance”) in these results, the blue-detuned cases would give negative contributions from these terms, which will lead to a ferromagnetic exchange interaction, unfavouring the spin-singlet d-wave. With an optimal choice of the driving field, γ attains a significantly large value ≃ 0.04t₀, which yields \({F}^{xy}({{{{{{{\boldsymbol{k}}}}}}}})\simeq 0.3{t}_{0}\times i{{\Delta }}\sin {k}_{x}\sin {k}_{y}\) (for δ = 0.2). This is remarkably large and comparable with the undriven \({F}^{{x}^{2}-{y}^{2}}({{{{{{{\boldsymbol{k}}}}}}}})\simeq 0.7{t}_{0}\times {{\Delta }}(\cos {k}_{x}-\cos {k}_{y})\), even though the coupling constant itself is much smaller than J. This is the first key result of the present work.

Phase diagram

Now we investigate the ground state of the effective static Bogoliubov-de Gennes Hamiltonian for several choices of the CPL parameters. Let us show the gap function and the density of states in Fig. 4. In the absence of the external field in Fig. 4a–c, the gap function has the \({d}_{{x}^{2}-{y}^{2}}\) symmetry with nodal lines that give the zero gap at the Fermi energy. We now switch on the CPL, with ω = 10.5t₀ = 0.88U, E = 16t₀/a for Fig. 4d–f, or with ω = 1.7t₀ = 0.14U, E = 3.5t₀/a for Fig. 4g–i. The former represents the case where J_χ is dominant (see Fig. 3e), while in the latter γ plays the central role. In both cases we have the \({d}_{{x}^{2}-{y}^{2}}+i{d}_{xy}\) pairing, where the nodal lines are gapped out due to the complex gap function. We can indeed see in Fig. 4f, i clear energy gaps in the density of states with a gap size comparable with the original \({d}_{{x}^{2}-{y}^{2}}\) superconducting gap (as measured by the energy spacing between the two peaks in the density of states). Note that the band width and the \({d}_{{x}^{2}-{y}^{2}}\) superconducting gap are renormalised there, due to the modified \(\tilde{t}\) and J.

To examine the robustness of the d + id superconductivity, we calculate ∣Δ∣ and \({{{\Delta }}}^{\prime}\equiv -i{{{\Delta }}}_{\pm (x+y)}=i{{{\Delta }}}_{\pm (x-y)}\) at finite temperatures, which gives a phase diagram against T and field amplitude E in Fig. 5. We can see that the critical temperature T_c, as delineated by the region for Δ ≠ 0, increases significantly as the laser intensity E is increased. The enhanced T_c originates from the fact that J [and \({{{{{{{\rm{Re}}}}}}}}\,{{\Gamma }}\); see Eq. (38) below] are enhanced when the laser is applied, as we have seen in Fig. 3c–e. The time-reversal breaking order parameter \({{{\Delta }}}^{\prime}\) emerges below T_c down to T = 0, with the field dependence emerging from those of γ and J_χ.

Transient dynamics

So far, we have investigated static properties of the effective Hamiltonian \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\). While the results clearly show the presence of Floquet topological superconductivity in a wide parameter region, an important question from a dynamical viewpoint is whether we can achieve the Floquet topological superconductivity within a short enough time scale, over which the description with the effective Hamiltonian Eq. (5) is valid (i.e., over which we can neglect the coupling to phonons, the heating process due to higher-order perturbations, etc.).

Thus let us study the dynamics of the superconducting gap, by solving a quench problem formulated as follows. We first prepare the initial state as the mean-field ground state of the equilibrium t–J model [i.e., Eq. (5) with E = 0], and then look at the transient dynamics when the CPL electric field is suddenly switched on, by changing the Hamiltonian at t = 0 to the effective Hamiltonian \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\) (Eq. 5) with E ≠ 0. In such a treatment, \(\langle {\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\rangle\) is a conserved quantity for t > 0, and the driven state is expected to thermalise to the equilibrium state of \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\) having a temperature that corresponds to the internal energy \(\langle {\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\rangle\) (which we call an effective temperature). Note that here we implicitly neglect the fact that \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\)s for E = 0 and E ≠ 0 are on different frames [specified by \(\hat{{{\Lambda }}}(t)\)].

Before exploring the dynamics, it is important to check whether the expected steady state of the quench problem remains superconducting. We can do this by looking at the above-defined effective temperature on the phase diagram as indicated by dotted lines in Fig. 5. While the effective temperature rises as we increase the field strength E when we start from temperatures below T_c, we can see that the effective temperature does not exceed T_c in a wide parameter region, which implies that the Floquet topological superconductivity should indeed appear as a steady state of the quench dynamics. A sudden increase of temperature around E ~ 4t₀/a in Fig. 5c, d is due to the band flipping (sign change of \(\tilde{t}\)), which makes the kinetic energy of the initial state quite large. As we further increase the field strength, the effective temperature eventually reaches infinity and will even become negative¹³.

Now, the question is whether the equilibration (known as the Floquet prethermalisation) occurs fast enough (e.g., faster than the neglected heating processes, within the experimentally accessible pulse duration, etc.). We compute the time evolution of the superconducting order parameter in the time-dependent mean-field approximation (with the Gutzwiller ansatz), and plot the time evolution of ∣Δ∣ and \(| {{{\Delta }}}^{\prime}|\) in Fig. 6a, b, where we set the initial temperature at k_BT_initial = 0 with E = 7t₀/a, ω = 10.5t₀ = 0.88U. Here, the unit of time ℏ/t₀ corresponds to ~1 fs for t₀ ≃ 0.4 eV. We can see that ∣Δ∣ rapidly evolves with an overshooting behaviour, and converges to a certain value with a damped oscillation. Similar behaviour can be found in a previous study of the quench problem for the d-wave (but within the \({d}_{{x}^{2}-{y}^{2}}\) pairing) superconductor³⁰. The time-reversal breaking order parameter \(| {{{\Delta }}}^{\prime}|\) also quickly converges to a nonzero value. This is the second key result in the present work. We note that the obtained steady state slightly deviates from the equilibrium state with respect to \({\hat{H}}_{{{{{{{{\rm{F}}}}}}}}}\), where the deviation arises because the pair-breaking scattering³⁰ is lacking in the mean-field treatment of the dynamics.

The curious oscillation in the amplitude of the gap function can be interpreted as an excitation of the Higgs modes in superconductors^31,32,33, and thus the typical time scale for the oscillation (and the emergence of \({{{\Delta }}}^{\prime}\)) can be roughly estimated as the inverse of the superconducting gap (~1/2∣F(k)∣ with an appropriate k-average). This should be much faster than the time scale for heating, although its nonempirical evaluation would be difficult. We can analyse the present quench dynamics in a linearised form³⁴ if the change in the coupling constant is small, which reveals that the appearance of \({{{\Delta }}}^{\prime}\) is described as the A_2g Higgs mode in the d-wave sector with an amplitude proportional to E⁴.

If we have a closer look at the effective temperature in Fig. 5, we find an intriguing phenomenon: the superconducting state can appear even when we start from an initial temperature that is above the T_c before laser illumination. Namely, some dotted lines that start from T above T_c at E = 0 do plunge into the superconducting region as E is increased. We show the time evolution of the order parameter for this “nonequilibrium-induced superconductivity” in Fig. 6c, d, where we set k_BT_initial = 0.2t₀ as an initial temperature, and choose E = 16t₀/a, ω = 10.5t₀. Since the homogeneous mean-field ansatz without fluctuations cannot describe the spontaneous symmetry breaking, we instead inspect the growth of a tiny (homogeneous) perturbation Δ = 10⁻⁴ on the initial state. After the quench, both of the gap functions Δ and \({{\Delta }}^{\prime}\) grow exponentially and converge respectively to nonzero values, although the damped oscillation is slower than the previous case, and the converged values are far below those (∣Δ∣ ≃ 0.18 and \(| {{{\Delta }}}^{\prime}| \simeq 0.03\)) expected for equilibrium with the effective temperature.