Combining experiments on luminescent centres in hexagonal boron nitride with the polaron model and ab initio methods towards the identification of their microscopic origin

The two-dimensional material hexagonal boron nitride (hBN) hosts luminescent centres with emission energies of ∼2 eV which exhibit pronounced phonon sidebands. We investigate the microscopic origin of these luminescent centres by combining ab initio calculations with non-perturbative open quantum system theory to study the emission and absorption properties of 26 defect transitions. Comparing the calculated line shapes with experiments we narrow down the microscopic origin to three carbon-based defects: C2CB, C2CN, and VNCB. The theoretical method developed enables us to calculate so-called photoluminescence excitation (PLE) maps, which show excellent agreement with our experiments. The latter resolves higher-order phonon transitions, thereby confirming both the vibronic structure of the optical transition and the phonon-assisted excitation mechanism with a phonon energy ∼170 meV. We believe that the presented experiments and polaron-based method accurately describe luminescent centres in hBN and will help to identify their microscopic origin.


Nomenclature of defect transitions
In the main text, we suppress the notations for the spin majority channel α as well as non-triplet spin states.For C 2 C B , this results in where ≡ stands for an equal representation.S2.Zero-phonon lines (ZPL), momentum displacements (∆Q) and Huang-Rhys factors (HR) of all studied defects.Here, α refers to the majority spin channel and β to the minority spin channel.The energy E in the second column gives the energy difference between the singlet and triplet states.It is only given for defects that show both a singlet and a triplet state.The PBE functional is used for all calculations except for the last three rows where the HSE06 functional was used.

III. FITTING AND COMPARISON WITH EXPERIMENT
In this section we outline the fitting procedure used to compare the polaron method to experiment.The fitting is done in two parts, in the first we extract the width of ZPL, the second fits the inhomogeneous broadening about the ZPL associated to acoustic phonons.To do both of these fits, we consider the emission lineshape driven by a excitation laser with detuning ∆ = 168 meV.We also assume that the phonons decay on a much faster timescale than light emission occurs, which allows us to ignore the excitation effects during the fitting procedure.This allows us to disregard the coherent driving term Ω R σ x /2 and the phonon dissipator K PH in Eq. 10, and assume that the emitter is initially in its excited state.

A. Fitting the zero-phonon line
To fit the ZPL, we restrict the experimental data to a frequency window ±0.008 eV around the measured ZPL.This allows us to ignore the phonon sideband contribution, fitting only the ZPL with only function S opt (ω).For an emitter initially in its excited state, the ZPL spectrum reduces to a simple Lorentzian lineshape: where ∆ω is the detuning from the emitter frequency.After first normalising the data and the above expression to the peak of the ZPL, we do a least min-squared fit.Note that we can fit only a single rate, since Γ and γ combine to give a single width of the Lorentzian, thus we set γ = 0 without loss of generality.We find an optimal Γ = 3.15 meV with residual of r = 0.098 for Emitter A.

B. Describing electron-phonon sidebands in photoluminescence
As discussed in the main tex, we divide the phonon spectral density into two components J Ph (ω) = J FP (ω) + J A (ω).The first contribution, J FP (ω) = k S k δ(ω − ω k ), captures the phonons calculated from first principles, where S k is the partial Huang-Rhys factor associated with a phonon with frequency ω k [2].To move to a continuum limit version, we account for the natural lifetime of phonons in a material by approximating the δ-functions in the definition of J Ph (ν) with a Gaussian function δ(ω) , where σ is the phonon broadening parameter.We choose the broadening parameter as σ = 5 meV, such that it phenomenologically reproduces features in the emission spectrum at low temperature.
The second contribution to the phonon spectral density, J A (ω), is associated to acoustic phonons, and leads to the inhomogeneous broadening around the ZPL, which is not captured by the first-principle calculations.We assume that this spectral density takes the form , where α and ω c are left as fitting parameters.This form of spectral density is commonly used in the semiconductor quantum dot literature, where it describes bulk acoustic phonons.Similarly to the ZPL fits we restrict ourselves to a frequency window of −0.06 eV < ∆ω < 0.01 eV.In contrast to the ZPL fits, however, we use the full expression for the emission spectrum, including the full phonon sideband given in Section XII C. The value of the optimal phonon parameters extracted from these fits changes depending on the defect.It was not possible to obtain a fit for the acoustic phonon sideband for V N C B due to the discrepancy between the sideband produced by first principles calculations and the measured PL.For illustrative reasons, we use the acoustic phonon parameters associated to C 2 C N when plotting the theoretical V N C B PL.
Fig. S3 shows the fitted PL lineshapes for the three defects considered in the manuscript.We compare the cases where acoustic phonons are neglected (dashed) and included (solid).For the C 2 C N and C 2 C B defect complexes, we see improved agreement about the ZPL when acoustic phonons are included.V N C B remains a are poor fit, with and without acoustic phonons.Notably, when acoustic phonons are included, we observe additional acoustic sidebands on the phonon replicas found from DFT calculations.

C. Theoretical photoluminescence emission spectroscopy
The model presented above, and the resulting fits, allow us to calculate the PLE directly from the full phonon model.We do this by varying the detuning of the continuous wave driving laser, ∆, and calculating the emission Comparison of the experiment (blue) and theory fits for three main defects considered in the main text.The dashed curves include only phonons calculated through ab initio calculations, while the solid orange curves show the fits with acoustic phonons.The spectra plotted in (b,d,f) are the same as (a,c,e) respectively, plotted over a larger range.
Temperature is set to T = 10 K.
spectrum for each detuning.For the driven case, we consider the steady-state emission spectrum, defined as: where the integral over t has been reduced to the steady state limit.

IV. DETAILED SCATTER PLOTS
The singlet system is shown on the left and the triplet system on the right.Emission from the triplet state is not possible by direct excitation because the ground state of the triplet state is not populated in equilibrium.However, optical excitation of the singlet system followed by an intersystem crossing (ISC) is possible.This requires laser detunings that are higher than 2.46 − 1.75 eV = 0.71 eV.This is significantly larger than our experimental laser detunings from ∼ 150 to 400 meV.Fig. S9a shows the photoluminescence spectra for several laser detunings.We obtain strong variations with laser detuning that are investigated in detail under continuous laser detuning, as shown in Fig. S9c.To investigate the absorption characteristic, we study the ZPL intensity which is defined as the area under the ZPL with a bandwidth of 10 meV (Fig. S9b).We observe a strong ZPL intensity at a detuning of 168 meV while an intermediate ZPL intensity is obtained at a detuning of 341 meV.These two important detunings are highlighted by horizontal lines in Fig. S9c and the corresponding photoluminescence spectra are shown in Fig. S9a.Another important detuning is exactly between the two aforementioned, i.e. at 255 meV, which yields a small ZPL intensity because it lays between the two states |e, 1⟩ and |e, 2⟩.

VIII. RESCALING OF PLE SPECTRA FOR EMITTER A
The PLE measurements shown in Fig. 2a are done by two laser sweeps: Sweep 1 and Sweep 2. We use two different laser settings to keep the power between 0.75 and 1.93 mW (see Fig. S11a).We first describe the rescaling within one sweep, followed by the rescaling between two sweeps.We finally describe the rescaling to one excitation power.
We measure the excitation power as a function of excitation energy prior to the PLE measurements (Fig. S11a).Within one sweep, the exposure time at power P i is set to where t 1 and P 1 are the exposure time and power at the beginning of the sweep while P i is the current power (at a specific laser energy or laser detuning).This adjustment in exposure time is meaningful, since the luminescent centres studied show a linear dependence of the intensity with excitation power in the studied power range, as shown in Fig. S11d.Furthermore, the detector operates in a constant-count regime which avoids saturation.
At each laser detuning in Sweep 2, the exposure time is adjusted as described by equation (4).To rescale the spectra of Sweep 2 to Sweep 1, we multiply each spectrum of Sweep 2 (after adjusting the exposure time) by where t * 1 and P * 1 are the exposure time and the excitation power at the beginning of Sweep 2. We will give an example why this multiplication is meaningful.At the end of Sweep 2, the excitation energy is identical to the beginning of Sweep 1.However, the excitation power at the end of Sweep 2 is P * end .Now we calculate the exposure time at the end of Sweep 2. Calculating the exposure time of Sweep 2, we multiply equation (4) with equation ( 5) and get which is identical to the exposure time given by equation ( 4) where P i is replaced by P * end .After carrying out all the rescaling mentioned above, we rescale the data to 1 mW at an excitation energy of 2.431 eV.This is the reason why we have shown absolute values in the photoluminescence (i.e.constant dose) in Fig. 2.

IX. PLE OF OPTICAL AND ACOUSTIC PSB OF EMITTER A
Fig. S12a shows how background-corrected PLE intensities are obtained, by using the ZPL intensity as an example.We use an unrelated photoluminescence range with the same width as the ZPL to define a background intensity (grey curve in the central panel).To obtain the background-corrected ZPL intensity, we subtract this background intensity.The resulting background-corrected ZPL intensity is shown in the right panel of Fig. S12a.Similarly, we obtain the background-corrected intensities for several PSB, as shown in Fig. S12b and in the main text.
Fig. S12b shows the optical PSB intensities.The right panel shows that the optical PSB intensities at one-phonon detuning (170 meV) and the ZPL intensity at two-phonon detuning (340 meV) show comparable strength.All these processes are two-phonon processes (see Methods of the main text).
Fig. S12c shows the low-energy PSB at detunings around 10 meV.The low-energy PSB intensity is lower than for the ZPL, since it is a two-phonon process (see Methods and discussion in the main text).For the backgroundcorrected intensities, we observe an almost identical qualitative behaviour of the low-energy PSB and ZPL intensity.This supports our interpretation that the low-energy PSB and the ZPL have the same excitation mechanism.S3.Zero-phonon line energies (in electron-volt) of Emitter M to Q.

XI. OPTICAL CHARACTERISATION
For optical characterisation, we use an objective (details below) to focus a laser onto the sample.For photoluminescence (PL), we use a continuous-wave (CW) laser with a photon energy of 2.37 eV.A helium flow cryostat is used for the low temperature measurements.
For PLE measurements at 10 K, we use a supercontinuum white light laser (78 MHz repetition rate) which is filtered down to a bandwidth of ∼ 1 nm with a tunable laser line filter.In order to compensate for intensity changes at different excitation wavelengths, we adjust the integration time of the spectrometer (see Supplementary Information VIII).For both photoluminescence (PL) and PLE measurements, the collected light is filtered by a longpass dichroic mirror (550 nm) and a 550 nm longpass filter.The filtered light is focused on the input slit of a spectrometer where the spectrum is obtained.The grating and detector information are given below.
For room temperature PL measurements, a 50X objective (NA = 0.6), a 150 lines/mm grating, and an EMCCD are used.The low-temperature PL and PLE experiments are performed using a home-made confocal microscope with a 50X objective (NA = 0.42), a 500 lines/mm grating, and a CCD detector.Furthermore, a spatial filter is used in front of the spectrometer for low-temperature measurements.

XII. MODELLING PHONON EFFECTS IN PHOTOLUMINESCENCE EMISSION FROM LUMINESCENCE CENTERS A. The model Hamiltonian
We model the studied defect complex as a two level system (TLS) with ground state |g⟩ and excited state |e⟩, with splitting ω e (ℏ = 1).The TLS is driven by a continuous wave laser with a frequency ω L and Rabi coupling Ω.The emitter couples to both a vibrational and optical environment, characterised by the Hamiltonian [3]: where we have defined a † l as the creation operator for a photon with energy ω l , and b † k as the creation operator for a phonon with energy ν k and wavevector k.We have also introduced the system operators σ x = σ † + σ and σ = |g⟩⟨e|.The coupling to the vibrational and electromagnetic environments are characterised by their respective spectral densities: for the optical environment we make the standard quantum optics approximation that the coupling constants h l are frequency independent, such that J EM (ω) ≈ Γ 2π , where Γ is the emission rate of the optical transition; the phonon spectral density takes the form are the partial Huang-Rhys parameters, and contain contributions from phonons calculated through ab initio methods and those fitted to the emission as detailed below.
We can simplify the above equation by making the rotating-wave approximation and moving to a frame rotating with respect to the laser frequency ω L , yielding the Hamiltonian: where ∆ = ω e − ω L is the detuning between the driving field and the exciton transition.This Hamiltonian forms the starting point of our analysis of the dynamical and optical properties of the defect.

B. Polaron theory for a driven emitter
In order to account for the strong coupling to the vibrational environment, we apply a polaron transformation to the global Hamiltonian, i.e. the unitary transformation U P = |e⟩⟨e| ⊗ B + + |g⟩⟨g|, where ) are displacement operators of the phonon environment [3,5,6].This transformation dresses the excitonic states with vibrational modes of the phonon environment.In the polaron frame, the Hamiltonian may be written as , where where we have introduced the phonon operators , with the Gibbs state of the phonon environment in the polaron frame given by ρ Notice that the polaron transformation has dressed the operators in the light-matter coupling Hamiltonian, and renormalised the system parameters, such that Ω R = ΩB and ∆ = ∆ − k ν −1 k g 2 k .The dynamics of the electronic states of the emitter are described by a 2 nd -order Born-Markov master equation [7], where both the electromagnetic field and the transformed vibrational environment are perturbatively eliminated.The polaron transformation means that the vibrational interaction is included non-perturbatively in the interaction strength [3].We can use the fact that only terms quadratic in field operators are non-zero when traced with a Gibbs state, so that there will be no cross-terms between the vibrational and electromagnetic dissipators, such that the master equation can be written in two parts: where ρ(t) is the reduced state of the emitter, L σ [ρ] = 2σρσ † − {σ † σ, ρ} is the Lindblad form dissipator for the electromagnetic environment, and similarly K PH describes the dissipation due to the vibrational environment.Note that we have also included a source of homogeneous broadening given by rate γ, to account for experimental imperfections and other dephasing mechanisms such as charge noise.For the optical dissipator, we refer the reader to standard texts in quantum optics for the derivation of the optical component [7], in the following section we outline the derivation of the phonon dissipator.

Deriving the phonon dissipator
Considering only the coupling to the vibrational modes, we move into the interaction picture with respect to the Hamiltonian H 0 , the interaction term becomes: where the interaction picture system operators are given by: and we have defined the generalised Rabi frequency η = ∆2 + Ω 2 R .The interaction picture bath operators are found as B x/y (t) = e iH0t B x/y e −iH0t , with B x/y defined above.In the Schrödinger picture and making the Born-Markov approximation in the polaron frame [7], the dissipator describing the electron-phonon interaction is given by: where we have introduced the rate operators: and defined the terms Γ a 0 = ∞ 0 Λ aa (τ )dτ , Γ a c = ∞ 0 Λ aa (τ ) cos(ητ )dτ , Γ a s = ∞ 0 Λ aa (τ ) sin(ητ )dτ , with a ∈ {x, y}.We may understand these terms as the rates at which transitions occur between the eigenstates of the system (i.e. the dressed states) induced by phonons.Qualitatively, Γ a s and Γ a c describe processes where a polaron is formed during excitation of the emitter, resulting in an exchange of energy with the environment.In contrast, Γ a 0 leads to no net exchange of energy, inducing a pure dephasing type process.
Since we are operating in the polaron frame, the above master equation is non-perturbative in the electron-phonon coupling strength, capturing strong-coupling and non-Markovian influences in the lab frame, despite having made a Born Markov approximation.This provides us with a simple, intuitive, and computationally straight-forward method for describing exciton dynamics in a regime where a standard weak coupling master equation would break down [3].
FIG. S1.Low-temperature photoluminescence of luminescent centres belonging to group I. (a) Photoluminescence line shapes at T = 10 K under 2.37 eV excitation, shifted vertically for clarity.The spectra from bottom to top correspond to Emitter A to Emitter L, given by the x axis in (c).(b) Optical PSB of the luminescent centres shown in (a) with corresponding colors, shifted vertically for clarity.The vertical dashed lines are at detunings of 160 and 190 meV, respectively.(c) Euclidean distance for Emitter A up to Emitter L, with corresponding spectra shown in (a) from bottom to top.(d) Histograms of Emitter A to Emitter L. The colors correspond to the defect transitions given by the legend.All line shapes are background corrected as outlined in Supplementary Information V and more details are given in the main text.

3 II.
DETAILS ON AB INITIO CALCULATIONS FIG. S2.Schematics, electronic levels and partial Huang-Rhys factors for C2C B , C2C N [β], and V N C B [S = 1].(a) Schematics and electronic levels of C2CB, C2CN[β], and VNCB[S = 1] where the latter is taken from Ref. [1].(b) The partial Huang-Rhys factors for C2CB, C2CN[β], and VNCB[S = 1] are shown from left to right.We highlight that the HSE06 functional was used for all three defects.
FIG. S4.Theoretical PLE maps fitted to Emitter A. The calculated PLE maps include electron-phonon interactions as described in the main text.

Fig. S5 shows
Fig. S5 shows the photoluminescence for C 2B V −1 N [S = 1] which does not conform nicely with the experiment.a

FIG. S5 .
FIG. S5.Comparison of the experimental photoluminescence (blue) with the theoretical spectrum of C2BV −1 N [S = 1] (red).The theoretical spectrum is calculated without acoustic phonons.The spectra in (a) and (b) are identical but plotted over different ranges.
FIG. S6.Detailed scatter plot for all studied defects.(a) Emitter A. The top left panel shows the experimental photoluminescence as well as the theoretical line shapes of C2CB[β] (red), C2CN (purple), and VNCN[S = 1] (dark green).The bottom panels show the scatter plot with a zoom-in on the right side, indicated by a square in the left plot.The top right panel shows the euclidean distance.The colors of the defects are given on the right.(b) Emitter B, labelling as in (a).(c) Emitter C, labelling as in (a).

Fig
Fig. S9 and S10 show the PLE data for Emitter A, B and C without any background correction, i.e. both photoluminescence and PLE intensities are shown without background correction.Fig.S9ashows the photoluminescence spectra for several laser detunings.We obtain strong variations with laser detuning that are investigated in detail under continuous laser detuning, as shown in Fig.S9c.To investigate the absorption characteristic, we study the ZPL intensity which is defined as the area under the ZPL with a bandwidth of 10 meV (Fig.S9b).We observe a strong ZPL intensity at a detuning of 168 meV while an intermediate ZPL intensity is obtained at a detuning of 341 meV.These two important detunings are highlighted by horizontal lines in Fig.S9cand the corresponding photoluminescence spectra are shown in Fig.S9a.Another important detuning is exactly between the two aforementioned, i.e. at 255 meV, which yields a small ZPL intensity because it lays between the two states |e, 1⟩ and |e, 2⟩.
FIG. S9.PLE of Emitter A at 10 K. (a) Photoluminescence spectra of the luminescent centre at several laser detunings, given by the legend.(b) ZPL intensity as a function of laser detuning.The spectral range for the ZPL intensity is shown by the vertical lines in (c).(c) PLE map of a luminescent centre in hBN.The horizontal lines correspond to the spectra shown in (b) and the dim upward line in the bottom right corner corresponds to the silicon Raman ∼ 520 cm −1 .
FIG. S11.Power dependence of excitation source and photoluminescence.(a) Photoluminesence spectra at several excitation powers, given by the legend.(b) The excitation power as a function of excitation energy.We performed two laser sweeps with different laser settings to obtain a similar power for all laser detunings.(c) Power-dependent photoluminescence.cps, counts per second.(d) Integrated zero-phonon line (ZPL) as a function of excitation power (not background corrected).The integrated spectral range (2.155 to 2.170 eV) is shown by vertical dashed lines in (c).We observe a linear dependence of the ZPL intensity as a function of excitation power.
FIG. S12.Background-correction of PLE intensities and acoustic PSB.(a) Background-correction of the ZPL intensity.The spectral ranges of ZPL and background are shown in the left panel, on top of the photoluminescence spectrum.The resulting PLE intensities are shown in the central panel and the background-corrected ZPL intensity in the right panel.The latter is defined as the green ZPL intensity minus the grey background intensity.The photoluminescence spectrum in (a) is taken at a detuning of 168 meV.(b) Optical PSB intensities, labelled as in (a).In the right panel, the flipped photoluminescence spectrum is shown in blue.(c) Low-energy PSB intensity, labelling as in (a).

TABLE S1 .
Zero-phonon line energies (in electron-volt) of Emitter A to L.
and the Frank-Condon factor of the phonon environment B = tr B