Prediction of Temperature Field and Thermal Damage of Multilayer Skin Tissues Subjected to Time-Varying Laser Heating and Fluid Cooling by a Semianalytical Method

A semianalytical method based on the separation of variables is presented to solve the problem of heat transfer in multilayer skin tissues. )is method is used to predict the temperature field and thermal damage of skin tissues, which are subjected to the timevarying laser heating and fluid cooling. )e effect of convective heat transfer coefficient of the skin surface on temperature field and thermal damage is studied. )e temperature field and thermal damage of skin tissue with different laser heating cycles are calculated. )e relationship between the laser heating period and the temperature fluctuation is obtained. )e effects of the peak time of laser heating on temperature and the thermal damage of skin tissues are investigated. )e relationship between the peak time of laser and the temperature fluctuation of skin tissue is drawn. )e results calculated in this paper can provide some guidance for hyperthermia of skin tissue. )e presented semianalytical method can be applied to solve the general heat transfer problems of multilayer structures.


Introduction
ermotherapy has been widely used in the treatment of skin tissue diseases. Laser, microwave, and radio frequency are often used to irradiate the surface of skin tissue during thermotherapy [1][2][3]. Accurate prediction of skin tissue temperature and thermal damage during hyperthermia plays an important role in thermotherapy [4]. More and more researchers have studied the heat transfer and thermal damage of the skin tissue. Important research can be summarized as follows.
Deng and Liu gave several closed-form analytical solutions to the bioheat transfer problems which are based on Green's function method [5]. Xu et al. studied the heat transfer of skin tissue by using pence model; they established a mathematical model of heat transfer of three-layer skin tissue and solved the model by using finitely different method [6]. Liu et al. established DPL model for heat transfer of multilayer skin tissue; the Laplace transform method was used to solve this model [7]. Tian et al. established a skin-fabric five-layer heat transfer model; a semianalytical method was employed to solve this model [8]. Ahmadikia et al. used analytical methods based on Laplace transform to study the distribution of temperature field in skin tissue under constant or transient heat flux conditions; parabolic and hyperbolic skin tissue heat transfer equations were studied [9]. In the study of Kumar et al., a dual-phaselag model of bioheat transfer had been studied by using Gaussian distribution source term with the most generalized boundary condition during hyperthermia treatment [10]. Lin solved the analytical solution of heat transfer of single skin layer which was based on the separation of variable methods; Fourier and non-Fourier cases were studied [11]. Hooshmand et al. developed an analytical solution method for solving the problem of heat transfer in biological tissue, generalized dual-phase-lagging (DPL) model was studied, and the calculation results were compared with the classical DPL model [12]. Kumar et al. dealt with the study of heat transfer and thermal damage in triple layer skin tissue using fractional bioheat model. An implicit finite difference scheme is obtained by a new method [13]. Kumar and Rai studied the thermal characteristics in skin tissue using time fractional dual-phase-lag bioheat transfer (DPLBHT) model subjected to Dirichlet boundary condition. ey solved this problem using the finite element Legendre wavelet Galerkin method. ey treated skin tissue as a single layer [14]. Li et al. predicted the temperature of different biological tissues using Fourier and non-Fourier heat transfer models. And the temperature response of different tissues was measured.
e experimental results were compared with those predicted by the model [15]. Hobiny and Abbas used the Laplace transform method to solve the problem of biological tissues subjected to moving heat source [16]. Ghanmi and Abbas derived an analytical solution for the fractional transient heating within the skin tissue during hyperthermia [17]. Sur et al. used Laplace transform to study the thermal behavior of skin tissue when it was subjected to moving heat sources [18]. Sharma and Kumar established a nonlinear dual-phase-lag heat transfer model of biological tissue. e finite element Runge-Kutta method was used to solve this model. e temperature field and thermal damage of skin tissue were studied using the numerical solution of the model [19].
Skin tissue is the largest organ of the human body, including epidermis, dermis, and subcutaneous tissue [20]. ermophysical parameters of different layers of skin tissue are different. So it is difficult to solve the heat transfer model of skin tissue by analytical method [9,11,21,22]. Numerical methods were often used to solve this problem [23]. Several papers used analytical methods, but all of them treated skin tissue as a single layer [12]. In this paper, a semianalytical method is developed to solve heat transfer in multilayer skin tissues. A special method of separation of variables is used to derive the analytical solution. In the derivation process, the time domain is firstly divided into several small time periods. en the temperature function at the interface of different layers of skin tissue is assumed to be linear in a small period of time. After that, the analytical solution of separated variables of skin tissue heat transfer is deduced in a small period of time. According to the condition that the temperature at the interface is equal and the heat flux at the interface is continuous, the temperature at the interface is calculated. Temperature field of skin tissue during laser heating and fluid cooling is obtained by the semianalytical solution. e effects of several parameters such as convective heat transfer coefficient, period, and peak time of laser heating on temperature field are studied. ese parameters can be used to control the temperature of the skin tissue and the thermal damage in the appropriate range.

Mathematical Model
Human skin tissue consists of epidermis, dermis, and subcutaneous. e geometric model of skin tissue is shown in Figure 1. e corresponding mathematical model can be established by combining the boundary conditions. e governing equation can be expressed as (1) [2,[24][25][26][27]: e boundary conditions can be expressed as follows: e initial condition is e corresponding dimensionless governing equation is e dimensionless boundary conditions are as follows: At ξ � 1, e dimensionless initial condition is

Epidermis
Dermis Subcutaneous Laser heating x Equations (7) and (9)- (11) show that the boundary conditions on both sides of each layer of skin tissue are known in a short time interval [0 Δτ]. A extended shifting variable method is developed by Lee and Lin, which can be used to solve the temperature field of each layer of skin tissue [29]. Assume Here, g i,m (ξ), i � 1, 2, 3; m � 1, 2 are the shifting functions to be specified and v i (ξ, τ) is the transformed function. Substituting (12) into (5), the two following subsystems are obtained. e first subsystem is expressed in terms of the transforming variables v i (ξ, τ) as follows; the transform governing equation is If i � 1, substitute (12) into (6) and (10): In order to homogenize the boundary conditions when separating variables, (17) and (18) need to be satisfied: So, Assume If i � 2, 3, substitute (12) into (10) and (11): In order to homogenize the boundary conditions when separating variables, (24) and (25) need to be satisfied: Assume e characteristic equation of (13) is Assume Substitute (30) into (29): If i � 1, substitute equation (30) into (17) and (18):

Mathematical Problems in Engineering
Equation (31) is a second-order homogeneous ordinary differential equation. Its characteristic function is as follows: where λ j is the positive solution of the following equation: If i � 2, 3, substitute (30) into (24) and (25): where λ j is the positive solution of the following equation: Obviously, ψ i,j (ξ) satisfies the following equation: e solution of (41) is

Thermal Damage
In addition to predicting the temperature field of skin tissue during thermotherapy, it is often necessary to predict its thermal damage. Formula (47) is developed by Henriques and Moritz which can be used to calculate thermal damage of skin tissue [1]: where Ω represents the thermal damage factor of skin tissue. E a represents the activation energy. T represents the temperature of the skin tissue. R represents the universal gas constant. According to the thermal damage factors, the thermal damage is divided into three levels in the literatures. Ω � 0.53, the first level of thermal damage; Ω � 1, the second level of thermal damage; Ω � 10 4 , the third level of thermal damage.

Numerical Results
When the surface of skin tissue is irradiated by pulsed laser, the variation of laser energy with time can be shown in Figure 2. Assume the initial temperature of the skin tissue T . e thermophysical parameters of skin tissue can be shown in Table 1 [4]. Figure 3 shows the temperature history at the skin surface, ED interface, and DS interface in a laser heating cycle. It is found that the temperature rises first and then decreases with time at x � 0, x � 0.1 mm. If h � 7 W/m 2°C , the temperature peaks appear at t � 2 s and t � 2.4 s, respectively. If h � 700 W/m 2°C , the temperature peaks appear at t � 1.4 s and t � 1.7 s, respectively. So we can speculate that the larger the value of x, the shorter the time required for temperature peak. If h � 7 W/m 2°C , at x � 1.6 mm, the temperature first remains constant and then continues to rise. When t � 10 s, the peak value of temperature still does not appear. If h � 700 W/m 2°C , at x � 1.6 mm, because this location is far away from the heating source, temperature varies a little with time. It can be inferred from the comparison of solid and dotted lines that the larger the value of h, the lower the temperature. Investigating its reason we can know that the larger the value of h is, the more heat is taken away by water cooling on the surface of skin tissue, so the lower the temperature of skin tissue is. Figure 4 shows the temperature distribution in the x direction at different times. As shown in the picture, if h � 7 W/m 2°C , t � 1 s, the temperature of the skin tissue surface is the highest, which is 65°C. e temperature distribution of skin tissue decreases in the x-axis. Comparing the temperature distribution curves of t � 0.5 s and t � 1 s, it is found that the temperature of T � 1 s is higher than that of t � 1 s. is is because the heat flux of laser heating keeps rising in 0 to 1 s. Because the heat flux of laser heating begins to decrease after t � 1 s, the temperature distribution at t � 5 s and 10 s is much lower than that at t � 1 s. Comparing the temperature distribution of h � 7 W/m 2°C and h � 700 W/ m 2°C , it is found that, at the same time, the temperature distribution of skin tissue at h � 700 W/m 2°C is significantly lower than that at h � 7 W/m 2°C . Because the cooling effect is greater than the heating effect, the temperature of skin tissue surface is lower than 37°C at t � 5 s and t � 10 s. Figure 5 demonstrates the relationship between heat flux and heating time at different locations of skin tissue. As can be seen from the figure, when h � 7 W/m 2°C , the heat flux density on the surface of skin tissue (x � 0) increases first and then decreases, which is in accordance with the trend of surface heating heat flux. At the skin-dermis interface (x � 0.1 mm), the heat flux is slightly lower than that at the skin surface. At the epidermal-subcutaneous interface (x � 1.6 mm), the heat flux increases first and then decreases, but the range of variation is very small. On the left side of the skin tissue (x � 6 mm), the heat flux is 0 and remains unchanged. When h � 700 W/m 2°C , the heat flux of heating is less than that of cooling in the initial period of time. e heat flux of skin surface (x � 0) and dermal-epidermal interface (x � 0.1 mm) changed from negative to positive, continued to increase, and then gradually decreased. Because the boundary condition we set is that the right side of the skin tissue is adiabatic, on the left side of the skin tissue (x � 6 mm), the heat flux is 0 and remains unchanged too. Figure 6 demonstrates the relationship between thermal damage and heating time at different locations of skin tissue.     As can be seen from the figure, thermal damage of skin tissue increases with time accumulation. e surface of skin tissue has the greatest thermal damage. If h � 7 W/m 2°C , the thermal damage value of skin tissue surface is close to 100, reaching the second-degree thermal damage, and the thermal damage value at the junction of epidermis and dermis exceeds 0.53, reaching the first-degree thermal damage. Because the greater the value of h, the better the heat dissipation and the smaller the thermal damage of skin tissue, the thermal damage at h � 700 W/m 2°C is six orders lower than that at h � 7 W/m 2°C . It can be seen from the figure that although the heat flux of laser heating increases first and then decreases, the heat damage increases continuously. As the heat flux of laser heating decreases, the trend of increasing skin tissue thermal damage slows down. In order to analyze the effect of R and H on thermal damage, a variance analysis was carried out. Table 2 shows the calculated results. By calculation, R 2 (x) � 0.163812, R 2 (h) � 0.318886, R 2 (interaction) � 0.318877. erefore, it can be inferred that the influence level of x on thermal damage is lower than the average and the influence level of h on thermal damage is above average. ese two factors and their interaction have a great influence on thermal damage. Figure 7 demonstrates the effect of heating cycle Tc on the temperature history at skin surface. As can be seen from the figure, because heating and cooling are the same in the first five seconds, so the temperature curve coincides completely. If T C � 5 s, the first pulse of laser heating is completed after 5 seconds, the heat flux of heating no longer decreases, and the second laser heating pulse starts. erefore, the temperature on the surface of skin tissue starts to rise when t � 5 s, and the rising process is fluctuant, and the fluctuating period is T C � 5 s. After five seconds, the shorter the heating cycle is, the faster the temperature rises. Temperature fluctuates with heating period. After 30 seconds, the amplitude of temperature fluctuation tends to be balanced. e shorter the heating period T C is, the greater the temperature fluctuation amplitude T s,max is. Figure 8 shows the effect of heating cycle T C on the temperature history at the interface between epidermis and dermis. Comparing with Figure 7, it can be seen that the temperature fluctuation at the interface of epidermis and dermis is smaller than that at the skin surface. It is the same as the result shown in Figure 6, e shorter the heating period T C is, the greater the temperature fluctuation amplitude T s,max is. When T C � 5 s, the amplitude of fluctuation reaches the maximum value; at the junction of epidermis and dermis, T ed1,max � 82.5°C, it is 9.5°C lower than the surface temperature of skin tissue. When T C � 15 s, the amplitude of fluctuation reaches minimum value; at the junction of epidermis and dermis, T ed3,max � 65.2°C, it is 7.9°C lower than the surface temperature of skin tissue. Figure 9 shows the effect of heating cycle T C on the temperature history at the interface between dermis and subcutaneous. Because the location is far from the heat source, the range of temperature fluctuation is much smaller than the skin surface and the interface between epidermis and dermis. It takes a long time for temperature fluctuation to reach a steady state. As can be seen from the figure, when T C values are different, the time required for temperature fluctuation to reach stability is 75 s, 100 s, and 150 s, respectively. It can be inferred that the greater the value of temperature T C is, the less time it takes for temperature fluctuation to reach stability. It can also be seen from the figure that the larger T C value is, the larger temperature fluctuation range is. Figure 10 shows the effect of heating cycles on the thermal damage at different locations of skin tissue. Because the surface of skin tissue is close to the heat source, the thermal damage on the surface of skin tissue (x � 0) is much greater than the interface between dermis and epidermis (x � 0.1). With the accumulation of heating time, the thermal damage of skin surface in all three cases exceeds 10 4 , reaching the third-degree burn. e thermal damage with T C � 5 s is four orders higher than the thermal damage with T C � 10 s. It can be inferred from the figure that the smaller the heating period is, the more serious the thermal damage of skin tissue is. e thermal damage at the interface between epidermis (x � 0.1 mm) and dermis exceeds 10 4 , reaching the third-degree burn when T C � 5 s or T C � 10 s. When T C � 15 s, the thermal damage value at the interface between epidermis and dermis exceeds 100 to reach second-degree burn. Because the dermis-subcutaneous interface (x � 1.6 mm) is far away from the heat source, no burn occurs with T C � 10 s and T C � 15 s, and the first-degree burn occurs with T C � 5 s. In order to study the influence of T C and x on thermal damage, variance analysis was carried out. Table 3 shows the calculation results. By calculation, R 2 (x) � 0.165575, R 2 (T C ) � 0.167528, and R 2 (Interaction) � 0.3308943. erefore, it can be inferred that the influence level of x and h on thermal damage is lower than the average, respectively. ese two factors and their interaction have a significant effect on thermal damage.  Figure 11 indicates the relationship between laser heating cycle and temperature fluctuation amplitude. T s,max and T s,min represent the maximum and minimum values of temperature fluctuation amplitude of skin tissue surface, respectively. T ed,max and T ed,min represent the maximum and minimum values of temperature fluctuation amplitude of the interface between epidermis and dermis, respectively. T ds,max and T ed,min represent the maximum and minimum values of temperature fluctuation amplitude at the interface between dermis and subcutaneous, respectively. It can be seen from the figure that the amplitude of temperature fluctuation increases first and then decreases with the increase of laser heating period. e range of temperature fluctuation increases with the increase of laser heating period. Skin tissue surface temperature fluctuation range is the largest, the epidermis and dermis interface is the second, and the dermal and subcutaneous interface is the smallest. With the increase of laser period, the minimum amplitude of temperature fluctuation becomes closer and closer; when T C > 16 s, the lowest temperature of skin tissue is equal at all locations. Figure 12 shows the effect of time t p on skin surface temperature history. It can be seen from the figure that the skin tissue surface temperature fluctuates and rises at the   Figure 7: Effect of heating cycle T C on the temperature history at skin surface (t p � 1 s, q peak � 3 W/cm 2 , T ∞ � 20°C, x � 0).

Mathematical Problems in Engineering 7
beginning time and then tends to fluctuate steadily. If t p � 3 s, the fluctuation range of skin surface temperature is 80-107°C. If t p � 1 s, the fluctuation range of skin surface temperature is 45-75°C. If t p � 0.5 s, the fluctuation range of skin surface temperature is 42-63°C. We can infer that the larger the t p value is, the larger the maximum value of skin tissue surface temperature fluctuation is and the smaller the minimum value of skin tissue surface temperature fluctuation is. Figure 13 shows the effect of peak time t p on the temperature history at the interface between epidermis and dermis. Comparing with Figure 12, we can see that the temperature of this interface is lower than that of the skin surface. e trend of fluctuation is consistent with Figure 12. Figure 14 shows the effect of peak time t p on the temperature history at the interface between dermis and subcutaneous. Because this location is far away from the heat source, the temperature of this location is much lower than that of the skin tissue surface and the interface between the epidermis and dermis. e temperature at this position fluctuates with the heat flux heated by laser. Figure 15 indicates the relation between peak time t p and temperature fluctuation Amplitude. e meaning of T s,max , T s,min , T ed,max , and T ed,min is the same as shown above. It can be seen from the figure that, with the increase of peak time t p , the amplitude of temperature fluctuation becomes larger and larger at all positions of skin tissue. e temperature fluctuation range of skin surface and epidermis-dermis interface increases with the increase of peak time t p . e temperature fluctuation range of dermal-subcutaneous interface increases first and then remains constant with the increase of t p . e larger the t p value is, the smaller the range of temperature fluctuation is. Figure 16 indicates the effect of peak time t p on the thermal damage history at different locations of skin tissue. It can be seen from the figure that the greater the t p value, the more severe the thermal damage of the skin tissue. With t p � 3 s thermal damage value of skin surface and epidermal-dermal interface exceeds 10 4 , reaching the third-degree burn, and thermal damage value of dermalsubcutaneous interface exceeds 1 to reach a second-degree      Figure 13: Effect of peak time t p on the temperature history at the interface between epidermis and dermis (T C � 10 s, q peak � 3 W/cm 2 , burn. If t p � 0.5, the thermal damage on the skin surface (x � 0) exceeds 100 to reach a second-degree burn. e thermal damage on the skin surface (x � 0.1 mm) exceeds 1 to reach a first-degree burn. Because the dermal-subcutaneous interface is far away from the laser heating heat source, there is no burn. If t p � 1, the thermal damage on the surface of skin tissue exceeds 10 6 , reaching the thirddegree burn. e thermal damage on the epidermisdermis interface of skin tissue exceeds 10 6 , reaching the second-degree burn. e thermal damage on the dermalsubcutaneous interface of skin tissue is less than 0.53, lower than the first-level burn standard. In order to study the influence of t p and x on thermal damage, variance analysis was carried out. Table 4 shows the calculation results. By calculation, R 2 (x) � 0.165721, R 2 (t p ) � 0.167244, and R 2 (Interaction) � 0.3325. erefore, it can be inferred that the influence level of x and t p on thermal damage is lower than the average, respectively. ese two factors and their interaction have a significant effect on thermal damage.

Discussion
It is of great significance to accurately predict the temperature field and control the thermal damage during skin tissue hyperthermia. In this paper, a special separation variable method is used to derive the semianalytical solution of temperature field in skin tissue hyperthermia. Compared with the traditional finite element method, finite difference method, Laplace transform method, and so on, this method has clearer physical meaning and higher calculation efficiency. e temperature field and thermal damage of skin tissue during laser hyperthermia are predicted by the analytical solution derived in this paper. Several parameters that can affect the temperature field and thermal damage of skin tissue during heating have been studied. rough analysis we can find that the convective heat transfer coefficient between the skin tissue surface and the cooling fluid has a greater effect on the temperature field and thermal damage. We can change the injection speed of the cooling fluid to change the convection heat transfer coefficient and then  Figure 14: Effect of peak time t p on the temperature history at the interface between dermis and subcutaneous (T C � 10 s, q peak � 3 W/cm 2 , T ∞ � 20°C, x � 0.1 mm).
control the temperature and thermal damage of the skin tissue in an appropriate range. In addition, the heating period and peak time of the laser also have a greater impact on the skin tissue temperature field and thermal damage during hyperthermia. Adjusting the parameters such as laser heating period or peak time can control the skin tissue