Next Article in Journal
Performance, Modeling, and Cost Analysis of Chemical Coagulation-Assisted Solar Powered Electrocoagulation Treatment System for Pharmaceutical Wastewater
Next Article in Special Issue
Modelling Plume Development with Annual Pulses of Contaminants Released from an Airport Runway to a Layered Aquifer, Evaluation of an In Situ Monitoring System
Previous Article in Journal
Spatial Distribution of Pine Pollen Grains Concentrations as a Source of Biologically Active Substances in Surface Waters of the Southern Baltic Sea
Previous Article in Special Issue
Combined Column Test for Characterization of Leaching and Transport of Trace Elements in Contaminated Soils
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

A Method for Calibrating the Transient Storage Model from the Early and Late-Time Behavior of Breakthrough Curves

1
Department of Industrial Engineering, University of Padova, 35131 Padova, Italy
2
Department of Land Environment Agriculture and Forestry, University of Padova, 35020 Padova, Italy
*
Author to whom correspondence should be addressed.
Water 2023, 15(5), 979; https://doi.org/10.3390/w15050979
Submission received: 25 January 2023 / Revised: 21 February 2023 / Accepted: 1 March 2023 / Published: 3 March 2023

Abstract

:
Solute transport in rivers is controlled by mixing processes that occur over a wide spectrum of spatial and temporal scales. Deviations from the classic advection–dispersion model observed in tracer test studies are known to be generated by the temporary trapping of solutes in storage zones where velocities and mixing rates are relatively small. In this work, the relation between the early and late-time behavior of solute breakthrough curves (BTCs) and the key parameters of the Transient Storage Model (TSM) is analyzed using non-asymptotic approximations of the model equations. Two main slopes are identified corresponding to the rising and decreasing limbs of the BTCs which are linked by specific relationships to transport and storage parameters. The validity of the proposed approximations is demonstrated with both synthetic and experimental data. Consistent with the TSM assumptions, the range of validity of the proposed approximations represents the limit of separability between surface dispersion and transient storage and can be expressed as a function of a nondimensional parameter. The results of this work can help environmental scientists identify solute transport and transient storage parameters and support the design of enhanced field tracer experiments.

1. Introduction

Tracer tests are a widespread technique for characterizing the fate and transport of solutes in streams and rivers. Experiments are typically conducted by injecting a tracer upstream and observing the evolution of the tracer concentration at one or more sections downstream. Reach-averaged transport parameters are then estimated by fitting a one-dimensional model to the observed breakthrough curves (BTCs).
Initial studies on solute transport in rivers were based on Taylor’s theory of longitudinal dispersion, which was found to be incapable of fully representing the skewed and long-tailed BTCs observed in tracer tests. Chatwin proposed that the longitudinal dispersion coefficient could be estimated from the rising limb of the BTCs and argued that the late-time behavior of the BTCs resulted from solute trapping in the viscous sublayer [1]. In rivers, the deviations from the advection–dispersion equation are exacerbated by the temporary storage of solutes in sediment [2,3,4], vegetated zones [5] and lateral pockets [6,7,8], characterized by small velocities and slow mixing rates compared to the transport in the main channel.
Several transport models have been proposed in the literature to represent deviations from Taylor’s dispersion theory, particularly to describe long-term retention phenomena due to exchanges between the surface water and the hyporheic zone [8]. These include the Aggregated Dead Zone model (ADZ) [9,10], general Multi-Rate Mass-Transfer (MRMT) formulations [11], Fractional advection–dispersion models [12], the STIR model [13,14], the Continuous Time Random Walk (CTRW) [15], data-driven models [16], and models representing hyporheic exchange as an exponentially attenuated diffusion process [3,17]. Further advancement in the field of tracer test and solute-transport modelling includes the use of reactive tracers, such as the resazurin-resorufin tracers, for quantifying the exchange with transient storage zones that are metabolically active [18,19,20]. These experimental techniques require the use of more sophisticated solute transport models [21,22,23] and quantification of deficiencies in the mass balance [24,25,26,27].
The Transient Storage Model (TSM) [28] has been widely used to analyze tracer test data in field and laboratory settings. The TSM combines a classical advection–dispersion equation with a first-order mass transfer into a storage zone of finite volume and has been shown to provide reasonably good approximations of tracer test breakthrough curves [28,29,30], although better fits can be obtained by assuming two storage components with different detention time [31,32,33]. A number of works have analyzed the capability of the TSM to consistently represent the dynamics of solute retention in streams [34,35,36].
Whilst these models are often capable of providing a good representation of the observed BTCs, uncertainties can arise in the interpretation of the model parameters. Nordin et al. [37], and more recently Gonzalez-Pinzon et al. [38], showed that the BTCs of tracer tests in rivers exhibit a persistent skewness that is inconsistent with the predictions of traditional 1D solute transport models. Recent studies [32,39] also observed that the transient storage parameters obtained from calibration of a 1D solute transport model with tracer test data can depend on the length of the study reach. Haggerty et al., found that distinct storage mechanisms acting at different timescales can combine to generate a power-law behavior of the tail of the concentration curves [40].
Even though more sophisticated models are available, the TSM is still widely used to investigate solute transport and exchange processes in rivers, thanks to its simplicity, accessibility, and ability to provide useful hints for understanding solute transport [41]. In the present work, a graphical approach is presented for estimating the parameters of the TSM. This is done by (i) deriving specific approximations for the initial rising limb and the final decreasing limb of the concentration curves, and (ii) identifying the conditions which allow separate determination of transport and storage parameters from the shape of the BTCs. The applicability of the newly derived approximations as a method to derive the parameters of the TSM from the independent fit of the rising and decreasing limb of a BTC is demonstrated with a numerical test and with experimental data. The method also provides insights into the limit of separability between surface dispersion and transient storage according to the Transient Storage Model.

2. Methods

In the TSM and ADZ models, the river channel is conceptually divided into two mutually interacting domains: a main flow channel with cross-sectional area A , in which the bulk flow occurs, and a storage domain of finite cross-sectional area A S adjacent to the main channel. The cross-sectional average concentration of solute is assumed to be homogeneous in the main channel and in the storage area and is denoted by C W and C S , respectively. The coupled differential equations governing the evolution of the concentration along the river are:
C W t + U C W x D W 2 C W x 2 = A S A C S t
C S t = α A A S C W C S
where U is the advective velocity, D W is the longitudinal dispersion coefficient in the bulk flow, and α is the exchange rate of solute at the storage zone-bulk flow interface. The model above implies an exponential residence time distribution (RTD) in the dead zones, φ t = ( 1 / T D ) e t T D , with T D = 1 α A S A   as the mean residence time in the storage zone and can be shown to be equivalent to other, more general residence time formulations under the assumption of a single storage domain with exponential RTD [13,14,40]. Following Davis et al. [42], the solution to Equations (1) and (2) for an instantaneous injection of mass M of conservative tracer at x = 0 and at t = 0 can be written as:
C W x , t = C T x , t e α t + e t T D 0 t C T x , τ e τ T D α T D e α τ τ t τ I 1 2 α T D τ t τ d τ
where I 1 is the first-order modified Bessel function of the first kind, and
C T x , t = M 2 A π D W t e x U t 2 4 D W t
is the solution to Equation (1) when transient storage is disregarded. In the alternative formulation proposed by Davis et al. [10], referred to as the Aggregated Dead Zone (ADZ) model, the longitudinal dispersion term in Equation (1) is dropped, D W = 0 , and the corresponding solution for an instantaneous mass injection at time t = 0 is:
C W x , t = M A U δ t x U e α t +
+ M A U H t x U e t x U T D e α x U α T D x U t x U I 1 2 α T D x U t x U
where δ t is the delta Dirac function and H t is the Heaviside function.

2.1. Approximation for the Rising Limb of the BTC

Equation (3) implies that C T e α t is the portion of mass M that never entered the storage domain. As a first approximation, this term describes the rising limb of a BTC, observed at a section x following an instantaneous injection of mass M at x = 0 :
C R δ x , t = M 2 A π D W t e x U t 2 4 D W t   e α t
An approximation for this limb of the curve in case of a step injection of tracer of duration T S at x = 0 can be derived from the convolution of C W , expressed by Equation (3), with the concentration at the boundary. Neglecting the contribution of the second term of Equation (3), the approximating expression for the rising limb of the BTC is found as
C R H x , t = Q C T x , t M × C 0 H t H t T S = 0 T S C 0 Q M C T x , t τ d τ = C 0   2   e α x U 3 U 2 α D W erf U 3 t U 2 x + 2 α D W x 2 U 3 D W x erf 2 α D x U 2 U T S + x U t 2 U 3 D W x
where C R H is the approximation of the rising limb of the BTC, and H(t) is the Heaviside function. The error function, erf y , can be approximated near y = 0 by its McLaurin series:
erf y = 2 π n = 0 1 n y 2 n + 1 2 n + 1 n !
By limiting the series expansion to the first two terms and substituting into Equation (7) under the frozen cloud approximation used by Chatwin [1], a linear approximation is found for the rising limb at a distance x from the injection point in the case of a finite-length step injection at x = 0 :
C R H x , t = m x t + q x
where
m x = C 0 U e α   x U 4 π D W x U
and
q x = C 0 e α   x U 2 1 2 π x 4 D W x U
Equation (9) represents the rising limb of a BTC at a section x resulting from a step concentration C 0 at x = 0 from t = 0 to the time t = T S . This approximation applies around the point t = x / U , where C W x , t C 0 / 2 .

2.2. Approximation for the Decreasing Limb of the BTC

Equation (5) can be used to obtain an approximation for the decreasing limb C D δ X , t of a BTC following an instantaneous injection of mass M at x = 0 . In the early stages of development, the decreasing limb is represented by Equation (5), provided that the term describing the concentration distribution in Equation (3), given by C T x , t e α t , is negligible for t 0 . The modified Bessel function of first order and first kind can be approximated by a linear equation I 1 c c / 2 in the range 0 < c < 1 [43], or, more accurately, by the series [44]:
I ν c c ν ν Γ ν 1 + c 2 4 ν + 1 + c 4 32 ν + 1 ν + 2 + c 8 6144 n = 1 4 ν + n
which is valid for c 12 or c ν . Using a second order approximation, thus limiting the series expansion of I 1 to the second term of Equation (12), the expression of I 1 in Equation (5) can be approximated as follows:
I 1 2 α T D x U t x U α T D x u t x U + 1 2 α T D 3 / 2 x U 3 / 2 t x U 3 / 2
Equation (13) is valid for x / U < t < x / U + T D / 2 α x / U . By substituting Equation (13) into Equation (5), the following approximating expression for the decreasing limb of the BTC is obtained:
C D δ x , t = M α   x U T D Q e α x U e t   x U T D 1 + α 2 T D x U t x U
For x / U < t < x / U + T D / 2 α x / U , another approximation can be derived,
1 + α 2 T D x U t x U e α 2 T D x U t x U
which leads to the final approximated expression for Equation (14):
C D δ x , t = a x e n x t
where:
a x = M α x U T D Q e α 2 T D x U 2 α T D 1 T D   x U
and:
n x = α 2 T D x U 1 T D
Equation (16) is convenient because it makes it possible to fit the observed BTC at a distance x from the injection point by linear regression in semi-log scale.
Equation (5) describes the behavior of the decreasing limb of a BTC generated by an instantaneous injection of tracer when the effect of the mass transfer into the dead zones is much greater than the effect of longitudinal dispersion at late times. For a step injection of tracer at x = 0 , the expression for the decreasing limb must be suitably modified. In this case, the concentration at the injection section can be expressed as C 0 H t H t T S . If a solute mass M is injected at x = 0 with a constant rate over the time T S starting at t = 0 , the approximated expression of the decreasing limb can be found from the convolution integral,
C D H x , t = 0 T S Q C 0 M a x e n x t τ d τ = Q C D δ x , t M × C 0 H t H t T S
The solution to Equation (19) is again an exponential function:
C D H x , t = b x e n x t
where:
b x = e α 2 T D   x U 1 T D T S + x U + e α 2 T D   x U 1 T D x / U 2 C 0 α α x U 2 x U e α x U
Given an observed BTC at a distance x, the coefficients m x , q x , b x and n x can be determined by performing a linear regression of the rising limb on a linear scale and a linear regression of the decreasing limb on a semi-log scale. Equations (10), (11), (18), and (21) then represent a system of nonlinear equations that can be solved numerically to obtain the unknown transport and storage parameters, U , D W , α , and T D .
To determine the limit of applicability of the proposed approximation, it is convenient to introduce the non-dimensional time, defined as t * = t U / x . Using these non-dimensional quantities, Equation (13) can be written as:
I 1 2 α T D x U t * 1 α T D x U t * 1 + 1 2 α T D 3 / 2 x U 3 t * 1 3 / 2
valid for 1 < t * < t L I M * , where t L I M * is a non-dimensional parameter defined as:
t L I M * = 1 + T D U 2 2 α x 2
The parameter t L I M * represents the upper time limit for which Equation (22) approximates the modified Bessel function. The comparison between the approximated and the original Bessel function is reported in Figure 1 for four different Bessel functions, characterized by different values of the parameter t L I M * . The lower limit, t * = 1 , corresponds to the non-dimensional advective time t A D * = t A D / x / U , where t A D = x / U . For the complete list of the symbols used in the equations, refer to Table A1 in Appendix A.

3. Application

A numerical test was performed to validate the approximations derived in the previous Section 2.1 and Section 2.2. Experimental data were used for calibrating the unknown TSM parameters using the proposed graphical method.

3.1. Numerical Experiment

The validity of the approximations derived above is tested on an ideal case, by (i) simulating a BTC from a TSM model using previously fixed A , A S , D W , and α , and (ii) plotting the approximation lines obtained from the Equations (10), (11), (18) and (21). The flow discharge is assumed to be Q = 0.4 m3 s−1 with stream width b = 5.0 and depth d = 0.4 m, respectively. The dispersion coefficient is estimated using Fischer’s formula, D W = 0.011 U 2 b 2 / U * d = 0.735 m2s−1 [45], where the shear velocity is U * = g d j = 0.037 m s−1 and j is the energy slope under the assumption of uniform flow and given a Manning coefficient n = 0.05 m−1/3s−1. The model output is generated using a single compartment TSM with storage area A S = 0.1 m2, and exchange rate α = 10 4 s−1 at three cross-sections located at X 1 = 500 m, X 2 = 1000 m and X 3 = 1500 m from the injection point. The boundary condition at x = 0 is a constant concentration C 0 injection of a solute mass M = 192 g, applied for a period T S = 480 s starting at t = 0 s.
In Figure 2, the dimensionless output concentration   C * x , t = C W x , t / C p x , t p , where C p x , t p is the peak concentration at time t p , is plotted at the three sections mentioned above in both semi-log and linear scales. The approximations for the rising and the decreasing parts of the BTC are plotted according to Equations (9) and (20) in linear and semi-log scale, respectively, where parameters m * x = m x / C p x , t p , q * x = q x / C p x , t p , b * x = b x / C p x , t p and n x are evaluated from Equations (10), (11), (18) and (21) using the same parameters used for the TSM. The values of A , A S , D W , α , m * , q * , b * , n and t L I M * are reported in Table 1 for each section.

3.2. Field Experiment

The applicability of the model for evaluating unknown TSM parameters from measured BTCs is demonstrated on experimental data obtained from a real tracer test. The data are taken from the experimental work of Bottacin-Busolin et al. [32], conducted on the Desturo canal, an irrigation channel located in Northern Italy. The studied reach is divided into two sections located at X 1 = 262 and X 2 = 567 m, respectively, from the injection point. The plateau injection of a mass of tracer M = 3.24 g is T S = 1470 s long and the flow discharge is Q = 0.045 m3 s−1.
At both section X 1 and X 2 , a linear and an exponential trendline are applied to the leading and the decreasing limb, respectively, considering the portion of data in the leading and decreasing limb which maximizes the coefficient of determination R2. The fitting lines are reported in Figure 3, showing the observed concentration data, normalized by the peak concentration, and the fitting trendlines in linear scale for the rising limb (panels a and b) and semi-log scale for the decreasing limb (panels c and d). The fitting lines provide the values of m * x and q * x of the linear trendline for the rising part, C R * x , t = m * x t + q * x , and the values of b * x and n of the exponential trendline for the decreasing limb of the curve, C D * x , t = b * x e n x t . The unknown transport parameters A , A S , D W and α are obtained by solving Equations (10), (11), (18), and (21) as a system of equations. Note that, here, the flow rate Q is known and therefore the cross-sectional velocity U can be determined from the cross-sectional area A as U = Q/A, whereas the parameter T D in Equations (18) and (21) depends on α, A, and AS, according to the relationship T D = 1 α A S A   . The resulting transport parameters and t L I M * are reported in the first two rows of Table 2.
To check the quality of the fitting, the obtained transport parameters are then used to generate a BTC using the TSM. A comparison between the simulated BTCs (dashed-dot blue line) and the observed data (circles) is reported in Figure 3 in linear (panel a and b) and semi-log scales (panel c and d).
The observed BTCs are then numerically fitted with the TSM by using differential evolution as an optimization algorithm. The computed parameters A , A S , D W , and α for the optimized fitted BTC are reported in Table 2 for comparison with the parameters obtained using the proposed approximations. The best fitting BTCs are reported in Figure 3 (dashed red line).

4. Results and Discussion

The results of the application presented in Figure 2 show that the goodness of the graphical estimation of the transient storage parameters and longitudinal dispersion coefficient varies at the three sections, suggesting that the method can be applied under certain conditions. In particular, it can be noted that the approximations are good at the closest section X 1 , while at X 2 , they start to deviate from the rising and decreasing part of the BTC. At the farthest section, X 3 , the decreasing part of the BTC is not well fitted by the trendline obtained via Equations (19) and (21). This indicates that this curve does not behave as the exponential function at late time of Equation (20). Similarly, the rising part is no longer described by Equation (9). There are two reasons for this behavior: firstly, the term C T = C W e α t in Equation (4) becomes negligible at longer times and for higher exchange rates α ; secondly, at late times, Equation (13) is not a good approximation of the Bessel function. At late times and large distances from the injection point, the contribution of transient storage to the mixing process becomes more important than longitudinal dispersion, and they become two overlapping processes. This means that the rising limb of the curve no longer represents the pure effect of longitudinal dispersion. At the same time, the decreasing limb of the curve does not behave as an exponential function (Equation (13)) anymore. The dimensionless upper time limit of validity for Equation (22), t L I M * , can represent a numerical index for identifying the conditions of applicability of the approximation method proposed in this work. For sections X 1 = 500 m, X 2 = 1000 m, and X 3 = 1500 m, the values of t L I M * are 1.40, 1.10, and 1.04, respectively. These values suggest that transient storage does not produce an exponential signature on the decreasing limb when the upper time limit of validity of Equation (22), t L I M * , is close to 1. Although there is no definite value for t L I M * , the results suggest that, in general, the exponential signature is not easily distinguishable from the longitudinal dispersion effect when t L I M * < 1.1 .
The application of this approach to the experimental data confirmed its validity for t L I M * > 1.1 . The results show that when the approximations for the leading and decreasing limbs of the observed BTCs are applied, the BTC generated by the TSM with the graphically estimated parameters for section X 1 = 262 m (dashed-dotted lines in Figure 3a,c) optimally fits the observed data. At this section, the upper nondimensional time limit of validity for Equation (22) is t L I M * = 1.33 > 1.10 . The set of parameters obtained with the proposed graphical method can provide a good first approximation of the parameters obtained by the numerically calibrated parameters. Their relative percentage errors are: +15% for D W , +1% for A , −21% for A S , and −29% for α . For the mean residence time in the storage zone, T D , evaluated according to Equation (3), the relative error between the graphical parameter estimation and the numerically calibrated model is +10%. In contrast, at section X 2 = 567 m, where t L I M * = 1.02 , the graphical approximation does not lead to a well-fitted BTC (dashed-dotted lines in Figure 3a,b): the time of the peak concentration is poorly caught and both the leading and decreasing limbs are not well approximated. At this section, the model calibration based on a global optimization algorithm yields a very different set of parameters, as shown in Table 2. The relative errors between the graphically obtained and the numerically optimized parameters are −89% for D W , −9.5% for A , −0% for A S , −230% for α and −67% for T D .
The graphical method presented in this work thus provides a simple preliminary way to evaluate transient storage parameters in field experiments. The empirical rule t L I M * > 1.1   can be also expressed in terms of the advective time t A D according to the following relations:
t L I M * = 1 + T D U 2 2 α X 2 = 1 + A s 2 α 2 A 1 t A D 2 > 1.1  
t A D < 1 0.1 A s 2 A α 2
which are more likely satisfied when the BTC is observed at relatively short distances from the injection.

5. Conclusions

Approximations were derived for the early and late time behavior of tracer test breakthrough curves for both slug and plateau injections. The approximations assumed that the transport of a solute in a river can be described by an advection-dispersion first-order mass transfer equation. Consistent with the model, the tracer BTCs measured at a section downstream of the injection point exhibited two main slopes, associated with the initial rising part of the concentration curve and the final decreasing part after the concentration peak. The parameters of the model could be inferred directly from these two slopes using the analytical relations presented in this work. The results showed that the graphical estimation of transient storage parameters and longitudinal dispersion coefficient could be made if two conditions were satisfied: first, when the two mechanisms did not produce an overlapping effect on the rising limb of the BTC and, secondly, when the condition of approximation of the decreasing limb as an exponential function was satisfied. However, a specific condition applies to the distance from the injection point for the validity of these approximations, which represents a limit of separability between longitudinal dispersion and first-order mass transfer in immobile domains. As a general empirical rule, this limit can be expressed by the condition t L I M * < 1.1 , where t L I M * is a nondimensional time given by Equation (23). When this condition is not satisfied, longitudinal dispersion and transient storage become two undistinguishable mixing processes that cannot be uniquely separated by fitting the model to a tracer BTC. The condition above also provides a lower limit for the timescales of solute retention that can be inversely modelled from a tracer experiment. This limit should therefore be considered in the design of tracer tests and subsequent analyses to provide sensitive parameter estimates and identify sources of uncertainty. The results of this work provide a twofold contribution to current modelling and experimental approaches for field tracer studies: (1) by providing a simple method for estimating the surface and subsurface storage parameters; and (2) by identifying a temporal limit for the separation of surface and transient storage processes which can be used for the design of improved field tracer studies.

Author Contributions

Conceptualization, M.Z., A.B.-B. and A.M.; methodology, E.D., M.Z., A.B.-B. and A.M.; data curation, E.D. and M.Z.; formal analysis, E.D., M.Z. and A.B.-B.; writing—original draft preparation, E.D., M.Z. and A.B.-B.; writing—review and editing, E.D., A.B.-B. and A.M.; visualization, E.D. and A.B.-B.; supervision, A.M. All authors have read and agreed to the published version of the manuscript.

Funding

This work is partly funded by the Italian National Recovery and Resilience Plan (PNRR), Extended Partnership 3 (PE3) “RETURN”.

Data Availability Statement

The sources of the field tracer data used for validating the proposed method were properly cited (and referred to in the reference list).

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A

Table A1. Notation list of all the symbols used in the manuscript.
Table A1. Notation list of all the symbols used in the manuscript.
SymbolUnitDescription
a (kg m−3) coefficient   of   the   trendline   for   C D δ
A (m2)mean flow area
A S (m2)transient storage area
b ,   b * (kg m−3), (-) coefficient   and   normalized   coefficient   of   the   trendline   for   C D H
C 0 (kg m−3)concentration at the injection section
C R H (kg m−3)approximation for the rising limb of the BTC
C R * (-)dimensionless linear approximation of the rising limb
C D H ,   C D δ (kg m−3)approximation for the decreasing limb of the BTC
C D * (-)dimensionless exponential approximation of the decreasing limb
C S (kg m−3)concentration in the storage area
C T (kg m−3)elementary solution of the advection–dispersion equation
C W (kg m−3)concentration in the main flow channel
D W (m2 s−1)longitudinal dispersion coefficient
H (-)Heaviside function
I 1 (-)modified Bessel function of the first order and first kind
m ,   m * (kg m−3 s−1), (s−1) slope   and   normalized   slope   of   the   trendline   for   C R H   and   C R *
n (s−1) exponent   of   the   trendline   for   C D H   and   C D *
q ,   q * (kg m−3), (-) intercept   and   normalized   intercept   of   the   trendline   for   C R H   and   for   C R *
Q (m3 s−1)discharge
t (s)time
t A D (s)advective time
T D (s)mean residence time in the storage area
T S (s)time duration of a plateau injection
t L I M * (-)dimensionless time limit
U (m s−1)flow velocity
x (m)longitudinal distance from injection point
α (s−1)exchange rate
δ (s−1)Dirac delta function
τ (s−1)dummy variable for time
φ (s−1)exponential RTD in the dead zones

References

  1. Chatwin, P.C. On the interpretation of some longitudinal dispersion experiments. J. Fluid Mech. 1971, 48, 689–702. [Google Scholar] [CrossRef]
  2. Boano, F.; Harvey, J.W.; Marion, A.; Packman, A.I.; Revelli, R.; Ridolfi, L.; Wörman, A. Hyporheic flow and transport processes: Mechanisms, models, and biogeochemical implications. Rev. Geophys. 2014, 52, 603–679. [Google Scholar] [CrossRef]
  3. Bottacin-Busolin, A. Non-Fickian dispersion in open-channel flow over a porous bed. Water Resour. Res. 2017, 53, 7426–7456. [Google Scholar] [CrossRef] [Green Version]
  4. Kwaw, A.K.; Dou, Z.; Wang, J.; Zhang, X.; Chen, Y. Advancing the knowledge of solute transport in the presence of bound water in mixed porous media based on low-field nuclear magnetic resonance. J. Hydrol. 2023, 617, 129059. [Google Scholar] [CrossRef]
  5. Salehin, M.; Packman, A.I.; Wörman, A. Comparison of transient storage in vegetated and unvegetated reaches of a small agricultural stream in Sweden: Seasonal variation and anthropogenic manipulation. Adv. Water Resour. 2003, 26, 951–964. [Google Scholar] [CrossRef]
  6. Magazine, M.K.; Pathak, S.K.; Pande, P.K. Effect of Bed and Side Roughness on Dispersion in Open Channels. J. Hydraul. Eng. 1988, 114, 766–782. [Google Scholar] [CrossRef]
  7. Jackson, T.R.; Haggerty, R.; Apte, S.V.; O’Connor, B.L. A mean residence time relationship for lateral cavities in gravel-bed rivers and streams: Incorporating streambed roughness and cavity shape. Water Resour. Res. 2013, 49, 3642–3650. [Google Scholar] [CrossRef]
  8. Chen, Z.; Tian, Z.; Zhan, H.; Huang, J.; Huang, Y.; Wei, Y.; Ma, X. The Effect of Roughness on the Nonlinear Flow in a Single Fracture with Sudden Aperture Change. Lithosphere 2022, 2022, 5775275. [Google Scholar] [CrossRef]
  9. Wallis, S.G.; Young, P.C.; Beven, K.J. Experimental investigation of the aggregated dead zone model for longitudinal solute transport in stream channels. Proc. Inst. Civ. Engin. 1989, 87, 1–22. [Google Scholar]
  10. Davis, P.M.; Atkinson, T.C. Longitudinal dispersion in natural channels: 3. An aggregated dead zone model applied to the River Severn, U.K. Hydrol. Earth Syst. Sci. 2000, 4, 373–381. [Google Scholar] [CrossRef]
  11. Haggerty, R.; McKenna, S.A.; Meigs, L.C. On the late-time behavior of tracer test breakthrough curves. Water Resour. Res. 2000, 36, 3467–3479. [Google Scholar] [CrossRef] [Green Version]
  12. Deng, Z.; Bengtsson, L.; Singh, V.P. Parameter estimation for fractional dispersion model for rivers. Environ. Fluid Mech. 2006, 6, 451–475. [Google Scholar] [CrossRef]
  13. Marion, A.; Zaramella, M. A residence time model for stream-subsurface exchange of contaminants. Acta Geophys. Pol. 2005, 53, 527. [Google Scholar]
  14. Marion, A.; Zaramella, M.; Bottacin-Busolin, A. Solute transport in rivers with multiple storage zones: The STIR model. Water Resour. Res. 2008, 44, W10406. [Google Scholar] [CrossRef] [Green Version]
  15. Boano, F.; Packman, A.I.; Cortis, A.; Revelli, R.; Ridolfi, L. A continuous time random walk approach to the stream transport of solutes. Water Resour. Res. 2007, 43. [Google Scholar] [CrossRef]
  16. Young, P.; Garnier, H. Identification and estimation of continuous-time, data-based mechanistic (DBM) models for environmental systems. Environ. Model. Softw. 2006, 21, 1055–1072. [Google Scholar] [CrossRef] [Green Version]
  17. Bottacin-Busolin, A. Modeling the effect of hyporheic mixing on stream solute transport. Water Resour. Res. 2019, 55, 9995–10011. [Google Scholar] [CrossRef]
  18. Haggerty, R.; Martí, E.; Argerich, A.; von Schiller, D.; Grimm, N. Resazurin as a “smart” tracer for quantifying metabolically active transient storage in stream ecosystems. J. Geophys. Res. Atmos. 2009, 114, G03014. [Google Scholar] [CrossRef] [Green Version]
  19. González-Pinzón, R.; Peipoch, M.; Haggerty, R.; Martí, E.; Fleckenstein, J.H. Nighttime and daytime respi-ration in a headwater stream. Ecohydrology 2016, 9, 93–100. [Google Scholar] [CrossRef]
  20. Knapp, J.L.A.; González-Pinzón, R.; Haggerty, R. The resazurin-resorufin system: Insights from a decade of “smart” tracer development for hydrologic applications. Water Resour. Res. 2018, 54, 6877–6889. [Google Scholar] [CrossRef]
  21. Argerich, A.; Haggerty, R.; Martí, E.; Sabater, F.; Zarnetske, J. Quantification of metabolically active transient storage (MATS) in two reaches with contrasting transient storage and ecosystem respiration. J. Geophys. Res. Atmos. 2011, 116, G03034. [Google Scholar] [CrossRef] [Green Version]
  22. Lemke, D.; Liao, Z.; Wöhling, T.; Osenbrück, K.; Cirpka, O.A. Concurrent conservative and reactive tracer tests in a stream undergoing hyporheic exchange. Water Resour. Res. 2013, 49, 3024–3037. [Google Scholar] [CrossRef]
  23. Bottacin-Busolin, A.; Dallan, E.; Marion, A. STIR-RST: A Software tool for reactive smart tracer studies. Environ. Model. Softw. 2020, 135, 104894. [Google Scholar] [CrossRef]
  24. Haggerty, R.; Argerich, A.; Martí, E. Development of a “smart” tracer for the assessment of microbiological activity and sediment-water interaction in natural waters: The resazurin-resorufin system. Water Resour. Res. 2008, 44, W00D01. [Google Scholar] [CrossRef]
  25. Lemke, D.; González-Pinzón, R.; Liao, Z.; Wöhling, T.; Osenbruck, K.; Haggerty, R.; Cirpka, O.A. Sorption and transformation of the reactive tracers resazurin and resorufin in natural river sediments. Hydrol. Earth Syst. Sci. 2014, 18, 3151–3163. [Google Scholar] [CrossRef] [Green Version]
  26. Chen, J.L.; Steele, T.W.J.; Stuckey, D.C. Metabolic reduction of resazurin; location within the cell for cytotoxicity assays. Biotechnol. Bioeng. 2017, 115, 351–358. [Google Scholar] [CrossRef] [PubMed]
  27. Dallan, E.; Regier, P.; Marion, A.; González-Pinzón, R. Does the Mass Balance of the Reactive Tracers Resazurin and Resorufin Close at the Microbial Scale? J. Geophys. Res. Biogeosci. 2020, 125, e2019JG005435. [Google Scholar] [CrossRef]
  28. Bencala, K.E.; Walters, R.A. Simulation of solute transport in a mountain pool-and-riffle stream: A transient storage model. Water Resour. Res. 1983, 19, 718–724. [Google Scholar] [CrossRef]
  29. Mulholland, P.J.; Marzolf, E.R.; Webster, J.R.; Hart, D.R.; Hendricks, S.P. Evidence of hyporheic retention of phosphorus in Walker Branch. Limnol. Oceanogr. 1997, 42, 443–451. [Google Scholar] [CrossRef]
  30. Valett, H.M.; Morrice, J.A.; Dahm, C.N.; Campana, M.E. Parent lithology, surface-groundwater exchange, and nitrate retention in headwater streams. Limnol. Oceanogr. 1996, 41, 333–345. [Google Scholar] [CrossRef]
  31. Briggs, M.A.; Gooseff, M.N.; Arp, C.D.; Baker, M.A. A method for estimating surface transient storage parameters for streams with concurrent hyporheic storage. Water Resour. Res. 2009, 45, W00D27. [Google Scholar] [CrossRef]
  32. Bottacin-Busolin, A.; Marion, A.; Musner, T.; Tregnaghi, M.; Zaramella, M. Evidence of distinct contaminant transport patterns in rivers using tracer tests and a multiple domain retention model. Adv. Water Resour. 2011, 34, 737–746. [Google Scholar] [CrossRef]
  33. Kerr, P.; Gooseff, M.; Bolster, D. The significance of model structure in one-dimensional stream solute transport models with multiple transient storage zones–competing vs. nested arrangements. J. Hydrol. 2013, 497, 133–144. [Google Scholar] [CrossRef]
  34. Marion, A.; Zaramella, M.; Packman, A.I. Parameter Estimation of the Transient Storage Model for Stream–Subsurface Exchange. J. Environ. Eng. 2003, 129, 456–463. [Google Scholar] [CrossRef] [Green Version]
  35. Zaramella, M.; Packman, A.I.; Marion, A. Application of the transient storage model to analyze advective hyporheic exchange with deep and shallow sediment beds. Water Resour. Res. 2003, 39, 1198. [Google Scholar] [CrossRef] [Green Version]
  36. Marion, A.; Zaramella, M. Diffusive Behavior of Bedform-Induced Hyporheic Exchange in Rivers. J. Environ. Eng. 2005, 131, 1260–1266. [Google Scholar] [CrossRef]
  37. Nordin, C.F.; Troutman, B.M. Longitudinal dispersion in rivers: The persistence of skewness in observed data. Water Resour. Res. 1980, 16, 123–128. [Google Scholar] [CrossRef]
  38. González-Pinzón, R.; Haggerty, R.; Dentz, M. Scaling and predicting solute transport processes in streams. Water Resour. Res. 2013, 49, 4071–4088. [Google Scholar] [CrossRef] [Green Version]
  39. Gooseff, M.N.; Briggs, M.A.; Bencala, K.E.; McGlynn, B.L.; Scott, D.T. Do transient storage parameters directly scale in longer, combined stream reaches? Reach length dependence of transient storage interpretations. J. Hydrol. 2013, 483, 16–25. [Google Scholar] [CrossRef]
  40. Haggerty, R.; Wondzell, S.M.; Johnson, M.A. Power-law residence time distribution in the hyporheic zone of a 2nd-order mountain stream. Geophys. Res. Lett. 2002, 29, 1640. [Google Scholar] [CrossRef] [Green Version]
  41. Knapp, J.L.; Kelleher, C. A Perspective on the Future of Transient Storage Modeling: Let’s Stop Chasing Our Tails. Water Resour. Res. 2020, 56, e2019WR026257. [Google Scholar] [CrossRef] [Green Version]
  42. Davis, P.M.; Atkinson, T.C.; Wigley, T.M.L. Longitudinal dispersion in natural channels: 2. The roles of shear flow dispersion and dead zones in the River Severn, U.K. Hydrol. Earth Syst. Sci. 2000, 4, 355–371. [Google Scholar] [CrossRef] [Green Version]
  43. Olver, F.W.J. Asymptotic and Special Functions; Academic Press: New York, NY, USA, 1974. [Google Scholar]
  44. Abramowitz, M.; Stegun, I.A. Handbook of Mathematical Functions; National Bureau of Standards Applied Mathematics Series-55; US Government Printing Office: Washington, DC, USA, 1972.
  45. Fischer, H.B. Discussion of “Simple method for predicting dispersion in stream”. J. Environ. Eng. Div. 1975, 101, 453–455. [Google Scholar] [CrossRef]
Figure 1. Comparison between the modified Bessel function I 1   and its approximation (Equation (22)) for four different values of t L I M * : (a) t L I M * = 1.5 , (b) t L I M * = 1.10 , (c) t L I M * = 1.05 , (d) t L I M * = 1.02 .
Figure 1. Comparison between the modified Bessel function I 1   and its approximation (Equation (22)) for four different values of t L I M * : (a) t L I M * = 1.5 , (b) t L I M * = 1.10 , (c) t L I M * = 1.05 , (d) t L I M * = 1.02 .
Water 15 00979 g001
Figure 2. BTCs generated by the TSM model at X 1 = 500 m, X 2 = 1000 m and X 3 = 1500 m and approximations for (a) the decreasing limb in semi-log scale and for (b) the rising limb in linear scale.
Figure 2. BTCs generated by the TSM model at X 1 = 500 m, X 2 = 1000 m and X 3 = 1500 m and approximations for (a) the decreasing limb in semi-log scale and for (b) the rising limb in linear scale.
Water 15 00979 g002
Figure 3. Normalized observed concentration data (“Cobs”, circles), the approximating trendlines (solid black lines), the simulated BTCs from the approximations (“Csim”, dashed-dot blue lines) and the optimized simulated BTCs (“Copt”, dashed red lines), for section X1 (panels (a) and (c) and section X2 (panels (b) and (d)). Panel (a) and (b) refer to the approximation of the rising limb of the BTC. Panel (c) and (d) are in semi-log scale and refer to the approximation of the decreasing part of the BTC.
Figure 3. Normalized observed concentration data (“Cobs”, circles), the approximating trendlines (solid black lines), the simulated BTCs from the approximations (“Csim”, dashed-dot blue lines) and the optimized simulated BTCs (“Copt”, dashed red lines), for section X1 (panels (a) and (c) and section X2 (panels (b) and (d)). Panel (a) and (b) refer to the approximation of the rising limb of the BTC. Panel (c) and (d) are in semi-log scale and refer to the approximation of the decreasing part of the BTC.
Water 15 00979 g003
Table 1. TSM parameters A , A S , D W , and α are used for the concentration curves in Figure 2, coefficients m * , q * , b * , and n of the approximation trendlines for the decreasing and rising parts at the three output sections X i , and peak concentration C p at each section used for normalizing the BTC in Figure 2.
Table 1. TSM parameters A , A S , D W , and α are used for the concentration curves in Figure 2, coefficients m * , q * , b * , and n of the approximation trendlines for the decreasing and rising parts at the three output sections X i , and peak concentration C p at each section used for normalizing the BTC in Figure 2.
SectionDistance D W A A s α t L I M * n b * m * q * C p
(m)(m2 s−1)(m2)(m2)(s−1)(-)(s−1)(-)(s−1)(-)(g m−3)
X 1 5000.7352.001.0 × 10−11.0 × 10−41.400−1.75 × 10−3 4.63 × 1012.04 × 10−3−4.320.502
X 2 10001.100−1.50 × 10−32.22 × 1031.63 × 10−3−7.270.347
X 3 15001.044−1.25 × 10−32.00 × 1041.31 × 10−3−8.960.274
Table 2. Coefficients m * , q * , b * , n of the approximation trendlines in Figure 3; TSM parameters A , A S , D W , α obtained from the trendline approximations (“Approx.”) and from a numerical optimized fit (“Optim.”) of the BTCs; estimation of t L I M * for each section and type of fitting, and peak concentration C p at each section used for normalizing the BTC in Figure 3.
Table 2. Coefficients m * , q * , b * , n of the approximation trendlines in Figure 3; TSM parameters A , A S , D W , α obtained from the trendline approximations (“Approx.”) and from a numerical optimized fit (“Optim.”) of the BTCs; estimation of t L I M * for each section and type of fitting, and peak concentration C p at each section used for normalizing the BTC in Figure 3.
SectionDistanceFit n b * m * q * D W A A s α t L I M * C p
(m) (s−1)(-)(s−1)(-)(m2 s−1)(m2)(m2)(s−1)(-)(mg m−3)
X 1 262Approx.−2.11 × 10−398.52.21 × 10−3−2.2780.2790.2022.93 × 10−24.00 × 10−41.3348.12
Optim.----0.2430.2003.72 × 10−25.65 × 10−41.21
X 2 567Approx.−1.80 × 10−34252.98.86 × 10−4−2.5820.0900.2441.86 × 10−24.55 × 10−41.0242.8
Optim.----0.5350.2701.86 × 10−21.37 × 10−41.16
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Dallan, E.; Bottacin-Busolin, A.; Zaramella, M.; Marion, A. A Method for Calibrating the Transient Storage Model from the Early and Late-Time Behavior of Breakthrough Curves. Water 2023, 15, 979. https://doi.org/10.3390/w15050979

AMA Style

Dallan E, Bottacin-Busolin A, Zaramella M, Marion A. A Method for Calibrating the Transient Storage Model from the Early and Late-Time Behavior of Breakthrough Curves. Water. 2023; 15(5):979. https://doi.org/10.3390/w15050979

Chicago/Turabian Style

Dallan, Eleonora, Andrea Bottacin-Busolin, Mattia Zaramella, and Andrea Marion. 2023. "A Method for Calibrating the Transient Storage Model from the Early and Late-Time Behavior of Breakthrough Curves" Water 15, no. 5: 979. https://doi.org/10.3390/w15050979

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop