Parametric Analysis of Effective Tissue Thermal Conductivity, Thermal Wave Characteristic, and Pulsatile Blood Flow on Temperature Distribution during Thermal Therapy.

Parametric Analysis of Effective Tissue Thermal

Conductivity, Thermal Wave Characteristic, and Pulsatile

Blood Flow on Temperature Distribution during Thermal

Therapy

Tzu-Ching Shih)a

Department of Biomedical Imaging and Radiological Science, China Medical University, Taichung, 40402, Taiwan; and Department of Biomedical Informatics, Asia University, Taichung, 41354, Taiwan

Huang-Wen Huang

Department of Innovative Information and Technology, Langyang Campus, Tamkang University, ILan County, 26247, Taiwan

Wei-Che Wei

Department of Biomedical Imaging and Radiological Science, China Medical University, Taichung, 40402, Taiwan

Tzyy-Leng Horng

Department of Applied Mathematics, Feng Chia University, Taichung, 40724, Taiwan Corresponding author:

Tzu-Ching Shih, Ph.D.

Department of Biomedical Imaging and Radiological Science, China Medical University, Taichung, Taiwan

Address: 91 Hsueh-Shih Road, Taichung, 40402, Taiwan Telephone: +886-4-2205-3366 # 7709

Fax: +886-4-2208-1447

ABSTRACT

This study examines the coupled effects of pulsatile blood flow in a thermally significant blood vessel, the effective thermal conductivity of tumor tissue, and the thermal relaxation time in solid tissues on the temperature distributions during thermal treatments. Due to the cyclic nature of blood flow as a result of the heartbeat, the blood pressure gradient along a blood vessel was modeled as a sinusoidal change to imitate a pulsatile blood flow. Considering the enhancement in the thermal conductivity of living tissues due to blood perfusion, the effective tissue thermal conductivity was investigated. Based on the finite propagation speed of heat transfer in solid tissues, a modified wave bio-heat transfer transport equation in cylindrical coordinates was used. The numerical results show that a larger relaxation time results in a higher peak temperature. In the rapid heating case I (i.e., heating power density of 100 W cm-3 _{and heating duration of 1 s) and a heartbeat frequency of 1 Hz,} the maximum temperatures were 62.587 and 63.107 °C for thermal relaxation times of 0.464 and 6.825 s, respectively. In contrast, the same total heated energy density of 100 J cm-3 _{in a} slow heating case (i.e., heating power density of 5 W cm-3 _{and heating duration of 20 s)} revealed maximum temperatures of 57.724 and 61.233 °C for thermal relaxation times of 0.464 and 6.825 s, respectively. In rapid heating cases, the occurrence of the peak

temperature exhibits a time lag due to the influence of the thermal relaxation time. In contrast, in slow heating cases, the peak temperature may occur prior to the end of the heating period. Moreover, the frequency of the pulsatile blood flow does not appear to affect the maximum temperature in solid tumor tissues.

Keywords: Effective tissue thermal conductivity, thermal wave characteristics, pulsatile blood flow, bio-heat transfer equation, thermal therapy

1. INTRODUCTION

Quantitative heat transfer in living tissues is an essential issue in many medical treatments, such as cryosurgery [13] and hyperthermia [46]. Heat transfer in living tissues is a complicated process that involves heat conduction through solid tissues, heat convection between moving fluids (e.g., blood flow, lymphatic fluids, and interstitial fluids) and solid tissues, heat generation by tissue metabolism, non-directional tissue blood perfusion (e.g., thermoregulatory mechanisms, i.e., the scalar blood perfusion term may be used as a heat sink for thermal therapy or a heat source for cryosurgery), which was first quantitatively proposed by Pennes (1948) [7], and heat deposition through an external heating source in thermal treatments. In the Pennes model for a blood-perfused tissue, there is an essential

assumption that energy exchange between blood vessels and the surrounding tissues occurs mainly across the vascular wall of capillaries, which have diameters less than approximately 200 µm [8]. During thermal treatments, it is important to have complete knowledge of the temperature distribution in living tissues. Thermal therapy using high temperatures can kill cancer cells [913]. High temperature induces the thermal denaturation of proteins and can cause thermal coagulation necrosis in biological tissues [1416].

Blood flow can significantly affect the temperature distributions during thermal treatments, particularly in large blood vessels [1720]. Previous studies have investigated the significance of thermally significant blood vessels (i.e., larger than 200 µm in diameter) in the absorbed power density and the temperature distributions during thermal therapies. In their models, however, these researchers did not consider the impact of pulsatile blood flow due to the periodic-in-time nature of the heart pumping. The velocity profile inside a straight blood vessel with pulsatile blood flow, which is driven by an oscillating pressure gradient, was first analyzed by Womersley [21]. He used a pressure gradient with a varying time period coupled with the frequency to describe the periodic phenomenon of pulsatile blood flow in blood vessels and obtained an exact solution of the equations of viscous fluid motion under a pressure gradient with a periodic function of time. Rohlf and Tenti used the techniques of dimensional analysis to investigate the meaning of the Womersley number for pulsatile blood flow in small vessels [22]. Moreover, Craciunescu and Clegg studied the

effect of pulsatile blood velocities on bio-heat transfer in a straight rigid blood vessel [23]. These researchers demonstrated that the velocity pulsations of thermally terminal arteries (0.04  1 mm) have a small influence on the temperature distribution. Nevertheless, these researchers only focused on the heat transfer inside a single blood vessel and did not consider the thermal interaction between the pulsatile blood flow and its surrounding perfused tissue. Furthermore, Horng and his colleagues investigated the effects of pulsatile blood flow on temperature distributions but ignored the effective thermal conductivity of the perfused tissue (i.e., the enhancement of thermal conduction in the tissue due to blood perfusion within living tissues) [24]. These researchers also found that pulsatile blood flow with large pulsation amplitudes may exhibit a downstream two-peak behavior in the thermal dose contour in middle-sized blood vessels with diameters between 0.6 and 1 mm. Thus, pulsatile blood flow is an important factor in thermal therapies.

The concept of the effective thermal conductivity of tissue, which is used to describe the enhancement in the thermal conductivity due to the microenvironment blood perfusion in living tissues, has been used by some investigators in the field of thermal physiology [2529]. These authors have noted that blood perfusion can affect the thermal conductivity of the tissue. Crezee and Lagendijk experimentally demonstrated the relationship among the effective thermal conductivity of tissue, the thermal conductivity of tissue, and the blood perfusion rate [30]. Here, we not only consider this influence of the effective tissue thermal

conductivity but also incorporate this factor into the wave bio-heat transfer model.

The thermal relaxation time of a biological tissue is widely discussed by several studies [3135]. The thermal relaxation time of biological tissue can describe the response between the heat flux and the temperature gradient. Furthermore, the thermal relaxation time illustrates the time lag between the heat flux and the temperature gradient. Wang and Fan suggested that the heterogeneous and non-isotropic nature of a biological tissue normally yields a strong non-instantaneous response between the heat flux and the temperature gradient in non-equilibrium heat transport [33]. Mitra and colleagues measured the thermal relaxation time τ of processed meat (Bologna) and reported that τ _{was approximately} 16 seconds [36]. Because the heat conduction term in the Pennes model (see Eq. (1)) is based on Fourier’s law of heat conduction (Eq. (2)), we used a modified unsteady heat conduction equation (Eq. (3)), which was formulated by Cattaneo [37] and Vernotte [38], to replace the heat conduction term in Fourier’s law as follows:

ρ_tc_t∂ T

∂ t =∇∙

[

]

(

)

Fourier’s law of heat conduction is the following:

The modified unsteady heat conduction is the following:

q ( ⃗r , t )+τ∂ q ( ⃗r , t )

∂ t =−kt∇ T ( ⃗r , t) (3)

where τ _{is the thermal relaxation time.}

Roetzel et al. [39] experimentally showed that the thermal relaxation time τ was approximately 1.77 s in inhomogeneous materials with a hyperbolic thermal propagation behavior. Zhang found that the dual-phase lag phenomenon is more pronounced in a large blood vessel [31]. Based on his study [31], we used a thermal relaxation time ranging from 0.464 to 6.825 s in our numerical simulations. In addition, Shih and coworkers demonstrated that the thermal wave characteristics (i.e., the thermal relaxation time of biological tissues) cause a delay in the appearance of the peak temperature during thermal therapies [34]. These studies show that the thermal relaxation time can affect the temperature distribution under different heating conditions. Therefore, it is important to understand the coupled effects of effective tissue thermal conductivity, the thermal wave characteristics, and the pulsatile blood flow during thermal treatments.

The physical model used in this study is illustrated in Fig. 1. As presented in Fig. 1, we consider a cylindrical cancer tumor tissue embedded in a healthy blood perfused tissue with an axisymmetric geometric configuration. The tumor tissue is localized close to a thermally significant blood vessel (i.e., larger than 200 µm in diameter). Heating a tumor tissue adjacent to a thermally significant blood vessel with pulsatile blood flow can be modeled by a conjugate problem containing a blood-perfused solid-volume (i.e., a tumor and its surrounding healthy tissue) and a liquid-volume (i.e., a blood vessel). This model contains three components: a solid tumor, the healthy blood perfused tissue, and a rigid thermally significant blood vessel with pulsatile blood flow. In our numerical simulations, we also consider the thermal properties of the tissue, such as the effective thermal conductivity of the tissue due to blood perfusion, the thermal relaxation time of the blood-perfused solid tissue due to the finite speed in the living tissue, and the pulsatile blood flow pattern in the thermally significant blood vessel due to the heart beat.

For the blood-perfused solid-volume tissues, we not only used the wave bio-heat transfer equation, i.e., the Pennes bio-heat equation modified with the Cattaneo-Vernotte formula, to simulate the heat transfer in a tumor tissue and its surrounding healthy solid tissue [34, 37, 38] but also incorporated the effective thermal conductivity of the tissue to replace the tissue conductivity due to the enhancement of the thermal conductivity induced by blood perfusion

in living tissues [25, 30]. For the liquid-volume blood vessel, the energy transport equation was employed, and a periodic pulsatile blood flow pattern was also considered.

2.1 Velocity profile of pulsatile blood flow in a circular rigid blood vessel

Pulsatile blood flow emanates from the heart and travels through the arteries. In this study, the pulsatile flow model involves the assumptions that the blood vessel segment is straight, the blood vessel wall is rigid and impermeable, and the pulsatile blood flow is an incompressible Newtonian fluid [24]. The pulsatile blood flow was described by the pulsating frequencies and amplitudes [24]. The axial HagenPoiseuille parabolic velocity profile of pulsatile blood flow passing through an axisymmetric rigid vessel with an inner radius r0 is as follows [40]: w (r , t )=−r0 2 −r2 4 μ dp dz (4)

In this equation, the symbol μ is the dynamic viscosity of blood. Moreover, the pulsatile blood flow that emanated from the heart is described by a variation in the time period.

which is the pressure drop along the blood vessel) in the blood flowing direction does not remain constant and is described by an additional sinusoidal component in time. For simplicity, the form of the pressure gradient was assumed to be a simple harmonic motion in this study. The corresponding pressure gradient along the z-axis of a blood vessel was modeled as follows:

∂ p

∂ z=c0+c1sin ωt=c0+c1sin (2 π f )t (5)

where ω is the angular frequency. Moreover, we defined the coefficient fac=c_c1 0 to

characterize the relative pulsation intensity in the blood flow. In other words, the parameter

fac was used to describe the magnitude of the pulsatile blood flow in blood vessels due to

the rhythmic nature of the heartbeat. The axial velocity profile W (r ,t ) can be rearranged to obtain the following equation:

W (r ,t )=r 2 −r0 2 4 μ c0+ r2−r0 2 4 μ c1sinωt (6)

The boundary conditions are axisymmetric at the center and no-slip on the blood vessel wall (i.e., at the radius r0 ):

∂W (r , t )

∂r =0 whenr=0, (7a) W (r ,t )=0 whenr=r_{0 (7b)}

Using the above velocity profile of pulsatile blood flow, the average volume flow rate Q in a blood vessel over the time period ´T and the average blood velocity ´W are shown in Eqs. (8) and (9), respectively.

Q=1 T

∫

∫

From Eq. (9), we obtained the coefficient c0=

−8 μ W

r02 and rewrote the fractional

coefficient fac= c1 c₀= c1 −8 μ W r₀2 =−c1r0 2

8 μ W , which represents the relative intensity of the

pulsatile blood flow. Note that the diameters of thermally significant blood vessels and their associated average velocities were employed and are listed in Table 1 [19]. Finally, the periodic velocity profile of pulsatile blood flow in a rigid blood vessel can be obtained as

follows [40]: W (r ,t )=2W

(

)

{

[

(

)

(

)

]

}

√

denotes the Womersley number, ρ is the density of blood, ℜ{}

represents the real part of a complex data structure, and J0 is the Bessel function of the first kind of order zero. The Womersley number α is a dimensionless expression of the pulsatile flow frequency in relation to viscous effects. As the Womersley number increases, the velocity profile may appear as two peaks [24, 40, 41]. At a Womersley number of approximately 2.568, a previous study found a two-peak velocity profile in a large blood vessel [24]. In this study, the diameter of thermally significant blood vessels was ranged from 1 to 2 mm, and the heart beat frequency was varied from 1 to 3 Hz [42].

2.2 Temperature governing equations

The governing equations for the temperature field are shown in equations (11) for solid tissue (i.e., for a solid tumor tissue and for a solid healthy perfused tissue) and (12) for blood

flow in cylindrical coordinates [24]. Note that the tissue metabolic heat production Qm

was ignored because it is markedly smaller than the heating power in this study.

ρ_sc_s∂Ts ∂ t =ks

[

(

)

]

(

)

[

]

[

(

)

]

where ρ , c , and k are the density, the specific heat, and the thermal conductivity of

solid blood-perfused tissues, Ts represents the temperature of solid tissues, wb _{is the}

blood perfusion rate in solid tissues, Ta is the arterial temperature, which was set to 37 °C,

Q (r , z ,t ) is the heating power, W (r ,t ) is the axial velocity of the pulsatile blood flow,

and the subscripts s and b represent the solid tissue and blood, respectively.

Considering the finite propagation speed in living solid tissue (i.e., the effect of the thermal relaxation time), we substituted the heat conduction term in Eq. (11) with the modified unsteady heat conduction term shown in Eq. (3), rearranged the equation, and acquired the following thermal wave bio-heat equation for solid tissues (Eq. (13)):

ρ_sc_s

(

)

[

(

)

]

[

]

₍

₎

2.3 Effective thermal conductivity equation

For solid blood-perfused tumor tissues, we considered the effective thermal conductivity of the tissue obtained due to the thermal enhancement in the thermal conductivity of solid tissues by the tissue blood perfusion. Based on to the experimental data reported by Crezee and Lagendijk [25, 30], we used the scalar effective thermal conductivity equation (Eq. (14)) proposed by Crezee and Lagendijk:

k_eff=k_s

₍

₎

We used the effective thermal conductivity of the solid tumor tissue to examine the influence of the tumor blood perfusion rate on the temperature distribution. By replacing the thermal conductivity term of Eq. (13) by the effective thermal conductivity term of solid tissue, Eq. (13) was rewritten as follows:

ρ_sc_s

(

)

[

(

)

]

[

]

₍

₎

where the effective thermal conductivity of blood-perfused solid tissue keff=ks

(

)

and the parameter β is equal to 0.02 kg1 _m3 _{s [25]. In this equation, the blood perfusion} rate of the solid tumor tissue ranges from 0.5 to 20 kg m3 _s1 _{[43], and the blood perfusion} rate of the surrounding normal solid tissue was 0.5 kg m3 _s1 _{in the numerical simulation.} Furthermore, the initial conditions for the blood vessel and the tissue are shown in Eq. (16). At the interface between the blood vessel and tissue, temperature and heat flux continuity conditions were imposed, as shown in Eqs. (17) and (18), respectively.

T_s(r , z ,0 )=T_b(r , z , 0)=37℃ (16) T_s=T_{b at  (17)} k_eff ∂ Ts ∂ n =kb ∂T_b ∂ n at  (18)

where  denotes the interface between the blood vessel, the solid tumor, and the solid healthy tissue, and n indicates the direction normal to . At r=0 , the r-axis-symmetry condition for a blood vessel was applied:

∂Tb

∂ r =0. (19)

The boundary conditions at r=rmax , z=0 , and z=zmax were all equal to 37 °C.

T_t=T_b=37℃ . ₍₂₀₎

In this study, the convective boundary condition of the blood flow at z=zmax was imposed

as follows:

∂T_b

∂t +W (r ,t ) ∂ T_b

∂ z =0 , z=zmax (21)

We previously described an approach to prescribe the boundary and interface conditions for simulations of pulsatile blood flow [24]. First, we solved Eqs. (12), (15), and (16) to (21) employing the method of lines (MOL) and then constructed a discrete form of the temperature governing Eqs. (12) and (15) using the multi-block Chebyshev pseudospectral method and the boundary and interface conditions shown by (16) to (21) in space into a semi-discrete system in time [24,44]. This coupled system consists of ordinary differential equations (ODEs) in time, which were mainly derived from Eqs. (12) and (15), and algebraic equations from the boundary and interface conditions (16) to (21). Using the implicit ODE

solver ODE15s in the MATLAB (MathWorks, Natick, Massachusetts, U.S.A.) mathematical computing software, this coupled system of differential-algebraic equations (DAEs) was solved [24, 44]. The blood vessel parameters used in the simulations are listed in Table 1, and the heating schemes used in this study are shown in Table 2. Moreover, for example, the delivery of the power density distribution used in heating case III is shown in Fig. 2.

3. RESULTS AND DISCUSSION

Fig. 3 demonstrates the development of the temperature distributions on the r-z plane in

heating case II (i.e., the heating power density Qt=Qb=50 W cm3, and the heating

duration th=2 s) with the thermal relaxation time τ =6 .825 s. Figs. 3(a) through 3(b)

show that the temperature increased during the heating time period of 0 to 2 s. Even after the heating power was turned off, the temperature continued to increase until it reached the peak temperature of approximately 62.609 °C at approximately t=5.208 s, as shown in Fig.

3(d). The finite propagation speed of heat transfer in living tissues explains why the peak

temperature occurred after the heating power was turned off. Due to heat dissipation by blood perfusion, heat conduction by tissues, and heat convection cooling by pulsatile blood flow, the temperature distribution then continuously decayed with time until it became flat, as shown in Figs. 3(e), 3(f), 5(g), and 3(h) for times 10, 20, 30, and 60 s, respectively. For

heating scheme III and a heartbeat frequency of 1 Hz, peak temperatures of approximately 62.268 °C, 62.401 °C, and 62.837 °C were obtained at times t=4.008 s, t=5.443 _s, and t=8.443 s for thermal relaxations τ =0.464 s, τ =1.756 s, and τ =6.825 s, respectively (as shown in Table 3 and Table 4). A larger thermal relaxation time resulted in a more postponed time at which a higher peak temperature was obtained because the thermal relaxation time of tissue leads to a finite propagation speed of heat transfer in tissues rather than an infinite speed [1, 31, 3335]. In other words, the heat flow in solid tissues does not start instantaneously; instead, heat travels gradually with a time lag after the application of a temperature gradient. The finite propagation speed of heat transfer in living tissues is

represented by the term

√

tctτ , where

k_t

is the thermal conductivity of the tissue,

ρ_{t is the density of the tissue,} c_t

is the specific heat of the tissue, and τ is the thermal relaxation time of solid tissues [1, 24, 34].

The peak temperature plays an important role in thermal treatments because the maximum (i.e., peak) temperature directly dominates the thermal dose levels [19, 24, 27, 34]. The maximum temperature decreases as the heating time period increased with a constant total heated power energy. For the same total deposited energy density of 100 J cm3_,

I and V were 63.109 and 61.232 °C, and their occurrence times were approximately

t=7.138 and 18.77 s, respectively. The different maximum temperatures obtained with

the different heating schemes are listed in Table 3. Due to the finite propagation speed of heat transfer in living tissues, the peak temperature increases significantly when the thermal relaxation increases. The data shown in Table 3 demonstrate that a larger thermal relaxation results in a higher peak temperature for the same total energy density. Moreover, the frequency of pulsatile blood flow does not appear to affect the maximum temperature in solid tumor tissues. Table 4 shows that the occurrence time of the peak temperature in heating cases I to III is not the end of the heating duration but rather occurs after heating due to the lagging response to the heating source for a finite propagation speed. For instance, the peak

temperature in heating scheme III (i.e., Qt=Qb=25 _{W cm}3, and the heating duration

t_h=4 s) occurred at approximately t=8.443 s for f =1 Hz, d=2 mm, and

τ =6.825 s. In this case, the peak temperature exhibits a time lag of 4.443 s. However,

considering the influence of the thermal relaxation time ( τ =0.464 s, τ =1.756 s, and

τ =6.825 s), the peak temperature occurred before the end of the heating period for the

slow heating scheme case V (i.e., Qt=Qb=5 _{W cm}3, and the heating duration th=20

influence on the peak temperature and its occurrence time, as shown in Tables 3 and 4.

The data shown in Table 5 demonstrate that the tumor blood perfusion rate and the heating scheme altered the level of the peak temperature and its occurrence time. For rapid heating (i.e., heating scheme I in Table 2), the time of occurrence of the peak temperature was retarded due to the thermal relaxation time (see cases #1, #6, and #11), and a higher tumor blood perfusion also shortened the retarded time, as shown in cases # 1 and #11. The occurrence time of the peak temperature was retarded by approximately 0.17 s due to the absence of tumor blood perfusion, as shown in cases #1 and #11 in Table 5. Furthermore, the heating scheme significantly affects the occurrence time of the peak temperature. A slower heating scheme results in a lower peak temperature. For instance, the peak temperatures in the tumor regions in cases # 6 and # 10 with heating scheme cases I and V were 62.558 °C and 53.082 °C, respectively. This phenomenon can be explained by the finding that for the same heating energy a slower heating results in a lower peak temperature due to a longer time for the cooling effect of heat conduction and heat sink (i.e., the tumor blood perfusion term). The occurrence time of the peak temperature was prior to the end of the heading period in slow heating cases, such as cases #5, #10, and #15.

The peak temperature in the tumor region decreases as the blood vessel increases in diameter, as shown in Table 6. In other words, the diameter of the blood vessel increases as the peak temperature decreases. For instance, the peak temperatures in the tumor region were

62.236 °C and 61.866 °C in blood vessels with diameters of 1 and 2 mm, respectively (cases #1 and #3 in Table 6). Furthermore, the peak temperature was significantly lower with a higher blood perfusion rate and a slower heating scheme. For example, the peak temperatures in cases #3 and #24 were 61.866 °C and 54.791 °C, respectively. The difference in the peak temperatures between cases # 3 and #24 was 7.075 °C. Moreover, a blood vessel with a larger diameter exhibits a lower peak temperature and an earlier occurrence time of the peak temperature. For example, peak temperatures of 62.236, 62.136 °C and 61.866 °C at occurrence times of 3.981, 3.951, and 3.891 s were obtained in cases #1, #2, and #3 (Table

6), respectively. The influence of the diameter of the blood vessel on the peak temperature is

apparent. The comparison of blood vessels with diameters of 1 mm and 2 mm under heating scheme III (Table 2) revealed differences in the peak temperatures of 0.370 °C and 0.567 °C between cases #1 and #3 and between cases #10 and #12, respectively. If the heating scheme is slow, the difference in the peak temperature in blood vessels with different diameters increases. For instance, the temperature differences between cases #13 and #15 and between cases #22 and #24 were 1.937 °C and 1.549 °C, respectively. This analysis of the influence of the effective thermal conductivity of tumor tissues (i.e., keff=ks

(

)

showed that the effective thermal conductivity of the tumor tissue significantly affects the peak temperature if the blood perfusion rate of the tumor tissue is high, the diameter of the blood vessel is large, and the heating scheme is slow, as shown in cases #15 and 24.

Figure 4 shows that the increase in the temperature profile was higher and steeper in a

smaller blood vessel. For instance, for heating scheme IV, the temperature profiles inside blood vessels with diameters of 1, 1.4, and 2 mm were 53.731, 47.821, and 41.356 °C, respectively. Furthermore, the peak temperature was shifted downstream due to the blood flow. The distance between the location of the peak temperature and the heating region (i.e., the heating range in the zaxis was located between 5 and 15 mm, as shown in Fig. 1) increased with an increase in the diameter. Moreover, the temperature approximately 45 mm downstream of the z-axis in a blood vessel with a dimeter of 2 mm was higher than that obtained in blood vessels with diameters of 1 and 1.4 mm.

4. CONCLUSIONS

The present work demonstrates a numerical analysis of the coupled effects of effective tissue thermal conductivity, thermal wave characteristics, and pulsatile blood flow on temperature distributions under thermal treatments. This coupled model can be used to predict a quantitative analysis of the temperature in blood-perfused tissue. The frequency of the pulsatile blood flow due to the heartbeat has a small influence on the temperature distribution during thermal therapy. For a slow heating scheme, the effective tissue thermal conductivity of the tumor tissue significantly affects the peak temperature, particularly for a

higher blood perfusion rate of tumor tissue and a larger blood vessel. In addition, a larger thermal relaxation time affects the temperature distribution and postpones the occurrence of the peak temperature during thermal treatments.

ACKNOWLEDGMENTS

This work was supported by the National Science Council of Taiwan under grant NSC 1002221E039002MY3 and by China Medical University under grant CMU 101-S-05(101514C*).

REFERENCES

[1]F. Xu, S. Moon, X. Zhang, L. Shao, Y.S. Song, U. Demirci, Multiscale heat and mass transfer modeling of cell and tissue cryopreservation, The Philosophical Transactions of The Royal Society A Mathematical, Physical & Engineering Sciences 368 (2010) 561583.

[2]J. Choi, J.C. Bischof, Review of biomaterial thermal property measurements in the cryogenic regime and their use for prediction of equilibrium and nonequilibrium freezing applications in cryobiology, Cryobiology 60 (2010) 5270.

[3]Y. Rabin, Key issues in bioheat transfer simulations for the application of cryosurgery planning, Cryosurgery 56 (2008) 248250.

[4]J. Fan, L. Wang, A general bioheat at macroscale, International Journal of Heat and Mass 54 (2011) 722726.

[5]P.K. Gupta, J. Singh, K.N. Rai, Numerical simulation for heat transfer in tissues during thermal therapy, Journal of Thermal Biology 35 (2010) 295301.

[6]Y. Zhang, Generalized dualphase lag bioheat equations based on nonequilibrium heat transfer in living biological tissues, International Journal of Heat and Mass Transfer 52 (2009) 48294834.

[7]H.H. Pennes, Analysis of tissue and arterial blood temperature in the resting human forearm, Journal of Applied Physiology 1 (1948) 93122.

[8]H. Arkin, and K.R. Holmes, Recent developments in modeling heat transfer in blood perfused tissues, IEEE Transactions on Biomedical Engineering 41 (1994) 97107.

[9]A.-R. A. Khaled, K. Vafai, The role of porous media in modeling flow and heat transfer in biological tissues, Int. J. Heat Mass Transfer 46 (2003) 49895003.

[10]R.B. Roemer, Engineering aspects of hyperthermia therapy, Annu. Rev. Biomed. Eng. 1 (1999) 347376.

[11]P. Wust, B. Hildebrandt, G. Sreenivasa, B. Rau, J. Gellermann, H. Riess, et al., Hyperthermia in combined treatment of cancer, Lancet Oncol. 3 (2002) 487497.

[12]W.C. Dewey, Arrhenius relationships from the molecule and cell to the clinic, Int. J. Hyperthermia 25 (2009) 320.

[13]P.K. Gupta, J. Singh, K.N. Rai, Numerical simulation for heat transfer in tissues during thermal therapy, Journal of Thermal Biology 35 (2010) 295301.

[14]B. Decadt, A.K. Siriwardena, Radiofrequency ablation of liver tumours: systematic review, Lancet Oncol. 5 (2004) 550560.

[15]N.T. Wright, J.D. Humphrey, Denaturation of collagen via heating: an irreversible rate process, Annu. Rev. Biomed. Eng. 4 (2002) 109128.

[16]C.J. Diederich, Thermal ablation and high-temperature thermal therapy: overview of technology and clinical implementation, Int. J. Hyperthermia 21 (2005) 745753.

[17]J. Crezee, J.J.W. Lagendijk, Temperature uniformity during hyperthermia: the impact of large vessels, Phys. Med. Biol. 37 (1992) 13211337.

[18]M.C. Kolios, M.D. Sherar, J.W. Hunt, Blood flow cooling and ultrasonic lesion formation, Med. Phys. 23 (1996) 12871298.

[19]T.C. Shih, H.S. Kou, W.L. Lin, The impact of thermally significant blood vessels in perfused tumor tissue on thermal distributions during thermal therapies, Int. Commun. Heat Mass Transfer 30 (2003) 975985.

[20]H.W. Huang, C.T. Liauh, T.C. Shih, T.L. Horng, W.L. Lin, Significance of blood vessels in optimization of absorbed power and temperature distributions during hyperthermia, International Journal of Heat and Mass Transfer 53 (2010) 56515662.

[21]J.R. Womersley, Method for the calculation of velocity, rate of flow and viscous drag in arteries when the pressure gradient is known, J. Physiol. 127 (1955) 553563.

[22]K. Rohlf, and G. Tenti, The role of the Womersley number in pulsatile blood flow: a theoretical study of the Casson model, Journal of Biomechanics 34 (2001) 141148. [23]O.I. Craciunescu, S.T. Clegg, Pulsatile blood flow effects on temperature distribution and

heat transfer in rigid vessels, ASME J. Biomech. Eng. 123 (2001) 500505.

[24]T.L. Horng, W.L. Lin, C.T. Liang, T.C. Shih, Effects of pulsatile blood flow in large vessels on thermal dose distribution during thermal therapy, Medical Physics 34 (2007) 13121320.

[25]J. Crezee, J.J.W. Lagendijk, Temperature uniformity during hyperthermia: the impact of large vessels, Phys. Med. Biol. 37 (1992) 13211337

Radiofrequency tumor ablation: Insight into improved efficacy using computer modeling, American Journal of Roentgenology 184 (2005) 13471352.

[27]T.C. Shih, H.S. Kou, and W.L. Lin, Effect of effective tissue conductivity on thermal dose distributions of living tissue with directional blood flow during thermal therapy, International Communications in Heat and Mass Transfer 29 (2002) 115126.

[28]S.I. Alekseev, A.A. Radzievsky, I. Szabo, M.C. Ziskin, Local heating of human skin by millimeter waves: effect of blood flow, Bioelectromagnetics 26 (2005) 489501.

[29]A. Shitzer, R.C. Eberhart, Heat Transfer in Medicine and Biology: Analysis and Applications Volume II, Plenum Press, New York, 1985.

[30]J. Crezee, J.J.W. Lagendijk, Experimental verification of bioheat transfer theories: measurement of temperature profiles around large artificial vessels in perfused tissue, Physics in Medicine and Biology 35 (1990) 905923.

[31]Y. Zhang, Generalized dual-phase lag bioheat equations based on nonequilibrium heat transfer in living biological tissues, International Journal of Heat and Mass Transfer 52 (2009) 48294834.

[32]A. Narasimhan, K.K. Jha, L. Gopal, Transient simulations of heat transfer in human eye undergoing laser surgery, International Journal of Heat and Mass Transfer 53 (2010) 482490.

[33]L. Wang, J. Fan, Modeling bioheat transport at macroscale, ASME Journal of Heat Transfer 133 (2011) 011010110.

[34]T.C. Shih, H.S. Kou, C.T. Liauh, W.L. Lin, The impact of thermal wave characteristics on thermal dose distribution during thermal therapy: A numerical study, Medical Physics 32 (2005) 30293036.

[35]J.A. López Molina, M.J. Rivera, M. Trujillo, Thermal modeling for pulsed radiofrequency ablation: Analytical study based on hyperbolic heat conduction, Med. Phys. 36 (2009) 11121119.

[36]K. Mitra, S. Kumar, A. Vedavarz, M.K. Moallemin, Experimental evidence of hyperbolic heat conduction in processed meat, ASME Trans. J. Heat Transfer 117 (1995) 568573. [37]C. Cattaneo, Sur une forme de l’équation de la chaleur éliminant le paradoxe d’une

propagation instantaneé, Comptes Rendus 247 (1958) 431433.

[38]P. Vernotte, Les paradoxes de la théorie continue de l’équation de la chaleur, Comptes Rendus 246 (1958) 31543155.

[39]W. Roetzel, N. Putra, S.K. Das, Experiment and analysis for non-Fourier conduction in materials with non-homogeneous inner structure, Int. J. Thermal Sci. 42 (2003) 541552. [40]Y.C. Fung, Biomechanics: Motion, Flow, Stress, and Growth (SpringerVerlag, New

York, 1996)

[42]Y. Huo, G.S. Kasab, Pulsatile blood flow in the entire coronary arterial tree: theory and experiment, Am. J. Physiol. Heart Circ. Physiol. 291 (2006) H1074H1087.

[43]D. Fuentes, R. Cardan, R.J. Stafford, J. Yung, G.D. Dodd III, Y. Feng, High-fidelity computer models for prospective treatment planning of radiofrequency ablation with in vitro experimental correlation, Journal of Vascular and Interventional Radiology 21 (2010) 17251732.

[44]T.C. Shih, T.L. Horng, H.W. Huang, K.C. Ju, T.C. Huang, P.Y. Chen, Y.J. Ho, and W.L. Lin, “Numerical analysis of coupled effects of pulsatile blood flow and thermal relaxation time during thermal therapy,” International Journal of Heat and Mass Transfer 55 (2012) 37633773.

Table and figure captions

Table 1

Blood vessel parameters used in the simulations [19].

Table 2

Parameters of the five different heating schemes used in the simulations.

Table 3

Maximum temperatures (in °C) for different thermal relaxation times and heating schemes at frequencies of 1 Hz and 2 Hz under the following conditions: wb=0.5 kg m3 s1,

fac=0.5 , and d=2 mm.

Table 4

relaxation times and heating schemes at frequencies of 1 and 2 Hz under the following conditions: wb=0.5 kg m3 s1, fac=0.5 , and d=2 mm.

Table 5

Effects of the tumor blood perfusion rate of solid tumor tissue on the maximum temperature and its occurrence time under the following conditions: f =1 Hz, wb=0.5 kg m3 s1,

α=1.482 _, d=1.4 _{mm, and} τ =1.756 _s.

Table 6

Effects of the blood vessel diameter, heating scheme, and tumor blood perfusion rate on the peak temperature and its occurrence time under the following conditions: f =1 Hz,

fac=0.5 , wb=0.5 kg m3 s1, and τ =1.756 s .

Figure 1

Geometric configuration of pulsatile blood flow in blood-perfused tissue solid tumor tissue. A velocity profile of the pulsatile blood flow in a blood vessel was indicated. The treatment target (i.e., the heating target) was specified as z1≤ z ≤ z2 , 0 ≤ r ≤ r1 , and z1=5 mm,

z₂=15 mm, and r₁=5 mm were considered here. The radius of thermally significant blood vessels was denoted r0 . In this numerical study, the diameters of thermally significant blood vessels were 1, 1.4, and 2 mm.

Figure 2

Heating delivery of the power density distribution for heating case III. A heating power

density of Qt=Qb=25 W cm3 and a heating duration of th=4 _{s were used. The total}

energy density delivered was 100 J cm3_.

Figure 3

Temperature distribution evolution and temperature contours for heating scheme II and a thermal relaxation time of τ =6.825 s. (a) t=1 s; (b) t=2 s; (c) t=4 s; (d)

t=5.208 s; (e) t=10 s; (f) t=20 s; (g) t=30 s; and (h) t=60 s. The peak

temperature of 62.609 °C occurred at time t=5.208 s.

Figure 4

Peak temperature inside blood vessels with diameters of 1, 1.4, and 2 mm for heating scheme IV and a tumor blood perfusion rate of wbt=¿ 10 kg.

Table 1

Blood vessel parameters used in the simulations [19].

Diameter (mm) Average blood velocity in tumor ( w¯ ) (mm s1 )

1.0 8

1.4 10.5

2.0 20

Table 2

Parameters of the five different heating schemes used in the simulations.

Heating scheme I II III IV V

Heating power density Q (W cm3_{) 100 50 25 10 5}

Heating duration th (s) 1 2 4 10 20

Total heated energy density (J cm3_{) 100 100 100 100 100}

Table 3

Maximum temperatures (in °C) for different thermal relaxation times and heating schemes at frequencies of 1 Hz and 2 Hz under the following conditions: wb=0.5 kg m3 s1, fac=0.5 , and d=2 mm.

Heating case I II III IV V

Frequency (Hz)

f =1 _{0 62.493 62.412 61.833 60.135 56.763} 0.464 62.587 62.347 62.268 60.865 57.724 1.756 62.815 62.669 62.401 61.827 58.557 6.825 63.107 63.043 62.837 62.366 61.233 f =2 _{0 62.094 62.058 61.573 59.263 55.382} 0.464 62.586 62.344 62.269 60.863 57.719 1.756 62.813 62.718 62.041 61.826 58.553 6.825 63.109 63.044 62.837 62.366 61.232 Table 4

Times (in seconds) at which the maximum temperatures occurred for different thermal relaxation times and heating schemes at frequencies of 1 and 2 Hz under the following conditions: 𝑤�=0.5 kg m3 _s1_{, 𝑓𝑎�=0.5, and �=2 mm.}

Heating case I II III IV V

Frequency (Hz)

Thermal relaxation time  (s) Occurrence time of the maximum temperature_(s)

f =1 _{0 1 2 4 10 20} 0.464 1.243 2.987 4.008 9.535 17.938 1.756 3.895 4.340 5.443 9.704 17.951 6.825 7.155 7.627 8.443 12.822 18.784 f =2 _{0 1 2 4 10 20} 0.464 1.233 2.993 4.003 9.560 17.936 1.756 2.064 4.337 5.430 9.733 17.952 6.825 7.138 7.616 8.476 12.791 18.770 Table 5

Effects of the tumor blood perfusion rate of solid tumor tissue on the maximum temperature and its occurrence time under the following conditions: f=1 Hz, 𝑤�=0.5 kg m3 _s1_{, �=1.482,} �=1.4 mm, and �=1.756 s.

Case # wbt (kg

m3_s1₎

Heating scheme Peak Temperature (°C) Occurrence time (s) 1 0 I 62.577 1.871 2 II 62.365 2.247 3 III 62.136 3.951 4 IV 57.809 8.803 5 V 53.149 16.578 6 0.5 I 62.558 1.862 7 II 62.349 2.226 8 III 62.108 3.947 9 IV 57.747 8.793 10 V 53.082 16.560 11 10 I 62.224 1.701 12 II 62.100 1.994 13 III 61.554 3.882 14 IV 56.657 8.625 15 V 51.930 16.220 Table 6

Effects of the blood vessel diameter, heating scheme, and tumor blood perfusion rate on the peak temperature and its occurrence time under the following conditions: f =1 Hz,

fac=0.5 _, w_b=0.5 _{kg m}3 _s1_{, and} τ =1.756 s . Case # Heating scheme w_bt (kg m3 _s1₎ Diameter (mm) Peak Temperature (°C) Occurrence time (s) 1 III 0 1 62.236 3.981 2 1.4 62.136 3.951 3 2 61.866 3.891 4 0.5 1 62.210 3.978 5 1.4 62.108 3.947

6 2 61.833 3.887 7 10 1 61.700 3.918 8 1.4 61.554 3.882 9 2 61.213 3.819 10 20 1 61.139 3.854 11 1.4 60.958 3.818 12 2 60.569 3.756 13 IV 0 1 58.697 8.951 14 1.4 57.809 8.803 15 2 56.760 8.686 16 0.5 1 58.628 8.940 17 1.4 57.747 8.793 18 2 56.702 8.679 19 10 1 57.425 8.745 20 1.4 56.657 8.625 21 2 55.694 8.566 22 20 1 56.340 8.580 23 1.4 55.669 8.486 24 2 54.791 8.468 Nomenclature c _{specific heat [J kg}1 _K1_] c_{0 coefficient in Eq. (5), c}0= −8 μ w r0 2 c_{1 coefficient in Eq. (5)}

d diameter of blood vessel [mm], d=2r0

f _{frequency in Eq. (5) [s}1_]

fac _{coefficient of the relative intensity of pulsation in a blood vessel,} fac=c1 c0

J₀

Bessel function of the first kind of order zero in Eq. (10) k _{thermal conductivity [W m}1 _K1_]

k_{b thermal conductivity of blood [W m}1 _K1_] k_{t thermal conductivity of tissue [W m}1 _K1_]

p _{pressure [kg m}2_] q _{heat flux [W m}2_]

r _{distance from the z-axis [mm]} r₀

radius of a blood vessel [mm]

r_{1 maximum radius of the heated target of tumor tissue [mm],} r₁=5 mm t _{time [s]}

t_{h heating duration [s]}

Q _{heating power density [W cm}3_]

Q _{average volume flow rate [m}3 _s1_{]; Eq. (8)} Q_m

rate of tissue metabolic heat generation in Eq. (1) [W cm3] T _{temperature [} K _]

T _{period of time [s] in Eq. (8)}

T_{a temperature of arterial blood in Eq. (1) [} K _], T_a=310 K T_{b temperature of blood [} K _]

T_{t temperature of solid tissue [} K _]

w _{axial Hagen-Poiseuille steady parabolic velocity [mm s}1_{] in Eq. (4)} W _{axial velocity [mm s}1_{] in Eq. (6)}

W averaged velocity [mm s1_{] in Eq. (9) and Table 1}

W_b

blood perfusion rate of solid tissue [kg m3 s1], Wb=0.5 kg m3 s1

w_{bt blood perfusion rate of solid tumor tissue [kg m}3 _s1_{], wb}

t=0 20 kg m3 s1

z _{distance along the z-axis [mm]} z₁

lower boundary limit of the heated solid tumor in the z direction [mm], z1=5 _mm z_{2 upper boundary limit of the heated solid tumor in the z direction [mm],} z₂=15

mm z_max

z_max=100 mm Greek symbols α _{Womersley number,} α= r₀

√

μ _{dynamic viscosity of blood [kg m}1 _s1_], μ=0.004 _{kg m}1 _s1 ρ _{density [kg m}3_]

ρ_b

density of blood [kg m3], ρb=1,050 kg m3 ρ_t

density of solid tissue/tumor tissue [kg m3], ρb=1,050 kg m3 τ _{thermal relaxation time [s]}

ω angular frequency of heart beating [s1_]

 interface between the blood vessel and tissue in Eqs. (17) and (18)

Subscripts

b _blood s _{solid tissue}

Fig. 1

05 15

50

100

0

1

5

10

0

5

10

15

20

25

z (mm)

r (mm)

P

ow

er

d

en

si

ty

(W

/c

m

)

(a

)

(b

)

(e

)

(f)

Fig. 3

(d