We study the energy spectrum of laser-induced conduction band (CB) electrons in water by multi-rate equations (MRE) with different impact ionization schemes. Rethfeld’s MRE model [Phys. Rev. Lett. 92, 187401(2004) Phys. Rev. B 79, 155424(2009)], but the corresponding rate equations are computationally very expensive. We introduce a simplified splitting scheme and corresponding rate equations that still agree with energy conservation but enable the derivation of an asymptotic SRE. This approach is well suited for the calculation of energy spectra at long pulse durations and high irradiance, and for combination with spatiotemporal beam propagation/plasma formation models. Using the energy-conserving MREs, we present the time-evolution of CB electron density and energy spectrum during femtosecond breakdown as well as the irradiance dependence of free-electron density, energy spectrum, volumetric energy density, and plasma temperature. These data are relevant for understanding photodamage pathways in nonlinear microscopy, free-electron-mediated modifications of biomolecules in laser surgery, and laser processing of transparent dielectrics in general.
© 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
Short-pulse lasers are widely used for plasma-mediated material processing [1–4], laser surgery [5–11], and photomodification of biomolecules [12–14]. The thresholds for phase transitions and ablation in transparent dielectrics depend on conduction band (CB) electron density and the resulting volumetric energy density. For a long time modeling efforts focused on the calculation of free electron density [2,15–18], and simple estimations of the average energy of free electrons were used to assess the volumetric plasma energy density [6,19,20].
However, exact knowledge of the average kinetic electron energy would create a more precise link between electron density and energy density, and the energy distribution of CB electrons determines which chemical changes may be induced by free electrons, e-. Thus, knowledge of the energy spectrum is pivotal for an understanding of nanosurgery by femtosecond (fs) pulse series that relies on the cumulative induction of molecular bond breaks [6,21], for analysis of fs laser-induced molecular modifications [12,14], and also for the mitigation of photodamage in nonlinear microscopy .
The time evolution of electron energy spectra was first analyzed for wide-bandgap solids such as SiO2 [23,24], and sophisticated simulations for SiO2 have been presented in recent years [25,26]. These models are based on solving the Boltzmann kinetic equations, which is computationally expensive and requires detailed knowledge on material parameters and various scattering cross sections. Such data are not available for water, and simpler approaches must be pursued.
In 2004, a multi-rate-equation (MRE) approach capable of tracing the time evolution of the energy distribution of CB electrons during fs breakdown was presented by Rethfeld  in order to describe the onset of impact ionization that influences the relative importance of strong-field ionization (SFI) and avalanche ionization (AI) [2,18,19]. The MRE model introduces different energy levels within the conduction band and traces the time evolution of the population density at each level. For longer pulse durations, Rethfeld derived an asymptotic solution, in which AI is described by a single-rate equation (SRE). However, the early work focused on the time evolution of the avalanche ionization rate and the number density of electrons available for impact ionization but did not explicitly address the CB electron energy spectrum and the volumetric energy density in the dielectric. Furthermore, it did not consider the excess energy of CB electrons remaining after impact ionization.
This problem was addressed by Christensen and Balling who introduced an energy splitting scheme to model the distribution of residual energy after impact ionization into the individual CB energy levels assumed in the MRE model . Unfortunately, this scheme is complex and no approximate solution is available. Therefore, it becomes computationally very expensive for long pulse durations and 3D modeling of plasma formation.
While for many dielectrics the band structure is fully characterized by the energy gap Egap between valance band (VB) and CB, water exhibits an intermediate energy level between VB and CB. Recent investigations of the wavelength dependence of the irradiance threshold, Ith, for nanosecond optical breakdown revealed steps in the Ith(λ) curve that are consistent with breakdown initiation by multiphoton ionization, with an initiation energy of about 6.6 eV . This value is considerably smaller than the excitation energy providing a significant rate of auto-ionization (≈9.5 eV), which defines the band gap relevant for avalanche ionization [28–31]. Thus, breakdown initiation is likely to occur via excitation of a VB electron into a solvated state, followed by rapid excitation into the conduction band .
In this paper, we adapt previous MRE models to the band structure and breakdown initiation dynamics of water. Furthermore, we design a simple impact ionization scheme that conserves energy and enables the derivation of an asymptotic SRE, like in Rethfeld’s model. The models are used to analyze the time evolution of the electron energy spectrum for different wavelengths under conditions typical for nonlinear imaging and nanosurgery by fs pulse series, and to calculate the irradiance dependence of electron spectra, mean CB electron energy, and volumetric energy density up to full ionization.
2. Ionization schemes for water
Figure 1 illustrates the three ionization schemes for water that are discussed in this paper. In all schemes, VB electrons are excited into the bottom of the CB via strong field ionization (SFI). SFI may occur either directly across the entire bandgap, Egap, or via an intermediate level, Esolv. The existence of this intermediate level is due to the network of weak hydrogen bonds between water molecules, in which thermal fluctuations can produce favorable constellations for electron abstraction from excited molecules. They are known as “pre-existing traps” [28,32] and are located around 6.6 eV above the top of the VB, i.e. at Esolv ≈6.6 eV. Excited electrons relax into the traps and solvate. The traps have a number density of about χtrap ≈1019 cm−3, corresponding to ≈3 traps in 1000 water molecules . The solvated electrons, eaq-, are then easily upconverted into the CB due to the fact that the remaining energy gap amounts to only 3 eV and that eaq- absorbs very strongly from ultraviolet (UV) to infrared (IR). Thus, the critical step for breakdown initiation is multiphoton excitation into the intermediate level Esolv. Subsequent upconversion into the CB is assumed to follow immediately and is not considered as a separate term of the rate equation [19,28].
The total energy required to excite electrons from the valence band into the conduction band is given by the effective ionization potentialthat accounts both for the band gap and the oscillation energy of free electrons in strong electromagnetic fields [15,19,33]:19,33]. The symbol ω denotes the circular frequency of the electric laser field, e is the electron charge, n0 is the refractive index of the medium at frequency ω, and ε0 and c are the vacuum dielectric permittivity and the vacuum speed of light. The reduced exciton mass is approximated by half of the mass mc of the conduction band electrons [15,18,27].
The SFI rate considering excitation through Esolv and via Egap is19].
Once in the CB, the electrons can gain kinetic energy via inverse Bremsstrahlung absorption (IBA). They generate further free electrons once their energy is larger than the critical energy εcrit required to cause impact ionization. To satisfy the conservation laws for energy and momentum, the energy of the impacting electron must exceed. For a parabolic bandgap, the minimum required energy is (3/2) [6,18,24,34].
Characteristic differences between the ionization schemes consist in the way in which the residual energy (1/2) remaining after impact ionization is treated. Rethfeld [Scheme in Fig. 1(a)] neglected the residual energy and assumed that both electrons start at zero energy . This way, their subsequent gain by IBA can be easily assigned to discrete energy levels. The energy loss associated with this approach leads to an underestimation of the speed of the avalanche. It can be tolerated for breakdown processes in which SFI dominates and impact ionization plays only a subordinate role because here the difference between absorbed laser energy and total energy of free electrons is small. The results in  supported such a scenario. However, they were obtained with a one-photon absorption rate, W1pt, for IBA that considered only electron-phonon scattering but neglected electron-ion collisions [18,24]. That W1pt value corresponds to τcoll = 14.3 fs in the Drude model, as has been shown in Ref . In later work, Rethfeld and associates included electron-ion collisions and assumed a collision time in the order of 1 fs , in accordance with experimental results on silicon and fused silica [36,37]. Our recent findings for water corroborate this value . The shorter collision time implies an avalanche-like behavior with a high rate of impact ionization events [19,35]. Therefore, the residual energy after impact ionization needs to be considered to obtain a correct energy balance of the optical breakdown process.
The scheme in Fig. 1(b) presents a simple approach for achieving energy conservation in the description of impact ionization. The residual energy (1/2)is split into three equal parts attributed to the three partners involved in the collision process (impacting electron, new CB electron, and hole in the VB [25,27,34]). The residual energy carried by the electrons contributes to AI, whereas the energy fraction imparted to the hole is thermalized. This approach results in a constant start energy for IBA, namely εstart = (1/6).
Since the start energy is determined by the effective band gap, it is usually not an integer multiple of the energy gain by one-photon absoption, . Therefore, the energy levels in Fig. 1(b) are shifted compared to the levels in Figs. 1(a) and 1(c) that start at εstart = 0. This is unproblematic for the description of avalanche ionization but must be considered in matching the seed electron generation by SFI to the description of IBA in avalanche ionization. We assume that the seed electrons rapidly acquire the start energy (1/6) either from the laser irradiation (the energy needed for multiphoton ionization will usually exceed ), or through collisions with other CB electrons. The necessary extra-energy is illustrated by the green arrow in Fig. 1(b).
The simplifying assumption of a finite start energy of electrons both after impact ionization and also for SFI initiation of AI creates a gap in the energy spectrum at 0 < ε < (1/6). Furthermore, the assumption of a finite start energy after SFI may imply a certain error in the energy balance. These shortcomings are negligible for material processing in which avalanche ionization dominates but become relevant for photomodifications of biomolecules in which SFI plays a significant role. Here, the energy splitting scheme of Fig. 1(c) suggested by Balling and associates is more appropriate [27,38]. Their splitting scheme for the residual energy remaining after impact ionization is compatible with the assumption that electrons excited via SFI are located at the bottom of CB possessing zero kinetic energy.
The gap in the energy spectrum at 0 < ε < (1/6)is a stronger drawback than the potential error in the energy balance for multiphoton ionization. An error may also be introduced by assuming zero start energy for multiphoton initiation because in most cases the photons will carry an energy larger than the exact value of the effective ionization potential. Thus, the assumption of zero start energy slightly underestimates the speed of the breakdown process, whereas the assumption of a finite start energy may lead to a slight overestimation.
For cases in which AI dominates over SFI, i.e. for high-intensities such as used for material processing, and for long pulse durations, the simple energy splitting scheme of Fig. 1(b) enables the derivation of a single rate equation (SRE) describing avalanche ionization as will be elaborated in section 3. The MRE for the more complex splitting scheme of Fig. 1(c) is described in section 4.
3. MRE with simple energy splitting scheme, and asymptotic SRE
We start with Rethfeld’s original MRE approach for treating the avalanche ionization process. Here it is assumed that electrons excited via SFI are initially located right at the bottom of the CB possessing zero kinetic energy, i.e. they start at the n0 energy level with ε0 = 0. Electrons then gain kinetic energy via inverse Bremsstrahlung absorption of photons and are excited into higher levels nj. In terms of the Drude model, the intraband one-photon excitation rate, W1pt, relates to the one photon absorption cross section σ1pt, laser intensity I, and photon energy (ħω) by [15,19,39]:
The stepwise gain of kinetic energy by the free electrons defines discrete energy levels εj = j × ħω in the CB. The critical energy, εcrit, at which impact ionization occurs requires a minimum number of IBA steps, . Here denotes the floor function. For ε > εcrit, the kinetic energy of CB electrons suffices to ionize a VB electron, which will result in two newly created CB electrons located at the bottom of CB. The impact ionization process occurs very fast; its rate, αimp, calculated for fused silica is about 1015 s−1 . No data are available for water, and we take the same value as for fused silica [18,40].
For discrete energy levels and ε0 = 0, the number density nj in each level can be fully described by the following set of equations :
When the individual rate equations are added up, we obtain an expression resembling a single-rate equation but with time dependent avalanche ionization rate:
The AI rate will grow during a laser pulse as long as nkcrit/ntotal increases but will approach a constant value when IBA is balanced by impact ionization. Rethfeld showed that the stationary regime with constant AI rate is reached after a transition time [18,40]18].
Equations (6) to (8) were derived under the assumption that the impact ionization rate is much larger than the one-photon absorption rate of CB electrons. Thus, there are two conditions for stationary avalanche ionization:
Numerical evaluation of Eq. (4) yields the temporal evolution of the electron number density in the discrete CB energy levels, i.e. the time-varying energy spectrum of CB electrons. The average energy of all CB electrons is given by:Eq. (11) is presented in the appendix.
Equations (10) and (11) corresponds to the impact ionization scheme of Fig. 1(a) and do not yet consider energy conservation upon impact ionization. Starting with Eq. (4) but based on the scheme of Fig. 1(b), we now derive a modified MRE and its approximate solution. Furthermore, we adapt the MRE to irradiation conditions relevant for material processing, in which and kcrit vary during the pulse and the valance band may be depleted.
In Fig. 1 (b), the residual energy remaining after collisional ionization is partitioned equally among the two CB electrons and the hole. This corresponds to a start energy εstart = (1/6) of the electrons. As mentioned before, the same start energy is also attributed to the electrons produced by SFI in order to obtain a discrete set of energy levels applicable to all CB electrons.
At large irradiance, the oscillation energy of CB electrons in the laser field becomes significant and varies during the laser pulse. Therefore, also the critical energy for impact ionization varies. It will be approximately equal to (3/2)Egap at the beginning and the end of the pulse but may significantly exceed this value at the peak of the pulse. Thus, the number of IBA events necessary to reach εcrit becomes a function of irradiance. Due to the finite start energy, this number is smaller than for the ionization scheme of Fig. 1(a). It is given byFig. 2.
To describe the occurrence of impact ionization when εcrit is exceeded, the Heaviside step function is used :
At large irradiance, the valance band may be depleted during the breakdown process. For water, the number density of bound electrons that can be ionized is nbound = 6.68 × 1023 cm−3, considering two 1b1 electrons per water molecule. During breakdown, the number density of valance band electrons still available for ionization is reduced to nval = nbound - ntotal. Therefore, a depletion factor, nval/nbound, is applied to both SFI and AI processes, following  and . Valence band depletion inhibits impact ionization, and the CB electrons may now reach energy levels above. The highest occupied level in a given breakdown process is represented by the symbol .
Based on the above considerations, the MRE with simple energy splitting scheme [denoted by MRE (1) in the rest of the paper] reads as follows:
The SRE corresponding to MRE (1) [denoted by SRE (1)] reads
Considering the finite start energy (1/6), the average kinetic energy of CB electrons computed from MRE (1) becomes
Note that in Eqs. (16) and (17) the kinetic energy includes the part of the electron energy, which is associated with its quiver motion in the electromagnetic field, i.e. its ponderomotive energy. The asymptotic average kinetic energy of CB electrons is a function of wavelength, as shown in Fig. 3. The fluctuations of the curve are due to the stepwise increase of k′crit with λ.
We can now assess the volumetric energy density U in the target that is composed of the energy deposited in electrons and holes. The kinetic energy of a hole produced by impact ionization is εh,imp = (1/6). Neglecting the energy of the holes produced by SFI, we obtain
4. Application range of SRE(1) in the (τL, λ) parameter space
SRE (1) for stationary avalanche ionization is applicable under the conditions summarized in Eq. (9). The ratio tasymp/τL is a measure for how fast the transient AI rate can reach its asymptotic value within the pulse duration, and W1pt/αimp measures how fast energy levels above εcrit are depleted by impact ionization. At very short pulse durations, when tasymp > τL, the free electron density is overestimated by the SRE. However, if the one-photon absorption rate exceeds the impact ionization rate, increases, and Eq. (17) for transiently underestimates the average free electron energy.
To study the validity range of the single-rate equation in the (τL, λ) parameter space, the normalized transition time tasymp/τL and the normalized one-photon excitation rate W1pt/αimp are plotted in Fig. 4 as a function of wavelength and pulse duration. The wavelength dependence of both parameters is shown for pulse durations ranging from 50 fs to 2 ns. In each case, the parameters values are given for the respective irradiance threshold for bubble formation, which for fs breakdown corresponds to a temperature of 441 K .
Generally, tasymp/τL decreases with increasing wavelength, whereas W1pt/αimp increases. This is because the one-photon absorption cross section σ1pt scales quadratically with wavelength [Eq. (3)]. Since W1pt = σ1pt × I/(ħω), this results in W1pt ∝ λ3 and tasymp ∝ λ−3. The steps appearing in the tasymp/τL curves are due to the change of k′crit with λ (see Fig. 2).
Many researchers argued that αimp is generally much larger than W1pt [6,15,18,27,43]. However, for a value αimp = 1015 s−1, we see that W1pt is close to or even larger than αimp for τL ≤ 250 fs and IR wavelengths.
The criterion tasymp < τL for stationary avalanche ionization is most easily violated under irradiation conditions featuring a slow, truncated ionization avalanche, i.e. for short wavelengths, short pulse durations, and at low irradiance. By contrast, the criterion W1pt < αimp is violated when the ionization avalanche becomes too fast, i.e. at high irradiance values, particularly at long wavelengths (since W1pt ∝ λ3). That may happen already under conditions in which the plasma density stays well below the critical electron density . Note that the data in Fig. 4 refer to the bubble threshold, which for fs pulses corresponds to nth = 1.8 × 1020 cm−3 . This is well below ncrit, which is ≈1021 cm−3 at 1050 nm.
Figure 5 compares simulation results for the time evolution of free electron density and average energy obtained by MRE (1) and SRE (1) for different values of τL, λ, and irradiance. We assumed a Gaussian temporal profile with peak irradiance I0. Both MRE (1) and SRE (1) were solved numerically using a Runge-Kutta method with adaptive step size control.
Predictions of SRE (1) and MRE (1) for the free electron density evolution, n(t), agree well for τL = 250 fs. Nevertheless, n(t) values are always slightly larger with the SRE because the finite time needed to reach a stationary AI rate is not considered. This problem can be fixed by using a slightly longer collision time for calculations with SRE(1) (τcoll = 1 fs) than with the MRE (τcoll = 0.9 fs) .
At τL = 50 fs, a larger deviation is observed, especially for the IR wavelength, Here the SRE (1) predictions for ntotal are, at the end of the pulse, 2.9 times higher than the results of the MRE (1) model. As seen in Fig. 4, the condition W1pt < αimp is not fulfilled, which leads to a slowing down of avalanche ionization. This is portrayed by MRE (1) but not by SRE (1) that neglects the time needed for impact ionization.
For all investigated parameters, SRE (1) and MRE (1) predictions for (t) and differ transiently during the laser pulse but at large irradiance values corresponding to the bubble threshold they converge towards the end of the laser pulse [Figs. 5(b) and 5(d)]. Here, transiently overshoots the average stationary value , partly because the quiver energy of the CB electrons is largest around the peak of the laser pulse [Eq. (16)], and partly because energy levels above are occupied when inverse Bremsstrahlung absorption is faster than impact ionization (W1pt > αimp). These effects are most pronounced for τL = 50 fs, where the irradiance is largest.
The kinks in the curve in Figs. 5(b) and 5(d) reflect changes of k′crit and εcrit during the laser pulse that are due to the irradiance dependence of the effective ionization potential [Fig. 2(b)]. The drop of the irradiance during the trailing edge of the pulse results in a stepwise decrease of k′crit and εcrit. At each step, large amounts of CB electrons can suddenly take part in the impact ionization process, resulting in an abrupt drop of . Finally, approaches Egap and converges against .
For τL = 250 fs and low irradiance values, such as used for nanosurgery by fs pulse series, the results of the SRE (1) and MRE (1) models do not converge [Fig. 5(f)]. Here, the values at the end of the laser pulse stay far below the values. SRE (1) is not applicable because at I0 = 1.0 × 1012 W/cm2 the transition to the stationary avalanche regime would take much longer than the laser pulse duration of 250 fs.
Overall, SRE (1) is useful for assessing the volumetric energy density of plasmas produced at large irradiance, i.e. in the context of laser material processing and surgery with single laser pulses. By contrast, MREs are needed to analyze molecular modifications by fs pulse series that rely on the electron energy available for molecular interactions.
5. Christensen and Balling’s MRE model applied to water
The simple energy splitting scheme of Fig. 1(b) is appropriate for conditions in which avalanche ionization dominates but for processes in which SFI plays a significant role, the scheme of Fig. 1(c) suggested by Christensen and Balling is more appropriate . As already mentioned, their scheme is compatible with the assumption that electrons excited via SFI are located at the bottom of CB possessing zero kinetic energy.
We saw from Fig. 5(b) that two or more energy levels may contribute to impact ionization, when kcrit drops at the trailing edge of a laser pulse. This is accounted for by considering several energy levels j′ with . While the simplified energy splitting scheme of Fig. 1(b) assumes a constant residual energy of (1/2), Christensen and Balling consider the exact values, . Now, the residual energy after impact ionization from the j′-th level that is carried by each of the two newly created CB electrons and the hole is:
While in the simple energy splitting scheme of MRE (1), a constant start energy (1/6) of the newly created electrons and hole is assumed, Christensen and Balling distributed the residual energy into the levels located around (1/6) such that the energy balance is optimized. In doing this, they take into account that electrons from more than one energy level with may contribute to impact ionization and that the residual energies for these levels may differ. This is achieved by a splitting function :
Like for MRE (1), the Heaviside step function [Eq. (13)] is used to express the probability for the occurrence of impact ionization at the j′-th level. Thus, the product of and describes both the probability for the occurrence of impact ionization for an electron located at the j′-th level and the allocation of electrons after impact ionization to the j-th level. That approach leads to the MRE :
The volumetric energy density, UMRE, and the part stored in holes, Uh,imp, are given byEq. (19)], and the valence band depletion.
Figure 6 shows simulation results of the MRE (2) model for the temporal evolution of energy deposition by a 100-fs, 1035 nm laser pulse at I0 = 2.0 × 1012 W/cm2. Such irradiance values are often used for nano-surgery with fs pulse series . Figures 6(a)-6(d) display the temporal evolution of electron density at each individual CB energy level, the total number density of electrons, their average kinetic energy, and finally the evolution of the focal temperature reached after thermalization.
Figure 6(a) shows the evolution of the number density of CB electrons at each individual energy level. The bottom level n0 is filled by SFI, and its occupation density thus follows the laser pulse shape. Seed electrons produced by SFI gain energy by IBA and successively occupy higher levels. The electron density is smallest at the level from where impact ionization occurs. Around the pulse peak, when irradiance is largest, kcrit increases from 12 to 13, and at t = 46 fs it drops from 13 to 12. Electrons from the well-occupied 12th level can now suddenly also contribute to impact ionization, and consequently n12 decreases sharply, accompanied by a sharp increase of n1 and n2. This results in a sudden increase of the total number density of CB electrons [Fig. 6(b)]. A slight drop of the average kinetic energy of CB electrons accompanies the loss of high-energy electrons from the 12th level [Fig. 6(c)].
The abrupt jump of occupation density associated with changes of kcrit does not reflect physical reality but is an artifact that arises from the assumption of discrete energy levels in MRE models. However, this artifact does not much affect integral parameters such as free electron density and total deposited energy, as seen in Figs. 6(b) and 6(d).
We see from Fig. 6(b) that at I0 = 2.0 × 1012 W/cm2 typical for nanosurgery by pulse series ≈70% of CB electrons is produced by SFI. The average kinetic energy of CB electrons approaches 6.6 eV at the end of the pulse [Fig. 6(c)]. The energy of holes is about 1.6 eV (Egap/6) at the beginning and end of the pulse, and reaches about 2.0 eV around the peak. The increased hole energy around t = 0 reflects the dependence of and kcrit on irradiance.
The increase of the focal temperature ΔT after thermalization is calculated by relating the energy density U carried by CB electrons and holes [Eq. (22)] to the mass density and heat capacity of water: ΔT = U/(ρwater × Cp). The rise of the thermodynamic temperature is very small (< 3 × 10−4 K) for the given laser parameters although the electron temperature (corresponding to the average kinetic energy ) is transiently large enough for the modification of biomolecules. Biomolecular changes induced by fs laser pulse series are, thus, chemical in nature and not caused by thermal effects .
6. Time evolution of CB electron energy spectra in low-density plasmas
Figure 7 shows the time evolution of CB electron energy spectra calculated with MRE (2) for the same laser parameters as in Fig. 6. At t = - 50 fs, a few seed electrons appear at the n0 level, and at t = - 25 fs, the bottom levels start to get populated. At the peak of the pulse (t = 0 fs), population of the n0 level is about to reach its maximum and the highest levels start to get occupied. At t = 25 fs, impact ionization begins to occur from the n13-th level. At t = 50 fs, impact ionization occurs from both n12 and n13 levels, and newly created electrons populate the n1 and n2 levels. From t = 50 fs to t = 100 fs, the energy spectrum assumes a nearly constant shape, which is characteristic for stationary avalanche ionization. Afterwards, dissipation mechanisms come into play. Thermalization of the electron energy through collisions and recombination occurs on a time scale of picoseconds. At the end of the laser pulse, most CB electrons have energies between 2 eV and 13 eV, with an average of 6.6 eV.
Figure 8 compares predictions of MRE (1) and MRE (2) for energy spectra at different laser wavelengths and irradiance values. At 1035 nm [Figs. 8(b) and 8(d)], the simulation results from both MREs are quite similar, except for the lowest energy levels. Here the assumption of a finite start energy (1/6) in the MRE (1) model results in a gap in the energy spectrum, whereas the MRE (2) model distributes the residual energy into two levels around (1/6) that belong to a ladder of energy levels starting at ε = 0. With both MREs, the level at ε = 0 is hardly occupied because the energy spectrum at λ = 1035 nm is characterized by avalanche ionization. That is different at λ = 515 nm [Figs. 8(a) and 8(c)], where multiphoton ionization plays a more prominent role and provides many electrons at the lower edge of the CB. This feature is correctly portrayed by MRE (2) but not by MRE (1). Thus, the MRE (2) model is the best choice for the calculation of energy spectra at short wavelengths, especially for low-density plasmas.
On the other hand, the MRE (1) model provides fairly accurate predictions for the energy spectra arising during the generation of high-density plasmas, when the stationary AI regime is reached. Here, the spectral gap in the range 0 < ε < (1/6) has little influence on the value that is of major interest for laser ablation. Thus, one can take the advantage of the smaller computational expense associated with the use of MRE (1) and SRE (1).
7. Irradiance dependence of total energy deposition and energy spectra
Figure 9 presents the irradiance dependence of volumetric energy density and temperature rise for 250 fs laser pulses at λ = 515 nm (a), and at λ = 1035 nm (b). It compares the results of MRE (1) model, MRE (2) model, and SRE (1). The energy density was calculated using Eqs. (18) and (22), respectively, and the temperature rise was obtained from ΔT = U/(ρwater × Cp). Irradiance values are normalized by a reference value of 1012 W/cm2, which serves as a benchmark for the onset of photomodifications by fs pulse series. Actual modification thresholds depend on wavelength and pulse duration  but the benchmark provides, nevertheless, some orientation in the laser parameter space. A second benchmark is defined through the threshold for bubble formation by a single laser pulse, i.e. the optical breakdown threshold. The respective irradiance values are marked by arrows in Figs. 9(a) and 9(b).
For both wavelengths, predictions of the two MRE models are very close to each other. This is because the volumetric energy density is largely determined by the effective band gap , and small differences in the assessment of have only a small influence on U. For the same reason, predictions of SRE (1) for U and ΔT are, at large irradiance, close to the predictions by the MRE models. Nevertheless, SRE (1) predictions are always slightly higher than those from MREs because the SRE assumes stationary avalanche ionization throughout the entire laser pulse. The largest relative difference between SRE and MREs is observed for λ = 515 nm and small irradiance, where it reaches 45.5% at I0 = 0.1 × 1012 W/cm2. However, as long as the conditions for stationary avalanche ionization are fulfilled [see Eq. (9) and Fig. 4], the SRE is a useful tool for assessing the total amount of deposited energy.
The irradiance dependence of average CB electron energy for different wavelengths and pulse durations is shown in Fig. 10. For these calculations, we used MRE (2) to obtain realistic values also for low irradiances, where SFI prevails and is small. With growing irradiance, increases continuously until it reaches a constant level in the stationary regime. At very large irradiance, the (I) curve exhibits a sharp increase when full ionization is reached and additional laser energy heats the CB electrons beyond the level in the stationary regime.
Since the AI rate increases with wavelength, the stationary regime is reached earliest for λ = 1035 nm and much later for λ = 350 nm. We also see that for shorter τL a higher irradiance is needed to reach the stationary AI regime. The dependence of (I) on wavelength and pulse duration offers interesting options for tuning the interactions between free electrons and biomolecules in aqueous media. Molecular changes with high activation energy can be avoided by using short wavelengths, or promoted via the use of IR laser pulses.
Figure 11 presents the irradiance dependence of the energy spectrum at the end of the laser pulse for 250-fs pulses of 1035 nm wavelength. The reference irradiance is Iref = 1012 W/cm2. Simulations are performed using MRE (1).
For I0 = 0.5 × Iref [Fig. 11(a)], electrons have already populated higher energy levels of the CB but the stationary avalanche regime featuring a constant shape of the post-pulse energy spectrum begins only at slightly higher irradiance. It is fully developed at I0 = Iref [Fig. 11(b)] and reaches up to I0 = 9.5 × Iref. The average energy at the end of the pulse amounts to = 6.8 eV in the entire stationary range [see insert in Fig. 11(b)], even though during the pulse (t) may transiently reach much higher values, as seen in Fig. 5. As mentioned above, the transient peak is due to the large ponderomotive electron energy at high irradiance, in conjunction with a fast IBA process exceeding the speed of impact ionization. The constant post-pulse value for electron densities between ≈1013 cm−3 and full ionization (6.7 × 1022 cm−3) can be explained by the fact that electrons with ε > εcrit quickly lose their kinetic energy through impact ionization of a valance band electron. This feature differs dramatically from laser-induced energy deposition into metals, where the free-electron density remains constant but the kinetic energy of carriers increases . By contrast, in dielectrics the electron energy distribution remains constant, whereas the number density increases by 9 orders of magnitude from I0 = Iref to I0 = 9.5 × Iref [see insert in Fig. 11(a)].
For I > 9.5 × Iref, full ionization is reached and extra laser energy will excite the CB electrons to higher energy levels [Figs. 11(d) and 11(e)]. This leads to a very rapid increase of . At I0 = 11 × Iref, it already reaches a value as high as 47 eV. Here, a division of CB into 70 discrete levels is required to calculate the spectrum, i.e. k′max = 70. It must be noted, however, that the validity of the model for such large irradiances is questionable because it does not consider the changes of optical plasma properties at n > ncrit  and alterations of the band structure in high-density plasma. The latter will firstly affect the initiation channel that relies on the existence of preexisting traps formed by relatively weak hydrogen bonds. Therefore, the MRE models and the asymptotic SRE model should be used only up to the onset of full ionization. Nevertheless, we may conclude that the electron temperature increases very rapidly when nonlinear energy deposition continues after full ionization, where a metal-like condition is reached.
8. Application fields for the different modeling approaches
We have seen in the previous section that the complex energy splitting approach provided by Christensen and Balling [MRE (2)] is needed to accurately portray the CB electron energy spectrum in fs laser generated low-density plasmas. By contrast, for large irradiance or long (ps and ns) laser pulse durations the simpler approach of MRE (1) provides sufficient accuracy. Moreover, for τL ≥ 250 fs, it enables to use a single-rate equation that reduces computation time. For a 250-fs, 1035-nm pulse at the bubble threshold, use of SRE (1) speeds up the computation by a factor of 8 compared to MRE (1), and by a factor of 8.4 compared to MRE (2).
The latter feature is of great interests for spatiotemporal simulations of energy deposition that combine modeling of nonlinear beam propagation and plasma formation. Spatiotemporal modeling of energy deposition plays a pivotal role in laser-material processing, since it relates the microscopic energy deposition process to macroscopic phase changes. This requires a spatial-temporal beam propagation model and a plasma formation model. Initially, simple geometric beam propagation models were combined with a single rate equation for plasma formation and estimates of the average kinetic energy of CB electrons [46,47]. More recently, researchers presented solutions of the full Maxwell equations in the bulk of dielectrics [20,48,49] or around a nanoparticle ; or of its approximations such as the nonlinear Schrödinger equation  or the vector Helmholtz equation . Nevertheless, most researchers still use a simple rate equation to model plasma formation [48,50,51], and few have incorporated an MRE model due to its complexity and large computational expense [27,43,49]. The SRE developed in this paper facilitates the combined modeling of beam propagation and plasma formation by providing accurate values of and volumetric energy density in the stationary avalanche regime. This will support the current trend towards multi-scale and multi-physics modeling of material modifications [20,50,51].
While determination of average kinetic energy is sufficient for the modeling of plasma-mediated phase transitions and ablation, detailed knowledge of the entire energy spectrum is needed for analyzing free-electron-mediated chemical modifications of biomolecules. Cumulative effects of this type may be relevant for photodamage in nonlinear microscopy [6,22], and are the mechanism of nanosurgery by fs pulse series [6,7,10].
Femtosecond laser induced populations of free electrons with known energy spectrum are an interesting tool for the investigation of e--mediated DNA damage and repair mechanisms in live cells [12,44,52]. DNA damage by low-energy electrons (LEEs) has already been studied in the context of radiation damage and radiation therapy [21,53,54]. In these investigations, mostly DNA films in vacuum were irradiated by an electron beam with known energy spectrum . However, investigations of the interaction of LEEs with DNA in aqueous environment under biological condition are still lacking. Recently, a new method for studying DNA damage and repairing mechanisms in live cells after laser micro-irradiation by tightly focused fs pulse series has been introduced [12,44,52]. An understanding of the damage mechanisms and their parameter dependence requires knowledge of both the number density and the energy distribution of laser-induced free electrons. MRE (2) provides information that will help deciphering the interaction mechanism of LEEs with biomolecules in radiation therapy as well as in nonlinear imaging under realistic biological conditions.
In this paper, we have studied the energy spectrum of laser-induced CB electrons in water by two multi-rate equation models that provide a wealth of information about the features and time course of nonlinear energy deposition. Both approaches consider the excess energy remaining after impact ionization in different ways, and are, thus, complementary. The approach with a simple energy splitting scheme [MRE (1)] provides sufficient accuracy for large irradiance and also for long (ps and ns) laser pulse durations. It enables the derivation of a single-rate equation [SRE (1)] that significantly reduces computation time and is valid for large irradiance and τL ≥ 250 fs. This will be useful for spatiotemporal simulations of energy deposition in dielectrics, which combine modeling of nonlinear beam propagation and plasma formation. Another approach based on a more complex energy splitting scheme [MRE (2)] is well suited for modeling energy spectra of CB electrons in low-density plasmas. It will be a valuable tool for analyzing free-electron-mediated chemical modifications of biomolecules. Potential application fields are the exploration of photodamage mechanisms in nonlinear microscopy and the investigation of free-electron mediated DNA damage and repair mechanisms in live cells. Moreover, understanding the dependence of electron energy spectra on laser pulse duration, wavelength, and irradiance opens pathways for inducing energy-specific molecular modifications.
The starting point for the derivation of in the stationary AI regime given in Eq. (11) is the electron number density at the critical level nkcrit and the total number density ntotal of CB electrons. They have been derived by Rethfeld for a rectangular pulse profile :
An expression for nj (t) can be derived in an iterative fashion: First, an analytical expression for n0(t) is obtained by inserting the above expression for nkcrit (t) into the rate equation for n0 in Eq. (4) and solving this linear ordinary differential equation. We then insert n0(t) into the rate equation for n1 and get the expression for n1(t). Continuing like this, we iteratively obtain a general expression for nj(t) with j = 0 to (kcrit −1)
The expression for the mean kinetic energy of CB electrons reads
Air Force Office of Scientific Research (FA9550-15-1-0326, FA9550-18-1-0521); National Natural Science Foundation of China (61335012, 61727823); China Scholarship Council.
We thank Elisa Ferrando-May, Michael Schmalz, and Norbert Linz for stimulating discussions.
1. R. R. Gattass and E. Mazur, “Femtosecond laser micromachining in transparent materials,” Nat. Photonics 2(4), 219–225 (2008). [CrossRef]
3. M. V. Shugaev, C. P. Wu, O. Armbruster, A. Naghilou, N. Brouwer, D. S. Ivanov, T. J. Y. Derrien, N. M. Bulgakova, W. Kautek, B. Rethfeld, and L. V. Zhigilei, “Fundamentals of ultrafast laser-material interaction,” MRS Bull. 41(12), 960–968 (2016). [CrossRef]
4. K. Sugioka and Y. Cheng, “Ultrafast lasers—reliable tools for advanced materials processing,” Light Sci. Appl. 3(4), e149 (2014). [CrossRef]
6. A. Vogel, J. Noack, G. Huettman, and G. Paltauf, “Mechanisms of femtosecond laser nanosurgery of cells and tissues,” Appl. Phys. B 81(8), 1015–1047 (2005). [CrossRef]
8. D. V. Palanker, M. S. Blumenkranz, D. Andersen, M. Wiltberger, G. Marcellino, P. Gooding, D. Angeley, G. Schuele, B. Woodley, M. Simoneau, N. J. Friedman, B. Seibel, J. Batlle, R. Feliz, J. Talamo, and W. Culbertson, “Femtosecond laser-assisted cataract surgery with integrated optical coherence tomography,” Sci. Transl. Med. 2(58), 58ra85 (2010). [CrossRef] [PubMed]
9. C. L. Hoy, O. Ferhanoglu, M. Yildirim, K. H. Kim, S. S. Karajanagi, K. M. C. Chan, J. B. Kobler, S. M. Zeitels, and A. Ben-Yakar, “Clinical ultrafast laser surgery: Recent advances and future directions,” IEEE J. Sel. Top. Quantum Electron. 20(2), 14 (2014). [CrossRef]
11. N. J. Friedman, D. V. Palanker, G. Schuele, D. Andersen, G. Marcellino, B. S. Seibel, J. Batlle, R. Feliz, J. H. Talamo, M. S. Blumenkranz, and W. W. Culbertson, “Femtosecond laser capsulotomy,” J. Cataract Refract. Surg. 37(7), 1189–1198 (2011). [CrossRef] [PubMed]
12. E. Ferrando-May, M. Tomas, P. Blumhardt, M. Stöckl, M. Fuchs, and A. Leitenstorfer, “Highlighting the DNA damage response with ultrashort laser pulses in the near infrared and kinetic modeling,” Front. Genet. 4, 135 (2013). [CrossRef] [PubMed]
13. J. S. D’Souza, J. A. Dharmadhikari, A. K. Dharmadhikari, B. J. Rao, and D. Mathur, “Effect of intense, ultrashort laser pulses on DNA plasmids in their native state: strand breakages induced by in situ electrons and radicals,” Phys. Rev. Lett. 106(11), 118101 (2011). [CrossRef] [PubMed]
14. C. Wang, M. Fomovsky, G. X. Miao, M. Zyablitskaya, and S. Vukelic, “Femtosecond laser crosslinking of the cornea for non-invasive vision correction,” Nat. Photonics 12(7), 416–422 (2018). [CrossRef]
15. P. K. Kennedy, “A first-order model for computation of laser-induced breakdown thresholds in ocular and aqueous media: Part I - Theory,” IEEE J. Quantum Electron. 31(12), 2241–2249 (1995). [CrossRef]
16. B. C. Stuart, M. D. Feit, S. Herman, A. M. Rubenchik, B. W. Shore, and M. D. Perry, “Nanosecond-to-femtosecond laser-induced breakdown in dielectrics,” Phys. Rev. B Condens. Matter 53(4), 1749–1761 (1996). [CrossRef] [PubMed]
17. J. Noack and A. Vogel, “Laser-induced plasma formation in water at nanosecond to femtosecond time scales: calculation of thresholds, absorption coefficients, and energy density,” IEEE J. Quantum Electron. 35(8), 1156–1167 (1999). [CrossRef]
19. N. Linz, S. Freidank, X. X. Liang, and A. Vogel, “Wavelength dependence of femtosecond laser-induced breakdown in water and implications for laser surgery,” Phys. Rev. B 94(2), 024113 (2016). [CrossRef]
20. N. M. Bulgakova, V. P. Zhukov, S. V. Sonina, and Y. P. Meshcheryakov, “Modification of transparent materials with ultrashort laser pulses: What is energetically and mechanically meaningful?” J. Appl. Phys. 118(23), 233108 (2015). [CrossRef]
21. B. Boudaïffa, P. Cloutier, D. Hunting, M. A. Huels, and L. Sanche, “Resonant formation of DNA strand breaks by low-energy (3 to 20 eV) electrons,” Science 287(5458), 1658–1660 (2000). [CrossRef] [PubMed]
22. D. Débarre, N. Olivier, W. Supatto, and E. Beaurepaire, “Mitigating phototoxicity during multiphoton microscopy of live Drosophila embryos in the 1.0-1.2 µm wavelength range,” PLoS One 9(8), e104250 (2014). [CrossRef] [PubMed]
23. D. Arnold and E. Cartier, “Theory of laser-induced free-electron heating and impact ionization in wide-band-gap solids,” Phys. Rev. B Condens. Matter 46(23), 15102–15115 (1992). [CrossRef] [PubMed]
24. A. Kaiser, B. Rethfeld, M. Vicanek, and G. Simon, “Microscopic processes in dielectrics under irradiation by subpicosecond laser pulses,” Phys. Rev. B Condens. Matter Mater. Phys. 61(17), 11437–11450 (2000). [CrossRef]
25. N. Brouwer and B. Rethfeld, “Excitation and relaxation dynamics in dielectrics irradiated by an intense ultrashort laser pulse,” J. Opt. Soc. Am. B 31(11), C28–C35 (2014). [CrossRef]
26. A. Ramer, O. Osmani, and B. Rethfeld, “Laser damage in silicon: Energy absorption, relaxation, and transport,” J. Appl. Phys. 116(5), 053508 (2014). [CrossRef]
27. B. H. Christensen and P. Balling, “Modeling ultrashort-pulse laser ablation of dielectric materials,” Phys. Rev. B Condens. Matter Mater. Phys. 79(15), 155424 (2009). [CrossRef]
28. N. Linz, S. Freidank, X. X. Liang, H. Vogelmann, T. Trickl, and A. Vogel, “Wavelength dependence of nanosecond infrared laser-induced breakdown in water: Evidence for multiphoton initiation via an intermediate state,” Phys. Rev. B Condens. Matter Mater. Phys. 91(13), 134114 (2015). [CrossRef]
29. M. U. Sander, M. S. Gudiksen, K. Luther, and J. Troe, “Liquid water ionization: mechanistic implications of the H/D isotope effect in the geminate recombination of hydrated electrons,” Chem. Phys. 258(2-3), 257–265 (2000). [CrossRef]
30. C. G. Elles, A. E. Jailaubekov, R. A. Crowell, and S. E. Bradforth, “Excitation-energy dependence of the mechanism for two-photon ionization of liquid H(2)O and D(2)O from 8.3 to 12.4 eV,” J. Chem. Phys. 125(4), 44515 (2006). [CrossRef] [PubMed]
31. C. G. Elles, I. A. Shkrob, R. A. Crowell, and S. E. Bradforth, “Excited state dynamics of liquid water: insight from the dissociation reaction following two-photon excitation,” J. Chem. Phys. 126(16), 164503 (2007). [CrossRef] [PubMed]
32. R. A. Crowell and D. M. Bartels, “Multiphoton ionization of liquid water with 3.0-5.0 eV photons,” J. Phys. Chem. 100(45), 17940–17949 (1996). [CrossRef]
33. L. V. Keldysh, “Ionization in the field of a strong electromagnetic wave,” Sov. Phys. JETP 20, 1307–1314 (1965).
34. B. K. Ridley, Quantum processes in semiconductors (Clarendon, 1999).
35. B. Rethfeld, O. Brenk, N. Medvedev, H. Krutsch, and D. H. H. Hoffmann, “Interaction of dielectrics with femtosecond laser pulses: Application of kinetic approach and multiple rate equation,” Appl. Phys., A Mater. Sci. Process. 101(1), 19–25 (2010). [CrossRef]
36. K. Sokolowski-Tinten and D. von der Linde, “Generation of dense electron-hole plasmas in silicon,” Phys. Rev. B Condens. Matter Mater. Phys. 61(4), 2643–2650 (2000). [CrossRef]
37. Q. Sun, H. Jiang, Y. Liu, Z. Wu, H. Yang, and Q. Gong, “Measurement of the collision time of dense electronic plasma induced by a femtosecond laser in fused silica,” Opt. Lett. 30(3), 320–322 (2005). [CrossRef] [PubMed]
38. K. Waedegaard, M. Frislev, and P. Balling, “Femtosecond laser excitation of dielectric materials: experiments and modeling of optical properties and ablation depths,” Appl. Phys., A Mater. Sci. Process. 110(3), 601–605 (2013). [CrossRef]
39. A. Ramer, L. Haahr-Lillevang, B. Rethfeld, and P. Balling, “Modeling the transient optical parameters in laser-excited band gap materials,” Opt. Eng. 56(1), 011015 (2016). [CrossRef]
40. B. Rethfeld, “Free-electron generation in laser-irradiated dielectrics,” Phys. Rev. B Condens. Matter Mater. Phys. 73(3), 035101 (2006). [CrossRef]
41. A. Vogel, N. Linz, S. Freidank, and G. Paltauf, “Femtosecond-laser-induced nanocavitation in water: implications for optical breakdown threshold and cell surgery,” Phys. Rev. Lett. 100(3), 038102 (2008). [CrossRef] [PubMed]
43. C. L. Arnold, A. Heisterkamp, W. Ertmer, and H. Lubatschowski, “Computational model for nonlinear plasma formation in high NA micromachining of transparent materials and biological cells,” Opt. Express 15(16), 10303–10317 (2007). [CrossRef] [PubMed]
44. M. F. Schmalz, I. Wieser, F. Schindler, C. Czada, A. Leitenstorfer, and E. Ferrando-May, “Highly standardized multicolor femtosecond fiber system for selective microphotomanipulation of deoxyribonucleic acid and chromatin,” Opt. Lett. 43(12), 2877–2880 (2018). [CrossRef] [PubMed]
45. L. Jiang and H. L. Tsai, “Improved two-temperature model and its application in ultrashort laser heating of metal films,” J. Heat Transfer 127(10), 1167–1173 (2005). [CrossRef]
46. J. Jiao and Z. Guo, “Modeling of ultrashort pulsed laser ablation in water and biological tissues in cylindrical coordinates,” Appl. Phys. B 103(1), 195–205 (2011). [CrossRef]
47. J. Zhou, J. K. Chen, and Y. Zhang, “Numerical modeling of transient progression of plasma formation in biological tissues induced by short laser pulses,” Appl. Phys. B 90(1), 141–148 (2008). [CrossRef]
48. A. Rudenko, J. P. Colombier, and T. E. Itina, “Graphics processing unit-based solution of nonlinear Maxwell’s equations for inhomogeneous dispersive media,” Int. J. Numer. Model 31(2), e2215 (2018). [CrossRef]
49. A. Rudenko, J. P. Colombier, and T. E. Itina, “From random inhomogeneities to periodic nanostructures induced in bulk silica by ultrashort laser,” Phys. Rev. B 93(7), 075427 (2016). [CrossRef]
50. A. Dagallier, E. Boulais, C. Boutopoulos, R. Lachaine, and M. Meunier, “Multiscale modeling of plasmonic enhanced energy transfer and cavitation around laser-excited nanoparticles,” Nanoscale 9(9), 3023–3032 (2017). [CrossRef] [PubMed]
51. V. Jukna, A. Jarnac, C. Milián, Y. Brelet, J. Carbonnel, Y. B. André, R. Guillermin, J. P. Sessarego, D. Fattaccioli, A. Mysyrowicz, A. Couairon, and A. Houard, “Underwater acoustic wave generation by filamentation of terawatt ultrashort laser pulses,” Phys. Rev. E 93(6), 063106 (2016). [CrossRef] [PubMed]
52. D. Träutlein, M. Deibler, A. Leitenstorfer, and E. Ferrando-May, “Specific local induction of DNA strand breaks by infrared multi-photon absorption,” Nucleic Acids Res. 38(3), e14 (2010). [CrossRef] [PubMed]
53. E. Alizadeh, T. M. Orlando, and L. Sanche, “Biomolecular damage induced by ionizing radiation: the direct and indirect effects of low-energy electrons on DNA,” Annu. Rev. Phys. Chem. 66(1), 379–398 (2015). [CrossRef] [PubMed]
54. L. Sanche, “Interaction of low energy electrons with DNA: Applications to cancer radiation therapy,” Radiat. Phys. Chem. 128, 36–43 (2016). [CrossRef]
55. L. Sanche, “Low energy electron-driven damage in biomolecules,” Eur. Phys. J. D 35(2), 367–390 (2005). [CrossRef]