Next Article in Journal
Identification of the Order of the Fractional Derivative for the Fractional Wave Equation
Next Article in Special Issue
An ε-Approximate Approach for Solving Variable-Order Fractional Differential Equations
Previous Article in Journal
Unique Solvability of the Initial-Value Problem for Fractional Functional Differential Equations—Pantograph-Type Model
Previous Article in Special Issue
Finite Element Approximations to Caputo–Hadamard Time-Fractional Diffusion Equation with Application in Parameter Identification
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Fractional Dual-Phase-Lag Model for Nonlinear Viscoelastic Soft Tissues

by
Mohamed Abdelsabour Fahmy
1,* and
Mohammed M. Almehmadi
2
1
Department of Basic Sciences, Adham University College, Umm Al-Qura University, Makkah 28653, Saudi Arabia
2
Department of Mathematical Sciences, Faculty of Applied Sciences, Umm Al-Qura University, Makkah 24381, Saudi Arabia
*
Author to whom correspondence should be addressed.
Fractal Fract. 2023, 7(1), 66; https://doi.org/10.3390/fractalfract7010066
Submission received: 4 November 2022 / Revised: 29 December 2022 / Accepted: 1 January 2023 / Published: 5 January 2023

Abstract

:
The primary goal of this paper is to create a new fractional boundary element method (BEM) model for bio-thermomechanical problems in functionally graded anisotropic (FGA) nonlinear viscoelastic soft tissues. The governing equations of bio-thermomechanical problems are briefly presented, including the fractional dual-phase-lag (DPL) bioheat model and Biot’s model. The more complex shapes of nonlinear viscoelastic soft tissues can be handled by the boundary element method, which also avoids the need for the interior domain to be discretized. The fractional dual-phase-lag bioheat equation was solved using the general boundary element method (GBEM) based on the local radial basis function collocation method (LRBFCM). The poroelastic fields are then calculated using the convolution quadrature boundary element method (CQBEM) The numerical findings show that our proposed numerical model is valid, efficient, and accurate.

1. Introduction

Arsene d’Arsonval and Claude Bernard demonstrated that the temperature difference between arterial and venous blood is due to blood oxygenation [1]. Pennes [2] was the first to use the blood flow effect to characterize temperature distribution in human tissue. To address a drawback of the linear bioheat transfer model, Cattaneo [3] and Vernotee [4] proposed a single-phase-lag (SPL) model associated with heat flow. When the SPL model is combined with the energy equation, it yields the thermal wave bioheat transfer (TWBHT) model. The dual-phase-lag (DPL) model was created by expanding the SPL model to include another phase lag caused by a temperature gradient [5,6]. Because analytical solutions to the current situation are extremely difficult to achieve, numerical methods have emerged as the primary tool for resolving these problems, as seen in research by Pei et al. [7], Ooi et al. [8], Zhou et al. [9], Ng et al. [10], Majchrzak and Turchan [11], Bottauscio et al. [12], Deng and Liu [13], and Partridge and Wrobel [14]. Dual-phase-lag thermoelastic models have been developed to study the effects of rotation [15], gravity [16], laser pulse [17], and temperature dependence [18]. Laplace domain fundamental solutions are generally easier to obtain than time domain fundamental solutions for scientific and engineering problems [19,20]. Engineering problems involving electrical machines, actuators, antennas, waveguides, and other components can be effectively solved with BEM. BEM simplifies the data needed to process and solve the problem by simply requiring discretization on a body’s surface. Because it builds surface mesh quickly and readily using discontinuous elements, BEM is a flexible numerical method. It also offers extremely accurate solutions derived from boundary elements. Fahmy [21,22,23] developed BEM models for solving dual-phase-lag bioheat problems. We refer the interested reader to BEM references in soft tissues studies [24,25,26] and engineering studies [27,28], and the references cited therein.
The primary goal of this paper is to propose a novel numerical BEM formulation for solving fractional bio-thermomechanical problems in FGA nonlinear viscoelastic soft tissues. To determine the temperature distribution, first, the bioheat equation is solved using the general boundary element method (GBEM) based on the local radial basis function collocation method (LRBFCM). The CQBEM is then used to solve the mechanical equation. The resulting linear systems were solved using a preconditioner called generalised modified shift-splitting (GMSS), which speeds up the convergence of the proposed technique. In comparison to the literature [20,21], this paper discusses the viscosity effects on the biothermal stress wave of FGA nonlinear viscoelastic soft tissues in the context of the classical Fourier, single-phase-lag, and dual-phase-lag bioheat models. Additionally, it investigates the viscosity effects on the biothermal stress wave of FGA nonlinear viscoelastic soft tissues for different values of fractional parameter, and functionally graded parameter. The present paper investigates the viscosity effects on the biothermal stress wave of isotropic, transversely isotropic, and anisotropic functionally graded nonlinear soft tissues. In this paper, the validity of the BEM is demonstrated through comparison with the generalized finite difference method (GFDM), and finite element method (FEM) under the special circumstances of the current problem. This paper also compares three efficient and powerful iterative methods used for solving the linear systems generated by the proposed technique. The generalized modified shift-splitting (GMSS) method outperforms the communication-avoiding Arnoldi (CA-Arnoldi) and Dehdezi and Karimi (D-K) methods. The numerical results demonstrate that the proposed model is correct, efficient, and precise.

2. Formulation of the Problem

Consider a region Ω bounded by a boundary Γ in the Cartesian coordinate system O x y . According to Biot’s model [29] and Darcy’s law [30], the fractional bio-thermomechanical governing equations in the context of FGA nonlinear viscoelastic soft tissues may be written as
T σ T + F = ρ u ¨ + ϕ ρ f u ¨ f - u ¨
ζ ˙ + T q = 0
σ = C a j l g   χ   t r   ϵ - A p I - B T
ϵ = 1 2 u T + u T T
ζ = A   t r   ϵ + ϕ 2 R P
q = - K p + ρ f u ¨ + ρ e + ϕ ρ f ϕ u ¨ f - u ¨
The fractional bioheat equation can be expressed without dual-phase-lag as
D τ a T x , τ = ξ K ˇ T ( x , τ ) + Q m x , τ , ξ = 1 c ρ s   a n d   ρ e = C ϕ ρ f
where we assumed that C = 0.66 is the shape factor as in Bonnet and Auriault [31].
Based on Bonnet [32], the governing equations may be expressed as
B ^ x ~ u ^ g x ~ = 0   f o r   x ~ Ω
where
B ^ x ~ = B x ~ e + s 2 ρ - β ρ f I α - β x ~ - B x ~ s α - β x ~ T - β s ρ f x ~ + s ϕ 2 R 0 , u ^ g x = G ^ D   f o r   x Γ 1
t ^ g x = T x e - α n x 0 s β n x T β s ρ f n x T x 0 u ^ x p ^ x θ x , β = ϕ 2 s K ρ f ϕ 2 + s K ( ρ e + ϕ ρ f ) , t ^ g x = G ^ N   f o r   x Γ 2

3. BEM Temperature Analysis

After discretization of the time interval 0 , F into F + 1 equal subintervals, we suppose that the temperature’s time derivative in τ f , τ f + 1 is given by
T ˙ x , τ = T f + 1 ( x ) - T f ( x ) τ + O τ , T f ( x ) = T ( x , τ f ) , τ f = f τ , f = 0 , 1 , 2 , . , F
where
D τ a T α x , τ = 1 Γ 1 - a 0 τ T α r , s s d s τ - s a , 0 < a < 1
Based on (11), we have
D τ a T α f + 1 + D τ a T α f j = 0 k W a , j T f + 1 - j x - T f - j x , f = 1 , 2 , . , F
in which
W a , 0 = τ - a Γ 2 - a , W a , j = W a , 0 j + 1 1 - a - j - 1 1 - a , j = 1 , 2 , . , F
The fractional bioheat equation without dual-phase-lag (7) can now be written using Equation (12) as
W a , 0 T f + 1 x - K ˇ x T , I I f + 1 x - K ˇ , I x T , I f + 1 x = W a , 0 T f x - K ˇ x T , I I f x - K ˇ , I x T , J f x - j = 1 f W a , j T f + 1 - j ( x ) - T f - j ( x ) + Q m f + 1 x , τ + Q m f x , τ , f = 0 , 1 , 2 , , F
The solution of (14) can be derived by Krahulec et al. [33] as
W a , 0 E - Λ Φ ˇ I I Φ ˇ - 1 - Λ I Φ ˇ I Φ ˇ - 1 T ˇ f + 1 = W a , 0 E + Λ Φ ˇ I I Φ ˇ - 1 + Λ I Φ ˇ I Φ ˇ - 1 T ˇ f - j = 1 f W a , j T f + 1 - j x - T f - j x + Q m f + 1 x , τ + Q m f x , τ , f = 0 , 1 , 2 , , F .
Based on [21], the dual-phase-lag bioheat equation without fractional derivative is as follows:
c ρ T τ + τ q 2 T τ 2 = K ˇ 2 T + K ˇ τ T τ 2 T + W b C b T ˘ b - T + Q m - W b C b τ q T τ
Now, at each transition τ f - 1 τ f we implement the boundary element approach to solve (see Appendix A)
2 U 1 x - B U 1 x + R T k - 1 f x = 0
According to [34,35], the following equation can be written
j = 1 N G i j W 1 x j = j = 1 N H i j U 1 x j + l = 1 L P i l R T k - 1 f x l
in which
G i j = 1 K ˇ Γ j T ξ i , x d Γ j , H i j = Γ j q ξ i , x d Γ j , i j - 0.5 , i = j , P i l = Ω j T ξ i , x d Ω j
Hence, the dual-phase-lag bioheat equation without fractional derivative has the following solution:
U 1 ξ i = j = 1 N H i j U 1 x j - j = 1 N G i j W 1 x j + j = 1 N P i l R T k - 1 f x l
The general solution of the non-linear fractional derivative dual-phase-lag bioheat equation T N L F D D P L has been demonstrated to be the sum of two solutions T N L F D + T N L D P L , where T N L F D denotes the solution of non-linear fractional derivative bioheat equation without dual-phase-lag (NLFD without DPL) which is represented by Equation (15) and T N L D P L denotes the solution of the non-linear dual-phase-lag bioheat equation without fractional derivative (NLDPL without FD) which is represented by Equation (29), i.e., T N L F D D P L = T N L F D + T N L D P L .

4. BEM Displacement Analysis

According to (8), the representation equation can be expressed as
u ^ g x ~ = V ^ t ^ g Ω x ~ - K ^ u ^ g Ω x ~   f o r   x ~ Ω
where the integral operators, fundamental solution and traction are [9]
V ^ t ^ g Ω x ~ = Γ . U ^ T y - x ~ t ^ g y d s y
K ^ u ^ g Ω x ~ = Γ . T ^ y U ^ T y - x ~ u ^ g y d s y
U ^ r = U ^ s r U ^ f r 0 P ^ s T ( r ) P ^ f ( r ) 0 , T ^ y = T y e s α n y 0 - β n y T β s ρ f n y T 0   w i t h   r y - x
where
U ^ s r = 1 4 π r ρ - β ρ f R 1 k 4 2 - k 2 2 k 1 2 - k 2 2 e - k 1 r - R 2 k 4 2 - k 1 2 k 1 2 - k 2 2 e - k 2 r + I k 3 2 - R 3 e - k 3 r
and
R j = 3 y r y T r - I r 2 + k j 3 y r y T r - I r + k j 2 y r y T r
Thus, the fundamental solution U ^ s r = U ^ s s r + U ^ r s r that can be divided into singular U ^ s s r and regular U ^ r s r parts may be expressed as [30]
U ^ s r = 1 4 π μ r λ + 2 μ λ + μ y r y T r + I λ + 3 μ + O ( r 0 )
Then, we have
U ^ s r = 1 μ I y - λ + μ λ + 2 μ y y T y x ^ r - 1 μ k 1 2 + k 2 2 y - k 1 2 k 2 2 I - k 1 2 + k 2 2 - k 4 2 - k 1 2 k 2 2 k 3 2 y y T x ^ r
U ^ f r = ρ f α - β y r 4 π r β λ + 2 μ k 1 2 - k 2 2 k 1 + 1 r e - k 1 r - k 2 + 1 r e - k 2 r = O r 0
P ^ s r = U ^ f r s = O r 0
P ^ f r = s ρ f 4 π r β k 1 2 - k 2 2 k 1 2 - k 4 2 e - k 1 r - k 2 2 - k 4 2 e - k 2 r = s ρ f 4 π r β + O r 0
in which
x ^ r = 1 4 π r e - k 1 r k 2 2 - k 1 2 k 3 2 - k 1 2 + e - k 2 r k 2 2 - k 1 2 k 2 2 - k 3 2 + e - k 3 r k 1 2 - k 3 2 k 2 2 - k 3 2 = - 1 k 1 2 - k 2 2 k 1 2 - k 3 2 k 3 2 - k 2 2 + O r 2
By using the limit x ~ Ω x Γ , we can write (21) as
lim x ~ Ω x Γ V ^ t ^ g Ω x ~ = V ^ x ^ g x Γ . U ^ T y - x t ^ g y d s y
Based on the limiting method [30], the following equation can be achieved
lim x ~ Ω x Γ K ^ u ^ g Ω x ~ = - I x + C x u ^ g x + K ^ u ^ g x
in which
C x = lim ε 0 y Ω : y - x = ε . T ^ y U ^ T y - x d s y
and
K ^ u ^ g x = lim ε 0 y - x ε . T ^ y U ^ T y - x u ^ g y d s y
By using Equations (32)–(35) we get
x u ^ g x = V ^ t ^ g x - K ^ u ^ g x
By using the inverse Laplace transformation, we have
C x u g x , t = V t g x , t - K u g x , t
Now, the fundamental poroelastodynamic solution can be written as
T ^ y U ^ T = T ^ y e s α n y - β n y T β s ρ 0 f n y T y U ^ s U ^ f P ^ s T P ^ f T = T ^ s T ^ f Q ^ s T Q ^ f T
where
Γ . y × a , n y d s y = - Γ , a , v d γ y = - ϕ , a , v d γ y = 0 , y Γ
and
Γ . n y × y , a d s y = 0
Now, we can write [36]
Γ . M y a d s y = 0
where M y = y y T T - y y T .
By applying (41) to a = v u we get [37]
Γ . M y v u d s y = - Γ , v M y u d s y
Γ . M y v T u d s y = Γ . v T M y u d s y
By using (27) and (38), we obtain
T ^ s T = T y e U ^ s i n g s + U ^ r e g s T + s α P ^ s n y T = T y e U ^ s i n g s T + O r 0
According to [36], we get
T ^ s T = λ + 2 μ n y y T U ^ s i n g s - μ n y × y × U ^ s i n g s + 2 μ M y U ^ s i n g s + o r 0
which can be expressed using (27) as
T ^ s T = M y y 2 X ^ + I n T y y 2 X ^ + 2 μ M y U ^ s i n g s T + o r 0
Substituting (22) into (46) to obtain
k ^ u ^ s x ~ = . M y y 2 X ^ u ^ + I n T y y 2 X ^ u ^ + 2 μ M y U ^ s i n g s T u ^ + O r 0 u ^ d s y
Based on [36], we get
K ^ u ^ s x ~ = . - y 2 X ^ M y u ^ + I n T y y 2 x ^ u ^ + 2 μ U ^ s s M y u ^ + O r 0 u ^ d s y
The second phrase on the right side of (48) has been changed to obtain
n T y y 2 x ^ r = n T y r 4 π r 2 + O r 0
in which
C s x = I x c x
where c x = ϕ x 4 π .
Based on [36], we can write
lim Ω x ~ x Γ K ^ u ^ Ω s x ~ = - I x - 1 + c x u ^ x + K ^ u ^ s x
By employing (43), Equation (48) can be expressed as
K ^ u ^ Ω s x ~ = Γ . - y 2 x ^ M y u ^ + I n T y y 2 x ^ u ^ + 2 μ U ^ s M y u ^ + O r 0 u ^ d s y
The convolution integral can be described as
f g t = 0 t f t - τ g τ d τ   f o r   t 0 , T
where
f g t n k = 0 n ω n - k t f ^ g t k , t n = n t ,   t > 0
Based on [38,39], we obtain
ω n t f ^ 1 2 π i z = R f ^ γ z t z - n + 1 d z
Now, by employing z = R e - i φ , the integral (55) can be approximated as follows
ω n t f ^ R - 1 L + 1 = o L f ^ s ζ n   w i t h   ζ = e 2 π i L + 1   a n d   s = γ R ζ - t
By plugging Equation (56) into Equation (54), we get
f g t n k = 0 N R - n - k N + 1 = 0 N f ^ s ζ n - k g t k R - n N + 1 = 0 N f ^ s g ^ s ζ n
with
g ^ s = k = 0 N R k g t k ζ - k .
Based on the procedure [38], and using (37), we have
C x u g x , t = v     t g x , t - k u g x , t
which can be written as
C x u ^ g x , s = v ^ t ^ g x , s - k ^ u ^ g x , s , = 0 . . N
The following discretization can be used:
Γ Γ h = e = 1 N e τ e
Also, we consider the following definitions
S h k Γ N , h s p a n φ i α k i = 1 𝕚 , α 1
S h k Γ D , h s p a n ψ j β k j = 1 𝕛 , β 0
The Neumann and Dirichlet datum can be calculated, respectively, as
u ^ g k x u ^ h g k x = i = 1 I u ^ h , i g k φ i α k x S h k Γ N , h , i = 1 , , I , k = 1 , 2 , 3 , 4
t ^ g k x t ^ h g k x = j = 1 J t ^ h , j g k ψ j β k x S h k Γ D , h , j = 1 , , J , k = 1 , 2 , 3 , 4
Thus, we obtain
V ^ D D V ^ N D - K ^ D N - C + K ^ N N t ^ D , h g u ^ N , h g = - V ^ D N - V ^ N N C + K ^ D D K ^ N D g ^ N , h g g ^ D , h g = 0 . . . N
where
S ^ N N V ^ N D V ^ D D - 1 K ^ D N - C + K ^ N N

5. Numerical Results and Discussion

The current study’s proposed fractional BEM method can be applied to a wide variety of bio-thermomechanical problems in FGA nonlinear viscoelastic soft tissues. The communication-avoiding Arnoldi (CA-Arnoldi) [40], Dehdezi and Karimi method (D-K) [41] and generalized modified shift-splitting (GMSS) [42] are the three efficient and powerful iterative methods used for solving the linear systems generated by the proposed technique. As in [43,44], we considered isotropic, transversely isotropic, and anisotropic soft tissues in the current paper. We also used the BEM discretization with 84 boundary elements and 404 internal nodes, as depicted in Figure 1. The calculations are carried out by Matlab in OS Windows 10 (64 bit) with AMD Core R7-4800H and 16G memory.
Figure 2 shows the biothermal stress σ 11 wave fluctuations along the x 1 axis in the absence and presence of a viscosity effect for different values of the fractional parameter a . It is noticed from present figure that the biothermal stress σ 11 begins from a negative value and satisfies the boundary condition at x 1 = 0 . In the context of different values of the fractional parameter and in the absence of viscosity effect, the biothermal stress σ 11 increases to a maximum value in the range 0 x 1 1.4 , then decreases in the range 1.4 x 1 7 .
Figure 3 displays the biothermal stress σ 11 wave fluctuations along the x axis in the absence and presence of a viscosity effect for anisotropic (AI), transversely isotropic (TI), and isotropic (I) soft tissues. It is shown from this figure that the biothermal stress σ 11 in the absence of a viscosity effect decreases to a minimum value in the distance 0 x 1 0.5 , then increases in the distance 0.5 x 1 2.25 , and then decreases in the distance 2.25 x 1 7 .
Figure 4 shows the biothermal stress σ 11 wave fluctuations along the x 1 axis in the absence and presence of a viscosity effect in the context of the classical Fourier F o u r i e r , τ q = τ T = 0 , single-phase-lag S P L , τ q = 0   a n d   τ T = 25 , and dual-phase-lag D P L , τ q = τ T = 25 models. It is seen from present figure that the biothermal stress σ 11 begins from negative values in the absence and presence of a viscosity effect. Additionally, it is noticed in the absence of a viscosity effect that the biothermal stress σ 11 increases to a maximum value in the range of 0 x 1 0.25 , then decreases in the range of 0.25 x 1 1.5 , and then increases in the range of 1.5 x 1 7 for the current model.
Figure 5 shows that the biothermal stress σ 11 wave fluctuations along the x 1 axis in the absence and presence of a viscosity effect, in the context of various functionally graded parameter values, begin from positive values. In the context of the higher values of the graded parameter m = 0.6 and m = 1.0 , the biothermal stress σ 11 decreases to a minimum value in the distance 0 x 1 1.1 , then increases in the distance 1.1 x 1 2.4 , and then decreases in the distance 2.4 x 1 7 . However, in the context of the lower value m = 0.2 , the biothermal stress σ 11 decreases to a minimum value in the distance 0 x 1 1.1 , then increases, and moves in wave fluctuations in the absence of a viscosity effect. In the context of the higher graded parameter values m = 0.6 and m = 1.0 , the biothermal stress σ 11 begins by decreasing to a minimum value, then increasing, and then decreasing in the presence of a viscosity effect. However, in the context of the lower graded parameter value m = 0.2 , the biothermal stress σ 11 decreases to a minimum value in the distance 0 x 1 1.25 , then increases, and moves in wave fluctuations in the presence of the viscosity effect.
To understand the physical meaning of the biothermal stress σ 11 negative and positive variations, we should know that the hypothalamus helps humans to self-regulate their body temperature by comparing their current internal temperature to their body’s “normal” temperature. When the body’s temperature falls below what is considered normal, it sends the signal to generate heat. Unless you have hypothermia, your core temperature should be stable. As cold air steals heat from your body, the temperature of your skin—how your fingers, toes, legs, arms, and forehead feel—may begin to drop. It is possible for your body to become extremely cold, affecting your core temperature. This is hazardous to your health and is a medical emergency. Continue reading to learn how to increase your body temperature. There were no published findings to back up the proposed technique’s findings. Some of the literature, however, can be viewed as a part of the considered study. Therefore, we considered a special case of our study and implemented generalized finite difference method (GFDM) [45], and finite element method (FEM) [46] to this special case, as shown in Figure 6.
Figure 6 depicts the biothermal stress σ 11 fluctuations along x 1 axis ( x 2 = 0.5 ) for BEM, GFDM, and FEM. This figure shows that the BEM results are in excellent agreement with the GFDM and FEM results.
The CPU time and number of iterations for the CA-Arnoldi, D-K, and GMSS iterative methods at each discretization level are shown in Table 1, with the equation numbers enclosed in parentheses. The GMSS method outperforms the CA-Arnoldi and D-K methods, as shown in the table.
Table 2 contrasts the computing resources for modeling nonlinear fractional bio-thermomechanical problems in FGA viscoelastic soft tissues using the generalized finite difference method (GFDM) of Turchan [45], the finite element method (FEM) of Marin et al. [46], and the boundary element method (BEM). The effectiveness of our suggested BEM method is displayed in the Table 2.

6. Conclusions

It was proposed to use a numerical BEM technique to describe fractional bio-thermomechanical interactions in FGA nonlinear viscoelastic soft tissues. The numerical procedure is made easier by boundary discretization. BEM makes mesh formation for the considered problem easier. Because BEM is a semi-analytical method, a high accuracy is obtained. Both open and moving boundary systems can benefit from this. BEM can be used in conjunction with analytical techniques and numerical techniques such as FEM. To calculate the temperature, the LRBFCM and GBEM were used to solve the bio-heat governing equation. To obtain displacements and stresses at each time step, the CQBEM was used to solve the poroelastic governing equation. The GMSS was developed to solve equations generated by boundary element discretization. The numerical outcomes are graphed to demonstrate the effects of the functionally graded parameter, fractional parameter, and anisotropic material property on biothermal stress. The numerical results also show how the Fourier, single-phase-lag, and dual-phase-lag bio-heat models differ. The primary benefit of the current technique is its generality and simplicity. The numerical results show that the proposed method outperforms other domain discretization techniques. The research presented in this paper contributes to the development of mathematical models for various biothermal engineering applications. We propose expanding the boundary element technique proposed in this paper for applications in three-dimensional bioheat transfer models involving multilayer problems, complex geometries, and the presence of convective terms in future work. The local radial basis function collocation method based on bioheat transfer modeling and biomechanics should be investigated because it may present complex modeling phenomena such as soft tissue movement and thermal expansion effect.

Author Contributions

Conceptualization, M.A.F. and M.M.A.; methodology, M.A.F.; software, M.A.F.; validation, M.A.F. and M.M.A.; formal analysis, M.A.F.; investigation, M.A.F.; resources, M.A.F.; data curation, M.A.F.; writing—original draft preparation, M.A.F.; writing—review and editing, M.A.F.; visualization M.A.F.; supervision, M.A.F.; project administration, M.A.F.; funding acquisition, M.A.F. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Deanship of Scientific Research at Umm Al-Qura University, grant number [22UQU4340548DSR11], and the APC was funded by the Deanship of Scientific Research at Umm Al-Qura University.

Data Availability Statement

All data generated or analysed during this study are included in this published article.

Acknowledgments

The authors would like to thank the Deanship of Scientific Research at Umm Al-Qura University for supporting this work through Grant Code: (22UQU4340548DSR11).

Conflicts of Interest

The authors declare no conflict of interest.

Nomenclature

* Convolution with respect to time
γ Internal angle at point source i
Γ Boundary
Γ D Dirichlet boundary condition
Γ N Neumann boundary condition
δ i j Kronecker delta ( i , j = 1 , 2 )
ε i j Strain tensor
ε i j k Alternate tensor
ϵ Linear strain tensor
μ 0 Shear moduli
χ Viscoelastic constant
ζ Fluid volume variation
ρ = ρ s 1 - ϕ + ϕ ρ f Bulk density
ρ s Soft tissue density
ρ f Fluid density
σ Total stress tensor
τ Time
τ ¯ Thermal relaxation time
τ q Phase lag of heat flux
τ T Phase lag of temperature gradient
φ i α k Continuous shape functions
ϕ = V f V Porosity
ψ j β k Discontinuous shape functions
Ω Region
A = ϕ 1 + Q / R Biot’s coefficient
a Fractional parameter
B Stress-temperature coefficients
B x ~ e Linear elastostatics operator
C Shape factor
c Specific heat of soft tissue
C b Blood specific heat
C i j k l Constant elastic moduli
E Young’s moduli
F Bulk body forces
g ^ D Dirichlet datum
g ^ N Neumann datum
K Iterations number
K ˇ Soft tissue thermal conductivity
m Iterative parameter
n Outward unit normal vector
p Fluid pressure in the vasculature
q Fluid specific flux
Q , R Fluid-soft tissue coupling parameters
r Euclidean distance between source point and field point
Q m Metabolic heat generation
T r , τ Soft tissue temperature
T ˘ b r , τ Blood temperature
T x e Traction derivative
t ^ g Generalized tractions
T r Trace of a matrix
u Soft tissue displacement
u f Fluid displacement
V Bulk volume
V f Fluid volume
V s Soft tissue volume
V Poisson’s ratio
W a , j Fractional derivative coefficient
W b Blood perfusion rate
x , y Space coordinates
x Source point
y Considered point
R , Z , G Matrices

Appendix A

Based on the finite difference scheme [47], Equation (16) yields the following formula
2 T f x - B T f x + C 2 T f - 1 x + D T f - 1 + E T f - 2 x + F = 0
where
B = c ρ + W b C b τ τ + τ q K ˇ τ τ + τ T , C = τ T τ + τ T , D = c ρ τ + 2 τ q + W b C b τ q τ K ˇ τ τ + τ T ,
E = - c ρ τ q K ˇ τ τ + τ T , F = τ W b C b T ˘ b + Q m K ˇ τ + τ T
and
T f x = T b f x   f o r   x Γ 1
Z f x = w b f x = - K ˇ T f x n   f o r   x Γ 2
where
w b f x = τ τ + τ T q b f x + τ q q b x , τ τ t = t f - K ˇ τ T τ + τ T T f - 1 x n
Based on [48,49], the following equation can be derived
1 - P L Ψ x ; P - U x = - P A Ψ x ; P
and
Ψ x ; P = P T b f x + 1 - P U x , x Γ 1
- K ˇ Ψ x ; P n = P w b f x + 1 - P - K ˇ U x n , x Γ 2
Then the linear operator is
L U = 2 u - B u = 2 u x 1 2 + 2 u x 2 2 + 2 u x 3 2 - B u
According to [28], the following differential formula can be derived
- L W x ; P - U x + 1 - P L W x ; P P - U x P = - Ζ W x ; P - P A W x ; P P
and
W x ; P P = T b f x + T b f x P - U x + 1 - P U x P   f o r   x Γ 1
- K ˇ n W x ; P P = w b f x + P w b f x P + K ˇ U n - 1 - P K ˇ n U x P   f o r   x Γ 2

References

  1. Cho, Y.I. Advances in Heat Transfer: Bioengineering Heat Transfer; Academic Press Inc.: San Diego, CA, USA, 1992. [Google Scholar]
  2. Pennes, H.H. Analysis of tissue and arterial blood temperatures in the resting human forearm. J. Appl. Physiol. 1948, 1, 93–122. [Google Scholar] [CrossRef] [PubMed]
  3. Cattaneo, C. Sur une forme de i’equation de la chaleur elinant le paradox d’une propagation instantance. Comptes Rendus L’académie Sci. 1958, 247, 431–433. [Google Scholar]
  4. Vernotee, M.P. Les paradoxes de la theorie continue de i equation de la chleur. Comptes Rendus L’académie Sci. 1958, 246, 3154–3155. [Google Scholar]
  5. Tzou, D.Y. A unified field approach for heat conduction from micro to macroscale. J. Heat Transf. 1995, 117, 8–16. [Google Scholar] [CrossRef]
  6. Tzou, D.Y. Macro- to Microscale Heat Transfer: The Lagging Behavior; Taylor & Francis: Washington, DC, USA, 1996. [Google Scholar]
  7. Pei, R.Z.; Jing, L.; Cheng, W.C.; Xue, J.P. Boundary element method (BEM) for solving normal or inverse bio-heat transfer problem of biological bodies with complex shape. J. Therm. Sci. 1995, 4, 117–124. [Google Scholar]
  8. Ooi, E.H.; Ang, W.T.; Ng, E.Y.K. A boundary element model for investigating the effects of eye tumor on the temperature dis-tribution inside the human eye. Comput. Biol. Med. 2009, 39, 667–677. [Google Scholar] [CrossRef]
  9. Zhou, J.; Chen, J.K.; Zhang, Y. Simulation of Laser-Induced Thermotherapy Using a Dual-Reciprocity Boundary Element Model with Dynamic Tissue Properties. IEEE Trans. Biomed. Eng. 2010, 57, 238–245. [Google Scholar] [CrossRef]
  10. Ng, E.; Tan, H.; Ooi, E. Boundary element method with bioheat equation for skin burn injury. Burns 2009, 35, 987–997. [Google Scholar] [CrossRef]
  11. Majchrzak, E.; Turchan, L. The general boundary element method for 3D dual-phase lag model of bioheat transfer. Eng. Anal. Bound. Elements 2015, 50, 76–82. [Google Scholar] [CrossRef]
  12. Bottauscio, O.; Chiampi, M.; Zilberti, L. Boundary Element Solution of Electromagnetic and Bioheat Equations for the Sim-ulation of SAR and Temperature Increase in Biological Tissues. IEEE Trans. Magn. 2012, 48, 691–694. [Google Scholar] [CrossRef]
  13. Deng, Z.-S.; Liu, J. Modeling of multidimensional freezing problem during cryosurgery by the dual reciprocity boundary element method. Eng. Anal. Bound. Elements 2004, 28, 97–108. [Google Scholar] [CrossRef]
  14. Partridge, P.W.; Wrobel, L.C. A coupled dual reciprocity BEM/genetic algorithm for identification of blood perfusion param-eters. Int. J. Numer. Methods Heat Fluid Flow 2009, 19, 25–38. [Google Scholar] [CrossRef] [Green Version]
  15. Abd-Elaziz, E.; Othman, M.I.A. Effect of rotation on a micropolar magnetothermo-elastic medium with dual-phase-lag model under gravitational field. Microsyst. Technol. 2017, 23, 4979–4987. [Google Scholar]
  16. Othman, M.I.; Eraki, E.E. Effect of gravity on generalized thermoelastic diffusion due to laser pulse using dual-phase-lag model. Multidiscip. Model. Mater. Struct. 2018, 14, 457–481. [Google Scholar] [CrossRef]
  17. Othman, M.I.A.; Abd-Elaziz, E.M. Dual-phase-lag model on micropolar thermoelastic rotating medium under the effect of thermal load due to laser pulse. Indian J. Phys. 2020, 94, 999–1008. [Google Scholar] [CrossRef]
  18. Othman, M.I.A.; Zidan, M.E.M.; Mohamed, I.E.A. Dual-phase-lag model on thermo-microstretch elastic solid under the effect of initial stress and temperature-dependent. Steel Compos. Struct. 2021, 38, 355–363. [Google Scholar]
  19. Hosseini, M.; Paparisabet, M.A. The effects of blood flow on blood vessel buckling embedded in surrounding soft tissues. Int. J. Appl. Mech. 2016, 8, 1650065. [Google Scholar] [CrossRef]
  20. Fahmy, M.A. A Computational model for nonlinear biomechanics problems of FGA biological soft tissues. Appl. Sci. 2022, 12, 7174. [Google Scholar] [CrossRef]
  21. Fahmy, M.A.; Almehmadi, M.M.; Al Subhi, F.M.; Sohail, A. Fractional boundary element solution of three-temperature thermoelectric problems. Sci. Rep. 2022, 12, 6760. [Google Scholar] [CrossRef]
  22. Fahmy, M.A. 3D Boundary Element Model for Ultrasonic Wave Propagation Fractional Order Boundary Value Problems of Functionally Graded Anisotropic Fiber-Reinforced Plates. Fractal Fract. 2022, 6, 247. [Google Scholar] [CrossRef]
  23. Fahmy, M.A. A new boundary element algorithm for a general solution of nonlinear space-time fractional dual-phase-lag bio-heat transfer problems during electromagnetic radiation. Case Stud. Therm. Eng. 2021, 25, 100918. [Google Scholar] [CrossRef]
  24. Fahmy, M.A. A new LRBFCM-GBEM modeling algorithm for general solution of time fractional-order dual phase lag bioheat transfer problems in functionally graded tissues. Numer. Heat Transfer Part A Appl. 2019, 75, 616–626. [Google Scholar] [CrossRef]
  25. Fahmy, M.A. Boundary Element Algorithm for Modeling and Simulation of Dual Phase Lag Bioheat Transfer and Biome-chanics of Anisotropic Soft Tissues. Int. J. Appl. Mech. 2018, 10, 1850108. [Google Scholar] [CrossRef]
  26. Chan, C.L. Boundary Element Method Analysis for the Bioheat Transfer Equation. J. Biomech. Eng. 1992, 114, 358–365. [Google Scholar] [CrossRef]
  27. Wrobel, L.; Kassab, A. Boundary Element Method, Volume 1: Applications in Thermo-Fluids and Acoustics. Appl. Mech. Rev. 2002, 56, B17. [Google Scholar] [CrossRef]
  28. Lua, W.Q.; Liub, J.; Zenga, Y. Simulation of the thermal wave propagation in biological tissues by the dual reciprocity boundary element method. Eng. Anal. Bound. Elem. 1998, 22, 167–174. [Google Scholar] [CrossRef]
  29. Biot, M.A. Theory of propagation of elastic waves in a fluid- saturated porous solid. II. Higher frequency range. J. Acoust. Soc. Am. 1956, 28, 179–191. [Google Scholar] [CrossRef]
  30. Schanz, M. Wave Propagation in Viscoelastic and Poroelastic Continua. In Lecture Notes in Applied Mechanics; Springer: New York, NY, USA, 2001; Volume 2. [Google Scholar]
  31. Bonnet, G.; Auriault, J.L. Dynamics of saturated and deform- able porous media: Homogenization theory and determination of the solid-liquid coupling coefficients. In Physics of Finely Divided Matter; Springer: Berlin/Heidelberg, Germany, 1985; pp. 306–316. [Google Scholar]
  32. Bonnet, G. Basic singular solutions for a poroelastic medium in the dynamic range. J. Acoust. Soc. Am. 1987, 82, 1758–1762. [Google Scholar] [CrossRef]
  33. Krahulec, S.; Sladek, J.; Sladek, V.; Hon, Y.-C. Meshless analyses for time-fractional heat diffusion in functionally graded materials. Eng. Anal. Bound. Elements 2016, 62, 57–64. [Google Scholar] [CrossRef]
  34. Steinbach, O. Numerical Approximation Methods for Elliptic Boundary Value Problems; Springer: New York, NY, USA, 2008. [Google Scholar]
  35. Wang, P.; Becker, A.A.; Jones, I.A.; Glover, A.T.; Benford, S.D.; Vloeberghs, M. Real-time surgical simulation for deformable soft-tissue objects with a tumour using Boundary Element techniques. J. Phy. Conf. Ser. 2009, 181, 012016. [Google Scholar] [CrossRef]
  36. Fahmy, M.A. Design Optimization for A Simulation of Rotating Anisotropic Viscoelastic Porous Structures Using Time-Domain OQBEM. Math. Comput. Simul. 2019, 66, 193–205. [Google Scholar] [CrossRef]
  37. Kielhorn, L. A time-domain symmetric Galerkin BEM for viscoelastodynamics. In Computation in Engineering and Science; Brenn, G., Holzapfel, G.A., Schanz, M., Steinbach, O., Eds.; Graz University of Technology, Institute of Applied Mechanics: Graz, Austria, 2009. [Google Scholar]
  38. Lubich, C. Convolution quadrature and discretized operational calculus. I. Numer. Math. 1988, 52, 129–145. [Google Scholar] [CrossRef]
  39. Lubich, C. Convolution quadrature and discretized operational calculus. II. Numer. Math. 1988, 52, 413–425. [Google Scholar] [CrossRef]
  40. Hoemmen, M. Communication-Avoiding Krylov Subspace Methods. Ph.D. Thesis, University of California, Berkeley, CA, USA, 2010. [Google Scholar]
  41. Dehdezi, E.K.; Karimi, S. A rapid and powerful iterative method for computing inverses of sparse tensors with applications. Appl. Math. Comput. 2022, 415, 126720. [Google Scholar]
  42. Huang, Z.-G.; Wang, L.-G.; Xu, Z.; Cui, J.-J. The generalized modified shift-splitting preconditioners for nonsymmetric saddle point problems. Appl. Math. Comput. 2017, 299, 95–118. [Google Scholar] [CrossRef]
  43. Gilchrist, M.; Murphy, J.; Parnell, W.; Pierrat, B. Modelling the slight compressibility of anisotropic soft tissue. Int. J. Solids Struct. 2014, 51, 3857–3865. [Google Scholar] [CrossRef]
  44. Morrow, D.A.; Donahue, T.L.H.; Odegard, G.M.; Kaufman, K.R. Transversely isotropic tensile material properties of skeletal muscle tissue. J. Mech. Behav. Biomed. Mater. 2010, 3, 124–129. [Google Scholar] [CrossRef] [Green Version]
  45. Turchan, Ł. Solving the dual-phase lag bioheat transfer equation by the generalized finite difference method. Arch. Mech. 2017, 69, 389–407. [Google Scholar]
  46. Marin, M.; Hobiny, A.; Abbas, I. Finite Element Analysis of Nonlinear Bioheat Model in Skin Tissue Due to External Thermal Sources. Mathematics 2021, 9, 1459. [Google Scholar] [CrossRef]
  47. Messner, M.; Schanz, M. A regularized collocation boundary element method for linear poroelasticity. Comput. Mech. 2011, 47, 669–680. [Google Scholar] [CrossRef]
  48. Liao, S. General boundary element method for non-linear heat transfer problems governed by hyperbolic heat conduction equation. Comput. Mech. 1997, 20, 397–406. [Google Scholar] [CrossRef]
  49. Liao, S.; Chwang, A. General boundary element method for unsteady nonlinear heat transfer problems. Numer. Heat Transf. Part B 1999, 35, 225–242. [Google Scholar] [CrossRef]
Figure 1. BEM modeling of the present problem.
Figure 1. BEM modeling of the present problem.
Fractalfract 07 00066 g001
Figure 2. Biothermal stress σ 11 wave fluctuations along x 1 axis for different values of the fractional parameter.
Figure 2. Biothermal stress σ 11 wave fluctuations along x 1 axis for different values of the fractional parameter.
Fractalfract 07 00066 g002
Figure 3. Biothermal stress σ 11 wave fluctuations along x 1 axis for different soft tissues.
Figure 3. Biothermal stress σ 11 wave fluctuations along x 1 axis for different soft tissues.
Fractalfract 07 00066 g003
Figure 4. Biothermal stress σ 11 wave fluctuations along x 1 axis for various bioheat models.
Figure 4. Biothermal stress σ 11 wave fluctuations along x 1 axis for various bioheat models.
Fractalfract 07 00066 g004
Figure 5. Biothermal stress σ 11 wave fluctuations along x 1 axis for different functionally graded parameter values.
Figure 5. Biothermal stress σ 11 wave fluctuations along x 1 axis for different functionally graded parameter values.
Fractalfract 07 00066 g005
Figure 6. Biothermal stress σ 11 fluctuations along the x 1 axis for various techniques.
Figure 6. Biothermal stress σ 11 fluctuations along the x 1 axis for various techniques.
Fractalfract 07 00066 g006
Table 1. CPU times and iterations number for CA-Arnoldi, D-K and GMSS.
Table 1. CPU times and iterations number for CA-Arnoldi, D-K and GMSS.
Discretization
Level
Preconditioning
Level
CA-ArnoldiD-KGMSS
CPU TimeIterations NumberCPU TimeIterations NumberCPU TimeIterations Number
1 (32)00.0770.0870.057
2 (56)0
1
0.20
0.16
9
7
0.24
0.20
9
7
0.16
0.12
8
6
3 (104)0
1
2
0.56
0.50
0.44
11
9
6
0.62
0.56
0.48
12
10
8
0.42
0.34
0.26
10
6
4
4 (200)0
1
2
3
2.44
1.90
1.67
1.42
14
12
8
7
2.54
2.12
1.84
1.54
18
16
11
8
1.88
1.58
1.42
1.38
12
8
6
4
5 (392)0
1
2
3
4
10.03
9.11
8.21
7.26
6.62
16
12
10
8
6
12.04
10.18
9.31
8.50
7.02
20
18
16
12
8
7.90
6.89
6.11
5.84
5.20
14
10
8
6
4
6 (776)0
1
2
3
4
5
42.2
37.1
34.7
29.8
27.6
25.7
20
18
16
12
10
8
48.5
43.7
41.8
37.9
33.6
29.8
22
20
18
14
12
10
36.1
34.1
30.3
25.9
21.9
19.8
16
12
10
8
6
3
Table 2. Calculating the computational power required to model fractional bio-thermomechanical problems in FGA nonlinear viscoelastic soft tissues.
Table 2. Calculating the computational power required to model fractional bio-thermomechanical problems in FGA nonlinear viscoelastic soft tissues.
GFDMFEMBEM
Number of Nodes56,00054,00050
Number of Elements11,50010,50015
CPU Time (min)1801902
Memory (MB)1301401
Disc Space (MB)2602800
Accuracy of Results (%)2.22.01.0
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Fahmy, M.A.; Almehmadi, M.M. Fractional Dual-Phase-Lag Model for Nonlinear Viscoelastic Soft Tissues. Fractal Fract. 2023, 7, 66. https://doi.org/10.3390/fractalfract7010066

AMA Style

Fahmy MA, Almehmadi MM. Fractional Dual-Phase-Lag Model for Nonlinear Viscoelastic Soft Tissues. Fractal and Fractional. 2023; 7(1):66. https://doi.org/10.3390/fractalfract7010066

Chicago/Turabian Style

Fahmy, Mohamed Abdelsabour, and Mohammed M. Almehmadi. 2023. "Fractional Dual-Phase-Lag Model for Nonlinear Viscoelastic Soft Tissues" Fractal and Fractional 7, no. 1: 66. https://doi.org/10.3390/fractalfract7010066

APA Style

Fahmy, M. A., & Almehmadi, M. M. (2023). Fractional Dual-Phase-Lag Model for Nonlinear Viscoelastic Soft Tissues. Fractal and Fractional, 7(1), 66. https://doi.org/10.3390/fractalfract7010066

Article Metrics

Back to TopTop