Origin and stability of generalized Wigner crystallinity in triangular moiré systems

Abstract

Generalized Wigner crystals (GWC) on triangular moiré superlattices, formed from stacking two layers of transition metal chalcogenides, have been observed at multiple fractional fillings [Nature 587, 214 - 218 (2020), Nat. Phys. 17, 715 - 719 (2021), Nature 597, 650 - 654 (2021)]. Motivated by these experiments, tied with the need for accurate microscopic descriptions of these materials, we explore the origins of GWC at $n=1/3$ and $2/3$ filling. We demonstrate the general limitations of theoretical descriptions relying on finite-range, versus long-range interactions, however, we clarify why some properties are captured by an effective nearest-neighbor model. We study both classical and quantum effects at zero and finite temperatures, discussing the role of charge frustration, identifying a ``pinball" phase, a partially quantum melted GWC, with no classical analog. Our work addresses several experimental observations and makes predictions for how many of the theoretical findings can be potentially realized in future experiments.

Figures (11)

• Figure 1: Hamiltonian, classical phases, and energetics of the extended Hubbard model on the triangular lattice. (a) A schematic of the extended Hubbard model as in Eq. \ref{['eq:model_ham']}, with nearest neighbor hopping $t$, on-site Hubbard interaction $U$, and nearest, next-nearest and next-next-nearest neighbor interactions denoted by $V_1,V_2,V_3$ respectively. (b) Locations of lowest energy states in phase space, restricted to $n=1/3$ charge configurations shown in panel c, in the $V_1-V_2-V_3$ model ($V_1>0, U \to \infty$, $t=0$).. (c) Charge configurations for $n=1/3$ showing the triangular, dimer, trimer, fourmer crystals, from left to right, and expressions for the corresponding total energy for $N$ sites. (d, e) The individual energy differences of the dimer and trimer configurations, each with respect to the triangular configuration, as a function of interaction cutoff $r_c/a$ for a fixed $d/a=10$ (panel d), and gate distance $d/a$ for $r_c \rightarrow \infty$ (panel e). The calculation is carried out on a $N=72\times72$ system with periodic boundary conditions. The y-axis label for panels d,e is the same.
• Figure 2: Charge ordered quantum phases for a series of extended Hubbard models. (a-d) show the expectation value of the local charge density $\langle n_i\rangle$ on every site $i$ of a $6 \times 6$ torus for $n=1/3$, computed in the quantum ground state (obtained from DMRG) for a series of models. These models were obtained by truncating the LR model (model 3 in Table \ref{['table:moiré_parameters']}) to various ranges (a) model 4 but with $V_2=V_3=0$ (b) model 4 with $V_3=0$ (c) model 4 (d) model 3. The $y$ axis is vertical and the $x$ axis is horizontal. (e-h) show the corresponding absolute value of the Fourier transform of the charge distribution in momentum space; $n(\vec{k})$ is defined as $|\sum_i \left( \langle n_i \rangle - \frac{1}{3} \right) e^{i \mathbf{k} \cdot \mathbf{r_i}}|$ and $N_p$ is the number of particles (i.e. $N/3$). The DMRG simulations utilized a maximum bond dimension of up to 10000 and were carried out in the $S_z=0$ sector.
• Figure 3: Specific heat per site ($C_v/N$) as a function of temperature $(T)$, for various models considered in this work. (a) $C_v(T)/N$ for the LR model (model 2 in Table \ref{['table:moiré_parameters']}) and the NN model (model 2 but truncated to have only $V_1$ non-zero) computed with classical Metropolis Monte Carlo for multiple system sizes. The classical model with $t=0$, $U \rightarrow \infty$ is particle-hole symmetric, so the results for $n=1/3$ and $n=2/3$ are identical. (b,c) show $C_v (T)/N$ for the quantum models obtained from exact diagonalization for $N=18$ sites (see Supplementary Note 3) for three densities $n=1/3,2/3$ and $5/3$. (b) $C_v(T)/N$ for the LR quantum model (model 3) and (c) $C_v(T)/N$ for a NN (model 1). Additionally, in panels (b) and (c) we show the "classical" result for $n=1/3$ for the same system size, obtained by setting $t=0$ but not changing the other parameters (i.e. keeping $U$ large, but finite).
• Figure 4: Dependence of effective ground state and finite temperature properties of GWCs on gate distance. (a) Order parameter (defined in terms of the charge structure factor, $S( {\bf k}) \equiv \frac{1}{N} \sum_{i,j} \langle n_i n_j \rangle e^{i \mathbf{k} \cdot (\mathbf{r}_i - \mathbf{r}_j)}$ computed at a representative momentum point $k_{\sqrt{3} \times \sqrt{3}}$ that shows a peak for the $\sqrt{3}\times\sqrt{3}$ triangular charge order) in the quantum ground state, as a function of $d/a$. Calculations were carried out in the momentum (0,0) sector for $N=27$ and $n=1/3$ (in the $S_z=1/2$ sector) and $n=2/3$ (in the maximally polarized $S_z$ sector).The parameters correspond to model 3 (Table \ref{['table:moiré_parameters']}) for variable $d/a$. The dotted line represents the reference value for a perfect $\sqrt{3}\times\sqrt{3}$ triangular charge order. (b) Effective NN $V_{1, \text{eff}}/t$ interaction as a function of $d/a$ obtained from: the bare interaction potential ($V_{1, \text{Coulomb}} \equiv V(r=a)$, see Eq. \ref{['eq:coulomb_long_range']}), one defect energy ($V_{1,\Delta E}$), and from fitting the specific heat curve of the LR model (model 3 with variable $d$) and NN model ($V_{1,\text{C,fit}}$). (c) The critical temperature $T_c$ corresponds to charge order melting transition for three cases: classical, quantum $n=1/3$, and $n=2/3$, as a function of $d/a.$ (Inset) Ratio of the $V_{1, \text{eff}}$ to $T_c$ using $V_{1,\Delta E}$ (classical) and $V_{1,\text{C,fit}}$ (quantum) from (b).
• Figure 5: Schematic of energetic contributions arising from a single defect in the triangular GWC. The left panel shows the $\sqrt{3} \times \sqrt{3}$ triangular GWC and the underlying original moiré triangular lattice. The right panel corresponds to the configuration where a single charge in this GWC is moved to a vacant nearest neighbor site, creating a "defect" in the GWC. In both cases some (not all) of the non-zero contributions to the Coulomb energy are highlighted, the expression shown corresponds to the energy difference between the two configurations, which is thus the energy of a single defect $2 V_{1,\Delta E}$ in the effective nearest neighbor model. $V_n$ corresponds to the potential energy arising from the $n^{\rm{th}}$ nearest neighbors.