Analysis of hyperbolic Pennes bioheat equation in perfused homogeneous biological tissue subject to the instantaneous moving heat source

The one-dimensional hyperbolic Pennes bioheat equation under instantaneous moving heat source is solved analytically based on the Eigenvalue method. Comparison with results of in vivo experiments performed earlier by other authors shows the excellent prediction of the presented closed-form solution. We present three examples for calculating the Arrhenius equation to predict the tissue thermal damage analysis with our solution, i.e., characteristics of skin, liver, and kidney are modeled by using their thermophysical properties. Furthermore, the effects of moving velocity and perfusion rate on temperature profiles and thermal tissue damage are investigated. Results illustrate that the perfusion rate plays the cooling role in the heating source moving path. Also, increasing the moving velocity leads to a decrease in absorbed heat and temperature profiles. The closed-form analytical solution could be applied to verify the numerical heating model and optimize surgery planning parameters.


Introduction
Increasing attention has arisen recently in the use of moving heat source, which has found various applications related to medical therapies, e.g., laser cutting tools [1], bone drilling [2], and neurosurgery bone grinding [3]. Since the definition of temperature distribution in biological living tissues and obtaining exact solutions are necessary, researchers have suggested several bioheat transfer models [4,5]. One of the most well-known models is the Pennes bioheat equation [6]: where x is the spatial coordinate, t is the time, t and C t refer to the tissue density and specific heat; k t ∇ 2 T models the thermal diffusion and k t is the thermal conductivity of tissue; b b C b T b − T expresses the heat caused by convection within the heat per unit mass of the tissue. The b , b ,C b and T b are the density, perfusion rate per unit volume of tissue, specific heat and the temperature of blood, respectively. The Q is the source heating rate per unit volume of tissue.
The above simple bioheat equation is based on Fourier theory, which assumes an infinite heat conduction speed. It is applicable for long-time heating processes such as radiofrequency [7], microwave [8], hyperthermia cancer therapy [9], and cryosurgical ablation [10]. Non-Fourier or hyperbolic heat transfer equation (HHTE) assumes a relaxation time between heat flux and temperature gradient. Thus, a finite speed and wave behavior are considered for the heat propagation. Hence, Cattaneo and Vernotte [11, (1) 12] introduced an improved constitutive heat flux model in the following form: where q is the heat flux vector, and q is the thermal relaxation time.
For instantaneous heating of biological tissue which short heating times are employed and have long thermal relaxation time, a non-Fourier model should be considered employing the hyperbolic Pennes bioheat equation [13]: Various researchers have studied the hyperbolic Pennes bioheat conduction of the biological tissues [14,15]. According to the literature, some studies were focused on the heat transfer of the biological tissues subjected to moving heat sources. Kamis et al. [16] worked on the DPL heat equation of skin under the moving heat source. Sur et al. [17] investigated the thermal behavior of skin subjected to the moving heat source. They modeled the moving source with the Dirac function and calculated the temperature and thermal damage. Ma et al. [18] solved the DPL heat conduction model of living tissue utilizing the Green function approach. Kabiri and Talaee [19] investigated the thermal field and damage of biological tissues under the Gaussian laser moving heat source. Also, they obtained the exact solution of the parabolic Pennes bioheat equation for the kidney, liver, and skin tissues under the Dirac moving heat source. They parametrically studied the tissue temperature distribution caused by the external heat source different moving speeds [20].
This article is explored a purely analytical solution of the hyperbolic Pennes bioheat transfer equation (Eq. (3)) in one-dimensional slab subjected to the instantaneous moving heat source using the Eigenvalue method. The solution is applied for skin, liver, and kidney as three examples for calculating the Arrhenius equation to predict the tissue thermal damage analysis. The closed-form solution is introduced, and the effects of moving speed and perfusion term are studied parametrically. Hence, this article answers the following questions: • How to achieve a closed-form solution without any additional numerical procedures such as inverse Laplaces against some recent researches? • What is the difference between parabolic and hyperbolic models results in bioheat modeling of the tissues under moving heat?
• What is the effect of moving velocity, perfusion rate, and thermophysical properties of tissues on temperature profiles and thermal tissue damage?
The rest of the manuscript is arranged as follows: mathematical modeling and analytical procedure are presented in Sect. 2. Thermal damage calculations are included in Sect. 3. In Sect. 4, comparison and validation are discussed. Results and discussion are presented in Sect. 5. Finally, a conclusion is illustrated in Sect. 6.

Mathematical modeling and analytical procedure
The homogenous isotropic finite biological tissue slab subjected to a point moving heat source with a constant speed was considered for modeling ( Fig. 1). The thermal parameters of the tissues were assumed to be independent of the temperature and were given constant values in the temperature range of interest. The moving heat source was assumed a moving point with constant velocity [18]: where Q 0 is a constant power and is the Dirac function. The Dirac function could be interpreted as limiting of the laser heating Gaussian distribution function when the laser focus increases. Furthermore, the intensity distribution over the point of the heat source was assumed to be uniform. The adiabatic conditions were considered for both sides of the tissue. The boundary and initial conditions were considered as (Temperature of body is assumed 37 • C): For simplicity in the consequent analysis, the following dimensionless variables were introduced: Thus, the Eq. (3) can be rewritten in terms of dimensionless variables: where is the non-dimensional moving heat source function: where 0 is the strength of the heating source term The boundary conditions of Eq. (5) and (6) in the dimensionless form are: And the dimensionless initial conditions of Eq. (7) and (8) are: All the analytical solutions of the hyperbolic heat conduction (Eq. (15)) with source term are done with the Laplace scheme, which most needs the numerical procedures. However, a close-form solution using the Eigenfunction series do not need any numerical calculation [20]. Hence, the below solution series was adopted as the solution of Eq. (15) considering boundary conditions (17) and (18): where A n ( ) is the time-dependent constant which can be determined from substitution of the Eq. (21) into the Eq. (15): whereV n ( ) is the Fourier expansion coefficient of the source term function 2 + 4 : Thus, V n ( ) can be determined considering the Fourier expansion: Also, A n ( ) can be expressed as: The coefficient of A n ( ) is obtained from Eq. (25) (Appendix). Hence, the closed-form solution of the problem is obtained from Eq. (21).

Thermal damage calculation
The surgical applications involve protein denaturation, coagulation, vaporization, and tissue ablation. Magnitude and time duration of heat source affect the tissue damage degree. Without including the non-Fourier effect, the calculation of tissue damage might not be creditable. In order to evaluate the thermal damage of biological tissues, the following well-known Arrhenius equation was used [21]: where T is the absolute tissue temperature, is the frequency factor, and E a (Jmol −1 ) is the activation energy of protein denaturation reaction. To model three different tissues, i.e., tissues of skin, liver, and kidney, we used values for A 0 and E a as presented in experimental studies by [22][23][24], cf. (Table 1). Average of the many reactions taking place in heated tissue and various situations may be occurred for the tissue due to the values of thermal damage (Ω). An increase of Ω over the value of 1 leads to burn a higher degree and complete necrosis of the tissue [25].

Comparison and validation
Museux et al. [26] studied and measured the thermal effect of laser irradiation on pig skin numerically and experimentally. Their in-vivo experimental results were used to validate our present mathematical model. Figure 2 shows the temperature history profiles of the hyperbolic bioheat model and experimental analysis. A strong correlation is found between the predicted and measured results. It is demonstrated that both of them have the same trends, but the predicted temperature is higher than the measured one due to the natural convection and evaporative cooling have not been considered in the present study. Also, due to the tissues thermal parameters were assumed independent of the temperature, over prediction of the temperature is happened in the analytical model result. However, this comparison verifies that the hyperbolic bioheat model is practicable for estimating biological tissue temperature distribution.

Results and discussion
It is vital to design and characterization of heating therapy procedures. The tissues temperature could be obtained at any time and point of these with the presented analytical solution. In this study, three different tissues (kidney, liver, and skin) with different perfusion rates were considered under the concentrated moving heat source. The moving path begins from l 6 to the point of 5l 6 , which is equal to non-dimensional length X = 112.5 to X = 562.5. The source term was considered to move at = 0 from starting point and reach to the end point at = 3.33 . The heat source was turned off when the end of the path was reached, but the temperature profiles for cooling are continued.
The temperature history of skin tissue at various points under moving source with speed v = 0.01(ms −1 ) is shown in Fig. 3. Since the tissue temperature profiles not raised before reaching the moving source, the non-Fourier behavior of temperature profiles can be seen in this figure. It is also seen that each point temperature increases fast due to the passage of the moving heat source and decreases slowly after turning it off at the end of the path.
Also, skin tissue temperature distribution at various times is shown in Fig. 4. The fluctuation of the temperature profiles, which is the other property of hyperbolic profiles, can be seen in this figure. Progress of heat source is realized with the edge of temperature profiles in the tissue layer affected zone. The tissue temperature profile is proportional to the absorbed heat from the moving source in its passage, which is a general rule in the hyperbolic temperature profiles. Due to thermal relaxation behavior, the non-equilibrium temperature profiles in hyperbolic models have many dominant features that cannot be described with parabolic models in equilibrium state [27]. Therefore, the non-Fourier bioheat equation can be better estimate temperature distribution in tissues rather than the parabolic Fourier bioheat.
The temperature profiles are increased to a constant value which is dependent on the model parameters such as adiabatic boundary conditions, the tissue cooling perfusion rate, and moving source velocity. The presented analytical solution that can estimate the effects of these parameters is stated in the following.
The effect of heat source velocity on the temperature profiles is investigated in Fig. 5. Due to the less time available to heat the tissue, an increase of the moving velocity leads to a decrease of tissue maximum temperature and increases linear regression slope (absolute value) in the cooling state. Wave propagation is observed in the temperature distributions when the moving source velocity decreased. The thermal wave occurred along with the peak temperature profile in Fig. 5.
The blood perfusion rate makes convective heat transfer between the blood and tissue. Hence, blood perfusion convective cooling has a vital role in thermal therapy. Therefore, blood flow affects the efficiency of thermal treatments. The comparison of temperature profiles and temperature history of three tissues of skin, liver, and kidney is shown in Figs. 6 and 7 to study the effects of perfusion rate. The effect of a higher perfusion rate to decrease temperature profiles is seen in Fig. 6. The kidney tissue get a lower temperature than the liver and liver compared with skin tissue under the same heating conditions due to their higher perfusion rates. Figure 7 compare these three tissues temperature profiles under the same moving heat source when they reached to the end of the moving path ( = 3.33 and v = 0.01(ms −1 ) ). At first, temperature increases and then decreases due to cooling by blood perfusion. The higher perfusion rate leads to the more fast reduction of the temperature after heat source passage and increases the linear regression slope (absolute value) in cooling state. The temperature growth is dependent on conductivity, diffusivity, and perfusion rate. As expected, the maximum temperature is achieved in skin tissue because of its lower perfusion rate. The passage of the instantaneous moving heat source in higher perfusion rate tissue (kidney) may produce approximately uniform temperature compared with skin or liver tissues. Thus, due to the blood perfusion cooling function, the production of hot spots in the temperature profiles under concentric moving heat source is more expected in lower perfusion rate tissues. Also, this figure shows the effect of the heating frequency on the temperature-fluctuation in the skin tissue. Tissue necrosis and treatment efficiency depend on the thermal damage. The tissues thermal damage history is shown in Fig. 8. The skin tissue thermal damage is occurred highly compared with liver and kidney tissues due to its higher temperature and lower perfusion rate. Also, skin tissue necrosis occurs faster than the liver, and the liver is faster than the kidney.

Conclusion
In this paper, a closed-form analytical solution was introduced for the hyperbolic Pennes bioheat transfer of the finite medium subjected to the instantaneous moving heat source. It was applied to model the characteristics of three tissues of skin, liver, and kidney. The effects of perfusion rate and moving velocity on the temperature profiles were investigated. The higher perfusion term leads to lower temperature amplitude and more uniform profiles under the concentric moving heat source motion. The presented analytical solution has been shown better estimation the thermal distribution in living tissue than the other thermotherapy models. Finally, this study can be used to verify the other heat therapies under the moving heat source, and the knowledge acquired for its development should be used in future works.

Compliance with ethical standards
Conflict of interest The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/. A n ( ) = 4 0 L e −(1+ ) n L 2 n − ( + 1) 1 + n + L 2 (2 + V ) + n 2 2 V 2 (4 + 2 + V ) Cosh n + n e (1+ ) L L 1 − n + 1 + n + L 2 (2 + V ) − n 2 2 V 2 (4 + 2 + V ) Cos n V L + n V L 2 + 2 + 2 V + 6h + 5 − 2 n − n 2 2 V 2 Sin n V L + n 4 4 V 4 − (1 + ) (1 + ) 2 − 2 n L 4 (2 + V ) + L 2 n 2 2 V 2 2 n + 2 − (1 + )V − 1 Sinh n ∕ n 2 n − (1 + ) 2 2 L 4 + 2 2 n + (1 + ) 2 L 2 n 2 2 V 2 + n 4 4 V 4