A temperature model for laser lithotripsy

Objective To derive and validate a mathematical model to predict laser-induced temperature changes in a kidney during kidney stone treatment. Methods A simplified mathematical model to predict temperature change in the kidney for any given renal volume, irrigation flow rate, irrigation fluid temperature, and laser power was derived. We validated our model with matched in vitro experiments. Results Excellent agreement between the mathematical model predictions and laboratory data was obtained. Conclusion The model obviates the need for repeated experimental validation. The model predicts scenarios where risk of renal tissue damage is high. With real-time knowledge of flow rate, irrigating fluid temperature and laser usage, safety warning levels could be predicted. Meanwhile, clinicians should be aware of the potential risk from thermal injury and take measures to reduce the risk, such as using room temperature irrigation fluid and judicious laser use.


Introduction
Over the last 30 years, Holmium lasers have been used to fragment stones within the urinary tract. Initially, lowpower lasers (20 Watts) were used to fragment the stones into pieces for manual basket extraction. Relatively few high-energy impulses were delivered. Over time, lasers for stone surgery have become more powerful (up to 120 Watts) to provide the ability to modulate the frequency and pulse energy across a greater range. This has made the technique of "dusting" feasible. In this technique, the laser is used at a higher frequency with a lower energy to gradually break off tiny pieces of stone that can be passed spontaneously in the urine and do not require time-consuming extraction. Other surgical techniques such as "pop-corning" rely on firing the laser at higher power-typically around 40 Watts-for a prolonged period of time to agitate the stone fragments in the kidney. As the fragments move around, they fragment into smaller pieces when they come in contact with the laser. In the last year, the new thulium-doped fibre laser has been launched for clinical use and offers the potential for delivering even higher power levels. As lasers have evolved, the power (energy/second) delivered can be increased. An unintended consequence is the risk of thermal tissue damage due to heating of the irrigation fluid within the ureter or renal pelvis.
The potential for thermal tissue damage during holmium laser lithotripsy, and its dependence on procedural parameters, is currently under investigation. Recent in vitro [1,2] and in vivo [3] experimental work has aimed to measure fluid temperatures resulting from holmium laser activation. These studies considered different laser settings-energy and pulse frequency-as well as irrigation flow rate; irrigation is the continuous delivery of saline solution to clear stone fragments resulting from lithotripsy [1][2][3]. All measured high fluid temperatures-particularly with high laser power and low irrigation flow-in their experimental setups and voiced a concern for thermal tissue damage as indicated by the commonly used t 43 metric [4].
Mathematical modelling provides a platform for predicting temperatures within the renal pelvis, and subsequent risk of thermal tissue damage, without the need to experimentally test every clinical condition. A computational model to predict the distribution of temperature over time and space due to holmium laser lithotripsy was developed in a commercial finite-element software, COMSOL Multiphysics, by [3]. The model includes equations for heat transfer through the fluid, the solid boundary, and couples temperature to fluid flow through temperature-dependent fluid properties. Simulations were performed to mimic their in vivo and in vitro experimental setups, and the calculated volume-averaged temperatures were in good agreement with their experimental results. However, implementation of this model requires some knowledge of COMSOL Multiphysics, and a simulation for a single set of parameter values requires a few minutes of computational time.
Therefore, it is of interest to investigate whether a significantly reduced mathematical model can still capture the necessary physics to provide an accurate estimate of volume-averaged fluid temperature as a function of relevant clinical parameters. To motivate our model, we conduct a series of in vitro experiments, similar to those performed by [1], although with different irrigation fluid temperatures and a larger experimental vessel. This provides two sets of experimental data against which to validate our mathematical model. The model itself comprises a single analytical expression, which can be evaluated with no computational expense for any clinical parameters on a standard calculator. We determine good agreement between our model and experimental data, along with the experimental data of [1]. We then discuss the potential for thermal tissue damage as predicted by the validated mathematical model, and predict clinically safe ranges of laser settings and irrigation flow rate such that fluid temperatures remain below a critical threshold.

Experimental setup
Wet-lab experiments consisted of a 38.3-mL cylindrical vessel, submerged in a 1-L container, both filled with room temperature saline. A Boston Scientific LithoVue™ ureteroscope was positioned with its tip flush with an 11/13 F hole at the top of the container. A 365-μm laser fibre (Flexiva™ Boston Scientific) was inserted through the ureteroscope and secured 10 mm distal to the ureteroscope end ( Fig. 1). Two T-type theromocouple wires were positioned proximal and distal to the laser (red circle and triangle, respectively, in Fig. 1) and measured temperature at 0.1 s intervals. Irrigation was switched as we commenced temperature recording. After 20 s, the laser was switched on for a total of 60 s. After the laser was switched off, irrigation was maintained and temperatures recorded for a further 20 s.
We considered six sets of experimental settings-two laser powers and three flow rates-summarised in Table 1. We performed five runs of each configuration.
Our experiments were similar to those performed by [1], although with different values for the vessel volume and material properties, initial fluid temperature, and irrigation temperature. The volume of our experimental vessel (38.3 mL) is representative of a typical kidney, whereas the volume of the experimental vessel used by [1] is more indicative of an isolated calyx (5.9 mL). We will consider both sets of experiments in this article to compare to our modelling predictions. We refer to experiments performed by [1] as Set A and our experiments as Set B. All experiment parameters are summarised in Table 2. Adjusting for the position of the laser fibre with respect to the scope had no notable effect on our findings (see Appendix B).

Theoretical model
A mathematical model for spatially averaged temperature as a function of time T(t), in a fluid-filled vessel of volume V, subject to irrigation at flow-rate Q and temperature T in and laser lithotripsy at power W, is derived in Appendix C. The model comprises an analytic formula predicting an exponential temperature rise during lasering, followed by an exponential decay when the laser is switched off. If the laser is switched on for sufficient time, a steady-state temperature is reached. The volume of the vessel sets only the time required to reach the steady-state temperature, and not its magnitude. The model takes, as input parameters, V, Q, W, T in , along with the initial temperature in the vessel T(0) = T 0 , the time the laser is switched on t 1 , the time the laser is switched off t 2 , and k, , and c which are the thermal conductivity, density, and specific heat capacity of the irrigation fluid, respectively. The final model parameter is = hs , where h is the heat transfer coefficient through the walls of the vessel and s is the cross-sectional area of the walls through which heat diffuses. As h is highly dependent on the material properties of the vessel under consideration, we treat this as an unknown parameter, and determine it through a fit to experimental data. The material properties of the vessel in Set A and Set B are potentially different, so we determine one value of for each data set.
We define the laser firing time as t f = t 2 − t 1 and the temperature change as

Fit for W
e fit for , which characterises the thermal conductivity through the walls of the experimental vessel. For both sets of experiments A and B, we fit for using the case where Q = 0 mL min −1 and the laser is switched on with power W = 40 Watts. The experimental data, along with the lines of best fit, are plotted in Fig. 2a, b for Set A and Set B, respectively.
The best-fit values of are found to be for the experiments from Set A and Set B, respectively. (1)

Results
Using the best-fit values for , Eq. (1), we have the required information to compare the experimental data from both sets to the predictions of the theoretical model.

Comparison with experimental data
In Fig. 3, we present the data extracted from Set A (symbols) along with model predictions (solid lines). Experimental parameters used in the model are given in Table 2.
In Fig. 3a, we consider fixed W = 40 Watts, and vary Q. We note that although was only obtained through a fit to the exponential rise for Q = 0 mL min −1 it accurately predicts temperatures for Q = 7.5 , 14.5, and 39 mL min −1 . In Fig. 3b, we fix Q = 14.5 mL min −1 and vary laser power W. The data and model prediction for W = 40 are repeated from the red data set and model prediction in Fig. 3a. Again, we observe good agreement between theory and experiment. We speculate that discrepancies between theory and data in Fig. 3b may be due to variability in flow rates which are given in [1] as ranges: Q = 0 , 7 − 8 , 14 − 15 , 38 − 40 mL min −1 -estimated from a measured driving pressure head.
In Fig. 4, we present experimental results from Set B (symbols) compared to the predictions of the mathematical model (solid lines). Experimental parameters are given in Table 2.  indicate the two different temperature probes with corresponding symbols in Fig. 1. We observe good agreement between theory and experiment in both setups. The overlapping triangle and circle data points in Fig. 4 validates the use of a spatially averaged mathematical model by indicating a fairly homogeneous temperature throughout the experimental vessel.
It is illuminating to compare Figs. 3a to 4b, which are both for W = 40 Watts, and observe the comparatively larger temperature increases observed in the experimental data of Set A. This is due to the volume disparity between the two experiment sets. The mathematical model described in Appendix C predicts that the time required to reach the steady-state temperature is inversely proportional to the volume of the vessel. Thus, within a smaller vessel, the temperature equilibrates more quickly, and thus higher temperatures are achieved for the same firing time than in a larger vessel. The inverse relationship between rise time and vessel volume is captured by the shallower slopes of the temperature curves in Fig. 4 when compared to the curves in Fig. 3.

Effects of volume, flow rate, and laser power
The separate effects of volume V, flow rate Q and laser power W, are illustrated as contour plots in Fig. 5 after firing for t f = 60 s.
For all plots in Fig. 5 the initial temperature T 0 = 37 • C . The top row of plots in Fig. 5 are for T in = 37 • C and the bottom row for T in = 23 • C . Colours denote temperature change ΔT , and white regions of the plots present for T in = 23 • C (bottom row) indicate where ΔT < 0 . We note an increase in ΔT with W, a decrease in ΔT with Q, and a decrease in ΔT with V. Figure 5 also illustrates, by comparing the bottom row with the top row, the advantage of irrigating at room temperature ( T in = 23 • C ) rather than body temperature ( T in = 37 • C ) to maintain low temperatures.

Thermal dose
Sustained temperatures above 43 °C are known to cause thermal damage to many biological tissues [4].
A common metric for evaluating thermal dose is to convert a given temperature curve to the equivalent time at constant 43 • C by the formula where Equation (2) indicates that if T(t) ≡ 43 • C then t 43 = t . Any increase in T above 43 • C results in an exponential increase in the thermal dose. This mathematical description of thermal dose is based on evidence from in vitro and in vivo systems [4]. The time threshold for damage to occur is often taken to be t 43 = 120 min [3]. Using our mathematical model for temperature as a function of system parameters, we can determine the firing time t f such that t 43 remains below the safe thermal dose threshold of 120 min. We will refer to this as t safe f . Details of the  Fig. 6c, d, where it can be seen that the curves plateaus below T = 37 • C , indicated by the dashed blue line. As T ⋆ is independent of V, the white region is nearly identical between Fig. 6a, b (any discrepancies due to the slightly different value for used for the two experiment setups). A second set of example curves are shown in Fig. 6c, d for T in = 28 • C and Q = 20 mL min −1 . The value of t safe f is indicated by the red dots. We see, by comparing the black curve in Fig. 6c to the black curve in Fig. 6d, the comparatively sharper temperature rise due to the smaller volume. Thus, t safe f for Set A conditions is less than t safe f for Set B conditions if T ⋆ > 37 • C.

Discussion
We conducted a set of in vitro experiments to measure fluid temperature over time as a result of holmium laser lithotripsy. Resulting temperature curves agreed qualitatively with [1,3] and [2], with temperatures rising after the initiation of laser activation to reach a stable, elevated value, before decreasing after the laser was switched-off. We subsequently derived a mathematical model from conservation of energy principles, neglecting spatial temperature variation and considering only the volume-averaged temperature over time. This produced a single equation with which to compare against experimental data, and we fit for the parameter representing the thermal conductivity through the walls of the experimental vessel. We obtained excellent quantitative agreement with our experimental data and the data from [1] for all laser settings and irrigation flow rates considered.
The model predicts an increase in fluid temperature within the kidney with laser power and a decrease with irrigation flow rate. Irrigation flow rate may be altered clinically by either the adjusting the inflow (e.g. increasing pressure on the irrigation fluid) or outflow (e.g. use of access sheath) from the kidney. Our results corroborate the key findings of [1][2][3]. We also determined the volume of the working space as a key parameter in controlling the temperatures achieved; a smaller volume reaches higher temperatures more quickly, although the equilibrium temperature, which will be reached if the laser is fired for sufficient time, is independent of volume. This raises a pertinent question: what is the relevant volume of the collecting system during laser lithotripsy? Large variation in pelvicalyceal volume are reported in the literature [5], and of course, the comprising features-the renal pelvis and multiple calyces-have different sizes.
Therefore, patient-specific anatomy, as well as the location of a stone, may have a significant effect on fluid temperature and the potential for thermal damage. In addition to irrigation flow rate-or equivalently driving pressure-the temperature of irrigation fluid before it enters the patient also affects renal temperatures during laser lithotripsy. As relatively small quantities of fluid are used for irrigation, there will only be a minimal effect on global body temperature. However, local heating of irrigation fluid by higher powered laser use has a risk of detrimental effect on renal tissue. Due to the specific heat capacity of water, clinicians should be advised to use room temperature fluids (rather than warmed fluids) for irrigation to minimise thermal damage to the kidney from laser use.
In our experiments, we recorded temperature at two locations within the fluid-filled cylinder, and found no significant difference between the readings. However, we realise that our simplified cylindrical geometry does not represent the complex anatomical structure of the renal collecting system, and that changes in geometry may lead to spatially heterogeneous temperatures. Thus, although our mathematical model has the benefit of computational simplicity, it only predicts a single volume-averaged temperature. This will always be a lower bound on the maximum temperature within the kidney and thus, model results must be applied with caution to not overlook potential hot spots proximal to the laser fibre. More intensive simulations, such as those performed by [3], are required to predict the spatial distribution of temperature. Another limitation of the mathematical model is the need to determine the thermal conductivity of the material surrounding the fluid-filled container experimentally; from an in vivo standpoint, this will be the material properties of the renal tissue. Thus, further characterisation of renal volume and tissue properties from in vivo experimental data will both contribute to improving the accuracy of our mathematical model.
It is also important to note that in our experiments and mathematical model, we consider a minute of uninterrupted lasering. However, it is unusual for a laser to be used continuously for long periods of time during a case due to practical factors such as repositioning of the laser fibre, movement of the stone and inadequate view of the stone. Further work to investigate how typical "operator duty cycle" impacts on temperature changes. This is within the scope of the model provided the on/off intervals and laser settings are known.

Conclusion
This validated mathematical model allows prediction of the change in temperature within a kidney for any given renal volume, irrigation flow rate, irrigation fluid temperature, and laser power. The model obviates the need for repeated experimental validation. The model predicts scenarios where risk of renal tissue damage are high. With real-time knowledge of flow rate, irrigating fluid temperature and laser usage, safety warning levels could be predicted. In the meanwhile, clinicians should be aware of the potential risk from thermal injury and take measures to reduce the risk such as using room temperature irrigation fluid and judicious laser use.

Error estimation
The thermocouples had a readability of 10 −3 • C . Therefore, the error in the temperature measurements were taken to be = 5 × 10 −4 • C . The error for the temperature difference, ΔT , is hence √ 2 . We assume that the standard deviation of the measured temperatures over each set of five runs and the measurement errors due to the readability of the thermocouples are independent and thus our total error estimation is which is taken as the vertical error bars on the experimental data in Figs. 4 and 8.

Appendix B: Laser positioned 2 cm distal to scope tip
We also considered the experimental setup pictured in Fig. 7, where the laser is positioned 20 mm (as opposed to 10 mm) distal to the scope tip. The experimental settings were the same as those summarised in Table 1.
Experimental results are plotted against model predictions in Fig. 8. The agreement is good between theory and experiment, although there is more variability in the triangle data points in Fig. 8; these correspond near the base of the vessel, positioned only 7 mm from the laser fibre (Fig. 7). Thus, a potential hypothesis for this variability is laser-thermocouple interference.  Fig. 7 The laser fibre is 20 mm distal to the scope tip. Two thermocouples measured temperature over time positioned at the scope tip and 9 mm from the base of the container. These are indicated by a circle and triangle, respectively. Diagram not to scale The considered domain advection by irrigation, the temperature will plateau at the temperature of the incoming irrigation fluid.
To emulate the experiments in "Experimental setup" section , we assume the laser is switched on at t = t 1 and off at t = t 2 and approximate this behaviour with a piecewise function The solution for all time, assuming a firing function of the form in Eq. (18) can be written in a piecewise manner We also define the quantities i.e., the temperature change from the initial temperature and the laser firing time, respectively.

Model solution (no flow)
The solution to Eq. (14) subject to (8) for Q = 0 (i.e., no irrigation flow) and W(t) = W is The associated steady state is thus  0, for t < t 1 , W, for t 1 ≤ t ≤ t 2 , 0, for t > t 2 .
T(t;T 0 , W) = T 0 + W − exp(− t∕ cV) (T 0 − T ic ) + W .  where C is an arbitrary constant. If T ⋆ > 43 • C , we split the integral into four regions: If T ⋆ < 43 • C , we simply consider: Finally, if T ⋆ < 37 • C , we determine For T ⋆ > 37 • C we use Matlab's fsolve function to solve for t f such that t 43 = 120 minutes.