Three-Phase-Lag Bio-Heat Transfer Model of Cardiac Ablation

Significant research efforts have been devoted in the past decades to accurately modelling the complex heat transfer phenomena within biological tissues. These modeling efforts and analysis have assisted in a better understanding of the intricacies of associated biological phenomena and factors that affect the treatment outcomes of hyperthermic therapeutic procedures. In this contribution, we report a three-dimensional non-Fourier bio-heat transfer model of cardiac ablation that accounts for the three-phase-lags (TPL) in the heat propagation, viz., lags due to heat flux, temperature gradient, and thermal displacement gradient. Finite element-based COMSOL Multiphysics software has been utilized to predict the temperature distributions and ablation volumes. A comparative analysis has been conducted to report the variation in the treatment outcomes of cardiac ablation considering different bio-heat transfer models. The effect of variations in the magnitude of different phase lags has been systematically investigated. The fidelity and integrity of the developed model have been evaluated by comparing the results of the developed model with the analytical results of the recent studies available in the literature. This study demonstrates the importance of considering non-Fourier lags within biological tissue for predicting more accurately the characteristics important for the efficient application of thermal therapies.


Introduction
Minimally invasive thermal therapies (viz., radiofrequency, laser, microwave, and cryo-and ultrasound ablations) are well-established as a promising localized treatment option for treating different types of cancers, cardiovascular diseases, neurodegenerative disorders, as well as chronic pain [1][2][3][4][5]. The goal of thermal ablation is to attain direct cellular injury by the application of extreme temperatures, either too high (hyperthermic ablation) or too low (cryoablation). During thermal therapies, the applicator/catheter is introduced percutaneously into the target tissue with the aid of image-guided modalities, such as computed tomography (CT), magnetic resonance imaging (MRI), or ultrasound [1,6]. The energy deposited within the biological tissue from the applicator in high-temperature thermal therapies leads to localized heating at the target site resulting in protein denaturation, mitochondrial dysfunction, halt in enzyme function, and cell membrane collapse [2]. Importantly, minimally invasive thermal therapies result in a potential reduction in the physical, emotional, and financial stress on patients as compared to surgery and are quite useful for treating patients on whom conventional surgery cannot be performed due to associated comorbidities.
Computational modeling and simulations have played an important role in better understanding the underlying physics of various thermal ablative procedures among

Materials and Methods
The three-dimensional computational domain of cardiac ablation considered in the present study comprising the cardiac tissue, blood chamber, and the applicator, is presented in Figure 1 [37]. Further, the heat transport utilizing Pennes' bio-heat transfer model is based on the classical Fourier's law, which assumes that heat propagation occurs at an infinite speed and relates to the heat flux (q) as [12]: (1) Fluids 2022, 7, 180 3 of 15 perature (K), k is the thermal conductivity (W/m/K), ρb is the density of blood (kg/m 3 ), cb is the specific heat capacity of blood (J/kg/K), ωb is the blood perfusion rate (1/s), Tb is the temperature of blood entering the tissue, the term [ρbcbωb(T − Tb)] represents the heat sink effect caused by the small capillary vasculature, Qm is the metabolic heat generation (W/m 3 ), which is normally ignored due to its minimal impact compared to other heat source terms during thermal ablation, and Qp (W/m 3 ) represents the volumetric heat generation during thermal ablative procedures (e.g., radiofrequency or laser heat sources). The thermal lag associated with the time delay between the heat flux and temperature gradient was independently introduced by Cattaneo [15] and Vernotte [16] in the Fourier's conduction law as: The generalized Fourier conduction of Pennes' bio-heat transfer equation based on energy conservation in the Cartesian coordinate system is represented as: where ρ is the density (kg/m 3 ), c is the specific heat capacity (J/kg/K), T is the tissue temperature (K), k is the thermal conductivity (W/m/K), ρ b is the density of blood (kg/m 3 ), c b is the specific heat capacity of blood (J/kg/K), ω b is the blood perfusion rate (1/s), T b is the temperature of blood entering the tissue, the term [ρ b c b ω b (T − T b )] represents the heat sink effect caused by the small capillary vasculature, Q m is the metabolic heat generation (W/m 3 ), which is normally ignored due to its minimal impact compared to other heat source terms during thermal ablation, and Q p (W/m 3 ) represents the volumetric heat generation during thermal ablative procedures (e.g., radiofrequency or laser heat sources). The thermal lag associated with the time delay between the heat flux and temperature gradient was independently introduced by Cattaneo [15] and Vernotte [16] in the Fourier's conduction law as: where τ q is the thermal relaxation time that represents the time delay between the heat flux vector and the temperature gradient. The constitutive relation obtained by incorporating the thermal relaxation time phase lag as the linear extension of Fourier law is known as the single-phase-lag (SPL) non-Fourier heat transfer model (or Cattaneo-Vernotte (C-V) heat transfer model). In this case, the resulting hyperbolic-type SPL heat equation is given by: Tzou [17] proposed another relaxation time to account for the effects of micro-structural interaction along with fast transient effects of thermal inertia, and accordingly, the conduction relation between heat flux and temperature gradient was modified as: where, as before, τ q is the phase lag for the heat flux, and τ t is the phase lag for the temperature gradient. For τ q < τ t , heat flux precedes the temperature gradient, and for τ q > τ t , the temperature gradient precedes the heat flux. It is noteworthy to mention that in order to obey the causality principle τ q ≥ τ t , as highlighted in [13]. By substituting Equation (2) in Equation (5), and after some mathematical manipulations, the dual-phaselag (DPL) non-Fourier heat transfer model can be obtained as: For τ t = 0, Equation (6) reduces to the SPL bio-heat transfer model (Equation (4)), and further reduces to Fourier's heat transfer model (Equation (2)) for τ q = 0.
Later, Choudhuri [30] proposed a three-phase-lag model by incorporating an additional phase lag of thermal displacement gradient in the DPL model based, in its turn, on the thermoelastic model developed by Green and Naghdi [31,32] and as given by: where τ q , τ t , and τ v are phase lags associated with heat flux, temperature gradient, and thermal displacement gradient, respectively, v is thermal displacement satisfying . v = T, and k * is the rate of thermal conductivity. After some mathematical manipulations of Equation (7) combined with Equation (2), the following constitutive equation for TPL can be obtained [33]: Importantly, the phase lags should obey the following inequality to obtain a stable solution, 0 ≤ τ v < τ t < τ q [38]. For k* = 0 and integrating w.r.t time, Equation (8) will reduce to the DPL bio-heat transfer model. Further with τ t = 0, this equation will reduce to the SPL model, and additionally, with τ q = 0, the equation will reduce to the Pennes' bio-heat transfer model. In the present study, we have simulated a cardiac ablative procedure at a constant power. As such, for a constant power source, the time-derivative of the volumetric heat generation (i.e., Q p ) would be zero in the TPL model (Equation (8)). To avoid this paradox, we have considered a constant temperature at the tip of the applicator, which is maintained at 95 • C. It is noteworthy to mention that the constant temperature can be easily maintained during the thermal ablative procedures utilizing a closed-loop feedback proportionalintegral-derivative (PID) controller, as demonstrated in previous studies [39][40][41][42][43]. A comparative analysis has been conducted by comparing the temperature distributions predicted within the myocardium tissue incorporating Pennes' (Equation (2)), SPL (Equation (4)), DPL (Equation (6)), and TPL (Equation (8)) bio-heat transfer models.
The blood flow within the cardiac chamber was modeled utilizing the laminar and incompressible flow as given by the coupled momentum and continuity, Navier-Stokes equations: where u is the velocity vector (m/s), ρ b is the density of blood (kg/m 3 ), p is the blood pressure (Pa), µ is the dynamic viscosity of blood (=2.1 × 10 −3 Pa.s), and F is the body force neglected in this study [43][44][45]. It is noteworthy to mention that Newtonian behavior has been used to model the blood, which is adequate for shear rates >100 s −1 [46,47]. The fluid-structure interaction (FSI) has received significant attention in research related to cardiovascular modeling [48][49][50][51][52][53][54][55]. In this study, FSI only pertains to the convective heat-sink (cooling) effect induced by the blood flow during cardiac ablation, as the cardiac tissue has been modeled as a rigid body. However, future studies are warranted to incorporate the mechanical deformation of the cardiac tissue caused by the pulsatile blood flow within the cardiac chamber, as has been recently reported in [56,57]. Further, the blood perfusion rate has been modeled utilizing a temperature-dependent piecewise model, whereby a constant value of blood perfusion rate prevails below the temperature of 50 • C and ceases beyond it [40,58], as given by: The ablation volume (V) attained during the cardiac ablation was obtained utilizing the isotherm of 50 • C, i.e., the volume of myocardium tissue having a temperature of ≥50 • C after the selected treatment time, as given by [37,40,59]: A finite-element method (FEM)-based COMSOL Multiphysics software [60] was used to solve different bio-heat transfer models for predicting the temperature distributions and ablation volumes during 60 s of cardiac ablation. The initial temperature of the entire computational domain presented in Figure 1b was 37 • C, and the tip of the applicator was maintained at a temperature of 95 • C throughout the simulation. Further, the blood flow inside the cardiac chamber was modeled by considering two values for the blood velocity at the inlet, viz., 8.5 cm/s and 3 cm/s representing the high and low blood flow rates, respectively [19,[43][44][45]58]. The outlet of the cardiac chamber was subjected to zero pressure conditions. No-slip boundary conditions were imposed on the wall of the applicator in the cardiac chamber as well as the tissue-fluid interface. The thermophysical properties of cardiac tissue considered in the present study were: ρ = 1060 kg/m 3 , c = 3111 J/kg/K, ρ b = 1000 kg/m 3 , c b = 4180 J/kg/K, ω b = 0.017 s −1 , k = 0.54 W/m/K, k* = 0.1 W/m/K/s, τ q = 16 s, τ T = 6 s and τ v = 2 s [19,33,[43][44][45]58]. The computational domain presented in Figure 1b was discretized using the heterogeneous tetrahedral mesh elements with additional refinements close to the applicator where higher thermal gradients were ex-pected. Furthermore, a grid independence study was performed to determine the optimal number of mesh elements that would result in a mesh-independent solution. The mesh refinements were carried out using the ablation volume convergence criterion, i.e., the mesh was progressively refined until the absolute error for ablation volume was less than 0.5 % compared to the previous mesh size. The final mesh comprises of 551,690 domain elements, 26,674 boundary elements, and 1288 edge elements. The constitutive equations for heat transfer and fluid flow, presented earlier, were solved using the multifrontal massively parallel sparse direct solver (MUMPS) [61], employing a segregated solution approach. All of the simulations were conducted on a Dell Precision 7920 Tower workstation with 96 GB RAM and 2.20 GHz Intel ® Xeon ® processors with a computational time of less than 4 h for each simulation.

Results and Discussion
The developed TPL model has been validated with the results of one of the recent studies [33] that provides the analytical solution for the bio-heat response of skin tissue based on the TPL model. In particular, the same parameters and boundary conditions of [33] have been utilized, viz., outer skin surface temperature = 80 • C, arterial blood temperature = 37 • C, thickness of the tissue = 9 mm, density of the tissue = 1190 kg/m 3 , specific heat capacity of tissue = 3600 J/kg/K, density of the blood = 1060 kg/m 3 , specific heat capacity of the blood = 3770 J/kg/K, the blood perfusion rate = 1.87 × 10 −3 s −1 , thermal conductivity of the tissue = 0.235 W/m/K, the rate of thermal conductivity = 0.1 W/m/K/s, phase lag associated with heat flux = 16 s, phase lag associated with temperature gradient = 6 s, and the phase lag associated with thermal displacement gradient = 2 s. The comparative analysis of the present results with that of [33] is presented in Figure 2. As evident from Figure 2, a good agreement has been obtained between the results predicted from the present model with those reported in [33]. Hence the TPL model used in the present study is efficient and lends great confidence to the results reported in the next sections.  Figure 3 presents the comparison between the predictions with the Pennes', SPL, DPL, and TPL bio-heat transfer models during cardiac ablation under low blood flow conditions within the cardiac chamber. The temporal evolutions of the temperature at a distance of 0.5 mm, 1 mm, and 2 mm below the applicator tip are presented for different bioheat transfer models in Figure 3a-c, respectively. It can be seen from Figure 3 that for all locations, the thermal propagation speed is fastest in the Pennes' model, which is based on the Fourier law of heat conduction that presumes an infinitely fast propagation of thermal signals within the medium. Furthermore, the thermal propagation speed is lowest in the SPL model, owing to the delayed response between the heat flux and temperature gradient that leads to the propagation of thermal disturbance at a finite speed. The thermal  Figure 3 presents the comparison between the predictions with the Pennes', SPL, DPL, and TPL bio-heat transfer models during cardiac ablation under low blood flow conditions within the cardiac chamber. The temporal evolutions of the temperature at a distance of 0.5 mm, 1 mm, and 2 mm below the applicator tip are presented for different bio-heat transfer models in Figure 3a-c, respectively. It can be seen from Figure 3 that for all locations, the thermal propagation speed is fastest in the Pennes' model, which is based on the Fourier law of heat conduction that presumes an infinitely fast propagation of thermal signals within the medium. Furthermore, the thermal propagation speed is lowest in the SPL model, owing to the delayed response between the heat flux and temperature gradient that leads to the propagation of thermal disturbance at a finite speed. The thermal propagation speed in the DPL and TPL models lies between those of the Pennes' and SPL, with the TPL model thermal response being lower than that of DPL. It is noteworthy to mention that the TPL model exhibits both the diffusion and wave characteristics of heat propagation, similar to the DPL model. Furthermore, as evident from Figure 3, the temperature prediction with the SPL model is always on the higher side for all three points. The deviations in the temperature prediction between the Pennes' and TPL models are significantly higher for the first 30 s or so and eventually diminish with the passage of time.  Figures 4a,b present the temperature distribution predicted with different bio-heat transfer models along the line below the tip of the applicator at 30 s and 60 s, respectively. The distance of 0 mm reflects the tip of the applicator that is maintained at the temperature of 95 °C. As evident from Figure 4a, among all the bio-heat transfer models, the temperature predicted by the SPL model is on the slightly higher side until it is close to 1.5 mm. After this location, the Pennes' model predicts a higher temperature compared to the other three models' predictions. This can be attributed to the fact that the location where the Pennes' model temperature prediction becomes higher than that of the SPL model shows the position where the thermal wave arrives in the SPL model. So, before this position, more heat was concentrated and confined due to the finite thermal propagation speed in the SPL model as compared to the infinite thermal propagation speed in the Pennes' bioheat transfer model, which leads to more thermal energy transfer to the unaffected zone away from the applicator. These findings are consistent with those recently reported in  The distance of 0 mm reflects the tip of the applicator that is maintained at the temperature of 95 • C. As evident from Figure 4a, among all the bio-heat transfer models, the temperature predicted by the SPL model is on the slightly higher side until it is close to 1.5 mm. After this location, the Pennes' model predicts a higher temperature compared to the other three models' predictions. This can be attributed to the fact that the location where the Pennes' model temperature prediction becomes higher than that of the SPL model shows the position where the thermal wave arrives in the SPL model. So, before this position, more heat was concentrated and confined due to the finite thermal propagation speed in the SPL model as compared to the infinite thermal propagation speed in the Pennes' bio-heat transfer model, which leads to more thermal energy transfer to the unaffected zone away from the applicator. These findings are consistent with those recently reported in [33]. Furthermore, the deviations between the temperature predictions of different bio-heat transfer models are more pronounced at 30 s and eventually diminish with the passage of treatment time of cardiac ablation, as can be seen from Figure 4a Figure 5 presents the comparison of predicted ablation volume with different bioheat transfer models after 10 s, 30 s, and 60 s of the cardiac ablative procedure. As depicted in Figure 5, there prevail significant deviations in the ablation volume predicted by Pennes', SPL, DPL, and TPL bio-heat transfer models post 10 s of cardiac ablation. In particular, the ablation volumes predicted by SPL, DPL, and TPL models are 44%, 57%, and 75%, respectively, lesser as compared to the Pennes' model predictions. This deviation decreases with an increase in the treatment time of cardiac ablation. For instance, at the end of 30 s and 60 s ablations, the differences between Pennes' and TPL model predictions are 37% and 21%, respectively. Similarly, at the end of 30 s and 60 s ablation, the differences between Pennes' and DPL model predictions are 11% and 2%, respectively. Thus, the importance of considering the non-Fourier heat transfer phenomena in the computational models of thermal ablations becomes quite essential in therapies for which treatment time ≤ 60 s. It is noteworthy to mention that the treatment time for cardiac ablation varies between 30-60 s, so it becomes extremely important to consider the various lags introduced due to the non-Fourier phenomena in the bio-heat transfer models to obtain accurate predictions. Similarly, for the thermal ablative procedures that utilize pulsed laser ablation (On-Off cycles), the non-Fourier phenomena will play a more pronounced role in the accurate treatment outcome predictions from computational studies.  Figure 5 presents the comparison of predicted ablation volume with different bio-heat transfer models after 10 s, 30 s, and 60 s of the cardiac ablative procedure. As depicted in Figure 5, there prevail significant deviations in the ablation volume predicted by Pennes', SPL, DPL, and TPL bio-heat transfer models post 10 s of cardiac ablation. In particular, the ablation volumes predicted by SPL, DPL, and TPL models are 44%, 57%, and 75%, respectively, lesser as compared to the Pennes' model predictions. This deviation decreases with an increase in the treatment time of cardiac ablation. For instance, at the end of 30 s and 60 s ablations, the differences between Pennes' and TPL model predictions are 37% and 21%, respectively. Similarly, at the end of 30 s and 60 s ablation, the differences between Pennes' and DPL model predictions are 11% and 2%, respectively. Thus, the importance of considering the non-Fourier heat transfer phenomena in the computational models of thermal ablations becomes quite essential in therapies for which treatment time ≤60 s. It is noteworthy to mention that the treatment time for cardiac ablation varies between 30-60 s, so it becomes extremely important to consider the various lags introduced due to the non-Fourier phenomena in the bio-heat transfer models to obtain accurate predictions. Similarly, for the thermal ablative procedures that utilize pulsed laser ablation (On-Off cycles), the non-Fourier phenomena will play a more pronounced role in the accurate treatment outcome predictions from computational studies. The influence of the blood perfusion rate on the TPL model predictions for cardiac ablation has been presented in Figure 6. Importantly, blood perfusion within the biological tissue accounts for microvascular tissue perfusion in absence of large blood vessels [9]. The temporal evolution of temperature with and without blood perfusion at a point 1 mm below the tip of the applicator has been presented in Figure 6a. As evident from Figure  6a, the temperature profile obtained with accounting for blood perfusion effects is significantly lower as compared to the temperature profile obtained without considering the blood perfusion. The lower temperature profile can be attributed to the heat-sink effect caused by the blood perfusion during the thermal ablative procedure. Furthermore, as depicted in Figure 6a, the deviation in the two curves of the temperature profile is lower till the first 40 s or so and increases afterward as the steady-state temperature is attained. Thus, the effect of the consideration of the blood perfusion within the bio-heat transfer model of cardiac ablation is to lower the steady-state temperature within the tissue. Figure  6b,c present the comparison of temperature distribution obtained with and without considering the effect of blood perfusion in the TPL bio-heat transfer model along the line below the tip of the applicator at 30 s and 60 s, respectively. The deviation between the temperature distribution obtained with and without blood perfusion at 30 s is lower as compared to that obtained at 60 s. Again, this highlights the fact that blood perfusion has more influence on the steady-state temperature during the thermal ablation, as at 60 s, the steady-state temperature has been reached, which is evident from Figures 3 and 6a. The temporal variation of the ablation volume predicted with the TPL bio-heat transfer model under low blood flow (3 cm/s), and high blood flow (8.5 cm/s) conditions within the cardiac chamber are presented in Figure 7. It can be seen from Figure 7 that initially, the deviation between the ablation volume is quite negligible, and it increases with the passage of time as the steady-state temperature is attained. Thus, both the macro-and microblood perfusions significantly affect the treatment outcomes of the cardiac ablation as predicted using the TPL bio-heat transfer model. The influence of the blood perfusion rate on the TPL model predictions for cardiac ablation has been presented in Figure 6. Importantly, blood perfusion within the biological tissue accounts for microvascular tissue perfusion in absence of large blood vessels [9]. The temporal evolution of temperature with and without blood perfusion at a point 1 mm below the tip of the applicator has been presented in Figure 6a. As evident from Figure 6a, the temperature profile obtained with accounting for blood perfusion effects is significantly lower as compared to the temperature profile obtained without considering the blood perfusion. The lower temperature profile can be attributed to the heat-sink effect caused by the blood perfusion during the thermal ablative procedure. Furthermore, as depicted in Figure 6a, the deviation in the two curves of the temperature profile is lower till the first 40 s or so and increases afterward as the steady-state temperature is attained. Thus, the effect of the consideration of the blood perfusion within the bio-heat transfer model of cardiac ablation is to lower the steady-state temperature within the tissue. Figure 6b,c present the comparison of temperature distribution obtained with and without considering the effect of blood perfusion in the TPL bio-heat transfer model along the line below the tip of the applicator at 30 s and 60 s, respectively. The deviation between the temperature distribution obtained with and without blood perfusion at 30 s is lower as compared to that obtained at 60 s. Again, this highlights the fact that blood perfusion has more influence on the steady-state temperature during the thermal ablation, as at 60 s, the steady-state temperature has been reached, which is evident from Figures 3 and 6a. The temporal variation of the ablation volume predicted with the TPL bio-heat transfer model under low blood flow (3 cm/s), and high blood flow (8.5 cm/s) conditions within the cardiac chamber are presented in Figure 7. It can be seen from Figure 7 that initially, the deviation between the ablation volume is quite negligible, and it increases with the passage of time as the steady-state temperature is attained. Thus, both the macro-and micro-blood perfusions significantly affect the treatment outcomes of the cardiac ablation as predicted using the TPL bio-heat transfer model.      The result of parametric studies of the influence of single-phase-lag (τ q ), dual-phaselag (τ t ), and three-phase-lag (τ v ) on the thermal behavior obtained with the TPL model is presented in Figures 8 and 9. The effect of τ q on temperature response at a point 1 mm below the tip of the applicator is presented in Figure 8a. As evident from this figure, the higher the value of τ q , the higher would be the lagging characteristic of thermal behavior that would lead to a steeper rise in temperature. Thus, the heat propagation within the tissue decreases with an increase in the value of τ q . The result of parametric studies of the influence of single-phase-lag (τq), dual-phaselag (τt), and three-phase-lag (τv) on the thermal behavior obtained with the TPL model is presented in Figures 8 and 9. The effect of τq on temperature response at a point 1 mm below the tip of the applicator is presented in Figure 8a. As evident from this figure, the higher the value of τq, the higher would be the lagging characteristic of thermal behavior that would lead to a steeper rise in temperature. Thus, the heat propagation within the tissue decreases with an increase in the value of τq. The effect of τt on the temperature response at a point 1 mm below the tip of the applicator is presented in Figure 8b. As evident from Figure 8b, the influence of variation in τt is smaller on the thermal propagation and lag characteristics as compared to τq. Figure 8c presents the influence of variation in the thermal displacement gradient τv on the temperature response at a point 1 mm below the tip of the applicator. As evident from Figure 8c, the increase in the value of τv results in a decrease in the predicted temperature during cardiac ablation. The effect of the influence of variation in τq, τ t , and τv on the ablation volumes predicted after 10 s, 30 s, and 60 s of cardiac ablation are presented in Figure  9. As evident from this figure, the maximum variation in the ablation volume for different values of lags occurs at 10 s and then decreases with the increase in treatment time. The effect of τ t on the temperature response at a point 1 mm below the tip of the applicator is presented in Figure 8b. As evident from Figure 8b, the influence of variation in τ t is smaller on the thermal propagation and lag characteristics as compared to τ q . Figure 8c presents the influence of variation in the thermal displacement gradient τ v on the temperature response at a point 1 mm below the tip of the applicator. As evident from Figure 8c, the increase in the value of τ v results in a decrease in the predicted temperature during cardiac ablation. The effect of the influence of variation in τ q , τ t , and τ v on the ablation volumes predicted after 10 s, 30 s, and 60 s of cardiac ablation are presented in Figure 9. As evident from this figure, the maximum variation in the ablation volume for different values of lags occurs at 10 s and then decreases with the increase in treatment time. The limitations of this numerical study are mainly related to the considerations of unidirectional and constant blood flow, which is highly complex and variable in the actual patient scenarios. Moreover, the flow rate and flow direction of blood in the cardiac chamber are highly dependent on the cardiac cycle. Further, we have considered the homogeneous and isotropic properties of cardiac tissue with a sharp insertion of the applicator without considering any mechanical deformation. In addition, the material properties of the cardiac tissue are considered to be constant, and in actuality, they vary with temperature. Albeit, this has been done purposely to quantify the variations only induced due to the phase lags and to suppress the variations induced due to other properties. Another limitation is associated with the validation of the numerical model with the results reported in the previous numerical study due to the lack of experimental in vivo data/findings available in the literature. Future studies are required to address these limitations to provide further improve the predicted temperature distribution and ablation volume during cardiac ablation. The extension of these models in terms of multiscale coupling at the tissue-cellular scale can be found in [62][63][64][65]. Despite these limitations, we were able to highlight the importance of considering the phase lags in the computational modeling of cardiac ablation for predicting accurate treatment outcomes. The limitations of this numerical study are mainly related to the considerations of unidirectional and constant blood flow, which is highly complex and variable in the actual patient scenarios. Moreover, the flow rate and flow direction of blood in the cardiac chamber are highly dependent on the cardiac cycle. Further, we have considered the homogeneous and isotropic properties of cardiac tissue with a sharp insertion of the applicator without considering any mechanical deformation. In addition, the material properties of the cardiac tissue are considered to be constant, and in actuality, they vary with temperature. Albeit, this has been done purposely to quantify the variations only induced due to the phase lags and to suppress the variations induced due to other properties. Another limitation is associated with the validation of the numerical model with the results reported in the previous numerical study due to the lack of experimental in vivo data/findings available in the literature. Future studies are required to address these limitations to provide further improve the predicted temperature distribution and ablation volume during cardiac ablation. The extension of these models in terms of multiscale coupling at the tissue-cellular scale can be found in [62][63][64][65]. Despite these limitations, we were able to highlight the importance of considering the phase lags in the computational modeling of cardiac ablation for predicting accurate treatment outcomes.

Conclusions
In the present study, we have reported a three-dimensional TPL non-Fourier bio-heat transfer model for cardiac ablation, considering the phase lags of heat flux, temperature gradient, and thermal displacement gradient. The model has also accounted for fluid-structure interactions. A comparative analysis has been conducted to highlight the deviations between the temperature distribution and ablation volume obtained with Pennes', SPL, DPL, and TPL bio-heat transfer models. Our findings suggest that there are significant prevailing variations in the predicted temperature distribution and ablation volumes among different bio-heat transfer models during cardiac ablation. These variations are very pronounced up to the first 30 s or so and decrease thereof, justifying the need for the consideration of these phase lags in the bio-heat transfer models for predicting more realistic estimates of ablation volume and temperature distribution, thus minimizing the deviations between the experimentally obtained and numerically predicted results. Importantly, these considerations would play a more pronounced role in thermal ablative modalities where the treatment time is 60 s or less (as in the case of cardiac ablation) or thermal ablative procedures that utilize pulsed laser ablation (On-Off cycles). Furthermore, parametric studies have also been conducted to quantify the effect of variations in the magnitude of different phase lags, blood perfusion rate, and blood flow on the treatment outcomes of cardiac ablation as predicted by the TPL bio-heat transfer model.