COBISS: 1.01
ON DENUDATION RATES IN KARST O HITROSTI DENUDACIJE NA KRASU
Franci GABROVŠEK1
Abstract                                           UDC 551.331.24:551.44
Franci Gabrovšek: On denudation rates in Karst
Paper presents a simple mathematical model, which enables study of denudation rates in karst. A vertical fow of water which is uniformly infltrated at the surface is assumed. Denudation rate is calculated from the time needed to remove certain thickness of rock. Tis is done concretely on a limestone block dissected by a vertical array of fractures. It is shown that denudation rate increases with the thickness of removed layer and approaches an upper limit which is defned by the maximum denudation equations, which are based on assumption that all dissolution potential is projected into a surface lowering. Keywords: karst, denudation rate, limestone dissolution, math-ematical model.
Izvleček                                            UDK 551.331.24:551.44
Franci Gabrovšek: O hitrosti denudacije na Krasu
V prispevku predstavim enostaven matematični model s katerim raziskujem dinamiko zniževanja kraškega površja. Predpostavim enakomerno napajanje s površja in vertikalno pronicanje vode. Denudacijsko stopnjo izračunam iz časa, ki je potreben za odstranitev določene debeline kamninskega sloja. Konkretno to naredim na primeru apnenca v katerem se voda pretaka v sistemu vertikalnih razpok. Hitrost denudacije narašča z debelino odstranjene plasti in doseže zgornjo mejo, ki je določena z enačbami, ki temeljijo na predpostavki, da se celoten korozivni potencial vode manifestira v zniževanju površja. Ključne besede: kras, denudacijska stopnja, raztapljanje apnenca, matematični model.
INTRODUCTION
Uniform lowering or surface denudation is a dominant karstifcation process (Dreybrodt, 1988; Ford & williams, 1989; white, 1988). Te denudation rate is defned as the rate (LT-1) of lowering of a karst surface due to the dissolution of bedrock. A common approach used to estimate the denudation rate is based on the presumed equilibrium concentration (or hardness) and the amount of water which infltrates into the subsurface. It is sum-marized in the famous Corbel’s equation (Corbel, 1959):
Dc(m/Ma) = (P   E)H f c                 1000-p
1
Te infltrated water in mm/y is the diference between precipitation P and evapotranspiration E. H is the equilib-
rium concentration (Hardness) in mg/L of dissolved rock, ? is the density of limestone in g/cm3, f denotes the portion of soluble mineral in the rock, which will be 1 in this paper. Te factor 1000 corrects for the mixture of units used in the equation.
Tere are more general equations of this kind like that of white (1984, this issue). For a Limestone terrain in a temperate climate all these equations give denudation rate of the order of several tens of meters per million years. Similar results are obtained from fow and con-centration measurements in rivers which drain a known catchment area. From the measured data the total rock volume removed from the area in a given time period can be calculated. Dividing the removed volume by the surface of the area and the time interval gives the denudation rate.
1 Karst Research Institute ZRC SAZU, Postojna, Slovenia, e-mail: gabrovsek@zrc-sazu.si Received/Prejeto: 01.02.2007
TIME in KARST, POSTOJNA 2007, 7–13
FRANCI GABROVŠEK
Eq.1 implies that all dissolution capacity of water is used in the rock column, i.e. the solution at the exit of rock block is close to saturation. Among the many assumptions behind such estimations of the denudation rate I will address two of which at least one must be valid:
1.  Most of the dissolution occurs close to the sur-face, i.e. within epikarst.
2.  In the long term, the dissolution at depth is inte-grated into a surface denudation.
It is the intention of this paper to theoretically vali-date “maximum denudation” approach.
SURFACE LOwERING AND THE VOLUME OF DISSOLVED ROCK
Dissolution of any rock is not instantaneous, but proceeds at some fnite rates. In conditions of difuse infltration through the karst surface and prevailing vertical fow, the concentration of dissolved rock in the infltrating water will normally increase with the depth as schematically shown by color intensity in Fig.1.
Fig. 1. Section of a terrain with a uniform surface infltration of aggressive solution and prevailing vertical fow. Color intensity denotes that the concentration of dissolved rock increases with depth.
Fig. 2 presents point at some depth z below the sur-face. Te volume ?V of rock dissolved per unit surface area S in time ?t between the surface and the point is given by
Fig. 2: Idealized profle through the rock column at time t = 0 (lef) and t > 0 (right). Te depth of the point which is at z decreases in time due to the surface lowering.
äV/S = c(z)-q-M I p
2
depth at t = 0 and D is the denudation rate (Fig. 2). Te volume of dissolved rock per surface area in time T above the point is then given by:
AV/S (T) = ^ fc(z(t)) dt = - fc(z0 - D • t) dt            3
Po                   Po
Introducing a new variable z=z0 - D·t into the right hand integral in Eq. 3 gives:
where c(z) is the concentration of dissolved rock [M/L3] at the depth z, q is the infltration rate at the surface [L3/ (L2T)] and ? is the density of the rock [M/L3].
Due to the surface lowering, the depth of the point is decreasing according to z(t) = z0 – D·t, where z0 is the
8     TIME in KARST – 2007
AV
S             D-PZ„JDR.T
4
ON DENUDATION RATES IN KARST
Te complete volume of rock initially above z0 is z0·S. To remove this volume a time TD is needed, where z0
D·TD. Using all this in Eq. 4, we obtain:
plete layer and D an avarege denudation rate. It is easy to see that if the solution quickly attains equilibrium Eq. 6 gives maximum denudation rates:
Sfc(z)dz
q-TD z0-p
fc(z)dz
a'C    %        a -c
D_^_3fdz-^—S-Dc VPo           P
D
VPo
fc(z)dz
As given, D is an average denudation rate, cal-culated from the time TD needed to remove a layer of thickness z0 from the rock column with initial a uniform porosity distribution in vertical direction. If a rock layer has a fnite thickness, z0 can be taken as the layer thickness, TD the time needed to remove the com-
If this is not the case D will be below Dc, since inte-
gral with ceq is maximal. In this case we
rewrite
D = ^-lfc(z)dz = fcJDc
C          7-J                                    P.
Eq. 6 as: 8
with increasing layer thickness an average concen-tration within the layer increases and average denudation rates approach maximal.
CALCULATION OF THE CONCENTRATION PROFILE
Te results given so far are valid for any “natural” c(z). To obtain some quantitative results we revert to a special case where the calcite aggressive water is infltrating into a vertical fracture network. Terefore we need to couple the rate equation for limestone and fow of laminar flm down a vertical fracture wall.
LIMESTONE DISSOLUTION RATES Recently Kaufmann & Dreybrodt (2007) published the corrected rate equation with two linear regions and a non-linear region of dissolution kinetics:
'«i(0.3c„-c) c<Q3ceq
R = ' a2(ceq -c) 0.3ceq<c< 0.9ceq           9
0(c„-c)"   c>0.9c
, *   v  eq        s                         eq
Te kinetic constants and rate orders are derived from theoretical and experimental results (Buhmann & Dreybrodt, 1985; Dreybrodt, 1988; Eisenlohr et al, 1999; Kaufmann & Dreybrodt, 2007). Values depend on the temperature, pCO2 and laminar layer thickness and are given in Kaufmann & Dreybrodt (2007). we will use ?1 = 3·104 cm/s and ?2 = 8·10 6 cm/s, values which are valid at 10°C for the open system dissolution (Kaufmann & Dreybrodt, 2007). Nonlinear kinetics will not be discussed here. It is valid close to equilibrium and does not change the results substantially. c   depends on the pCO2, temper-
ature, the presence of the foreign ions and the nature of the system where dissolution proceeds (open, closed, in-termediate) (Appelo & Postma, 1993; Dreybrodt, 1988). Te calcium equilibrium concentration normally takes values between 0.5 mmol/l – 3 mmol/l, which in terms of dissolved calcite means 50—300 mg/l .
LAMINAR FLOw DOwN A SMOOTH VERTICAL wALL Only rough assumptions can be made about the fow regime of infltrated water. we assume a laminar flm fow down the walls of vertical fractures. Te velocity of such flm is given by (Bird et al., 2002):
3v
,
10
where ? is the flm thickness, g gravitational acceleration and ? the kinematic viscosity. More suitable master variable is a fow density along the fracture walls qf (cm2/s). Applying qf = vö we get:
3qfu
and v
*"(L)'
11
TIME in KARST – 2007
5
7
6
9
FRANCI GABROVŠEK
THE COUPLING OF FLOw AND DISSOLUTION To get the evolution of a concentration in a falling flm, the mass balance for dissolved calcite must be coupled to the rate laws given in Eq. 9. In a small volume of water flm with thickness ?, width b, length dz and concentration c, the mass balance requires:
Aax(c   -c)dt=Vdc=Aödc;c<03ceq Aa2(ceq- c)dt=Vdc=A8dc; 0.3 cn< c < 0.9^
12
where A is the surface of the water rock contact (b·dz). Integration of Eq. 12 gives the evolution of concentration in time as the flm fows down the fracture wall:
c(t)
U"'Ceqli     e        J>C<«"'C«1                     13
c n-0.7e-^);0.3c   <c<0.9c
where  xi = 8/cij . Note that ?1 is more than an order of magnitude smaller than ?2. The time domain can easily be converted into the space domain using z = v ¦ t and qf =v-8
Fig.3: Concentration profle in a flm fowing down a smooth vertical fracture and dissolving limestone walls. values of ?2 are calculated from the fracture fow density obtained if the fracture spacing is 1 m and infltration intensity is 0.114 mm/h, 10 mm/h, 20 mm/h, 30 mm/h and 40 mm/h for curves 1-5 respectively.
c(z)
[ceq(l-0.7e-^);0.3ceq<c<0.9ceq
where X{ =qf/<Xj . Fig. 3 presents the evolution of saturation ratio c(z)/c   for diferent ?2. For most reasonable
 eq
scenarios, the frst linear kinetics is active only close to the surface. Terefore, it will be integrated directly into the surface lowering (i.e. c(z= 0) = 0.3 c ).
 eq
SATURATION LENGTH ? AND THE FRACTURE FLOw DENSITy
Saturation length ?2 controls the vertical evolution of concentration profle. It depends on the kinetic constant and the fracture fow density. To estimate the latter we assume that the rain falling to the surface with an inten-sity q is evenly infltrated into a regular grid of fractures as shown on Fig. 4. Te fow density in each fracture is proportional to the ratio between the surface of the infl-tration area and the total breadth of the fractures drain-ing the area. In a regular grid of fractures with fracture spacing d we obtain:
qf = N-q-d ,
15
Fig. 4: Rain falling with intensity q [Lt-1] is uniformly distributed into the fractures with fow density qf [Lt-1] according to Eq.15.
10     TIME in KARST – 2007
ON DENUDATION RATES IN KARST
where N depends on the geometry of fracture grid (e.g. 1 for a series of parallel fractures and 1 for a square grid). For thin flms (? < 0.005 cm) the rates are controlled by conversion of CO2 into H+ and HCO3- and therefore increase linearly with flm thickness. Infltration intensity below 3 mm/h into a series of fractures with d = 100 cm produces flm thicknesses below 0.005 cm. Extremely low infltration intensities (< 1 mm/h) and fracture spacing in the range of few centimeters result in flm thicknesses in the order of 0.001 cm. Tis reduces the kinetic con-
stant ? approximately by a factor of 5. For thin flms, the flm thickness and ? decrease with qf1/3. Consequently ?2 is proportional to the qf2/3. Tis has limited consequences for the dissolved volume and denudation rate discussed in the next section.
In the early stages, the fracture fow is not expected to be in the form of a free surface flm, but full fracture fow instead. Te saturation lengths in that case would be smaller than those derived here. Te evolution of such fractures is given in Dreybrodt et al. (2005).
RESULTS AND DISCUSSION
Inserting the concentration profle from Eq. 14 (second linear region only) into Eq. 7 gives:
D = ^^f(l-0.7e-z'^) dz = P   Jo
P                P      Jo
= Dc-0.7-Dc-^-(l-e-z°"12)
z„
16
Although it is a matter of a defnition, the average denudation rate given in Eqs. 8 and 16 are not exactly what we are afer. what we look for is the actual lowering of karst surface, which is given by dz0 / dtd. we will not go into mathematical details of derivation, but instead discuss its consequences on a plot of z0(TD). Note that the z0(TD) has no explicit form, but its inverse function does:
TD(z0) =
Dc-0.7¦T>c-y^{\-s-1^1)
17
we will demonstrate the results on a characteristic data for a moderate climate with I=1000 mm/y and rela-tively bare karst area with ceq= 1mmol/l or H = 100 mg/L. For ? = 2.5 g/cm3. DC for this case is 40 m/Ma. we assume that the rain infltrates into a parallel set of fractures with spacing d = 1 m.
Fig. 5a shows z0(TD) for four different saturation lengths arising from different infiltration intensities. yearly infiltration is 1000 mm/y for all curves. Therefore, the time period of dissolution is inversely proportional to the infiltration intensity. Dashed line shows the uniform lowering by DC. we wee that all lines become practically parallel to maximum denudation line for zg > 2?2.The actual denudation rate becomes “maximal 0 when the removed thickness is larger than 2?2. This is about the depth where the concentration reaches 90% of saturation. The slope of the dotted lines presents the averaged denudation rates for curve with ?2 = 70 m.
Fig. 5b shows the averaged rate for the same scenarios as Fig. 5a. Red dashed red curve clearly shows the fast approach of the actual rate to maximal for ?2 = 52.5 m.
Another interesting conclusion can be made from Fig. 5a. Diferent saturation lengths ? can also arise from diferent fracture spacing (see Eq. 15 for qf) . If we imagine a region with high fracture density within a region of low fracture density, the frst will initially be denuded faster, but latter on both actual rates will become the same. Terefore the diference made at the onset will stay projected in the surface. Tis is shown by the double arrow between lines 3 and 4.
TIME in KARST – 2007
11
FRANCI GABROVŠEK
Fig. 5: a) Te time dependence of removed thickness for several infltration intensities. I=1000 mm/y, h = 100 mg/L, ? = 2.5 g/cm3, d = 100 cm, N = 2. dashed line show the “maximum denudation” rate which is 40 m/ma. dotted lines present the time averaged denudation rates (Eq.16). double arrow demonstrates the diference between the denuded thicknesses which is kept in time due to the initial rate diferences.
b) dependence of average denudation rates on the removed thickness for the same scenarios as in Fig. 5a. dashed line presents the actual surface lowering for ?2 = 52.5 m.
12     TIME in KARST – 2007
ON DENUDATION RATES IN KARST
CONCLUSION
Denudation rate in a block with initially uniform porosity increases as the denudation proceeds and becomes maximum denudation (Eq.1), when the thickness of removed layer is about twice the typical saturation length. Initial diferences arising from diferent saturation lengths remain imprinted in the surface.
If a soluble layer has a fnite thickness, the average denudation rate increases with the thickness, i.e. denudation is more efective on thick rock layers.
Te presented results are based on many assump-tions which might not be valid. Nevertheless, it gives some theoretical validation of maximum denudation formulae and suggest some mechanisms that can cause irregularities in karst surface.
REFERENCES
Appelo, C. A. J. & D. Postma, 1993: Geochemistry, ground-water and pollution. A.A. Balkema, xvi, 536 pp, Rotterdam; Brookfeld, VT.
Bird, R. B., Stewart, w. E. & E.N. Lightfoot, 2002: transport phenomena. John wiley & Sons, Inc., xii, 895 p. pp, New york, Chichester.
Buhmann, D. & w. Dreybrodt, 1985: Te kinetics of cal-cite dissolution and precipitation in geologically relevant situations of karst areas.1. Open system.-Chemical geology, 48, 189-211.
Corbel, J., 1959: Vitesse de l’erosion.- Zeitschrif fur Geomorphologie, 3, 1-2.
Dreybrodt, w. , Gabrovšek, F. & D. Romanov, 2005: Pro-cesses of speleogenesis: A modeling approach. Vol. 4, Carsologica, Založba ZRC, 375 pp, Ljubljana.
Dreybrodt, w. , 1988: Processes in karst systems: physics, chemistry, and geology. Springer-Verlag, xii, 288 p. pp, Berlin; New york.
Eisenlohr, L., Meteva, K., Gabrovšek, F. & w. Dreybrodt, 1999: Te inhibiting action of intrinsic impurities in natural calcium carbonate minerals to their dissolu-tion kinetics in aqueous H2O-CO2 solutions.- Geo-chimica Et Cosmochimica Acta, 63, 989-1001.
Ford, D.C. & P. williams, 1989: Karst geomorphology and hydrology. Unwin Hyman, 601 pp, London.
Kaufmann, G. & w. Dreybrodt, 2007: Calcite dissolutio n kinetics in the system CaCO3-H2O-CaCO3 at high undersaturation.- Geochimica Et Cosmochimica Acta, In Press.
white, w.B., 1984: Rate processes: chemical kinetics and karst landform development. In: La Fleur (Ed.): Groundwater as a geomorphic agent. Allen and Un-win, 227-248.
white, w. B., 1988: Geomorphology and hydrology of karst terrains. Oxford University Press, ix, 464 pp, New york.
TIME in KARST – 2007
13