Binding energies of trions and biexcitons in two-dimensional semiconductors from diffusion quantum Monte Carlo calculations

TL;DR

This paper provides highly accurate diffusion quantum Monte Carlo data on the binding energies of excitons, trions, and biexcitons in 2D semiconductors, aiding interpretation of experimental spectra and understanding of charge interactions.

Contribution

It offers the first statistically exact binding energy data for 2D semiconductor models, along with interpolation formulas for various material parameters.

Findings

Trion binding energies are generally larger than biexciton binding energies.

Contact pair densities enable perturbation theory descriptions of exchange interactions.

Interpolation formulas relate binding energies to effective masses and polarizability.

Abstract

Excitonic effects play a particularly important role in the optoelectronic behavior of two-dimensional (2D) semiconductors. To facilitate the interpretation of experimental photoabsorption and photoluminescence spectra we provide statistically exact diffusion quantum Monte Carlo binding-energy data for Mott-Wannier models of excitons, trions, and biexcitons in 2D semiconductors. We also provide contact pair densities to allow a description of contact (exchange) interactions between charge carriers using first-order perturbation theory. Our data indicate that the binding energy of a trion is generally larger than that of a biexciton in 2D semiconductors. We provide interpolation formulas giving the binding energy and contact density of 2D semiconductors as functions of the electron and hole effective masses and the in-plane polarizability.

Tables

Table 1. Table 1: Trion and biexciton binding energies for an important group of TMDCs. We compare our fitting formulas with path integral Monte Carlo (PIMC) Kylanpaa_2015 and DMC data Mayers_2015 in the literature. Effective masses (in units of the bare electron mass) are taken from G ​ W 𝐺 𝑊 GW calculations in the literature Cheiwchanchamnangij_2012 ; Shi_2013 .

Negative-trion binding energy ${ℰ}_T$ (meV)

Biexciton binding energy ${ℰ}_{\mathrm{XX}}$ (meV)

Material

$r_{\ast }$ (Å)

$m_e$

$m_h$

PIMC

DMC

Eq. (LABEL:eq:Trion_binding_fit)

Exp.

PIMC

DMC

Eq. (LABEL:eq:Biex_binding_fit)

Exp.

$\bullet$

MoS2

$39$ Cheiwchanchamnangij_2012

$0.35$

$0.43$

$32.0$

$33.8$

$32$

$34$ Trions2 , $35$ ZhangCleo

$22.7\left(3\right)$

$22.7\left(5\right)$

$24$

$60$ Sie_2015 , $70$ Biexcitons1

$■$

MoSe2

$40$ Kumar_2012

$0.38$

$0.44$

$27.7$

$28.4$

$31$

$30$ Ross2013 ; Singh2014

$19.3\left(5\right)$

$17.7\left(3\right)$

$23$

N/A

$▲$

WS2

$38$ Berkelbach_2013

$0.27$

$0.32$

$33.1$

$34.0$

$31$

$34$ Zhu2015 , $36$ Heinz1

$23.9\left(5\right)$

$23.3\left(3\right)$

$23$

$65$ Plechinger_2015

$▼$

WSe2

$45$ Berkelbach_2013

$0.29$

$0.34$

$28.5$

$29.5$

$27$

$30$ Trions3 ; Trions6

$23.9\left(5\right)$

$23.3\left(3\right)$

$20$

$52$ Biexcitons3

Equations

$$\left\{\sum _{k={\mathrm{e}}_1,{\mathrm{e}}_2}\left[-\frac{{\hslash }^2{\nabla }_{{\mathbf{r}}_k}^2}{2m_{\mathrm{e}}}-U\left({\mathbf{r}}_{k\mathrm{h}}\right)\right]-\frac{{\hslash }^2{\nabla }_{{\mathbf{r}}_{\mathrm{h}}}^2}{2m_{\mathrm{h}}}+U({\mathbf{r}}_{{\mathrm{e}}_1{\mathrm{e}}_2})\right\}\Psi$$

$$=\left[-\mathcal{U}-{\mathcal{E}}_{\mathrm{T}}\right]\Psi ,$$

$$U(\mathbf{r})$$

$$\mathcal{U}/R_{\mathrm{y}}^{\ast }$$

$$g_{\mathrm{eh}}^{\mathrm{X}}$$

$$\frac{{\mathcal{E}}_{\mathrm{T}}}{R_{\mathrm{y}}^{\ast }}$$

$$g_{\mathrm{eh}}^{\mathrm{T}}$$

$$g_{\mathrm{ee}}^{\mathrm{T}}$$

$$\frac{{\mathcal{E}}_{\mathrm{XX}}}{R_{\mathrm{y}}^{\ast }}$$

$$g_{\mathrm{eh}}^{\mathrm{XX}}$$

$$g_{\mathrm{ee}}^{\mathrm{XX}}$$

Full text

The optical properties of two-dimensional (2D) semiconductors such as monolayer MoS2, MoSe2, WS2, WSe2, InSe, and phosphorene have recently attracted a great deal of interest 2Doptical1 ; 2Doptical2 ; 2Doptical3 ; 2Doptical4 ; 2Doptical5 ; 2Doptical6 ; 2Doptical7 ; InSepaper . Numerous observations have been made of the rich structure of the luminescence spectra of these 2D materials, in which the most pronounced features have been interpreted in terms of neutral excitons Excitons1 ; Excitons2 ; Excitons3 ; Excitons4 ; Excitons5 ; Excitons6 , charged excitons (trions) Trions1 ; Trions2 ; Trions3 ; Trions4 ; Trions5 ; Trions6 , and biexcitons Biexcitons1 ; Biexcitons2 ; Biexcitons3 , while recent experiments on higher-quality monolayer transition-metal dichalcogenide (TMDC) samples have revealed additional structure in their spectra New1 ; New2 ; New3 ; New4 .

In this work we study a Mott-Wannier model of excitons and excitonic complexes in monolayer 2D semiconductors, taking into account the polarizability of the 2D crystal Heinz1 ; Heinz2 ; Bogdan and providing data to allow for a perturbative treatment of contact interactions between carriers. We use the diffusion quantum Monte Carlo (DMC) approach Ceperley1980 ; Foulkes2001 ; Needs2010 to find the energies of trions and biexcitons, and we provide approximate formulas for the exciton ($\cal U$), trion (${\cal E}_{\rm T}$), and biexciton (${\cal E}_{\rm XX}$) binding energies as functions of the in-plane polarizability and the electron and hole effective masses, which fit the DMC data to within 5%. We calculate and report contact pair densities, enabling the evaluation of perturbative corrections to the energies of charge-carrier complexes, as well as intervalley scattering, due to contact (exchange) interactions between charge carriers. The strength of the contact interactions could in principle be determined from first-principles calculations for different 2D semiconductors; alternatively, the strengths of the contact interactions can be regarded as parameters to be determined using experimental data in conjunction with the contact pair densities reported here.

The energy $-{\cal U}-{\cal E}_{\rm T}$ of a trion with one hole (h) and two electrons (e1 and e2) can be found by solving the Schrödinger equation (in Gaussian units)

[TABLE]

where $m_{\rm e}$ and $m_{\rm h}$ are the electron and hole effective masses and ${\bf r}_{ij}={\bf r}_{i}-{\bf r}_{j}$ is the position of particle $i$ relative to particle $j$. The Keldysh potential $U$ describes the Coulomb interaction screened by the polarization of the electron orbitals in the 2D lattice Heinz1 ; Heinz2 ; Bogdan ; Keldysh ,

[TABLE]

where $r_{*}=2\pi\kappa_{\bot}$ is a parameter directly related to the in-plane susceptibility $\kappa_{\bot}$ of the material, which has dimensions of length, and $\epsilon$ is the average permittivity of the media on either side of the 2D semiconductor. The potential is expressed in terms of a Struve function $H_{0}$ and a Bessel function of the second kind $Y_{0}$. Equation (1) determines the main contribution ${\cal E}_{\rm T}$ towards the trion binding energy, which is counted from the exciton binding energy $\cal U$. Similar Schrödinger equations can be written for an exciton and a biexciton.

Numerical solution of the Mott-Wannier Schrödinger equation for an exciton yields the $r_{*}$-dependent binding energy ${\cal U}$ (see the inset in Fig. 1), which agrees with the asymptotic limits Excitonlog1 ; Excitonlog2 ; Heinz1 ; Heinz2 ; Bogdan ${\cal U}(r_{*}\rightarrow\infty)\sim\frac{e^{2}}{\epsilon r_{*}}\ln\frac{r_{*}}{a_{\rm B}^{*}}$ and ${\cal U}(r_{*}\rightarrow 0)=-4\ R_{\rm y}^{*}$, as well as the contact pair density $g^{\rm X}_{\rm eh}=\left<\delta({\bf r}_{\rm eh})\right>$. Their $r_{*}$ dependence was fitted by

[TABLE]

where $x=r_{*}/(a_{\rm B}^{*}+r_{*})$, $a_{\rm B}^{*}=\epsilon\hbar^{2}/(\mu e^{2})$ is the excitonic Bohr radius, $\mu=(m_{\rm e}m_{\rm h})/(m_{\rm e}+m_{\rm h})$ is the reduced mass, and $R_{\rm y}^{*}=\mu e^{4}/(2\epsilon^{2}\hbar^{2})$ is the excitonic Rydberg. Here, the $\sqrt{x}$ term is not a physical singularity; rather, it reflects the enhanced contact density of a 2D hydrogen atom. $g^{\rm X}_{\rm eh}$ is plotted in Fig. 2.

The ground-state solution to Eq. (1) for a trion was obtained using the DMC approach Ceperley1980 ; Foulkes2001 , with the trial wave function being optimized using variational Monte Carlo (VMC). The trial wave function was of the Jastrow form $\Psi=\exp[J({\bf R})]$, where the Jastrow exponent $J({\bf R})$ consisted of a pairwise sum of terms of the form $u_{0}(r)=\left[c_{1}r^{2}\log(r)+c_{2}r^{2}+c_{3}r^{3}\right]/\left(1+c_{4}r^{2}\right)$, where $c_{1}$, $c_{2}$, $c_{3}$, and $c_{4}$ are optimizable parameters (different for each particle-pair type), together with two-body and three-body polynomial terms that are truncated at finite range Neil2004 ; Lopez2012 . The short-range behavior of $u_{0}(r)$ is such that the analogs of the Kato cusp conditions for the logarithmic interaction are satisfied (see the Supplemental Material in Ref. Bogdan ). Trial wave functions were optimized by unreweighted variance minimization Umrigar_1988 ; Drummond_2005 and energy minimization Umrigar_2007 . The ground-state wave functions for these systems are nodeless; hence the fixed-node DMC algorithm is exact. The DMC calculations were performed using the casino code Needs2010 with time steps in the ratio $1:4$ and the corresponding target configuration populations in the ratio $4:1$. Afterwards, the energies were extrapolated linearly to zero time step and hence, simultaneously, to infinite population. The resulting DMC trion binding energies, shown in Fig. 1, agree with the asymptotic binding energies found earlier in the limits of $r_{*}\rightarrow\infty$ Bogdan and $r_{*}\rightarrow 0$ Spink2015 . The trion contact pair densities $g^{\rm T}_{\rm ee}=\left<\delta({\bf r}_{{\rm e}_{1}{\rm e}_{2}})\right>$ and $g^{\rm T}_{\rm eh}=\left<\delta({\bf r}_{{\rm e}_{1}{\rm h}})+\delta({\bf r}_{{\rm e}_{2}{\rm h}})\right>$ were obtained by binning the radial distances sampled in the VMC and DMC calculations, evaluating the extrapolated estimate of the pair density Foulkes2001 , and then extrapolating the pair density to zero radial separation. The resulting contact pair densities are shown in Fig. 2. The trion binding-energy data are fitted (to an accuracy within 5%: see Fig. 1) by the formula

[TABLE]

while the contact pair densities are fitted by

[TABLE]

where $y=\mu/m_{\rm h}$. The term proportional to $\sqrt{1-y}$ in Eq. (LABEL:eq:Trion_binding_fit) describes the contribution to the ground-state energy due to the harmonic zero-point vibration of two heavy electrons treated using the Born-Oppenheimer approximation Spink2015 .

Similarly, the binding energies ${\cal E}_{\rm XX}$ of biexcitons were calculated using DMC and the results are presented in Fig. 3. ${\cal E}_{\rm XX}$ is the energy required to dissociate a biexciton into two separate excitons. A fitting formula with up to 5% accuracy,

[TABLE]

incorporates the fact that the biexciton binding energy is symmetric under the exchange of electrons and holes and includes the correct behavior in the Born-Oppenheimer/harmonic-approximation limit of extreme mass ratio. The biexciton electron-hole and electron-electron contact pair densities can be approximated as

[TABLE]

The ratio of the negative-trion to the biexciton binding energy is plotted against $x=r_{*}/(a_{\rm B}^{*}+r_{*})$ and $y=\mu/m_{\rm h}=m_{\rm e}/(m_{\rm e}+m_{\rm h})$ in Fig. 4. Although the biexciton binding energy is larger than the trion binding energy for the Coulomb interaction ($x=0$), the situation is generally reversed when the interaction is of logarithmic form ($x=1$). However, at extreme mass ratios, especially where the hole is heavy, the biexciton is stabilized with respect to the negative trion. In practice 2D materials typically have $x>0.9$ and $y\approx 0.5$ (see Table 1), and hence are strongly in the regime in which the trion binding energy exceeds the biexciton binding energy. The qualitative form of our predicted trion spectrum is shown in Fig. 5. The trion peak occurs at lower energy than the biexciton peak, in stark contradiction to the classification of experimental peaks reported in Ref. Plechinger_2015 . In fact several experimental works Biexcitons1 ; Biexcitons2 ; Biexcitons3 have reported biexciton binding energies of TMDCs that are about twice as large as the reported trion binding energies Trions1 ; Trions2 ; Trions3 ; Trions4 ; Trions5 ; Trions6 . However the physical origins of experimentally observed peaks in optical spectra are not always clear. Our conclusion that the trion binding energy is larger than the biexciton binding energy is robust against large changes in the values of the effective masses and the susceptibility and, taken at face value, suggests that the experimental “trion” and “biexciton” peaks may be misclassified.

In Table 1 we compare the trion and biexciton binding energies obtained using Eqs. (LABEL:eq:Trion_binding_fit) and (LABEL:eq:Biex_binding_fit) with previous theoretical calculations and experimental results in the literature for molybdenum and tungsten dichalcogenide materials. The theoretical results are in good agreement with each other, and also with experimental results for the trion. However, for the biexciton, there is a major disagreement between theory and experiment: the experimental binding energies are around three times larger than the theoretical biexciton binding energies. Since our DMC solution of the Mott-Wannier model is exact, the quantitative disagreement between the positions of the theoretical and experimental trion and biexciton peaks must indicate either a serious inaccuracy in the Keldysh interaction between charge carriers or a misinterpretation of experimental spectra. One possibility is that contact (exchange) interactions between charge carriers as well as intervalley scattering effects could play a significant role in charge-carrier complexes. Using our contact pair density data together with ab initio calculations of contact interaction parameters could provide a promising avenue for improving the quantitative description of the measured photoemission spectra.

In summary we present exact numerical data for the ground-state solutions of Mott-Wannier models of trions and biexcitons in 2D semiconductors in which the charge carriers interact via the Keldysh interaction. We have evaluated the contact pair density between charge carriers, to permit subsequent perturbative evaluations of the energy contribution due to contact exchange interactions. Our results suggest that experimental spectra have been misclassified, because the trion binding energy should exceed the biexciton binding energy, but they also indicate that the Keldysh interaction fails to give a quantitative description of the observed excitonic properties of 2D TMDCs. The contact pair density data that we provide will enable the theoretical and experimental exploration of the role played by contact interactions between charge carriers and intervalley scattering in 2D semiconductors.