Analytical optimization of the cutting efficiency for generic cavitation bubbles

: A theoretical method to determine the optimum laser parameters for maximizing the cutting efficiency for different materials (in particular human cornea) is proposed. The model is simple and reduced to laser beam characteristics and cavitation properties. The model further provides a method to convert energy fluctuations during the cutting process to equivalent deviations in the cavitation bubbles. The proposed model can be used for calibration, verification and validation purposes of laser systems used for cutting processes at relatively low cost and may improve the quality of the results.


Introduction
The available methods that allow for the correction of refractive defects have been expanded with the advent of ultrashort pulses laser systems [1]. Intrastromal dissection is necessary for femtosecond laser-assisted refractive surgery, in particular for the creation of a flap for the laser assisted in-situ keratomiluesis, or for the creation of a intrastromal lenticule of corneal tissue for lenticule extractions techniques [2]. It evolved from picosecond laser systems to ones operating in the femtosecond regime [3], with recent attempts on sub-nanosecond regimes [4], but the underlying mechanism remains the same [5]: tightly focusing an ultrashort laser pulse, nonlinear absorption processes within the focal volume, generation of dense free electron plasma. If a critical value is exceeded, laser induced optical breakdown occurs with very fast increase of temperature and pressure, which leads to a rapid plasma expansion. This results in a shock wave propagating into the surrounding medium, causing the formation of a cavitation bubble, which may undergo a series of oscillations before ending in a small persistent gas bubble after some microseconds. The temporally and spatially scanned process eventually leads to a confluent (or coalescent) bubble (conforming a cleavage plane), enabling the separation of the tissue. The high temperature and pressure of the gas in the bubble then leads to the bubble oscillation. However, this picture is valid for water and may not be transferred completely to corneal tissue, in which during bubble expansion, work must be done against the restoring forces of the lamellar structure of the cornea. This has strong influences on the intrastromal bubble dynamics, eventually leading to a smaller bubble size (and less oscillations) in cornea than in water for the same pulse energy.
As presented in several research works, the (asymptotic) size of the bubble grows with the cubic root of the applied (suprathreshold) pulse energy [6][7][8]. For energies closer to the bubblethreshold, bubble size rapidly increases for slightly increased energies. Achieving accurate clinical outcomes and reducing the likelihood of a retreatment procedure are major goals of refractive surgery. For that, accurately calibrated lasers with high stability are required.
The understanding of the underlying processes may help improving systems. It is well known that a successful surgery depends on the appropriate design of the corrective profile, precise delivery of laser energy to the corneal position, and reliable understanding of the corneal tissue response. A large number of factors influence the outcome. Among them, laser energy delivery technique [9], ablation decentration and registration [10], eye tracking [11], wound-healing [12], and biomechanics of the cornea [13], have been explored to predict or explain the clinically observed discrepancy between intended and actual outcomes. The quantification of influence of these factors is important for providing the optimal outcome for the refractive surgeries.
With the introduction of the femtosecond laser technologies for refractive surgery, corneal cuts are more accurate and precise than ever before [14]. The three-dimensional geometry of the cutting pattern must be converted into a series of pulse positions, often requiring several million pulses for the surgery. The procedure is nowadays a successful technique, due to its micrometric precision and the high predictability and repeatability of the process accompanied by minimal side effects [15]. Although the procedure and techniques are quite mature and distributed, still the topic "Optimization of the cutting efficiency" is worth to be analyzed and considered, because its clinical implications are not yet deeply explored. The real impact of resolution in laser corneal refractive surgery is still discussed in a controversial way. There may be essential differences between the coalescence of one single, small bubble with another single, small bubble and the coalescence of one single, small bubble with a residual bubble layer which has already been accumulated by thousands of bubbles before. Every new single, small bubble contributes to the cutting of the cornea, but they may be governed by different underlying mechanisms. For instance, very different regimes have been observed for the interaction mechanisms of cavitation bubbles induced by spatially and temporally separated fs-laser pulses (at least for pulse energies well above the breakdown threshold) [16][17][18][19].
In a previous work [20] we presented an analytical model to optimize the laser settings of ablation processes governed by linear absorption. For this work, we have taken similar expressions (and a conceptually similar approach) adapted to processes governed by non-linear (multi-photon) absorption. The main purpose of this paper is to present a simple theoretical framework (in the form of a reduced model) to determine the optimum laser parameters for maximizing the cutting efficiency. The latter is defined here as the completion of a corneal cut ensuring the generation of coalescent bubbles, while imparting the minimum amount of energy to the tissue (i.e. minimum dose).

Methods
The interaction of near infrared femtosecond laser radiation and corneal tissue is a complex process [21]. The generation of cavitation bubbles for a given wavelength is mainly governed by the relationship between the duration of the laser pulse and the numerical aperture of the optical system (i.e. the available photon density per space and time in the focal volume) [22,23].
Corneal cutting is essentially a form of micro-machining. Although the single pulse energy is low, given the small size of the focal diameter (beam waist) when focusing with moderately-high numerical apertures (NA), energy density can be high. Additionally, given the short pulse duration, the peak power provided can be high.
Many parameters have to be considered in designing an efficient laser cutting process. One of which is the selection of the appropriate wavelength (typically in the 1000-1100 nm window for corneal tissue) with very low linear absorption in the corneal tissue while tightly focused, which results in a high-energy deposition in a small volume for a speedy and complete plasma mediated ablation [24]. A second important parameter is a short pulse duration τ to maximize peak power and minimize thermal conductivity to the adjacent tissue (femtosecond based τ<1 ps).
It has been reported that for surface ablation, the plateau-like region observed between 100 fs and 1 ps for the corneal layers indicates that for use in laser surgery, laser pulse durations chosen within this range should be practically equivalent [25]. Further to that, previous works show that while pulse duration plays an important role in the determination of breakdown thresholds, the threshold grows slower than the statistically expected square root of the pulse duration [26].
The affected area of a single spot is much smaller than the total treatment zone. Due to this, multiple laser pulses are sequentially delivered on to the cornea. Each laser pulse locally creates a small cavitation bubble separating the corneal tissue. The global process is an integral effect of the local process of each individual laser pulse. A larger bubble gives a faster treatment, but a lower resolution. On the other hand, a smaller bubble increases the resolution at the cost of increasing the treatment time due to increased number of laser pulses invested in cutting the same area. Additionally, using laser pulses close to the threshold of the material (for generating cavitation bubbles) would mean only imparting thermal effects instead of generating stable bubbles. All these factors make the energy selection a sensitive criterion. While selecting the optimum energy, a delicate balance needs to be accounted for, between the material thresholds, resolution of the laser pulses, thermal and dose effects of the material and the total time required. In the following steps, based on the size of the cavitation bubble, we calculate the local and overall dose per treatment, further optimizing this dose for a minimum energy per area while also accounting for fluctuations in the laser tissue interaction from fluctuations potentially appearing in laser output energy.

Calculation of the bubble radius per laser pulse
The bubble radius R B can be calculated (in good agreement with the experience [27]) as shown in [28]: where E P is the energy of a single pulse, E Th is the threshold for the formation of cavitation bubbles for the irradiated tissue or material below which no bubbles occur, and K is a coefficient of the irradiated tissue or material. For human corneal tissue irradiated with ultrashort near infrared pulses, the threshold for bubble formation takes values in the order of nanojoules [26].
Despite of its simplicity, the method remains general. Different materials would have thresholds, and specific laser beam characteristics would lead to a particular threshold. Parameters like pulse duration, tissue properties are considered using previous works. In our simulations, we can determine the threshold as being proportional to with the corresponding pulse duration τ, the wavelength λ, the beam quality M 2 and the numerical aperture NA.

Calculation of the local dose per pulse
The following equation applies for the calculation of the local dose per pulse, i.e. per cavitation bubble D P (A B represents the cross-sectional area of the cavitation bubble): (3)

Calculation of the overall dose per treatment
The following equation applies for the calculation of the overall dose per treatment D T : where E T is the total energy deposited during a treatment, A T is the area of the treatment and R T its radius. The number n represents the amount of pulses assumed that all pulses carry the same energy, and pulses are approximately evenly distributed in the treatment area. The total number of pulses (for a single pass over the treatment area) can be calculated as: where ∆ x and ∆ y are the distances between deposited neighboring spots along two essentially perpendicular axes (e.g. x and y but also r and θ). Therewith, the overall dose simplifies to Since usually, the goal is to reduce the number of remaining tissue bridges, the interspot distance ∆ is usually (theoretically or empirically) related to the bubble radius: where F is a scaling factor, i.e. dimensionless, and G is an additive value i.e. with length dimensions (in particular µm). In general, the overlap factor will be F ≥ 1 and the reduction in distance G ≥ 0 (i.e. the interspot distance is reduced compared to the bubble size). So that the overall dose takes the form:

Calculation of the linear dose per treatment
Another engineering simplification implies setting one of the interspot distances (∆ 1 ) very small (setting very dense pulse trains in a given direction) and independent of the applied energy, i.e. this interspot distance is smaller than the smallest bubble for the considered range of pulse energies. Therefore, the linear dose D L can be defined as follows:

Optimization process
Optimization involves minimizing the energy per area, i.e. the dose D in order to minimize the thermal effects on the tissue as well as the effects on the neighboring laser spots, while maintaining the delicate balance between thermal effects and the resolution of the treatment. With regard to the optimum energy value E opt , the minimization can be expressed as The value of E P calculated from this expression is the optimum pulse energy for the considered metric. Following the optimization process, in the following expressions, the optimum pulse energy for the considered metric, is represented by E opt , and is calculated for the minima of the considered metrics.

Determination of the fluctuations in laser-tissue interaction from the fluctuations observed in the laser energy output
The deviations in laser-tissue interaction derived from the fluctuations observed in the laser energy output (∆E P ) with regard to a certain metric M can be equated as:  . At threshold energy, bubble radius collapses to zero, they all approach the asymptotic value for suprathreshold pulse energies with different rates.

Optimum local dose per pulse
Minimizing Eq. (3) following Eq. (10), the optimum energy can be calculated as which leads to optimum local dose per laser pulse D opt P when inserted again in Eq. (3):

Optimum overall dose per treatment
The optimum overall dose per treatment can be calculated through the following minimization: while the differential as derived from Eq. (8) is not a very comfortable equation to find the root, yet a number of simplifications may be applied.
Since the cross-sectional area of the bubbles is approximately circular, one can set: so that Since the size of the bubbles increases for larger pulse energies, it is more common (and logical) to set the interspot distance in relative relation (ratio) to the bubble size, so that G can be set to zero.
Analogous to the optimum local dose, the optimum pulse energy and therewith the optimum overall dose per laser pulse can be calculated as E opt P = 3E Th (see Eq. (12)) and respectively.

Optimum linear dose per treatment
The corresponding calculations for the optimum linear dose per laser pulse lead to and respectively. Figure 2 represents the corresponding dose for each metric (and the minimum dose corresponds to the optimum energy). Figure 3 aims to show the impact of F and G in the determination of the optimum. The optimum for each metric is clearly represented with a minimum valley in the curve. Before the optimum energy, the dose decreases (with maximum cutting efficiency at minimum dose). Beyond the optimum value the dose starts to increase (representing a decline in cutting efficiency). The rate of increment in efficiency below optimum is observed to be higher compared to the rate of decrement in efficiency beyond optimum.

Determination of the fluctuations in laser-tissue interaction from the fluctuations observed in the laser energy output
The deviations in dose ∆D derived from the fluctuations observed in the laser energy ∆E P output can be equated as but we have defined the optimum as the energy for which So that in the surroundings of the optimum pulse energy, fluctuations in the pulse energy have little impact in the applied dose (local, overall, or linear). And as seen in Figs. 2 and 3, the fluctuations are higher for pulse energies below optimum than for pulse energies above the optimum. However, in reality, it may be a non-linear effect: Not only the dose is affected (directly) by potential energy fluctuations but also the spacing (indirectly). In a real laser system, internal structures may control the spacing in a way, that changes in energies lead to corrected spacings. The deviations in bubble radius ∆R B derived from the fluctuations observed in the laser energy output can be equated as or, in relative terms  (19)). From the fluctuations observed in energy delivery, the deviations in bubble radius can be estimated, which is shown in Fig. 4 for the case E Th = 150 nJ.
The same dependency as for the bubble radius R B can similarly be calculated for the cross- and or, in relative terms For a given energy fluctuation, as the energy increases, the relative bubble size deviation decreases. This was also predicted by the smaller slope for high-above threshold energies in Fig. 1. The relative deviation is the corresponding relative change in the metric value in percentage. For energies beyond optimum energy, deviations of the bubble size are smaller than energy fluctuations. The effects of the energy fluctuation are magnified at energies lower than the optimum energy.

Discussion
This study provides analytical expressions for calculation of some optimum energy parameters for maximizing cutting efficiency and for determination of the fluctuations in corneal cuts from the fluctuations observed in energy delivery. Separate analysis of the effect of each parameter was performed. Several metrics relevant to laser refractive surgery have been proposed and optimized (local dose D P , overall dose D T , linear dose D L ). In a sense, these are the metrics characterizing the safety of the irradiation.
The optimization of intrastromal dissection first requires an understanding of the underlying mechanisms, thus the simplistic modeling approach presented in this work may not cover all the details required to provide a rigorous analytical solution. Yet, since the cornea is composed of ≈ 85% water, we believe that the model may still hold as a valid approximation. Human cornea is a viscoelastic material and the dynamics of cavitation bubbles is quite different from the case of water. However, even in the case of viscoelastic materials Eq. (1) holds [29,30]. Both references seem to support our simplistic model. In a viscoelastic material and for large values of the laser pulse energies, the shock wave emitted during optical breakdown has a bipolar component (a compressive wave followed by a tensile wave); but this is likely beyond the range 3× Threshold determined with our simplistic approach. In that regime of high pulse energies (probably clearly in the photodisruption regime), the tensile wave may produce negative collateral effects which are not included in the present model and, thus, the proposed model may have limitations for energies exceeding 10× the threshold.
The differences between bubble dynamics in water and cornea may be actually more significant than assumed in this work. Even though the intrastromal collagen tissue has ≈ 85% water content, the lamellar structure of the tissue leads to much more complex bubble dynamics than stated here, with correspondingly large implications for the cutting process. However, this manuscript does not aim to investigate or model intrastromal corneal dissection, but rather to present a simple model for dose optimizations. This complexity may limit the extent this approach can really contribute to the optimization of cutting parameters for ultrashort pulse laser systems.
At the core of our calculation is the cubic root relationship between energy (exceeding the threshold) and bubble size suggested by Vogel and associates [31]. We present a very simple model based on the above mentioned and well-known formula of bubble radius R B depending on laser pulse energy E P . We consider the simplicity of the model as one of the beauties of this work. All the optimization problem can be reduced to the value of the energy threshold to generate cavitation bubbles.
The size of the cavitation bubble increases with pulse energy. In our simulations, we considered K as constant coefficient depending only on the tissue. This is in agreement with the findings by Vogel et al. [32]. This can be interpreted in the sense, that the initial amount of energy (subthreshold) is invested to satisfy the condition to generate the plasma; and only (a fraction of) the energy above threshold (the amount of energy exceeding the threshold) can be invested in generating the shock wave which ultimately leads to the cavitation bubble. So that the threshold is the only characterizing part of the bubble dimensions, i.e. all curves of bubble size vs energy asymptotically trend to the same curve (the one for a hypothetical zero-threshold). Hence, for very high energies, the bubble size is irrespective of laser parameters such as wavelength, numerical aperture, or pulse duration (or even M 2 ).
Even in the case that K may not be an absolute constant, and may depend upon wavelength, numerical aperture, repetition rate, or pulse duration (or even M 2 ), our model remains valid as long as K is not a function of pulse energy. It may be argued that Eq. (1) is valid only for energies far above the threshold value for bubble formation, and that closer to the threshold value the coefficient K (i.e. in a sense the conversion efficiency) may be a function of the laser pulse energy. Previous works [29] actually, seem to support our simplistic model. For energies far above the threshold for bubble formation (the threshold becomes negligible compared to pulse energy), Eq. (1) takes the form of E  3 . This can be explained without the need for a dynamic behavior of K or the conversion efficiency close to the threshold. It suffices that K (i.e. in a sense the conversion efficiency) affects (E − E Th ) 1 3 instead of E 1 3 . Actually, taking 12% as the asymptotic non-linear conversion efficiency η nl (as estimated from [29] for a viscoelastic material with 85% water), the overall ratio from bubble energy to deposited pulse energy equates: where α nl is the non-linear absorption coefficient. This leads to the following values: perfectly matching the fact that the conversion efficiency of incident laser energy into bubble energy strongly depends on the nonlinear absorption process. At threshold it is virtually 0 and increases to an approximately constant value of 10%-15% at energies larger than ten times the threshold energy. Therefore, the proposed simple model is capable of reasonably reproducing the conversion efficiency for the intermediate energy range, which is absolutely relevant for refractive surgery.
The understanding of the underlying processes may help improving systems. Our model aims to a simple approach for describing the global effects, without the need to enter the details of the physical processes. We suggest a simple way to optimize the settings to drive the basic physical process for such structure formation within specific range of parameters (mainly fs, near-IR pulses).
In essence, laser focus, cavitation bubbles and the disruption process are three dimensional processes. All the used calculations, however, were done in two dimensions since the aim of generation of cleavage planes was approached as a dose optimization process, i.e. energy deposited per unit of area.
In this work, it is assumed that the equatorial plane of the bubble has a circular shape, whereas the axial plane will generally represent a vertically elongated ellipse; thus the bubble can be represented as a vertically elongated ellipsoid of revolution/spheroid.
The presented results are valid as long as effects of the repetition rate of the laser remain negligible. The repetition rate could play a crucial role. If the lifetime of the cavitation bubble (or the remaining gas bubble) is longer than the inverse of the pulse repetition rate (i.e. the interpulse time), the subsequent pulse may target the cavitation and will be distracted so that no breakdown can occur. Such effects can be inferred from previous works on the interaction mechanisms of cavitation bubbles induced by spatially and temporally separated fs-laser pulses (at least for pulse energies well above the breakdown threshold) [16][17][18][19]. This effect (as far as we can judge, not deeply characterized yet) may/would open a new degree of freedom in the optimization process.
The interplay between bubble size, spot spacing and repetition rate is very complex. In the work by Juhasz et al. [5] photographs of the bubble dynamics in bovine cornea showed a bubble lifetime of up to 1 ms. The maximum bubble size of 45 µm diameter is too large for today's application in the refractive surgery field. More recent work by Tinne et. al. [19] showed that bubbles in water and 1% gelatin have lifetimes above 10 µs (for energies above 1 µJ and relative energies above 6× the threshold). In the same work Tinne et al. [19] found in water that below 3.6× the threshold (≈ 500 nJ in their setup) no influence between the cavitation bubble oscillations was observed, and below 7.2× the threshold (≈ 1000 nJ in their setup) an asymmetric bubble collapse can be observed. For the 2% and 5% gelatin an asymmetric bubble collapse was observable during the interaction oscillation for energies as high as 8.4× the threshold (≈ 1200 nJ in their setup) due to the decrease in maximum bubble radius and lifetime with increasing gelatin concentration [19], whereas the temporal overlap parameter was neglected for there was no influence on the interaction effects.
In fact, the repetition rate is crucial for intrastromal dissection in refractive surgery. The bubble remains within the cornea for a few microseconds (for low pulse energies slightly exceeding the breakdown threshold, only one expansion and collapse of the bubble occurs). So that below some 10 kHz Repetition rate would probably not play a relevant role, and beyond some 100 kHz Repetition rate would probably play a crucial role and interactions with the bubble of the previous laser pulse have to be considered, which is additionally dependent on spot spacing and bubble size/laser energy. Systems on the market tend to operate at repetition rates above 100 kHz, and would be "exactly at the border" in the several 100 kHz range. We acknowledge this as a limitation of the current model.
In addition, closely neighboring pulses can strongly influence each other: A subsequent pulse may have to interact with material optically changed by the previous pulse, thus causing varying conditions for the laser-tissue interaction [16,19,33].
Also, as reported in previous works [7], we only consider here the threshold for bubble formation. This threshold usually is lower than the threshold for dielectric breakdown [34], and lies below the actual photodisruption regime (high density plasma), maybe in the region of low density plasma [35] (i.e. plasma mediated ablation [36]).
We acknowledge that for the typical surgeon/user of a system determining the threshold empirically may be difficult (or actually impossible given the limitations imposed by the manufacturers of commercial systems), yet for researchers and manufacturers this information may well be available. Further to that, by applying comprehensive physical modeling [26], simple expressions may be derived to determine (or at least reasonably estimate) the breakdown threshold for specific laser configurations.
Pulse energy plays the most important role in the determination of bubble size. The range for energies of the femtosecond systems for refractive surgery available in the market runs from about 50 nJ to about 1 µJ, whereas the threshold for generating bubbles runs from about 30 nJ to about 300 nJ.
Equation (4) shall be understood as the energy reaching the focus (i.e. after reducing the component of the linear absorption) and not the energy delivered by the laser (or the energy reaching the corneal surface), but linear absorption in this regime is so weak, that this reduction may represent only about 1% of the pulse energy. This 99% transmission (in general constant for a constant depth) would not (or only minorly) change the obtained results (but may marginally affect 3-dimensional cuts across different depths in the cornea).
The calculation of the local dose per pulse by calculating the laser pulse energy per bubble cross-sectional area does not correspond to the actual dose imparted by the individual laser pulses in focus. From a purely physical point of view, the energy of the laser pulse must be normalized to the focal area. We defined the local dose per laser pulse as the energy per area of the cavitation bubble, and not per area of the spot size. Spot size will not increase with energy, so dose per pulse over spot area would actually be linear with energy, so that trivially zero energy or E P = E Th would lead to minimum spot dose. The pulse energy required for the optimum local dose per laser pulse was found to be 3× the threshold. The approach of taking the laser pulse energy per bubble cross-sectional area as the "effective" dose of the individual pulses is practical from a modeling point of view (especially because as the pulse energy increases, the increase in "effective" dose is attenuated by the increasing bubble size), but physically close to the focus, the increase in pulse energy is associated with an actual higher local dose per pulse.
For the overall dose per treatment we used two approaches, one in which two perpendicular spot distances both depending on pulse energy are used, and another one in which one spot distance is fixed at a very small value, and only the spot distance in the perpendicular direction depends on pulse energy. The former may represent well the conditions for pulse-to-pulse scanning systems using reticular-grids as scan patterns (e.g. Johnson & Johnson Vision iFS, or Alcon FS200), whereas the latter may represent the conditions for pulse-train scanning systems using fast and slow tracks as scan patterns (e.g. Ziemer Z8).
In order to generate coalescent bubbles (cleavage plane), the interspot distances have to be set smaller than the bubble diameter, and in the pulse-train scanning systems the fixed interspot distance should be smaller than the smallest bubble for the considered range of pulse energies (and for circular bubble areas the factor √ 2 would ensure tissue-bridge free dissection). This implies that the linear dose (actually a misnomer, but to reflect the fact that only the interspot distance along one direction is varied) is always higher than the overall treatment dose, and both higher than the local dose per laser pulse.
In general, the interspot distance is set as a function of the pulse energy (ideally through Eq. (7)). For ∆ = 2 · R B , adjacent bubbles would be "just touching" in the two principal directions; but free, uncut areas (tissue bridges) will remain on the diagonals. For ∆ = 2·R B √ 2 , adjacent bubbles would strongly overlap in the two principal directions; but just touching along the diagonals. Thus, for ∆ = 2·R B 2 0.25 , adjacent bubbles would slightly overlap in the two principal directions and close to overlap in the two diagonal directions (leaving minimum tissue bridges). For that reason, we consider the values F = 2 0.25 and G = 0.1 µm (since the bubble are a few µm in size) as a reasonable compromise. This is depicted in Fig. 5 for a cartesian arrangement. , adjacent bubbles strongly overlap in the two principal directions; and are tangential along the diagonals. In this ideal case, no uncut areas are left and a tissue-bridge-free dissection can be expected. C: For ∆ = 2·R B 2 0.25 , adjacent bubbles slightly overlap in the two principal directions and are close to overlap in the two diagonal directions (leaving minimum tissue bridges between the 4 bubbles). The size of the diamond-like area is 1 4 of the size in A. For that reason, we consider the values F = 2 0.25 and G = 0.1 µm (since the bubble are a few µm in size, not depicted) as a reasonable compromise.
The pulse energy required for the optimum overall dose per treatment was found to be also 3× the threshold, whereas the pulse energy required for the optimum linear dose per treatment resulted in 1.5× the threshold. These optimum values do not change for different overlap factors (F), but shift to higher optimum energies for absolute reductions in distance (G).
Linking this model with experimental data or cutting efficiency data from the systems in the market would strengthen its value. However, those systems offer the freedom to select spacing and pulse energy with a relatively high degree of independence. First determining the theoretical optimum, and then "strongly" modifying the settings in all directions (increasing energy shall not (largely) improve separation, while reducing energy shall be clearly detrimental; decreasing spacing shall not (largely) improve separation, while increasing spacing shall be clearly detrimental) may help at this aim. Considering experiments on some ex-vivo model to further stress the efficiency curve, while qualifying the cutting quality and the easiness of dissection remains a subjective matter and may mislead the clean model. It may actually be interpreted that the determination of threshold shall be improved.
Perhaps the strongest criticism of this work is the lack of comparison with experimental results. Despite the absence of verification or confirmation by experimental methods, recent reports [37][38][39] work in energy regimes of less than 2× the threshold associated with tight/dense spacings and repetition rates in the some MHz range, well below the 10 µs interpulse time.
For pulse energies above 1.5× the threshold, energy fluctuations are attenuated in smaller relative fluctuations in bubble radius, whereas for energies above 3× the threshold, energy fluctuations are attenuated in smaller relative fluctuations in bubble area. This optimum is irrespective of the actual energy fluctuations.
The effect of energy deviation in bubble size accounts for 100% (bubble or no bubble) close to the threshold. This implies that working at energies above and beyond the optimum value can reduce the effects of energy fluctuations on the outcomes by reducing the deviations in bubble size. Observing the fluctuations in the energy delivery, the deviations in the size of the cavitation bubble can be estimated. For a given energy fluctuation, as the energy increases, the relative deviation in the size of the cavitation bubble decreases. Beyond the optimum energy, the effect of energy fluctuations on the size of the bubble is lesser compared to energies below optimum. This has also been recently reported for other types of femtosecond micromachining [40].
There is another way of interpreting optimum efficiency in terms of rate of change as shown in Fig. 4. The rate of change of the selected metric with respect to the energy can be evaluated. At the optimum fluence for a given metric, the rate of change vs. energy equals unity (so relative energy fluctuations correspond 1:1 to relative deviations in bubble size). This can be explained with the effect of energy fluctuations on different metrics with respect to the energy (Fig. 4). Below optimum, relative energy fluctuations are "amplified" in larger deviations in the bubble size, i.e. investing a bit more energy, drastically increases bubble size. Beyond optimum, relative energy fluctuations are "attenuated" in smaller relative deviations in the size of the bubble, i.e. to slightly increase bubble size, much more energy shall be invested. Additionally, the rate of increment in efficiency below optimum is observed to be higher compared to the rate of decrement in efficiency beyond optimum (Fig. 2). In Fig. 4, for explanatory purposes, the energy fluctuation was set at 1% to 9% level, corresponding to a realistic value, which helps to understand the key concept of amplifying versus attenuating the energy fluctuation.
For the selected valuable metric to be optimized (slightly different versions of dose) the final determined value shall be equal or higher than the theoretical optimum in order to have a safety margin for achieving the highest cutting efficiency for this energy, but attenuated relative deviations.
The calculation of the final optimum energy depends mainly on two factors. First, the total number and kind of metrics being considered for optimization of the cutting efficiency and second the energy fluctuations of the system. Therefore, the final energy can be calculated in two steps in order to account for the aforementioned factors. In the first step the optimum energy for the considered metrics is calculated. In case one considers multiple metrics an optimum energy between lowest and highest values amongst the considered metrics should be chosen. In the second step the calculated optimum fluence can be added with a factor that accounts for the fluctuations in energy. Adding this factor shall provide additional stability from the energy fluctuations of the system. For femtosecond laser systems, energy deviation typically lies below ±3% Root-Mean-Square (RMS). Considering the RMS as one standard deviation, we arrive at around ±5% energy deviation for a typical femtosecond laser system. Hence, we propose to add +5% to the analytically calculated optimum energy. This factor will ensure that the energy level will remain optimum even if energy instability is introduced because of the system itself.
We consider dose (in different variants) as most valuable metric for two reasons, one to reduce the total amount of energy delivered and second since it is the base for the determination of the maximum permissible exposure as described in the standards for laser safety [41]. The cutting resolution is also made more stable for the optimum criteria as the bubble size saturates quickly after optimum. Furthermore, the energy optimized for minimum dose shall reduce the thermal effects induced in the tissue, as well as the incidence of clinical complications as e.g. opaque bubble layers [42].
Being in the energy regime above optimum (i.e. for which the relative energy fluctuations are "attenuated" in smaller relative deviations in the size of the bubble) also sets the energy sufficiently above threshold as to reduce the incidence of clinical complications as e.g. black spots [43].
As a conclusive example, we would choose an optimum energy of which lies between 3/2× and 3× the threshold (theoretical optima for dose and for attenuating energy fluctuations) +5% (from typical energy fluctuations of femtosecond laser systems). This model can be easily generalized to any material for which the bubble coefficient and the threshold for bubble formation for the specific wavelength and laser characteristics are known. With that, it may complement previous analytical approaches to the efficiency problem and may sustain the observations reported by others.
Clinical evaluations on human eyes are needed to confirm the preliminary simulated results presented herein. Accurate knowledge for the coefficient of the human cornea is imperative in validating the presented results.

Conclusions
The model introduced in this study provides an analytical expression for optimizing cutting parameters for ultrashort pulse laser systems. Furthermore, due to its analytical approach, it is valid for different laser devices used in refractive surgery, as well as for any materials for which the coefficients and the bubble characteristics for the specific laser system are known.
The development of more accurate models to improve emmetropization and the correction of ocular aberrations is an important issue. We hope that this model will be an interesting and useful contribution to refractive surgery and will take us one-step closer to this goal.