\chapter{Measurement of Ultrafast Excitonic Dynamics of Few-Layer MoS$\mathbf{_2}$ using
  State-Selective Coherent Multidimensional Spectroscopy} \label{cha:mx2}

\textit{This Chapter borrows extensively from \textcite{CzechKyleJonathan2015a}. The authors are:
  \begin{denumerate}
    \item Kyle J. Czech
    \item Blaise J. Thompson
    \item Schuyler Kain
    \item Qi Ding
    \item Melinda J. Shearer
    \item Robert J. Hamers
    \item Song Jin
    \item John C. Wright
  \end{denumerate}
}

\clearpage

We report the first coherent multidimensional spectroscopy study of a MoS\textsubscript{2} film.  %
A four-layer sample of MoS\textsubscript{2} was synthesized on a silica substrate by a simplified
sulfidation reaction and characterized by absorption and Raman spectroscopy, atomic force
microscopy, and transmission electron microscopy.  %
State-selective coherent multidimensional spectroscopy (CMDS) on the as-prepared
MoS\textsubscript{2} film resolved the dynamics of a series of diagonal and cross-peak features
involving the spin---orbit split A and B excitonic states and continuum states.  %
The spectra are characterized by striped features that are similar to those observed in CMDS
studies of quantum wells where the continuum states contribute strongly to the initial excitation of
both the diagonal and cross-peak features, while the A and B excitonic states contributed strongly
to the final output signal.  %
The strong contribution from the continuum states to the initial excitation of both the diagonal
and cross-peak features, while the A and B excitonic states contributed strongly to the final
output signal.  %
The strong contribution from the continuum states to the initial excitation shows that the continuum
states are coupled to the A and B excitonic states and that fast intraband relaxation is occurring
on a sub-70 fs time scale.  %
A comparison of the CMDS excitation signal and the absorption spectrum shows that the relative
importance of the continuum states is determined primarily by their absorption strength.  %
Diagonal and cross-peak features decay with a 680 fs time constant characteristic of exciton
recombination and/or trapping.  %
The short time dynamics are complicated by coherent and partially coherent pathways that become
important when the excitation pulses are temporally overlapped.  %
In this region, the coherent dynamics create diagonal features involving both the excitonic states
and continuum states, while the partially coherent pathways contribute to cross-peak features.  %

\section{Introduction}  % =========================================================================

Transition metal dichalcogenides (TMDCs), such as MoS\textsubscript{2}, are layered semiconductors
with strong spin-orbit coupling, high charge mobility, and an indirect band gap that becomes direct
for monolayers. \cite{WangQingHua2012a, MakKinFai2010a}  %
The optical properties are dominated by the A and B excitonic transitions between two HOMO
spin-orbit split valence bands and the lowest state of the conduction band at the $K$ and $K^\prime$
valleys of the two-dimensional hexagonal Brillouin zone. \cite{MolinaSanchezAlejandro2013a}  %
The spin and valley degrees of freedom are coupled in individual TMDC layers as a result of the
strong spin-orbit coupling and the loss of inversion symmetry.  %
The coupling suppresses spin and valley relaxation since both spin and valley must change in a
transition.  %
These unusual properties have motivated the development of TMDC monolayers for next-generation
nano/optoelectronic devices as well as model systems for spintronics and valleytronics
applications. \cite{MakKinFai2010a, XuXiaodong2014a, XiaoDi2012a}  %

Ultrafast dynamics of the MoS\textsubscript{2} A and B electronic states have been measured by
pump-probe, transient absorption, and transient reflection spectroscopy. \cite{FangHui2014a,
  KumarNardeep2013a, NieZhaogang2014a, SunDezheng2014a, SimSangwan2013a}  %
The spectra contain A and B excitonic features that result from ground-state bleaching (GSB),
stimulated emission (SE), and excited-state absorption (ESA) pathways.  %
The excitons exhibit biexponential relaxation times of $\approx$10--20 and $\approx$350--650 fs,
depending on the fluence and temperature.  %
The dependence on excitation frequency has not been explored in previous ultrafast experiments on
MoS\textsubscript{2}, but it has played a central role in understanding exciton cooling dynamics
and exciton-phonon coupling in studies of quantum dots. \cite{KambhampatiPatanjali2011a}  %

Coherent multidimensional spectroscopy (CMDS) is a complementary four wave mixing (FWM) methodology
that differs from pump-probe, transient absorption, and transient reflection methods.
\cite{XiaoDi2012a, FangHui2014a, KumarNardeep2013a, NieZhaogang2014a, SimSangwan2013a,
  MakKinFai2012a, SunDezheng2014a}  %
Rather than measuring the intensity change of a probe beam caused by the state population changes
induced by a pump beam, CMDS measures the intensity of  a coherent output beam created by
interactions with three excitation pulses.  %
The interest in CMDS methods arise from their ability to remove inhomogeneous broadening, define
interstate coupling, and resolve coherent and incoherent dynamics. \cite{CundiffStevenT2008a,
  TurnerDanielB2009a, KohlerDanielDavid2014a, YursLenaA2011a, GriffinGrahamB2013a, HarelElad2012a,
  CundiffStevenT1996a, BirkedaDl1996a, WehnerMU1996a}  %
CMDS typically requires interferometric phase stability between excitation pulses, so CMDS has been
limited to materials with electronic states within the excitation-pulse bandwidth.  %
Multiresonant CMDS is a particularly attractive method for the broader range of complex materials
because it does not require interferometric stability and is able to use independently tunable
excitation pulses over wide frequency ranges.  %

The multiresonant CMDS used in this work employs two independently tunable excitation beams with
frequencies $\omega_1$ and $\omega_2$.  %
The $\omega_2$ beam is split into two beams, denoted by $\omega_2$ and $\omega_2^\prime$.  %
These three beams are focused onto the MoS\textsubscript{2} thin film at angles, creating an output
beam in the phase-matched direction
$\mathbf{k}_{\mathrm{out}}=\mathbf{k}_1-\mathbf{k}_2+\mathbf{k}_{2^\prime}$ where $\mathbf{k}$ is
the wave vector for each beam and the subscripts label the excitation frequencies.  %
Multidimensional spectra result from measuring the output intensity dependence on frequency and
delay times.  %

\autoref{fig:Czech01} introduces our conventions for representing multidimensional spectra.  %
\autoref{fig:Czech01}b,d are simulated data.  %
\autoref{fig:Czech01}a shows one of the six time orderings of the three excitation pulses where
$\tau_{22^\prime}\equiv t_2-t_{2^\prime}>0$ and $\tau_{21}\equiv t_2-t_1<0$; that is, the
$\omega_{2^\prime}$ pulse interacts first and the $\omega_1$ pulse interacts last.  %
\autoref{fig:Czech01}b illustrates the 2D delay-delay spectrum for all six time orderings when
$\omega_1$ and $\omega_2$ are both resonant with the same state.  %
The color denotes the output amplitude.  %
Along the negative ordinate where $\tau_{22^\prime}=0$, interactions with the $\omega_2$ and
$\omega_{2^\prime}$ pulses create a population that is probed by $\omega_1$.  %
Similarly, along the negative absicissa where $\tau_{21}=0$, interactions with the $\omega_2$ and
$\omega_1$ pulses create a population that is probed by $\omega_{2\prime}$.  %
The decay along these axes measures the population relaxation dynamics.  %
Note that these delay representations differ from previous publications by our group.
\cite{PakoulevAndreiV2006a}  %
This paper specifically explores the dynamics along the ordinate where $\tau_{22^\prime}$ is zero
and the $\tau_{21}$ delay is changed.  %

\autoref{fig:Czech01}c depicts the A and B excitonic transitions between the spin-orbit split
valence bands and the degenerate conduction band states of MoS\textsubscript{2}.  %
\autoref{fig:Czech01}d illustrates the 2D frequency-frequency spectrum when $\omega_1$ and
$\omega_2$ are scanned over two narrow resonances.  %
The spectrum contains diagonal and cross-peaks that we label according to the excitonic resonances
AA, AB BA, and BB for illustrative purposes.  %
The dynamics of the individual quantum states are best visualized by 2D frequency-delay plots,
which combine the features seen in \autoref{fig:Czech01}b,d.  %

This works reports the first multiresonant CMDS spectra of MoS\textsubscript{2}.  %
It includes the excitation frequency dependence of the A and B excitonic-state dynamics.  %
These experiments provide a fundamental understanding of the multidimensional MoS\textsubscript{2}
spectra and a foundation for interpreting CMDS experiments on more complex TMDC
heterostructures.  %
The experimental spectra differ from the simple 2D spectrum shown in \autoref{fig:Czech01}d and
those of earlier CMDS experiments with model systems.  %
The line shape of the CMDS excitation spectrum closely matches the absorption spectrum, but the
line shape of the output coherence is dominated by the A and B excitonic features.  %
The difference arises from fast, $<70$ fs intraband relaxation from the hot A and B excitons of the
continuum to the band edge.  %
A longer, 680 fs relaxation occurs because of trapping and/or exciton dynamics.
\cite{SimSangwan2013a}  %
The intensity of the cross-peaks depends on the importance of state filling and intraband
relaxation of hot A excitons as well as the presence of interband population trnasfer of the A and
B exciton states.  %

\begin{figure}
  \includegraphics[width=0.5\textwidth]{MX2/01}
  \caption[Introduction to CMDS.]{
    (a) Example delays of the $\omega_1$, $\omega_2$, and $\omega_{2^\prime}$ excitation pulses.
    (b) Dependence of the output intensity on the $\tau_{22^\prime}$ and $\tau_{21}$ time delays
    for $\omega_1=\omega_2$.
    The solid lines define the regions for the six different time orderings of the $\omega_1$,
    $\omega_2$, and $\omega_{2^\prime}$ excitation pulses.
    We have developed a convention for numbering these time orderings, as shown.
    (c) Diagram of the band structure of MoS\textsubscript{2} at the $K$ point.
    The A and B exciton transitions are shown. (d) Two dimensional frequency-frequency plot
    labeling two diagonal and cross-peak features for the A and B excitons.}
  \label{fig:Czech01}
\end{figure}

\section{Methods}  % ==============================================================================

MoS\textsubscript{2} thin films were prepared \textit{via} a Mo film sulfidation reaction, similar
to methods reported by \textcite{LaskarMasihhurR2013a}.  %
A 1 nm amount of Mo (Kurt J. Lesker, 99.95\%) metal was electron-beam evaporated onto a fused
silica substrate at a rate of 0.05 \AA/s.  %
The prepared Mo thin films were quickly transferred to the center of a 1-inch fused silica tub
furnace equipped with gas flow controllers (see \autoref{fig:CzechS1}) and purged with Ar.  %
The temperature of the Mo substrate was increased to 900 $^\circ$C over the course of 15 min, after
which 200 mg of sulfur was evaporated into the reaction chamber.  %
Sulfidation was carried out for 30 min, and the furnace was subsequently cooled to room
temperature; then the reactor tube was returned to atmospheric pressure, and the
MoS\textsubscript{2} thin film samples were collected.  %
The MoS\textsubscript{2} samples were characterized and used for CMDS experiments with no further
preparation.  %

MoS\textsubscript{2} thin film absorption spectra were collected by a Shimadzu 2401PC
ultraviolet-visible spectraphotometer.  %
Raman and photoluminescence experiments were carried out in parallel using a Thermo DXR Raman
microscope with a 100x 0.9 NA focusing objective and a 2.0 mW 532 nm excitation source.  %
Raman/PL measurements were intentionally performed at an excitation power of $<$8.0 mW to prevent
sample damage. \cite{CastellanosGomezA2012a}  %
Contact-mode atomic foce microscopy was performed with an Agilent 5500 AFM.  %
MoS\textsubscript{2} film thickness was determined by scratching the sample to provide a clean
step-edge between the MoS\textsubscript{2} film and the fused silica substrate.  %
TEM samples were prepared following the method outlined by \textcite{ShanmugamMariyappan2012a}
using concentrated KOH in a 35 $^\circ$C oil bath for 20 minutes.  %
The delaminated MoS\textsubscript{2} sample was removed from the basic solution, rinsed five times
with DI water, and transferred to a Cu-mesh TEM grid.  %
TEM experiments were performed on a FEI Titan aberration corrected (S)TEM under 200 kV accelerating
voltage.  %

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/S1}
  \caption[Schematic of sulfidation reaction chamber.]{
    Schemiatic of the synthetic setup used for Mo thin film sulfidation reactions.
  }
  \label{fig:CzechS1}
\end{figure}

The coherent multidimensional spectroscopy system used a 35 fs seed pulse, centered at 800 nm and
generated by a 1 kHz Tsunami Ti-sapphire oscillator.  %
The seed was amplified by a Spitfire-Pro regenerative amplifier.  %
The amplified output was split to pump two TOPAS-C collinear optical parametric amplifiers.  %
OPA signal output was immediately frequency doubled with BBO crystals, providing two $\approx$50 fs
independently tunable pulses denoted $\omega_1$ and $\omega_2$ with frequencies ranging from 1.62
to 2.12 eV.  %
Signal and idler were not filtered out, but played no role due to their low photon energy.  %
Pulse $\omega_2$ was split into pulses labeled $\omega_2$ and $\omega_{2^\prime}$ to create a total
of three excitation pulses.  %

In this experiment we use motorized OPAs which allow us to set the output color in software.  %
OPA1 and OPA2 were used to create the $\omega_1$ and $\omega_2$ frequencies, respectively.  %
In \autoref{fig:CzechS4} we compare the spectral envelope generated by the OPA at each set
color.  %
Negative detuning values correspond to regions of the envelope lower in energy than the
corresponding set color.  %
The colorbar allows for comparison between set color intensities.  %
The fluence values reported correspond to the brightest set color for each OPA.  %
A single trace of OPA2 output at set color = 1.95 eV can be found in \autoref{fig:Czech02}.

After passing through automated delay stages (Newport SMC100 actuators), all three beams were
focused onto the sample surface by a 1 meter focal length spherical mirror in a distorted BOXCARS
geometry to form a 630, 580, and 580 $\mu$m FWHM spot sizes for $\omega_1$, $\omega_2$, and
$\omega_{2^\prime}$, respectively.  %

\autoref{fig:CzechS5} represents delay corrections applied for each OPA.  %
The corrections were experimentally determined using driven FWM output from fused silica.  %
Corrections were approximately linear against photon energy, in agreement with the normal dispersion
of transmissive optics inside our OPAs and between the OPAs and the sample.  %
OPA2 required a relatively small correction along $\tau_{22^\prime}$ (middle subplot) to account
for any dispersion experienced differently between the two split beams.  %
OPA1 was not split and therefore needed no such correction.  %

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/S4}
  \caption{OPA outputs at each color explored.}
  \label{fig:CzechS4}
\end{figure}

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/S5}
  \caption{Spectral delay correction.}
  \label{fig:CzechS5}
\end{figure}

\autoref{fig:Czech10}a represents the to-scale mask that defines our distorted BOXCARS
configuration.  %
Relative to the center of the BOXCARS mask (small black dot), $\omega_1$, $\omega_2$, and
$\omega_{2^\prime}$ enter the sample at angles of 5.0, 1.5, and 1.0 degrees.  %
Each is angled only along the vertical or horizontal dimension, as indicated in
\autoref{fig:Czech10}a.  %
This distortion allowed us to remove a large amount of unwanted $\omega_2$ and $\omega_{2^\prime}$
photons from our signal path (\autoref{fig:Czech10}a red star).  %
$\omega_1$ photons were less efficiently rejected, as we show below.  %
The center of the BOXCARS mask was brought into the sample at $\approx$45 degrees.  %
All three beams had S polarization.  %
After reflection, the output beam was isolated using a series of apertures, spectrally resolved
with a monochromator (spectral resolution 9 meV). and detected using a photomultiplier (RCA
C31034A).  %

Our experimental setup allowed for the collection of both transmissive and reflective
(epi-directional) FWM signal.  %
The 2D delay spectra in \autoref{fig:Czech10}b show the presence of a large nonresonant
contribution at the origin for the transmissive FWM signal and weaker signals from the
MoS\textsubscript{2} thin film at negative values of $\tau_{21}$ and $\tau_{22^\prime}$.  %
The nonresonant contribution is much weaker than the signals from the film for the reflective
signal id is the geometry chosen for this experiment.  %
This discrimination between a film and the substrate was also seen in reflective and transmissive
CARS microscopy experiments. \cite{VolkmerAndreas2001a}  %


\begin{figure}
  \includegraphics[width=\textwidth]{MX2/10}
  \caption[Mask and epi vs transmissive.]{
    (a) Mask.
    (b) 2D delay spectra at the BB diagonal ($\omega_1=\omega_2\approx1.95$ eV) for transmissive
    and reflective geometries.
    Transmissive signal is a mixture of MoS\textsubscript{2} signal and a large amount of driven
    signal from the substrate that only appears in the pulse overlap region. Reflective signal is
    representative of the pure MoS\textsubscript{2} response.}
  \label{fig:Czech10}
\end{figure}

Once measured, the FWM signal was sent through a four-stage workup process to create the data set
shown here.  %
This workup procedure is visualized in \autoref{fig:Czech11}.  %
We use a chopper and boxcar in active background subtraction mode (averaging 100 laser shots) to
extract the FWM signal from $\omega_1$ and $\omega_2$ scatter.  %
We collect this differential signal (\autoref{fig:Czech11}b) in software with an additional 50
shots of averaging.  %
In post-process we subtract $\omega_2$ scatter and smooth the data using a 2D Kaiser window.  %
Finally, we represent the homodyne collected data as (sig)$^{1/2}$ to make the dynamics and line
widths comparable to heterodyne-collected techniques like absorbance and pump-probe spectra.  %
Throughout this work, zero signal on the color bar is set to agree with the average rather than the
minimum of noise.  %
Values below zero due to measurement uncertainty underflow the color bar and are plotted in
white.  %
This is especially evident in lots such as +120 fs in \autoref{fig:Czech08}, where there is no real
signal.  %
IPython \cite{PerezFernando2007a} and matplotlib \cite{HunterJohnD2007a} were important for data
processing and plotting in this work.

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/11}
  \caption[MoS\textsubscript{2} post processing.]{
    Visualization of data collection and processing.
    With the exception of (c), each subsequent pane represents an additional processing step on top
    of previous processing.
    The color bar of each image is separate.
    (a) Voltages read by the detector at teach color combination.
    The large vertical feature is $\omega_1$ scatter; the shape is indicative of the power curve of
    the OPA.
    MoS\textsubscript{2} response can be barely seen above this scatter.
    (b) Data after chopping and active background subtraction at the boxcar (100 shots).
    (c) The portion of chopped signal that is not material response.
    This portion is extracted by averaging several collections at very positive $\tau_{21}$ values,
    where no material response is present due to the short coherence times of MoS\textsubscript{2}
    electronic states.
    The largest feature is $\omega_2$ scatter.
    Cross-talk between digital-to-analog channels can also be seen as the negative portion that
    goes as $\omega_1$ intensity.
    (d) Signal after (c) is subtracted.
    (e) Smoothed data.
    (f) Amplitude level (square root) data.
    This spectrum corresponds to that at 0 delay in \autoref{fig:Czech03}.
    Note that the color bar's range is different than in \autoref{fig:Czech03}.}
  \label{fig:Czech11}
\end{figure}

\section{Results and discussion}  % ===============================================================

The few-layer MoS\textsubscript{2} thin film sample studied in this work was prepared on a
transparent fused silica substrate by a simple sufidation reaction of a Mo thin film using a
procedure modified from a recent report. \cite{LaskarMasihhurR2013a}  %
\autoref{fig:Czech02}a and b show the homogeneous deposition and surface smoothness of the sample
over the centimeter-sized fused silica substrate, respectively.  %
The Raman spectrum shows the $E_{2g}^1$ and $A_{1g}$ vibrational modes (\autoref{fig:Czech02}c)
that are characteristic of MoS\textsubscript{2}. \c ite{LiSongLin2012a}  %
The transmission electron micrograph (TEM) in \autoref{fig:Czech02}d shows the lattice fringes of
the film with an inset fast Fourier transform (FFT) of the TEM image indicative of the hexagonal
crystal structure of the film corresponding to the 0001 plane of MoS\textsubscript{2}.
\cite{LukowskiMarkA2013a}  %
The MoS\textsubscript{2} film thickness was determined to be 2.66 nm by atomic force microscopy and
corresponds to approximately four monolayers.  %
\autoref{fig:Czech02}3 shows the absorption and fluorescence spectrum of the film along with the A
and B excitonic line shapes that were extracted from the absorption spectrum. A representative
excitation pulse profile is also shown in red for comparison.  %

Extracting the exciton absorbance spectrum is complicated by the large ``rising background'' signal
from other MoS\textsubscript{2} bands.  %
With this in mind, we fit the second derivative absorption spectrum to a sum of two second
derivative Gaussians, as seen in \autoref{fig:CzechS3}.  %
Conceptually, this method can be thought of as maximizing the smoothness (as opposed to minimizing
the amplitude) of the remainder between the fit and the absorption spectrum.  %
The fit parameters can be found in the inset table in \autoref{fig:CzechS3}.  %
The Gaussians themselves and the remainder can be found in \autoref{fig:CzechS3}.  %

\begin{figure}
  \includegraphics[width=0.75\textwidth]{MX2/02}
  \caption[Few-layer MoS\textsubscript{2} thin film characterization.]{
    Characterization of the few-layer MoS\textsubscript{2} film studied in this work.
    Optical images of the
    MoS\textsubscript{2} thin film on fused silica substrate in (a) transmission and (b)
    reflection.
    (c) Raman spectrum of the $E_{2g}^1$ and $A_{1g}$ vibrational modes.
    (d) High-resolution TEM image and its corresponding FFT shown in the inset.
    (e) Absorption (blue), photoluminescence (green), Gaussian fits to the A and B excitons, along
    with the residules betwen the fits and absorbance (dotted), A and B exciton centers (dotted)
    and representative excitation pulse shape (red).}
  \label{fig:Czech02}
\end{figure}

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/S3}
  \caption[MoS\textsubscript{2} absorbance.]{Extraction of excitonic features from absorbance
    spectrum. (a) Second derivative spectra of absorbance (black) and fit second derivative
    spectrum (green). Gaussian fit parameters are shown in the inset table. (b) Absorption curve
    (black), Gaussian fits (blue and red), and remainder (black dotted).}
  \label{fig:CzechS3}
\end{figure}

The multiresonant CMDS experiment uses $\approx$70 fs excitation pulses created by two
independently tunable optical parametric amplifiers (OPAs).  %
Automated delay stages and neutral density filters set the excitation time delays over all values
of $\tau_{21}$ with $\tau_{22^\prime}=0$ and the pulse fluence to 90 $\mu$J/cm$^2$ (114
$\mu$J/cm$^2$) for the $\omega_1$ ($\omega_2$ and $\omega_{2^\prime}$) beam(s).  %
Each pulse was focused onto the sample using a distorted BOXCARS configuration.
\cite{EckbrethAlanC1978a}  %
The FWM signal was spatially isolated and detected with a monochromator that tracks the output
frequency so $\omega_m = \omega_1$.  %
In order to compare the FWM spectra with the absorption spectrum, the signal has been defined as
the square root of the measured FWM signal since FWM depends quadratically on the sample
concentration and path length.  %

The main set of data presented in this work is an $\omega_1\omega_2\tau_{21}$ ``movie'' with
$\tau_{22\prime}=0$.
\autoref{fig:Czech03} shows representative 2D frequency-frequency slices from this movie at
increasingly negative $\tau_{21}$ times.  %
Each 2D frequency spectrum contains side plots along both axes that compare the absorbance spectrum
(black) to the projection of the integrated signal onto the axis (blue).  %
Along $\omega_1$ (which for negative $\tau_{21}$ times acts as the ``probe'') we observe two peaks
corresponding to the A and B excitons.  %
In contrast, we see no well-defined excitonic peaks along the $\omega_2$ ``pump'' axis.  %
Instead, the signal amplitude increases toward bluer $\omega_2$ values.  %
The decrease in FWM above 2.05 eV is caused by a drop in the $\omega_2$ OPA power.

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/03}
  \caption[MoS\textsubscript{2} frequency-frequency slices.]{2D frequency-frequency spectra of the
    MoS\textsubscript{2} sample in the epi configuration. In all spectra $\tau_{22^\prime}=0$ fs,
    while $\tau_{21}$ is designated in the bottom-right corner of each spectral panel. The color
    bar defines the square root of the intensity normalized to the most intense feature in the
    series of spectra. The integration of the signal onto the $\hbar\omega_1=\hbar\omega_m$ and
    $\hbar\omega_2$ axes are represented ans the blue curves in the top and right side plots,
    respectively. The side plots also contain the absorbance spectrum (black line) to aid
    intepretation of the dynamics of the integrated 2D signals. The dashed lines mark the centers
    of the A and B excitons, as designated from the absorption spectrum.}
  \label{fig:Czech03}
\end{figure}


Figures \ref{fig:Czech04} and \ref{fig:Czech05} show representative 2D frequency-delay slices from
this movie, where the absicissa is the $\omega_1$ or $\omega_2$ frequency, respectively, the
ordinate is the $\tau_{21}$ delay time, and the solid bold lines represent five different
$\omega_2$ or $\omega_1$ frequencies.  %
The color bar is normalized to the brightest feature in each subplot.  %
This normalization allows comparison of the time dependence of the line shapes, positions, and
relative signal amplitudes along the $\omega_1$ or $\omega_2$ axis directly.  %

Each subplot in \autoref{fig:Czech04} is similar to published pump-probe, transient absorption,
and transient reflection experiments that have measured the electronc dynamics of the A and B
excitons. \cite{XiaoDi2012a, FangHui2014a, KumarNardeep2013a, NieZhaogang2014a, SunDezheng2014a,
  SimSangwan2013a, MakKinFai2012a, ThomallaMarkus2006a}  %
These previous experiments measure relaxation dynamics on the same $\approx$400-600 fs time scale
that is characteristic of \autoref{fig:Czech04}.  %

Our experiments also show how the spectral features change as a function of the $\omega_2$
excitation frequency.  %
The top to subplots of \autoref{fig:Czech04} reflect the changes in the AA and BA features, while
the bottom two subplots reflect the changes in the AB and BB features.  %
The figure highlights the changes in the relative amplitude of the A and B features as a function
of excitation frequency.  %
Both the line shapes and the dynamics of the spectral features are very similar.  %
\autoref{fig:Czech05} is an excitation spectrum that shows that the dynamics of the spectral
features do not depend strongly on the $\omega_1$ frequency.

\begin{figure}
  \includegraphics[width=0.75\textwidth]{MX2/04}
  \caption[MoS\textsubscript{2} $\omega_1$ Wigner progression.]{Mixed $\omega_1$---$\tau_{21}$
    time---frequency representations of the 3D data set at five ascending $\omega_2$ excitation
    frequencies (solid black lines) showing the impact of the $\omega_2$ excitation frequency on
    the $\omega_1$ spectral line shape as a function of time. The A and B exciton energies are
    marked as dashed lines within each spectrum.}
  \label{fig:Czech04}
\end{figure}

\begin{figure}
  \includegraphics[width=0.75\textwidth]{MX2/05}
  \caption[MoS\textsubscript{2} $\omega_2$ Wigner progression.]{Mixed $\omega_2$---$\tau_{21}$
    time---frequency representations of the 3D data set at five ascending $\omega_1$ probe
    frequencies (solid black lines) showing the impact of the $\omega_1$ excitation frequency on
    the $\omega_2$ spectral line shape as a function of time. The A and B exciton energies are
    marked as dashed lines within each spectrum.}
  \label{fig:Czech05}
\end{figure}


The spectral features in Figures \ref{fig:Czech03}, \ref{fig:Czech04} and \ref{fig:Czech05} depend
on the quantum mechanical interference effects caused by the different pathways.  %
features. \cite{PakoulevAndreiV2010a, PakoulevAndreiV2009a}  %
These pathways correspond to the time orderings labeled V and VI in \autoref{fig:Czech02}b. The
letters denote the density matrix elements, $\rho_{ij}$, where g representest the ground state and
e, e$^\prime$ represent any excitonic state.  %
Interaction with the temporally overlapped $\omega_2$ and $\omega_{2^\prime}$ pulses excites the ee
excited-state population and bleaches the ground-state population.  %
Subsequent interaction with $\omega_1$ creates the output coherences for the diagonal spectral
features when $\omega_1=\omega_2$ or the cross-peak features when $\omega_1\neq\omega_2$.  %
The stimulated emission (SE) and ground-state bleaching (GSB) pathways create the eg or e$^\prime$g
output coherences from the ee and gg populations, respectively, while the excited-state absorption
pathway creates the 2e,e or e$^\prime$+e,e biexcitonic output coherences.  %
There is also a population transfer pathway from the ee excited-state population to an
e$^\prime$e$^\prime$ population from which similar SE and ESA pathways occur.  %
Since the ESA pathways destructively interfere with the SE and GSB pathways, the output singal
depends on the differences between the pathways.  %
Factors that change the biexcitonic output coherences such as the transition moments, state filling
(Pauli blocking), frequency shifts, or dephasing rate changes will control the output signal.  %
State filling and ground-state depletion are important factors for MoS\textsubscript{2} since the
transitions excite specific electron and hole spin and valley states in individual layers.  %

The state-filling and ground-state bleaching effects on the diagonal and cross-peak features in
\autoref{fig:Czech03} depend on the spin and valley states in the output coherence.
\cite{XuXiaodong2014a}  %
The effects of these spins will disappear as the spin and valley states return to equilibrium.
\cite{ZengHualing2012a, ZhuBairen2014a, MaiCong2013a}
If we assume spin relaxation is negligible, the FWM transitions that create the diagonal features
involve either the A or B ESA transitions, so the resulting 2e,e state includes two spin-aligned
conduction band electrons and valence band holes.  %
Similarly, the cross-peak regions denoted by AB or BA in \autoref{fig:Czech01}d will have two
transitions involving an A exciton, so the initial e$^\prime$+e,e state will include antialigned
spins.  %
A quantitative treatment of the cancellation effects between the GSB, SE, and ESA pathways requries
knowledge of the transition moments and state degeneracies and is beyond the scope of this paper.
\cite{WongCathyY2011a}  %

The most important characteristic of the experimental spectra is the contrast between the absence
of well-resolved excitonic features that depend on $\omega_2$ in Figures \ref{fig:Czech03} and
\ref{fig:Czech05} and the well-defined excitonic features that depend on $\omega_1$.  %
It is also important to note that the projections of the signal amplitude onto the $\omega_2$ axis
in \autoref{fig:Czech03} match closely with the continuum features in the absorption spectrum and
that the line shapes of the features along the $\omega_2$ axis in \autoref{fig:Czech05} do not
change appreciably for different delay times or $\omega_1$ values. %
It shoudl be noted that the excitation pulse bandwidth (see \autoref{fig:Czech02}e) contributes to
the absence of well-resolved A and B excitonic features along $\omega_2$.  %
The similarity to the continuum states in the absorption spectrum and the absence of a strong
dependence on $\omega_1$ show that the continuum states observed at higher $\omega_2$ frequencies
participate directly in creating the final output coherences and that their increasing importance
reflects the increasing absorption strength of higher energy continuum states.  %
In contrast, the features dependent on $\omega_1$ in Figures \ref{fig:Czech03} and
\ref{fig:Czech04} match the line shapes of the A and B excitonic resonances.  %
Although the relative amplitudes of the spectral features depend on the $\omega_2$ frequency, they
do not depend strongly on the delay times.  %
These characteristics show that hot A and B excitonic states undergo rapid intraband population
relaxation over a $<$70 fs time scale set by excitation pulses to the A and B excitonic states
excited by $\omega_1$.

A central feature of Figures \ref{fig:Czech03} and \ref{fig:Czech04} is that the AB region is much
brighter than the BA region.  %
This difference is suprising because simple models predict cross-peaks of equal amplitude, as
depicted in \autoref{fig:Czech01}d.  %
The symmetry in simple models arises because the AB and BA cross-peaks involve the same four
transitions.  %
The symmetry between AB and BA may be broken by material processes such as population relaxation
and transfer, the output coherence dephasing rates, and the bleaching and state-filling effects of
the valence and conduction band states involved in the transitions.  %

We believe that ultrafast intraband population transfer breaks the symmetry of AB and BA
cross-peaks.  %
For the BA peak, the interactions with $\omega_2$ and $\omega_{2^\prime}$ generate only A excitons
that do not relax on $<$70 fs time scales.  %
For the AB peak, $\omega_2$ and $\omega_{2^\prime}$ generate two kinds of excitons: (1) B excitons
and (2) hot excitons in the A band.  %
Relaxation to A may occur by interband transitions of B excitons or intraband transitions of hot A
excitons.  %
Either will lead to the GSB, SE, and ESA of the ee $\rightarrow$ e$^\prime$e$^\prime$ population
transfer pathways.  %
This relaxation must occur on the time scale of our pulse-width since the cross-peak asymmetry is
observed even during temporal overlap.  %
We believe that intraband relaxation of hot A excitons is the main factor in breaking the symmetry
between AB and BA cross-peaks.  %
\autoref{fig:Czech04} shows that B $\rightarrow$ A interband relaxation occurs on a longer time
scale.  %
The B/A ratio is higher when $\omega_2$ is resonant with the B excitonic transition than when
$\omega_2$ is lower than the A exciton frequency (the top subplot).  %
If population transfer of holes from the B to A valence bands occurred during temporal overlap, the
B/A ratio would be independent of pump frequency at $\tau_21<0$.
\autoref{fig:Czech07} shows the delay transients at the different frequencies shown in the 2D
spectrum.  %
The colors of the dots on the 2D frequency-frequency spectrum match the colors of the
transients.  %
The transients were taken with a smaller step size and a longer time scale than the delay space
explored in the 3D data set.  %
The transients are quite similar.  %
Our data are consistent with both monomolecular biexponential and bimolecular kinetic models and
cannot discriminate between them.  %
We have fit the decay kinetics to a single exponential model with a time constant of 680 fs and an
offset that represents the long time decay.  %
The 680 fs decay is similar to previously published pump-probe and transient absorption
experiments. \cite{NieZhaogang2014a, SunDezheng2014a, DochertyCallumJ2014a}  %

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/07}
  \caption[MoS\textsubscript{2} transients.]{
    Transients taken at the different $\omega_1$ and $\omega_2$ frequencies indicated by the
    colored markers on the 2D spectrum.
    The dynamics are assigned to a 680 fs fast time constant (black solid line) and a slow time
    constant represented as an unchanging offset over this timescale (black dashed line).}
  \label{fig:Czech07}
\end{figure}

The spectral features change quantitatively for delay times near temporal overlap.  %
\autoref{fig:Czech08} shows a series of 2D spectra for both positive and negative $\tau_{21}$ delay
times with $\tau_{22^\prime}=0$.  %
Each spectrum is normalized to its brightest feature.  %
The spectra at positive $\tau_{21}$ delay times become rapidly weaker as the delay times become
more positive until the features vanish into the noise at +120 fs.  %
The spectra also develop more diagonal character as the delay time moves from negative to positive
values.  %
The AB cross-peak is also a strong feature in the spectrum at early times.  %

The pulse overlap region is complicated by the multiple Liouville pathways that must be
considered.  %
Additionally, interference between scattered light from the $\omega_1$ excitation beam and the
output signal becomes a larger factor as the FWM signal decreases.  %
Since $\tau_{22^\prime}=0$, the initial $\omega_1$ pulse creates an excited-state coherence, while
the subsequent $\omega_2$ and $\omega_{2^\prime}$ pulses create the output coherence.  %
The output signal is only important at short $\tau_{21}>0$ values because the initially excited
coherence dephases very rapidly.  %
When $\omega_1\neq\omega_2$, the first two interactions create an e$^\prime$e zero quantum
coherence that also dephases rapidly.  %
However, when $\omega_1=\omega_2$, the first two interactions create an ee, gg population
difference that relaxes on linger time scales.  %
The resulting signal will therefore appear as the diagonal feature in \autoref{fig:Czech08}
(\textit{e.g.}, see the +40 fs 2D spectrum).  %
In addition to the diagonal feature in \autoref{fig:Czech08}, there is also a vertical feature when
$\omega_1$ is resonant with the A excitonic transition as well as the AB cross-peak.  %
Although these pathways are depressed when $\tau_{21}>0$, there is sufficient temporal overlap
between the $\omega_2$, $\omega_{2^\prime}$, and $\omega_1$ pulses to make their contribution
comparable to V and VI.  %
More positie values of $\tau_{21}$ emphasize the I, III pathways over the V, VI pathways,
accounting for the increasing percentage of diagonal character at increasingly positive delays.  %

\begin{figure}
  \includegraphics[width=\textwidth]{MX2/08}
  \caption[MoS\textsubscript{2} frequency-frequency slices near pulse overlap.]{2D
    frequency-frquency spectra near zero $\tau_{21}$ delay times. The signal amplitude is
    normalized to the brightest features in each spectrum.}
  \label{fig:Czech08}
\end{figure}

\section{Conclusions}  % ==========================================================================

This paper presents the first coherent multidimensional spectroscopy of MoS\textsubscript{2} thin
films.  %
CMDS methods are related to the earlier ultrafast pump-probe and transient absorption methods since
they all share bleaching, stimulated emission, Pauli blocking, and excited-state absorption
pathways, but they differ in how these pathways define the spectra.  %
In addition, CMDS methods have many additional pathways that become important when the coherence
dephasing times are longer than the excitation pulse widths.  %
In this work, the dephasing times are short so the pathways are identical to transient
absorption.  %
This work reports the first frequency-frequency-delay spectra of MX\textsubscript{2} samples.  %
These spectra are complementary to previous work because they allow a direct comparison between the
initially excited excitonic states and the states creating the final output coherence.  %
The spectra show that the same hot A and B exciton continuum states that are observed in the
absorption spectrum also dominate the CMDS excitation spectra.  %
They also show that rapid, $<70$ fs intraband relaxation occurs to create the band-edge A and B
excitonic features observed in the CMSD spectrum.  %
The relative intensity of the diagonal peak features depends on the relative absorption strength of
the A and B excitons.  %
The relative intensity of cross-peak features in the 2D spectra depends on the excitation
frequency.  %
Excitation at or above the B exciton feature creates strong cross-peaks associated with hot A and B
excitons that undergo ultrafast intraband population transfer.  %
Excitation below the B excitonic feature creates a weak cross-peak indicating A-induced B-state
bleaching but at a lower signal level corresponding to the lower optical density at this energy.  %
Population relaxation occurs over $\approx$680 fs, either by transfer to traps or by bimolecular
charge recombination.  %

These experiments provide the understanding of MoS\textsubscript{2} coherent multidimensional
spectra that will form the foundation required to measure the dynamical processes occurring in more
complex MoS\textsubscript{2} and other TMDC heterostructures with quantum-state resolution.  %
The frequency domain based multiresonant CMDS methods described in this paper will play a central
role in these measurements.  %
They use longer, independently tunable pulses that provide state-selective excitation over a wide
spectral range without the requirement for interferometric stability.  %

\section{Errata}  % ===============================================================================

The following is an errata for the of this chapter, published in November 2015.
\cite{CzechKyleJonathan2015a}  %

\begin{ditemize}
  \item Reference 13 is identical to reference 9.
  \item In the last paragraph of the introduction the sentence ``The experimental spectra differ
    from the simple 2D spectrum shown in Figure 1d and those of earlier CMDS experiments with model
    systems'' appears. This sentence cites references 6 through 10. Instead, it should cite
    references 15 through 20.
  \item In the last paragraph begining on page 12148, the text ``Automated delay stages and neutral
    density filters set the excitation time delays over all values of $\tau_{21}$ with
    $\tau_{22}=0$''    appears. For the second $\tau$, the subscript should read $22^\prime$, not
    $22$.
  \item Caption of Figure 5 reads, in part: "showing the impact of the $\omega_1$ excitation
    frequency on the $\omega_1$ spectral line shape". This should instead read "showing the impact
    of the $\omega_1$ excitation frequency on the $\omega_2$ spectral line shape". The subscript on
    the last $\omega$ should be a 2 and not a 1.
  \item Figure 6 e$^\prime$+e,e$^\prime$ should read e$^\prime$+e,e and vice versa.
\end{ditemize}