Fractional thermoelasticity problem for an infinite solid with a cylindrical hole under harmonic heat flux boundary condition

The time-fractional heat conduction equation with the Caputo derivative results from the law of conservation of energy and time-nonlocal generalization of the Fourier law with the “long-tail” power kernel. In this paper, we consider an infinite solid with a cylindrical cavity under harmonic heat flux boundary condition. The Laplace transform with respect to time and the Weber transform with respect to the spatial coordinate are used. The solutions are obtained in terms of integrals with integrands being the Mittag-Leffler functions. The numerical results are illustrated graphically.


Introduction
The theory of heat conduction is an integral part of thermoelasticity. The standard parabolic heat conduction equation ∂ T ∂t = aΔT (1) results from the law of conservation of energy, and the phenomenological Fourier law which states the proportionality between the heat flux vector q and the temperature gradient, where C is the heat capacity, k is the thermal conductivity of a solid, and a = k/C denotes the thermal diffusivity coefficient.
In media with complex internal structure, the standard Fourier law (3) and the classical heat conduction equation (1) are no longer accurate enough, and nonclassical theories, in which these equations are replaced by more general equations, are formulated.
For example, the time-nonlocal dependence between the heat flux vector q and the temperature gradient with the "long-tail" power kernel [1], where Γ (α) is the gamma function, can be interpreted in terms of fractional integrals and derivatives, In combination with the law of conservation of energy (3), the constitutive equations (6) and (7) result in the time-fractional heat conduction equation with Caputo derivative, The details of obtaining Eq. (8) from the constitutive equations (6) and (7) can be found in [2].
Recall that the Riemann-Liouville fractional integral is introduced as a natural generalization of the repeated integral written in a convolution type form [3], The Riemann-Liouville derivative of the fractional order α is interpreted as left-inverse to the fractional integral I α f (t) [3]: The Caputo fractional derivative is defined as [4] D α Equation (8) is known as the diffusion-wave equation. For 0 < α ≤ 1, it interpolates the Helmholtz equation and the classical diffusion (heat conduction) equation; for 1 < α ≤ 2, it interpolates the parabolic heat conduction equation and the wave equation for temperature. The case α = 2 is known as ballistic heat conduction and corresponds to the time-nonlocal generalization of the Fourier law with "total memory" [5,6].
Each generalization of the heat conduction equation leads to the formulation of the corresponding extended theory of thermal stresses. Thermoelasticity without energy dissipation [6] starts from the wave equation for temperature. The generalized thermoelasticity of Lord and Shulman [7] is based on the Cattaneo telegraph equation for temperature. The first paper on fractional thermoelasticity associated with the time-fractional diffusion-wave equation (8) was published in 2005 (see [8]).
Later on, this subject has attracted much attention from researchers. Thermoelasticity based on the fractional telegraph equation was proposed in [9][10][11][12]. The theory of generalized thermoelasticity with memorydependent derivatives was considered in [13]. The book [14] sums up investigations in the field of fractional thermoelasticity.
It should be emphasized that the term "diffusion-wave" is also used in another context. Ångström [15,16] was the first to consider the classical diffusion (heat conduction) equation (1) under harmonic (wave) impact. Sometimes, the term "oscillatory diffusion" is also used. For detailed discussion, see [17][18][19] and references therein. Nowacki [20,21] considered the standard heat conduction equation with a source term varying harmonically as a function of time and studied the associated thermal stresses. In the works devoted to "oscillatory diffusion," the quasi-steady-state oscillations were considered under the assumption that the solution can be represented as a product of a function of the spatial coordinate and the time-harmonic term without taking into account the initial conditions. The Caputo derivative of the exponential function has a more complicated form than the corresponding derivative of integer order [22]. For this reason, for equations with fractional derivative, the assumption that the solution is represented as a product of a function of the spatial coordinate and the time-harmonic term cannot be used.
A plane with a cylindrical cavity in the framework of thermoelasticity based on the time-fractional heat conduction equation (8) under different boundary conditions was studied in [23,24]. A similar problem in the framework of thermoelasticity based on the time-fractional telegraph equation was investigated in [25], where inversion of the Laplace transform was carried out numerically. The control problem of axisymmetric thermal stresses associated with the time-fractional heat conduction equation (8) in an infinite cylindrical domain was studied in [26]. The heat source was used as a control function which guaranteed the distribution of stresses at a prescribed level. Fractional heat conduction in a space with a source varying harmonically in time and the corresponding thermal stresses were investigated in [22]. In the present paper, we study the axisymmetric time-fractional heat conduction equation in a solid with a cylindrical hole under harmonic heat flux boundary condition as well as associated thermal stresses in this domain.

Statement of the problem
In the framework of fractional thermoelasticity proposed in [8], the stressed-strained state of a solid is governed by the equilibrium equation in terms of displacements, the Duhamel-Neumann equation the geometrical relations and the time-fractional heat conduction equation Here, u is the displacement vector, σ is the stress tensor, e denotes the linear strain tensor, λ and μ are Lamé constants, K T = λ+2μ/3, β T is the thermal coefficient of volumetric expansion, and I denotes the unit tensor.
We have restricted ourselves to the quasi-static statement of the thermoelasticity problem when the relaxation time of mechanical oscillations is considerably less than the relaxation time of the heat conduction process. Neglecting in the heat conduction equation the term which describes coupling of thermal and mechanical effects which is valid for problems when thermoelastic energy dissipation is not of special interest.
In the considered case, Eq. (15) is written as We assume zero initial conditions, and harmonic boundary conditions for the heat flux vector, The zero condition for temperature at infinity is also adopted, In the case of axial symmetry, the equation of equilibrium (12) is rewritten as where ν is Poisson's ratio. The cavity surface is traction free, r = R : σ rr = 0.
We also have the zero condition for displacement at infinity:

Equation (22) has the solution [27]
The integration constants C 1 and C 2 are found from the boundary conditions (23) and (24) and in the considered case are equal to zero: The components of the stress tensor are calculated as integrals of the temperature field [27], In what follows, we will use the integral transform technique. Recall the Laplace transform rules for fractional integrals and derivatives [1], where the asterisk denotes the transform, and s is the Laplace transform variable. The Weber integral transform with respect to the spatial coordinate r is defined as [1,28,29] and has the inverse The specific expression of the kernel K 0 (r, R, ξ) depends on the boundary conditions at r = R. For the Neumann boundary condition with the given boundary value of the normal derivative, where J n (r ) and Y n (r ) are the Bessel functions of the first and second kind, respectively. For the Laplace operator in polar coordinates in the case of a Neumann boundary condition at the surface of a cylindrical hole r = R, we have Applying to the heat conduction equation (16) the Laplace transform with respect to time t and the Weber transform with respect to the radial coordinate r , we get .
The inverse integral transforms give the solution where E α (z) is the Mittag-Leffler function in one parameter α [30], and the following formula for the inverse Laplace transform: has been used. Equations (27) and (28) lead to expressions for components of the stress tensor The results of numerical calculations are presented in Figs. 1, 2, 3, and 4. In the calculations, we have used the following nondimensional quantities: All the calculations have been carried out for the real part of the solution. Figures 1 and 2 show the dependence of temperature and stresses on the distance. Figure 3 presents the dependence of temperature at the surface of a hole on the frequency ω. It follows from this Figure that a decrease in the order of the fractional derivative α results in an increase in attenuation. The time dependence of temperature at the surface of the hole is depicted in Fig. 4.

Conclusions
We have considered the axisymmetric time-fractional heat conduction equation with the Caputo derivative of the order 0 < α ≤ 2 in an infinite solid with a cylindrical cavity under harmonic heat flux boundary condition. It should be emphasized that for the fractional heat conduction equation, two types of Neumann boundary condition (the boundary condition of the second kind) can be considered: the mathematical condition with the prescribed boundary value of the normal derivative and the physical condition with the prescribed boundary value of the heat flux according to time-nonlocal generalization of the Fourier law (6) and (7). To calculate the Mittag-Leffler function E α (z) appearing in the solutions, we have used the algorithms suggested in [31]. (The interested reader is also referred to the MATLAB programs that implement these algorithms [32].) According to (41) and the boundary condition (23)