Influence of CO2 dynamics on the longitudinal variation of incision rates in soluble bedrock channels: Feedback mechanisms

Influence of CO2 dynamics on the longitudinal variation of incision rates in soluble bedrock channels: Feedback mechanisms

Geomorphology 186 (2013) 85–95 Contents lists available at SciVerse ScienceDirect Geomorphology journal homepage: www.elsevier.com/locate/geomorph ...

2MB Sizes 1 Downloads 7 Views

Geomorphology 186 (2013) 85–95

Contents lists available at SciVerse ScienceDirect

Geomorphology journal homepage: www.elsevier.com/locate/geomorph

Influence of CO2 dynamics on the longitudinal variation of incision rates in soluble bedrock channels: Feedback mechanisms Matthew D. Covington a, b,⁎, Mitja Prelovšek a, Franci Gabrovšek a a b

Karst Research Institute ZRC SAZU, Postojna, Slovenia Department of Geosciences, University of Arkansas, Fayetteville, AR 72701, USA

a r t i c l e

i n f o

Article history: Received 10 March 2012 Received in revised form 11 December 2012 Accepted 12 December 2012 Available online 27 December 2012 Keywords: Karst Bedrock channel Speleogenesis Carbon dioxide Dissolution

a b s t r a c t We use measurements of dissolution rates of limestone tablets placed along a cave stream to estimate rates of modern incision. Dissolution rates within the stream display a systematic decrease with downstream distance. We discuss a variety of mechanisms that could be responsible for the longitudinal decrease in dissolution rates and develop simple mathematical models for each. The dissolutional length scales that arise from each model allow a first-order estimate of the plausibility of each mechanism and motivate further field studies to test each possibility. Water chemistry and other field data suggest that a decrease in the concentration of CO2 along the cave stream is responsible for the observed decrease in dissolution rates. We propose two potential mechanisms that could trigger this reduction in dissolved CO2 and discuss the plausibility of each mechanism in light of the field data collected. Either of these mechanisms introduces a feedback loop whereby the stream profile of a channel in soluble bedrock indirectly influences CO2 concentrations in the water, via either microbial or hydraulic processes. The CO2 concentration in turn effects the incision rates and therefore the future stream profile. This study illustrates the importance of CO2 dynamics in determining incision rates in a soluble channel and points to further modeling and field work that are needed in order to enable the development of realistic stream incision models in soluble strata. © 2013 Elsevier B.V. All rights reserved.

1. Introduction The initial formation stages of caves have been explored with a number of physically based mathematical models that couple dissolution and flow (Dreybrodt, 1990; Palmer, 1991; Groves and Howard, 1994; Dreybrodt, 1996). Beyond these initial one-dimensional models (1D), a variety of two-dimensional (2D) (e.g., Gabrovšek and Dreybrodt, 2000; Birk et al., 2005; Dreybrodt et al., 2005) and three-dimensional (3D; Kaufmann, 2009) models have been developed to study conduit network formation. While the dynamics of the early stages of cave formation have been extensively explored, the available models are typically inapplicable to the later stages of cave evolution. Similarly, mechanistic models are lacking for the incision of surface bedrock channels in soluble strata. In order to understand the development of mature cave systems or the incision of soluble bedrock channels, a number of additional processes and factors need to be included that are not typically accounted for in speleogenetic models, such as free surface flow (Ford and Ewers, 1978), sediment transport (Farrant and Smart, 2011), variable

⁎ Corresponding author at: Department of Geosciences, University of Arkansas, Fayetteville, AR 72701, USA. Tel.: + 1 479 575 3876; fax: + 1 479 575 3469. E-mail address: [email protected] (M.D. Covington). 0169-555X/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.geomorph.2012.12.025

discharge and chemistry (Groves and Meiman, 2005; Palmer, 2007), and carbon dioxide (CO2) dynamics (Anthony et al., 2003), among others. One reason for the lack of mathematical models for the late stages of cave formation is that many of the relevant processes are complex, relatively poorly understood, and difficult to study within a laboratory setting. One approach for improving knowledge of these processes is using in situ observation of incision. While this approach has its own limitations, because of potentially large variations in rates over time, it still provides useful data for the development and constraint of quantitative process models. A method that has been extensively used to quantify dissolution rates in karst terrain is weight-differencing of limestone tablets that are placed in the field. This approach has largely been used to determine the relative rates of denudation within different climates and different types of settings within a given location (e.g., Gams, 1981; Plan, 2005). Little of this work has focused specifically on stream incision. Prelovšek (2009) undertook a study of dissolution rates in cave streams across Slovenia. At many sites, tablets were distributed throughout a cave in order to explore variations in rates along a cave stream or at different elevations within a given channel cross section. Because of the coupling between flow and the processes controlling dissolution, longitudinal variations in dissolution rates provide a particularly powerful tool for constraining mechanistic models of incision (Covington et al., 2012). In some cases, characteristic length

86

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

scales will emerge from the combination of typical flow velocities and the characteristic timescale for a given process. In this case, observed length scales can be compared against the theoretically predicted length scales of different models. Here we use this approach to analyze the longitudinal decrease in dissolution rates observed in an allogenically recharged cave stream in Slovenia. Initial theoretical results motivate further observations of cave water chemistry that are used to test possible models for decreases in dissolution rates along the stream. In closing, we discuss the geomorphological implications of our results and note potential feedback mechanisms that result from the coupling of CO2 dynamics and channel incision.

2. Description of field site Lekinka is a 1-km-long cave situated at the NE edge of the Pivka basin in Slovenia (Fig. 1). It is a stream cave with an entrance at an elevation of 510 m at the lithological contact between Eocene flysch and overthrust upper Cretaceous limestones. A surface stream, Crni potok, flows into the entrance of Lekinka and drains a surface basin of noncarbonate sediments with an area of 1.05 km 2. The downstream end of Lekinka flows into an underwater tunnel that is too narrow to pass but emerges after only a few tens of meters into the underground Pivka River in the Postojna cave system. Along the 300 m of the cave study reach, one primary infeeder and two smaller infeeders are visible. The typical low flow discharge of the stream at the entrance of Lekinka is around 5 l s −1. During large floods, discharges in the Lekinka stream reach several cubic meters per second. The water that sinks into Lekinka is typically undersaturated with respect to calcite (Prelovšek, 2009).

3. Dissolution experiments 3.1. Methods Cumulative dissolution rates in the stream were measured with limestone dissolution tablets similar to those used by Gams (1981). The methodology is based on measuring the weight difference of limestone tablets before and after exposure to the stream (Fig. 2). Weight change was converted to incision rate using a limestone density of 2688 kg m −3 and the surface area of each tablet. The tablets were attached to the cave wall using a stainless steel screw through a central 8-mm hole in the tablet. Abrasion between metallic parts and the limestone tablet was minimized using felt washers. After removal from the field, tablets were dried for 15 d in the laboratory. A reference set of tablets was kept within the laboratory and weighed along with the field tablets. This allowed correction of errors resulting from the influence of relative humidity, which was found to be strongly correlated to weight variations (of about 0.01%) in the reference tablets. For more details on the methodology, see Prelovšek (2009). On the basis of past experiments with limestone tablets, the minimum observable change in surface height is ∼ 0.4 μm for 15–30 g limestone tablets (Prelovšek, 2009). This provides much higher accuracy than is possible using a microerosion meter (Spate et al., 1985). The accuracy of the weight-differencing method, in combination with the highly undersaturated water sinking into Lekinka, enabled measurements in 15-d intervals. 3.2. Results Dissolution was measured during the period from 15 April 2007 until 4 April 2009. Limestone tablets were placed at six locations

Fig. 1. Geological and hydrological setting of Lekinka (Prelovšek, 2009).

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

87

Fig. 2. The limestone tablets used in the dissolution experiments and the manner of affixing them to the cave wall (Prelovšek, 2009).

along the first 250 m of the stream in Lekinka (Fig. 3). Tablets were measured and replaced with fresh tablets every 15 d during the measurement period. During the 2-year measurement period, weather conditions were typical for the field site, with wet autumn and spring periods, a dry summer period, and one large flood. The dissolution rates summed over the entire measurement period are displayed in Fig. 4. The cumulative dissolution rates show an exponential decrease with downstream distance, where the best fit exponential function has a length scale of λ = 255 m (Fig. 4).

where F is the rate of limestone dissolution, α is the kinetic rate constant, C is the concentration of dissolved calcite, and Ceq is the equilibrium concentration of dissolved calcite (Palmer, 1991; Liu and Dreybrodt, 1997). Linear kinetics result in an exponential decrease in dissolution rates with distance along a cave stream; therefore, a first hypothesis might be that the dissolution rates are decreasing because the water is dissolving limestone and approaching saturation via this process. The e-folding length scale for dissolution is

4. Potential mechanisms for the longitudinal decrease in dissolution rates

λ¼

A number of potential mechanisms could produce the observed longitudinal decrease in dissolution rates. Here we use theoretical models of these mechanisms to assess their ability to explain the observed behavior and to illuminate potential field data that could constrain the mechanism.

(Eq. (15) in Covington et al., 2012), where Q is the volumetric discharge and Pw is the conduit wetted perimeter. Using a typical wetted perimeter for Lekinka of Pw =1.5 m and a kinetic constant α=2 ×10−7 m s−1, one can estimate that the discharge required to produce the observed e-folding length, λ=300 m, is Q=0.45 l s−1. This discharge is one to two orders of magnitude below the typical discharges observed in Lekinka (Table 1), suggesting that introduction of additional calcite into solution from the dissolution process itself is unlikely to be responsible for the observed decrease in rates. Expected dissolutional length scales under the flow conditions within Lekinka are in the range of 104 m≲λ ≲106 m. This conclusion can be tested in the field by determining the upstream and downstream dissolved loads, as well as the dissolved load that is being introduced via any significant infeeding streams.

4.1. Saturation caused by dissolution of calcite Calcite dissolution kinetics are linear in the highly undersaturated conditions that typically apply at the insurgence of Lekinka and can be approximated using   F ðC Þ ¼ α C eq −C

ð1Þ

Q Pw α

ð2Þ

Map of the Study Reach Q1

Primary infeeder (Q2, S2) S3 Cave reach (Lekinka)

continues to siphon

Q3,S4

Downstream Entrance Sinkhole - Location of limestone tables - Surface stream bank - Cave wall

S1

Surface reach (Crni potok)

N 0

15

30 m

Upstream Fig. 3. Map of the study reach showing 100 m of the surface stream, Črni potok, followed by the first 300 m of stream passage in Lekinka. The sites for the dissolution experiments, discharge measurements, (Q1-Q4) and water sampling (S1-S4 ) are depicted.

Fig. 4. Cumulative dissolution rates measured using limestone tablets along the stream in Lekinka over a two year period with increasing downstream distance from the entrance. The dashed line represents a best fit using the exponential equation and parameter values in the top right corner.

88

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

Table 1 Discharges measured within the main cave stream (Q1 and Q3) and in the largest infeeder (Q2).a Date

Q1

20 Sept 20 4 Nov 2011 5 Mar 2012 18 Apr 2012 13 June 2012

3.0 6.8 16.5 28.2 36.9

Q2 l l l l l

s−1 s−1 s−1 s−1 s−1

0.5 1.5 1.3 5.5 4.2

Q3 l l l l l

s−1 s−1 s−1 s−1 s−1

4.5 8.9 21.1 46.6 42.2

Other infeeders l l l l l

s−1 s−1 s−1 s−1 s−1

1.0 0.6 3.3 12.9 1.1

l l l l l

s−1 s−1 s−1 s−1 s−1

Percent diffuse 34% 24% 22% 39% 12%

a The discharge from other infeeders is calculated by subtracting Q2 from the difference between Q3 and Q1. The percent diffuse is the percentage of the total flow at the downstream end (Q3) that is not measured within the main stream at the upstream end of the reach (Q1).

4.2. Inflow of saturated water Another possible explanation for the decrease in dissolution rates is the addition of saturated water from more diffuse flow paths in the fractured and porous matrix. In fact, Lekinka contains a number of small tributary streams that carry much higher dissolved loads than the main stream (Prelovšek, 2009). Because the dissolutional length scale is orders of magnitude longer than the physical length, L, of the conduit (i.e., the process number Λ ≡ L/λ ≪ 1), we can treat the addition of diffuse water independently from the dissolution process (Covington et al., 2012). Consequently, if we assume that the added matrix water is at Ceq then, to explain the reduction in dissolution rates, the fraction of the water that has come from the matrix at the downstream station, fm is f m ¼ 1−

F s;downstream F s;upstream

ð3Þ

where Fs,downstream and Fs,upstream are the downstream and upstream dissolution rates, respectively. Therefore, to explain the observations, with a reduction in dissolution rates of Fs,downstream/Fs,upstream ≈ 0.36, ∼ 60 % of the water must have arrived from the matrix. This calculation assumes that the diffuse inputs are coming in at the equilibrium concentration for the main stream (Ceq). More likely the infeeding streams, which are entering through the soil zone, have a somewhat higher partial pressure of CO2 than the main stream and therefore may contain larger dissolved loads than Ceq for the main stream. As a result, the quantity of diffuse input required to produce the observed decrease in rates is probably somewhat less than the 60% figure above. Regardless of this complexity, the potential reduction of dissolution rates from diffuse input can be tested by measuring upstream and downstream discharges, comparing the upstream and downstream dissolved loads in the main stream, and quantifying the saturation state and dissolved load of the tributary streams. 4.3. Degassing of CO2 The final mechanism we consider is the degassing of CO2 from the water. In this case, the lowering of the concentration of CO2 within the water reduces the dissolution rates by decreasing the value of Ceq along the length of the stream (Eq. (1)). In general, the dissolved concentration of CO2 in the water is set by three factors: (i) the partial pressure of CO2 (pCO2) in the air above the water, (ii) the rate of production of CO2 in the water, and (iii) the reaeration rate. The flux of CO2 into the air is typically modeled using a linear rate law,   F CO2 ¼ k C w −α CO2 C a

ð4Þ

(Eq. (2) in (Wanninkhof et al., 2009)) where F CO2 (mol s −1 m −2) is the flux of CO2 from the water, k (m s −1) is the gas transfer velocity, Cw (mol m −3) is the bulk concentration of CO2 in the water, Ca

(mol m −3) is the concentration of CO2 in the air, and α CO2 is the dimensionless Ostwald solubility coefficient for CO2, which can be calculated using CO2 solubility and Henry's law and is on the order of one for the relevant temperatures (Wanninkhof et al., 2009). Because degassing is a linear process, a sudden change in the equilibrium concentration of dissolved CO2 as the stream enters the cave would result in an exponential longitudinal profile of CO2 concentration within the water, approaching the new equilibrium value over a characteristic length scale (Covington et al., 2012). To examine to first order whether CO2 degassing is a reasonable candidate mechanism, we calculate the degassing length scale and compare it against the observed length scale over which dissolution rates are changing. Typical degassing timescales of small streams are on the order of 10−1 to 10−2 d (Wanninkhof et al., 1990). Consequently, to produce the observed length scale of 300 m, flow velocities in the range of 3 cm s−1 to 0.3 m s−1 would be required. These flow velocities are close enough to expected flow velocities within the Lekinka stream to suggest that CO2 degassing is a plausible mechanism. The degassing of CO2 reduces dissolution rates by decreasing the equilibrium concentration of calcium ions, Ceq. The effect of CO2 on equilibrium for limestone dissolution can be approximated using "

C eq

K1Kc ¼ Cw 4K 2 γ Ca γ 2HCO3

#1=3 ;

ð5Þ

where K1, K2, and Kc are mass action constants and γCa and γ HCO3 are the activity coefficients for calcium and carbonate ions (Eq. (2.35c) in Dreybrodt, 1988). Therefore, Ceq varies with the one-third power of Cw, and the calcite dissolution rate varies nonlinearly with changes in Cw according to the combination of Eqs. (1) and (5). Whether CO2 degassing is responsible for the observed decrease in dissolution rates can be tested using measurements of stream pH and water chemistry along the length of the stream to estimate the concentration of CO2 in the water at each location, or via direct measurement of dissolved CO2 concentrations. 5. Field tests of mechanisms To test the potential mechanisms for the longitudinal decrease in dissolution rates, we collected a variety of field data including stream discharges, water quality parameters, and a survey of the stream longitudinal profile. Additionally, water samples were collected under a range of flow conditions, direct measurements of stream CO2 were made, and CO2 production rates were measured from samples of stream sediment. 5.1. Methods 5.1.1. Discharge Discharges were measured at three locations (Fig. 3) in order to quantify the addition of diffuse flow into the cave stream. These locations were near the upstream end of the cave (Q1), in the largest infeeder within the study reach (Q2), and at the downstream end of the study reach (Q3). Discharges were measured primarily during periods of relatively low flow, as this enables an estimation of the maximum contribution from diffuse recharge. All discharge measurements were conducted using the salt tracing method and a Schlumberger CTD Diver instrument that logged specific conductivity (SpC) values once per second. The CTD Diver measures SpC with an accuracy of ± 1 %. Salt traces were conducted using 100 to 400 g of salt depending on the discharge. The SpC values were converted to concentrations of dissolved NaCl using a linear relationship set to the values measured in three standard NaCl solutions with known concentrations of NaCl.

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

5.1.2. Field measurements of water quality Spot measurements of pH, temperature, and SpC were taken under a variety of flow conditions at locations in the surface stream (S1), near the upstream end of the cave stream (S2), in the largest infeeding stream along the studied reach (S3), and at the downstream end of the studied reach (S4). These locations are shown in Fig. 3. The measurements were acquired using a WTW multimeter that was calibrated prior to each field excursion. 5.1.3. Longitudinal profile and cave survey The studied surface and cave reaches were surveyed using a DistoX cave survey instrument, which has an angular accuracy of ± 0.5 ∘ for inclination and azimuth and a distance accuracy of ±1.5 mm. At each survey station, numerous splayed survey measurements allowed precise location of the cave walls and other features. To improve the vertical control of the survey, a water-leveling device was used. 5.1.4. Water samples Water samples were collected during a variety of flow conditions at four different locations: in the surface stream (S1), in the largest infeeder along the studied reach (S2), near the upstream end of the cave stream (S3), and at the downstream end of the studied reach (S4). All samples were analyzed in the laboratory within several hours of collection. Alkalinity was determined via titration with 0.02 N HCl; Ca and Mg concentrations were determined using titration with 0.01 M EDTA. The saturation index for calcite (SIcalc) and the dissolved CO2 were calculated by using PHREEQC (Parkhurst and Appelo, 1999) with the phreeqc.dat database and by employing the pH and temperature values measured in the field. In order to estimate errors in the calculations, 500 random initializations of the model inputs were run for each sample. For each of these random initializations, the pH and chemical parameters were chosen from a Gaussian distribution with a given variance about the measured value. Variances for alkalinity, Ca, and Mg were assumed to be 10% of the measured value; and the variance for pH measurements was assumed to be 0.05. Additionally, the influence of other ions (N, S, Cl, and P) was approximated using average values from historical water quality analyses at the site. Concentrations for these ions were also drawn randomly from Gaussian distributions with an average from the historical values and a variance of 50 %. In practice, we find that the concentrations of these ions are small enough to have negligible effects on SIcalc. This Monte Carlo procedure allows a relatively conservative estimate of error in the chemical analysis. 5.1.5. Direct measurements of dissolved CO2 During the later stages of this study, when it became apparent that CO2 dynamics played an important role in determining the longitudinal changes in dissolution rate, we equipped a Vaisala GMT series non-dispersive infrared (NDIR) CO2 sensor with a waterproof breathable membrane to allow for direct measurement of dissolved CO2. This technique is described in detail by Johnson et al. (2010). The sensor was used to measure dissolved CO2 at a higher spatial frequency along the stream than the water sampling sites. This sensor also allowed measurement of the dissolved concentration of CO2 at a higher accuracy than possible with the water chemistry analysis. For each measurement, the sensor was placed in the water and allowed to equilibrate to the dissolved CO2, which usually required 5–15 min. 5.1.6. CO2 production experiments Production rates of CO2 were measured in sediment samples from the field site. Sediment was taken from the main channel outside the cave, the main channel inside the cave, and from backwater locations inside the cave where finer sediments are deposited. These sampling locations covered the range of sediment types observed at the field site. Samples were taken back to the laboratory within an hour of

89

collection, and the volume of sediment in each sample was measured, with sediment volumes ranging from 50 to 300 ml. Each sample was placed in a bottle, which was filled to the top with degassed tap water and sealed. The dissolved CO2 was measured continuously using the NDIR sensor and a datalogger. Samples were stirred every 5 min to homogenize CO2 within the water column. Each experiment was run for an initial equilibration stage of 1 h, after which samples displayed roughly linear production curves. Experiments were then run for an additional hour to allow estimation of the CO2 production rate within each sample using a best fit linear relation for the postequilibration production curve. Initial experiments utilized water from the field site rather than tap water degassed to atmospheric pCO2. However, for samples obtained in the outside reach, the relatively high initial CO2 concentrations and high rate of production led to almost immediate saturation of the sensor, making it impossible to obtain accurate production rates. For the same reason, we used smaller volumes of sediment (typically 50 ml) for the outside samples, so that longer duration production curves could be obtained. The main channel cave sediment samples could not be similarly reduced in volume, because the larger grain sizes made it impossible to obtain representative samples with less than 200–300 ml of volume. To account for this variation in sample volume, production rates were normalized by sample volume. We normalize using volume instead of weight, because volume provides a rough proxy for stream bed surface area. 5.2. Results 5.2.1. Discharge Discharge values were measured at Q1, Q2, and Q3 on 20 September 2011 (during low flow conditions), on 4 November 2011 (during moderately low flow conditions), on 5 March 2012 (during moderate flow conditions), on 18 April 2012 (during moderate flow conditions), and on 13 June 2012 (during the recession following a small storm) (Table 1). The diffuse component of the flow varied from 12% to 39%. 5.2.2. Cave survey The planview of the studied reach can be seen in Fig. 3 along with the locations of the measurement, sampling, and dissolution experiment sites. The stream profile for the study reach is depicted in Fig. 5. Upon contact with the limestone and entry into the cave, the stream significantly steepens. The peak in stream bed elevation just inside the entrance corresponds to the accumulation of cobbles and boulders in the entrance collapse zone. 5.2.3. Water chemistry Water samples were collected on 8 June 2011 during high flow after a rain storm (on this day a sample was not collected at S3), on 6 July 2011 during a base flow period, on 20 September 2011 during a period of low flow during the recession following a small rain event, on 18 April 2012 during moderate flow, and on 13 June 2012 during a recession following a small rain event. A summary of these water chemistry data and field measurements is provided in Table 2. Alkalinity, pH, SI, and dissolved CO2 are depicted in Fig. 6 as a function of longitudinal distance in the stream. This allows visualization of the evolution of water chemistry along the studied reach. For all sampled dates, alkalinity remains relatively constant along the stream length, with the largest increase in alkalinity occurring during the recession measurement of 20 September 2011. Alkalinity within the main stream is lowest during the highest observed flows and highest during the lowest observed flows. For all dates, pH gradually increases along the length of the stream, with the lowest values observed on the surface and the highest values observed at the downstream end of the cave stream. This pH trend was also observed on numerous other dates where water quality parameters were measured without water sampling and during long-term datalogger

90

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

main stream. The pH is lower during higher discharges. The SI is also more stable within the infeeding stream than the main stream, with the infeeder displaying values that are slightly undersaturated during high flows and slightly oversaturated during lower flows. Dissolved CO2 within the infeeding stream is typically higher than in the main stream, but no clear correlation was seen with discharge or season for the dates measured.

Fig. 5. Elevation profile of the study reach. The dashed line at 0 indicates the entrance of the cave.

measurements of pH at S2 and S4 (data not shown). Additional spot measurements of pH along a 100-m reach of the surface stream upstream of S1 showed constant values of pH and SpC, suggesting relatively stable conditions in comparison to the cave reach. The SI increases monotonically with distance downstream for all samples and is inversely correlated with discharge, such that high discharges result in low (undersaturated) values of SI. The estimated dissolved CO2 decreases monotonically with downstream distance for all sampling periods. However, CO2 does not correlate tightly with discharge. Instead, values of CO2 are similar for the June and July 2011 sampling dates and somewhat lower for the September 2011 and spring 2012 sampling dates. For every sampling date, the alkalinity values of the infeeding stream were higher than in the main stream; however, the magnitude of the difference between the two values is positively correlated with discharge, such that the biggest differences are observed during the highest flows. The infeeder displays a more stable pH than the

5.2.4. Direct measurement of dissolved CO2 On 13 June and 12 July of 2012, dissolved CO2 was measured directly in the stream along the study reach. Measurements were made at 7 locations in June and 9 locations in July. Both dates show patterns of higher CO2 upstream (outside) and decreasing CO2 downstream (Fig. 7). Higher values of CO2 were observed in July. On both dates, the infeeding stream contained high values of CO2 that raise the values of CO2 downstream of the infeeder for several tens of meters. The steepest decreases in CO2 coincide with the steepest sections of the study reach, at a distance of roughly 50 m inside the cave. In July, CO2 measurements were taken at a higher spatial resolution, and a decrease of 500 ppm was observed in a short distance of only 5 m. Between these two measurement points the stream falls through two short (10–20 cm) cascades. 5.2.5. CO2 production experiments The sediment outside of the cave is characterized by sand and silt particles with significant particulate organic debris. Inside the cave, sediments are primarily gravels with some sands deposited among the gravel banks (Fig. 8). Finer sediments and particulate organic matter are also found in a few backwater locations in the cave, which make up only a few percent of the streambed surface area within the cave section of the study reach. In order to bracket the sediment types within the study reach, samples were taken from the outside, from typical gravel banks in the main stream channel in the cave, and from backwater deposits of fine sediment within the cave. Experiments with these sediment types showed that CO2 production was highest in the outside sediments, lowest in the main channel

Table 2 Field measured water quality parameters and results of water sample analysis. Location

CO2a (mmol kg−1)

SIa

Temp (°C)

SpC (μS cm−1)

pH

Alkalinity (as CaCO3 mg L−1)

Ca (mg L−1)

Mg (mg L−1)

8 June 2011 Surface (S1) Infeeder (S2) Downstream cave (S4)

14.7 8.6 14.1

131. 433. 163.

6.69 7.30 7.11

44.2 166.4 57.1

22.1 72.1 26.0

2.4 2.6 2.9

0.36 ± 0.075 0.40 ± 0.072 0.16 ± 0.043

−1.87 ± 0.10 −0.30 ± 0.08 −1.37 ± 0.14

6 July 2011 Surface (S1) Infeeder (S2) Upstream cave (S3) Downstream cave (S4)

18.3 8.7 10.8 9.9

397. 496. 444. 436.

7.31 7.58 7.73 7.80

193.6 216.1 201.1 198.1

71.8 89.3 76.6 78.1

5.1 4.3 6.1 5.3

0.37 ± 0.070 0.27 ± 0.047 0.17 ± 0.031 0.14 ± 0.027

−0.10 ± 0.09 0.17 ± 0.09 0.24 ± 0.09 0.30 ± 0.10

20 September 2011 Surface (S1) Infeeder (S2) Upstream cave (S3) Downstream cave (S4)

12.6 9.3 11.4 10.9

367. 582. 435. 432.

7.31 7.45 7.80 7.95

95.4 231.4 129.3 128.4

50.2 86.7 65.1 64.7

2.5 4.2 1.8 1.4

0.18 ± 0.041 0.40 ± 0.071 0.082 ± 0.020 0.057 ± 0.013

−0.67 ± 0.11 0.07 ± 0.08 0.02 ± 0.12 0.16 ± 0.13

18 April 2012 Surface (S1) Infeeder (S2) Upstream cave (S3) Downstream cave (S4)

8.4 8.5 8.1 8.1

175. 389. 217. 237.

7.33 7.34 7.42 7.55

64.2 157.3 80.9 88.2

27.3 65.4 33.9 37.3

1.5 1.9 1.9 1.6

0.11 ± 0.033 0.34 ± 0.060 0.13 ± 0.031 0.10 ± 0.025

−1.2 ± 0.17 −0.33 ± 0.09 −0.89 ± 0.14 −0.69 ± 0.14

13 June 2012 Surface (S1) Infeeder (S2) Upstream cave (S3) Downstream cave (S4)

14.3 8.9 12.9 12.6

162. 452. 199. 205.

7.23 7.50 7.46 7.50

69.3 186.7 83.1 85.0

28.6 79.6 34.6 35.9

1.3 0.9 1.5 1.7

0.15 ± 0.037 0.27 ± 0.047 0.10 ± 0.027 0.098 ± 0.025

−1.1 ± 0.13 −0.017 ± 0.09 −0.75 ± 0.14 −0.70 ± 0.14

a

Calculated using PHREEQC. Errors estimated using Monte Carlo sampling of input values for PHREEQC.

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

91

Fig. 6. Evolution of water chemistry along the studied reach: (A) alkalinity, (B) pH, (C) saturation index, and (D) dissolved CO2. The lines connect values along the main stream measured on the same day (8 June 2011=solid/squares; 6 July 2011=dash-dot/diamonds; 20 Sept 2011=dashed/circles; 18 Apr 2012=dotted/x's; 13 June 2012=solid/triangles). The symbols clustered at a stream distance of 70 m that are separated from the lines represent the values for the primary infeeding stream.

cave sediments, and intermediate in the backwater cave sediments. CO2 production occurred roughly at a constant rate in each sample over the one-hour experiment period following a one-hour equilibration period (Fig. 9). Production rates were normalized to sediment volumes to allow a more direct comparison from one sample to the

next despite the variation in sediment volume (see Section 5.1.6). Experiments were run on three samples from the outside channel, four samples from the in-cave backwater sediments, and three samples from the main channel in the cave. Average rates of production of CO2 per ml of sediment for the three sample types were 1.56 × 10 −10 mol s−1 ml−1 (outside sample), 2.40× 10−11 mol s −1 ml−1 (inside backwater sample), and 6.38× 10−12 mol s −1 ml−1 (inside main channel sample). The observed production rates displayed significant variability within a given sediment type, particularly for the outside samples. However, in all experiments sediment from the outside stream channel produced CO2 at a higher rate than sediment from inside the cave, and all backwater cave sediments produced CO2 at a higher rate than samples from the main channel in the cave (Fig. 10). 6. Discussion The field data collected allows discrimination between the possible mechanisms for the longitudinal decrease in dissolution rates, and the results of this analysis indicate processes that are of potentially broad importance in determining the rates of speleogenesis. We begin by examining each potential mechanism in detail and then discuss the broader implications of this work. 6.1. Discrimination between mechanisms for the decrease in dissolution rates

Fig. 7. Directly measured dissolved CO2 as a function of longitudinal distance along the stream on two dates. The square and circular symbols that are disconnected from the lines show the CO2 values measured in the infeeding stream in June and July, respectively. The dashed line depicts the stream gradient calculated from differencing of the surveyed stream profile. The sharpest drop in dissolved CO2 coincides with the steepest portion of the channel.

The decrease of dissolution rates because of dissolution of calcite along the flow path was rendered implausible as an explanatory mechanism using theoretical considerations alone. The typical length scales of the dissolution process are much longer than the observed

92

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

Fig. 10. Production rates for each sediment sample (triangles) divided along the x-axis into sediment types. Rates are normalized to sediment volume (ml). The squares connected by a solid line indicate mean production rates for each sediment type, and the gray region depicts the sample estimate of the population standard deviation for each sediment type.

Fig. 8. Typical settings in the surface reach (A) and cave reach (B). The surface is characterized by fine sediments and organic debris, whereas the cave channel contains a thin layer of sediment, primarily comprised of gravels and cobbles, overlying bedrock.

scale of 300 m. However, this conclusion is further supported by the water chemistry data. During all sampled flow conditions, the increases in alkalinity (dissolved load) along the studied reach are small (Fig. 6A) in comparison to the drop in dissolution rates of more than a factor of two (Fig. 4). The small observed increase in alkalinity could plausibly result from the input of diffuse flow with a higher dissolved load rather than dissolution of limestone along the reach. This conclusion can be

Fig. 9. Production curves from one experiment on each sediment type. The curves begin at the end of the equilibration period and run for 1 h. The concentration of CO2 at the end of the equilibration phase has been subtracted so that the slopes (rates) of the three experiments can be more easily compared. Production rates were estimated using the slope of the best fit line for each experiment.

tested more quantitatively by examining the date of the largest increase in alkalinity (20 Sept 2011). On this date, discharge was also measured and the diffuse component was determined to make up 34%. If one assumes that the diffuse input has a chemistry comparable to that measured in the primary infeeding stream then a simple mixing model with no additional dissolution would produce a downstream alkalinity of 141 mg l −1. This is actually slightly higher than the observed alkalinity but within the error of the measurement values and assumptions of the simple mixing model. Therefore, the data confirms that any increase in dissolved load because of dissolution along the study reach is minimal. The inflow of saturated water was deemed a more plausible explanation based on theory alone. However, the relatively small increases in dissolved load along the reach (Fig. 6A) in comparison to the large changes in dissolution rates suggest that this mechanism is not particularly important at the observed site. This conclusion is further supported when one considers the relatively small fraction of diffuse input observed (Table 1), even during low flow and recession conditions when diffuse input is expected to be maximal. The highest observed percentage of diffuse flow is 39 %, whereas the theory suggested a value of around 60% was needed in order to explain the observed decrease in dissolution rates. Furthermore, during low flow levels when the diffuse infeeder provides a larger percentage of the total flow, the stream is relatively close to saturation, or even slightly oversaturated (Fig. 6C), such that little dissolution occurs at any point along the stream. During floods, when the dissolution rates are highest, the diffuse input component is expected to be smallest. Taken together, the data support the conclusion that diffuse input of saturated water is not a strong factor in producing the downstream decrease in dissolution rates at the study site. The final proposed mechanism for the decrease in dissolution rates along the reach is the degassing of CO2. Under all observed conditions, the pH increases with distance along the stream (Fig. 6B). This, in combination with the small change in dissolved load, suggests that the dissolved CO2 is decreasing along the reach. The speciation calculations of PHREEQC display this directly, as for all sample dates the calculated concentration of CO2 decreases systematically with distance downstream. This was additionally confirmed via direct measurements of dissolved CO2 that showed systematic downstream decreases (Fig. 7). Therefore, CO2 degassing is responsible for the gradual increase of SI with stream distance (Fig. 6C) and stands out as the most likely mechanism for the decrease in dissolution rates observed using limestone tablets. 6.2. Potential causes of the decrease in dissolved CO2 As noted above, the equilibrium concentration of dissolved CO2 is set by three factors: (i) the pCO2 in the air above the water, (ii) the

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

rate of production of CO2 in the water, and (iii) the reaeration rate. The drop in dissolved CO2 could be triggered by a sudden shift in one or more of these parameters. The pCO2 of cave air is typically somewhat higher than that on the surface because of high levels within the soil zone (e.g., Troester and White, 1984; Ek and Gewelt, 1985; Baldini et al., 2006). In fact, CO2 values measured within the cave air showed values near atmospheric during the winter and up to several times atmospheric during the summer (Prelovšek, 2009). Therefore, a change in the pCO2 of the air above the stream cannot explain the reduction in CO2 within the cave stream. However, either of the other mechanisms could be active. At equilibrium, the production rate of CO2 in the sediment and water (RCO2 ) is equal to the flux of CO2 from the water to the air given by Eq. (4). Therefore,   RCO2 ¼ k C w −α CO2 C a :

ð6Þ

In order to make estimates of the change in gas transfer coefficient, k, or production rate, RCO2 , that would be required to produce the observed downstream change in pCO2, we can assume that the stream is near equilibrium in the furthest upstream and downstream reaches. This is reasonable since pH and conductivity were found to be relatively constant along the surface reach and CO2 is changing relatively slowly by the time the stream reaches the downstream end of the study reach. With this assumption, and assuming that the concentration of CO2 in the air is similar in the cave and outside, RCO2 ;out kin RCO2 ;in kout

¼

C w;out −α CO2 C a C w;in −α CO2 C a

:

ð7Þ

Using the date when observed dissolved CO2 changes the most (12 July 2011), the ratio on the right hand side of Eq. (7) is approximately equal to 4. Therefore, on that date, the production rate or gas transfer coefficient must change sufficiently downstream such that RCO2 ;out kin RCO2 ;in kout

≈4:

ð8Þ

If a change in the production rate or gas transfer coefficient alone was responsible then that parameter would need to change by a factor of roughly four. The CO2 production rate may significantly decrease within the cave. The microbial production rate of CO2 is significantly correlated with the percent of fine sediment and organic debris (Baker, 1986; Hedin, 1990) and higher temperatures also typically result in higher production rates (Naiman, 1983; Bott and Kaplan, 1985). The surface stream is surrounded by a larger soil zone, where CO2 is being produced; and the surface channel itself contains abundant fine sediment and organic debris (Fig. 8A). On the other hand, the cave channel is surrounded by a relatively impermeable bedrock matrix, and the channel contains primarily gravels and cobbles with only limited patches of fine sediment and organic debris (Fig. 8B). The possibility of lower CO2 production rates in the cave is further evidenced by the CO2 production experiments, which showed that typical cave sediments produced CO2 at a rate significantly less than that seen in outside sediments. Many factors present within the natural setting, such as differences in temperature and the structure of the sediment deposits, were not accounted for in the experiments. Consequently, the experimental rates cannot be used to quantitatively predict the production rates at the field site. However, the experiments do provide strong evidence that a contrast in production rates exists at the site and that this contrast explains some of the observed downstream reduction in dissolved CO2. A second plausible explanation for the decrease in dissolved CO2 is an increase in the gas transfer velocity, k. The value of k is largely

93

controlled by the hydrodynamics within the top layer of the water (Jähne and Haußecker, 1998; Wanninkhof et al., 2009), and streams with shallower and more turbulent flow typically have higher k (Wanninkhof et al., 1990). Upon entry into the cave, the stream significantly steepens, such that flow is both shallower and more turbulent (Fig. 5). In fact, this increase in steepness is largely responsible for the change in average sediment sizes noted above. Given the increase in steepness, the gas transfer velocity inside the entrance zone of the cave must be higher than outside. We use an estimation of the scaling of k with channel steepness to determine whether observed changes in slope could explain the decrease in dissolved CO2. Melching and Flores (1999) (Eq. (10)) determined an approximate relationship for reaeration rates in small pool and riffle streams as a function of hydraulic parameters. For streams of the same depth, reaeration rate is proportional to the gas transfer coefficient. Since stream depths at the far upstream and downstream ends of the reach are comparable, we can use the scaling relation given by Melching and Flores (1999) and 0:524

k∝ðVSÞ

Q

−0:242

;

ð9Þ

where V is the flow velocity, S is the stream slope, and Q is the discharge. The increase in discharge is small enough, and the dependence weak enough, that discharge effects k by only about 10 %, and can be neglected in an order of magnitude estimation. Using the Darcy–Weisbach equation, V ∝ S 0.5, giving an approximate relation for the scaling of k as a function of slope, k∝S

0:75

:

ð10Þ

Using this relation, slope would have to increase by a factor of roughly six in order to explain the observed decrease in dissolved CO2 (Eq. (8)). On the contrary, the slope at the downstream end is only slightly steeper than at the upstream end (Fig. 7). Consequently, changes in k do not provide as plausible of an explanation for the observed differences in equilibrium CO2 concentrations at the upstream and downstream ends of the study reach. However, a short distance downstream of the cave entrance, the stream gradient reaches values nearly 10 times the slopes at the far ends of the study reach. Therefore, it is likely that CO2 dynamics within the middle portion of the study reach are influenced by significantly higher values of k. This is further evidenced in that directly measured values of stream CO2 drop steeply as the stream flows through the steepest portion of the study reach at the location of a series of small cascades (Fig. 7). Therefore, within the study reach, both proposed mechanisms for the reduction in CO2 appear to act, only at different scales. A short-term increase in k near the cave entrance causes a relatively rapid decrease in dissolved CO2 downstream, but, ultimately, a reduced production rate of CO2 sets the new downstream equilibrium value. 6.3. Feedback mechanisms Theoretical considerations in combination with field data demonstrate that the observed longitudinal decrease in dissolution rates primarily results from the downstream decrease in CO2. Two plausible mechanisms can trigger this change, and both appear to play some role at the study site. Additionally, the considerations in Section 6.2 suggest that both of these mechanisms are likely to be active under conditions that are relatively common in karst streams that contain significant changes in steepness and have a supply of organic material or other source of CO2. These two mechanisms introduce geomorphic feedback loops between channel profile, CO2 dynamics, and incision rates (Fig. 11). Dissolution processes influence the stream profile. The stream profile in turn affects the CO2 dynamics via the sediment and/or hydraulic characteristics of the channel. The CO2 dynamics feed back into the

94

M.D. Covington et al. / Geomorphology 186 (2013) 85–95

dissolution processes and stream profile evolution. While the geomorphological implications of such feedback loops over long timescales are not clear without a more solid understanding of the CO2 dynamics, the existence of such feedbacks is well-motivated by both theoretical considerations and the observations at the study site. The presence of similar feedback loops that drive the precipitation rates of calcite have been shown to drive the morphology of calcite speleothems, such as stalactites and rimstone dams (Short et al., 2005; Goldenfeld et al., 2006; Chan and Goldenfeld, 2007; Hammer et al., 2007; Meakin and Jamtveit, 2010). Previous work has demonstrated that production of CO2 within an incipient cave system can significantly influence the timescales and patterns of evolution (Gabrovšek et al., 2000). However, no speleogenesis or stream incision models have yet considered coupling of production rates with channel characteristics or degassing in open channel settings. Focused field and modeling studies may elucidate the implications of these processes and the resulting feedback loops. 7. Conclusions The studied allogenically fed cave stream exhibited a strong downstream decrease in dissolution rates over the observed 2-year period. We proposed a number of potential mechanisms to explain this decrease and used theoretical models to examine the plausibility of each mechanism. An analysis of the characteristic length scales of each process and comparison with the observed length scale provides a straightforward test of each model. Mechanisms were further tested using a variety of field data, the results of which suggest that a decrease in dissolved CO2 is primarily responsible for the decrease in dissolution rates. Two mechanisms for the reduction in CO2 are proposed: (i) a reduction in CO2 production due to a change in sediment characteristics and microbial activity, and (ii) an increase in the gas transfer velocity resulting from steepening of the channel. Field data suggest that both of these mechanisms play a role in the CO2 dynamics at the study site. An increase in the gas transfer coefficient in the steepest section of the study reach coincides with a sudden decrease in CO2, and the difference in CO2 production rates in sediments from the surface and cave can explain the lower equilibrium value of CO2 at the downstream end of the cave. Both mechanisms introduce a feedback loop between channel profile and CO2 dynamics, via sediment transport and biological activity, or via stream hydraulics alone. Current numerical models of speleogenesis do not consider any of these processes, and further

Fig. 11. Feedback loops introduced because of the coupling of stream profiles and CO2 dynamics. Arrows indicate the direction of causality. Field data indicates that stream profile indirectly influences both CO2 production rate and gas transfer velocity at the study site. These connections introduce feedback loops in the incision process. The feedback loops should occur commonly in channels incising soluble strata, particularly near significant changes in steepness and when organic material is present in the stream.

modeling and field work are warranted to better constrain the controls on CO2 in a number of environments and the potential long-term geomorphic implications. Acknowledgments MDC is supported by a National Science Foundation International Research Fellowship (OISE-0754495). Any opinions, findings, and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the NSF. References Anthony, D., Groves, C.G., Meiman, J., 2003. Preliminary investigations of seasonal changes in the geochemical evolution of the Logdson River, Mammoth Cave, Kentucky. Speleogenesis and Evolution of Karst Aquifers 1, 1–6. Baker, J., 1986. Relationship between microbial activity of stream sediments, determined by three different methods, and abiotic variables. Microbial Ecology 12, 193–203. Baldini, J., Baldini, L., McDermott, F., Clipson, N., 2006. Carbon dioxide sources, sinks, and spatial variability in shallow temperate zone caves: evidence from Ballynamintra Cave, Ireland. Journal of Cave and Karst Studies 68, 4–11. Birk, S., Liedl, R., Sauter, M., Teutsch, G., 2005. Simulation of the development of gypsum maze caves. Environmental Geology 48, 296–306. Bott, T., Kaplan, L., 1985. Bacterial biomass, metabolic state, and activity in stream sediments: relation to environmental variables and multiple assay comparisons. Applied and Environmental Microbiology 50, 508–522. Chan, P., Goldenfeld, N., 2007. Steady states and linear stability analysis of precipitation pattern formation at geothermal hot springs. Physical Review E 76, 1–11. Covington, M.D., Luhmann, A.J., Wicks, C.M., Saar, M.O., 2012. Process length scales and longitudinal damping in karst conduits. Journal of Geophysical Research 117, 1–19. Dreybrodt, W., 1988. Processes in Karst Systems: Physics, Chemistry, and Geology. Springer, New York, USA. Dreybrodt, W., 1990. The role of dissolution kinetics in the development of karst aquifers in limestone: a model simulation of karst evolution. Journal of Geology 98, 639–655. Dreybrodt, W., 1996. Principles of early development of karst conduits under natural and man-made conditions revealed by mathematical analysis of numerical models. Water Resources Research 32, 2923. Dreybrodt, W., Gabrovsek, F., Romanov, D., 2005. Processes of Speleogenesis: A Modeling Approach. ZRC Publishing, Ljubljana, Slovenia. Ek, C., Gewelt, M., 1985. Carbon dioxide in cave atmospheres. New results in Belgium and comparison with some other countries. Earth Surface Processes and Landforms 10, 173–187. Farrant, A.R., Smart, P.L., 2011. Role of sediment in speleogenesis; sedimentation and paragenesis. Geomorphology 134, 79–93. Ford, D., Ewers, R., 1978. The development of limestone cave systems in the dimensions of length and depth. Canadian Journal Earth Science 15, 1783–1798. Gabrovšek, F., Dreybrodt, W., 2000. Role of mixing corrosion in calcite-aggressive H2O– CO2–CaCO3 solutions in the early evolution of karst aquifers in limestone. Water Resources Research 36, 1179–1188. Gabrovšek, F., Menne, B., Dreybrodt, W., 2000. A model of early evolution of karst conduits affected by subterranean CO2 sources. Environmental Geology 39, 531–543. Gams, I., 1981. International comparative measurement of surface solution by means of standard limestone tablets. Proceedings of the 8th International Congress of Speleology, Bowling Green, KY. National Speleological Society, Huntsville, AL, USA, pp. 273–275. Goldenfeld, N., Chan, P., Veysey, J., 2006. Dynamics of precipitation pattern formation at geothermal hot springs. Physical Review Letters 96, 1–4. Groves, C.G., Howard, A.D., 1994. Minimum hydrochemical conditions allowing limestone cave development. Water Resources Research 30, 607. Groves, C.G., Meiman, J., 2005. Weathering, geomorphic work, and karst landscape evolution in the Cave City groundwater basin, Mammoth Cave, Kentucky. Geomorphology 67, 115–126. Hammer, O., Dysthe, D., Jamtveit, B., 2007. The dynamics of travertine dams. Earth and Planetary Science Letters 256, 258–263. Hedin, L., 1990. Factors controlling sediment community respiration in woodland stream ecosystems. Oikos 57, 94–105. Jähne, B., Haußecker, H., 1998. Air-water gas exchange. Annual Review of Fluid Mechanics 30, 443–468. Johnson, M., Billett, M., Dinsmore, K., Wallin, M., Dyson, K., Jassal, R., 2010. Direct and continuous measurement of dissolved carbon dioxide in freshwater aquatic systems — method and applications. Ecohydrology 3, 68–78. Kaufmann, G., 2009. Modelling karst geomorphology on different time scales. Geomorphology 106, 62–77. Liu, Z., Dreybrodt, W., 1997. Dissolution kinetics of calcium carbonate minerals in H20−CO2 solutions in turbulent flow: the role of the diffusion boundary layer and the slow reaction H2O+CO2 ⇌H+ +HCO3−. Geochimica et Cosmochimica Acta 61, 2879–2889. Meakin, P., Jamtveit, B., 2010. Geological pattern formation by growth and dissolution in aqueous systems. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 466, 659–694.

M.D. Covington et al. / Geomorphology 186 (2013) 85–95 Melching, C., Flores, H., 1999. Reaeration equations derived from US Geological Survey database. Journal of Environmental Engineering 125, 407–414. Naiman, R., 1983. The annual pattern and spatial distribution of aquatic oxygen metabolism in boreal forest watersheds. Ecological Monographs 53, 73–94. Palmer, A.N., 1991. Origin and morphology of limestone caves. Bulletin of the Geological Society of America 103, 1. Palmer, A.N., 2007. Variation in rates of karst processes. Acta Carsologica 36, 15–24. Parkhurst, D.L., Appelo, C.A.J., 1999. User's guide to PHREEQC (Version 2) — a computer program for speciation, batch-reaction, one-dimensional transport, and inverse geochemical calculations. Water-Resources Investigations Report 99-4259, p. 326. Plan, L., 2005. Factors controlling carbonate dissolution rates quantified in a field test in the Austrian Alps. Geomorphology 68, 201–212. Prelovšek, M., 2009. Present-day speleogenetic processes, factors and features in the epiphreatic zone. Ph.D. thesis. University of Nova Gorica. Nova Gorica, Slovenia.

95

Short, M., Baygents, J., Beck, J., Stone, D., Toomey, R., Goldstein, R., 2005. Stalactite growth as a free-boundary problem: a geometric law and its platonic ideal. Physical Review Letters 94, 14–17. Spate, A., Jennings, J., Smith, D., 1985. The microerosion meter: use and limitations. Earth Surface Processes and Landforms 10, 427–440. Troester, J., White, W.B., 1984. Seasonal fluctuations in the carbon dioxide partial pressure in a cave atmosphere. Water Resources Research 20, 153–156. Wanninkhof, R., Mulholland, P., Elwood, J., 1990. Gas exchange rates for a first-order stream determined with deliberate and natural tracers. Water Resources Research 26, 1621–1630. Wanninkhof, R., Asher, W.E., Ho, D.T., Sweeney, C., McGillis, W.R., 2009. Advances in quantifying air-sea gas exchange and environmental forcing. Annual Review of Marine Science 1, 213–244.