Axisymmetric Solutions to Time-Fractional Heat Conduction Equation in a Half-Space under Robin Boundary Conditions

Abstract

The time-fractional heat conduction equation with the Caputo derivative of the order 0<𝛼<2 is considered in a half-space in axisymmetric case under two types of Robin boundary condition: the mathematical one with the prescribed linear combination of the values of temperature and the values of its normal derivative and the physical condition with the prescribed linear combination of the values of temperature and the values of the heat flux at the boundary.

1. Introduction

The generalized Fourier law, the time-nonlocal dependence between the heat flux vector 𝐪 and the temperature gradient grad𝑇 with the “long-tail memory” power kernel [1, 2] (see also [3]) as 𝐪(𝑡)=−𝑘Γ(𝛼)𝜕𝜕𝑡^𝑡0𝐪(𝑡)=−𝑘(𝑡−𝜏)^𝛼−1grad𝑇(𝜏)𝑑𝜏,0<𝛼≤1;Γ(𝛼−1)^𝑡0(𝑡−𝜏)^𝛼−2grad𝑇(𝜏)𝑑𝜏,1<𝛼≤2,(1.1) where Γ(𝛼) is the gamma function, can be interpreted in terms of the fractional calculus: 𝐪𝐪((𝑡𝑡))==−−𝑘𝑘𝐷𝐼^{1−𝛼𝑅𝐿}_𝛼−1grad𝑇(𝑡),0<𝛼≤1,grad𝑇(𝑡),1<𝛼≤2,(1.2) and in combination with the law of conservation of energy as 𝜌_𝑚𝐶𝑑𝑇𝑑𝑡=−div𝐪(1.3) results in the time-fractional heat conduction equation with the Caputo fractional derivative 𝜕^𝛼𝑇𝜕𝑡^𝛼=𝑎Δ𝑇,0<𝛼≤2.(1.4) Here, 𝜌_𝑚 is the mass density, 𝐶 denotes the specific heat capacity, and 𝑎=𝑘/(𝜌_𝑚𝐶) is the thermal diffusivity coefficient.

Recall that the Riemann-Liouville fractional integral 𝐼^𝛼𝑓(𝑡) and derivative 𝐷^𝛼𝑅𝐿𝑓(𝑡) are defined as follows (see [4–6]): 𝐼^𝛼𝑓(𝑡)=1Γ(𝛼)^𝑡0𝐷(𝑡−𝜏)^𝛼−1𝑓(𝜏)𝑑𝜏,𝛼>0,(1.5)_𝛼𝑅𝐿𝑓(𝑡)=𝑑^𝑚𝑑𝑡^𝑚Γ(𝑚1−𝛼)^𝑡0(𝑡−𝜏)^{𝑚−𝛼−1}𝑓(𝜏)𝑑𝜏,𝑚−1<𝛼<𝑚,(1.6) whereas the Caputo fractional derivative has the following form [5–7]: 𝑑^𝛼𝑓(𝑡)𝑑𝑡^𝛼=1Γ(𝑚−𝛼)^𝑡0(𝑡−𝜏)^{𝑚−𝛼−1}𝑑^𝑚𝑓(𝜏)𝑑𝜏^𝑚𝑑𝜏,𝑚−1<𝛼<𝑚.(1.7)

A detailed explanation of derivation of time-fractional heat conduction equation (1.4) from the constitutive equations (1.2) and the law of conservation of energy (1.3) can be found in [8]. Here we briefly present the main idea. In the case 0<𝛼≤1, as a consequence of (1.2), (1.3), and (1.6), we have 𝜕𝑇𝜕𝑡=𝑎𝜕𝜕𝑡Γ1(𝛼)^𝑡0(𝑡−𝜏)^𝛼−1Δ𝑇(𝜏)𝑑𝜏,(1.8) or after integration with respect to time as 𝑇(𝑡)−𝑇(0)=𝑎𝐼^𝛼Δ𝑇.(1.9) Applying to both sides of (1.9) the Caputo derivative 𝜕^𝛼/𝜕𝑡^𝛼, and taking into account that for 𝛼>0 [6], 𝜕^𝛼𝜕𝑡^𝛼𝐼^𝛼𝑇(𝑡)=𝑇(𝑡),(1.10) we obtain the time-fractional heat conduction equation (1.4) for 0<𝛼≤1.

Similarly, for 1<𝛼≤2, we get 𝜕𝑇𝜕𝑡=𝑎𝐼^𝛼−1Δ𝑇.(1.11) Applying 𝜕^𝛼−1/𝜕𝑡^𝛼−1 to both sides of (1.11) gives (1.4) for 1<𝛼≤2 as the following equality fulfills [5] 𝜕^𝛼−1𝜕𝑡^𝛼−1𝜕𝑇𝜕𝑡=𝜕^𝛼𝑇𝜕𝑡^𝛼.(1.12)

If the heat conduction equation is investigated in a bounded domain, the boundary conditions should be prescribed. The mathematical Robin boundary condition is a specification of a linear combination of the values of temperature and the values of its normal derivative at the boundary of the considered domain 𝑐¹𝑇+𝑐₂𝜕𝑇𝜕𝑛|||_𝑆=𝐹₀𝐱_𝑆,𝑡,(1.13) with some nonzero constants 𝑐₁ and 𝑐₂, while the physical Robin boundary condition specifies a linear combination of the values of temperature and the values of the heat flux at the boundary of the domain. For example, the condition of convective heat exchange between a body and the environment with the temperature 𝑇_𝑒𝐪⋅𝐧||_𝑆=ℎ𝑇||_𝑆−𝑇_𝑒,(1.14) where ℎ is the convective heat transfer coefficient, leads to ℎ𝑇+𝑘𝐷^{1−𝛼𝑅𝐿}𝜕𝑇𝜕𝑛|||_𝑆ℎ𝑇+𝑘𝐼_𝛼−1=ℎ𝑇_𝑒𝐱_𝑆,𝑡,0<𝛼≤1,𝜕𝑇𝜕𝑛|||_𝑆=ℎ𝑇_𝑒𝐱_𝑆,𝑡,1<𝛼≤2.(1.15)

The literature on mathematical aspects concerning correctness of initial-boundary-value problems for time-fractional diffusion equation and form and properties of its solutions is quite extensive (see, e.g., [9–16], among others). Geometrical explanation of fractional calculus is given in [17–19].

There are only a few papers [20, 21] in which the fractional diffusion is investigated under the mathematical Robin boundary condition. In previous publications, problems for a cylinder [22] and a sphere [23] under mathematical and physical Neumann boundary conditions were considered. In the present paper, for the first time, the solutions to time-fractional heat conduction equation in a half-space are studied under both the mathematical and physical Robin boundary conditions. The Laplace integral transform with respect to time 𝑡, the Hankel transform with respect to the spatial coordinate 𝑟, and the sin-cos-Fourier transforms with respect to spatial coordinate 𝑧 are used. The solutions under the mathematical and physical Neumann boundary conditions are obtained as particular cases.

2. Mathematical Preliminaries

2.1. Laplace Transform

The Laplace transform is defined as ℒ{𝑓(𝑡)}=𝑓^∗(𝑠)=^∞0𝑓(𝑡)𝑒^−𝑠𝑡𝑑𝑡,(2.1) where 𝑠 is the transform variable.

The inverse Laplace transfrom is carried out according to the Fourier-Mellin formula: ℒ⁻¹𝑓^∗(𝑠)=𝑓(𝑡)=12𝜋𝑖^{𝑐+𝑖∞𝑐−𝑖∞}𝑓^∗(𝑠)𝑒^𝑠𝑡𝑑𝑠,𝑡>0,(2.2) where 𝑐 is a positive fixed number.

The Laplace transform of the Riemann-Liouville fractional integral of the order 𝛼 is carried out according to the formula similar to the Laplace transform of 𝑚-fold primitive of a function as ℒ{𝐼^𝛼𝑓(𝑡)}=1𝑠^𝛼𝑓^∗(𝑠).(2.3)

The Caputo derivative for its Laplace transform requires the knowledge of the initial values of the function 𝑓(0⁺) and its integer derivatives 𝑓^(𝑘)(0⁺) of the order 𝑘=1,2,…,𝑚−1ℒ𝑑^𝛼𝑓(𝑡)𝑑𝑡^𝛼=𝑠^𝛼𝑓^∗(𝑠)−^𝑚−1_𝑘=0𝑓^(𝑘)0⁺𝑠^{𝛼−1−𝑘},𝑚−1<𝛼<𝑚,(2.4) whereas the Riemann-Liouville derivative for its Laplace transform rule requires the knowledge of the initial values of the fractional integral 𝐼^𝑚−𝛼𝑓(0⁺) and its derivatives of the order 𝑘=1,2,…,𝑚−1ℒ𝐷^𝛼𝑅𝐿𝑓(𝑡)=𝑠^𝛼𝑓^∗(𝑠)−^𝑚−1_𝑘=0𝐷^𝑘𝐼^𝑚−𝛼𝑓0⁺𝑠^{𝑚−1−𝑘},𝑚−1<𝛼<𝑚.(2.5) The reader interested in applications of integral transforms in fractional calculus is referred to [4–7, 24].

2.2. Hankel Transform

The Hankel transform is used to solve problem in cylindrical coordinates in the domain 0≤𝑟<∞ and is defined as ℋℋ{𝑓(𝑟)}=𝑓^∗(𝜂)=^∞0𝑓(𝑟)𝐽_𝜈(𝜂𝑟)𝑟𝑑𝑟,₋₁𝑓^∗(𝜂)=𝑓(𝑟)=^∞0𝑓^∗(𝜂)𝐽_𝜈(𝜂𝑟)𝜂𝑑𝜂,(2.6) where 𝐽_𝜈(𝑟) is the Bessel function of the order 𝜈.

The following formula is fulfilled: ℋ𝑑²𝑓(𝑟)𝑑𝑟²+1𝑟𝑑𝑓𝑑𝑟(𝑟)−𝜈𝑟²²𝑓(𝑟)=−𝜂²𝑓^∗(𝜂).(2.7)

2.3. Sin-Cos-Fourier Transform

In the case of boundary condition of the third kind with the prescribed boundary value of linear combination of a function and its normal derivative 𝑧=00𝑥00899∶−𝑑𝑓𝑑𝑧+𝐻𝑓=𝜑⁰,(2.8) the following sin-cos-Fourier transform [25] is employed: ℱℱ{𝑓(𝑧)}=𝑓^∗(𝜉)=^∞0𝐾(𝑧,𝜉)𝑓(𝑧)𝑑𝑧,₋₁𝑓^∗(𝜉)=𝑓(𝑧)=𝜋2^∞0𝐾(𝑧,𝜉)𝑓^∗(𝜉)𝑑𝜉,(2.9) with the kernel 𝐾(𝑧,𝜉)=𝜉cos(𝑧𝜉)+𝐻sin(𝑧𝜉)√𝜉²+𝐻².(2.10)

Application of the sin-cos-Fourier transform to the second derivative of a function gives ℱ𝑑²𝑓(𝑧)𝑑𝑧²=−𝜉²𝑓^∗(𝜉)−𝜉√𝜉²+𝐻²𝑑𝑓𝑑𝑧(𝑧)−𝐻𝑓(𝑧)|||||_𝑧=0.(2.11)

3. Solution to the Problem under Mathematical Robin Boundary Condition

Consider the axisymmetric time-fractional heat conduction equation in cylindrical coordinates 𝜕^𝛼𝑇𝜕𝑡^𝛼=𝑎𝜕²𝑇𝜕𝑟²+1𝑟𝜕𝜕𝑇𝑟+𝜕²𝑇0≤𝑟<∞,0<𝑧<∞,𝜕𝑧0²<𝑡,<∞,0<𝛼≤2,(3.1) with zero initial conditions 𝑡𝑡==0000𝑥𝑥0000889999∶∶𝜕𝑇𝑇=0,0<𝛼≤2,𝜕𝑡=0,1<𝛼≤2,(3.2) and the mathematical Robin boundary condition 𝑧=00𝑥00899∶𝐻𝑇−𝜕𝑇𝜕𝑧=𝑓(𝑟,𝑡).(3.3)

The zero conditions at infinity are also assumed lim_𝑟→∞𝑇(𝑟,𝑧,𝑡)=0,lim_𝑧→∞𝑇(𝑟,𝑧,𝑡)=0.(3.4)

The solution to the initial-boundary-value problem (3.1)–(3.4) can be written as 𝑇=^𝑡0^∞0𝑓(𝜌,𝜏)𝒢_𝑚(𝑟,𝑧,𝜌,𝑡−𝜏)𝜌𝑑𝜌𝑑𝜏,(3.5) where 𝒢_𝑚(𝑟,𝑧,𝜌,𝑡) is the fundamental solution being the solution of the following problem: 𝜕^𝛼𝒢_𝑚𝜕𝑡^𝛼=𝑎𝜕²𝒢_𝑚𝜕𝑟²+1𝑟𝜕𝜕𝒢𝑟^𝑚+𝜕²𝒢_𝑚0≤𝑟𝑡<=∞0,0𝑥000<89𝑧9<∶∞𝒢,_𝑚0<𝜕𝑡𝑧<²∞,,0<𝛼≤2,𝑡=00𝑥00899∶𝜕𝒢_𝑚=0,0<𝛼≤2,𝑧=00𝑥00899∶𝐻𝒢_𝑚𝜕𝑡=0,1<𝛼≤2,−𝜕𝒢^𝑚𝜕𝑧=1𝑟𝛿(𝑟−𝜌)𝛿(𝑡).(3.6) Here, 𝛿(𝑟) is the Dirac delta function.

The Laplace transform with respect to time 𝑡, the Hankel transform with respect to the spatial coordinate 𝑟, and the sin-cos-Fourier transform with respect to the spatial coordinate 𝑧 result in 𝒢^{𝑚∗∗∗}=𝑎𝜉𝐽⁰(𝜌𝜂)√𝜉²+𝐻²1𝑠^𝛼+𝑎𝜉²+𝜂²,(3.7) where each of the integral transforms is denoted by the asterisk.

Invertion of the integral transforms leads to 𝒢_𝑚(𝑟,𝑧,𝜌,𝑡)=2𝑎𝑡^𝛼−1𝜋^∞0𝜉²cos(𝑧𝜉)+𝜉𝐻sin(𝑧𝜉)𝜉²+𝐻²×𝐸_𝛼,𝛼−𝑎𝜉²+𝜂²𝑡^𝛼𝐽₀(𝑟𝜂)𝐽₀(𝜌𝜂)𝜂𝑑𝜉𝑑𝜂.(3.8) Here, 𝐸_𝛼,𝛽 is the Mittag-Leffler function in two parameters 𝛼 and 𝛽 defined by the series representation: 𝐸_𝛼,𝛽(𝑧)=^∞_𝑚₌₀𝑧^𝑚Γ(𝛼𝑚+𝛽),𝛼>0,𝛽>0,𝑧∈𝐶.(3.9) The essential role of the Mittag-Leffler function in fractional calculus results from the following formula for the inverse Laplace transform ℒ⁻¹𝑠𝑠^𝛼+^𝛼−𝛽𝑏=𝑡^𝛽−1𝐸_𝛼,𝛽(−𝑏𝑡^𝛼).(3.10) Consider several particular cases of solution (3.8).

The case 𝐻=0 corresponds to the mathematical Neumann boundary condition with the prescribed boundary value of the normal derivative of temperature, and the solution reads 𝒢_𝑚(𝑟,𝑧,𝜌,𝑡)=2𝑎𝑡^𝛼−1𝜋^∞0×cos(𝑧𝜉)𝐸𝐽_𝛼,𝛼0−𝑎𝜉²+𝜂²𝑡^𝛼(𝑟𝜂)𝐽₀(𝜌𝜂)𝜂𝑑𝜉𝑑𝜂.(3.11)

In the case of classical heat conduction (𝛼=1), the Mittag-Leffler function: 𝐸_1,1−𝑥²=𝑒^−𝑥2,(3.12) and taking into account the following integrals [26, 27]: ^∞0𝑥𝑒^−𝑎𝑥2𝐽₀(𝑏𝑥)𝐽₀(𝑐𝑥)𝑑𝑥=12𝑎exp−𝑏²+𝑐²4𝑎𝐼⁰2𝑏𝑎𝑐,𝑎>0,_∞0𝑒^−𝑎𝑥2cos(𝑐𝑥)𝑑𝑥=√𝜋2√𝑎exp−4𝑐𝑎²,𝑎>0,_∞01𝑥²+𝑏²=𝜋𝑒^−𝑎𝑥2cos(𝑐𝑥)𝑑𝑥4𝑏𝑒^𝑎𝑏2𝑒^−𝑏𝑐erfc√𝑎𝑏−𝑐2√𝑎+𝑒^𝑏𝑐erfc√𝑎𝑏+𝑐2√𝑎,𝑎>0,𝑏>0,_∞0𝑥𝑥²+𝑏²=𝜋𝑒^−𝑎𝑥2sin(𝑐𝑥)𝑑𝑥4𝑒^𝑎𝑏2𝑒^−𝑏𝑐erfc√𝑎𝑏−𝑐2√𝑎−𝑒^𝑏𝑐erfc√𝑎𝑏+𝑐2√𝑎,𝑎>0,𝑏>0,(3.13) where 𝐼₀(𝑥) is the modified Bessel function, erfc𝑥 is the complementary error function, we get 𝒢_𝑚(𝑟,𝑧,𝜌,𝑡)=12𝑡exp−𝑟²+𝜌²+𝑧²×⎧⎪⎨⎪⎩14𝑎𝑡𝐼⁰2𝑟𝑎𝜌𝑡√𝜋𝑎𝑡−𝐻exp⎡⎢⎢⎣√𝑎𝑡𝐻+𝑧2√𝑎𝑡²⎤⎥⎥⎦erfc√𝑎𝑡𝐻+𝑧2√𝑎𝑡⎫⎪⎬⎪⎭.(3.14)

For the wave equation (𝛼=2), 𝐸_2,2−𝑥²=sin𝑥𝑥,(3.15) and after evaluation of the necessary integrals appearing in solution (3.8) (under assumption 𝜌=0) [26, 27] ^∞0𝑥√𝑥²+𝑝²sin𝑏√𝑥²+𝑝²𝐽⁰(𝑐𝑥)𝑑𝑥=⎧⎪⎨⎪⎩√𝑏²1−𝑐²0,cos𝑝√𝑏²−𝑐²,00<<𝑐𝑏<<𝑏𝑐,._∞0cos(𝑝𝑥)cos(𝑞𝑥)𝑥²+𝑏²𝑑𝑥=⎧⎪⎨⎪⎩2𝜋𝑏𝑒^−𝑏𝑝𝜋cosh(𝑏𝑞),0<𝑞<𝑝,2𝑏𝑒^−𝑏𝑞cosh(𝑏𝑝),0<𝑝<𝑞,𝑏>0,_∞0𝑥cos(𝑝𝑥)sin(𝑞𝑥)𝑥²+𝑏²𝑑𝑥=⎧⎪⎨⎪⎩−𝜋2𝑒^−𝑏𝑝𝜋sinh(𝑏𝑞),0<𝑞<𝑝,2𝑒^−𝑏𝑞cosh(𝑏𝑝),0<𝑝<𝑞,𝑏>0,(3.16) we obtain the solution 𝒢_𝑚(𝑟=,𝑧⎧⎪⎪⎨⎪⎪⎩,0,𝑡)√𝑎√𝑎𝑡²−𝑟²𝛿√𝑎𝑡²−𝑟²−𝑧−𝐻𝑒^{√−𝐻(𝑎𝑡2−𝑟2−𝑧)}00<<𝑧𝑟<<√√𝑎𝑡,0√𝑎𝑡²−𝑟²,𝑎𝑡<𝑟<∞,√𝑎𝑡²−𝑟²<𝑧<∞(3.17)

Of particular interest is also the case 𝛼=1/2 for which 𝐸_1/2,1/2(−𝑥)=1√𝜋−𝑥𝑒^𝑥2erfc𝑥=2√𝜋^∞0𝒢𝑒^{−𝑢2−2𝑥𝑢}𝑢𝑑𝑢,𝑥≥0,_𝑚(𝑟,𝑧,𝜌,𝑡)=12√𝜋𝑡^∞0exp−𝑢²−𝑟²+𝜌²+𝑧²×⎧⎪⎨⎪⎩18𝑎√𝑡𝑢𝐼⁰4𝑎𝑟√𝜌𝑡𝑢√2𝜋𝑎𝑢𝑡^1/4−𝐻exp⎡⎢⎢⎣√2𝑎𝑢𝐻𝑡^1/4+𝑧2√2𝑎𝑢𝑡^1/4²×erfc√⎤⎥⎥⎦2𝑎𝑢𝐻𝑡^1/4+𝑧2√2𝑎𝑢𝑡^1/4𝑑𝑢.(3.18)

4. Solution to the Problem under Physical Robin Boundary Condition

Consider the following axisymmetric time-fractional heat conduction equation: 𝜕^𝛼𝑇𝜕𝑡^𝛼=𝑎𝜕²𝑇𝜕𝑟²+1𝑟𝜕𝜕𝑇𝑟+𝜕²𝑇0≤𝑟<∞,0<𝑧<∞,𝜕0𝑧²<𝑡,<∞,0<𝛼≤2,(4.1) under zero initial conditions 𝑡𝑡==0000𝑥𝑥0000889999∶∶𝜕𝑇𝑇=0,0<𝛼≤2,𝜕𝑡=0,1<𝛼≤2,(4.2) and the physical Robin boundary condition 𝑧=00𝑥00899∶𝐻𝑇−𝐷^{1−𝛼𝑅𝐿}𝜕𝑇𝑧=00𝑥00899∶𝐻𝑇−𝐼_𝛼−1𝜕𝑧=𝑓(𝑟,𝑡),0<𝛼≤1,𝜕𝑇𝜕𝑧=𝑓(𝑟,𝑡),1<𝛼≤2,(4.3) where 𝐻=ℎ/𝑘.

The solution to the initial-boundary-value problem (4.1)–(4.3) has the form 𝑇=^𝑡0^∞0𝑓(𝜌,𝜏)𝒢_𝑝(𝑟,𝑧,𝜌,𝑡−𝜏)𝜌𝑑𝜌𝑑𝜏,(4.4) where the fundamental solution 𝒢_𝑝(𝑟,𝑧,𝜌,𝑡) fulfills the following equation 𝜕^𝛼𝒢_𝑝𝜕𝑡^𝛼=𝑎𝜕²𝒢_𝑝𝜕𝑟²+1𝑟𝜕𝜕𝒢𝑟^𝑝+𝜕²𝒢_𝑝0≤𝑟<∞,0<𝑧<∞,0<𝜕𝑡𝑧<²∞,,0<𝛼≤2,(4.5) under the following conditions: 𝑡𝑡==0000𝑥𝑥0000889999∶∶𝜕𝒢𝒢_𝑝𝑝=0,0<𝛼≤2,𝑧=0∶𝐻𝒢_𝑝𝜕𝑡=0,1<𝛼≤2,−𝐷^{1−𝛼𝑅𝐿}𝜕𝒢_𝑝𝑧=0∶𝐻𝒢_𝑝𝜕𝑧=1𝑟𝛿(𝑟−𝜌)𝛿(𝑡),0<𝛼≤1,−𝐼^𝛼−1𝜕𝒢_𝑝𝜕𝑧=1𝑟𝛿(𝑟−𝜌)𝛿(𝑡),1<𝛼≤2.(4.6)

The Laplace transform with respect to time 𝑡 leads to the following equation: 𝑠^𝛼𝒢^∗𝑝=𝑎𝜕²𝒢^∗𝑝𝜕𝑟²+1𝑟𝜕𝜕𝒢𝑟^∗𝑝+𝜕²𝒢^∗𝑝𝜕𝑧²,0≤𝑟<∞,0<𝑧<∞,0<𝛼≤2,(4.7) and the following boundary condition: 𝑧=0∶𝑠^𝛼−1𝐻𝒢^∗𝑝−𝜕𝒢^∗𝑝𝜕𝑧=1𝑟𝛿(𝑟−𝜌)𝑠^𝛼−1.(4.8)

In this case, the kernel of the sin-cos-Fourier transform with respect to the spatial coordinate 𝑧 depends on the Laplace transform variable 𝑠𝐾(𝑧,𝜉,𝑠)=𝜉cos(𝑧𝜉)+𝑠^𝛼−1𝐻sin(𝑧𝜉)𝜉²+𝑠^𝛼−1𝐻²,(4.9) and in the transform domain we obtain 𝒢^{𝑝∗∗∗}=𝑎𝜉𝐽⁰(𝜌𝜂)𝜉²+𝑠^𝛼−1𝐻²𝑠^𝛼−1𝑠^𝛼+𝑎𝜉²+𝜂².(4.10)

Inversion of the Laplace transform in (4.10) depends on the value of 𝛼. For 0<𝛼≤1, we have 𝒢_𝑝(𝑟,𝑧,𝜌,𝑡)=2𝑎𝜋^∞0^𝑡0×(𝑡−𝜏𝜂)𝐽_0−𝛼(𝑟𝜂)𝐽₀(𝜌𝜂)𝜏^𝛼−1𝐸_𝛼,𝛼−𝑎𝜉²+𝜂²𝜏^𝛼𝐸_{2−2𝛼,1−𝛼}−𝐻𝜉²²(𝑡−𝜏)^2−2𝛼+𝐻cos(𝑧𝜉)𝜉(𝑡−𝜏)^1−2𝛼𝐸_{2−2𝛼,2−2𝛼}−𝐻𝜉²²(𝑡−𝜏)^2−2𝛼sin(𝑧𝜉)𝑑𝜏𝑑𝜉𝑑𝜂,(4.11) whereas for 1<𝛼≤2, we get 𝒢_𝑝(𝑟,𝑧,𝜌,𝑡)=2𝑎𝜋^∞0^𝑡0×𝜉₂𝐽₀(𝑟𝜂)𝐽₀(𝜌𝜂)𝐸_𝛼−𝑎𝜉²+𝜂²𝜏^𝛼𝐻²(𝑡−𝜏)^2𝛼−3𝐸_{2𝛼−2,2𝛼−2}−𝐻𝜉²²(𝑡−𝜏)^2𝛼−2+𝜉cos(𝑧𝜉)𝐻(𝑡−𝜏)^𝛼−2𝐸_{2𝛼−2,𝛼−1}−𝐻𝜉²²(𝑡−𝜏)^2−2𝛼sin(𝑧𝜉)𝑑𝜏𝑑𝜉𝑑𝜂.(4.12) The fundamental solution under physical Neumann boundary condition is obtained for𝐻=0 as 𝒢_𝑝(𝑟,𝑧,𝜌,𝑡)=2𝑎𝜋^∞0𝐸_𝛼−𝑎𝜉²+𝜂²𝑡^𝛼𝐽₀(𝑟𝜂)𝐽₀(𝜌𝜂)𝜂cos(𝑧𝜉)𝑑𝜉𝑑𝜂.(4.13)

Formulae (4.11) and (4.12) simplify for 𝛼=1/2 and 𝛼=3/2, respectively, taking into account (3.12) and that 𝐸_1,1/2−𝑥²=1√𝜋[1−2𝑥𝐷(𝑥)],(4.14) where 𝐷(𝑥)=𝑒^−𝑥2^𝑥0𝑒^𝑢2𝑑𝑢(4.15) is the Dawson integral.

5. Conclusion

We have derived the analytical solutions to time-fractional heat conduction equation in a half-space under mathematical and physical Robin boundary conditions. The integral transform technique has been used. It should be emphasized that in the case of physical Robin boundary condition, the order of integral transforms is important as the kernel of the sin-cos-Fourier transform depends on the Laplace transform variable. The limiting case 𝐻=0 corresponds to solutions of problems under mathematical and physical Neumann boundary conditions with the prescribed boundary value of the normal derivative and with the prescribed boundary value of the heat flux, respectively. The difference between mathematical and physical boundary conditions (as well as the difference between the solutions) disappears in the case of standard heat conduction equation (𝛼=1).