AN ADVECTION-DILUTION MODEL TO ESTIMATE CONDUIT
GEOMETRY AND FLOW
OCENA GEOMETRIJE IN TOKA V KRAŠKIH PREVODNIKIH Z ADVEKCIJSKO RAZREDČITVENIM MODELOM
Guangquan LI1* & Hong LIU1
Abstract	UDC 556.3:551.44
Guangquan Li & Hong Liu: An advection-dilution model to estimate conduit geometry and flow
As conduit water carries pollutants introduced from a point recharge downgradient, a contaminant-free seepage from the surrounding limestone is added into the conduit water and actively dilutes the pollutants. In this study, the transport with advection and this active dilution but no dispersion, is solved using the method of characteristics. The complete solution considering initial condition, boundary condition, and production is presented. Scale analysis reveals that the model is applicable when injection of pollutants at sinkholes persists long enough but unsuitable for analysis of dye tracing experiments. An approach combining the model with dye tracing experiment is used to quickly estimate the geometry and flow of a major conduit from Ames Sink to Wakulla Spring, Northwest Florida.
Key words: karst conduit, dilution, method of characteristics, breakthrough curve.
Izvleček	UDK 556.3:551.44
Guangquan Li & Hong Liu: Ocena geometrije in toka v kraških prevodnikih z advekcijsko razredčitvenim modelom
V članku predstavimo model advekcijskega prenosa točkovno injiciranega sledila skozi kraški prevodnik, ob upoštevanju razredčenja zaradi dotoka sveže vode iz okoliškega masiva. Disperzije ne upoštevamo. Diferencialno enačbo z robnimi in začetnimi pogoji, rešimo analitično z metodo metodo karakteristik. Analiza velikostnih redov pokaže, da je model uporaben pri dolgotrajnem vnosu sledila, manj primeren pa pri analizi sledilnih poskusov. V članku prikažemo uporabo modela v kombinaciji s sledilnim poskusom pri oceni geometrije in pretoka v kraškem prevodniku med ponorom Ames in izvirom Wakulla na severovzhodnem delu Floride. Ključne besede: kraški prevodnik, razredčenje, metoda karakteristik, prebojna krivulja.
INTRODUCTION
A karst aquifer is featured with interconnected solution caves or conduits which provide preferential pathways for groundwater flow and contaminant migration (Shuster & White 1971; Kiraly 1998). Transport of contaminants in karst aquifers is driven by a fast turbulent flow in cavernous conduits as well as a very slow laminar flow in the limestone matrix (including small pores, fissures and fractures). These dual flows
are the most distinct feature of karst aquifers, and thus contaminant migration in the aquifers has two completely different time scales. Contaminants entering the aquifers from sinkholes with point recharges flush through the conduits after days to weeks, while pollutants seeping downward with non-point recharges via the surface reside in the aquifers for tens of years (Even et al. 1986).
1 International Joint Research Center for Karstology, yunnan University, 5 xueyun Rd., Kunming, yunnan 650223, China,
e-maü: guangquanli74@gmaü.com (Guangquan LI), hongliu@ynu.edu.cn (Hong LIU) Received/Prejeto: 18.04.2013
During a storm event, the head in a conduit may exceed the head in the matrix, such that water and contaminants are emplaced from the conduit into the matrix, which later get released back into the conduit when the pressure difference between the conduit and the matrix becomes reversed. That scenario has been discussed by Field (1993) and Li et al. (2008). In this paper, we shall only discuss an ordinary hydrogeologic condition in which contaminant-free water is driven from the matrix into a conduit and actively dilutes pollutants introduced from a sinkhole or doline. Water seeping from the matrix into conduits is typically slow and pollutant-free, as the potential contaminants have enough time to experience adsorption by rock minerals as well as decomposition by bacteria. In contrast, pollutants entering from sinkholes, dolines, swallow holes, grikes, etc. have much less time and available surface area for such physical adsorption or biochemical reactions (Field & Pinsky 2000), and can therefore contaminate the aquifer to a large extent in a short period.
Decrease of contaminant concentration in a conduit is caused by three primary mechanisms: a passive decrease induced by hydrodynamic dispersion in the conduit (Taylor 1954), a passive retention/release between mobile- and immobile-water regions, and an active dilution by the matrix seepage. There is significant prior work regarding transport in a karst conduit. An advection-dispersion equation with a supplemental equation accounting for retention in immobile-water regions has been developed by Toride et al. (1993) to describe transport in soils and was first introduced by Field & Pinsky (2000) to transport in a single conduit. The model is called 2RNE (two-region non-equilibrium). The weakness of that model includes -1) it requires specifying fitting parameters such as an exchange coefficient between mobile- and immobile-water regions; 2) water exchange between the conduit and the surrounding matrix is not considered such that conduit flow was assumed to have a constant velocity. Despite these disadvantages, the 2RNE is widely used for analysis
of transport in karst conduits, because retention/release was well modeled and it can successfully reproduce the often-observed strong skewness in spring breakthrough curves, e.g., Field & Li (2011). Birk et al. (2006) used a hybrid method (with MODFLOW and the Darcy-Weis-bach equation to simulate matrix seepage and conduit flow, respectively) to model the discharge response and the breakthrough curve at a spring. In that scenario, water entering the sinkhole is solute-free, while water released from the matrix is rich in Calcium. Using the method of characteristics, Li (2009) derived an analytical solution for the case where the seepage flow from the matrix carries solutes into a conduit, being diluted by a solute-free sinkhole flow. Li & Loper (2011) developed a new model in which advection, dilution and dispersion were included, sought an approximate analytical solution for the initial-value problem, and successfully simulated a dye tracing experiment performed by Davies et al. (2004) between Ames Sink and Indian Spring, Northwest Florida. Later, Li (2011) used transform of variables to obtain the solution for the initial-value problem as well as the solution for the boundary-value problem of that model. In this paper, we are interested with a simple transport model that can be used to quickly estimate geometry and flow of a conduit.
In our conceptual model, contaminant enters a conduit with the sinkhole recharge, is diluted by the uncon-taminated seepage released from the surrounding limestone, and eventually is transported to a spring. (Dissolved minerals and decomposed organic matters from the matrix into the conduit will be considered afterward in this paper). We firstly construct the mathematical model and present its solution by intentionally neglecting conduit dispersion. Secondly, an approach is used to estimate the geometry and flow of a major conduit in the Woodville Karst Plain, Northwest Florida. Thirdly, the role of conduit dispersion is analyzed and a sufficient condition for applying the advection-dilution model is offered.
ADVECTION-DILUTION IN A CONDUIT
THE GOVERNING EQUATION The essential features of the advection-dilution in a conduit may be described using a relatively simple one-dimensional equation. Ignoring longitudinal dispersion and assuming zero production (i.e., no contaminant enters the conduit from the matrix), solute mass conservation requires:

(1)
where C [M/L3] is the solute concentration in the conduit (averaged over the cross-sectional area), W [L/T] is the speed of conduit water averaged over the conduit cross-section, z [L] is the distance downstream from the
sinkhole, and t [T] is time. Assuming a circular conduit with constant radius a [L] and a constant specific flux q [L/T] of water entering the conduit from the matrix, the velocity of conduit flow is increased in the downstream direction according to water mass conservation, i.e.,

(2)
where W0 is the averaged speed at the sinkhole (z = 0) and
T =
2q-
(3)
The conduit domain is 0 < z < Z, where Z is the downstream location of the spring. Equation differs from the traditional advection-dispersion equation in that it neglects dispersion but incorporates the specific water flux q, which causes active dilution. The linear increase of velocity with downstream position z is due to the addition of matrix water uniformly (assumed) into the conduit. Equation (2) does not contradict the Dar-cy-Weisbach equation. Combining Equation (2) with the Darcy-Weisbach equation would yield a non-linear pressure distribution along the conduit that is possible. If the seepage flux q is set to zero, the flow velocity along the conduit would be a constant and the pressure would degenerate into a linear distribution.
Equation (1) is a first-order partial differential equation, the solution of which requires specifying an initial condition in the conduit, i.e.,
C(2;,0)=Cj(2:),for2:a0,	(4)
and a boundary condition at the conduit entrance, i.e.,
C{0,t)=Cß)Jort2:0.	(5)
ANALYTICAL SOLUTION USING METHOD OF CHARACTERISTICS Substituting Equation (2) into (1) yields a first-order wave equation with the characteristic equation (Strauss 1992) being
dz dt
Z+Wj T
C ■
(6)
The characteristic curves, plotted in Fig. 1, are obtained by integrating the first two terms of Equation (6) and expressed as
z+W„r=W„Te«-i>\
(7)
where n is a constant of integration that represents the intersection of the curves at the t axis. From the Lagrangian viewpoint, the characteristic is the displacement-time curve for fluid particle, and thus
according to (7), ^^
dt
^W^, which is correct because
fluid particle at the conduit entrance has the speed W0.
The characteristics fill the first quadrant of the z-t plane as n varies from to The characteristic for n = 0 passes through the points (0, 0) and (Z, T) where
r=rlnfj-+l]	(8)
is the travel time for fluid particle to reach the spring from the sinkhole. This characteristic divides the first quadrant of the z-t plane into two domains. Characteristics in domain B, having n < 0, intersect the line z = 0, t > 0 and solutions in this domain satisfy the boundary condition, while characteristics in domain I, having n < 0, intersect the line z > 0, t = 0 and solutions in this
Fig. 1: Characteristic curves, z + W0T = W0re<t-n)/T. As n changes, the curves fill the first quadrant of the plane. I and B denote the domains determined by the initial condition and the boundary condition, respectively.
domain satisfy the initial condition. The characteristics have positive slope, so that as time goes on, progressively more solute concentration in the conduit is determined by the solution having characteristics in domain B.
The two families of solutions applicable in the two domains are easily obtained from Equation (6). Integrating the second and third terms of Equation (6) yields a family of solutions satisfying the initial condition, while the second family of solutions is obtained by integrating the first and third terms of that equation and satisfying the boundary condition. The combined solution is expressed as (see Appendix A for detailed derivation):
C(z,0=Cj((z+in domain I. (9a)
Ciz,t)=C^{t-Tln (^+l)) ^^ , in domain B. (9b)
For t < t, the solute concentration at the spring (z = Z) is determined by the initial condition, while for t > T, the concentration is determined by the boundary condition.
It is interesting to explore the physical meaning of Equation (9). Similar to Equation (8), the travel time for fluid particle to reach the downstream location z from
the conduit entrance is
MwA
Therefore,
+ljj in Equation (9b) represents the solute
concentration of the sinkhole flow at time t-r\a.

fl
For Equation (9a), supposing the fluid particle reaching z at time t starts from downstream location z*, the traveling
time from z* to z will be	This time is equal
to t, from which we get z* = (z + W0r)e'^T-W0r that appears in Equation (9a). This is the starting position of fluid particle in the conduit that will reach z at time t. To facilitate understanding the exponential term in Equation (9a), we need to rewrite the formula for z* into
—e""'^. Using Equation (2), this equation can be
Z+WqT
further rewritten as / = er^'', the left side of which W{z)
represents exactly the ratio of the contaminated water at z* diluted by the side seepage from z* through z.
The solution presented in Equation (9) assumes that water entering the conduit from the matrix contains no solute. The case where the matrix water contains solute is a production-value problem, which has been analytically solved by Li (2009). The general problem in which solute can preexist in the conduit (the initial-value problem), can enter the conduit from the sinkhole (the boundary-value problem), and can release from the matrix into the conduit (the production-value problem), is a superposition of these three problems, because the governing equation is linear with respect to concentration. The complete solution is listed in Appendix B.
examples
This section contains two simple solutions illustrating the behavior of the breakthrough curves generated at the spring (z = Z) from two classes of initial and boundary conditions.
STEP INITIAL CONDITION AND ZERO BOUNDARY CONDITION In this case, the initial concentration of solute in the conduit is assumed to be a step function, i.e.,
C,(z)=
0
0
0<z<zj z<z<z^, z<z<Z
(10)
and water entering the conduit at the sinkhole contains no solute: CB(t) = 0. Equation (10) is an idealization of contaminant hotspot. Substituting Equation (10) into (9) and evaluating the result at the spring yields
0
az,t)=
where f=rln(
Q<t<t
C^e-"^ t<t<t^, 0

Z+WqT

(11)
(12)
for i = 1, 2. This breakthrough curve is plotted in Fig. 2. Water containing solute that is initially at z = z2 arrives at the spring first and is relatively undiluted, while water that travels farther experiences greater dilution, so that the concentration at the spring decays with time. Since longitudinal dispersion has been ignored, there is no spreading in the breakthrough curve.
Fig. 2: Spring breakthrough curve generated by solute hotspot preexisting in a conduit (no solute enters from the sinkhole).
C(z = 0)'			
	-1		
			
C{z = Z)			
Q	-1		
		W.T + Z	
			


fig. 3: Spring breakthrough curve generated by a rectangular pulse of solute introduced from the sinkhole (the conduit is initially solute-free).
ZERO INITIAL CONDITION AND STEP BOUNDARY CONDITION In this case, the conduit initially contains no solute, i.e., CI(z) = 0, and a pulse of solute enters the conduit at the sinkhole uniformly over a prescribed interval of time:

0
0
o<^<^J t<t<t^.
(13)
Substituting Equation (13) into (9) yields
C{t,Z)=
0
t<T<+t^ T+t<t<T+t^, T+t<t
(14)
where t is given by (8). As can be seen in Fig. 3, this breakthrough curve is the result of a simple binary mixing between conduit water and matrix water. In contrast to the first example, the solute concentration at the spring is constant.
APPLICATION AND LIMITATION
This section focuses on the extraction of information regarding geometry and flow of a conduit from analysis of the breakthrough curve, employing the solution developed above.
Referring to the solution in Equation (14) and to Fig. 3, the peak concentration of solute at the spring, Cp, relative to that entering at the sinkhole is given by

Qo

Z+W,r C„ Qs'
(15)
where Qp and Qs are the water fluxes at the sinkhole and the spring, respectively.
Substituting Equation (15) into (8), it follows that
T =

(16)
Substituting Equation (16) into (3) gives
Qs

(17)
which can be used to estimate the specific seepage from the matrix into the conduit. On a related note, it is
necessary to know the travel time T. Noting
Equations (15) and (16) may be combined to estimate the conduit length:
(Q,-Q„)r
^ 7tann(Q/Q„)"
(18)
There are two ways to get the ratio Cp/C0, or Q0/QS. The first is to chemically measure the peak concentration of a long-persisting tracer rich in surface water but
devoid in matrix water, while the second is to physically measure the water fluxes at the sinkhole and spring. The travel time t from sinkhole to spring can be measured by a dye tracing experiment. Thus these two quantities can be entered into the model to calculate the seepage flux and the conduit length via Equations (17) and (18), respectively. This approach does not require dealing with dispersion.
A REAL WORLD EXAMPLE The results given above can be applied to extract information about conduit geometry and flow. To demonstrate this, a dye tracing experiment conducted by Davies et al. (2004) in the Woodville Karst Plain, Northwest Florida, USA, is chosen. Conduits in this plain are phreatic, as the plain is very close to the coast of the Gulf of Mexico. In the dye tracing experiment, dye traveled from Ames Sink to Wakulla Spring, for a straight-line distance of 9.3 km after 528 hours (Dav-ies et al. 2004). Allowing for an estimated tortuosity 1.3 (Li 2012), the actual traveling distance is approximately 12 km. The conduit radius ranges from 1 to 15 m and the spring discharge ranges from 5 to 20 m3/s (Davies et al. 2004). The recharge to Ames Sink is estimated to be 0.01 m3/s (Kincaid 2004). Thus the conduit radius and matrix seepage are calculated through Equations (17) and (18); see Tab. 1.
Essentially, this approach to estimate the parameters of a conduit inputs the time of the peak concentration
Tab. 1: Parameters of geometry and flow of the major conduit from Ames Sink to Wakulla Spring, Northwest Florida. Numbers with * present calculated values.
we use the transport time in the whole conduit t, the total length of the conduit Z, the sinkhole recharge Q0 and the spring discharge QS, to calculate the radius of the upstream segment and the velocity of the seepage from matrix into conduit (refer to Appendix C). The result is listed also in Tab. 1. With the radius ratio being 0.7, the cross-sectional area of the downstream segment is double that of the upstream segment. Though the cross-section varies significantly, the resulting seepage velocity only changes a little. Therefore, the seepage velocity is not sensitive to the variation of conduit cross-section and the uniform-radius model is reliable.
THE ROLE OF DISPERSION Releasing a tracer with a step function of time into a sinkhole, the spreading length can be estimated as Ax=8-'IdJ (see Figure 10.6.1, Bear 1972, p. 628), wherein Dc is dispersion in the conduit. Assuming Dc~aW (i.e., conduit dispersivity is scale-variant and of the same order as conduit radius.) and inserting T~Z/W, it follows that Ax=8^(äZ. Thus the spreading time in the breakthrough curve is	This spreading can be neglected when and only when
At^ SVÖZ^SVÖZ«! T. Wr, Lp '
(19)
Parameter	Value		Units
	Uniform conduit	Two segments	
Conduit length, Z	12	12	km
Travel time, T	528	528	hour
Spring discharge, QS	10	10	m3/s
Water flux at sinkhole, Q0	0.01	0.01	m3/s
Conduit radius, a	8.54*	7.54*, 10.8*	m
Specific water flux across wall, q	0.0155*	0.0145*	mm/s
on the breakthrough curve from dye tracing experiment into our advection-dilution model. (The peak amplitude on the breakthrough curve from dye tracing experiment can never be used into our model, for conduit dispersion invariably decreases it significantly). This approach is based upon the fact that the time of the peak concentration is not sensitive to dispersion; see Fig. 2, Li & Loper (2011).
To test validation of the uniform-conduit assumption, we now divide the conduit into two segments with different radius. Assume the radius radio between the upstream segment and the downstream segment, k = 0.7,
where Tc is the duration of injection, and Lp = Wtc is the length scale of the plume. Therefore, as the conduit length Z increases, dispersion shall not be ignorable.
Only for release of a long solute pulse (large Tc) at the sinkhole, or very small a or Z, conduit dispersion may be neglected.
The Peclet number is defined as
A a '
being widely used to quantify the ratio between advection and dispersion at a certain distance (Bear 1972). It is worth noting that a large Peclet number only states that dispersion may be neglected at that distance, but in no way necessarily assures that dispersion can be ignored in the transport process from the sinkhole to the spring. Taking dye tracing experiments as example, it is not recommended to neglect dispersion. Inequality (19) suggests how long the solute should last (or how short the domain should be) such that the above dispersion-free model can behave well.
Inequality (19) can be rewritten as 8Z
M
Tc

Supposing Lp < Z, it is evident that only at very large Peclet number (i.e., Pe > 10000), dispersion can be neglected with respect to advection. This value is far larger than any Peclet numbers ever observed for conduit flows in karstic aquifers. Therefore, it is concluded that our model is generally not applicable to transport in a real conduit. The only hope for application of this model is under an extreme condition that the duration of solute
injection must be very long, or equivalently Lp >> Z. As a reminder, the above condition is a very conservative (sufficient) one, because it is based upon a rectangular tracer pulse injected at sinkhole. In practice, a natural solute at sinkhole may have a much more gradual shape, rather than a sharp rectangle, which has been analyzed in Li et al. (2010).
DISCUSSION
Conduit flow in a karst aquifer is often a combination of water entering at a sinkhole or doline and water released from the surrounding limestone matrix. Contaminants can enter a conduit with the sinkhole recharge and will be diluted by the uncontaminated seepage from the matrix. This advection-dilution process has been solved analytically using the method of characteristics and the solution is presented in Equation (9). Without conduit dispersion, there is no spreading in the shape of spring breakthrough curves; see Fig. 2 and 3.
The issue whether dispersion can be neglected or not depends on whether it is evaluated locally or globally. The former refers to solute transporting over a length scale of the solute plume, while the latter refers to solute transporting over the conduit length. The latter is more rigorous than the former. Li et al. (2010) analyzed the former and got «1, as a sufficient and necessary con-
Lp
dition for neglecting dispersion locally. That inequality can be rewritten as ^ -Lj«l, being more accurate
than the traditional condition, -i^«!. In this paper we are further concerned with the global ignorability of dispersion (a stricter condition), and inequality (19) is offered as a sufficient criteria for assessment of the applicability of the model in this paper.
Numerical simulation often suffers from a notorious problem that it cannot accurately simulate the transport of contaminant plume near the sharp edge of the plume. This is because when the Peclet number is large, dispersion is small with respect to advection, such that undesired numerical dispersion will become dominant to decrease the accuracy of numerical solutions, even causing unphysical numerical oscillations. In that case, very fine grids would have to be adopted to guarantee a small cell Reynolds number (Thomas 1995) such that a realistic solution can be obtained. Sun (1989) has drawn a similar conclusion that when dispersion is very small, it is better to adopt particle-tracking method (ignoring dispersion) to simulate transport.
CONCLUSIONS
In this paper, an analytical solution to dilution of a conservative contaminant in a conduit due to water infiltrating from the surrounding limestone (referred to as seepage) is presented. The main assumptions involved in the model include -1) there is negligible dispersive transport in the conduit; 2) the conduit radius is uniform along the conduit; and 3) the seepage rate across the wall is uniform along the conduit. The first assumption limits the applicability of the model (e.g., unsuitable for dye tracing experiments in which dye is injected instantaneously so that Tc is too small to satisfy inequality (19)). Nevertheless, the model may be applicable to the case that contaminant is naturally introduced into a
conduit by long-duration rainstorm events. The second assumption has been tested to be valid. The effect brought by the third assumption is unclear and may be a good topic for future study.
The solutions presented by Li & Loper (2011) and Li (2011) considered advection, dilution and dispersion, but were approximate solutions. This paper presents the exact solution to a transport model only considering advection and dilution. Therefore, this model may be used as a benchmark. Scale analysis shows that for a rectangular injection at a sinkhole, inequality (19) provides a sufficient condition for neglecting conduit dispersion. Nonetheless, this inequality is conservative
(too restricting) and is not necessary for application of the model to a natural/gradual injection.
The model in this paper is useful mainly in two aspects. 1) The solution provides a useful benchmark to check numerical solutions. That is, numerical solutions of transport in a conduit at the limit of small dispersion (large Peclet number) should converge to the analytical solution expressed by Equation (9). 2) The model can be applied to field study with dye tracer test. As demonstrat-
ed, if travel time from sinkhole to spring is measured by a dye tracing experiment and the conduit length, sinkhole recharge and spring discharge are estimated beforehand, then they can be inputted into the model to quickly extract the average conduit radius and seepage velocity, without any needs to deal with dispersion. In the future, we will investigate how the variation of specific seepage affects transport in conduits.
ACKNOWLEDGEMENTS
This research was supported by the National Science Foundation of China under grant 41162008 and grant 41371040. The authors are deeply grateful to an anony-
mous reviewer and Editor-in-Chief Franci Gabrovšek for their insightful comments and constructive suggestions.
REFERENCES
Bear, J., 1972: Dynamics of Fluids in Porous Medium. New York, Dover, 764 pp.
Birk, S., Liedl, R. & M. Sauter, 2006: Karst spring responses examined by process-based modeling.-Ground Water, 44(6), 832-836, doi: 10.1111/j.1745-6584.2006.00175.x.
Davies, G., Kincaid, T., Hazlett, T., Loper, D., DeHan, R. & C. McKinlay, 2004: Why do quantitative groundwater tracing? Lessons and examples from the Woodville Karst Plain of North Florida. GSA 2004 Denver Annual Meeting, 18, November 7-10, Paper No. 49-22.
Even, H., Magaritz, M. & R. Gerson, 1986: Timing the transport of water through the upper vadose zone in a karstic system above a cave in Israel.- Earth Surface Processes and Landforms, 11, 181-191.
Field, M., 1993: Karst hydrology and chemical contamination.- Journal of Environmental Systems, 22(1), 1-26.
Field, M. & G. Li, 2011: Inversion for the input history of a dye tracing experiment.- Journal of Cave and Karst Studies, 73(1), 16-20, doi: 10.4311/jcks2010es0143.
Field, M. & P. Pinsky, 2000: A two-region nonequilib-rium model for solute transport in solution conduits in karstic aquifers.- Journal of Contaminant Hydrology, 44(3-4), 329-351.
Kincaid, T., 2004: Ames Sink Tracer Test. http://www. globalunderwaterexplorers.org/content/ames-sink-tracer-test (accessed November 1, 2012).
Kiraly, L., 1998: Modeling karst aquifers by the combined discrete channel and continuum approach.- Bulletin du. Centre d'Hydrogeologie, 16, 77-98.
Li, G., 2009: Analytical solution of advective mixing in a conduit.- Ground Water, 47(5), 714-722, doi: 10.1111/j.1745-6584.2009.00575.x.
Li, G., 2011: Spatially varying dispersion to model breakthrough curves.- Ground Water, 49(4), 584592, doi: 10.1111/j.1745-6584.2010.00777.x.
Li, G., 2012: Calculation of karst conduit flow using dye tracing experiments.- Transport in Porous Media, 95(3), 551-562, doi: 10.1007/s11242-012-0061-6.
Li, G., Loper, D. & R. Kung, 2008: Contaminant sequestration in karstic aquifers: Experiments and quantification.- Water Resources Research, 44, W02429, doi:10.1029/2006WR005797.
Li, G. & D. Loper, 2011: Transport, dilution, and dispersion of contaminant in a leaky karst conduit.- Transport in Porous Media, 88(1), 31-43, doi:10.1007/s11242-011-9721-1.
Li, G., Shang, Y. & J. Gao, 2010: Scale analysis of the significance of dispersion in mixing-transport in conduits.- Journal of Cave and Karst Studies, 72(3), 150-155, doi: 10.4311/jcks2009es0106.
Shuster, E. & W. White, 1971: Seasonal fluctuations in the chemistry of limestone springs: A possible means for characterizing carbonate aquifers.- Journal of Hydrology, 14, 93-128.
Strauss, W, 1992: Partial Differential Equations: An Introduction. New York, John Wiley & Sons, 425 pp.
Sun, N., 1989: Mathematical Modeling of Groundwater Pollution. Beijing, The Geological Press of China, 361 pp (in Chinese).
Taylor, G., 1954: tte dispersion of matter in turbulent flow through a pipe.- Proceedings of the Royal Society (London) A, 223, 446-468, doi:10.1098/ rspa.1954.0130. ttomas, J., 1995: Numerical Partial Differential
Equations. New York, Springer, 436 pp. Toride, N., Leij, F. & M. van Genuchten, 1993: A comprehensive set of analytical solutions for nonequi-librium solute transport with first-order decay and zero-order production.- Water Resources Research, 29(7), 2167-2182.
APPENDIX A
ttis appendix contains a detailed derivation of solution (9). From the second and third terms of Equation (6), along the characteristic (refer to Fig. 1) we have
dC(z.t) _ -dt C	T ■
Integrating the two sides of Equation (A1) and using initial condition (4), we get Inserting Equation (7) into (A2), it follows that
Please note that the initial condition requires (see Fig. 1)
ri=t-rln(^+l)<0
On the other hand, from the first and third terms of Equation (6), along the characteristic we have
dC(z,t) ^ -dz C WoT+z'
Integrating the two sides of Equation (A5) and using boundary condition (5) yield
Inserting Equation (7) into (A6), it follows that
C(z.t)=CJt-rH^+l))
tte boundary condition requires (see Fig. 1)
^WoT "WoT+Z
ri>0.
Equations (A3), (A4), (A7) and (A8) can be summarized as follows:
C{z,t)=
C,((2:+in domain I
CJt-rln
+1
W,r jlW^r+z
^^^ , in domain B
(Al)
(A2)
(A3)
(A4)
(A5)
(A6)
(A7)
(A8)
(A9)
APPENDIX B
ttis appendix presents a generalization of solution (9) to the case that water entering the conduit from the matrix contains solute (a production-value problem). tte solution for this case has been presented by Li (2009). Since the mathematical problem is linear with respect to concentration, that solution can be superposed on solution (9); the complete solution is
t
C{z.t)=
+	in domain I.
0
(Bl)
^qT
in domain B.
where CM is the solute concentration in the matrix, j(z, t) is the dimensionless specific flux of solute at the conduit wall, n is defined in Equation (7) and ZD is defined by
(B2)
APPENDIx C
This appendix seeks how to calculate the radius and seepage velocity of a conduit consisting of two segments with different radius. First we suppose the conjunction is at the middle of the conduit and the radius ratio between the upstream segment and the downstream segment is k. tte conduit discharge at the conjunction is denoted as Qm. Water mass conservation yields
Dividing (C1) by (C2) yields
O Qo+kQs
1+k
Now we apply Equation (18) to the upstream segment and the downstream segment, which yields
7t«^ln(Q^/Q„) 2
Qs-Qm(t-T,) ^ Z nal\n(Q/QJ 2
where T^ is the transport time within the upstream segment. Diving (C4) by (C5) yields


Substituting (C3) into the above equation and solving for T^ yields
\Qo+kQs
fcrin

(i+fc)Qo
In
(1+k)Q,
Q0+kQ.

Qo+fcQs
(i+fc)Qo
(Cl) (C2)
(C3)
(C4) (C5)
(C6)
(C7)
With T1 known, we substitute (C3) into (C4) to solve for the radius of the upstream segment.

(l+fc)Zln
. (i+fc)Q„
With a1 known, substituting (C3) into (C1) yields the seepage velocity
^ ^(1+fc)Zfl1
In summary, (C7), (C8), and (C9) can yield the transport time within the upstream segment, the radius of the upstream segment, and the velocity of seepage from matrix into conduit.