Assessment of Specific Yield in Karstified Fractured Rock through the Water-Budget Method

Abstract

In this note, the Water Budget Method (WBM) is applied to estimate local values of the specific yield of the deep karst aquifer of Salento peninsula. A selection in a period of two years of relevant short precipitation events has been considered and the related localized recharges have been compared to the water table fluctuations measured at two selected wells. The recharge amounts have been corrected by using data of evapotranspiration and soil water storage available from a micrometeorological base. The results are very similar for both the wells and more consistent when the corrections are applied. A discussion involving frequency and apertures of the fractures in the rock mass of the aquifer suggests the effect of the karst dissolution to be dominant in determining these values of the specific yield.

1. Introduction

The assessment of hydraulic parameters of aquifers in karstified fractured rocks is always fraught with problems deriving from the inconsistency of well test data and the whimsical response of the rock formation. Among these parameters the specific yield S_y of an unconfined aquifer is of importance, especially in order to determine water availability [1,2]. It is defined as the volume of water released from storage by the aquifer per unit surface area and per unit decline of the water table [3].

In Salento peninsula (Apulian region, Southern Italy), the main groundwater reservoir is a deep aquifer, generally unconfined, consisting of Cretaceous limestones. When the limestones outcrop, meteoric water rapidly infiltrates from the outcrops into the aquifer for the high permeability of the vadose zone, apart from the volume of water collected into the surface soil [4]. Sometimes, upper and relatively thick layers of calcarenites may interfere, but the infiltration may be still rapid in case of favorable geological conditions. Apart from the moderately thick layer of soil, there is no retention in the vadose zone. Furthermore, cretaceous limestones do not retain water, therefore only fractures and karstic voids contribute to the storage when the water table raises after a rainfall. The occurrence of these features in limestones is random, and so is the value of the specific yield, therefore any direct measurement of this parameter strictly refers to the limited area of influence of the selected test. The simplest and most reliable method for this measurement is the Water-Budget Method (WBM), based on the equation obtained by combining the link between the change in groundwater storage ΔS (or recharge) and the water table fluctuation Δh, i.e. the equation ΔS = S_yΔh, and the aquifer budget Equation [5]. Given the specific features of the Salento Cretaceous aquifer, the budget equation may be written as follows:

ΔS = P − E − (θ_f − θ_i) t_s

(1)

where P is precipitation, E is evapotranspiration, θ_i and θ_f are respectively initial and final volumetric depth-averaged moisture content of the soil and t_s is the thickness of the soil. Note that in the last equation the contribution of irrigation or additional sources, as well as the withdrawal from pumping, is neglected. In addition, given the depth of the water table, the evaporation from the water table is negligible.

It is worth to note that, in soils or porous rocks, by using WBM one measures the fillable porosity rather than the effective porosity, coincident with S_y. The fillable porosity specifically represents the amount of water that an unconfined aquifer can store per unit rise in water table and per unit area [6], and can be smaller than the specific yield for the occurrence of air trapped in the pores below the water table when it quickly raises. In fractured rocks the chance of having trapped air is low, therefore in this context there is no ambiguity in the measurement of S_y.

In what follows local estimations of the specific yield of the deep aquifer of Salento peninsula by using WBM are reported and the steps of the method described. Data of water table fluctuations are extracted from a five-years governmental campaign of aquifer monitoring (Tiziano project), whereas data concerning the natural recharge derive from a micrometeorological base of CNR-ISAC (Consiglio Nazionale delle Ricerche, Istituto di Scienze dell’Atmosfera e del Clima, Italy) and rain gauges of the Civil Protection Agency. The value of S_y is then compared to semi-quantitative estimations of the porosity of the Cretaceous limestones by using geological and structural data.

In Salento, in the context of the 2007–2011 Tiziano project, a monitoring network was set up to assess both quantity and quality of the groundwater resource. Raw data are available online (http://www.sit.puglia.it/portal/portale_cis). This project was the first large management action of the regional authority to perform a general check on the availability of the local water resources. Such an issue started from a decennal decrease of the water reserves stored into huge karst aquifers, including the one located at the Salento peninsula [4], for which the depletion of the water availability is of a particular concern because it provides 75% of the local water supply and is subjected to saltwater intrusion. The monitoring network of the Tiziano project would allow the qualitative classification of groundwater bodies, the evaluation over time of the evolutionary dynamics of the qualifying state and the assistance in the formulation of the necessary actions to protect and safeguard the resources in the context of the regional protection plan [7].

The CNR-ISAC micrometeorological base, located in the campus of Salento University in Lecce and active since 2003, collects soil moisture θ, precipitation P, and evapotranspiration E on half-an-hour basis (Figure 1).

The analysis of the effect of the surface water balance for the application of the WBM in the deep aquifer recharge requires medium-long term estimations of the water budget, corresponding to several months/years available data of all the surface balance components, and evapotranspiration in particular. The eddy-covariance (e.c.) method for the evapotranspiration measurements allows a direct estimation of the whole evapotranspiration flux in the atmosphere (soil evaporation together with canopy transpiration) over a significant area (footprint) depending on the height of the measurement point and the local meteorological conditions. Since the development of the e.c. method that uses fast response devices to measure turbulent fluctuations of wind and scalars such as water vapor and other gases [8], the increasing assessment of the methodology allowed first the spreading of e.c. measurements for short campaigns and then the creation of networks aimed to long term monitoring of mainly water and carbon fluxes, mainly in Europe and North America (http://ameriflux.lbl.gov, http://www.europe-fluxdata.eu). More recently they have been integrated in the global flux data set (http://fluxnet.fluxdata.org/), containing also the increasing contribution from South America, Asia and Africa to the global network. Most of the stations are mainly devoted to surface-atmosphere exchanges, with focus to the carbon cycle and its implication on climate change, and the use of these data for ground water studies, and karst aquifers in particular, has not been greatly exploited. The ISAC micrometeorological base is part of the Hydrological cycle in the Mediterranean Experiment (HyMEx) (www.hymex.org), and is at present the only station providing long term micrometeorological measurements for a complete observational estimation of the surface water balance in Salento peninsula [9].

The content of this note is organized as follows: In Section 2 the geological setting of the Salento peninsula is shown; in Section 3 the data from the Tiziano Project and the CNR-ISAC base are reported; in Section 4 the derived results are illustrated; in Section 5 the interpretation of the data is given aimed at the direct and indirect estimation of S_y; finally, in Section 6, considerations arising from this interpretation are provided.

2. Geological Setting

The Salento peninsula is a crustal domain with a horst-and-graben tectonic arrangement. The surface is a plateau, gently dipping from SE (mean altitude 100 m a.s.l.) to NW (mean altitude 35 m a.s.l.). Limestones, calcarenites and clayey-marls are the prevailing rock formations cropping out (Figure 2). The horst (raised blocks) is constituted by Cretaceous limestone, while Tertiary and Quaternary stratigraphic units cover the graben (lowered blocks). The surface is rich in well-drained humus-enriched soils (luvisol type), whereas both weak-developed mineral soils (regosol) and weak-differentiated soils (cambisol) are less [10].

As previously mentioned, the main groundwater resource is a deep and extended unconfined aquifer, laterally bounded by semi-confined parts, according to the horst-and-graben tectonic structure. It mainly consists of Cretaceous limestones and contains a large lens of freshwater floating on saltwater of continental intrusion with a thick transition zone of variable saline content (see section A-A′ in Figure 2). The altitude of the water table is about of 3–4 m in the center of Salento peninsula and irregularly decreases toward the coastline, depending on the geological setting and the amount of well exploitation.

Rainfall rapidly infiltrates from the Cretaceous outcrops into the deep aquifer due to the high permeability of the vadose zone, where, as well as inside the deep aquifer, water flows in tectonic fractures and sedimentary planes sometimes affected by karst dissolution. Nine sets of fractures occur: eight are tectonic, secondary to the main NW-SE and NNW-SSE faults mapped in Figure 2, and one is bedding. All the tectonic sets dip about 85 degrees and average strikes are 0, 15, 45, 60, 90, 115, 135, 150 degrees, respectively. Bedding is sub-horizontal [11]. It is to remark that 4 out of 8 tectonic sets always show on the outcrops, while the others superpose randomly. The persistence of the fractures ranges from tens of meters up to hundreds of meters, whereas the spacing ranges from some centimeters up to several decimeters. As far as the apertures are concerned, specific data for the limestones of the deep Salento aquifer are not available. In technical literature, data are instead reported for other rock formations, not only calcareous. Rocha and Barroso [12] measured fracture apertures in granite using a special coring procedure, obtaining an average value of 400 µm at 2.5–4.5 m depth into the bedrock, decreasing to 100 µm at 13.5–14.5 m. Snow [13] back-calculated fracture apertures in granite, gneiss and metavolcanic rocks from packer tests, obtaining apertures 50–100 μm at 15–60 m depth. Similarly, Gale et al. [14] estimated fracture apertures between 5 and 20 μm at depths of several hundred metres in granites at Stripa, Sweden. Motyka and Wilk [15] observed in the carbonate Triassic rocks of the Olkusz mining district (Poland) apertures of around 200 µm (most frequent value). Morgan and Murray [16], by using a hand-held air permeameter, measured values in the range <0.3 µm to 701.5 µm, with an average of of 4.4 µm, in 20 fractures detected in cores (depths between 850 m and 1050 m) of the Arbuckle Group in Oklahoma, a very fine-grained limestone. In consideration of these results, one may consider 200 µm as a generic reference value for the aperture of non-evolutive fractures and in what follows it is assumed for comparisons.

Due to the geological age of the deep aquifer and the evolutive history of the Salento peninsula, the Cretaceous limestones present variable degrees of karstification [11]. Karst systems composed by large cavities and galleries are next to fractured limestone masses seemingly not affected by carbonate dissolution, thus the effects of the speleogenesis processes within the Cretaceous limestones appear quite random.

3. Data Collection

Data set of water table fluctuations was extracted from the Tiziano project database (www.sit.puglia.it/portal/portale_cis). Piezometric levels, recorded by automatic instruments (half hour interval reading), are available for several wells at 1 cm resolution as daily averages. Twenty monitoring wells on the peninsula operated among which LE-PS17 and BR-NM2 are considered herein, for the reasons described below. Both these wells are located within large agriculture areas, where the exploitation of the groundwater resource for irrigation is high, especially during the warm-dry spring-summer season (Figure 2 and Figure 3).

The second data set concerns the ISAC-CNR micrometeorological base. It is located at 40°20′11.7′′ N, 18°7′16.6′′ E (Figure 2), at an altitude of 21 m above sea level and, as previously mentioned, is active since 2003 providing data of surface water and energy balance in the area of the Salento University Campus. The station includes an eddy covariance system on a 15-m mast measuring turbulent vertical sensible heat flux, water vapor flux (evapotranspiration), wind speed, air temperature and humidity. It is completed by an ancillary automatic meteorological station collecting long and short-wave radiation, precipitation and soil data at two depths (temperature at 2 and 5 cm, soil heat flux at 2 and 5 cm, soil moisture at 2 and 35 cm). The e.c. system is composed by a sonic anemometer (Gill R2) coupled with a fast response krypton hygrometer (Campbell Kh20), both located over the telescopic mast at about 12 m over the surface (Figure 1). This allows for a 90% footprint for the scalar fluxes of some km around the mast (fetch along where 90% of the collected flux originated) [4]. The raw data are processed in real time to obtain averaged half hour vertical fluxes from the calculated covariance between wind velocity and water vapor fluctuations. The post processing and validation of the raw e.c. data is relatively complex and its description is out of the scope of this presentation and can be found elsewhere [8,9,17]. The soil moisture data are collected at two depths by a couple of capacitive sensors (Decagon Ec5) previously calibrated using plates filled with local soil and weighted at variable water content. All data are automatically sent to the main frame and stored in the web database (www.basesperimentale.le.isac.cnr.it), and periodically also in the HyMeX database. Figure 4a–c show some seasonal averages and totals for the period of study from the ISAC-CNR database. The air temperature, soil moisture and air humidity data show a clear separation between the warm-dry season (spring–summer) and the cool-wet season (autumn–winter). The total precipitation of the two hydrological years of study (from spring 2009 to winter 2011) amounts to 1691 mm while the total evapotranspiration amounts to 650 mm at the ISAC CNR station.

For the application of WBM, a relatively long period of monitoring has been considered, in which both the mentioned data are continuous and of good quality. The period spans from 2009 to 2011, including 11 large precipitation events that have been selected for the analysis. These events constitute the basis of the comparison between the net infiltration, coinciding with the groundwater storage of Equation (1), and as previously mentioned, the piezometric heights recorded in wells BR-NM2 and LE-PS17 of the Tiziano project. These two wells have been selected for the following reasons:

(1)

the CNR-ISAC micrometeorological station is located approximately in the middle between the wells (Figure 2);

(2)

they show quite continuous and reasonable quality data during the monitoring period;

(3)

the rock of the aquifer is supposed to structurally similar at the two sites; regarding LE-PS17, the upper portion of the vadose zone is made up by fine calcarenites (Figure 2);

(4)

two meteorological stations providing precipitation data are close to the two sites (Manduria and Corigliano rain gauges in Figure 2 and Figure 3); they belong to the network of the Civil Protection Agency of the Apulian region (data available at http://www.protezionecivile.puglia.it/centro-funzionale). The altitudes of the stations are 80 m and 95 m, respectively. They have been used for the local estimation of the precipitation over the two wells respectively.

The rationale underpinning the assessment of the specific yield, in considering short-time relevant precipitation events only, is based on the following:

(1)

a short-time signal in the cross-correlation between surface water infiltration and variations of the deep aquifer level is apparent (see Delle Rose and Martano [4]); according to these results the recharge events studied here have been considered as generally extinguished after three dry days;

(2)

the selection of strong short-time rainy events only is expected to minimize uncontrolled disturbances in the piezometric level variation due to longer time discharge effects and to the possible use of the selected wells for agricultural irrigation;

(3)

a similar procedure has been already used elsewhere [6], although for a shallow aquifer.

In Figure 5 the water levels measured at the two wells are plotted against time. It must be noted that 11 out of the 10 considered precipitation events occurred during Autumn/Winter semesters, i.e., the recharge periods of the deep aquifer [4].

With reference to Equation (1), for the selected events the precipitation amount P (measured at the meteorological stations Manduria e Corigliano), the corresponding values of evapotranspiration E, soil moisture change θ_f−θ_i (measured at the CNR-ISAC station) and the piezometric level change Δh at the two wells are shown in

Table 1. Synopsis of the data used for the estimation of S_y: final day of the event day#, precipitation P, evapotranspiration E, piezometric level change Δh; superscript 1 and 2 refer to well LE-PS17 and well BR-NM2, respectively, 3 refers to the CNR-ISAC base, * indicates that consecutive events have been aggregated; all quantities in mm except (θ_f−θ_I) (volumetric depth-averaged moisture content difference, non-dimensional).

Year	Day#	P1	P2	P3	E	θf−θi	Δh1	Δh2
2009	69	72.4	42.8	51.6	11.1	0.019	120	20
2009	121	31.0	19.0	48.6	15.7	0.006	40	40
2009	283	50.2	40.8	33.8	8.0	0.005	40	50
2009 *	298	78.2	59.0	76.0	3.0	0.017	120	70
2009 *	315	60.8	65.4	135.6	1.5	0.029	200	120
2009	353	75.0	39.4	19.0	2.2	0.010	80	40
2010	50	90.2	74.6	63.4	7.1	0.009	180	150
2010	72	136.4	58.8	65.8	9.1	0.021	120	50
2010 *	293	18.0	116.6	155.0	9.2	0.027	60	120
2010 *	314	74.6	49.4	139.0	15.2	0.013	20	160
2011	34	50.8	54.6	61.0	1.2	0.028	10	10

. Note that θ_f−θ_i values are the average of the values collected at 2 cm and 35 cm. The criterion for the selection of such events is the following: at least one day with more than 30 mm precipitation (note that the piezometric height resolution is 10 mm) and more than 1mm precipitation per consecutive day of event (conventionally ‘rainy’ day). As expected, although the precipitation events are common to both the sites, they differ for the total amount. Considering this data availability and variability and the short-term storage, the closest precipitation event data for each of the two wells have been considered as the main short-term contribution to the storage; then, the precipitation data have been corrected with the evapotranspiration and the soil storage measured at the CNR-ISAC station. It is to remark that in two months (October 2009 and October 2010) two couples of strong events occurred one very close to the other (in October 2010 even with a relevant precipitation in between); they have been included in a first analysis as four separate events and in a second analysis as two compounded and large events (and considering also the contributions in between).

In terms of magnitude, and in particular for short rainy events, E and (θ_f−θ_i)t_s represent smaller corrections to the main contribution term P in Equation (1), and the thickness t_s of the surficial soil cover is fixed to 0.5 m [4].

4. Results

In Figure 6 and Figure 7 the values of ΔS of Equation (1) are plotted against the corresponding change Δh (in mm) of piezometric level for the two wells, respectively.

The solid lines are least-squares regression through the origin, whose slope is the specific yield for the following assumptions:

• no corrections applied in estimating the storage, say ΔS = P (△);
• ΔS calculated by Equation (1) (○);
• as above, but with disaggregation of the two largest 2009 and 2010 events, say 11 total events instead of 9 (+).

Table 2. S_y values resulting from least-square regression of Figure 3 and Figure 4 (slopes of the regression lines) and residual standard deviation values; (△) no corrections applied, (○) with correction of Equation (1) and disaggregated values, (+) with correction of Equation (1) and aggregated values.

Well	Sy(△)	S.d. (△)	Sy(○)	S.d. (○)	Sy(+)	S.d. (+)
LE-PS17	0.59	39	0.48	34	0.48	37
BR-NM2	0.55	34	0.45	32	0.49	33

resumes the results of the least-square regression in terms of slope and residual standard deviation in mm. It shows that the application of the (although approximated) corrections has a general tendency to increase the confidence of the regression, decreasing the residual standard deviation, and that the obtained result for the specific yield is quite robust in terms of aggregating/disaggregating precipitation events, showing almost no significant variation between the two sites. However, the application of the corrections has a more significant effect on the estimation of S_y.

5. Interpretation of the Results

As previously mentioned, one intent is to compare the value of S_y as obtained with WBM to a semi-quantitative estimation of the porosity of the rock mass of the aquifer by means of geological and structural data.

In order to properly define S_y in fractured and karstic rocks, the following quantities are introduced with reference to a volume V_t of rock: V_k, volume of all the openings due to the karst dissolution, V_f, volume of all the pre-existing fractures, including sedimentation joints and tectonical fractures, V_m, volume of all the pores in the matrix blocks. Consequently, three porosities are defined: n_k = V_k/V_t, porosity of the karstic features, n_f = V_f/(V_t − V_k), fracture (secondary) porosity, n_m = V_m/(V_t − V_k − V_f), matrix porosity. The specific yield S_y is therefore [18]:

S_y = β_kn_k + β_f (1 − n_k) n_f + β_m (1 − n_k − n_f) n_m

(2)

where β_k, β_f, β_m are the fractions of the void volumes V_k, V_f, V_m, respectively, containing the water that is released when the water table drops. The accurate measurement of n_k, n_f, β_k, β_f is impracticable, the values depending upon extent, spacing and pattern of fractures and karst openings, features that are in turn intrinsically random. At least for n_f, approximated values may be provided by analyzing the fracturing state of the formation. In addition, β_k, β_f are close to unity and β_m is usually very small in limestones, therefore Equation (2) may be reduced to the following:

S_y − n_k + (1 – n_k) n_f

(3)

With reference to Section 2, by assuming the occurrence of four to eight tectonic sets and the bedding, and assuming a spacing of 10 cm, in a unit cubic meter of the rock of the aquifer approximately 50 to 80 highly persistent fractures could be found with a total planar extension of 50 to 80 m² (as reported above, the fracture length is at least tens of meters). A specific yield S_y of almost 0.5, corresponding to the secondary porosity n_f, would be consistent with apertures of 5.5 to 10 mm, consistently larger than the upper bound of 200 µm for non-evolutive fractures (see Section 2). Even by assuming larger fracture intensities the value would be still too high, therefore it can be stated that the fracturing state in the portion of the rock where fluctuations occur is strongly evolutive, with a pronounced dissolution. In fact, the seasonal fluctuations of the water table are accompanied by upward/downward fluid flows enhancing dissolution along the zone of fluctuation. Over time more open pathways and large voids around the water table generate [19]. Karst dissolution processes close to the water table of an unconfined limestone aquifer are expected by the speleogenesis theory proposed by Swinnerton [20] and further refined by Ford and Ewers [21]. Due to the geological history of the Salento Peninsula, the water table of the deep aquifer reached the current depth 6000 years ago, at the end of the so-called Flandrian Transgression. Since then, an epiphreatic zone, characterized by fractures enlarged by dissolution of calcium carbonate, developed in the zone of the fluctuations. Especially for high fracture frequencies, the karst processes may have produced the development of water-table caves [22]. Evidences of water-table caves are observed along the sea border of the Salento peninsula [23]. As a consequence, by recalling Equation (2) and by considering 200 µm for the apertures, the porosity n_k of the karstic features would result equal to 0.48. This pronounced value seems compatible with the opening of large features in thousands of years when first-order kinetic dissolution of limestones takes place [24].

One should weigh up the high value of S_y (0.45–0.49) of Table 2. Similar values were estimated through direct porosity measurements for the unconfined Floridan aquifer and range between 0.1 and 0.45 [25]. Bolster et al. [26] found large values of S_y, up to 0.57, for the Byscaine aquifer by using canal-drawdown tests. It is opinion of the authors that the value estimated herein are plausible, however, in consideration of the limited amount of available wells for the application of WBM, further studies necessitate in order to verify such a statement, also resorting to the available analytical/numerical models of the karstic dissolution of fractures in limestones [22].

6. Summary and Final Considerations

In this note the Water Budget Method (WBM) is applied to a number of relevant precipitation events at two selected wells in the unconfined karst aquifer of Salento peninsula, resulting in a specific yield of 0.45–0.49 at both the sites, more consistently when the precipitation recharge is corrected by estimates of the evapotranspiration and soil storage. The application of WBM over short-time recharge events is supported by previous observation of a short-time cross correlation between surface infiltration and piezometric response for selected wells of the same area, and by the assumption of negligible agricultural pumping withdrawal during rainfalls. However, this range is way far higher than a (rough) estimate obtained by considering sound values of frequency and apertures of the fractures when non-evolutive, thus requiring a dominant effect of the karstic dissolution to explain the discrepancy.