# Time-resolved Broadband Spectroscopy of Coherent Lattice Dynamics in a
Rashba Material

Peter Fischer  Julian Bär  Moritz Cimander  Volker Wiechert  Alfred
Leitenstorfer  Davide Bossini  Department of Physics and Center for Applied
Photonics, University of Konstanz, D-78457 Konstanz, Germany

###### Abstract

Spintronics aims at establishing an information technology relying on the
spins of electrons. Key to this vision is the spin-to-charge conversion (and
vice versa) enabled by spin-orbit coupling, which represents the intimate bond
between spins and charges. Here, we set the foundations to manipulate the
spin-orbit coupling in Rashba systems in a coherent way at THz frequencies
without parasitic electronic dissipation. Our approach relies on the
generation of coherent lattice dynamics in BiTeI. Using ultrashort laser
pulses ranging from the visible to the mid-infrared spectral range we
establish quantitatively the conditions for an effective generation of
coherent THz phonons, with the required symmetry for a manipulation of the
Rashba coupling. Furthermore, we conclude an open discussion on the
possibility of photoinducing of a spin current in BiTeI, proving that this is
not the case. We envision a direct application of our pumping scheme in a
photoemission experiment, which can directly disclose the THz dynamics of the
Rashba splitting.

††preprint: APS/123-QED

##  I Introduction

The central role of spin-orbit coupling in spintronics arguably makes it the
most relevant physical interaction and process in this research field. This
includes that in extended systems spin-orbit interaction can lift the spin-
degeneracy of the electron dispersion [1]. This phenomenon is the Rashba
effect and it stems from both the atomic spin-orbit coupling and an effective
electric field $E_{z}$ that breaks the crystal inversion symmetry. In the one-
particle Hamiltonian describing the electron motion within the crystal
potential, the Rashba effect is taken into account by the interaction term

$H_{\textnormal{R}}=-\alpha_{\textnormal{R}}(\bm{k}\times\bm{\sigma})\cdot\bm{\hat{z}},$
(1)

where $\alpha_{\textnormal{R}}$, $\bm{k}$, $\bm{\sigma}$ and $\bm{\hat{z}}$
represent the Rashba constant, the electron wave vector, the Pauli matrices
and the unit vector along the direction of the inversion symmetry breaking,
respectively. The Rashba effect has been widely employed in the investigation
of magneto-transport effects [2, 3]. Recently, even photoinduced modifications
of the entire Rashba energy in Eq. (1) (i.e. dynamics of $\bm{k}$,
$\bm{\sigma}$ and $\alpha_{\textnormal{R}}$) have been reported [4]. This
impressive amount of experimental and theoretical work has so far addressed
heterostructures, where the interface inherently breaks the inversion
symmetry. However, heterostructures require non-trivial engineering and
fabrication technology, as the quality of the interface determines the
functionality of a device [2, 3]. Moreover, from the experimental point of
view, heterostructures can be studied almost exclusively by surface-sensitive
techniques. However, the Rashba effect was originally predicted for bulk
materials. As a matter of fact, bulk semiconductors exhibiting a huge Rashba
effect have been discovered, such as BiTeI [5]. The research field of
spintronics would benefit greatly from a coherent control of the Rashba
coupling in bulk systems on the ultrashort timescale. This scenario would
allow arbitrary suppression and enhancement of the Rashba coupling and thus of
the spin-to-charge conversion capability. In fact, it would directly open the
doors to bulk-compatible spintronic devices able to perform tasks at THz
operational frequencies and thus faster than ever before. Furthermore, a fully
coherent approach would also not be restricted by parasitic thermal energy
dissipation. In relation to the expression of the Rashba Hamiltonian (Eq.
(1)), this concept implies a pure modification of $\alpha_{\textnormal{R}}$
whereby the electronic and magnetic degrees of freedom remain unperturbed.
Hence, an urgent fundamental question arises: how can a fully coherent control
of the Rashba coupling on the ultrashort timescale be realised?

We propose that an excitation of coherent longitudinal optical phonons with an
ionic displacement along the symmetry-breaking axis, can modulate the
effective electric field generating the Rashba coupling. As a first step
towards this unprecedented concept, here we experimentally and quantitatively
establish the conditions that allow efficient excitation of coherent lattice
dynamics with adequate symmetry in BiTeI. In addition , we investigate the
optical generation of a spin polarisation in BiTeI, as this concept is also
highly relevant for spintronics. The literature reports contrasting claims
about this issue [6, 7, 8]. Our experimental results definitely set the
debate, demonstrating that no spin polarisation can be photoinduced in BiTeI
on viable timescale.

##  II Methods and Materials

Figure 1:  (a) Atomic displacement of the 2.7-THz A1g phonon mode in BiTeI.
(b) Sketch of the pump-probe setup. The pump pulse (pink) excites coherent
lattice dynamics in BiTeI. The probe pulse (cyan) is reflected off the surface
of the sample with a variable time delay, allowing to track the dynamics of
the reflection coefficient $\Delta R/R_{0}$.

For our experiments we choose the nonmagnetic polar semiconductor BiTeI. Our
specimen was grown by a modified Bridgmen method using a rotating heat field.
Its layered structure is built from covalently bound trilayers stacked by van
der Waals forces. The opposite partial charges of iodine and tellurium ion
cores break the inversion symmetry along the stacking direction, generating
the electric field $E_{z}$ responsible for the Rashba effect (see Fig. 1).
Hence, not only the surface states exhibit a giant Rashba splitting but also
the bulk states [5]. The excitation of a specific A1g mode, i.e., a coherent
displacement of iodine and tellurium ion cores in opposite directions but
parallel to the stacking, as shown in Fig. 1 (a), is expected to modulate
$E_{z}$ and thus the Rashba splitting for the bulk states. The eigenfrequency
of this specific phonon has already been identified via spontaneous Raman
scattering [9] and it amounts to 2.7 THz. Since the A1g mode in BiTeI
classifies as longitudinal optical phonon, a direct excitation with THz
radiation is forbidden by selection rules. However, Raman-active modes can be
coherently excited via laser pulses, provided that their duration is shorter
than the period of the phonon [10]. Aiming at establishing the most favourable
conditions to realise such an excitation, we set out to perform time-resolved
pump-probe reflection spectroscopy on BiTeI with a broad spectral range of
pump wavelengths. We thus employ two different light sources, both amplified
laser systems emitting fs laser pulses at kHz repetition rates (details in the
Supplementary Materials). Having access to both of these light sources grants
us a broadband tunability of the central photon energy of the excitation beam,
ranging from the visible to the mid-infrared. Thus, we are able to employ
photon energies both above and below the fundamental band gap of BiTeI. In
both setups, the central wavelength of the probe beam is 1.20 µm. However, the
duration of the probe pulses differs significantly between the two systems and
is 20 fs for the mid-infrared measurements, while it is 90 fs for the other
measurements. A quantitative analysis of our data must take into account this
pivotal difference (see Supplementary Materials).

The pump-probe technique employed to investigate the coherent lattice dynamics
in BiTeI is sketched in Fig. 1 (b). To track the dynamics of the normalised
reflection coefficient $\Delta R/R_{0}$, the probe pulses are reflected off
the sample surface with a variable time delay with respect to the arrival of
the pump pulses. The intensity of the reflected probe beam is measured with a
photo diode in a lock-in scheme. The dynamics of $\Delta R/R_{0}$ were
measured both at room temperature and 10 K. We have also performed time-
resolved magneto-optical experiments to investigate a possible photoinduction
of a spin polarization in BiTeI. These methods are known to be sensitive to
such excitations in semiconductors [11]. The rotation of the polarisation is
measured by installing a quarter-wave plate and a Wollaston prism in the
optical path of the probe pulse after the reflection.

##  III Time-resolved broadband optical spectroscopy

Figure 2:  (a) Transmissivity measurement performed with an FTIR on a µm-thin
film exfoliated from a bulk specimen of BiTeI. The inset shows the band
structure of BiTeI obtained by first-principles calculations [12]. The
transition $\beta$ (red dashed line at 0.37 eV) represents the highest energy
for which a transition between the Rashba-split bands is possible and $\alpha$
(green) the lowest. The transition $\gamma$ (blue dashed line at 0.52 eV)
represents the lowest energy for which electrons can be excited from the
valence band to a free state above the Fermi energy $E_{\textnormal{F}}$
(dashed black line). In the narrow window between $\beta$ and $\gamma$ lies a
range of increased transmission. (b) Spectral absorption coefficient of BiTeI
digitized from [13]. The triangular markers indicate the values of the
wavelength of the pump beam employed in our experiment. Figure 3: (a) Dynamics
of the normalised reflection coefficient of BiTeI photoinduced by 0.52 µm
laser pulses. The fluence was set to 0.77 mJ/cm2. (b) Coherent component of
the transient reflection coefficient isolated from the incoherent
background.(c) Spectrum of the time trace shown in panel (b), obtained via
Fourier transform. The solid black line represents a fit with a Lorentzian
function.

We first measure the transmissivity of a thin film of BiTeI exfoliated from
the bulk sample. As shown in Fig. 2 (a) the spectrum is dominated by two
transitions. The transmissivity is smaller than 5$\%$ in the energy range
below $\beta=0.37$ eV (red) from which it sharply increases. However, near
$\gamma=0.52$ eV (blue) a sharp drop of the transmissivity is observed. The
spectrum in Fig. 2 is consistent with the data reported in the literature
[14]. Although the highest valence and the lowest conduction bands in BiTeI
are all symmetrically of the same character, optical transition are
nevertheless allowed due to the presence of Rashba effect [12] even between
the Rashba-split conduction bands. In fact, BiTeI is a degenerate n-type
semiconductor, which implies that the Fermi energy $E_{\textnormal{F}}$ is
above the conduction-band minimum [5]. Thus, $\beta$ represents the highest
energy at which such a transition is possible, while $\gamma$ is the lowest
energy at which electrons can be promoted from the highest valence band to the
lowest conduction band. The latter is therefore the effective band-gap energy
shifted from the fundamental band gap by $\approx 0.17$ eV due to the Moss-
Burstein effect [15]. Note that the energy corresponding to $\gamma$ is
strongly dependent on the density of carriers in the conduction band, which
determines $E_{\textnormal{F}}$. Comparing to the literature [14] we estimate
the intrinsic conduction-band carrier density in our sample to be
$n_{\textnormal{i}}=6\cdot 10^{19}$ cm-3. On the other hand, due to the linear
dispersion of the Rashba-split conduction bands for $|k|>0.1$ Å$\beta$ is far
less sensitive to $E_{\textnormal{F}}$ and should only change significantly if
$E_{\textnormal{F}}$ is lowered near the conduction band minimum.

Fig. 2 (b) shows the absorption coefficient $\alpha_{\textnormal{abs}}$ of
BiTeI digitized from data in the literature [13]. We note that the intrinsic
carrier density corresponding to the data shown in Fig. 2 (b) was estimated to
be 6 times smaller than in our sample [13]. However, this difference affects
only the exact value of the energy corresponding to $\gamma$. The values of
the wavelength chosen for the pump beam in our experiments are marked in Fig.
2 (b). For the interband transitions (0.52 µm to 2.08 µm)
$\alpha_{\textnormal{abs}}$ is of the order of $10^{5}$ cm-1. These
considerations about the absorption coefficient of BiTeI are relevant to our
objective of impulsive excitation of coherent lattice dynamics. In particular,
this process can be successfully induced in two different light-matter
interaction regimes. If dissipation is negligible (i.e. negligible absorption
of the pump beam), impulsive stimulated Raman scattering (ISRS) of elementary
excitations can be realised [10, 16]. On the other hand, the absorption of
light by electrons can also lead to an impulsive force triggering Raman-active
phonons. In this case a displacive excitation of coherent phonons (DECP) is
achieved [10]. We note that the ISRS scenario seems favourable in view of our
goal, as it allows to exclusively excite the quasiparticles of interest,
without parasitic electronic dissipation. As these two mechanisms were
established in dielectrics and semiconductors with a much more simple band
structure than BiTeI, a spectrally broad tunability of the pump-photon energy
is key for our purposes. In particular, it can be expected that the interband
transitions in BiTeI may allow a DECP generation [17, 10] of the A1g lattice
mode. However, the role of the incoherent electronic background must be
quantitatively assessed by an experimental procedure, to realise the
conditions under which such background does not hinder the coherent lattice
dynamics [18]. The situation for the intraband transitions is less obvious, as
it has not yet been explored. The intraband absorption coefficient is an order
of magnitude smaller than in the case of interband processes (33 times lower
at 7.5 µm than at $0.52$ µm). Importantly, it is not straightforward to
predict how the vastly different natures of intraband and interband processes
affect the photoinduced coherent lattice dynamics.

Fig. 3 (a) shows the dynamics of the normalised reflection coefficient $\Delta
R/R_{0}$ photoinduced by laser pulses with a wavelength of 0.52 µm. The pulse
duration is 198 fs and the fluence is set to $\Phi=0.77$ mJ/cm2. The signal
consists of two components, (i) coherent oscillations with a frequency of 2.7
THz superimposed to (ii) an exponentially decaying background with an
amplitude of $\approx 15\cdot 10^{-3}$. We ascribe the incoherent component of
the observed signal to the relaxation dynamics and thermalisation processes of
the photoinduced electrons reported in the literature [7]. As we are
interested in the coherent response of the lattice, we subtract the incoherent
background (details of data analysis in Supplementary Materials) and analyse
the oscillatory pattern depicted in Fig. 3 (b). The peak-to-peak amplitude of
the oscillations amounts to $\approx 1\cdot 10^{-3}$ and decays in
approximately 7 ps. A Fourier transform of the time trace in Fig. 3 (b)
reveals the spectrum of the oscillations displayed in Fig. 3 (c), which we fit
with a Lorentzian function (see Supplementary Materials). As a result we
estimate the center frequency to be $f_{0}=2.76$ THz and the lifetime to be
$\tau=1.82$ ps. These values are consistent with the particular A1g mode, in
which the atomic displacement has the correct symmetry to modulate $E_{z}$, as
previously determined by spontaneous Raman scattering [9].

We investigate the dependence of the electronic background and the oscillation
amplitude on the pump wavelength. We perform measurements of $\Delta R/R_{0}$
employing the pump wavelengths indicated in Fig. 2 (b). For each pump
wavelength we analyze the dependence of electronic-background and oscillation
amplitude, $A_{\textnormal{e}}$ and $A_{\textnormal{osci}}$, on the excitation
fluence $\Phi$. The results are reported in Fig. 4 (a) and (b), in which
intraband and interband processes are color coded. We calculate the values of
the absorbed fluence for every pump wavelength, relying on the literature [13]
(see Supplementary Materials).

First, we focus on the interband excitation, highlighted in the figure by the
blue background. For each pump wavelength we find a linear dependence of
$A_{\textnormal{e}}$ on $\Phi$, indicating a strictly linear excitation of the
electronic system. As a general observation, Fig. 4 (a) shows that the slope
of the linear fluence dependence of $A_{\textnormal{e}}$ decreases as the
excitation wavelength increases. This observation is consistent with the
spectral absorption coefficient shown in Fig. 2, with the longest wavelength
employed in our experiment (2.08 µm) approaching the transparency range of
BiTeI ($\gamma$ transition, Fig. 2). The data set obtained for 0.75 µm
deviates from this trend. We observe that the duration of the pump pulses for
this wavelength is 4 times shorter than in all the other cases.

In Fig. 4 (b) the dependence of $A_{\textnormal{osci}}$ on $\Phi$ is depicted,
in which $A_{\textnormal{osci}}$ is obtained as in Fig. 3 (c) for each
measurement. We find linear fluence dependencies of $A_{\textnormal{osci}}$
with different slopes for different excitation wavelengths. These results are
fully consistent with the DECP mechanism [10]. Also in the case of
$A_{\textnormal{osci}}$, the data obtained by pumping BiTeI with 0.75 µm laser
pulses reveal a stronger signal. The already mentioned shorter pulse duration
for this excitation wavelength plays a major role in the discussion of
$A_{\textnormal{osci}}$. In fact both ISRS and DECP are impulsive processes
and thus tend to be more efficient as the excitation pulses approach the ideal
temporal profile of a delta function [10]. Additional to this consideration,
the spectrum of the Raman cross section of the A1g mode directly determines
the spectral dependence of $A_{\textnormal{osci}}$. Although we do not observe
sharp resonances in Fig. 4 (b), we ascribe the weakest response observed at
1.04 µm to a minimum of the cross section.

Figure 4:  Analysis of the dependence of the amplitudes of electronic
background $A_{\textnormal{e}}$ and 2.7-THz oscillations
$A_{\textnormal{osci}}$ on the fluence $\Phi$ in the $\Delta R/R_{0}$
measurement signal. Pump wavelengths (pulse durations) are indicated in the
top panel. For each measurement series, $A_{\textnormal{e}}$ (a) and
$A_{\textnormal{osci}}$ (b) show a linear dependence on the fluence. The
dotted lines serve as guides to the eye.

Next we turn to the light-matter interaction based on intraband processes
highlighted in Fig. 4 by the red background. The laser source for the mid-
infrared excitation (central wavelength 7.5 µm) is completely self-built and
has recently been characterised [19]. For the dependence of
$A_{\textnormal{e}}$ on $\Phi$ we find a linear dependence. Considering the
high value of the fluence employed (up to 10 mJ/cm2), it is rather surprising
that no traces of saturation are observed. In fact, the absorbed fluence
corresponds to an estimated excited carrier density of the order of $\approx
10^{21}$ cm-3 (see Supplementary Materials), which is two orders of magnitudes
larger than the estimated conduction-band charge-carrier density. However, in
the literature [6, 7] it has been reported that photoexcited excited electrons
scatter on a time scale shorter than 80 fs. As the pulse duration of the
7.5-µm pump pulses is 180 fs, the conduction-band electrons are rapidly
redistributed in momentum space within the duration of the optical excitation,
thus preventing the absorption from saturating. In addition, we note that at
the lowest fluence of 3.7 mJ/cm2 the electronic background possesses an
amplitude of $5\cdot 10^{-3}$. Although this value is much smaller than for
interband excitation at comparable fluences, it is not negligible. It follows
that also in this case significant electronic dissipation takes place,
implying that the light-matter interaction regime falls in the DECP
description. We thus conclude that the pure ISRS scenario cannot be realised
in BiTeI, in view of the photoinduced dynamics of the electrons already
present in the conduction bands (quantified by $A_{\textnormal{e}}$).

It is now crucial to compare the coherent lattice dynamics observed pumping
BiTeI in the mid-infrared spectral range, with the results obtained for
interband transitions. Also in the case of intraband processes
$A_{\textnormal{osci}}$ exhibits a linear dependence on $\Phi$, consistent
with the DECP picture [17]. The increase of $A_{\textnormal{osci}}$ with
$\Phi$ is significantly flatter as in the 2.08-µm measurement. We observe that
in the time traces of $\Delta R/R_{0}$ obtained pumping with 7.5 µm laser
pulses, a second eigenmode is visible. The Fourier analysis of this additional
feature allows to estimate its frequency to be 4.5 THz (see Supplementary
Materials). The value matches another A1g lattice mode [9], whose atomic
displacement cannot modulate $E_{z}$. The frequency of this mode demands
sub-100-fs resolution for its successful detection. This condition is
fulfilled by our mid-infrared optical setup while it is not by the visible and
near-infrared system. Comparing all data sets presented in Fig. 4 (a) and (b)
we conclude that $\Delta R/R_{0}$ shows a coherent 2.7-THz A1g mode induced
via DECP, even for pump photon energies smaller than the band gap, with
$A_{\textnormal{e}}\gg A_{\textnormal{osci}}$. Although we have already ruled
out the possibility of pure ISRS excitation, our broadband time-resolved
spectroscopy and its quantitative analysis suggest that it is far more
favorable and efficient to excite coherent phonons in BiTeI in the visible
spectral region. This statement challenges the commonly accepted wisdom, based
on the idea to minimize the optical absorption as much as possible to induce
only the elementary excitations of interest, while leaving unperturbed the
other degrees of freedom. Our conclusion, while counter-intuitive, is based on
the quantitative analysis of our results: the same oscillatory amplitude
$A_{\textnormal{osci}}$ induced by mid-infrared laser pulses can also be
observed employing visible pump pulses (0.52 µm) at an 3 times lower
excitation fluence. Our experiments provide this conclusion, despite the
striking observation that in BiTeI $\alpha_{\textnormal{abs}}$ for 0.52 µm is
more than 33 times higher than for 7.5 µm (Fig.2 (b)).

Throughout all of our datasets the frequency of the A1g phonon is constant.
The lifetime shows minor deviations for $\Phi>1.26$ mJ/cm2, ascribed to
electron-phonon scattering (see Supplementary Materials). The photoinduced
coherent lattice dynamics in BiTeI lies thus completely in a linear regime. In
addition, we have also performed measurements for visible and near-infrared
pump pulses at a sample temperature of 10 K. The analysis of these data,
presented in the Supplementary Materials does not reveal significant
differences compared to the room-temperature measurements.

Figure 5:  Probe-pulse polarization rotation upon the excitation with left
($\sigma_{-}$, red) and right ($\sigma_{+}$, red) circularly polarized light
as well as linearly ($\sigma_{0}$, black) polarized light. No asymmetry is
observed that would indicate a photoinduced spin polarization.

At last we investigate the possibility of generating a photoinduced spin
polarization by exciting BiTeI with circularly polarized light. In the
literature, this concept is controversial. In fact, on the one hand time-
resolved photo-emission spectroscopy (tr-ARPES) investigations [7, 6] have
shown that no spin polarization could be detected. On the other hand, optical
spectroscopy has revealed traces of a spin polarisation on a 100 fs timescale
[8]. As different methods report different results, we set out to settle this
debate by measuring time-resolved magneto-optical effects, which are sensitive
to the spin polarisation [11]. The measurements of the rotation of the probe
polarisation $\Delta\Theta$ induced by different polarisation states of the
pump beam are shown in Fig. 5. We compare the excitation of left
($\sigma_{-}$, red) and right ($\sigma_{+}$, red) circularly polarized light
with linearly ($\sigma_{0}$, black) polarized light and focus in the 100-fs
time scale. No polarization-dependent excitation occurs. We observe that the
report of a photoinduced spin polarisation relied on purely optical
spectroscopy [8]. However, because of the selection rules of optical
transitions this method is not sensitive to detect a net accumulation of
spins.

##  IV Conclusion

We have successfully demonstrated the coherent excitation of the 2.7-THz A1g
mode in BiTeI in a wide spectral range reaching from the visible to the mid-
infrared. The symmetry of this particular mode, i.e., the displacement of
iodine and tellurium atoms parallel to the stacking but in different
directions, has the potential to manipulate the Rashba spin splitting in the
band structure. Thus, the main, albeit counter-intuitive, conclusion of our
work is that excitation in the mid-infrared with photon energies smaller than
the band gap cannot compete with excitation in the visible spectral range,
although it is associated with the lowest values of the absorption
coefficient. In addition, working with visible laser pulses is both
technically much simpler and economically more favorable, since, in contrast
to the still demanding mid-infrared range, standard components can be used. To
observe the coherent modulation of the Rashba splitting in BiTeI, we propose a
combination of our pumping scheme with tr-ARPES. The electronic background
appearing in our data would not hinder the tr-ARPES detection of the dynamics
of the Rashba splitting, owing to the momentum-selectivity of this technique.

## References

  * Rashba [2015] E. Rashba, Semiconductors with a loop of extrema, Journal of Electron Spectroscopy and Related Phenomena 201, 4–5 (2015).
  * Sinova  _et al._ [2015] J. Sinova, S. O. Valenzuela, J. Wunderlich, C. H. Back, and T. Jungwirth, Spin hall effects, Rev. Mod. Phys. 87, 1213 (2015).
  * Manchon  _et al._ [2019] A. Manchon, J. Železný, I. M. Miron, T. Jungwirth, J. Sinova, A. Thiaville, K. Garello, and P. Gambardella, Current-induced spin-orbit torques in ferromagnetic and antiferromagnetic systems, Rev. Mod. Phys. 91, 035004 (2019).
  * Michiardi  _et al._ [2022] M. Michiardi, F. Boschini, H.-H. Kung, M. X. Na, S. K. Y. Dufresne, A. Currie, G. Levy, S. Zhdanovich, A. K. Mills, D. J. Jones, J. L. Mi, B. B. Iversen, P. Hofmann, and A. Damascelli, Optical manipulation of rashba-split 2-dimensional electron gas, Nature Communications 13, 10.1038/s41467-022-30742-5 (2022).
  * Ishizaka  _et al._ [2011] K. Ishizaka, M. S. Bahramy, H. Murakawa, M. Sakano, T. Shimojima, T. Sonobe, K. Koizumi, S. Shin, H. Miyahara, A. Kimura, K. Miyamoto, T. Okuda, H. Namatame, M. Taniguchi, R. Arita, N. Nagaosa, K. Kobayashi, Y. Murakami, R. Kumai, Y. Kaneko, Y. Onose, and Y. Tokura, Giant rashba-type spin splitting in bulk bitei, Nature Materials 10, 521 (2011).
  * Mauchain  _et al._ [2013] J. Mauchain, Y. Ohtsubo, M. Hajlaoui, E. Papalazarou, M. Marsi, A. Taleb-Ibrahimi, J. Faure, K. A. Kokh, O. E. Tereshchenko, S. V. Eremeev, E. V. Chulkov, and L. Perfetti, Circular dichroism and superdiffusive transport at the surface of bitei, Phys. Rev. Lett. 111, 126603 (2013).
  * Ketterl  _et al._ [2021] A. S. Ketterl, B. Andres, M. Polverigiani, V. Voroshnin, C. Gahl, K. A. Kokh, O. E. Tereshchenko, E. V. Chulkov, A. Shikin, and M. Weinelt, Effect of rashba splitting on ultrafast carrier dynamics in bitei, Phys. Rev. B 103, 085406 (2021).
  * Deng  _et al._ [2022] H. Deng, C. Zhang, W. Liang, X.-X. Zhang, and S.-N. Luo, Hot carrier dynamics of bitei with large rashba spin splitting, RSC Adv. 12, 16479 (2022).
  * Sklyadneva  _et al._ [2012] I. Y. Sklyadneva, R. Heid, K.-P. Bohnen, V. Chis, V. A. Volodin, K. A. Kokh, O. E. Tereshchenko, P. M. Echenique, and E. V. Chulkov, Lattice dynamics of bismuth tellurohalides, Phys. Rev. B 86, 094302 (2012).
  * Merlin [1997] R. Merlin, Generating coherent thz phonons with light pulses, Solid State Communications 102, 207 (1997), highlights in Condensed Matter Physics and Materials Science.
  * Gupta  _et al._ [2001] J. A. Gupta, R. Knobel, N. Samarth, and D. D. Awschalom, Ultrafast manipulation of electron spin coherence, Science 292, 2458 (2001), https://www.science.org/doi/pdf/10.1126/science.1061169 .
  * Bahramy  _et al._ [2011] M. S. Bahramy, R. Arita, and N. Nagaosa, Origin of giant bulk rashba splitting: Application to bitei, Phys. Rev. B 84, 041202 (2011).
  * Makhnev  _et al._ [2014] A. A. Makhnev, L. V. Nomerovannaya, T. V. Kuznetsova, O. E. Tereshchenko, and K. A. Kokh, Optical properties of bitei semiconductor with a strong rashba spin-orbit interaction, Optics and Spectroscopy 117, 764 (2014).
  * Lee  _et al._ [2011] J. S. Lee, G. A. H. Schober, M. S. Bahramy, H. Murakawa, Y. Onose, R. Arita, N. Nagaosa, and Y. Tokura, Optical response of relativistic electrons in the polar bitei semiconductor, Phys. Rev. Lett. 107, 117401 (2011).
  * Burstein [1954] E. Burstein, Anomalous optical absorption limit in insb, Phys. Rev. 93, 632 (1954).
  * Bossini  _et al._ [2016] D. Bossini, V. I. Belotelov, A. K. Zvezdin, A. N. Kalish, and A. V. Kimel, Magnetoplasmonics and femtosecond optomagnetism at the nanoscale, ACS Photonics 3, 1385 (2016), https://doi.org/10.1021/acsphotonics.6b00107 .
  * Zeiger  _et al._ [1992] H. J. Zeiger, J. Vidal, T. K. Cheng, E. P. Ippen, G. Dresselhaus, and M. S. Dresselhaus, Theory for displacive excitation of coherent phonons, Phys. Rev. B 45, 768 (1992).
  * Cheng  _et al._ [1990] T. K. Cheng, S. D. Brorson, A. S. Kazeroonian, J. S. Moodera, G. Dresselhaus, M. S. Dresselhaus, and E. P. Ippen, Impulsive excitation of coherent phonons observed in reflection in bismuth and antimony, Applied Physics Letters 57, 1004–1006 (1990).
  * Schoenfeld  _et al._ [2023] C. Schoenfeld, L. Feuerer, A.-C. Heinrich, A. Leitenstorfer, and D. Bossini, Nonlinear generation, compression and spatio-temporal analysis of gv/cm-class femtosecond mid-infrared transients (2023), arXiv:2310.06619 [physics.optics] .


††thanks: Current address: Beijing Jiaotong University, Beijing, China

# Non-thermal evolution of dense plasmas driven by intense x-ray fields

Shenyuan Ren<EMAIL_ADDRESS>Department of Physics, Clarendon
Laboratory, University of Oxford, Parks Road, Oxford OX1 3PU, UK  Yuanfeng Shi
Department of Physics, Clarendon Laboratory, University of Oxford, Parks Road,
Oxford OX1 3PU, UK  Quincy Y. van den Berg  Department of Physics, Clarendon
Laboratory, University of Oxford, Parks Road, Oxford OX1 3PU, UK  Muhammad
Firmansyah  Department of Physics, Clarendon Laboratory, University of Oxford,
Parks Road, Oxford OX1 3PU, UK  Hyun-Kyung Chung  Korea Institute of Fusion
Energy, Daejeon, 34133, Republic of Korea  Elisa V. Fernandez-Tello  Instituto
de Fusión Nuclear, Universidad Politécnica de Madrid, José Gutiérrez Abascal
2, 28006 Madrid, Spain  Pedro Velarde  Instituto de Fusión Nuclear,
Universidad Politécnica de Madrid, José Gutiérrez Abascal 2, 28006 Madrid,
Spain  Justin S. Wark  Department of Physics, Clarendon Laboratory, University
of Oxford, Parks Road, Oxford OX1 3PU, UK  Sam M. Vinko
<EMAIL_ADDRESS>Department of Physics, Clarendon Laboratory,
University of Oxford, Parks Road, Oxford OX1 3PU, UK  Central Laser Facility,
STFC Rutherford Appleton Laboratory, Didcot OX11 0QX, UK

###### Abstract

The advent of x-ray free-electron lasers (XFELs) has enabled a range of new
experimental investigations into the properties of matter driven to extreme
conditions via intense x-ray-matter interactions. The femtosecond timescales
of these interactions lead to the creation of transient high-energy-density
plasmas, where both the electrons and the ions may be far from local
thermodynamic equilibrium (LTE). Predictive modelling of such systems remains
challenging because of the substantially different timescales on which
electrons and ions thermalize, and because of the vast number of atomic
configurations that are required to describe the resulting highly-ionized
plasmas. Here we explore the evolution of systems driven to high energy
densities using CCFLY, a non-LTE, Fokker-Planck collisional-radiative code. We
use CCFLY to investigate the evolution dynamics of a solid-density plasma
driven by an XFEL, and explore the relaxation of the plasma to local
thermodynamic equilibrium on femtosecond timescales in terms of the charge
state distribution, electron density, and temperature.

###### pacs:

Valid PACS appear here

##  I Introduction

Experiments carried out at x-ray free-electron lasers (FELs) over the past
decade have shown how intense x-ray fields can produce well-characterized
plasmas at high temperatures and densities via isochoric x-ray photoabsortion
Vinko  _et al._ (2012); Fletcher  _et al._ (2015); Levy  _et al._ (2015);
Ciricosta  _et al._ (2016a). This approach has proven valuable in providing
access to physical properties and atomic processes in warm dense matter, of
importance to applications across plasma physics Ciricosta  _et al._ (2012);
Vinko  _et al._ (2014), astrophysics and planetary studies Eggert  _et al._
(2010); Smith  _et al._ (2018), and fusion energy science Gaffney  _et al._
(2018); Vinko  _et al._ (2020). Experimentally, the systems produced are
typically studied via x-ray emission or inelastic scattering Toleikis  _et
al._ (2010) spectroscopy, while the theoretical description is provided by
time-dependent atomic kinetics calculations Chung  _et al._ (2017); Ciricosta
_et al._ (2016a).

Because of the short pulse duration of FEL x-ray sources, on the order of 100
fs, the time evolution of the ionization and ionic charge state populations
must be modelled directly via collisional-radiative rate equations, and local
thermodynamic equilibrium (LTE) cannot be assumed – such models are thus
called non-LTE. In contrast, these same timescales are typically long for the
purposes of electron equilibration, and the free electron distribution is
often modelled using a Maxwell-Boltzmann distribution function at some
evolving (but instantly equilibrating) electron temperature and density. Being
able to replace the time evolution of the full electron distribution function
by a time-dependent temperature and density is convenient computationally as
it leads to a substantial simplification in modelling complexity, but it is
also carries significant physical implications. For one, the shape of the
electron distribution determines multiple properties that govern the evolution
of the atomic kinetics, including photoemission and collisional ionization
rates; these can be altered significantly by non-thermal electron
distributions.

The importance of properly accounting for electron equilibration dynamics was
demonstrated in early ultrafast experiments at the FLASH FEL with XUV pulses a
few tens of fs in duration. There, the electron relaxation time was comparable
to the typical lifetime of radiative decay, leading to fluorescence emission
from a seemingly cold (non-equilibrated) distribution, alongside a
bremsstrahlung emission spectrum arising from the non-thermal electrons
consistent with a 20-times higher effective temperature Vinko  _et al._
(2010); Medvedev  _et al._ (2011). The contribution from non-thermal electrons
to the total emission spectrum is amplified in such experiments due to the
combination of low energy XUV photoexcitation ($<100$ eV) and short pulse
lengths. In contrast, similar effects are expected to be smaller at harder
x-ray photon energies ($>$1 keV) and for typical SASE pulse durations of 60-80
fs Vinko  _et al._ (2012); Preston  _et al._ (2017). Nevertheless,
experimental access to both high x-ray intensities and ultra-short pulse
lengths ($<$10 fs) will undoubtedly lead to novel x-ray matter interaction
regimes where non-thermal electron dynamics dominate the plasma evolution,
even on x-ray emission and scattering timescales. This poses a challenge to
established modelling approaches, but also provides new opportunities to study
ultrafast electron equilibration in a hot-dense plasmas experimentally in
situ.

The effect of non-thermal electrons in extreme conditions can be modelled in
various ways, and efforts have been made in both x-ray Thomson scattering
Milder  _et al._ (2019) and x-ray emission spectroscopy investigations via
non-thermal atomic kinetics de la Varga  _et al._ (2013). A kinetic modelling
addition to the particle-in-cell (PIC) simulation suite PICLS Sentoku and Kemp
(2008); Sentoku  _et al._ (2014) allowed for the inclusion of self-consistent
radiation transport, photoionization and Auger electron transport, enabling a
full kinetic simulation of XFEL-driven plasmas Royle  _et al._ (2017). While
PIC-based models are very capable in capturing the evolution of the electron
distribution in time and in space, they lack the full predictive capability
for detailed x-ray emission spectroscopy. In contrast, non-thermal
collisional-radiative codes deploying a Vlasov-Boltzmann-Fokker-Planck
equation to treat non-Maxwellian electrons Le  _et al._ (2019) are able to
treat a vast number of ionic charge configurations, essential for predicting
x-ray spectroscopy, and allow us to investigate the effect of non-thermal
distributions on spectroscopic experiments at FEL facilities.

In this context we present CCFLY, a non-thermal, non-LTE collisional radiative
code, designed to model the time-dependent evolution of of both electron
distributions and ion states interacting with intense x-ray fields on ultra-
short timescales. The code follows on from the SCFLY suite Chung  _et al._
(2017, 2007), itself based on the FLYCHK code Chung  _et al._ (2005), which
has proven successful in modelling and interpreting x-ray FEL experiments over
the past decade Ciricosta  _et al._ (2016a).

CCFLY allows us to model the plasma evolution by tracking arbitrary electron
distribution functions during the atomic-kinetics calculation, foregoing the
need to assume a Maxwellian distribution for the free electrons. More
precisely, the Fokker-Planck numerical framework is used to describe the
evolution of the electron distribution during the elastic electron-electron
collisions and account for the electrons released or absorbed during inelastic
electron-ion interactions. A self-consistent method is achieved to compute the
density and energy balance of the free-electron population together with the
ion population dynamics in a partially ionized warm dense plasma. Thus, it
becomes possible for us to investigate ionization dynamics in both
equilibrated and non-equilibrated systems at high density. We note here that
while CCFLY tracks the time-evolution of the electron distribution and of the
ionic charge states, it does not explicitly take into account the motion of
ions. As is normally the case in atomic-kinetics modelling, there is no
molecular dynamics, and the ions have no associated kinetic energy (and thus
no temperature). This is adequate to describe the evolution of x-ray-
irradiated samples on femtosecond timescales, but starts to breakdown as the
ions are heated. However, ion heating takes place on longer picosecond
timescales, normally dictated by electron-phonon coupling dynamics Cho  _et
al._ (2011); Vinko  _et al._ (2020); Lee  _et al._ (2021); Zhang  _et al._
(2022).

In this article use CCFLY to study the effect of non-thermal electrons in Al
and Mg plasmas irradiated by x-ray pulses with characteristics achievable at
current x-ray FEL facilities. We study both normal isochoric heating
experimental setups, where the x-ray photon energy is chosen to be higher than
the ionization threshold of several ionic charge states, and resonant pumping
setups where the x-ray pulse is tuned to a specific bound-bound atomic
resonance. In both cases we find that non-thermal electron processes become
important for short pulse durations and at high intensities. We quantify the
regimes for which instantaneous thermalization assumptions hold, and provide a
more rigorous justification for the modelling assumptions employed in previous
x-ray spectroscopy experiments at the Linac Coherent Light Source (LCLS)
Andreasson  _et al._ (2011); Vinko  _et al._ (2012); Ciricosta  _et al._
(2012); Cho  _et al._ (2012); Brown  _et al._ (2014); Rackstraw  _et al._
(2014); Vinko  _et al._ (2015); Ciricosta  _et al._ (2016b, a); Preston  _et
al._ (2017); Kasim  _et al._ (2018); van den Berg  _et al._ (2018); Vinko  _et
al._ (2020).

##  II Theory

The atomic kinetics modelling approach aims to track the evolution of the
charge states in the plasma as it interacts with an intense x-ray pulse. The
model assumes we have ions in various ionization states, interacting with a
classical free electron gas (the ionized free electrons), the density of which
is constrained by the requirement for the overall plasma to be charge-neutral.
We use a zero dimensional collisional-radiative model, where the transitions
between various ion states are driven by electron collisional processes, or
are mediated by the photon field. The evolution of the system is described in
terms of the changes in the relative population density $N_{i}(t)$ of each
ionic charge state (labeled by the index $i$), and of the free electron
distribution $f_{E}(E,t)$. A collisionally-induced transition between two ion
states given by some cross section $\sigma_{\pm}(E)$ is described in the model
by a transition rate per unit volume $S^{\pm}(t)/V$ given by

$S^{\pm}(t)/V=\int\sqrt{2E/m_{e}}\;\sigma_{\pm}(E)\;f_{E}(E,t)\;dE,$ (1)

where $E$ is the kinetic energy of the electron and $m_{e}$ the electron mass.
By computing the transition rates between all the ion states, we can construct
a transition matrix $S$. For each available ion state $i$, its population
density $N_{i}$ can then be determined by the following set of coupled rate
equations

$\frac{N_{i}(t+dt)-N_{i}(t)}{dt}=N_{i}(t)\left(\sum_{j\neq
i}S_{ji}-\sum_{k\neq i}S_{ik}\right).$ (2)

Because our system is globally neural, the total electron density $N_{e}$ can
be calculated from the ion density and the ionization $z$ as:

$N_{e}(t+dt)=\sum_{i}N_{i}(t+dt)z_{i},$ (3)

where we have summed across all ion states and multiplied them by their
ionization $z_{i}$. From the perspective of energy balance, the ion
transitions can either add energy or draw energy from the free electron
population. Thus, due to the requirement of energy conservation, the energy
change in the free electron distribution is given by:

$\frac{dE_{\rm{kin}}(t)}{dt}=-\sum_{i\neq
j}S_{ij}E_{ij}N_{i}(t)+\frac{dE_{rad}(t)}{dt},$ (4)

where $E_{ij}$ is the energy difference in between ionic charge states $i$ and
$j$ involved in the transition $i\rightarrow j$. This expression sums up all
the energies that are released into, or drawn from, the free electron
distribution due to photo-absorption, Auger decay, and collisional
interactions. This includes bound-bound radiative transitions, bound-bound
collisional transitions, bound-free radiative ionisation, bound-free radiative
recombination, bound-free collisional ionisation, bound-free collisional
recombination, autoionization, and electron capture.

We will call the model thermal if we assume that the electron distribution
function is Maxwellian, uniquely described by an electron density $N_{e}$ and
a temperature $T$. For situations where energy is entering the system, for
example via x-ray photoabsorption, the temperature and density are
instantaneously equilibrated at each time step. Then,

$\frac{3}{2}\,N_{e}(t)\,k_{\rm B}T(t)=E_{\rm kin}(t).$ (5)

In this case, Eqs. (3) and (4) completely determine the free-electron
distribution in terms of the changing temperature and density of the ionized
electrons.

In contrast, if a thermal electron distribution cannot be assumed, then the
change in the electron distribution needs to be modelled explicitly for all
electron energies. This can be accomplished using the Fokker-Plank equation:

$\displaystyle\frac{\partial f_{E}(E)}{\partial t}|_{ee}=$
$\displaystyle-\frac{\partial}{\partial E}[a(E)\;f_{E}(E)]$
$\displaystyle+\frac{1}{2}\frac{\partial^{2}}{\partial E^{2}}[D(E)\;f_{E}(E)]$
$\displaystyle+I(f_{E}(E))+S(E).$ (6)

The change in the electron distribution as a function of time is governed by
four terms on the right hand side of this equation. The first two terms,
containing $a(E)$ and $D(E)$, describe the elastic part of the evolution. The
third term $I(f_{E}(E))$ denotes all contributions of inelastic collisional
processes within the electron distribution, where the free electrons either
gain or lose energy to the ionic charge state. The last term $S(E)$ denotes
the collection of all source interactions that do not depend on the electron
distribution. By capturing the effects on the electron distribution in an
inelastic collision operator that considers the photo-ionisation and Auger
decay, a self-consistent approach for the ion populations and electrons can be
achieved. Then the ionic charge state populations, the electron distribution
and its time-derivatives can be fed into an ordinary differential equation
(ODE) solver. We refer to this situation as the non-thermal case. Note that
the ion population is never assumed to be thermal, and the time-dependent
modelling is purely non-LTE for both the thermal and non-thermal simulation
approaches.

###  II.1 Atomic data

CCFLY supports two levels of detail in the atomic data used in the atomic
kinetics modelling: superconfigurations and relativistic configurations. In
the super-configurational approach, we simply count how many electrons are
bound with each quantum number. States are denoted in the
$K_{a}L_{b}M_{c}N_{d}$ notation which means that there are $a$ electrons in
the K-shell (n = 1), $b$ electrons in the L-shell (n = 2), $c$ electrons in
the M-shell (n = 3) and $d$ electrons in the N-shell (n = 4). The advantage of
using superconfigurations is that the dataset can be kept relatively small
compared with a more detailed configurational dataset. In contrast,
relativistic configurations provide access to more detailed atomic physics
modelling, but at a substantially greater computational cost. A relativistic
configurations dataset can be generated using the Flexible Atomic Code Gu
(2008) Gu (2004), via the CFACDB library that provides access to the cFAC-
generated databases.

For the study presented in this work we will use the superconfigurational
dataset for our atomic kinetics modelling. In contrast, a more detailed atomic
model using relativistic configurations calculated by the DHS code Chen  _et
al._ (1980) was used for the spectral synthesis where energy levels,
statistical weights, and radiative data are applied to produce frequency-
dependent emissivities and opacities. Population distributions of relativistic
configurations belonging to a super-configuration (or to a super-shell) were
determined by Boltzmann statistics as a function of electron temperature, as
described in more detail in ref. Chung  _et al._ (2005, 2016); Ciricosta  _et
al._ (2016a).

###  II.2 Continuum lowering

The effect of continuum lowering, where the binding energies of atomic states
are lowered due to high electron and ion density, is accounted for in the
atomic kinetic modelling via an ionisation potential depression (IPD) model.
The CCFLY code supports two main IPD models: the modified Ecker-Kroll IPD
model (EK) Preston  _et al._ (2013) and the Stewart-Pyatt IPD model Stewart
and Pyatt (1966). All simulations presented here were carried out with the EK
model, to make the results directly comparable to previously published related
work.

###  II.3 Elastic Electron Collisions

The elastic part of the evolution of the electron distribution is given by the
first two terms in the Fokker Planck equation in Eq. II. The elastic collision
terms can alter the shape of the distribution, while preserving particle
density and total kinetic energy. Additionally, the distribution function
should equilibrate to the correct Maxwellian distribution. In CCFLY, the
electron distribution is initialised and evolved in energy space but in the
Fokker-Planck operator, the electrons are evolved in velocity space.

Figure 1: Elastic collisional thermalization of the electron distribution
function: we initialze a nonthermal distribution consisting of a Maxwellian at
a temperature of 40 eV combined with a hot-electron spike at 333 eV. The
distribution equilibrates to a hotter Maxwellian distribution on fs
timescales.

First derived by Bobylev and Chuyanov Bobylev and Chuyanov (1976), and cast
into a non-relativistic collision operator by Tzoufras et al. Tzoufras  _et
al._ (2011)Tzoufras  _et al._ (2013), this formulation is given in integro-
differential operator form as:

$\left(\frac{\partial f_{v}(v)}{\partial
t}\right)_{ee}=\frac{4\pi\Gamma_{ee}}{3}\frac{1}{v^{2}}\frac{\partial}{\partial
v}\left[\frac{1}{v}\frac{\partial W\left(v\right)}{\partial v}\right],$ (7)

where the function $W(v)$ is given by:

$\displaystyle W\left(v\right)=$ $\displaystyle
f_{v}\int_{0}^{v}f_{u}u^{4}du\;+\;v^{3}f_{v}\int_{v}^{\infty}f_{u}udu\;$
$\displaystyle-3\int_{v}^{\infty}f_{u}udu\int_{0}^{v}f_{u}u^{2}du.$ (8)

The function $W(v)$ contains the time evolution of the electron distribution
due to to electron-electron collisions that conserve particle number and
kinetic energy. The $\Gamma_{ee}$ is defined as:

$\Gamma_{ee}=\frac{e^{4}\ln\Lambda}{4\pi\varepsilon_{0}^{2}m_{e}^{2}}$ (9)

with $e$ the electron charge, and $\ln\Lambda$ the Coulomb logarithm. The
Coulomb logarithm for electron-electron interactions in hot plasmas is
computed with a semi-empirical formula constructed by Huba Huba (2013):

$\ln\Lambda=23.5-\ln\left(N_{e}^{1/2}T^{-5/4}\right)-\left[10^{-5}+(\ln
T-2)^{2}/16\right]^{1/2},$ (10)

where the density and temperature have units of cm-3 and eV, respectively.

For large velocities, Eq. (8) can be evaluated through standard discretisation
techniques such as the chained trapezoid method. For velocities that are small
compared with the characteristic velocity $v_{c}$ of the distribution, the
second and third term in Eq. (8) cancel up to $\mathcal{O}(v^{5})$. It is
necessary to modify the numerical scheme to capture this, in order for the
distribution to equilibrate to a Maxwellian. This can be done by writing
$W(v_{n})$ as a Taylor series $f_{n}\approx
f(0)+\left[f_{1}-f(0)\right]\left(v_{n}^{2}/v_{1}^{2}\right)$ and
$f(0)=\left(f_{0}-f_{1}v_{0}^{2}/v_{1}^{2}\right)/\left(1-v_{0}^{2}/v_{1}^{2}\right)$.
Setting $W\left(v_{n}<v_{c}\right)\equiv W_{s}(v)$, this results in:

$\displaystyle W_{s}(v)=$ $\displaystyle
f_{n}\times\left[f(0)\frac{v_{n}^{5}}{5}+\left(f_{1}-f(0)\right)\frac{v_{n}^{7}}{7v_{1}^{2}}\right]$
(11)
$\displaystyle+\left[\left(f_{n}-f(0)\right)v_{n}^{3}-\left(f_{1}-f(0)\right)\frac{3v_{n}^{5}}{5v_{1}^{2}}\right]$
$\displaystyle\times\sum_{k=n+1}^{N-1}\frac{1}{2}\left(f_{k}v_{k}+f_{k-1}v_{k-1}\right)\delta_{k-1/2}.$

Similarly, setting $W\left(v_{n}\geq v_{c}\right)\equiv W_{l}(v)$, we have

$\displaystyle W_{l}(v)=$ $\displaystyle
f_{n}\sum_{k=0}^{n}\frac{1}{2}\left(f_{k}v_{k}^{4}+f_{k-1}v_{k-1}^{4}\right)\delta_{k-1/2}$
(12)
$\displaystyle+v_{n}^{3}f_{n}\sum_{k=n+1}^{N-1}\frac{1}{2}\left(f_{k}v_{k}+f_{k-1}v_{k-1}\right)\delta_{k-1/2}$
$\displaystyle-3\left[\sum_{k=n+1}^{N-1}\frac{1}{2}\left(f_{k}v_{k}+f_{k-1}v_{k-1}\right)\delta_{k-1/2}\right]$
$\displaystyle\times\left[\sum_{k=1}^{n}\frac{1}{2}\left(f_{k}v_{k}^{2}+f_{k-1}v_{k-1}^{2}\right)\delta_{k-1/2}\right].$

We show an example of how this elastic collisional operator relaxes a non-
thermal electron distribution to a Maxwellian in Fig. 1. For simplicity, we
assume here that only free electrons are involved in the collisional
thermalization process, and all interactions with ions (e.g. collisional
ionization and recombination) are ignored. We initialize the electrons to a
thermal distribution at 40 eV, and add an additional non-thermal population of
electrons with energies around 333 eV. The collisional operator equilibrated
the distribution to a new overall temperature of 42 eV, the majority of which
takes place on a timescale of a few femtoseconds.

##  III Non-thermal electrons in x-ray isochoric heating

The capabilities of the CCFLY code allow us to look into the detailed
evolution of the electron distribution function in an isochoric heating
experiment at an x-ray FEL. In particular, we are interested in exploring what
role the non-thermal electrons play in the atomic kinetics of the system, and
wish to assess the range of validity of the instantaneous electron
thermalization assumption commonly used in the modelling and in the
interpretation of experimental results.

Figure 2: Time evolution of the electron distribution in Al irradiated by
x-rays at 1590 eV at various peak intensities. The non-thermal electrons
generated from photoionization and Auger processes at higher energies are
absent in the thermal calculations.  Figure 3: Time evolution of the charge
state populations in Al irradiated by x-rays at 1590 eV at various peak
intensities. Thermal and non-thermal calculations are compared.

We start by modelling the heating process of a thin Al foil irradiated by
x-rays at a photon energy of 1590 eV. This energy was deliberately chosen to
be above the Al K-edge (1560 eV) so that the dominant absorption mechanism is
K-shell (1s) photoionization. The x-ray pulse is chosen to be a Gaussian in
time with a FWHM duration of 40 fs, and has a spectral bandwidth of 0.4%,
typical of a SASE FEL pulse. Four different intensities are compared, between
$2.5\times 10^{15}$ Wcm-2 and $2.5\times 10^{18}$ Wcm-2, depositing an ever
increasing amount of energy into the target. Experimentally, these intensities
would be realized by focusing the x-rays to various spot sizes: focusing the
beam to a spot of a few micron in diameter typically yields intensities in
region of $10^{17}$ Wcm-2 and can be readily achieved on most facilities.

We show the evolution of the electron density in Fig. 2, where we compare the
results of the thermal and non-thermal models side by side. The x-ray pulse
peaks in intensity in the middle of the simulation window, at 40 fs. The
thermal simulations simply show a Maxwellian distribution that gets
increasingly hotter in time over the duration of the excitation pulse. More
deposited energy leads to increasing ionization, and thus the free electron
density also increases over time. The temperatures and densities are higher
for higher irradiation intensities, since more energy is deposited into
smaller volumes, as expected. The electron distribution function for the non-
thermal simulation is visibly different, and can be divided into two
components. The first is the bulk part of the free electron distribution,
lying at lower energies. This assumes a near-Maxwellian shape because the
electron thermalization time is short even on the 100 fs timescales used here
(see Fig. 1 and related discussion). This bulk component is largely
equilibrated, and grows in prominence as the x-rays are absorbed and the
sample heats. At higher energies, however, we can see a large energy spread
over which non-thermal electrons can be found. These are mostly created by the
equilibrating hot electrons generated in the x-ray-matter interaction process,
and subsequent atomic kinetics. Two main sources of quasi-mono-energetic hot
electrons can be identified: a feature at due to the photoionization of the Al
L-shell electrons (at $\approx$1520 eV), and a feature due to the Auger KLL
process (at $\approx$1420 eV) that efficiently fills a K-shell core hole. Both
features can be readily seen in the figure, but are clearest at the lowest
intensities. They are harder to discern at high intensity because their
relative populations, compared with the bulk free-electron distribution,
decreases with intensity.

Figure 2 shows the density distribution on a logarithmic scale, and thus
amplifies the non-thermal signature. But for the underlying physics, what
matters most is the relative contribution of the quasi-thermalized component
of the electron distribution in driving the atomic kinetics, compared with the
non-thermal hot electron population. To examine this, we plot the populations
for ionic charge states $3+$ to $10+$ in Fig. 3, in both the thermal and non-
thermal case. We immediately note that the results are comparable: the charge
state populations do not change dramatically under the instantaneous
thermalization approximation, validating this common assumption. Nevertheless,
we do see a lag in the time at which certain charge states peak in the
population distribution between the thermal and non-thermal calculations. In
particular, we note that lower charge states in the more detailed non-thermal
calculations tend to appear later in time compared with the thermal
prediction, whereas higher charge states appear earlier in time. This is
primarily due to the mechanism through which these charge states are created.
The lower charge states are more efficiently generated via collisional
ionization processes driven by the bulk free electrons; these electrons heat
up faster in the thermal calculation as all energy is placed into the thermal
distribution instantaneously, whereas in the non-thermal case the bulk
temperature rises more slowly as the energy needs to be transferred from the
highly energetic photoionization and Auger electrons. In contrast, the higher
charge states require substantially higher temperatures to be collisionally
ionized, and are thus more efficiently created via collisional interactions
with the non-thermal electrons. Even if there are fewer of them available,
they have sufficient energy to play an important, and possibly dominant role,
in the dynamics of collisional ionization.

In an experiment one cannot normally measure individual plasma charge state
populations. However, these populations can be inferred via x-ray emission
spectroscopy of Kα satellite lines, when the sample is excited by x-rays above
the relevant ionization thresholds Vinko  _et al._ (2012); Vinko (2015). For
charge states photo-pumped below their ionization threshold, the appearance of
emission lines can be instead indicative of collisional ionization processes,
and thus used to measure collisional rates Vinko  _et al._ (2015). For this
measurement to be reliable, however, the bulk free electrons must be the
dominant contributor to collisional ionization, and the contributions from
non-thermal electrons should be small.

We investigate the contribution of non-thermal electrons to the total Kα
emission spectrum in Fig. 4. The x-ray photon energy is set just above the
K-edges of the first two charge states, but below the ionization threshold of
higher charge states. Thus only the first two charge states can be revealed
directly via Kα emission spectroscopy. While higher charge states are still
present in the plasma, as illustrated in Fig. 3, their $1s$ states cannot be
photo-ionized and they therefore cannot relax emitting a Kα photon. Despite
this, emission from higher charge states can still be observed, because of
electron collisional ionization of the L-shell in the presence of a K-shell
hole van den Berg  _et al._ (2018). Non-thermal electrons affect the emission
spectrum in two key ways. Firstly, because the equilibration of
photoionization and Auger electrons is the main heating channel of the bulk
free-electron distribution, non-thermal electrons help determine the
equilibration timescales and drive the evolution of the charge state
population towards LTE. These evolving populations ultimately determine the
emission intensity of charge states pumped above their K-edge thresholds. The
second mechanism applies to charge states that are pumped below their K-shell
ionization threshold. These are primarily driven by collisional ionization of
the L-shell electrons, which can be either due to the bulk electrons at some
average temperature $k_{\rm B}T\approx 100$ eV, or to the non-thermal
electrons at a much higher energy (typically $\approx 10\,k_{\rm B}T$). Higher
electron energies are more efficient at L-shell ionization, and so a system
with a large non-thermal electron population and relatively cold bulk
distribution will lead to more collisional ionization, and consequently to
more intense emission from higher charge states. We see this play out in the
panels of Fig. 4. Because the x-ray pulse length is relatively long compared
with the electron relaxation timescale, the thermal and non-thermal modelling
predict very similar intensities for the first two emission lines. The charge
state population is very close to LTE for intensities up to around 1018 Wcm-2,
and the intensity of the first two lines is simply a reflection of the
population of those charge states modulated by the temporal intensity of the
x-ray pulse: the second charge state is most prominent at irradiation
intensities of around 1016 Wcm-2, and decreases for lower (insufficient
heating) and higher (too much heating) intensities. We note a slight
difference in the intensity of the second line at the highest intensities:
here the charge state is swept through very quickly, and the instantaneous
thermalization assumption starts to breakdown as it underestimates the time
needed to establish LTE in the plasma distribution. The emission intensity of
higher lines (VI and above) shows a considerably stronger dependence on the
x-ray intensity. For intensities up to several times 1016 Wcm-2 the thermal
and non-thermal predictions match, indicating that the L-shell ionization
process described above is indeed dominated by quasi-thermal bulk free
electrons. Above this, the number of non-thermal electrons remains
sufficiently high to more than double the effective ionization rate,
indicating a stark departure from LTE collisional dynamics.

Figure 4: Emission spectrum of Al irradiated by x-rays at 1590 eV at various
peak intensities, comparing thermal and non-thermal simulations. At this
photon energy the emission lines IV and V can be pumped directly via K-shell
photoionization, while lines VI to VIII emit because of collisional processes.

The results shown in Fig. 4 are for a single x-ray intensity, and are thus not
directly comparable with experimental results that are measured over a broad
range of x-ray intensities. However, the average experimental intensity
(rather than the often quoted peak intensity) is indicative of the overall
irradiation regime, and can be linked to the discussion above. Average
intensities typically achieved at the LCLS FEL, such as those fielded for the
work on collisional ionization presented in Refs. Vinko  _et al._ (2015); van
den Berg  _et al._ (2018), is on the order of $10^{16}$ Wcm-2.

###  III.1 The effect of x-ray pulse length

We have seen that for a given x-ray pulse duration, the effect of focusing the
x-ray beam to achieve increasing intensities is to increase the relative
contribution of non-thermal electrons to the evolution of the plasma, and thus
to the observed emission spectrum. We now turn our attention to the effect of
pulse length. Changing the pulse length is another way to change the intensity
on target, but it adds an additional layer of complexity because the pulse
duration also gates the observation timescale of the system. This is primarily
because we only observe x-ray emission in the Kα window when $1s$ electrons
can be ionized, and the overall temperatures reached in a typical
microfocusing experiment ($\approx 100$ eV) are insufficient to drive
collisional K-shell ionization to any notable extent. Thus, emission is only
observed during the x-ray pulse, and for a short time after depending on the
femtosecond core-hole lifetime. By changing the duration of the x-ray pulse,
while maintaining a constant pulse energy, we thus not only modify the x-ray
intensity seen by the plasma, but also change the period over which we view
the evolution of the charge state distribution. We show the simulated results
for a Al plasma in Figs. 5 and 6. For these simulations we have assumed a
pulse energy of 1 mJ, and a spotsize on-target of 12 $\mu$m2, so that a pulse
length of 80 fs corresponds to an x-ray intensity of about 1017 Wcm-2. These
are typical values achievable at FEL facilities today Vinko  _et al._ (2012).
The pulse energy is then kept constant, and the pulse duration lowered down to
8 fs.

Figure 5: Evolution of the charge state populations in Al irradiated at
constant energy for various pulse lengths. Thermal and non-thermal
calculations are compared.

For the longer pulse lengths of 40 and 80 fs we see little difference in the
computed emission spectrum for the thermal and non-thermal cases, and a lag in
the creation of various charge states analogous with that described in the
previous section. In these cases, the observation window is sufficiently long
for the evolution of the change state populations to agree well with the
assumption of instantaneous electron thermalization. The overall population of
the L-shell states is collisionally driven by the equilibrated bulk free
electron distribution, and the charge states evolve as if in near-LTE for each
point in time. For pulse durations of 20 fs and less, these assumptions start
to breakdown, and the calculated spectrum assuming instantaneous electron
thermalization underestimates the emission intensity.

Figure 6: Emission spectra of Al irradiated at constant energy for various
pulse lengths. Differences in the emission intensity can be seen for pulse
durations around 20 fs and below. Figure 7: Evolution of the charge state
populations in Mg irradiated on-resonance at 1304 eV, at constant energy for
various pulse lengths. Thermal and non-thermal calculations are compared.

###  III.2 Non-thermal electrons in resonant spectroscopy

Figure 8: Emission spectra of Mg irradiated on-resonance at 1304 eV, at
constant energy for various pulse lengths. The emission lines on two sides of
the resonant VII K$\alpha$ line are due to collisional ionization (right of
the K$\alpha$) and three-body recombination (left of the K$\alpha$). Non-
thermal electrons start to affect the spectrum for pulse lengths below 30 fs.

Recently, van den Berg et al. proposed driving a resonant $1s\rightarrow 2p$
transition in Mg to measure collisional ionization rates in dense plasmas van
den Berg  _et al._ (2018). The principle is straightforward: an x-ray pulse
form an FEL is tuned to resonantly drive a specific inner-shell Kα transition,
and a spectroscopic measurement is performed to look for emission from
neighbouring Kα emission lines. Because of the narrow bandwidth of the FEL
beam, the x-ray pulse can only excite a single Kα line, and so the observation
of multiple lines is indicative of L-shell collisional processes taking place
during the lifetime of a K-shell core hole. This method thus allows for a
measurement of an ultra-fast collisional ionization rate. However, the real
objective is not so much to infer a rate, but rather to measure a collisional
cross section. This can be accomplished using the expression in Eq. (1), under
the caveat that the collisions are predominantly due to an interaction of the
ion with the bulk free-electron distribution at some known temperature and
density. But if the electron distribution is non-thermal, then it is easy to
see that the cross section cannot uniquely be extracted from the integral
expression of Eq. (1). Estimating the effect of non-thermal electrons is then
important to guide the reliable interpretation of such experiments, especially
since we have shown that the contribution of non-thermal electrons can be
significant with increasing x-ray intensity or decreasing x-ray pulse length.

For this investigation we set out simulation parameters to match the
conditions fo Ref. van den Berg  _et al._ (2018): an optically thin Mg plasma
is pumped at a photon energy of 1304 eV, corresponding to the resonance energy
of the $1s\rightarrow 2p$ transition energy in a 7+ Mg ion. At the nominal
pulse duration of 60 fs the intensity on target is $10^{17}$ Wcm-2. We show
the evolution of the ionic populations as a function of time, using the
thermal and a non-thermal models, in Fig. 7. The related emission spectra are
shown in Fig. 8. Here we see that the assumption of instantaneous
thermalization holds well for nominal pulse durations, above $\sim$30 fs,
providing a more rigorous justification of the assumptions made in the work on
collisional ionization rates. For shorter pulse durations, however, the
contribution of non-thermal electrons starts to become significant, and even
dominant for pulse lengths below 10 fs. At these short pulse lengths,
modelling the full evolution of the electron distribution function becomes
essential to produce atomic kinetics simulations with predictive capability.

##  IV Conclusions

We have presented the atomic kinetics simulation package CCFLY, which can
model the non-thermal (electron) and non-LTE (ion) evolution of matter
interacting with intense x-ray radiation on ultra-short timescales. We have
used these capabilities to investigate the evolution of solid-density systems,
heated isochorically to high temperatures via femtosecond FEL pulses, and
assess the role of non-thermal electron distributions. We find that in mid-Z
systems such as Al and Mg, non-thermal electrons become important in
determining the short-time plasma kinetics and electron dynamics for
intensities above $10^{17}$ Wcm-2, and for pulse lengths below 10 fs.

To date, the use of Fokker-Planck models within collisional-radiative
simulations to interpret experiments at free-electron lasers has been limited,
largely due to the substantial computational cost the method requires. Whilst
we have shown that this is often justifiable for experiments using a typical
SASE FEL pulse on the timescales of order 100 fs, ongoing developments
enabling the routine generation of energetic few-fs and sub-fs pulses, and of
ultra-high intensities via narrow bandwidth self-seeding and x-ray
nanofocusing, are increasingly accessing x-ray regimes where the direct
modelling of non-thermal electron relaxation will be essential to retain
predictive power in the modelling of extreme x-ray matter interactions.

##  V Acknowledgments

S.R., M.F.K., J.S.W. and S.M.V. acknowledge support from the UK EPSRC grant
EP/P015794/1 and the Royal Society. S.M.V. is a Royal Society University
Research Fellow.

## References

  * Vinko  _et al._ (2012) S. M. Vinko _et al._ , Nature 482, 59 (2012).
  * Fletcher  _et al._ (2015) L. B. Fletcher _et al._ , Nature Photonics 9, 274 (2015).
  * Levy  _et al._ (2015) A. Levy _et al._ , Physics of Plasmas 22, 030703 (2015).
  * Ciricosta  _et al._ (2016a) O. Ciricosta _et al._ , Physics of Plasmas 23, 022707 (2016a).
  * Ciricosta  _et al._ (2012) O. Ciricosta _et al._ , Phys. Rev. Lett. 109, 065002 (2012).
  * Vinko  _et al._ (2014) S. M. Vinko, O. Ciricosta, and J. S. Wark, Nat. Commun. 5, 3533 (2014).
  * Eggert  _et al._ (2010) J. H. Eggert _et al._ , Nature Physics 6, 40 (2010).
  * Smith  _et al._ (2018) R. F. Smith _et al._ , Nature Astronomy 2, 452 (2018).
  * Gaffney  _et al._ (2018) J. Gaffney _et al._ , High Energy Density Physics 28, 7 (2018).
  * Vinko  _et al._ (2020) S. M. Vinko _et al._ , Phys. Rev. Lett. 124, 225002 (2020).
  * Toleikis  _et al._ (2010) S. Toleikis _et al._ , Journal of Physics B: Atomic, Molecular and Optical Physics 43, 194017 (2010).
  * Chung  _et al._ (2017) H.-K. Chung _et al._ , AIP Conference Proceedings 1811, 020001 (2017).
  * Vinko  _et al._ (2010) S. M. Vinko _et al._ , Phys. Rev. Lett. 104, 225001 (2010).
  * Medvedev  _et al._ (2011) N. Medvedev, U. Zastrau, E. Förster, D. O. Gericke, and B. Rethfeld, Phys. Rev. Lett. 107, 165003 (2011).
  * Preston  _et al._ (2017) T. R. Preston _et al._ , Phys. Rev. Lett. 119, 085001 (2017).
  * Milder  _et al._ (2019) A. L. Milder, S. T. Ivancic, J. P. Palastro, and D. H. Froula, Physics of Plasmas 26, 022711 (2019).
  * de la Varga  _et al._ (2013) A. G. de la Varga, P. Velarde, F. de Gaufridy, D. Portillo, M. Cotelo, A. Barbas, A. González, and P. Zeitoun, High Energy Density Physics 9, 542 (2013).
  * Sentoku and Kemp (2008) Y. Sentoku and A. Kemp, Journal of Computational Physics 227, 6846 (2008).
  * Sentoku  _et al._ (2014) Y. Sentoku, I. Paraschiv, R. Royle, R. C. Mancini, and T. Johzaki, Phys. Rev. E 90, 051102 (2014).
  * Royle  _et al._ (2017) R. Royle, Y. Sentoku, R. C. Mancini, I. Paraschiv, and T. Johzaki, Phys. Rev. E 95, 063203 (2017).
  * Le  _et al._ (2019) H. P. Le, M. Sherlock, and H. A. Scott, Phys. Rev. E 100, 013202 (2019).
  * Chung  _et al._ (2007) H.-K. Chung, M. Chen, and R. Lee, High Energy Density Physics 3, 57 (2007).
  * Chung  _et al._ (2005) H.-K. Chung, M. Chen, W. Morgan, Y. Ralchenko, and R. Lee, High Energy Density Physics 1, 3 (2005).
  * Cho  _et al._ (2011) B. I. Cho _et al._ , Phys. Rev. Lett. 106, 167601 (2011).
  * Lee  _et al._ (2021) J.-W. Lee _et al._ , Phys. Rev. Lett. 127, 175003 (2021).
  * Zhang  _et al._ (2022) J. Zhang, R. Qin, W. Zhu, and J. Vorberger, Materials 15 (2022), 10.3390/ma15051902.
  * Andreasson  _et al._ (2011) J. Andreasson _et al._ , Phys. Rev. E 83, 016403 (2011).
  * Cho  _et al._ (2012) B. I. Cho _et al._ , Phys. Rev. Lett. 109, 245003 (2012).
  * Brown  _et al._ (2014) C. R. D. Brown _et al._ , Scientific Reports 4, 5214 (2014).
  * Rackstraw  _et al._ (2014) D. Rackstraw _et al._ , High Energy Density Physics 11, 59 (2014).
  * Vinko  _et al._ (2015) S. M. Vinko _et al._ , Nat Commun 6, 6397 (2015).
  * Ciricosta  _et al._ (2016b) O. Ciricosta _et al._ , Nat. Commun. 7, 11713 (2016b).
  * Kasim  _et al._ (2018) M. F. Kasim, J. S. Wark, and S. M. Vinko, Scientific Reports 8, 6276 (2018).
  * van den Berg  _et al._ (2018) Q. Y. van den Berg _et al._ , Phys. Rev. Lett. 120, 055002 (2018).
  * Gu (2008) M. F. Gu, Canadian Journal of Physics 86, 675 (2008).
  * Gu (2004) M. F. Gu, AIP Conference Proceedings 730, 127 (2004).
  * Chen  _et al._ (1980) M. H. Chen, B. Crasemann, and H. Mark, Phys. Rev. A 21, 436 (1980).
  * Chung  _et al._ (2016) H.-K. Chung, S. B. Hansen, and H. A. Scott, “Generalized collisional radiative model using screened hydrogenic levels,” in  _Modern Methods in Collisional-Radiative Modeling of Plasmas_, edited by Y. Ralchenko (Springer International Publishing, Cham, 2016) pp. 51–79.
  * Preston  _et al._ (2013) T. R. Preston, S. M. Vinko, O. Ciricosta, H.-K. Chung, R. W. Lee, and J. S. Wark, High Energy Density Physics 9, 258 (2013).
  * Stewart and Pyatt (1966) J. C. Stewart and J. Pyatt, Kedar D., Astrophys. J. 144, 1203 (1966).
  * Bobylev and Chuyanov (1976) A. Bobylev and V. Chuyanov, USSR Computational Mathematics and Mathematical Physics 16, 121 (1976).
  * Tzoufras  _et al._ (2011) M. Tzoufras, A. R. Bell, P. A. Norreys, and F. S. Tsung, J. Comput. Phys. 230, 6475 (2011).
  * Tzoufras  _et al._ (2013) M. Tzoufras, A. Tableman, F. S. Tsung, W. B. Mori, and A. R. Bell, Physics of Plasmas 20, 056303 (2013).
  * Huba (2013) J. D. Huba,  _Plasma Physics_ (Naval Research Laboratory, Washington, DC, 2013) pp. 1–71.
  * Vinko (2015) S. M. Vinko, Journal of Plasma Physics 81, 365810501 (2015).