The solution of Pennes' bio-heat equation with a convection term and nonlinear specific heat capacity using Adomian decomposition

Pennes' bio-heat equation is the most widely used equation to analyze the heat transfer phenomenon associated with hyperthermia and cryoablation treatments of cancer. In this study, the semi-analytical and numerical solutions of Pennes' equation in a highly nonlinear form derived from renal cell carcinoma tissue's nonlinear specific heat capacity along with a freezing convection term were obtained and analyzed for the first time. Here, the governing equation was reduced to a lumped capacity form for simplification and exerted on a solid spherical renal tumor. In the following, two semi-analytical techniques, the adomian decomposition method (ADM) and the Akbari–Ganji's method (AGM) were evaluated in solving the governing ODE. The comparison revealed full conformity between ADM and AGM, in addition to an excellent agreement between the semi-analytical and the numerical results before the phase transition. The analysis highlighted a deviation between the semi-analytical and numerical results for the limited convergence of power-series-based semi-analytical methods throughout the phase change and beyond. For the investigated case with Dt=0.025,Biot=0.075\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$D_{t} = 0.025,\quad {\text{Biot}} = 0.075$$\end{document}, the convergence occurred while τ∈[0,0.7]\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau \in [0,0.7]$$\end{document} for both methods. Consequently, both semi-analytical techniques ADM and AGM, are applicable to find the solution of Pennes' bio-heat equation before the phase change, and there is no superiority in favor of one in accuracy. In contrast, numerical methods are reliable during the phase transition and after that. The analysis of the numerical solution showed that with the growth of the tumor, achieving the necrosis of the malignant tissue takes longer, and large tumors' temperature may not decrease to necrosis temperature. The interpretation of the numerical results indicated that cryoablation could be considered an effective thermal treatment for renal tumors with a diameter lower than 2.0 cm.


Introduction
Pennes' bio-heat equation is the most common mathematical model to formulate biological systems with heat transfer [1]. Due to the extensive application of this equation in the mathematical modeling of the bio-heat transfer phenomenon, many attempts have been made to solve this well-known equation in different forms using analytical and numerical methods. The exact classical linear Pennes' equation solutions were obtained in one-dimensional multiregion Cartesian and spherical geometries with temperature-invariant physiological parameters by Durkee et al. In contrast, Liu and Tu investigated the bio-heat transfer phenomenon based on Pennes' equation considering transient blood temperature [2,3]. In the following, Marin et al. did a finite element analysis (FEA) of the nonlinear Pennes' equation for skin tissue subjected to external thermal sources. Valente et al. implemented a computer simulation of hyperthermia with nanoparticles using the finite volume method (FVM). Ghazanfarian et al. discretized nonlinear Pennes' bio-heat equation, including the bio-heat transfer thermal wave model and the dual-phase-lag model using a mesh-free smoothed-particle hydrodynamics (SPH). Dehghani and Sabouri offered a spectral element method (SEM) for solving the Pennes' bio-heat equation through triangular and quadrilateral elements as a different numerical approach. Kalateh Bojdi and Askari Hemmat applied wavelet collection methods (WCMs) for solving Pennes' bio-heat 1 3 equation [4]- [8]. Hatami [20]- [22]. Zhang and Chauhan presented a cellular neural network (CNN) methodology as a real-time Pennes' equation to apply during the thermal ablation process [23]. A modified form of Pennes' equation to study the temperature distribution induced by any electrode array in an anisotropic tissue containing several nodules, including primary or metastatic with arbitrary shape, was proposed by Roca Oria et al. [24]. Zhao et al. created a two-level finite difference schema for one-dimensional Pennes's bio-heat equation. In the subsequent attempts, Luitel et al. extended the one-dimensional Pennes' equation by adding a protective clothing layer and derived the temperature profile at different time steps by employing the fully implicit finite difference method (FDM) [25,26]. The solution of nonlinear bio-heat transfer equation in living tissues under periodic heat flux in tissue surface using the dual-phase lagging (DPL) non-Fourier heat conduction model and the adomian decomposition method (ADM) was obtained by Ghasemi et al. [27,28]. Akbari et al. applied the Akbari-Ganji's method (AGM) to find the nonlinear governing equations solution of chemical reactors and also the arched beams with nonlinear vibration [29,30]. Finding the semi-analytical solutions of Pennes' bio-heat equation and validation by comparing them with the numerical outputs will help adapt and improve these techniques to obtain more accurate closed-form answers according to the nonlinear characteristics of malignant tissues. On the other hand, analyzing the numerical solutions of Pennes' bio-heat equation gives an excellent perspective of tumors' behavior during thermal treatments. This study evaluates the capability of the two wellknown semi-analytical techniques, the adomian decomposition method and the AGM, to obtain the solution of Pennes' bio-heat equation with a strong nonlinearity relevant to the nonlinear specific heat capacity of renal cell carcinoma tissue, in addition to finding the numerical solution and analyzing the results for the first time.

Material and methods
Hyperthermia treatment and cryoablation are two primary thermal methods to treat malignant tumors. Most developments of these two techniques are based on finding the solution of the Pennes' bio-heat transfer equation in different situations. The general form of Pennes' equation is [31]: T, t, t , c t , k, b , c b , p , T a , q m are the tissue temperature, time, the density of the tissue, the specific heat capacity of the tissue, the thermal conductivity of the tissue, the density of blood, the specific heat capacity of blood, the blood perfusion rate, core body temperature and the volumetric metabolic heat generation rate, respectively.
Here the equation is reduced to a lumped capacity form for simplification purposes [32], and a cooling convection term is also added, as shown below: V t , h, A, T ∞ are the volume of tissue, convection coefficient, the area subjected to convection and the temperature of the convective medium.
shows the exit of energy via blood flow passing through the tumor. −hA(T(t) − T ∞ ) presents the thermal energy exchange between the renal tumor and the surrounding environment with the temperature −160 • C,which ablates the malignant cells. Q m denotes the generated thermal energy by the metabolic reactions of a kidney. The purpose is to obtain T(t) , which represents the temperature of the tumor subjected to the freezing environment over time. A schema of the problem is shown in Fig. 1.
By applying these new nondimensional parameters: Equation (2) converts into the following nondimensional equation: Considering core body temperature as the initial temperature ( T 0 = T a (, Eq. (4) includes the following dimensionless initial condition: The exact solution of Eq. (4) regarding the initial condition is: Ponomarev and Pushkarev's experimental research demonstrated the temperature dependency of malignant renal tissue's specific heat capacity during cryoablation [33]. A two-term normal distribution with the following form was utilized to formulate this dependence and also cover the phase transition during this procedure: which c ∞ indicates the specific heat capacity of the renal tumor at T ∞ and a 1 , a 2 , b 1 , b 2 , c 1 , c 2 are constants. The numeric values of these constants are available in Table 1. The normal distribution associated with the specific capacity of renal cell carcinoma tumor is depicted in Fig. 2. It can be seen that near zero degrees Celsius, an intense increase in the value of specific heat capacity occurs, and it follows a nonlinear form. During a phase change, the temperature remains nearly constant at the freezing point, and all the exchanged thermal energy is consumed for the phase .
transition. After the phase transition termination, the specific heat capacity follows an almost linear function. By introducing the dimensionless groups below, the new form of (7) can be presented as (9).
Thus, the nondimensionalized form of the Pennes's bioheat transfer equation with a convective term and nonlinear temperature-dependent specific heat capacity associated with malignant renal tissue will be in the form: In sulation with the dimensionless initial condition. (0) = 0 = 0.
The numeric values of the physical parameters used to solve Pennes' equation exist in Table 2.
Here renal tumor is considered a solid sphere with the diameter D t and the characteristic length L c where: Because: and since the lumped analysis of a system is valid if only Biot < 0.1 [38], then: for the nondimensionalization of the governing Eq. (10), at: Therefore, the values of ( ) vary from 0 to 1 throughout the freezing process.

Adomian decomposition method (ADM)
Adomian decomposition method (ADM) is a practical semianalytical method based on infinite series to find the solution to nonlinear problems [39]. Exerting ADM on Eq. (10) as the governing equation of the problem leads [40]: operator on both sides of (10) and then applying the initial condition results: Assuming the solution in the form of: and expanding N( ( )) as an infinite series in the form of adomian polynomials leads: where: Substituting (17), (18) into (16) ends in finding the approximate solution of (10).

Akbari-Ganji's method (AGM)
Akbari-Ganji's method (AGM) is a semi-analytical technique to find the solution of linear and nonlinear problems based on finite series. In this method, a finite series with unknown coefficients is considered to approximate the solution of the governing differential equation [41]. Applying the initial condition and perhaps the boundary conditions to the supposed solution, to the governing equation, and its derivatives after substituting the solution into the differential equation, an algebraic equation system appears which determines the unknown coefficients [42]- [43].
In AGM, the solution of (10) is assumed in the form: Exerting the initial condition on (20) leads: By inserting (20) into (10) and then applying the initial value to it, the following equation will emerge: Substituting (21) into (22) will result 1 as: Through inserting (20) into (10) and deriving the outcome with respect to , and applying the initial condition and regarding the values of 0 and 1 , 2 will be obtained as: Repeating this procedure n − 2 times will compute the other unknown coefficients.

Numerical solution
The numerical methods that use the information at more than the last mesh point are referred to as multistep methods. The four-step Adams-Bashforth method and the modified Euler method as the starter technique were used to solve the governing Eq. (10) numerically. The solution was obtained based on the original variables (t, T(t)) for the various values of h and Biot number.

Adomian decomposition method
The first sixteen terms of the solution for D t = 0.025(m) and Biot = 0.075 using ADM are expressed as follows:  (25) in comparison to the numerical solution. As it can be seen, there is a perfect agreement between ADM and the numerical solution before the phase transition. The reason for the deviation after the phase change is the limited convergence domain of power-series-based methods, which is [0,0.7] for the investigated case here.

Numerical method
The first sixteen terms of the solution for D t = 0.025(m) and Biot = 0.075 applying AGM are denoted as below: such as ADM, and the convergence domain is [0,0.7] for the studied case.
As observed, in all the investigated cases, the temperature of the tumor falls from 37 • C(core body temperature) and approaches the convection medium's temperature −160 • C during cryoablation. However, the duration of this process increases with the growth of tumors' diameter. Consequently, larger tumors need more time to ensure their complete necrosis. Moreover, with the growth of the tumor's diameter, its temperature becomes stable at a higher value, and it does not reach the essential magnitude to guarantee necrosis. Therefore, it appears that the acceptable diameter to treat a renal cell carcinoma tumor through cryoablation is D t ≤ 0.02(m) . The numerical method was implemented for the tumor's diameters = 0.01, 0.015, 0.02, 0.025, 0.03 (m) and the Biot numbers = 0.025, 0.05, 0.075, 0.1. These numerical solutions are illustrated in Fig. 5.

Conclusions
Considering the solutions of the two semi-analytical methods, ADM and AGM, it is evident that the outputs of these two methods are entirely coincident. As observed in Figs. 3 and 4, the convergence domain of both methods is from the starting time to the phase transition. For D t = 0.025, Biot = 0.075 this interval is nearly ∈ [0, 0.7] . The semi-analytical methods and the numerical technique show excellent agreement in this domain. It approves that ADM and AGM are both suitable methods to solve the nonlinear form of the Pennes' bio-heat transfer equation before phase transition, while numerical techniques are still the appropriate procedures to find the solution during phase change and beyond. This is because of the limited convergence domain of power-series-based semi-analytical methods, which is a common characteristic. The convergence domain is limited to the phase transition for these two techniques. Comparing the numerical solutions determines that with the increase in the diameter of renal carcinoma tumors, approaching the suitable temperature for cryoablation takes longer. The results also confirm that cryoablation can effectively treat tumors with a diameter lower than 0.02(m). 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/.