1. Introduction

In recent decades, researchers have paid great attention to fractional calculus and its applications. It has been used successfully for modelling many phenomena in different areas of sciences and engineering as quantum physics, continuum mechanics, viscoelastic and viscoplastic flow, electrical circuits, control theory, image processing, viscoelasticity, biology and hydrodynamics (see, for example, [¹-¹³]). Historically, the emergence of this type of calculus was at the end of the seventeenth century, when Gottfried Leibniz sent a letter to L’Hospital which he raised a question about the possibility of the meaning of derivatives with integer order in which be generalized to derivatives with non-integer orders [¹⁴]. He then developed this calculus through several approaches, including Riemann, Liouville [¹⁵, ¹⁶], Caputo [⁸, ¹⁷, ¹⁸], and others, see references [¹⁹-²¹]. In several papers, authors have used the Caputo or Riemann-Liouville (RL) derivatives to describe the models. However, these derivative has an weakness because of singular kernel. To handle this weakness, Caputo and Fabrizio [²², ²³] has proposed a non-local derivative. This indicates that the Caputo-Fabrizio derivative has a non-singular exponential kernel. One of the advantages of using non-singular kernels in fractional calculus is that they provide a more well-behaved behavior than singular kernels, which can lead to better convergence properties and more precise solutions. Additionally, non-singular kernels can provide a more natural and physically relevant description of certain physical phenomena, such as diffusion. This is what motivated many researchers to use it to solving many equations and modelling many phenomena in various branches of science, for instance, analysis of logistic equation [²⁴], Korteweg-de Vries-Burgers Equation [²⁵], nonlinear Fisher’s reaction diffusion equation [²⁶], Heat transfer analysis of Walters’-B fluid with Newtonian heating through an oscillating vertical plate [²⁷], Bose-Einstein condensation [²⁸] and other works see [²⁹-³²].

For a function c(q) (q > 0) belongs to the Sobolev space, Caputo and Fabrizio [²², ²³] proposed the following fractional derivative:

where M(ν) is a normalization constant depending on ν(0 < ν < 1) and c ^’ (p) is the derivative of c(p).

Nieto and Losada [²³], suggested a particular method that enables to find the normalized function. This method depends on that the fractional integral in Eq. (1) is the average of the function and its antiderivative for 0 < ν < 1. Therefore the normalized function by using this way, takes the following form:

So that the fractional Caputo-Fabrizio derivative of a function c(q) can redefined as:

From the Eq. (3) the solution, by using new fractional derivative, of the following fractional differential equation,

is given by Losada as follows [²³]:

The main objective of our work is to solve space-time fractional diffusion equation by using Caputo-Fabrizio fractional derivative in a semi-infinite medium with an application exaample of this solution, it is represented by the diffusion of arsenic in germanium under intrinsic conditions, where the diffusion coefficient is concentration-independent and remains constant.

So, our paper is organized as the following: in the next section we present the theoretical model and the associated mathematical calculation of the fractional partial differential equation. In particular, we have introduced the Laplace transform to solve this type of equations. In Sec. 3, we take the diffusion of arsenic in germanium as an application of our obtained solution and we will give our results summarizing some numerical simulations and discuss the different cases corresponding to different values of the fractionality degrees α and β. The conclusion is reported in Sec. 4.

2. Theoretical implementation and calculations

The transport phenomena in the material, can be modeled by the space-time fractional diffusion equation, this equation is referred to as:

where the positive‚ λ _α,β is the diffusion coefficient, 0 < α < 1 and 0 < β < 1and x and t are dimensionless variables, whereas λ _α,β can be assimilated to a diffusion coefficient and is also dimensionless parameter. At the same time, the modeling parameter λ _α,β is related to the flow-flow correlation function [³³] and describes the mechanisms of particle diffusion processes. The parameters λ _α,β can also be associated with the nature of the interaction of particles and their distribution. For example ”x” can be as x = X=L where X is the spatial location and L a characteristic length, like a small depth in the diffusion medium, whereas ”t” can be as t = T / τ where T is the real time and τ a characteristic time, as a short period in the temporal range of the studied phenomenon.

We consider diffusion in a semi-infinite medium, x > 0, when the surface is kept at a constant concentration C _s . We need a solution of the space-time fractional diffusion equation satisfying the following initial condition and the boundary condition:

By multiplying both sides of Eq. (6) by exp(-st) and integrating, we obtain:

then

From the initial condition we have c(x, 0) = 0 therefore, the last equation write as:

We put

therefore

The solutions of the Eqs. (13) and (14) are given by (Losada) as [²²]:

The first and second derivatives of the last two equations give:

and

From the Eqs. (17), (18), (19) and (20) we have the following ordinary differential equation:

Then, the solution of the Eq. (21) which is finite as x → ∞ whatever α and β is:

As mentioned before the Eq. (1) and becquse the factor in front of x is negative in the last expresion, the spatial variable x must be positive, else, if x goes to + ∞, the concentration c(x, s) diverges what has no physical sense.

Using the boundary condition we observe that:

From the Equation (22), we get:

Then, the solution is giving by the following expression:

therefore

where L ^-1 denoted the inverse Laplace transform. In the case of α = 1 and β = 1(corresponding to the standard case), we get:

Reference to a table of Laplace transforms [³⁴] shows that the function whose transform is given by Eq. (28) is the complementary error function:

which is the solution of the standard diffusion equation Dt1c(x,t)=λDx2c(x,t) in a semi-infinite medium, with the given boundary and intial conditions see Ref. [³⁵]. Erfc(.) is the complementary error function.

3. Application to intrinsic arsenic diffusion in germanium

Recently germanium (Ge) has emerged as a promising candidate in order to improve the complementary metal oxide semiconductors (CMOS) devices, because of its distinctive physical characteristics such as high intrinsic carrier mobility, small band gap and possible monolithic integration with silicon (Si) based devices [³⁶-³⁸]. To create (Ge) based devices, P-type and N-type germanium are brought into direct contact with each other. P-type germanium is doped with elements such as boron and gallium, while N-type germanium is doped with elements such as phosphorus and arsenic. Germanium can be doped using various methods, including dopant diffusion and ion implantation. The basis of diffusion doping way is the injection of a dopants by thermal diffusion [³⁸, ³⁹]. We are concerned with the intrinsic arsenic diffusion in germanium, which we present as an example applied to our mathematical treatment of diffusion equation. In this section, we simulate an profile of arsenic diffusion in germanium taken from experimental data of Brotzmann and Bracht [⁴⁰] througt the obtained solution. The diffusion satisfies same the initial condition and the boundary condition shown in the two Eqs. (7) and (8) above.

The arsenic intrinsic diffusion in germanium, can be modelled by the Eq. (6), then the simulated profiles we can get from the solution that we got in the Eq. (26).

We have accomplished a program by FORTRAN language, where we relied on the numerical method of Gaver-Stehfest to find the inverse Laplace transform in the Eq. (26).

This method was developed in the late 1960s. It is very simple numerical inverse Laplace transform method which has been used in such diverse areas as chemistry, economics, mathematics, computational physics and engineering. For a function f: (0 ,∞) → R , such that its Laplace transform:

the Gaver-Stehfest method, transformed the Laplace space into the time domain as follows [⁴¹-⁴³]:

where F (.) is the Laplace transform of f(t):

Applying this method to Eq. (26), we get

where C(x, n ln 2/t) is given by Eq. (25). The coefficients K _n depends only on the (necessarily even) number of expansion terms, M, given by:

The parameter M is the number of terms used in Eq. (31). M must be an even integer and should be choosen by trial and error method.

The parameter M plays an important role in the accuracy of the obtained results. In our study, we evaluated it with 10 on the basis of a convergence test that we conducted by comparing an analytical curve obtained from the standard solution in Eq. (28). and the solution obtained by using the Gaver-Stehfest method for α = 1, β = 1 in Eq. (31).

The Fig. 1 shows a comparison of the curve obtained through the numerical solution with the curve of the analytical solution shown in Eq. (28) for (α = 1, β = 1).

Figure 1 Comparing the numerical solution with the analytical solution, where (α = 1, β = 1), t = 18000, (λ _1,1 = 3.35 × 10 ^-11 ). 

We choose the values of λ _α,β for different orders of fractional derivation (α, β), indicated in the Table 1. which bring the coincidence of simulated profile and experimental profile.

The way we fit the data is to give a certain expected value of λ _α,β and compare the obtained data with the experimental data. We change the value each time until we reach the best fit between them. All profiles, we have plotted with dimensionless parameters. We took the dimensionless concentration c(x, t) = C(x, t)/C where C(x, t) is concentration and C = 1 cm ^-3 . For dimensionless variables, x = X/L and t = T/τ and we took L = 1 μm and τ = 1 s.

Figures 2 and 3 show the simulated profiles (line) we have obtained, at the indicated parameters λ _α,β and for different space and time Fractional-order and compare them with the experimental profile(symbol) of arsenic diffusion in Ge. The simulation profile (line) shown in Fig. 2, represents the different simulation profiles, which were completely identical according to the different derivation orders and according to the corresponding diffusion coefficient value that we choose, which leads to the best concordance with the experimental profile.

Figure 2 Simulated profiles C (x, t = 18000) (solid line) for α = {0.2, 0.4, 0.6, 0.8, 1}, β = 1 and experimental profile (symbol) of arsenic diffusion in Ge. 

Figure 3 Simulated profiles C (x, t = 18000) (solid line) for α = 1, β = {0.2, 0.4, 0.6, 0.8, 1} and experimental profile (symbol) of arsenic diffusion in Ge. 

The simulated profiles and the experimental profile concordance, for different space and time fractional-order, in Figs. 2 and 3 confirm that the fractional derivative orders affect on the estimation of the diffusion coefficient.

We see from Figs. 2 and 3 for every binary combination (α, β), there exists a diffusion coefficient that shows a satisfactory agreement with the experimental diffusion profile.

Figures 4 and 5 show that the fractional derivative order α and β, respectively, affects the estimation of the diffusion coefficient.

Figure 4 The effect of fractional derivative order α on the estimation of the diffusion coefficient for β = 1. 

Figure 5 The effect of fractional derivative order β on the estimation of the diffusion coefficient for α = 1. 

It is clear that increasing α reduces the value of the diffusion coefficient, while increasing the space fractional derivative order β increases the value of the diffusion coefficient.

4. Conclusion

In this work, we have treated the space-time fractional diffusion equation in a semi-infinite medium by using the fractional derivative of Caputo-Fabrizio. We have introduced the Laplace transform to solve this type of equations. Based on the obtained solution, we have simulated an profile of arsenic diffusion in germanium. Accurate simulations have been achieved showing that the fractional derivative orders affect on the estimation of the diffusion coefficient, where increasing the time fractional derivative order α reduces the value of the diffusion coefficient, while increasing the space fractional derivative order β increases the value of the diffusion coefficient. This result opens research perspectives to find a model for calculating the diffusion coefficient if fractional derivation is used.