Analytical Solutions of the Space-Time Fractional Derivative of Advection Dispersion Equation

Abstract

Fractional advection-dispersion equations are used in groundwater hydrology to model the transport of passive tracers carried by fluid flow in porous medium. A space-time fractional advection-dispersion equation (FADE) is a generalization of the classical ADE in which the first-order space derivative is replaced with Caputo or Riemann-Liouville derivative of order 0<𝛽≤1, and the second-order space derivative is replaced with the Caputo or the Riemann-Liouville fractional derivative of order 1<𝛼≤2. We derive the solution of the new equation in terms of Mittag-Leffler functions using Laplace transfrom. Some examples are given. The results from comparison let no doubt that the FADE is better in prediction than ADE.

1. Introduction

The description of transport is closely related to the terms convection, diffusion, dispersion, and retardation as well as decomposition. First, it is assumed that there are no interactions between the species dissolved in water and the surrounding solid phase [1]. The primary mechanism for the transport of improperly discarded hazardous waste through the environment is by the movement of water through the subsurface and surface waterways. Studying this movement requires that one must be able to measure the quantity of waste present at a particular point in space time. The measure universally for chemical pollution is the concentration. Analytical methods that handle solute transport in porous media are relatively easy to use [1]. However, because of complexity of the equations involved, the analytical solutions are generally available restricted to either radial flow problems or to cases where velocity is uniform over the area of interest. Numerous analytical solutions are available for time-dependent solute transport within media having steady state and uniform flow. This work is devoted to the discussion underpinning the derivation of the analytical solution of space-time fractional derivative of advection-dispersion equation.

2. Governing Equations

A relatively complete set of one-dimensional analytical solutions for convective-dispersive solute equations has been recently published by Van Genuchten and Alves in 1982 [2]. Here we shall review a case having a practical application.

Let us consider a one-dimensional model consisting of infinitely ling homogenous isotropic porous media with steady state uniform flow with seepage velocity 𝑣. We inject a particular chemical from one end of the model for a period of time 𝑡₀ such that the input concentration varies as an exponential function of time [3]. The value of that chemical concentration at any time 𝑡 and at a distance 𝑥 from the injection boundary, allowing for the decay and adsorption, may be obtained from the solution of the following set of equations [3]: 𝐷𝜕²𝑐(𝑥,𝑡)𝜕𝑥²()−𝑣𝜕𝑐𝑥,𝑡()𝜕𝑥−𝜆𝑅𝑐=𝑅𝜕𝑐𝑥,𝑡,𝜕𝑡 (1) where 𝐷 is the dispersion coefficient and 𝑅 the retardation factor, subject to the initial condition: 𝑐()𝑐()𝑥,𝑡=0𝑡=0,0,𝑡=𝑐₀exp−𝛾𝑡0<𝑡≤𝑡₀, (2) which means that the system is initially free of that chemical, 𝛾 and 𝑐₀ are constants and boundaries conditions 𝜕𝑐(𝑥,𝑡)𝜕𝑥=0𝑥⟶∞. (3) This indicates that the concentration of the gradient at the other end of the model remains unchanged. Note that the standard version of advection-dispersion equation does not allow for predicting the mass transform through the geological formation accurately; it is then important to investigate a possible analytical partial differential equation that can describe better this problem. In this work this possibility is further investigated for a rectangular symmetric form of (1), by replacing the classical first-order derivative of the concentration by a fractional derivative. Because the concepts of fractional (or noninteger) order derivatives may not be widely known, the concept is first briefly discussed below.

3. Fractional Calculus

Fractional calculus has been used to model physical and engineering processes, which are found to be best described by fractional differential equations. It is worth noting that the standard mathematical models of integer-order derivatives, including nonlinear models, do not work adequately in many cases. In the recent years, fractional calculus has played a very important role in various fields such as mechanics, electricity, chemistry, biology, economics, notably control theory, and signal and image processing. In the past several decades, the investigation of travelling-wave solutions for nonlinear equations has played an important role in the study of nonlinear physical phenomena [4–12]. The concept of fractional order derivatives for a function, say 𝑓(𝑥), is based on a generalization of the Abel integral: 𝐷^−𝑛𝑓(𝑥)=𝑓(𝑥)𝑑𝑥_𝑛=1Γ(𝑛)^𝑥₀(𝑥−𝑡)^𝑛−1𝑓(𝑡)𝑑𝑡, (4) where 𝑛 is a nonzero positive integer and Γ(⋅) is the Gamma function [13].

This represents an integral of order 𝑛 for the continuous function 𝑓(𝑥), whenever 𝑓 and all its derivatives vanish at the origin, 𝑥=0. This result can be extended to the concept of an integral of arbitrary order 𝑐, defined as 𝐷^−𝑐𝑓(𝑥)=𝐷^{−𝑗−𝑠}𝑓(𝑥)=1Γ(𝑐)^𝑥₀(𝑥−𝑡)^𝑐−1𝑓(𝑡)𝑑𝑡, (5) where 𝑐 is a positive real number, 𝑗 an integer such that 0<𝑠≤1.

Let 𝑝 now be the least positive integer larger than 𝛼 such that 𝛼=𝑚−𝜌; 0<𝜌≤1. Equation (4) can then be used to define the derivative of (positive) fractional order, say 𝛼, of a function 𝑓(𝑥) as 𝐷^𝑐𝑓(𝑥)=𝐷^𝑝−𝜌𝑓(𝑥)=1Γ𝜌^𝑥₀(𝑥−𝑡)^𝜌−1𝑑^𝑝𝑓(𝑡)𝑑𝑡^𝑝𝑑𝑡. (6) Note that these results, like Abel’s integral, are only valid subject to the condition that 𝑓^(𝑘)(𝑥)∣𝑥=0 for 𝑘=0,1,2,…,𝑝.

3.1. Properties

Properties of the operator can be found in [14, 15]; we mention only the following.

For 𝑓∈𝐶_𝜇,  𝜇≥−1,𝛼,𝛽≥0 and 𝛾>−1: 𝐷^−𝛼𝐷^−𝛽𝑓(𝑥)=𝐷^{−𝛼−𝛽}𝑓(𝑥),𝐷^−𝛼𝐷^−𝛽𝑓(𝑥)=𝐷^−𝛽𝐷^−𝛼𝑓(𝑥),𝐷^−𝛼𝑥^𝛾=Γ𝛾+1Γ𝑥𝛼+𝛾+1^𝛼+𝛾. (7)

3.2. Formulation of Space-Time Fractional Derivative of Hydrodynamic Advection-Dispersion Equation

In order to include explicitly the possible effect of flow geometry into the mathematical model, the Cartesian component of the gradient of concentration, 𝜕_𝑥𝑐(𝑥,𝑡) is replaced by the Riemann-Liouville fractional derivatives of order 𝛽,  𝜕^𝛽_𝑥𝑐(𝑥,𝑡) and 𝜕²_𝑥𝑐(𝑥,𝑡) is replaced by 𝜕^𝛼_𝑥𝑐(𝑥,𝑡) with 0<𝛽≤1<𝛼<𝑙, as follows: 𝐷𝜕^𝛼𝑐(𝑥,𝑡)𝜕𝑥^𝛼𝜕−𝑣^𝛽𝑐(𝑥,𝑡)𝜕𝑥^𝛽()−𝜆𝑅𝑐=𝑅𝜕𝑐𝑥,𝑡.𝜕𝑡 (8) This provides a generalized form of the classical equation governing the transport of the solute (1): this integrodifferential equation does contain the additional parameter 𝛼 and 𝛽, which can be viewed as new physical parameters that characterize the transport through the geological formations. The same transformation generates also a more general form for the boundary condition at the other end of the model: 𝜕^𝛽𝑐(𝑥,𝑡)𝜕𝑥^𝛽=0𝑥⟶∞. (9) Relations (8) and (9), together with the initial condition described in (2), represent a complete set of equations for which a solution exists. The integrodifferential character of the relations makes the search for analytical solutions for the problem very difficult however. In this work an analytical solution in series form will be discussed in the next section.

4. Analytical Solution

4.1. The Riemann-Liouville Derivative

The method here consists of applying the Laplace transform on both sides of (8) to have 𝐷𝜕^𝛼𝑐(𝑥,𝑠)𝜕𝑥^𝛼𝜕−𝑣^𝛽𝑐(𝑥,𝑠)𝜕𝑥^𝛽−𝑅(𝜆+𝑠)𝑐(𝑥,𝑠)=𝑅𝑐(𝑥,0) (10) with the initial condition (2) and further transformation; the above equation can then be in the following form: 𝜕^𝛼𝑐(𝑥,𝑠)𝜕𝑥^𝛼𝜕−𝜇^𝛽𝑐(𝑥,𝑠)𝜕𝑥^𝛽−𝜏𝑐(𝑥,𝑠), (11) where 𝑠 is the variable of Laplace for the time-component, 𝜇=𝑣/𝐷, and 𝜏=𝑅(𝜆+𝑠)/𝐷. Let 𝑐(𝑥,𝑠)=𝑦(𝑥); then (11) becomes 𝜕^𝛼𝑦(𝑥)𝜕𝑥^𝛼𝜕−𝜇^𝛽𝑦(𝑥)𝜕𝑥^𝛽−𝜏𝑦(𝑥)=0. (12) Applying the Laplace operator on both sides of (11), on the space component and replacing, we have the following equation [14]: ℒ𝑦𝑝=^𝑙_𝑖=1ℎ_𝑖𝑝^𝑖−1𝑝^𝛼−𝜇𝑝^𝛽,−𝜏 (13) where 𝑝 is the Laplace variable for the space component and ℎ_𝑖=𝜕_𝑥^𝛼−𝑖𝑐(0⁺).

For 𝑝∈ℂ and |𝜏𝑝^−𝛽/(𝑝^𝛼−𝛽−𝜇)|<1, we have the following expression 1/(𝑝^𝛼−𝜇𝑝^𝛽−𝜏) which can be written in the form of series as follows [14]: 𝑝^𝑖−1𝑝^𝛼−𝜇𝑝^𝛽−𝜏=𝑝^∞𝑖−1_𝑗=0𝜏^𝑛𝑝^{−𝛽−𝛽𝑛}𝑝^𝛼−𝛽−𝜇^𝑛+1. (14) And hence replacing the above expression in (13) yields the following representation: ℒ𝑦𝑝=²_𝑖=1ℎ_𝑖^∞_𝑗=0𝜏^𝑛𝑝^{−𝛽−𝛽𝑛+𝑖−1}𝑝^𝛼−𝛽−𝜇^𝑛+1. (15) The above expression is then simplified further, for 𝑝∈ℂ and |𝜇𝑝^𝛽−𝛼|<1, we have first 𝜏^𝑛𝑝^{−𝛽−𝛽𝑛+𝑖−1}𝑝^𝛼−𝛽−𝜇^𝑛+1=𝜏^𝑛𝑝^{(𝛼−𝛽)−(𝛼+𝛽𝑛−𝑖+1)}𝑝^𝛼−𝛽−𝜇^𝑛+1. (16) And secondly the above equation can now be expressed as follows: =1ℒ𝑥𝑛!^{𝛼𝑛+𝛼−𝑖}𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛼+𝑛𝛽+1−𝑖}𝜇𝑥^𝛼−𝛽, (17) where 𝜕𝜕𝑥^𝑛𝐸_𝛼,𝛽(𝑥)=^∞_𝑗=0Γ𝑛+𝑗+1Γ𝑥𝑛𝛼+𝛽+𝛼𝑗^𝑗.𝑗! (18) Hence the solutions of (12) can be given as follows: 𝑦_𝑖(𝑥)=^∞_𝑛=0𝜏^𝑛𝑥𝑛!^{𝛼𝑛+𝛼−𝑖}𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛼+𝛽𝑛+1−𝑖}𝜇𝑥^𝛼−𝛽. (19) Thus it follows that the solution of (12) is given as 𝑦(𝑥)=²_𝑖=1ℎ_𝑖𝑦_𝑖(𝑥), (20) so that 𝑐(𝑥,𝑠)=²_𝑖=1ℎ_𝑖𝑐_𝑖(𝑥,𝑠). (21) Thus the series solution of (8) can be now given by applying the inverse Laplace operator on 𝑐(𝑥,𝑠) to have 𝑐_𝑖(𝑥,𝑡)=ℒ^−1∞_𝑛=0𝜏^𝑛𝑥𝑛!^{𝛼𝑛+𝛼−𝑖}𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛼+𝛽𝑛+1−𝑖}𝜇𝑥^𝛼−𝛽. (22) Since the inverse Laplace operator is a linear operator, it follows that 𝑐_𝑖()=𝑥,𝑡_∞𝑛=0ℒ⁻¹𝜏^𝑛𝑥𝑛!^{𝛼𝑛+𝛼−𝑖}𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛼+𝛽𝑛+1−𝑖}𝜇𝑥^𝛼−𝛽. (23) Replacing 𝜏^𝑛=(𝑅(𝜆+𝑠)/𝐷)^𝑛=(𝑅/𝐷)^𝑛(𝜆+𝑠)^𝑛 so that ℒ⁻¹𝜏^𝑛=𝑅𝐷^𝑛ℒ⁻¹(𝜆+𝑠)^𝑛=𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ(),𝑐−𝑛_𝑖(𝑥,𝑡)=^∞_𝑛=0(𝑅/𝐷)^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛/Γ(−𝑛)𝑥𝑛!^{𝛼𝑛+𝛼−𝑖}×𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛼+𝛽𝑛+1−𝑖}𝜇𝑥^𝛼−𝛽,𝑐(𝑥,𝑡)=²_𝑖=1ℎ_𝑖𝑐_𝑖(),𝑐𝑥,𝑡₁(𝑥,𝑡)=^∞_𝑛=0(𝑅/𝐷)^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛/Γ(−𝑛)𝑥𝑛!^{𝛼𝑛+𝛼−1}×𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛼+𝛽𝑛}𝜇𝑥^𝛼−𝛽,𝑐₂(𝑥,𝑡)=^∞_𝑛=0(𝑅/𝐷)^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛/Γ(−𝑛)𝑥𝑛!^{𝛼𝑛+𝛼−2}×𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛼+𝛼𝑛−1}𝜇𝑥^𝛼−𝛽. (24) To find the coefficient ℎ_𝑖,  𝑖=1,2, we need to apply the boundaries and initial condition on 𝑐(𝑥,𝑡) which yields to ℎ_𝑖=𝑐₀2. (25)

Example 1. Our concern here is to consider (8) when 𝛼=2 and 0<𝛽≤1. Following the discursion presented earlier, the analytical solution of space-time fractional derivative of hydrodynamic advection-dispersion equation has its two solutions given by 𝑐₁(𝑥,𝑡)=^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^2𝑛+1×𝜕𝜕𝑥^𝑛𝐸_{2−𝛽,𝑛𝛽+1}𝑅𝐷𝑥^2−𝛽,𝑐₂(𝑥,𝑡)=^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^2𝑛×𝜕𝜕𝑥^𝑛𝐸_{2−𝛽,𝑛𝛽+1}𝑅𝐷𝑥^2−𝛽. (26) The above solutions form the fundamental system of solution when 𝛽<1.

4.2. The Caputo Derivative

The Riemann-Liouville derivative has certain disadvantages when trying to model realworld phenomena with fractional differential equations [16–19]. Therefore, we investigate the solution of space-time Caputo fractional derivative of hydrodynamic advection-dispersion equation.

For the Caputo derivative, the Laplace transform is based on the formula ℒ𝑐𝐷^𝛼𝑦(𝑠)=𝑠^𝛼ℒ𝑦(𝑠)−¹_𝑖=0ℎ_𝑖𝑠^{𝛼−𝑖−1} (27) with ℎ_𝑖=𝑦^𝑖(0)(𝑖=0,1). (28) Thus applying the Laplace transform in both sides of (8) on the component of time and applying again the Laplace transform on the component of space yield ℒ𝑦𝑝=²⁻¹_𝑖=0ℎ_𝑖𝑝^{𝛼−𝑖−1}𝑝^𝛼−𝜇𝑝^𝛽−𝜏−𝜇¹⁻¹_𝑖=0ℎ_𝑖𝑝^{𝛽−𝑖−1}𝑝^𝛼−𝜇𝑝^𝛽.−𝜏 (29) For 𝑝∈ℂ and |𝜏𝑝^−𝛽/(𝑝^𝛼−𝛽−𝜇)|<1, in analogy of the discursion presented earlier for the case of Riemann-Liouville, we have the following: ℒ𝑦𝑝=²⁻¹_𝑖=0ℎ_𝑖^∞_𝑛=0𝜏^𝑛𝑝^{(𝛼−𝛽)−(𝛽𝑛+𝑖+1)}𝑝^𝛼−𝛽−𝜇^𝑛+1−𝜇¹⁻¹_𝑖=0𝜏^𝑛𝑝^{(𝛼−𝛽)−(𝛽𝑛+𝑖+1+𝛼−𝛽)}𝑝^𝛼−𝛽−𝜇^𝑛+1. (30) Hence for 𝑝∈ℂ and |𝜇𝑝^𝛽−𝛼|<1, we have that 𝑝^{(𝛼−𝛽)−(𝛽𝑛+𝑗+1)}𝑝^𝛼−𝛽−𝜇^𝑛+1=1ℒ𝑥𝑛!^{𝑛𝛼+𝑖}𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛽𝑛+𝑖+1}𝜇𝑥^𝛼−𝛽,𝑝^{(𝛼−𝛽)−(𝛽𝑛+𝑗+1+𝛼−𝛽)}𝑝^𝛼−𝛽−𝜇^𝑛+1=1ℒ𝑥𝑛!^{𝑛𝛼+𝑖+𝛼−𝛽}𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛽𝑛+𝑖+1+𝛼−𝛽}𝜇𝑥^𝛼−𝛽. (31) Thus from the above expression we derive the following solution to the space-time Caputo fractional derivative of hydrodynamic advection-dispersion equation (8): 𝑐(𝑥,𝑡)=²⁻¹_𝑖=0ℎ_𝑖𝑐_𝑖(𝑥,𝑡)−𝜇¹⁻¹_𝑖=0ℎ_𝑖𝑐_𝑖(𝑥,𝑡), (32) where for 𝑖=0𝑐_𝑖(𝑥,𝑡)=^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^{𝑛𝛼+𝑖}×𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛽𝑛+𝑖+1}𝜇𝑥^𝛼−𝛽−𝜇^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^{𝑛𝛼+𝑖+𝛼−𝛽}×𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛽𝑛+𝑖+1+𝛼−𝛽}𝜇𝑥^𝛼−𝛽 (33) and for 𝑖=1𝑐_𝑖(𝑥,𝑡)=^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^{𝑛𝛼+𝑖}×𝜕𝜕𝑥^𝑛𝐸_{𝛼−𝛽,𝛽𝑛+𝑖+1}𝜇𝑥^𝛼−𝛽. (34) And the coefficients ℎ_𝑖 are found by applying the initial and boundary conditions on 𝑐(𝑥,𝑡).

Example 2. Our concern here is to consider (8) when 𝛼=2 and 0<𝛽≤1. Following the discursion presented earlier, the analytical solution of space-time fractional derivative of hydrodynamic advection-dispersion equation has its two solutions given by 𝑐₁(𝑥,𝑡)=^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^2𝑛×𝜕𝜕𝑥^𝑛𝐸_{2−𝛽,𝛽𝑛+1}𝜇𝑥^2−𝛽−𝜇^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^{2𝑛+2−𝛽}×𝜕𝜕𝑥^𝑛𝐸_{2−𝛽,𝛽𝑛+3−𝛽}𝜇𝑥^2−𝛽,𝑐₂(𝑥,𝑡)=^∞_𝑛=0𝑅𝐷^𝑛exp[−𝜆𝑡]𝑡^−1−𝑛Γ()𝑥−𝑛𝑛!^2𝑛+1×𝜕𝜕𝑥^𝑛𝐸_{2−𝛽,𝛽𝑛+2}𝜇𝑥^2−𝛽. (35) These solutions are linearly independent and they provide the fundamental system of solutions to space-time Caputo fractional derivative of hydrodynamic advection-dispersion equation. An approximation of this series is given below for possible simulation. Some other analytical methods and their recent development for solving nonlinear fractional partial differential equation can found in the work done by [15–26] and the excellent book for analytical and numerical methods.

5. Numerical Simulation

Up to this section we expressed the solution of the fractional advection-dispersion equation in terms of Mittag-Leffler function. This function is cumbersome to be used in real world problem, especially when the users of this solution are from the field of geohydrology. Since the solution is in series form, one will need first to know how many terms of the series expansion can be used to simulate real world problem. Therefore to accommodate the users of this solution, we propose the approximate solution of the fractional advection-dispersion equation to be in the form of 𝑐()=𝑐𝑥,𝑡⁰exp−𝛾𝑡2×𝑥exp^𝛼𝑞_𝑟−𝑢_𝑟2𝐷_𝑟𝑥erfc^𝛼−𝑢_𝑟𝑡2𝐷_𝑟𝑡^1/𝛽𝑥+exp^𝛼𝑞_𝑟+𝑢_𝑟2𝐷_𝑟𝑥erfc^𝛼+𝑢_𝑟𝑡2𝐷_𝑟𝑡^1/𝛽. (36) The above solutions take into account the effect of the fractional derivative order. Now notice that if one set 𝛼=1 and 𝛽=2, we recover the solution of the advection-dispersion equation 𝑐()=𝑐𝑥,𝑡⁰exp−𝛾𝑡2×𝑥𝑞exp^𝑟−𝑢_𝑟2𝐷_𝑟erfc𝑥−𝑢_𝑟𝑡2𝐷_𝑟𝑡𝑥𝑞+exp^𝑟+𝑢_𝑟2𝐷_𝑟erfc𝑥+𝑢_𝑟𝑡2𝐷_𝑟𝑡. (37) To access the effect of the fractional order derivative into the solution of the advection-dispersion equation, we compare both solutions (36) and (37) with the theoretical values firstly and secondly we compare both solutions with experimental data obtained from one of the experimental sites of the Institute for Groundwater Studies (IGS). We shall start with the theoretical values.

Figures 1 and 2 show the numerical simulation of the plume first as a function of time and second as a function of time and space with the FADE and ADE for theoretical values. Figures 1 and 2 show the comparison of the approximate and exact solutions of FADE and ADE, respectively, as function of time for a fixed distance 𝑥. These figures are plotted via Mathematica.

Figure 1 

Comparison of FADE and ADE for 𝛼=0.95 and 𝛽=1.9.

Figure 2 

Comparison of FADE and ADE for 𝛼=1 and 𝛽=1.98.

Figures 3–7 show the density plots of the theoretical simulation of the plume by the FADE as function of time and fractional order derivative. The figures are simulated via MATLAB. Here, Figure 3 is the simulation of the concentration for FADE through the geological formation, for 𝑐₀=100,  𝛼=0.45,  𝛽=2,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1. Figure 4 is the simulation of FADE for 𝑐₀=100,  𝛼=1,  𝛽=2,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1. Figure 5 is the simulation of FADE for 𝑐₀=100,  𝛼=0.55,  𝛽=1.55,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1. And finally Figure 6 is the simulation of FADE for 𝑐₀=100,  𝛼=0.25,  𝛽=1.55,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1. Figure 7 is the simulation of FADE for 𝑐₀=100,  𝛼=0.25,  𝛽=1.95,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1.

Figure 3 

Simulation of the FADE (𝑐0=100,  𝛼=0.45,  𝛽=2,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1).

Figure 4 

Simulation of the ADE (𝑐0=100,  𝛼=1,  𝛽=2,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1).

Figure 5 

Simulation of the FADE  (𝑐0=100,  𝛼=0.55,  𝛽=1.55,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1).

Figure 6 

Simulation of the FADE (𝑐0=100,  𝛼=0.25,  𝛽=1.55,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1).

Figure 7 

Simulation of the FADE (𝑐0=100,  𝛼=0.55,  𝛽=1.95,  𝐷=2;  𝑞=1,  𝛾=0.25, and 𝜆=1).

From Figures 3–7 one can see that the solutions of FADE are not only a function of time and space but also a function of the order of the derivative. If these orders are integer, we recover the standard ADE. Figures 3 and 5 show that the order of the derivative can be used to simulate the real-world problem and this makes the fractional version of ADE better than the ADE.

To test the accuracy and efficiency of FADE, we compare the solution of FADE, ADE, and the experimental data from field observation. Figures 8 and 9 show the comparison between FADE, ADE, and measured data for different values of 𝛼 and 𝜇.

Figure 8 

Comparison of FADE, ADE, and experimental data from real world; 𝐷𝑟=4.5,  𝛽=1.95,  𝛼=0.99, and 𝑞𝑟=0.51.

Figure 9 

Comparison of FADE, ADE, and experimental data from real world; 𝐷𝑟=2.5, 𝛽=1.36, 𝛼=0.3,  𝑞𝑟=0.4, and 𝑐0=150.

The numerical simulation in Figures 8 and 9 lead us to believe that the order of the derivative plays an important role while simulation of the plume of the pollution in the aquifer. The comparison revealed that the fractional advection-dispersion equation is compatible with observations of the plume in the laboratory and the field. It predicts power law, faster than the apparent plume variance. It is shown that the traditional 2nd-order advection equation does not adequately describe the movement of solute trace in the aquifer. On the basis of this assertion we conclude that the fractional advection-dispersion equation is better than the classical version of advection-dispersion equation. Plumes observed in natural systems are used in this section to distinguish further the time and space nonlocalities. Realistic data from natural systems provide the most important criteria for distinguishing the space- and time-nonlocal processes and evaluating the applicability of the FADE models. Analysis, comparison, and application of various FADEs in this study are intended to provide a general guidance for model selection. Natural geological deposits with highly contrasting permeability may form mobile and relatively immobile zones, where the potential mass exchange between mobile and immobile zones results in a wide time distribution for solute “trapping”. The transport process, combined with the distinct particle status, can be characterized efficiently by the time-nonlocal model, including the time FADE. If the high-permeable material tends to form preferential flow paths, such as the interconnected paleochannels observed in alluvial depositional systems, then the solute transport may show a heavy leading edge, which can be described by the space FADE with maximally positive skewness as shown in Figures 3, 5, 6, 7, and 8. Development of partial differential equations such as the advection-dispersion equation (ADE) begins with assumptions about the random behavior of a single particle: possible velocities it may experience in a flow field and the length of time it may be immobilized. When assumptions underlying the ADE are relaxed, a fractional ADE (FADE) can arise, with a noninteger-order derivative on time or space terms. Fractional ADEs are nonlocal; they describe transport affected by hydraulic conditions at a distance. Space fractional ADEs arise when velocity variations are heavy tailed and describe particle motion that accounts for variation in the flow field over the entire system. Time fractional ADEs arise as a result of power law particle residence time distributions and describe particle motion with memory in time. As shown, the best fitting curve from the classical radial flow model greatly underestimated early arrival. We also remark that this solute flow model and its numerical solution match the test data closely only up to the peak.

An excellent literature review revealed that the fractional advection-dispersion equation has proven to be useful in modeling contaminant flow in heterogeneous porous media [26–34]. The fractional advection-dispersion equation is known to be a special case of a general transport equation with convolution flux [28] and a limit case of the continuous time random walk with power-law particle jumps [26–34]. It is a simple matter to derive the fractional advection-dispersion equation from the fractional conservation of mass equation using a moving coordinate system at the plume centre of mass, in exactly the same way that the usual advection-dispersion equation follows from the traditional conservation of mass equation. This approach validates the utility of the fractional advection-dispersion equation and related theories, by highlighting the scaling factor that renders the fractional equation scale invariant. We believe that this scaling captures the fractal nature of the porous medium.

6. Conclusion

The classical hydrodynamic advection-dispersion equation has been generalized using the concept of fractional order derivatives. This leads to the formulation of a new (generalized) form of the hydrodynamic advection-dispersion equation. A general solution of the new equation was given in terms of Mittag-Leffler functions for two general cases including Riemann-Liouville fractional derivative and the Caputo fractional derivative. The solutions of FADE are not only function of time and space but also a function of the order of the derivative. If these orders are integer, we recover the standard ADE. Figures 3–8 show that the order of the derivative can be used to simulate the real-world problem and this makes the fractional version of ADE better than the ADE. The comparison of FADE, ADE, and experimental data shows that the FADE is better in prediction than ADE.