Vectorial nonlinear coherent response of a strongly confined exciton–biexciton system

The vectorial four-wave mixing response of an individual strongly confined exciton–biexciton system with fine-structure splitting in a GaAs/AlGaAs quantum dot is measured by dual-polarization heterodyne spectral interferometry. The results are compared with theoretical predictions based on the optical Bloch equations. The system is described by a four-level scheme, which is a model system of the nonlinear excitonic response in low-dimensional semiconductors. We measure its coherence properties and determine the underlying dephasing mechanisms. An impact of the inhomogeneous broadening by spectral wandering on the coherent response is investigated. We further discuss the different four-wave mixing pathways, polarization selection rules, the time-resolved polarization state, the vectorial response in two-dimensional four-wave mixing and ensemble properties.


Introduction
An individual exciton in a strongly confined quantum dot (QD) is a model system for quantum optics, demonstrating single-photon emission [1], two-photon indistinguishability [2], coherent control [3,4], entanglement [5], optical dressing [6][7][8] and slow light [9]. Excitons are also candidates for ultrafast information processing schemes in solids [10], requiring accurate measurement, manipulation and modeling of their coherent response. The coherent response of excitons can be optically inferred via nonlinear techniques [11][12][13]. In particular, fourwave mixing (FWM) and spectral hole burning spectroscopy are efficient approaches to extract coherence of excitonic ensembles in bulk semiconductors and quantum well structures. Singlechannel heterodyne detection of the FWM enabled probing coherence of excitons in QD ensembles [14]. Combining multichannel optical heterodyning with spectral interferometry (HSI) enabled the observation of FWM [15] and two-dimensional FWM (2DFWM) [16] from individual excitons. Up to now, these measurements were restricted to either cavity-enhanced responses [17][18][19] or individual excitons with a high optical dipole moment µ [15,16,[20][21][22][23], since the FWM amplitude scales with the fourth power of µ. Such localized excitons gain their large dipole moment by a weak confinement, providing a large coherence volume. However, the weak confinement energy of the order of 10 meV is comparable with their Coulomb-correlation energies, such that the related localized biexcitons are quite different from a product state of the excitons. Furthermore, the nearby continuum states allow for biexciton-mediated inter-exciton coupling [22]. Their physics is therefore more complex than that of an exciton-biexciton fourlevel system (4LS).
In this work, we measure the coherent nonlinear response of an exciton-biexciton system in an individual strongly confined QD. Owing to their well-defined level structure, excitons in this model system are better suited for quantum information processing applications than weakly confined excitons. Conversely, measuring on individual strongly confined QDs is more challenging due to their smaller transition dipole moment µ. The excitation intensity needed to reach a given third-order coherence scales as µ 6 , and constitutes a background for the measurement of the nonlinear emission. Moreover, their dephasing is typically not radiatively limited, as they exhibit inhomogeneous broadening due to spectral wandering, and hence their time-integrated coherent emission is reduced. To enable the present measurements, we employed a detector of higher saturation, readout rate and quantum efficiency (Princeton Instruments PIXIS:100BR), and increased the transmission in the detection path by optimized anti-reflection coatings. The paper is organized as follows. In section 2 the investigated sample and the experimental method are described. In section 3 the experimental results are given, containing the response of an individual exciton, starting with the photoluminescence (PL) in section 3.1. The delay time, real time and spectral dynamics of the FWM are investigated to infer the underlying exciton's dephasing mechanisms in section 3.2.1. The polarization state dynamics of an exciton with fine-structure splitting (FSS) is retrieved in section 3.2.2. Next, we provide textbook examples of the vectorial coherent response in a 4LS with FSS. In section 3.2.3 we measure and simulate its response in 2DFWM spectroscopy. Then, in section 3.2.4, by controlling the arrival sequence and polarization of the driving fields and resolving the polarization of the emitted fields [24], we select different FWM creation pathways involving different second-order stages, such as density grating (DG), two-photon coherence (TPC) [25,26] or Raman coherence (RC) [27,28]. These findings show quantitative agreement with the predictions of the optical Bloch equations in a 4LS, confirming the textbook behavior. In section 3.3 a characterization of the coherent response in a small ensemble of QDs employing FWM hyperspectral imaging and 2DFWM is reported.

Sample and experimental methods
The investigated GaAs QD sample was fabricated by hierarchical self-assembly [29,30], and has a layer sequence shown in figure 1(a). The sample growth starts with an undoped [001] GaAs substrate onto which a GaAs buffer layer and self-assembled InAs QDs are grown, which are then partially covered with GaAs and removed by in situ etching, providing a GaAs surface with voids. Next, an Al 0.45 Ga 0. 55 As barrier layer is grown with a thickness of 5 or 7 nm, partially retaining the voids of about 3 nm depth and tens of nanometers width. The subsequent growth of a thin GaAs layer is filling the voids and leaving a thin wetting layer. The resulting QD structure is finally overgrown with an Al 0.35 Ga 0.65 As barrier, followed by an Al 0.45 Ga 0.55 As barrier and a 10 nm GaAs cap layer. The structure was anti-reflection coated with a 98 nm thick layer of hafnium dioxide, reducing its reflectivity to the 0.1% range over the spectral range 1.65-1.7 eV.
The FWM is retrieved by the HSI technique [16] and its dual-polarization extension DHSI [24]. We use two trains of laser pulses E 1 and E 2 with a repetition rate of τ −1 r = 76 MHz. We up-shift their frequencies by acousto-optic modulators driven at 1 /(2π) = 79 MHz and 2 /(2π) = 80.77 MHz. This introduces the phase shifts e i 1 nτ r and e i 2 nτ r for the nth pulse in E 1 and E 2 , respectively. The induced degenerate FWM signal R −1,2 ∝ E * 1 E 2 2 (also see equation (A.8) in the appendix) is retrieved using the time invariance in the repetitive experiment by selecting the field with the phase shift e i(2 2 − 1 )nτ r from light reflected from the sample. This is realized by measuring the spectrally resolved interference J (2 2 − 1 , ω) between a frequency unshifted reference field and the reflection from the sample detected at the frequency of the FWM heterodyne beat 2 2 − 1 . The delay between the pulse arrival times at the samples is denoted by τ , which is positive for E 1 arriving before E 2 , as depicted in figure 5(a). The amplitude and phase of the FWM field are obtained via spectral interferometry [31,32]. More details of the experimental methodology and the corresponding notation are given in [16,22,24]. We indicate the polarization states of E 1 and E 2 and detection by a triplet ( , , +), where the individual symbols are lines along the direction of linear polarization (horizontal is along the [110] crystal axis of the sample) and circles ( ) for right-handed (left-handed) circular polarization, respectively. Dual polarization detection is indicated by +.  [30]. The data on the right of the dashed line are multiplied by 15.
All presented FWM data are taken in the χ (3) regime, with a fifth-order contribution of less than 10%, as determined by the power dependence. Single QD spectroscopy was performed on a sample with 7 nm bottom barrier with a low density (10 6 cm −2 ) of QDs. Specifically, we present results obtained on two spectrally similar QDs A and B showing different FSS of A = 61 µeV and B = 27 µeV, respectively, which are typical of these samples [33]. The statistical properties of the QDs were studied on a sample with a 5 nm bottom barrier layer and a high density (10 9 cm −2 ) of QDs. The excitation pulses are of 230 fs duration unless specified otherwise. The experiments were performed at a temperature T = 14 K. Throughout the paper we give frequency values in energy units by multiplying withh.

Photoluminescence of single quantum dots (QDs)
The investigated QDs are strongly oblate, i.e. they are significantly larger in the in-plane direction than in the growth direction (001), resulting in a large quantization energy difference between heavy-and light-hole states, with the heavy hole in the (001) direction (z) forming the lowest hole states, with two spin-orbit projections J z = ±3/2. Together with the lowest electron states of spin projection S z = ±1/2, they form two bright exciton states S z + J z = ±1 and two dark exciton states S z + J z = ±2. The two bright exciton states are mixed by non-circular confinement and the zincblende crystal structure, creating two linearly polarized eigenstates (GX and GY), spectrally separated by the FSS. Its value, , can reach up to a few hundreds of µeV.
The lowest energy two-exciton state occupies the two available electron and hole spin states, forming a spin-singlet state, which is called a biexciton B. By Coulomb and exchange interactions, it carries an energy renormalization, typically in the meV range.
Hence, the lowest excitonic and biexcitonic transitions in such QDs can be represented as a 4LS comprising the ground state G, the two linearly polarized fine-structure split exciton states X and Y, and the biexciton state B, as shown in the inset of figure 3 The excitonic transition energies in such a 4LS in a QD can be measured using nonresonantly excited PL. We discuss here the PL of QD A in the low-density sample under non-resonant excitation by a CW laser at a wavelength λ = 532 nm. At low excitation power (see figure 1(b)), the (−) polarized emission is dominated by a single line around 1697 meV attributed to the exciton transition GX. With increasing excitation power (see figure 1(c)) multiple additional emission lines are observed, which are attributed to charged excitons (below the GX energy) and charged and neutral excited state excitons some 20 meV above the GX energy, respectively. From the exciton energy we estimate [30] a confinement energy of about 200 meV. Upon a further increase of the excitation power (see insets in figure 1(b)) we identify at approximately 1691.4 meV the bound biexciton transitions XB and YB, by their complementary FSS XB-YB with respect to GX-GY. The biexciton binding energy is around 5.6 meV, which is confirmed by FWM in sections 3.2.2 and 3.3.
The trion formation, observed in figure 1(c), is due to the charging of the QDs by random capture of the non-resonantly excited electrons and holes [30]. At low excitation intensities, the capture rate is lower than the tunneling rate through the barrier to the GaAs substrate of the capping layer, and no charges build up in the QD. With increasing excitation intensities, random charging with electrons or hole occurs, with a preference for holes due to their longer tunneling times, and the time-averaged probability of the QD being in the neutral state reduces. We have assessed the temporal stability of the QD and its environment by measuring time sequences of PL spectra, as shown in figure 2. The majority of the QDs investigated on the sample with a 5 nm barrier showed substantial intensity and spectral fluctuations (within a range of a meV), as well as long-lived charging. Such features were found to virtually disappear for an excitation power below 1 µW. Conversely, the QD emission from the sample with a 7 nm barrier was stable, as presented in figure 2(b). For resonant excitation as used in the FWM measurements, the rate with which non-resonant carriers are created (e.g. by two-photon absorption or interband excitation of residual dopants) is low. The QDs are thus typically neutral. A minority of QDs showed a blinking in the seconds to minutes scale, possibly due to nearby defect states. These were excluded from further studies. The FWM experiments presented here were done on QDs which were stable during the timescale of the experiments.

Four-wave mixing of single QDs
To investigate the coherent response of a single QD, we apply FWM spectroscopy. To excite only one of the fine-structure states (X), we use a (−, −, −) polarization configuration. The excitation pulses, with a spectrum given by the green line in figure 3(b), provide a total fluence on the sample surface of 18 µJ cm −2 . We chose |E 2 | 2 = 2|E 1 | 2 , maximizing the FWM in third order for a given total fluence. The fluence is about an order of magnitude higher than that used for weakly confined QDs [22,34].
An example of a measured spectral interferogram J (2 2 − 1 , ω, τ ) is given in figure 3(a). For τ = 0 ps, wide fringes are observed, created by the quasi-instantaneous FWM of the GaAs continuum of the substrate (see also figure 13(d)). The fringe width is given by the time delay between E 2 and the reference pulse of 2.00 ± 0.05 ps. For τ = 0.4 ps, the excitation pulse overlap is negligible and thus non-resonant FWM is suppressed. The remaining FWM is dominated by the long-lived excitonic coherence. The features at 1697.2 and 1691.7 meV are due to the GX and XB transitions, respectively, consistent with the PL results shown in figure 1 The FWM responseR −1,2 − − − (ω, τ ) is retrieved in amplitude and phase from the interferogram by spectral interferometry, using the prior knowledge that the reference field is arriving earlier in time than the FWM field. The resulting response is given in figure 3(b). The amplitude shows peaks at the GX and XB transitions, as well as π-shifts of the FWM phase across the resonances. Furthermore, a mutual phase difference between GX and XB of about π is measured, as expected from an exciton-biexciton pair [22,23].
From the measured FWM amplitudes A GX and A XB (spectrally integrated over GX and XB, respectively) we retrieve the dipole moment ratio α = µ GX /µ XB between the two transitions, as detailed in equation (A.21). Note that the measured linewidths of both transitions are equal within error (the experimental error is approximated by the spectral resolution of around 15 µeV), and are dominated by inhomogeneous broadening, as shown in section 3.2.1. The linewidth is therefore not affecting the analysis of α, different from the weakly confined exciton-biexciton systems measured in [23] which showed less inhomogeneous broadening. The retrieved value of α = 1.03 is close to the strong-confinement limit of 1, in which B is a product state of two X. This is in contrast to previously studied weakly confined QDs for which this ratio was typically several times bigger [22,23].

3.2.1.
Dephasing. The dephasing of the two-level transitions can be separated into two parts. The homogeneous part, with a correlation time much shorter than the coherence time, is described by an exponential coherence decay with a decay rate of γ and decay time T 2 = 1/γ resulting in a Lorentzian lineshape with a full-width at half-maximum (FWHM) of 2γ . The inhomogeneous part, which has a correlation time much longer than the coherence time, is described by a Gaussian distribution with standard deviation σ , related to the inhomogeneous dephasing time T * 2 = 1/σ . The combination of both leads to a Voigt lineshape in the spectral domain, and a real-time dynamics of FWM described for τ > 0 by where E 2 arrives at t = 0, and ω 0 denotes its center frequency. The measured temporal transients |R −1,2 (t, τ )| 2 of GX are given in figure 4(a). The data have been corrected for the amplitude response r S (t) of the setup due to the finite spectrometer resolution (see figure 8(a)), yielding a diagram shown in figure 4(a). A shift of the signal maximum in real time according to the delay time τ is observed, in agreement with equation (1), showing the echo formation. Note that the photon echo is developing from a single transition exhibiting spectral wandering over the measurement time of the FWM experiment in the 1-1000 s range [20]. The observed FWHM of the echo amplitude is 43 ± 5 ps, which in the model equation (1) is given by √ 8 ln 2/σ , from which we deduce σ = 36 ± 4 µeV ≡ σ B . This inhomogeneous broadening is larger than the GX homogeneous broadening γ of around 8 µeV (see figure 4(d)). Using the above parameters and equation (1) the temporal transients for different delays were simulated as given in figure 4(b). The photon echo formation has been recorded with a more dynamic range for a third QD C as shown in figure 4(c). The temporal FWHM of the echo amplitude is here 48 ± 2 ps, corresponding to σ = 32 ± 1 µeV ≡ σ C , similar to σ B . The spectral wandering, generating the inhomogeneous broadening, is attributed to fluctuating charge carriers in the vicinity of the QDs, shifting the exciton energy by the quantum-confined Stark effect. Even though the FWM   (1) and prediction (red, dotdashed line) using equation (2). |R −1,2 − − − (ω XB , τ )| 2 is shown by green squares. The exponential decay with γ XB = 56 ± 8 ps is marked with a green dashed line. The noise level is shown with open circles. Both data were corrected for the spectrometer response r S (t).
excitation is resonant, the GaAs of the cap layer and the substrate is excited in its continuum, enabling free carrier dynamics and the observed spectral diffusion.
The spectrally and delay-resolved FWM of an inhomogeneously broadened two-level system can be modeled as [35] The measured |R −1,2 −−− (ω, τ )| 2 integrated in ω over the GX resonance (we used an integration range extending a FWHM to either side of the resonance in the data shown), denoted by is shown as black circles in figure 4(d). The initial rise of the FWM signal is a result of the inhomogeneous broadening σ and the related photon echo formation in real time [35], as depicted in figure 5(b). For τ σ 1 the real time response corresponds to a half-Gaussian (orange trace). For larger delays, such that τ σ 1, the photon echo is fully developed (blue trace) and the integrated FWM power shows an exponential decay rate versus τ of 4γ .
The simulation of |R −1,2 −−− (ω GX , τ )| 2 is given as a dash-dotted line in figure 4(d), and yields the parameters (γ , σ ) = (8 ± 1, 36 ± 4) µeV, corresponding to T 2 = 82 ± 3 ps and T * 2 = 18 ± 2 ps. The T 2 time of the investigated QD is about five times shorter than the typical PL lifetime [36] of T PL 400 ps. However, these PL lifetime measurements under non-resonant excitation can be influenced by dark states and can be significantly longer than the radiative lifetime of the bright state. Nevertheless, a radiative lifetime given by the T 1 = T 2 /2 41 ps is unlikely. Therefore, we assume that there is significant pure dephasing present, possibly by scattering with free carriers accumulated in the GaAs layer at 7 nm below the QD.  In figure 4(d) we also present the FWM dynamics |R −1,2 − − − (ω XB , τ )| 2 of the XB transition, denoted by green squares. The quantitative analysis is less reliable, due to the limited dynamic range of the data. Still, from the exponential decay at larger delays, we can estimate its homogeneous dephasing time T 2, XB = 58 ± 8 ps. The shorter T 2 time of XB with respect to GX indicates that the spectral fluctuations of the exciton and biexciton are not fully correlated [37].
The inhomogeneous broadening and the associated photon echo formation create a delaydependent FWM spectrum |R −1,2 − − − (ω, τ )|. For τ = 38 ps (figure 5(d)) the Gaussian echo is fully developed in the time domain, creating a Gaussian in the frequency domain (blue line), which is also measured. For τ = 0.4 ps (figure 5(c)) instead, the FWM transient has the form of a half-Gaussian (t) exp(iω GX t − σ 2 t 2 /2). The FWM spectrum is then given by its Fourier transform which reads where D is the Dawson function. The resulting FWM amplitude |L(ω)| for σ = σ B shown in figure 5(c) has Lorentzian tails (∝ |ω − ω GX | −1 ), in agreement with the measured spectrum. A smooth evolution between the two limiting spectral shapes is observed for intermediate delays (not shown). To determine the excitation fluence range for which the signal is dominated by the thirdorder response we have measured |R −1,2 − (ω GX , 2 ps)| 2 as a function of the fluence of E 1 , for a fixed fluence of E 2 of 12 µJ cm −2 , as shown in figure 6. The co-circular polarization configuration was chosen to suppress biexcitonic contributions. To reduce the effect of nonresonant two-photon absorption for a given pulse area of the excitation, we have spectrally centered the excitation to the exciton transition and narrowed it to 2.4 meV, increasing the pulse duration to 760 fs. |R −1,2 − (ω GX , 2 ps)| 2 peaks at a E 1 fluence of around 15 µJ cm −2 . The peak is expected [38] for a pulse area of π/2. Further increasing the fluence and thus the pulse area of E 1 reduces the signal, with a zero expected at a pulse area of π . The predicted fluence dependence given by equation (A.25) is shown as a violet line in figure 6(a). From such a Rabi oscillation, we estimate a transition dipole moment [39,40] of 10 ± 3 Debye, which is lower than that for single excitons weakly confined by interface fluctuations in a quantum well [38]. Most of the data in this work were taken at pulse areas of about π/5 corresponding to a fluence of about 2 µJ cm −2 for 760 fs pulses shown as a vertical dashed line in figure 6(a). For this fluence, χ (5) contributions are expected to be more than an order of magnitude below the χ (3) response. Their influence can be still observed in |R −1,2 | (ω GY , τ )| 2 and |R −1,2 | (ω YB , τ )| 2 , shown in figures 6(b) and (c) (see also figures 10 and 12). A beating with a period of about 0.75 ps corresponding to a spectral splitting of 5.5 meV is due to exciton-biexciton beating occurring in fifth order [41], and is absent in third order, as shown in the appendix.

Vectorial polarization dynamics.
To resolve the vectorial polarization dynamics of the emitted FWM we use DHSI [24]. By setting the heterodyne frequency to 1 and 2 and measuring the spectral interferences of E 1 and E 2 in J ( 1 , ω, τ ) and J ( 2 , ω, τ ), respectively, we determine the polarization state of the excitation fields, which is controlled by combinations of quarter-and half-wave plates for each excitation pulse.
The results of figure 1. To excite all transitions in the 4LS, we use the configuration ( , , +), which creates a superposition of GX and GY transitions introducing a dynamics of the polarization state [24] observable both in the emitted FWM polarization and in a beat versus τ , as presented in figure 7(b) for (−) detection. The quantum beat (black dots), which is observed on both GX and GY (latter not shown) with a period of T = 2π/ = 65 ± 5 ps, which yields = 63 ± 5 µeV, in agreement with the splitting measured in the FWM and PL spectra. Such an agreement of FSS and beat period was also observed for QD B, exhibiting B = 27 ± 7µeV and, respectively, T = 153 ± 39 ps.
The delay time dynamics of the FWM at XB instead shows only the dynamics of the corresponding excitonic resonances GX (see blue triangles in figure 7(b)) since it is created by the DG on X. The equivalent holds for the GY and YB transitions (not shown). The coupling of X and Y via G is not influencing the FWM signal at the biexciton transition in third order.
For co-circular configuration ( , , +) a coherence dynamics similar to ( , , +) is found for the excitonic transitions (not shown), while biexcitonic transitions do not emit, as shown in figure 7(a) (bottom), independent of τ . Even though the FSS transfers first-order coherence to the polarization, at the arrival of the second pulse, it cannot be converted into a DG by the polarized second pulse, excluding FWM at the biexcitonic transitions (see table 1).
In DHSI the relative phase between | and − FWM components can be measured such that the time-dependent polarization state can be reconstructed [24]. Here we determine this dynamics for the excitons of QDs A and B exhibiting considerably different FSSs. We use both ( , , +) and ( , , +) configurations, exciting different initial relative phases of GX and GY transitions. The real-time dynamics of the polarization-resolved FWM amplitudes |R −1,2 + (t, 0.5 ps)| and |R −1,2 + (t, 0.5 ps)| is shown in figures 8(a) and (b). The finite spectrometer resolution yields the instrument response r S (t) shown by the dashed-dotted line. To reduce noise and eliminate the emission of XB and YB, the data were filtered in the spectral domain with a Gaussian (FWHM of 0.62 meV) centered on the GX transition, resulting in a temporal resolution of 3 ps, as depicted by the horizontal bar in figure 8(a). Table 1. Relative third-order FWM amplitudes of an exciton-biexciton system with a FSS (4LS), as displayed in figure 3(b), for various excitation and detection polarizations for positive (>) and negative (<) delays τ . GX and GY are transitions from the ground state to (−) polarized and (|) polarized excitons, respectively. Transitions from these excitons toward the biexciton are assigned as XB and YB, respectively. We consider a dipole moment ratio between excitonic and biexcitonic transitions α = µ GX /µ XB = µ GY /µ YB with µ XB µ GX = µ YB µ GY = µ 2 and the ratio between the excitation field amplitudes at biexcitonic and excitonic transitions, β = |E(ω XB )/E(ω GX )|. The table is normalized with respect to equation (A. 19) in the appendix.
It is interesting to note that the phase slope (figures 8(c) and (d)) is constant over approximately 70 ps, in spite of the inhomogeneous broadening σ by spectral wandering. For QD B, this amounts to σ √ 8 ln 2 = 100 µeV 4 . This shows that the spectral wandering is affecting GX and GY similarly, such that their relative phase is conserved over a timescale longer than σ −1 . It is expected for the spectral wandering mechanism due to locally fluctuating electric fields, as the two transitions GX and GY involve exciton states which are different only in their spin component, but have the same spatial exciton wave function. The latter determines the static dipole moment, which governs the spectral wandering by its interaction energy with an electrostatic field. Conversely, the spatial wave functions of exciton and biexciton are typically different [42]. Therefore electric fluctuations influence them differently, such that spectral wandering of both is not completely correlated, as previously concluded from |R −1,2 − − − (ω XB , τ )| 2 shown in figure 4(d).
The time-resolved polarization state can be represented [24] by the angles (ψ, χ) of the polarization ellipse, describing azimuth angle and ellipticity. These are plotted in figures 8(e) and (f) and show agreement with the theoretical prediction for equal dipole moments of the transitions. A sketch of the time evolution of the polarization state is given in figure 8(g). For QD A and ( , , +) (upper row) we observe the evolution from the excited ( ) polarization

Dual polarization 2DFWM.
The exciton-biexciton 4LS gives rise to a characteristic 2DFWM, as was shown for quantum wells [43]. The 2DFWM representation of the data is obtained by Fourier-transformingR −1,2 (ω, τ ) for τ > 0 from τ to ω 1 yielding |Ř −1,2 (ω, ω 1 )| 2 . We use here the data set analyzed in figure 7(b). For τ > 0 the biexciton FWM at XB (YB) is driven by GX (GY), respectively. We can therefore use the FWM on the XB (YB) transitions as phase reference [22]. The resulting polarization-resolved 2DFWM |Ř −1,2 − (ω, ω 1 )| 2 and |Ř −1,2 | (ω, ω 1 )| 2 is shown in figure 9(a). We observe that the FWM on GX and GY is created by both GX and GY (peaks at (ω, . This coherent coupling of GX and GY is due to their common ground state. The FWM on XB (YB) is only created by the first-order coherence on GX (GY), respectively. The experimental data are reproduced by the theoretical simulations shown in figure 9(b). As discussed in the appendix, the modeling is based on solving time evolution of the density matrix for the exciton-biexciton system, with delta-like excitation pulses [44]. In addition, the polarization control of E 1 , E 2 was introduced and the resulting field is resolved in both linear polarizations. Specifically, in simulations shown in figure 9(b) we have assumed equal dipole moments and equal driving of all four transitions, = 61 µeV, a biexciton binding energy of 5.6 meV and a spectral resolution of 15 µeV FWHM. We have set the pulse area of π/50 for E 1 and π/25 for E 2 . Note a reduction of biexcitons' intensity with respect to excitons' ones, as expected (see table 1).

Polarization selection rules.
The strong exciton confinement in the QDs studied here creates a 4LS without a nearby biexcitonic continuum, with defined polarization selection rules as derived in the appendix. We demonstrate here on this textbook four-level scheme the different FWM creation pathways. They can be separated by polarization-selective excitation and detection for positive and negative delays [41,43,45,46], realized by τ = ±0.5 ps in the measurements. The polarization selection rules in a 4LS in the framework of optical Bloch equations truncated to third order in the pulse areas are given in table 1. We have introduced an exciton-biexciton dipole moment ratio α and an excitation field amplitudes ratio β between the biexciton and exciton. In the following discussion of the relative FWM amplitudes, we assume α = 1 and β = 1.
Using linear polarization (−) of E 1 (see figure 10) we excite transitions GX and XB. For τ > 0, E 1 comes first, creating a GX coherence. If E 2 is co-linearly polarized, it creates a DG in G and X, and drives the FWM on both GX and XB transitions, with the amplitude of XB half of GX, as the former is created only by the X grating. For τ < 0, E 2 arrives first and creates a TPC [ 25,26] between G and B, which is converted by E 1 into FWM on GX and XB of equal amplitudes. Using a cross-linearly (|) polarized E 2 , the FWM is unchanged for τ < 0, as the TPC remains the same. For τ > 0, in contrast, E 2 excites only the GY transition and no DG is formed. Instead, a second-order RC [27,28] between GX and GY is created, and converted by E 2 , driving the YB transition to FWM on the XB transition. Exchanging (−) and (|) polarizations for all fields results in an equivalent finding in experiment and theory (not shown), with the role of X and Y exchanged.
In the third-order limit, FWM of a 4LS can therefore originate from three different, secondorder components: DG, RC and TPC. The corresponding FWM creation pathways are illustrated at the right sides of figures 10 and 12. Alternatively, they can be represented by double-sided Feynman diagrams [47], as shown in figure 11. In such diagrams, density matrix elements carrying DG, RC and TPC can be readily identified.
Moving to circular polarization ( ) of E 1 (see figure 12), we now excite both excitonic transitions with a π/2 phase shift. The delay τ = ±0.5 ps is sufficiently small to neglect the FS dynamics during the time between the arrival of the two pulses, so that one can discuss the transitions in terms of excitons with spin ±1, |P = (|X + i|Y )/ √ 2 and |M = (|X − i|Y )/ √ 2. For co-circular polarization ( ) of E 2 , no biexcitonic FWM is observed (see figure 7(a), bottom left), due to the alternate circular polarization of the GP ( ) and PB ( ) transitions. For the same reason, no TPC can be created by E 2 and no FWM is observed for negative delay. The excitonic transition is initially ( ) polarized, and shows a polarization beat due to the FS splitting (i.e. coupling to M). For cross-circular polarization ( ) of E 2 , no FWM is created, as the RC cannot be converted into FWM via the B.  For circular polarization ( ) of E 1 and linear polarization (−) of E 2 , different FWM pathways are open. For positive delay, we form a DG on G and X, which results in FWM on GX and XB twice larger than that on GX. We also form an RC between X and Y, resulting in a FWM signal on YB. For negative delay, TPC on GB is created by E 2 , which is converted into FWM on all four transitions with equal amplitude by the circularly polarized E 1 .
A selection of experimental results shown in figures 10 and 12 shows qualitative agreement with the predictions. For a quantitative comparison, we use previously retrieved α = 1.03 and β = 0.87 taken from the measured excitation spectrum (see the green line in figure 10(a)). Small deviations from the third-order response are observed. Specifically, in figure 10(b) we observe that for (−, |, −) the amplitude of XB is about 1.3 times larger for τ = +0.5 ps than for τ = −0.5 ps, while it is predicted to be equal. Additionally, a weak FWM is observed on GX for τ = +0.5 ps. The latter is absent in the FWM in the third-order response. Such a signal is, however, allowed by taking into account fifth-order contributions. The FWM analytical modeling of the 4LS, which is described in the appendix, includes all orders. The outcome is presented by red traces in figures 10 and 12, using for all eight spectra one common absolute scaling as the only adjustable parameter. Quantitative agreement between the simulations and the experimental results is achieved.

Four-wave mixing imaging and two-dimensional four-wave mixing
To characterize the spatial distribution and the statistical properties of the investigated QDs, we use FWM hyperspectral imaging [22,34], measuring the FWM power |R −1,2 −−− (ω, τ ; r)| 2 as a function of the spatial position r = (x, y). Results on the high-density sample at τ = 0 are shown in figures 13(a)-(c). Individual excitons are identified by their localized response in spatial and spectral domains. In the χ (3) regime, the imaging has a spatial resolution of about 300 nm. At higher excitation power beyond the third-order response, the spatial resolution gradually deteriorates [34]. The FWM power, spatially averaged over a 10 µm 2 area, |R −1,2 −−− (ω, τ ; r)| 2 r , is shown in figure 13(d). Some 70 individual FWM resonances are detected within the energy range of 15 meV. For τ = 0 we observe a spectrally broad background (see the orange, dashed line) attributed to FWM from the GaAs substrate and cap layer, which is excited 100 meV above its band gap. This signal decays on a 100 fs timescale [48] via LO-phonon and carrier-carrier scattering of the underlying interband electron-hole excitations. This quasi-instantaneous signal is removed from the measured FWM by time-filtering the FWM, retaining only the signal for t > 0.5 ps (black line in figure 13(d)). Such a posteriori time or spectral filtering is enabled by the retrieval of the FWM field in amplitude and phase and provides a powerful tool to discriminate between different components of the signal [17].
The observed spatial and spectral vicinity of QDs indicates that coherent coupling might occur between individual X. To measure coherent coupling, we perform 2DFWM spectroscopy that efficiently retrieves such mutual couplings [16,19,22,23]. The resulting |Ř −1,2 −−− (ω, ω 1 ; r)| 2 at r = (0, 0), as marked by the white crosses in figures 13(a)-(c), is shown in figure 14. We observe peaks from exciton transitions GX 1..9 . They are aligned along the diagonal ω = ω 1 (orange dashed line and magenta dashed line). The multiple diagonals are due to the discrete Fourier transform along the equidistant delay step of τ = 0.5 ps, resulting in a determination of ω 1 only up to integer multiples 2π/ τ 8.27 meV. A second set of red-shifted peaks (dotted, orange line) is observed for each exciton at ω = ω 1 − δ B , which are biexciton transitions XB 1..4 with a binding energy of δ B 5.5 meV.
Coherent coupling between excitons creates off-diagonal signals connecting the corresponding diagonal exciton signals. No such features are observed in figure 14, despite the fact that all observed excitons are located within the excitation focus of 0.5 µm diameter. To determine spatial positions of the GX transitions we fitted the data cube |R −1,2 (ω, 0; r)| 2 shown in figure 13 with a Lorentzian in the spectral domain and a 2D Gaussian in real space, yielding a position with an accuracy of about 50 nm, an order of magnitude better than the spatial resolution. The positions are shown in the inset of figure 14, and show distances in the 200 nm range. The lack of coherent coupling over such distances is consistent with the predicted strength of the radiative coupling between QDs [49] in the µeV range. We conclude that the coherent response of these excitons is uncoupled. Absence of coupling among these strongly confined excitons is in contrast to the results on weakly confined QDs [22], where coherent coupling over distances up to a micrometer was observed. We note that the strength and range of the radiative coupling between QDs is predicted to be enhanced by up to two orders of magnitude in photonic crystal structures [50]. Coherent coupling of different excitonic states within a single QD instead is strong, as shown by the biexciton renormalization, and the coupling of fine-structure split states shown in section 3.2.3.

Conclusions
Improving the sensitivity of the heterodyne spectral interferometry technique and its dualpolarization extension enabled a comprehensive investigation of the vectorial, coherent, nonlinear response in an exciton-biexciton model system in individual GaAs/AlGaAs strongly confined QDs. Using FWM measurements we have determined the exciton coherence time and assessed underlying homogeneous and inhomogeneous dephasing mechanisms. We have studied the coherent dynamics of the exciton-biexciton 4LS and identified fifthorder contributions to the FWM coherent dynamics and vectorial response. We have shown, both experimentally and theoretically, different FWM creation pathways, and their polarization selection rules. The results are promising for the application of heterodyne spectral interferometry to a large range of QD systems, like site-controlled QDs [51][52][53][54] or QD molecules [55,56] and-in a broader context-other individual optical transitions in solid state.
Following [44] we model the time dependence of the laser field amplitude by a δ function. This approximation holds in the limit of pulses with the same temporal shape for both polarization components of the jth pulse and shorter than the characteristic electron-phonon coupling time. However, figures 3(b) and 10(a) show that the spectral width of the pulses is finite, driving biexciton and exciton transitions differently. In order to take this into account we introduce the ratio β between the excitation field amplitudes at biexcitonic and excitonic transitions. Using this notation we have E ( j) ω XB = βE ( j) ω GX = β(E j /2)e iφ j , the label i has been dropped for brevity and j has been lowered to lighten the notation.
This analytical model considers a single repetition of the experiment, assigning a relative phase relation between the two exciting pulses with two phases φ j . In practice, this phase relation between the pulses is set by optical heterodyning. The experimental implementation (see section 2) employs pulse trains, with a repetition rate τ −1 r , which are phase shifted by φ j = j τ r by acousto-optic modulators operating at j . We define the arrival time of the first pulse as t = −τ and that of the second pulse as t = 0, as illustrated in figure 5(a). We denote a delay as positive (negative) when the pulse with φ 1 arrives first (second), but always use τ > 0 in the following notation.
We use the homogeneous broadenings γ GX , γ GY , γ XB , γ YB of the transitions given by the imaginary part of the corresponding resonance frequency. Even though inhomogeneous broadenings and phonon-induced dephasing are observed in the experiment, here they are disregarded for simplicity.
In order to model circular and linear polarizations for the exciting fields and to compare quantitatively the spectra, we express the electric fields in circular polarization and take the circular-linear transformation explicitly. Being contravariant components, the change of basis gives the following relations for the field components: With this choice, x-polarized light (−) will be modeled with E j = E j = E j / √ 2 and φ j = φ j = φ j , where for brevity we dropped the x label in E and we lowered the j related to the time ordering of the pulse. The y-polarization (|) corresponds to E j / √ 2 = E j = E j and φ j = φ j = φ j + π (dropping the corresponding y label). In general, the optical polarization is a sum over contributions from all transitions between states with electron-hole pair numbers differing by one [57]. For instance, the k-component of the electronic polarization reads where the density matrix elements are defined over all four states of the 4LS and the sum is understood for those ν to ν transitions parallel to the k-direction (i.e. if k = x, then the sum will be over GX and XB). In the present context the two-time FWM response can be extracted from the quantum mechanical electronic polarization as where ≷ stands for positive and negative delay, respectively. Its Fourier transformation with respect to the real time t yields the calculated spectra shown in figures 10 and 12 as This model admits rather compact analytical solutions at all orders in the exciting electric fields, giving an access to specific microscopic parameters (homogeneous T 2 and inhomogeneous T * 2 broadening of the transitions, and-to some extent-electric dipole moments) by a fitting procedure. The inclusion of higher-order nonlinear terms (specifically the fifth order) was necessary to attain quantitative agreement with experimental results. For instance, in (−, |, −) and τ > 0 configuration shown in figure 10(b) (top), a weak yet nonnegligible FWM signal on GX transition is explained by a fifth-order contribution and is not predicted in the third-order regime. We note that having found good parameters for the best fit in figure 10(a) (and of figure 4(d) for the analysis of GX dephasing), without any further adjustment (with only one common absolute scaling) our simulations and data superposed in all other spectra in figures 10 and 12.
In particular, the density matrixρ ν,ν = |ν ν | in equation (A.7) evolves according to the equation ih∂ tρ = [H,ρ]. We follow the step indicated in [44] in the form used in [22], by assuming delta-like pulses and a simple Lorentzian dephasing model, the time evolution of the density matrix and the FWM response can be calculated analytically. Explicitly, the free time evolution between two successive pulses is given by where ν,ν = ω ν − ω ν + iγ ν,ν and γ ν,ν are phenomenological damping rates, matrix elements of  The arrival of the jth pulse changes the density matrix to inside the brackets the matrix product is assumed and In general, the FWM response is extracted as shown in equation (A.8); in the light of the exact analytic solution the terms in the density matrix evolution contributing to the FWM signal are those proportional to exp{i(2φ 2 − φ 1 )}.
In the following, we summarize the most relevant results, including third-order limits in the excited fields E 1 and E 2 given in table 1. We first derive the two-time solutions in (−, −, −) polarization. We then present the solutions for spectrally resolved FWM signals at specific delays shown in figures 10 and 12. For simplicity we take µ GX = µ GY = µ and µ XB = µ YB = µ B , and define M = µ 2 + β 2 µ 2 B , ω XB = ω B − ω X , ω YB = ω B − ω Y , ω GX = ω X − ω G and ω GY = ω Y − ω G .
The real time and delay FWM response for positive delays and upon (−, −, −) polarization configuration reads R −1,2 − − − (t, τ ) > = θ (t) e −i(ω GX −iγ GX )t −e i(ω GX +iγ GX )τ C X GX + e i(ω XB +iγ XB )τ C X XB +e −i(ω XB −iγ XB )t e i(ω GX +iγ GX )τ C B GX − e i(ω XB +iγ XB )τ C B XB (A.14) with C X GX = i µ 4 M 5 sin The pulse E 1 arrives first at t = −τ and in the lowest order drives the GX coherence; this is the sin( 1 2 E 1 M) term in C X GX and C B GX . Between the two pulses (−τ, 0 − ) the system undergoes a free evolution giving the exponentials in τ at the transition that E 1 drove, at its corresponding order. In the degenerate FWM configuration employed here, the pulse E 2 arriving at t = 0 triggers a two-fold action: it creates a density grating (DG) in the GX transition and drives the FWM response at both resonances. The terms proportional to exp(−iω GX t) and exp(−iω XB t) correspond to the FWM on GX and XB, respectively. This solution is shown by Feynman diagrams in figure 11(a).
For negative delay, see diagrams in figure 11(c), one can readily identify the TPC that E 2 creates, the sin( 1 4 E 2 M) 2 term, with the lowest order of two in E 2 , and notice how E 1 , arriving second, converts this coherence into a FWM on GX and XB, which has up to the factor β equal amplitude in lowest order: The dipole moment ratio α can be retrieved from the ratio between the FWM amplitudes of the exciton and biexciton in third order in E 1 and E 2 , for which equation (A.17) reads We note that the normalization factor in table 1 reads 1 √ 2π Assuming an equal broadening for the two transitions, the ratio of the FWM spectral amplitudes integrated around exciton and biexciton transitions