Soil Science Society of America Journal 67:1079-1084 (2003)
© 2003 Soil Science Society of America
DIVISION S-1—SOIL PHYSICS
Application of the Fractional Advection-Dispersion Equation in Porous Media
Liuzong Zhou and
H. M. Selim*
Agronomy Dep., Sturgis Hall, Louisiana State Univ., Baton Rouge, LA 70803
* Corresponding author (mselim{at}agctr.lsu.edu)
 |
ABSTRACT
|
|---|
According to the classical advection-dispersion equation (ADE), the variance of travel distance increases linearly with the time elapsed after the release of a solute tracer. However, several field studies showed nonlinear relationships exist between the variance of travel distance for solute tracers and time. Therefore, the transport at the field scale often shows non-Gaussion property. To describe the non-Gaussian transport process, a fractional advection-dispersion equation (FADE) is a possibility. Fractional ADE has been applied to the transport process of non-reactive solutes in a laboratory sandbox as well as that for the Cape Cod field site. In this paper, we discuss several issues regarding the application of FADE. Because of the complexity of the FADE, its application often needs justification. We propose an F-test to judge whether a FADE is necessary. For the Cape Cod experiment, the FADE of order 1.82 is necessary according to our F-test results. In addition, we present a method to estimate the dispersion coefficient and fractional order simultaneously. The obtained parameters from joint estimation will recover the variance-time relationship. Finally, we pointed out that the dispersivity-time or dispersivity-mean travel distance relationship is needed to estimate the fractional order 
. Use of dispersivity-distance relationship for estimating the fractional order is theoretically unjustified and may result in unreasonable values.
Abbreviations: ADE, advection-dispersion equation • FADE, fractional advection-dispersion equation
|
INTRODUCTION
|
|---|
THE CLASSICAL ADE is used to describe Fickian diffusion. The characteristic of Fickian diffusion is that the variance of travel distance increases linearly with time. However, with increasing evidence from field experiments, this relationship may not always hold true. For example, data from the Borden site (Freyberg, 1986) and the Columbus Air Force Base site (Adams and Gelhar, 1992) displayed nonlinear (power-law) increase of the variance of travel distance with time. Obviously, the traditional ADE associated with Fickian diffusion is no longer applicable to the anomalous diffusion in heterogeneous media. Accordingly, the normal distribution was replaced by Lévy distribution to describe anomalous transport phenomenon. Chaves (1998) proposed a fractional diffusion equation to describe Lévy flights. Schumer et al. (2001) also showed that the FADE could be obtained by replacing the classic Fick's law in a Eulerian evaluation of solute transport in porous media with a fractional Fick's law. Benson et al. (2000a) developed a FADE to describe anomalous transport phenomena in aquifers. Cushman and Ginn (2000) showed that the FADE is a special case of the convolution-Fickian non-local ADE proposed earlier by Cushman and Ginn (1993).
The one-dimensional symmetrical FADE reads,
 | [1] |
Where D is the fractional dispersion coefficient (L T-1) and the superscript is the order of fractional differentiation, 0 < ≤ 2. It is worth noting that the classic ADE is recovered when the fractional order of FADE is set to 2. Fractional derivatives are integro-differential operators defined as (Podlubny, 1999),
 | [2] |
 | [3] |
Where > 0, 
is the gamma function, and k is the smallest integer number larger than . If the value of is a whole number, fractional derivatives reduce to ordinary derivatives. Some properties of fractional derivative are given in the appendix of Benson et al. (2000a).
Benson et al. (2000b) applied the FADE to two cases: a laboratory sandbox transport experiment and the Cape Cod field scale transport experiment. The fractional order for the laboratory data was found to be 1.55 whereas the Cape Cod bromide plumes could be modeled using a FADE with an order of 1.65 to 1.8. Pachepsky et al. (2000) simulated scale-dependent solute transport in soils using the FADE. They also presented a comparison between FADE and ADE based on statistical analysis. They found that FADE could describe transport processes of a solute tracer better. Benson et al. (2000b) also provided a method to estimate the fractional order and the fractional dispersion coefficient D separately. However, in our opinion, a simultaneous estimation of both and D is more appropriate. Specific reasons will be discussed later in this paper. Benson et al. (2000b) estimated the fractional order based on the relationship between measured apparent dispersivity and distance from the sandbox experiment. Pachepsky et al. (2000) justified the application of FADE based on the scale-dependent transport phenomenon. We believe that the order of FADE cannot be associated with scale-dependent transport directly. In fact, the relationship between and scale-dependent dispersivity is not theoretically supported. All the above issues about the application of FADE need to be clarified. Therefore, the objectives of this paper are:
- (i) to propose a method for simultaneously or jointly estimating the fractional order and dispersion coefficient,
- (ii) to propose a more appropriate method to test the applicability of FADE, and
- (iii) to point out a pitfall in the estimation of the fractional order parameter .
|
THEORY
|
|---|
Comparison between ADE and FADE
As indicated above, the classical ADE predicts a linear increase of the variance of travel distance with time or mean travel distance whereas FADE predicts a nonlinear increase of the variance of travel distance. Therefore, comparing linear models of variance pattern with nonlinear ones is equivalent to direct comparison between ADE and FADE. The comparison between these two models is achieved by conducting the appropriate F-test (Green and Caroll, 1978).
The relationship of the measured variance of travel distance to the mean travel distance or time (i.e., 
2x vs 
or 2x vs t) can thus be fitted using linear or nonlinear models through linear and nonlinear least square fitting techniques. Because of the simplicity of the linear least square method, nonlinear models are often converted to linear format by taking logarithm transformation. Therefore, there are three options to model 2x - or 2x - t relationship:
- Linear coordinate: Linear model (Model I):
 | [4] |
- Linear coordinate: Nonlinear model (Model II) (nonlinear least square fitting):
 | [5] |
- Log-log coordinate: Nonlinear model (Model III)(linear least square fitting):
 | [6] |
where Y refers to the variance of travel distance 2x, X stands for the mean travel distance or time t, and A, B, 10a, and b are regression coefficients that can be used to determine the transport parameters. Models II and III are equivalent. For a linear model, the regression coefficient A is twice the longitudinal dispersivity or dispersion coefficient depending on whether X refers to the mean travel distance or time. For nonlinear models, the power B or b can be used to determine the fractional order .
The F statistics for the comparison of different models was used by Green and Caroll (1978) and is defined as,
 | [7] |
where SSEf and SSEr are, respectively, the sums of squared errors of a full model, and a "restricted" model; dfef and dfer are the respective degrees of freedom of errors. Here the full model refers to the FADE while the restricted model refers to the classical ADE. Accordingly, for the variance models, the nonlinear model is the full model whereas the linear one is the restricted model. The SSE is calculated as,
 | [8] |
where 2i,f and 2i,o are fitted and observed variance of travel distance 2x, j = r (ADE, linear) or f (FADE, nonlinear), and n is the total number of data points. If the calculated F-statistics does not exceed the critical Fv1, v2 with degrees of freedom v1 = dfer - dfef, v2 = dfef at 0.05 significance level, the two models are not significantly different in describing the variance pattern of the transport process. For the opposite case, the full model, that is, FADE is superior to the classical ADE. It should be pointed out that the FADE and ADE with the parameters given by the regression coefficients of the fitted variance-time relationship actually mimic the variance pattern of the observed data. The direct comparison of FADE and ADE could also be achieved using Eq. [7] and [8]. In this case, the variance of travel distance must be replaced by solute concentration from breakthrough curves or snapshots.
Simultaneous Estimation of and D
Because the ADE predicts a linear increase of the variance of travel distance 2x, the measured relationship between 2x and time t or mean travel distance is often used to estimate the dispersion coefficient or dispersivity. Similarly, this method is also applicable to the parameter estimation of FADE. According to Benson et al. (2000b), the variance of travel distance and the scale of Lévy distribution have roughly the following relationship:
 | [9] |
For = 2, both sides of Eq. [9] are exactly equal for the classical ADE. By comparison of Eq. [9] with Eq. [5] and [6], one can find a relationship between the regression coefficients and fractional parameters. Based on Model II and the relationship between 2x and t, we have,
 | [10] |
 | [11] |
Simple algebraic calculation gives the estimations of parameters and D as,
 | [12] |
 | [13] |
For Model III and the relationship between 2x and t, similar equations can be obtained as,
 | [14] |
 | [15] |
In addition, if one assumes that the molecular diffusion can be ignored, the fractional dispersion coefficient can be expressed as,
 | [16] |
where v is the pore water velocity (L T-1) and 
is the fractional dispersivity (L - 1). Accordingly, Eq. [9] becomes,
Where is mean travel distance. In this case, we may also estimate according to the measured relationship between 2x and . The fractional dispersivity is thus given by,
Model II:
 | [17] |
Model III:
 | [18] |
The estimations of D or for the classical ADE will be recovered if B or b is reduced to 1. It is worth noting that the parameters in Eq. [12] and [13] or [14] and [15] are simultaneous estimation. Use of the parameters estimated jointly would more accurately reproduce the variance pattern.
|
RESULTS AND DISCUSSION
|
|---|
Justification for Use of FADE
The goal of a comparison between different models is to determine whether a FADE is necessary and whether the FADE is superior to the classical ADE in description of the given data set. In this study, the data from the Cape Cod site are used as an example. We fitted the above three models with the data from Garabedian et al. (1991). The functional relationships between 2x versus as well as 2x versus t were obtained. The autofluorescence data set on Day 349 (Garabedian et al., 1991) was not included in either variance-mean travel distance fitting or variance-time fitting. Comparisons among fitted curves and measured data are shown in Fig. 1 and 2
. The optimized regression coefficients are summarized in Tables 1 and 2
. Double asterisk indicates that the estimated is significantly <2 at the 0.05 significance level. If the estimated is not significantly different from 2, it is not necessary to use a fractional order other than 2. A quick way to test the significance of is to check whether = 2 is included in the 95% confidence interval of the estimated . According to Eq. [12] and [14], an significantly <2 is equivalent to an exponent B or b significantly >1. The obtained fractional order values based on the relationship between variance and mean travel distances are 1.82 and 1.65 for Model II and Model III, respectively. Both values are significantly <2.0. However, if we estimate the fractional order based on the variance-time relationship, we obtained different values of , that is, 1.91 and 1.79 for Model II and Model III, respectively. Between the two values, only the second ( = 1.79 based on Model III) is significantly <2.0. These values for the fractional order are almost the same as those reported in Benson et al. (2000b). Therefore, different estimated values for were obtained when the same data set was fitted using the two methods, that is, direct nonlinear fitting of the original data versus linear fitting of log-transformed data. The linear fitting of log-transformed data overestimates the slope (Model II versus Model III in Tables 1 and 2), which resulted in an underestimated value of . A smaller based on Model III is thus an artifact of the fitting technique. After we assess the sum of squared error (SSE) based on back-calculated parameters from Model III, we obtain values even greater than those for a linear model (last column in Tables 1 and 2). In other words, a FADE with parameters estimated from Model III is inferior to an ADE with parameters estimated from Model I with regard to the description of the variance pattern of the Cape Cod data. From this point of view, use of FADE with parameters from Model III is not recommended.
Results of F-test for the variance-mean travel distance relationship and the variance-time relationship are given in Tables 3 and 4
, respectively. The values in brackets are the critical F value; double asterisk indicates a significant difference at the 0.05 probability level. As shown in Tables 3 and 4, the calculated F values for comparison between Model I and III are negative. This fact implies that no extra sum of squares is gained through shifting from Model I to Model III. In other words, the classical ADE with parameters from Model I can give a better description of the variance pattern than the FADE with parameters from Model III. Therefore, the use of the fractional order derived from Model III, that is, of 1.65 and 1.79, is questionable. Similarly, the use of of 1.91 based on the variance-time relationship does not give an extra sum of squares significantly different from zero and thus is unjustified. Based on the variance-mean travel distance relationship, the F-test result showed that a FADE with of 1.82 does provide some extra sum of squares of error, however. Therefore, it is justified to use a FADE with an order of 1.82. Pachepsky et al. (2000) presented a different F-test to compare the performance of ADE and FADE. Their F statistics is computed as the ratio of the mean square error of the ADE to that of the FADE. According to our understanding, the objective of an F-test is to examine whether the mean square error in the numerator is significantly different from zero. When two related models are compared, the mean square error in the numerator of an F statistics should be a factor that indicates the difference of the two models. Therefore, the F statistics defined in Eq. [7] is more appropriate than a direct ratio.
The above conclusion is a weak one because both the variance of travel distance and the mean travel distance were estimated values based on moments analysis; none of these two variables were derived from direct measurements. Therefore, to strengthen the above conclusion, one needs to compare the predicted solute concentration profiles based on Models I and II with measured data. Thus the SSE in terms of solute concentration for each model must be calculated. Based on the SSE data, the F-statistics as discussed above should be computed. If the calculated F value is greater than the critical F, then a FADE is indeed necessary.
Jointly Estimated Parameters
The estimated fractional dispersivity and dispersion coefficient D, based on different the models, are given in Tables 1 and 2. The obtained dispersivity value for the classical ADE is 0.91 m according to our analysis. This value is very close to 0.96 m reported by Garabedian et al. (1991) based on moment analysis. Based on the measured variance-mean travel distance relationship, we obtained , doublets of 1.82, 1.43 m0.82 and 1.65, 1.10 m0.65. Using a pore-water velocity of 0.42 m d-1 (Garabedian et al., 1991), we obtained corresponding D of 0.25 m1.82 d-1 and 0.19 m1.65 d-1. Similarly, Based on the measured variance-time relationship, the estimated , D pairs are 1.91, 0.31 m1.91 d-1, and 1.79, 0.25 m1.79 d-1, respectively. The values of estimated fractional dispersion coefficient are comparable with the average values reported by Benson et al. (2000b). Parameters obtained by Benson et al. (2000b) were based on separate estimation. First, the fractional order was estimated. After that, the dispersion coefficient was calculated based on known . They showed that the estimated D changes slightly with time. Thus, this method to estimate the dispersion coefficient D is questionable. The estimation of was based on the time evolution of the variance of travel distance. In other words, the change of the variance of travel distance over a range of time shows that the transport is anomalous. The obtained is an optimal value to describe the entire pattern of the growth of the variance of travel distance over time. This value does not necessarily mean that the estimated is exactly the fractional order of the transport process at any specific time. One has to fit the snapshot at a time t to the Lévy distribution to find out the fractional order and related parameters for that specific time. Therefore, to calculate the fractional dispersion coefficient based on the estimated fractional order is problematic. On the other hand, each estimated D is based on one single data point. That means what one obtains is a point estimator. Because the total degrees of freedom are equal to the number of parameter to be estimated, one has no way to figure out the confidence interval of the estimated parameter. From statistical point of view, this type of estimation is unreliable. Another point we want to make is that whether the use of separately estimated and D would recover the variance pattern remains unknown. This is because the estimation of D is not related to the variance-time or the variance-mean travel distance relationship. Thus to mimic the variance pattern, one has to estimate both parameters and D simultaneously or jointly.
A Pitfall in the Estimation of
Several studies show that the following relationship holds true for anomalous diffusion (Weeks et al., 1995),
 | [19] |
where 2x(t) is the measured variance of travel distance at time t, and 
is a real number greater than zero. For a Fickian process, equals 1. However, for the case where > 1.0, the process is referred to as super-diffusion. Simple comparison between Eq. [19] and [9] gives a relationship between and as = 2/. Because the traditional ADE can only explain linear growth of 2x with time, the nonlinear increase can only be accounted for by using a time-dependent dispersion coefficient or dispersivity. The time-dependent longitudinal dispersivity L is obtained by taking the first derivative of 2x(t) with respect to time (Gelhar, 1993, p. 200–202.). With Eq. [19], we have the following power function of L,
 | [20] |
For a constant pore-water velocity 
, one obtains,
 | [21] |
where is the mean travel distance at time t, and = vt. Equations [20] and [21] describe a time-dependent dispersivity. Usually, one can obtain the apparent dispersivity 
L at different times by fitting analytical solutions to the solute concentration profiles. It is easy to show that the L increases at the same rate, that is, proportional to t - 1 or -1 if we assume the apparent dispersivity at time t is the same as the average dispersivity over the period of 0 to t. Based on the relationship between the apparent dispersivity and time or mean travel distance, one can estimate the fractional order . If the measured L is found to grow proportionally to tm or m, where m = - 1, then based on the above discussion, can be estimated as = 2/(1 + m). For a constant L over time (m = 0), the classical ADE is recovered, that is, = 2. Similar relation has been obtained by Benson et al. (2000b). The exponent m is simply the slope of the log-log plot of L versus t or . There is an important point that one must bear in mind when using this method to estimate the fractional order. That is, the relationship between apparent dispersivity and t or must be used.
Pachepsky et al. (2000) summarized several case studies in the literature regarding the relationship between apparent dispersivity and distance. They justified the use of FADE based on an apparent linear trend of dispersivity with distance in a log-log plot. However, the relationship between dispersivity and distance is not consistent in the literature. The apparent dispersivity may hold constant, increase linearly or nonlinearly, decrease linearly or nonlinearly, or randomly fluctuate with distance. A similar conclusion has been reached by Ellsworth et al. (1996). The slope m' of the log-log plot of dispersivity versus distance cannot be used to compute . We can use a simple counter example to support our argument. If the measured apparent dispersivity actually decreases with distance, one will obtain a negative slope m' in a log-log plot. As a result, an larger than 2 will be obtained. Obviously, such a value for is not reasonable. On the other hand, the slope m' for a log-log plot of dispersivity versus distance x cannot be associated with in the same way as that for a log-log plot of dispersivity versus t or . To verify our conclusion, we examined the data sets from Toride et al. (1995). Pachepsky et al. (2000) gave several optimized values for the data sets from Toride et al. (1995), for different soil depths and for both unsaturated and saturated soil column experiments. Based on the data in Table 2 of Pachepsky et al. (2000), we computed the apparent dispersivity L (commonly expressed as 
L) at different depths as,
 | [22] |
The obtained apparent dispersivity is then plotted against distance in our Fig. 3
. The fractional order based on m' is given in Table 5
together with results from Pachepsky et al. (2000). Figure 3 shows that the dispersivity increases with distance for the unsaturated soil column whereas it decreases with distance for the saturated soil column. The fractional order based on Fig. 3 is 1.538 for unsaturated soil condition. This value is not among the three optimized values given by Pachepsky et al. (2000). Though an order of 1.538 is very close to 1.574 for the length of 23 cm, one cannot thus conclude that a fractional order could be obtained through the slope of log-log graph of dispersivity versus distance. For saturated soil conditions, we obtained a value of 3.75 which is quite different from those obtained through optimization. The reason why computed according to the dispersivity-distance relationship is unreasonable is because mean travel distance and distance x are two different concepts and they are not interchangeable. The fundamental difference between and x lies in that x is an independent variable whereas is a dependent variable which depends linearly on time. On the other hand, the cornerstone to calculate based on dispersivity at different t or is the nonlinear growth of 2x. The estimation of 2x can only be achieved by conducting moment analysis on solute concentration distribution at a certain time. It is not possible to back-calculate and obtain the growth of 2x with time based on estimated dispersivity at different distances.
|
SUMMARY AND CONCLUSIONS
|
|---|
We discussed several issues related to the application of a FADE in porous media. Use of the FADE often needs justification because its solution is more complex than that of the classical ADE. If the order of a FADE is not significantly smaller than 2 or no improvement in describing the transport process is achieved, the use of FADE is unjustified. We proposed an F-test on the basis of extra sum of squared errors to examine whether a FADE is superior to the classical ADE. For the Cape Cod field experiment, because the variance-mean travel distance relationship shows significant non-linearity, the FADE of order 1.82 is necessary according to our F-test results. Benson et al. (2000b) presented a method to estimate and D separately. However, we found that separate estimation of both parameters is problematic. One major problem lies in that it is difficult to ascertain that the use of two separately estimated parameters would actually mimic the variance pattern. To overcome this, we presented a method to estimate both parameters simultaneously. Theoretically, the obtained parameters from joint estimation would more accurately recover the variance-time relationship. Finally, we pointed out that the dispersivity-time or dispersivity-mean travel distance relationship is needed to estimate the fractional order . The use of dispersivity-distance relationship to estimate the fractional order is theoretically unjustified and may result in unreasonable values.
|
NOTES
|
|---|
Louisiana Agric. Exp. Stat. Manuscript no. 02-14-0773.
Received for publication April 29, 2002.
|
REFERENCES
|
|---|
- Adams, E.E., and L.W. Gelhar. 1992. Field study of dispersion in a heterogeneous aquifer 2. Spatial moments analysis. Water Resour. Res.
28:3293–3307.
- Benson, D.A., S.W. Wheatcraft, and M.M. Meerschaert. 2000a. The Fractional order governing equation of Levy motion. Water Resour. Res.
36:1413–1423.
- Benson, D.A., S.W. Wheatcraft, and M.M. Meerschaert. 2000b. Application of a fractional advection-dispersion equation. Water Resour. Res.
36:1403–1412.
- Chaves, A.S. 1998. A fractional diffusion equation to describe Levy flights. Phys. Lett. A
239:13–16.
- Cushman, J.H., and T.R. Ginn. 1993. Nonlocal dispersion in media with continuously evolving scales of heterogeneity. Transp. Porous Media
13:123–138.
- Cushman, J.H., and T.R. Ginn. 2000. Fractional advection-dispersion equation: A classical mass balance with convolution-Fickian flux. Water Resour. Res.
36:3763–3766.
- Ellsworth, T.R., P.J. Shouse, T.H. Skaggs, J.A. Jobes, and J. Fargerlund. 1996. Solute transport in unsaturated soil: Experimental design, parameter estimation, and model discrimination. Soil Sci. Soc. Am. J.
60:397–407.[Abstract/Free Full Text]
- Freyberg, D.L. 1986. A natural gradient experiment on solute transport in a sandy aquifer 2. Spatial moments and the advection and dispersion of nonreactive tracers. Water Resour. Res.
22:2031–2046.
- Garabedian, S.P., D.R. LeBlanc, L.W. Gelhar, and M.A. Celia, 1991. Large-scale natural gradient tracer test in sand and gravel, Cape Cod, Massachusetts 2. Analysis of spatial moments for a nonreactive tracer. Water Resour. Res.
27:911–924.
- Gelhar, L.W. 1993. Stochastic subsurface hydrology. Prentice Hall, Englewood Cliffs, NJ.
- Green, P.E., and J.D. Caroll. 1978. Analyzing multivariate data. John Wiley & Sons, New York.
- Pachepsky, Y., D. Benson, and W. Rawls. 2000. Simulating scale-dependent solute transport in soils with the fractional advective-dispersive equation. Soil Sci. Soc. Am. J.
64:1234–1243.[Abstract/Free Full Text]
- Podlubny, I. 1999. Fractional differential equations. Academic Press, San Diego, CA.
- Schumer, R., D.A. Benson, M.M. Meerschaert, and S.W. Wheatcraft. 2001. Eulerian derivation of the fractional advection-dispersion equation. J. Contam. Hydrol.
48(1–2):69–88.[Medline]
- Toride, N., F. Leij, and M.Th. van Genuchten. 1995. The CXTFIT code for estimating transport parameters from laboratory or field tracer experiments. Version 2.0. Research Rep. 137. U.S. Salinity Lab., Riverside, CA.
- Weeks. E.R., T.H. Solomon, J.S. Urbach, and H.L. Swinney, 1995. Observation of anomalous diffusion and Levy flights. p. 51–71. In M.F. Shlesinger et al. (ed.) Levy flight and related topics in Physics. Springer, New York.
This article has been cited by other articles:
|
|
|
|
 
A. Cortis and B. Berkowitz
Anomalous Transport in "Classical" Soil and Sand Columns
Soil Sci. Soc. Am. J.,
September 1, 2004;
68(5):
1539 - 1548.
[Abstract]
[Full Text]
[PDF]
|
|
|
TOP
ABSTRACT
INTRODUCTION
