The conventional advection-dispersion equation (ADE) has been widely used to describe the solute transport in porous media. However, it cannot interpret the phenomena of the early arrival and long tailing in breakthrough curves (BTCs). In this study, we aim to experimentally investigate the behaviors of the solute transport in both homogeneous and heterogeneous porous media. The linear-asymptotic model (LAF solution) with scale-dependent dispersivity was used to fit the BTCs, which was compared with the results of the ADE model and the conventional truncated power-law (TPL) model. Results indicate that (1) the LAF model with linear scale-dependent dispersivity could better capture the evolution of BTCs than the ADE model; (2) dispersivity initially increases linearly with the travel distance and is stable at some limited value over a large distance, and a threshold value of the travel distance is provided to reflect the constant dispersivity; and (3) compared with the TPL model, both the LAF and ADE models can capture the behavior of solute transport as a whole. For fitting the early arrival, the LAF model is less than the TPL; however, the LAF model is more concise in mathematics and its application will be studied in the future.

Accurately predicting the solute transport in porous or fractured media is essential for groundwater quality control and pollution remediation. The advection-dispersion equation (ADE) with constant dispersivity, typically adopted to model and explain the solute transport process in groundwater system, does not satisfactorily describe the underground solute transport in many cases [13]. A growing number of laboratory and field experiments have shown that ADE could not capture the early arrival and long tailing of breakthrough curves (BTCs) [4, 5]. Dispersivity is widely thought to be related to scale dependence, and the scale-dependent dispersivity is largely due to the heterogeneity of porous media at different scales [68].

In order to predict solute transport in groundwater system and solve the problem of the scale-dependent dispersivity, many numerical models were developed. Previous studies showed that distance-dependent dispersivity can be described by five types of functions, namely, linear, linear-asymptotic, parabolic, exponential, and hyperbolic functions. Of these, the linear distance-dependent dispersivity is often obtained from analytical or semianalytical solutions for solute transport in infinite or semi-infinite domains. Pickens and Grisak [9] used a finite element (linear triangular) solute transport model to handle scale-dependent dispersion, which provided a continuous description of the dispersive process over all travel distances and times. The model was also successfully applied to previous tracer test results and exhibited a scale-dependent dispersion in heterogeneous porous media. Moradi and Mehdinejadiani [10] compared the fractional dispersion coefficient (Df) in the space fractional advection-dispersion equation (s-FADE) with the dispersion coefficient (D) in the advection-dispersion equation (ADE) by six solute transport experiments at three flow rates, and they concluded that the s-FADE model could efficiently explain the solute transport process in natural porous media because the scale dependency of the s-FADE model is less than that of the ADE model. Many studies confirmed that analytical solutions of various models have one-dimensional (1D) scale dependence dispersivity. Yates [11, 12] developed analytical solutions for dissolved substances in heterogeneous porous media using 1D generalized ADE with linear or exponential distance-dependent dispersion coefficients and found a distance-dependent dispersion relationship. Huang et al. [13] obtained an analytical solution for 1D solute transport with asymptotic scale-dependent dispersion by assuming that dispersivity linearly increased with transport distance within a certain distance and subsequently reached an asymptotic value. Pang and Hunt [14] also verified the applicability of the 1D scale-dependent dispersion model (SDM) using an 8 m long soil column experiment. Guerrero and Skaggs [15] developed a general 1D analytical solution for solute transport in a finite domain with the scale-dependent dispersivity, and they presented the solution in detail for the linear distance-dependent dispersivity. Besides, a mobile-immobile model was also proposed for the solute transport with scale-dependent dispersion in a heterogeneous porous medium, and a semianalytical solution with linear and exponential distance-dependent dispersivity was obtained [16]. The mathematical model considering the nonlinear change of the mass transfer coefficient can better fit the case of solute transport in porous media [17]. You and Zhan [18] proposed two new semianalytical solutions, namely, the linear-asymptotic (LAF solution) and the exponential distance-dependent (EF solution) dispersivities for solute transport in a finite column, and they found that breakthrough curves change faster and possess higher peak values than previous studies, and the differences increased with the dispersed growth rate and the asymptotic dispersivity. The LAF model, in which the dispersion is treated as a linear function to quantify the non-Fickian migration in a porous glass column, had a higher simulating accuracy than the ADE model, and the maximum error was only 0.62 g/L [19]. Swami et al. [20] found that the multiprocess nonequilibrium transport model better fitted the observed data of a reactive solute transport than the linear and exponential distance-dependent dispersion model in stratified porous media. The effects of distance-dependent dispersivity and tailing on the observed BTCs are visible.

From the above results of the analytical solution, the linear-asymptotic (LAF solution) model is more suitable for the experimental BTCs of solute transport in theory and can describe the solute transport accurately and interpret the distance-dependent dispersivities under certain boundary conditions in the subsurface, but it has not been verified by the experimental results.

In this study, we aim to verify the LAF model using the experimental results for 1D solute transport, assuming that the dispersivity is a nonlinear function of the migration distance, with the scale-dependent dispersivity also taken into account. The solute transport process in 1D homogeneous and heterogeneous porous media with a 600 cm long glass chamber was performed and simulated by the LAF and ADE models, and the fitting results were compared. The fitting of the LAF and TPL models which consider the variable dispersivities are also compared.

The ADE model for one-dimensional solute transport with the steady-state flow can be described as follows:
(1)RCt=xDCxVCx,
where c is the solute concentration (ML-3), t is time (T), D is the dispersion coefficient (L2T-1), V is the flow velocity (LT-1), x is the distance from the inlet boundary (L), R is the retardation factor.
Considering that D is spatially dependent rather than constant, the LAF model can be used for a solute transport in a one-dimensional finite domain [18], and the dispersivity could be written as follows:
(2)αx=λx,0xL0,α0,L0<xL,
where αx is the dispersivity (L), λ is the growing rate of the dispersivity near the source, α0 is the asymptotic dispersivity (L), L0 is the thresholds of distance where the dispersivity reaches the asymptotic value (L). This model implies that the dispersivity increases linearly with the transport distance in the region near the inlet boundary, and it becomes constant in the region far away from a large distance.
The piecewise function of equation (2) consists of two segments; i.e., the dispersivity increases linearly with distance as x is between 0 and L0, and when xis over L0, the dispersivity is constant. Therefore, the governing equations for solute transport in steady-state flow can be described as follows:
(3)RC1t=xD1xC1xVC1x,0xL0(4)RC2t=D22C2x2VC2x,L0<xL,
where C1 and C2 are the solute resident concentrations (ML-3) in the two segments, respectively; R is the retardation factor; D1 and D2 are the dispersion coefficients (L2T-1) in the two segments, respectively. According to equation (2), D1 and D2 could be written as follows:
(5)D1x=λxV+D0,(6)D2=α0V+D0,
where D0 is the molecular diffusion coefficient (L2T-1), D1x is the scale-dependent dispersion coefficient (L2T-1). D0 is generally neglected because it is much smaller than the mechanical dispersion coefficient in the field experiments [14]. However, it is not negligible for small-scale laboratory experiments [18]. Moreover, scale-dependent dispersivity increases from zero at the source and is very small at the distance near the source, and in this sense, molecular diffusion should play an important role [20]. In this study, molecular diffusion will be considered in the model.
The LAF solution in the Laplace domain is as follows:
(7)C¯1D=AεrKr2rsε+BεrIr2rsε,(8)C¯2D=Eea1xD+Fea2xD,r=1λ,ε=λxD+D0D,a1=1+4λβ+D0Ds2λβ+D0D,a2=11+4λβ+D0Ds2λβ+D0D,
where s is the dimensionless time of the dimensionless Laplace transform parameter; Ir and Kr are the rth-order modified Bessel function of the first kind and second kind, respectively. The variables in equation (7) and equation (8) can be calculated in equation (9) and equation (10).
The analytical solution for the one-dimensional ADE model with constant dispersivity is
(9)Cx,t=c02erfcxuxt2DLt+expuxxDLerfcx+uxt2DLt,
where erfc is the complementary error function, exp is the exponent, DL is the dispersion coefficient, ux is the flow velocity.
The governing equations of TPL in Laplace space are as follows:
(10)uc~x,uc0x=M~uvψxc~x,uDψ2x2c~x,u,
where
(11)M~utcharuψ~u1ψ~u,
where C0 is the initial concentration, u is the Laplace transform parameter, C¯x,u is the normalized concentration, M¯u is a memory function, ψ¯u is the Laplace transform form of ψt, vψ is the solute particle transport speed, Dψ is the dispersion coefficient which can be used to describe local dispersion of the solute transport. vψ and Dψ depend on ψt, which are different from the dispersion coefficient D and average flow velocity v in the classical ADE model.

3.1. Experimental Apparatus

The schematic setup for solute transport experiments in porous media is presented in Figure 1. The experimental apparatus consists of water flumes, a model body, and piezometer tubes. The water flumes include a recharge flume and a discharge flume, which could create a stable gradient of hydraulic head. The model body is a thin glass chamber with 600 cm length, 10 cm width, and 2 cm height and filled with different types of sands as the porous media. Piezometer tubes are separately connected with the recharge and discharge flume to measure the inlet and outlet water heads under different flow rates.

3.2. Experimental Materials

White sand, quartz sand, and glass beads are used as the porous medium materials and filled into the glass chamber, respectively. The particle size of white sand is in the range of 0.5 to 1.0 mm, the particle size of quartz sand is in the range of 1.0 to 2.0 mm, and the diameter of glass beads is 3.0 mm. Brilliant Blue (BBF) solution is used as the conservative tracer at a concentration of 0.5 g/L.

3.3. Experimental Design

The experiment was carried out in three groups according to the type of sand filled or flow direction. The chamber was filled with quartz sand in Experiment One. In Experiment Two, the chamber was filled with white sand, quartz sand, and glass beads, which were loaded into the chamber in sequence with the length ratio of 1 : 1 : 1. In Experiment Three, the chamber was filled with the same types of sand as in Experiment Two, but the water flow was in the opposite direction. Photos were taken and analyzed at 200, 400, and 600 cm along the flow direction under the flow rate of 0.7, 1.0, and 1.5 mL/s, respectively.

3.4. Experimental Procedure

The experimental procedure is composed of three major steps as follows:

  • (1)

    Step 1: the experimental apparatus is assembled and the glass chamber is filled with the specified types of sand and packed homogeneously, and then the chamber was filled with water to check the seal

  • (2)

    Step 2: solute transport experiments were performed. The step contains the following substeps:

    • (a)

      The hydraulic gradient is adjusted to obtain the predesigned flow rate which was determined by a cylinder and a chronograph. Timing is starded when BBF is observed flowing into the model entry

    • (b)

      Photos were taken by a camera, and the time was recorded during the experiment at three sampling locations. The experiment ended when the water in the device was completely replaced by the BBF solution; that is, the BBF concentration increased to 0.5 g/L at the outlet

    • (c)

      After the experiment, the solution in the chamber was pumped out, and the chamber was washed thoroughly using tap water

    • (d)

      Changing the flow rate, the above substeps (1)–(3) are repeated

    • (e)

      Changing the sand type in the chamber, the above substeps (1)–(3) are repeated

    • (f)

      Changing the flow direction, the above substeps (1)–(3) are repeated

  • (3)

    Step 3: the G value of photos was transformed into solute concentration, and the breakthrough curves (BTCs) were obtained

Each experiment was performed in triplicate and averaged.

4.1. Experimental Analysis of Hydraulic Conditions

The relationships between the hydraulic gradient (J) and the flow velocity (v) in the three groups of experiments were obtained and shown in Figure 2. From Figure 2, the correlation coefficients R2 were greater than 0.99 in the three experiments, showing each had a good linear relationship. Solid lines denote the best-fit lines following Darcy’s law. Therefore, the water flow in the chambers follows Darcy’s law.

The corresponding value of G was obtained from the standard concentration curve by the software processing based on image analysis. By comparing the three standard concentration curves, the best linear correlation is selected. So a mathematical model with a G value and concentration is established, which could describe the relationship between the G value and the concentration of the image accurately and quickly.

As shown in Figure 3, the value of G and concentration standard curves of the three media were established in a specified color space, and the G level and concentration relationships were linear.

4.2. Comparison of Results of Numerical Simulation of Solute Transport

Observed BTCs and model results simulated by the ADE and LAF models at different flow rates and different scales for Experiment One are compared in Figure 4. Fitting parameters of the LAF model and the ADE model for Experiment One are listed in Table 1. From Figure 4, the LAF model fitted the observed concentrations well. The precision of the two models were compared by the correlation coefficient (r2) and root mean square error (RMSE). The values of r2 and RMSE of the ADE and LAF models for Experiment One are shown in Table 1. The values of r2 were over 0.957 at different scales no matter how the flow velocity changes, and the corresponding RMSE was below 0.043. The LAF model is much better than the ADE model over a range of L values. The prediction value of the ADE model was always lower than that of LAF and was observed especially at lower flow velocity. The difference among the three BTCs becomes smaller as the velocity increases at the same scale. At the same flow rate, the BTCs observed and simulated by the LAF and ADE models in the homogeneous quartz chamber almost satisfied a linear-asymptote trend with time at distances of 200, 400, and 600 cm.

Figures 5 and 6 are BTCs obtained in heterogeneous porous media with the opposite flow direction for Experiments Two and Three. From Figure 5, comparing the three types of BTCs, there are obvious discrepancies between the observed and the simulated BTCs by the ADE model, while the fitted values of the LAF model are very close to the observed values. The overall shapes of BTCs simulated by the ADE model are different from the observed BTCs, and the predicted values are much smaller than those observed values after the intermediate stages of BTCs. The precision of the two models was also compared by r2 and RMSE. The values of r2 and RMSE of the ADE and LAF models for Experiments Two and Three are shown in Tables 2 and 3. The r2 value obtained by the LAF model is larger and the corresponding RMSE is smaller than that obtained by the ADE model. Most of the r2 values are over 0.95, and the corresponding RMSE values are below 0.040, illustrating that the LAF model is suitable for simulating and predicting solute transport in porous media and has a better fitting effect. The result is consistent with Abgaze and Sharma [21] who found that the asymptotic dispersion model simulated the observed BTCs was better than constant and linear distance-dependent dispersion models. Sabahi et al. [22] also achieved the solution of the 1D advection equation with the scale-dependent dispersion equation (ASDE) using the practical finite analytic (PFA) method in porous media under advection-dominated conditions, and they found that the linear function of the numerical solution had greater accuracy and stability in comparison with that of the analytical solution.

Solute transport in porous media does not satisfy the standard Fick transport. The simulated BTCs in Figures 4, 5, and 6 all show the phenomenon of early arrival to varying degrees, indicating that the fluid flow belongs to non-Fickian (or anomalous) behavior, and the ADE and LAF models do not capture the early arrival. There is no obvious trailing phenomenon of the BTCs in Figures 4, 5, and 6. The LAF model can fit well solute transport in porous media, and it has the best behavior in comparison with the ADE model for distance dependence of dispersivity.

4.3. Characteristics of L0 Value

The LAF solutions for one-dimensional solute transport are linear-asymptotic distance-dependent in a finite domain. A LAF solution includes two main parameters of L0 and k. The parameter L0 represents the maximum dispersivity that linearly changes with distance, and parameter k represents the growth rate of the dispersivity. L0 was fixed as 1.8 m for Experiment One in a homogeneous quartz chamber, and the fitted BTCs by the LAF model were in excellent agreement with measured values at distances of 200, 400, and 600 cm, respectively. In the heterogeneous porous media, when L0 was fixed as 7.4 m for Experiment Two and L0 was 5.7 m for Experiment Three, the LAF model could better fit the evolution of BTCs. The two groups of experiments had the same media, but the solute flow direction was opposite, which resulted in a L0 value with a difference of 1.7 meters. This might be attributed to the difference of porosity. The porosity gradually increased in the chamber with irregular sand particles, and the growth rate of the dispersivity k began to have deviating decreases. The estimated value of L0 is 1.8 m in homogeneous media, while L0 values are 7.4 and 5.7 m in heterogeneous media. The L0 value in homogeneous media is much smaller than that in heterogeneous media, indicating that the linear-asymptotic distance-dependent dispersivity reached a stable peak within a shorter travel distance. A larger L0 makes BTC smaller, which due to a larger L0 could result in a larger dispersivity [23, 24].

4.4. Comparison with TPL Model

The continuous time random walk (CTRW) model with a truncated power-law (TPL) memory function has been widely applied to characterize non-Fickian transport. In this study, the effects of flow travel distances (L=200,400,600cm) at the same flow rate (Q=0.7mL/s) in homogeneous porous media on the fitting curves of the LAF and TPL models were compared and analyzed. As shown in Figure 7, predictions of the LAF and TPL models have excellent performance in solute transport at flow travel distances of 200 and 400 cm, while the two models deviated from the observed value at 600 cm. Neither model captures the phenomenon of early arrival. From Figure 7, the TPL model better fits the experimental results. However, for the homogeneous porous media, good linear-asymptotic fitting results of LAF can be obtained when the distance from the critical value L0 is less than 1.8 meters. Therefore, the LAF model better characterizes the activity trend at the peak solute concentration stage.

This study mainly focuses on the investigation of 1D solute transport in porous media. Experiments were performed in a glass chamber of 600cmlength×10cmwidth×2cmheight and filled with different types of sand to verify the solution of the models. Three experiments were conducted to determine the changes in solute concentrations with distance downstream and flow rates. LAF and ADE models are used for the fitting of solute transport in the finite domain porous media and compared with the TPL model. The results showed that the LAF model fits better than the ADE model for the solute transport in the porous media, which is consistent with the results of r2 and RMSE. The thresholds of L0 value were obtained under the three conditions. When the porous media is homogeneous quartz sand, the value of L0 is 1.8 m. When the porous media is composed of the same heterogeneous media but the flow direction was opposite, the L0 values are 7.4 m and 5.7 m, respectively, indicating that the flow direction had a great effect on the L0 value in the same porous heterogeneous media. The inhomogeneity of the medium seriously affects the L0 value, and the L0 value in a homogeneous medium is significantly smaller than that in heterogeneous media, which may be due to the porosity of the medium. In the process of solute transport in the three experiments, the LAF and TPL models can capture the behavior of solute transport as a whole. However, for fitting the early arrival, the LAF model is less than the TPL. This takes into consideration that the LAF model is more concise in mathematics, which will have a greater advantage in site applications. Further study will be needed in the future.

The [DATA TYPE] data used to support the findings of this study are included within the article.

The authors declare that they have no conflicts of interest.

This study was supported by the National Natural Science Foundation of China (grant numbers 41831289 and 41877191) and Key R&D Program of Anhui Province (grant numbers 1804a0802198 and 201904a07020071). Authors thank Yang Zhou’s work in experiments.

Exclusive Licensee GeoScienceWorld. Distributed under a Creative Commons Attribution License (CC BY 4.0).