Efficient calculation of trion energies in monolayer transition metal dichalcogenides

The reduced dielectric screening in atomically thin semiconductors leads to remarkably strong electron interactions. As a result, bound electron-hole pairs (excitons) and charged excitons (trions), which have binding energies in the hundreds and tens of meV, respectively, typically dominate the optical properties of these materials. However, the long-range nature of the interactions between charges represents a significant challenge to the exact calculation of binding energies of complexes larger than the exciton. Here, we demonstrate that the trion binding energy can be efficiently calculated directly from the three-body Schrödinger equation in momentum space. Key to this result is a highly accurate way of treating the pole of the electronic interactions at small momentum exchange (i.e., large separation between charges) via the Landé subtraction method. Our results are in excellent agreement with quantum Monte Carlo calculations, while yielding a substantially larger ratio of the trion to exciton binding energies than obtained in recent variational calculations. Our numerical approach may be extended to a host of different few-body problems in 2D semiconductors, and even potentially to the description of exciton polarons.

💡 Research Summary

The paper addresses the longstanding challenge of accurately computing trion binding energies in atomically thin two‑dimensional (2D) semiconductors, especially transition‑metal dichalcogenide (TMD) monolayers where reduced dielectric screening yields exceptionally strong Coulomb interactions. While exciton (electron‑hole pair) binding energies in these materials are on the order of hundreds of meV, trions (charged excitons) exhibit binding energies of tens of meV, making them dominant contributors to optical spectra at room temperature. However, the long‑range 1/r Coulomb potential leads to a divergence (pole) at zero momentum transfer in momentum‑space formulations, rendering direct three‑body Schrödinger‑equation solutions numerically unstable.

The authors adopt a technique originally developed in nuclear physics—Landé subtraction—to regularize the singularity. They introduce a compensating kernel (g(p,k)=\frac{2p^{2}}{p^{2}+k^{2}}V(p-k)) that reproduces the pole of the interaction potential (V(q)). By subtracting and adding this term inside the integral equation, the pole is analytically cancelled, leaving a well‑behaved effective kernel. Crucially, the remaining term (K(p)=\sum_{k}g(p,k)) can be evaluated analytically for both the bare 2D Coulomb potential and the screened Rytova‑Keldysh potential, which incorporates a material‑specific screening length (r_{0}). This analytic treatment eliminates the most problematic part of the integral, dramatically accelerating convergence.

First, the method is demonstrated on the exciton problem. Using Gaussian quadrature with a modest grid (as few as 16 radial and 8 angular points), the authors recover the full 2D Rydberg series and its modification by finite (r_{0}). The exciton energies agree with previous high‑precision calculations, confirming that the subtraction scheme yields accurate results with minimal computational effort.

The core contribution is extending the same subtraction to the three‑body trion Schrödinger equation written in momentum space. The trion state is expressed as (|X^{-}\rangle = \sum_{k,p}\chi_{k,p}, e^{\dagger}{\uparrow,k} e^{\dagger}{\downarrow,p} h^{\dagger}_{-k-p}|0\rangle), leading to an integral equation containing two electron‑hole attraction terms and one electron‑electron repulsion term. After rearranging, each interaction kernel is regularized with its own (g) function, and the corresponding (K) terms are added back analytically. The resulting equation resembles the exciton case but includes off‑grid momentum combinations ((k+p-k’) and (k+p-p’)). To handle these efficiently, the authors employ a Lanczos‑type iterative algorithm that builds a Krylov subspace directly from the integral operator, allowing rapid convergence to the lowest eigenvalue (the trion ground‑state energy).

Numerical results are presented for the symmetric case of equal electron and hole masses and for a range of screening lengths. In the unscreened limit ((r_{0}=0)), the trion binding energy is found to be (\varepsilon_{T}=0.122,\varepsilon_{X}), where (\varepsilon_{X}) is the exciton binding energy. This ratio is noticeably larger than that obtained from recent variational wave‑function approaches ((\sim0.11)) but matches diffusion quantum Monte Carlo (QMC) calculations ((\sim0.12)) within statistical error. When the Rytova‑Keldysh screening is introduced, the trion‑to‑exciton binding‑energy ratio remains in excellent agreement with QMC across the entire experimentally relevant range of (r_{0}/a_{0}) (typically 10–100 for monolayer TMDs). The authors also compare with exact diagonalization data from the literature and find consistent agreement, confirming the robustness of their method.

Beyond the immediate results, the paper emphasizes the generality of the approach. The Landé subtraction can be applied to any long‑range interaction that exhibits a simple pole, including modified potentials arising from heterogeneous dielectric environments. Combined with the Lanczos iterative solver, the technique scales favorably with the number of particles, suggesting straightforward extensions to four‑body complexes (biexcitons), exciton‑polaron problems, and even many‑body Green’s‑function calculations in 2D materials. The authors argue that this framework bridges the gap between computationally cheap but approximate variational methods and the highly accurate but expensive QMC, offering a practical tool for theorists and experimentalists seeking quantitative predictions of few‑body bound states in emerging 2D semiconductors.