1. Introduction

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 $\overline{\epsilon }$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 [22].

The time evolution of electron energy spectra was first analyzed for wide-bandgap solids such as SiO₂ [23,24], and sophisticated simulations for SiO₂ 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 [18] 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 [27]. 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 E_gap 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, I_th, for nanosecond optical breakdown revealed steps in the I_th(λ) curve that are consistent with breakdown initiation by multiphoton ionization, with an initiation energy of about 6.6 eV [28]. 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 [28].

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, E_gap, or via an intermediate level, E_solv. 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 E_solv ≈6.6 eV. Excited electrons relax into the traps and solvate. The traps have a number density of about χ_trap ≈10¹⁹ cm⁻³, corresponding to ≈3 traps in 1000 water molecules [28]. The solvated electrons, e_aq^-, are then easily upconverted into the CB due to the fact that the remaining energy gap amounts to only 3 eV and that e_aq^- absorbs very strongly from ultraviolet (UV) to infrared (IR). Thus, the critical step for breakdown initiation is multiphoton excitation into the intermediate level E_solv. Subsequent upconversion into the CB is assumed to follow immediately and is not considered as a separate term of the rate equation [19,28].

 

Fig. 1 Ionization schemes for laser energy deposition in water with different impact ionization schemes. They include Rethfeld’s scheme with zero start energy after impact ionization (a), a simplified scheme with constant start energy of 1/6 of the effective bandgap $\tilde{\Delta }$ (b), and Christensen and Balling’s scheme, in which the residual energy remaining after impact ionization is split into different energy levels that depend on photon energy (c). Red arrows mark transitions involving inverse Bremsstrahlung absorption (IBA) and impact ionization. The interband energy level E_solv of water is relevant only for strong field excitation process, which are marked by blue arrows.

Download Full Size | PDF

The total energy required to excite electrons from the valence band into the conduction band is given by the effective ionization potential$\tilde{\Delta }$that accounts both for the band gap and the oscillation energy of free electrons in strong electromagnetic fields [15,19,33]:

The SFI rate considering excitation through E_solv and via E_gap is

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$\tilde{\Delta }$. For a parabolic bandgap, the minimum required energy is (3/2)$\tilde{\Delta }$ [6,18,24,34].

Characteristic differences between the ionization schemes consist in the way in which the residual energy (1/2)$\tilde{\Delta }$ 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 [18]. 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 [18] supported such a scenario. However, they were obtained with a one-photon absorption rate, W_1pt, for IBA that considered only electron-phonon scattering but neglected electron-ion collisions [18,24]. That W_1pt value corresponds to τ_coll = 14.3 fs in the Drude model, as has been shown in Ref [19]. In later work, Rethfeld and associates included electron-ion collisions and assumed a collision time in the order of 1 fs [35], in accordance with experimental results on silicon and fused silica [36,37]. Our recent findings for water corroborate this value [19]. 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)$\tilde{\Delta }$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)$\tilde{\Delta }$.

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, $\hslash \omega$. 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)$\tilde{\Delta }$ either from the laser irradiation (the energy needed for multiphoton ionization will usually exceed $\tilde{\Delta }$), 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)$\tilde{\Delta }$. 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)$\tilde{\Delta }$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 n₀ energy level with ε₀ = 0. Electrons then gain kinetic energy via inverse Bremsstrahlung absorption of photons and are excited into higher levels n_j. In terms of the Drude model, the intraband one-photon excitation rate, W_1pt, 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, $k_{\text{crit}}=\lfloor \frac{\left(3/2\right)\tilde{\Delta }}{\hslash \omega }+1\rfloor$. Here $\lfloor \rfloor$ 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 10¹⁵ s⁻¹ [18]. No data are available for water, and we take the same value as for fused silica [18,40].

For discrete energy levels and ε₀ = 0, the number density n_j in each level can be fully described by the following set of equations [18]:

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 n_kcrit/n_total 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]

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:

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 $\tilde{\Delta }$ and k_crit 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)$\tilde{\Delta }$ 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 $\tilde{\Delta }$ varies during the laser pulse. Therefore, also the critical energy for impact ionization varies. It will be approximately equal to (3/2)E_gap 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 by

 

Fig. 2 Temporal evolution of the effective bandgap and k′_crit for fs pulses of 1050 nm with Gaussian temporal shape at different irradiance values.

Download Full Size | PDF

To describe the occurrence of impact ionization when ε_crit is exceeded, the Heaviside step function is used [27]:

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 n_bound = 6.68 × 10²³ cm⁻³, considering two 1b¹ electrons per water molecule. During breakdown, the number density of valance band electrons still available for ionization is reduced to n_val = n_bound - n_total. Therefore, a depletion factor, n_val/n_bound, is applied to both SFI and AI processes, following [6] and [27]. Valence band depletion inhibits impact ionization, and the CB electrons may now reach energy levels above$k_{\text{crit}}^{\text{'}}$. The highest occupied level in a given breakdown process is represented by the symbol $k_{\max }^{\text{'}}$.

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)$\tilde{\Delta }$, 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 ${\overline{\epsilon }}_{\text{asymp}}\left(\lambda \right)$curve are due to the stepwise increase of k′_crit with λ.

 

Fig. 3 Asymptotic average kinetic energy of CB electrons as a function of wavelength, calculated from Eq. (17).

Download Full Size | PDF

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)$\tilde{\Delta }$. 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 t_asymp/τ_L is a measure for how fast the transient AI rate can reach its asymptotic value within the pulse duration, and W_1pt/α_imp measures how fast energy levels above ε_crit are depleted by impact ionization. At very short pulse durations, when t_asymp > τ_L, the free electron density is overestimated by the SRE. However, if the one-photon absorption rate exceeds the impact ionization rate, $\overline{\epsilon }$ increases, and Eq. (17) for ${\overline{\epsilon }}_{\text{asymp}}$ 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 t_asymp/τ_L and the normalized one-photon excitation rate W_1pt/α_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 [41].

 

Fig. 4 Normalized transition time t_asymp/τ_L (solid lines) and normalized one-photon excitation rate W_1pt/α_imp (dashed lines) at bubble formation thresholds, plotted as a function of wavelength and pulse duration. The respective peak irradiance values at threshold are: I₀ = 2.90 × 10¹³ W/cm² for τ_L = 50 fs; I₀ = 8.3 × 10¹² W/cm² for τ_L = 250 fs; I₀ = 1.63 × 10¹² W/cm² for τ_L = 3 ps; I₀ = 6.9 × 10¹¹ W/cm² for τ_L = 30 ps; and I₀ = 3.8 × 10¹¹ W/cm² for τ_L = 2 ns. These I₀ values refer to the average bubble threshold in the displayed wavelength range. The value for τ_L = 50 fs is calculated by the MRE (1) model, the value for τ_L = 250 fs is taken from [19], and the value for τ_L = 2 ns is taken from [28]. The others are taken from Fig. 5 in [42]. The collision time is set as 1.0 fs.

Download Full Size | PDF

Generally, t_asymp/τ_L decreases with increasing wavelength, whereas W_1pt/α_imp increases. This is because the one-photon absorption cross section σ_1pt scales quadratically with wavelength [Eq. (3)]. Since W_1pt = σ_1pt × I/(ħω), this results in W_1pt ∝ λ³ and t_asymp ∝ λ⁻³. The steps appearing in the t_asymp/τ_L curves are due to the change of k′_crit with λ (see Fig. 2).

Many researchers argued that α_imp is generally much larger than W_1pt [6,15,18,27,43]. However, for a value α_imp = 10¹⁵ s⁻¹, we see that W_1pt is close to or even larger than α_imp for τ_L ≤ 250 fs and IR wavelengths.

The criterion t_asymp < τ_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 W_1pt < α_imp is violated when the ionization avalanche becomes too fast, i.e. at high irradiance values, particularly at long wavelengths (since W_1pt ∝ λ³). That may happen already under conditions in which the plasma density stays well below the critical electron density $n_{\text{crit}}={\omega }^2{m}_c{\epsilon }_0/e^2$. Note that the data in Fig. 4 refer to the bubble threshold, which for fs pulses corresponds to n_th = 1.8 × 10²⁰ cm⁻³ [19]. This is well below n_crit, which is ≈10²¹ cm⁻³ 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 I₀. Both MRE (1) and SRE (1) were solved numerically using a Runge-Kutta method with adaptive step size control.

 

Fig. 5 Temporal evolution of free electron density (left column) and average kinetic energy of CB electrons (right column) at various laser parameters simulated by the MRE (1) model (solid lines) and by the SRE (1) model (dashed lines). First and second rows present results for irradiation close to the bubble threshold with pulse duration with 50 fs and 250 fs. The respective peak irradiance values are I₀ = 2.9 × 10¹³ W/cm² for τ_L = 50 fs, and I₀ = 8.3 × 10¹² W/cm² for τ_L = 250 fs. The same I₀ values are used for all wavelengths (350 nm, 500 nm and 1050 nm) although the bubble thresholds slight vary with wavelength. The third row contains results for laser pulse τ_L = 250 fs at a lower irradiance values I₀ = 1.0 × 10¹² W/cm², which is slightly above the threshold for nanosurgery by fs pulse series [6,44]. In all calculations τ_coll = 1.0 fs.

Download Full Size | PDF

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) [19].

At τ_L = 50 fs, a larger deviation is observed, especially for the IR wavelength, Here the SRE (1) predictions for n_total 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 W_1pt < α_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 $\overline{\epsilon }$(t) and ${\overline{\epsilon }}_{\text{asymp}}$ 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, $\overline{\epsilon }$ transiently overshoots the average stationary value ${\overline{\epsilon }}_{\text{asymp}}$, 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 $k_{\text{crit}}^{\text{'}}$are occupied when inverse Bremsstrahlung absorption is faster than impact ionization (W_1pt > α_imp). These effects are most pronounced for τ_L = 50 fs, where the irradiance is largest.

The kinks in the $\overline{\epsilon }\left(t\right)$ 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 $\tilde{\Delta }$ [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 $\overline{\epsilon }$. Finally, $\tilde{\Delta }$ approaches E_gap and $\overline{\epsilon }$ converges against ${\overline{\epsilon }}_{\text{asymp}}$.

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 $\overline{\epsilon }$ values at the end of the laser pulse stay far below the ${\overline{\epsilon }}_{\text{asymp}}$ values. SRE (1) is not applicable because at I₀ = 1.0 × 10¹² W/cm² 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 [27]. 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 k_crit drops at the trailing edge of a laser pulse. This is accounted for by considering several energy levels j′ with ${\epsilon }_{j\text{'}}\ge {\epsilon }_{\text{crit}}$. While the simplified energy splitting scheme of Fig. 1(b) assumes a constant residual energy of (1/2)$\tilde{\Delta }$, Christensen and Balling consider the exact values, ${\epsilon }_{j\text{'}}-\tilde{\Delta }$. 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)$\tilde{\Delta }$ of the newly created electrons and hole is assumed, Christensen and Balling distributed the residual energy into the levels located around (1/6)$\tilde{\Delta }$ such that the energy balance is optimized. In doing this, they take into account that electrons from more than one energy level with ${\epsilon }_{j\text{'}}\ge {\epsilon }_{\text{crit}}$ may contribute to impact ionization and that the residual energies for these levels may differ. This is achieved by a splitting function [27]:

Like for MRE (1), the Heaviside step function $\Theta$ [Eq. (13)] is used to express the probability for the occurrence of impact ionization at the j′-th level. Thus, the product of $\Theta$ and $\Upsilon$ 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 [27]:

The volumetric energy density, U_MRE, and the part stored in holes, U_h,imp, are given by

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 I₀ = 2.0 × 10¹² W/cm². Such irradiance values are often used for nano-surgery with fs pulse series [6]. 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.

 

Fig. 6 Time evolution of energy deposition simulated by the MRE (2) model. (a) Change of number density of CB electrons at each individual energy level as a function of time; (b) evolution of total number density of CB electrons, n_total, and of the contribution coming from SFI, n_SFI; (c) evolution of the mean kinetic energy of CB electrons and holes; (d) evolution of energy deposition expressed by the resulting temperature rise. Besides the total temperature rise, ΔT_total, also the contributions from electrons and holes, ΔT_e- and ΔT_h are displayed. Simulation parameters are I₀ = 2.0 × 10¹² W/cm², τ_L = 100 fs, λ = 1035 nm, and τ_coll = 0.9 fs.

Download Full Size | PDF

Figure 6(a) shows the evolution of the number density of CB electrons at each individual energy level. The bottom level n₀ 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, k_crit 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 n₁₂ decreases sharply, accompanied by a sharp increase of n₁ and n₂. 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 k_crit 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 I₀ = 2.0 × 10¹² W/cm² 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 (E_gap/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 $\tilde{\Delta }$ and k_crit 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 × C_p). The rise of the thermodynamic temperature is very small (< 3 × 10⁻⁴ K) for the given laser parameters although the electron temperature (corresponding to the average kinetic energy $\overline{\epsilon }$) 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].

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 n₀ 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 n₀ 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 n₁₃-th level. At t = 50 fs, impact ionization occurs from both n₁₂ and n₁₃ levels, and newly created electrons populate the n₁ and n₂ 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.

 

Fig. 7 Time evolution of CB electron energy spectra calculated using the MRE (2) model for water. Simu-lation parameters are the same as in Fig. 5: I₀ = 2.0 × 10¹² W/cm², τ_L = 100 fs, λ = 1035 nm, and τ_coll = 0.9 fs.

Download Full Size | PDF

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)$\tilde{\Delta }$ 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)$\tilde{\Delta }$ 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.

 

Fig. 8 Comparison of electron energy spectra at the end of a 100-fs pulse predicted by MRE (1) (dashed lines) and MRE (2) (solid lines) for two different wavelengths (515 nm and 1035 nm) and irradiance levels. Irradiance values in the top row are typical for nanosurgery by pulse series, and values in the bottom row represent the bubble threshold.

Download Full Size | PDF

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)$\tilde{\Delta }$ has little influence on the $\overline{\epsilon }$ 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 × C_p). Irradiance values are normalized by a reference value of 10¹² W/cm², which serves as a benchmark for the onset of photomodifications by fs pulse series. Actual modification thresholds depend on wavelength and pulse duration [6] 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).

 

Fig. 9 Irradiance dependence of energy deposition simulated by the MRE (1) model, the SRE (1) model, and the MRE (2) model for 250 fs laser pulses at λ = 515 nm (a) and at λ = 1035 nm (b). The respective irradiance ranges are chosen such that full ionization is not reached because otherwise the computational efforts becomes too large for calculations using MRE(2). The collision time τ_coll is 1.0 fs for the SRE model and 0.9 fs for the MRE models. The deposited energy is represented in temperature rise ΔT on the left axis and in volumetric energy density U on the right axis. Scaling of both axes is related by U = 4.2 × 10⁻³ (kJ cm⁻³ K⁻¹) × ΔT. Arrows indicate the bubble formation threshold.

Download Full Size | PDF

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 $\tilde{\Delta }$, and small differences in the assessment of $\overline{\epsilon }$ 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 I₀ = 0.1 × 10¹² W/cm². 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 $\overline{\epsilon }$ values also for low irradiances, where SFI prevails and $\overline{\epsilon }$ is small. With growing irradiance,_{$\overline{\epsilon }$} increases continuously until it reaches a constant level in the stationary regime. At very large irradiance, the $\overline{\epsilon }$(I) curve exhibits a sharp increase when full ionization is reached and additional laser energy heats the CB electrons beyond the $\overline{\epsilon }$level in the stationary regime.

 

Fig. 10 Average kinetic energy at the end of the pulse as a function of peak irradiance I₀, calculated by MRE (2) for 100 fs and 250 fs pulses at 350 nm (a), 515 nm (b), and 1035 nm (c). The sudden jump of $\overline{\epsilon }$for 250 fs pulses indicates that full ionization is reached; see Fig. 11.

Download Full Size | PDF

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 $\overline{\epsilon }$(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 I_ref = 10¹² W/cm². Simulations are performed using MRE (1).

 

Fig. 11 Irradiance dependence of electron energy spectra at the end of the pulse simulated by the MRE (1) model. Spectra are shown for: I₀ = 0.5 × I_ref (a); I₀ = 1.0 × I_ref (b); I₀ = 9.5 × I_ref (c); I₀ = 10 × I_ref (d); I₀ = 11 × I_ref (e), with I_ref = 1.0 × 10¹² W/cm². Simulations were performed for laser pulse of τ_L = 250 fs and λ = 1035 nm, using a collision time τ_coll = 0.9 fs. The insert in (a) shows the irradiance dependence of total number density of CB electrons, while the insert in (b) shows the irradiance dependence of the average kinetic energy of CB electrons at the end of the pulse, calculated using Eq. (16).

Download Full Size | PDF

For I₀ = 0.5 × I_ref [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 I₀ = I_ref [Fig. 11(b)] and reaches up to I₀ = 9.5 × I_ref. The average energy at the end of the pulse amounts to $\overline{\epsilon }$ = 6.8 eV in the entire stationary range [see insert in Fig. 11(b)], even though during the pulse $\overline{\epsilon }$(t) may transiently reach much higher values, as seen in Fig. 5. As mentioned above, the transient $\overline{\epsilon }$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 $\overline{\epsilon }$value for electron densities between ≈10¹³ cm⁻³ and full ionization (6.7 × 10²² cm⁻³) 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 [45]. By contrast, in dielectrics the electron energy distribution remains constant, whereas the number density increases by 9 orders of magnitude from I₀ = I_ref to I₀ = 9.5 × I_ref [see insert in Fig. 11(a)].

For I > 9.5 × I_ref, 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 $\overline{\epsilon }$. At I₀ = 11 × I_ref, 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 > n_crit [2] 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 [50]; or of its approximations such as the nonlinear Schrödinger equation [51] or the vector Helmholtz equation [43]. 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 ${\overline{\epsilon }}_{\text{asymp}}$ 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 [55]. 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.

9. Conclusions

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.