Variational Monte Carlo approach for the Hubbard model applied to twisted bilayer ${\mathrm{WSe}}_2$ at half-filling

We consider an effective Hubbard model with spin- and direction-dependent complex hoppings $t$, applied to twisted homobilayer ${\mathrm{WSe}}_2$ using a variational Monte Carlo approach. The electronic correlations are taken into account by applying the Gutzwiller on-site correlator as well as long-range Jastrow correlators subjected to noninteracting part being of Pfaffian form. Our analysis shows the emergence of the Mott insulating state at the critical value of Hubbard interaction $U_{c1}\approx 6.5\left|t\right|\sim 7\left|t\right|$ estimated by extrapolating the density-density equal-time two-particle Green's functions. The signatures of an intermediate insulating phase between $U_{c1}$ and $U_{c2}\approx 9.5\left|t\right|\sim 10\left|t\right|$ are also discussed. Furthermore, we report the formation of the ${120}^{\circ }$ in-plane Néel state indicated by the detailed analysis of the spin-spin correlation functions. As shown, switching between antiferromagnetic phases characterized by opposite chirality could be experimentally realized by the change of perpendicular electric field. In a proper range of electric fields, also a transition to the in-plane ferromagnetic state appears.

Figure 1

Schematic representation of hopping phases according to spin and direction.

Figure 2

Structure of variational parameters used for the formulation of $|{\mathrm{\Psi }}_T\rangle$.

Figure 3

Ground-state energy $E_G$ obtained for the wide range of $U/\left|t\right|$ by means of the VMC approach (symbols and dotted line) and unrestricted Hartree-Fock (dashed black line). The two vertical dashed regions indicate values of $U$ for which an anomalous behavior of the energy curve appears. This becomes more evident in the plot of $\partial E_G/\partial U$ provided in the inset. The estimated statistical error of energy per lattice site is $\sim {10}^{-4}\mathrm{meV}$ thus is much smaller the symbol size.

Figure 4

Double occupancy $\langle {\widehat{d}}^2\rangle$ as a function of $U$. An abrupt decrease of double occupancies is clearly visible in vicinity of $U_{c1}$, whereas for $U_{c2}$, it is much less pronounced. As stated in main text, the affection of data by statistical error is negligble since it is of order of ${10}^{-5}$ for $\langle {\widehat{d}}^2\rangle$.

Figure 5

The Fermi surface (solid lines) for the spin-splitted bare band ($U=0$). The colors distinguish spin-up and spin-down subbands. The red dashed arrows indicate the trajectory, containing the high-symmetry points, along which the momentum resolved electron occupation number has been calculated and is provided in Fig. 6.

Figure 6

[(a)–(d)] Momentum resolved occupation number for both spin directions $\sigma$, along the path defined in Fig. 5 for different values of $U$. (e) Fourier transform of electron occupancy at the high-symmetry point $\mathbf{q}=\mathbf{M}$.

Figure 7

Values of $q^2/\mathcal{N}\left(\mathbf{q}\right)\equiv {\mathrm{\Lambda }}_U\left(\mathbf{q}\right)$ as a function of $\left|\mathbf{q}\right|\in \mathrm{\Gamma }\text{−}\mathbf{K}$ for the selected values of $U/\left|t\right|$. The dashed lines refer to the polynomial (fifth order) fits. The hollow symbols refer to the calculations performed for the $24\times 6$ supercell (see main text).

Figure 8

The estimate of ${lim}_{\mathbf{q}\to 0}{q}^2/\mathcal{N}\left(\mathbf{q}\right)$ as a function of interaction magnitude $U$. This quantity is proportional to the magnitude of the correlation-induced gap [38].

Figure 9

Real space spin-spin correlation functions $S_{A,j}$ for different values of the interaction amplitude $U$. (a) Sketch of the path along which the correlation functions have been collected; for example, in the A-B segment, index $j$ refers to the lattice sites moving along the $R_{BA}$ vector. [(b)–(d)] Spatial dependence of the spin-spin correlation functions for representative values of $U$.

Figure 10

The value of ${\zeta }_{i,j}$ along the path defined in Fig. 9. Spatial dependence of spin correlation function $S_{i,j}$ divided by the corresponding value of ${\zeta }_{i,j}$ collected for the sites included in the loop defined by vectors ${\mathbf{R}}_{BA}$, ${\mathbf{R}}_{CB}$, and ${\mathbf{R}}_{AC}$ (b).

Figure 11

The expectation value of total spin squared per site as a function of $U$. The insets present the zoom of $S^2$ in the vicinity of the critical values of $U$.

Figure 12

The maps presenting amplitudes of $\mathcal{S}\left(\mathbf{q}\right)$ for the four representative values of $U$ [(a)–(d)]. The increase in intensity of the red color indicates a higher value, and the color scale is the same for all the plots. Note, the we applied interpolation based smoothing for the presentation purposes. Peaks (red circles) are located exactly at ${\mathbf{q}}_{peak}\in \{\mathbf{K},{\mathbf{K}}^{\prime }\}$. In (d), we present $\mathcal{S}\left(\mathbf{q}\right)$ for the $U/\left|t\right|$ shown in (a)–(d) along path defined in Fig. 5, the peaks located at $\mathbf{K}\left({\mathbf{K}}^{\prime }\right)$ are clearly visible.

Figure 13

Momentum-resolved spin-spin correlation function at $\mathbf{q}=\mathbf{K}$, i.e., peak value of $\mathcal{S}\left(\mathbf{K}\right)$. The insets show values in the vicinity of $U_{c1}$ and $U_{c2}$.

Figure 14

The sketch of magnetic phase diagram provided in Ref. [10] which can be deduced from the analysis of the effective anisotropic Heisenberg model supplied with the Dzyaloshinskii-Moriya term. The plus(minus) sign in superscript of AF indicates clock(anticlock)wise rotation of spin in the given frame of reference. In the upper part of the diagram, we mark the range $\varphi \in (\pi /2,7\pi /6)$ which we have examined. Red arrows indicate the range of $\varphi$ that has been scanned experimentally in Ref. [5] by applying the displacement field in the range $(0,0.9)$ V/nm.

Figure 15

The peak values of ${\mathcal{S}}^{x-y}\left(\mathbf{q}\right)$ functions (see main text) for the selected range of phase $\varphi$. A solid vertical black lines indicate $\varphi =2\pi /3$ and $\varphi =\pi$.

Figure 16

Fourier transform of the $\mathbf{z}$-component of spin-spin correlation functions at $\mathbf{q}\in \{\mathrm{\Gamma },\mathbf{K}\}$ for the considered range of phase $\varphi$. Vertical solid black line indicates $\varphi =2\pi /3$ and $\varphi =\pi$.