A Double Scale Methodology to Investigate Flow in Karst Fractured Media via Numerical Analysis: The Cassino Plain Case Study (Central Apennine, Italy)

Lancia, M.; Saroli, M.; Petitta, M.

doi:https://doi.org/10.1155/2018/2937105

Geofluids

On this page

Abstract Introduction Materials and Methods Results and Discussion Conclusions Conflicts of Interest References Copyright Related Articles

Research Article | Open Access

Volume 2018 | Article ID 2937105 | https://doi.org/10.1155/2018/2937105

A Double Scale Methodology to Investigate Flow in Karst Fractured Media via Numerical Analysis: The Cassino Plain Case Study (Central Apennine, Italy)

M. Lancia,¹M. Saroli,¹and M. Petitta²

Academic Editor: Mauro Cacace

Received20 Jun 2017

Revised08 Nov 2017

Accepted23 Nov 2017

Published18 Jan 2018

Abstract

A methodology to evaluate the hydraulic conductivity of the karst media at a regional scale has been proposed, combining pumping tests and the hydrostructural approach, evaluating the hydraulic conductivity of fractured rocks at the block scale. Obtaining hydraulic conductivity values, calibrated at a regional scale, a numerical flow model of the Cassino area has been developed, to validate the methodology and investigate the ambiguity, related to a nonunique hydrogeological conceptual model. The Cassino plain is an intermontane basin with outstanding groundwater resources. The plain is surrounded by karst hydrostructures that feed the Gari Springs and Peccia Springs. Since the 1970s, the study area was the object of detailed investigations with an exceptional density of water-wells and piezometers, representing one of the most important karst study-sites in central-southern Italy. Application of the proposed methodology investigates the hydraulic conductivity tensor at local and regional scales, reawakening geological and hydrogeological issues of a crucial area and tackling the limits of the continuum modelling in karst media.

1. Introduction

Carbonate rocks have a heterogeneous and anisotropic hydraulic conductivity difficult to characterize, with a variation in space and in geological time [1]. The hydraulic parameters depend on the hierarchal rock structure with a facilitated circulation in the downgradient direction [2]. With a negligible permeability matrix, interruptions in the continuity of the material (i.e., faults, fractures, dissolution conduits, and coarse-grained channel fills) constitute the preferential path for the groundwater discharge and contaminant transport. For their economic importance [3], karst media have been fully investigated, including the definition of water balance [e.g., [4]], use of tracers [5], correlation between the rainfall rate and the spring hydrographs [6], hydrochemical samplings [7], or mixing hydrogeological and geophysical data [8]. Usually, followed approaches take into account mean hydraulic parameters, valid at regional scales. Modern methods investigated the fractured and karst media by a discrete approach: Atkinson [9] described the groundwater flow in a Carboniferous karst aquifer combining turbulent conduit flow and Darcian flow in fine fractures; Andersson and Dverstorp [10] predicted the groundwater flows through a network of discrete fractures statistically generated, via a Monte-Carlo simulation; Berkowitz [11] analysed open questions of flow and transport in fractured geological media; Maramathas and Boudouvis [12] described a deterministic mathematical method by the characterization of the fractal dimension of the network; Pardo-Igúzquiza et al. [13] adopted a method to generate conduit geometries via a statistical approach, based on speleological surveys; Borghi et al. [14] proposed a stochastic method to generate a karst network; Bauer et al. [15] studied the karstification process of a single conduit via a hybrid continuum-discrete approach; Liedl et al. [16] used a continuum pipe-model describing the evolution of a karst media; Langevin [17] calibrated a flow model describing the fracture zone; Masciopinto et al. [18] focused on the fractures and a laminar/nonlaminar flow. However, improvements of the porous equivalent approach have been adopted in fractured or karst media, via numerical simulation of the flow [19, 20] or transport model [21]. Karst media offer several solutions that can be applied at the time, involving a discrete approach, a porous equivalent, or a combination of both [22].

In this work, we investigated the fractured karst rocks of the Cassino plain, located in Central Italy (Figure 1) and fully studied since the 1970s, due to its outstanding groundwater resources. Celico [23–25] and CASMEZ [26] focused on the Mesozoic and Cenozoic karst layers with a high density of pumping tests and other indirect investigations. Nevertheless, geological and hydrogeological setting of the Cassino plain is debated with two different conceptual models (Figure 2). Area is appropriate for a continuum approach in the study of the fractured karst rocks at the regional scale. In the area, investigated karst layers were not affected by the huge karst processes related to the Mediterranean salinity crisis event [27], due to the recent uplift of the Apennine orogen [28]. The subsequent normal faulting event deposited thick continental sequences in the low-stand sectors, with a function of permeability boundary. These sequences affected the groundwater circulation, hindering the development of a karst network in the spring areas. The combination between tectonic activity and sedimentation rate varied the base level of fractured carbonate aquifers [29]. As a consequence, basal springs are outstandingly high and steady, without impulsive response to seasonal recharge. In accordance with other hydrogeological studies in Central Apennine [30], groundwater flows through fracture networks, instead of karst conduits. Thus, fractured karst masses of the studied area cannot be compared with the other Mediterranean karst basins or older Paleozoic carbonate rocks.

Figure 2

(a) Geological and hydrogeological map of the Cassino plain, from Gari to Peccia Springs. Key to the legend: (1) Quaternary unit (QU) made of clastic aquifers and aquitards; (2) terrigenous unit (TU) made of clastic aquiclude and aquitards from Miocene foredeep clayey deposits; (3) limestone aquifer unit (CU) made of continuous layers of Mesozoic to Cenozoic limestones; (4) dolostone unit (DU) made of Lower Mesozoic massive dolostones; (5) main tectonic lineament; (6) karst punctual spring; (7) karst streambeds; (8) water-well exploiting the carbonate aquifer; (9) locations of the geomechanical analysis; (10) town; (11) mountain peak; (b) examples of boreholes stratigraphies from CASMEZ [26]. No karst caves have been reported in the original documents; (c) redrawn hydrogeological model from Boni and Bono [42]. The authors [4, 42] schematize the carbonate bedrock via a horst and graben setting and a flux from Gari Springs (G1, G2, and G3) to Peccia Springs. (d) modified hydrogeological section from Celico [23]. Celico [23–25] supposed the interruption of the carbonate bedrock and no underground connection between Gari Spring and Peccia Springs.

To characterize the groundwater flow through the carbonate layers of the studied area and verify the importance of the fracture networks over the karst conduits, a methodology has been proposed to estimate the hydraulic conductivity tensor (Figure 3). This methodology starts with geomechanical surveys to define average fracturing setting parameter, further implemented by a FEM (Finite Element Model) analysis at a 1 m scale. This first step is in accordance with the hydrostructural approach [31–33]. Estimated dataset has been calibrated with pumping tests to derive a reliable hydraulic conductivity tensor and successively validated with a Darcian numerical analysis at a regional scale. Methodology, validated for the Cassino plain, can be applied across the Central Apennine, thanks to the similar hydrogeological setting. Based on the new collected data and the review of the literature, this study aims to(i)define a hydraulic conductivity tensor of the investigated fractured karst media, representative at regional scale;(ii)verify the adoption of a downscaling procedure, able to shift from 1 m scale model to a regional one;(iii)offer an in-depth analysis related to the highly debated conceptual models of the study area.

Figure 3

Flow diagram for the proposed methodology. Geomechanical analysis including the evaluation of orientation, spacing, and aperture (a) allows constructing discrete-fractured cubes (b) representing the average fracturing condition. Numerical models at 1 m scale (c) have been calibrated via pumping tests (d) estimating the hydraulic conductivity tensor at regional scale (e). The latter has been tested via numerical simulation at regional scale and compared with conceptual models from literature (f).

2. Materials and Methods

2.1. Geological and Hydrogeological Setting

The Cassino plain is an intermontane basin, located downstream of the Latina Valley in the Central Apennine (Figure 1). Mesozoic dolostones constitute the base of the stratigraphic sequence, followed by Mesozoic and Cenozoic limestones. The entire carbonate sequence has a thickness of 2000 m [34]. Upwards in the sequence, Tortonian sandstone and clay levels represent the foredeep deposits of the Apennine chain, with a thickness of 800 m [35]. Quaternary continental deposits follow [36]. The spatial distribution of the stratigraphic sequence is strongly controlled by tectonics. Since the Messinian [28, 37], the Mesozoic and Cenozoic carbonate deposits experienced contractional tectonics and were structured in a thrust and fold geometry. The latter was completely masked by strike-slip tectonics along the Val Roveto-Atina-Caserta line [38, 39] and a subsequent normal faulting event. The last tectonic event dissected the Mesozoic and Cenozoic carbonate sequences by normal faults with NE-SW and NW-SE trends, creating a horst and graben geometry. Indeed, the flat alluvial morphology of the plain is interrupted by two main monoclines (Trocchio Mt. and Porchio Mt.) made of Mesozoic and Cenozoic limestones and bounded along three sides by normal faults (Figure 2(a)). Numerous minor carbonate structures outcrop from the alluvial plain. Unfortunately, in literature a unique structural setting has not been agreed: Zalaffi [40] assumed Trocchio and Porchio Mts. as horsts but, two years later, Accordi interpreted the latter as klippe [41]. The problem related to the structural setting was accentuated in case of hydrogeological modelling [4, 23–25, 42]. All the conceptual models are in accordance with defining the Gari Spring (elevation of 41–30 m a.s.l., discharge of 18 m³/s) and Peccia Spring (elevation of 29–25 m a.s.l., discharge of 5 m³/s) as related to the karst Mesozoic-Cenozoic aquifers. From a hydrogeological point of view, the sequence can be schematized in four units. The dolostone unit (DU) represents the basal aquifers (Figure 2(a)). Though it has several fractures and karst elements, the hydrogeological characteristics are similar to terrains with interstitial hydraulic conductivity and high storage capacity [4]. The Cinquina water-well (CWW in Figure 2(a)) suggests a hydraulic conductivity of [m/s]. Conversely, a continuous network of fractures affects the younger limestone unit (CU), contributing to a higher conductivity value. The limestone unit was investigated via 36 boreholes furnished with pumping tests [26], during the construction of the western Campania aqueduct (WCAWW in Figure 2(a)). The thickness of investigated limestones varies between 50 and 200 m (Figure 2(b)), well representing the fracturing average condition. The pumping tests show a different hydraulic conductivity value, fluctuating from to [m/s]. Upward in the sequence, the upper Tortonian clayey unit has no aquifer capacity and a reduced hydraulic conductivity, so it can be assumed as aquiclude if compared to the limestone and dolostone sequences (Figure 2(a)). In the end, the quaternary unit (QU) has a variable importance related to their extension (Figure 2(a)). At a regional scale, these deposits act as an aquitard of the limestone and dolostone sequence. From the reconstructed hydrogeological sequence the importance of the limestone unit (CU) and its spatial characterization is evident. Boni and Bono [42] linked the Gari Springs to the limestone sequences. They individuated in the Gari Spring area (Figure 2(a)) several carbonate horsts, delimitated by the quaternary unit. This setting gives origin to localized springs with different elevation: the Cassino urban area springs (G1 in Figure 2(a); elevation of 41 m a.s.l.), the Terme Varroniane Springs (G2 in Figure 2(a), elevation of 33-32 m a.s.l.), and the Borgo Mastronardi Springs (G3 in Figure 2(a); elevation of 31-30 m a.s.l.). Centamore et al. [37] proposed a horst and graben model from Gari to Peccia Springs, based on the Zalaffi [40] evidence (Figure 2(c)). Consequently, Boni et al. [4] supposed a groundwater flow of 4 m³/s from Gari Spring to Peccia Springs (Figure 1), assuming a single hydrostructure (Simbruini-Ernici-Cairo Venafro Mts.). Conversely, Celico and Stanganelli [43] assumed the Peccia Springs fed by the independent Venafro Mts. hydrostructure (Figure 1). Based on the Accordi [41] model and new collected data, Celico [23–25] schematized Trocchio and Porchio Mts. as klippe (Figure 2(d)).

2.2. Hydraulic Conductivity Evaluation

The study of the groundwater flow through the karst media requires the definition of the hydrodynamic parameters of the carbonate layers, their spatial variability, and the scale of the analysis. At regional scale, the hydraulic conductivity is one of the most important parameters and its definition allows characterizing the groundwater flows as well as the contaminant transport processes. Because groundwater flows along preferential ways (e.g., fractures and conduits), the magnitude of the hydraulic conductivity should be represented along the three main directions, for full characterization. The Cassino plain (Figure 1) represents an ideal case for the study of the fracturing setting; the carbonate layers outcrop along the whole plain and the area was fully investigated from a geological and hydrogeological point of view. To estimate the hydraulic conductivity in fractured karst media, we calibrated a methodology combining the hydrostructural approach and the pumping test data (Figure 3). The hydrostructural approach [31–33] requires geomechanical surveys to recognize and characterize the joint sets (i.e., a group of subparallel joints) (Figure 3(a)). The approach starts with geomechanical measures, whereas the karst layers outcrop. A total number of 40 geomechanical stations have been performed. Every geomechanical site is located where the bedrock is well exposed, without detritus or vegetation. If possible, measurement sites have been chosen in correspondence with excavations or quarries. Indeed, two exposed rock-walls with different orientations better describe the three-dimensionality of the rock masses. In addition, geomechanical stations were also located nearby roads or other linear features. According to the principles of the geomechanical standards [45], the joint sets have been identified, taking note of their spacing and opening. For our purposes, the joints were assumed smooth and plain; indeed the roughness of the joint wall along the topographic surface is not representative of the rock-mass condition in depth. Because the fracturing condition has a strong spatial variability difficult to characterize, in terms of density and characteristic of the joints, the definition of a mean fractured setting of the karst media was necessary. After the definition of the recognized joint sets, different values of spacing and opening have been considered, constructing a discrete fracture network (Figure 3(b)). Thus, the different fracturing cases have been considered. Fractures have been bounded in cubes with a side length of 1 m, oriented with the following annotation: positive = east, positive = north, and positive = up.

For every fracturing setting, numerical simulations by a FEM code have been performed to calculate the flow through the faces of the cubes (Figure 3(c)). A FEM gives the opportunity to easily modify the opening of the fractures, not varying the geometry of the model. Due to the high heterogeneities of the limestone layers, the spacing dataset has been separated in three intervals as well as the opening dataset in four intervals. The combination of the two parameters (spacing and opening) gave back a final dataset of twelve different fracturing cases. COMSOL® software and its Fracture Flow Interface in the Porous Media and Subsurface Package were used. Inside the single fractures, the following equations have been assumed:where is gradient along the fracture extent; is average fracture opening; is flow velocity; is liquid density; is outflow; is fracture permeability; is fluid viscosity; is pressure.

According to the cubic law [46],where represents the hydraulic head variation and a form factor; the permeability along a single fracture can be written in the following form [47]: Though several studies outlining the limit of the cubic law [48, 49], it is a point of reference for the rock mechanics and it is usually used for flow in fractured systems [50–52].

During the construction of the numerical model, fractures were represented as planes subdivided via a free-triangular mesh. The number of the triangles varies in function of the considered spacing set as well as the number of fractures inside the single cube; usually it is around 60,000. The dimension of the triangles strongly diminishes near the intersection between two planes. The flow through the rocky matrix was neglected, assuming that water cannot flow outside the fracture discrete network (Neumann condition type). As supplementary boundary condition, an average hydraulic gradient of 0.5%, typical for the limestones of the Central Apennine [4, 25] has been applied (Dirichlet condition type). The gradient was applied firstly along the -axis. Simulations were repeated applying the gradient along the -axis and the -axis. Measuring the simulated discharge, a porous equivalent hydraulic conductivity has been calculated for every cube and along the three directions. Therefore, the hydraulic conductivity dataset was computed, varying the fracturing condition. Dataset has been calibrated with pumping tests (Figure 3(d)), in order to individuate a permeability tensor (Figure 3(e)) valid at regional scale. Numerical analysis of the flow model at regional scale (Figure 3(f)) allows validating the proposed methodology for the study of karst media.

2.3. Geometry and Boundary Condition at Regional Scale

The calibrated hydraulic conductivity tensor has been validated via numerical analyses at regional scale. A geological section for numerical purposes has been constructed (Figure 4) along the A-B section trace of Figure 2(a), assuming the horst and graben conceptual model [4, 42], with some adjustments due to new surveys and the reinterpretations of borehole stratigraphies [25, 44]. Because the Cassino plain has a nonunique geological and hydrogeological setting, the most recent conceptual model [4] has been preferred. The model has a strong horizontal variability due to the presence of normal faults; therefore a FEM analysis has been preferred, using COMSOL® software and its Porous and Subsurface Package. According to the literature review, 40 domains have been classified in four hydrogeological units: the dolostone unit (DU), the limestone unit (CU), the terrigenous unit (TU), and the Quaternary unit (QU). At regional scale, every unit is considered homogeneous and porous and the Darcy Law valid [53]. The entire geometry was discretized via 39.457 elements, in a free-triangular mesh, with refinements in correspondence with the domain vertices. Because the numerical model just verifies the final part of the hydrostructures, a steady-state hydraulic head, from the inner Cairo Mt. to Gari Spring and from Camino Mt. (Venafro Mts.) to Peccia Spring, was applied (Dirichlet condition type). The hydraulic gradient of 0.5% has been chosen in accordance with Boni et al. [4] and Celico [25]. Along the vertical borders and at the bottom, no flow conditions (Neumann condition type) have been applied. To simulate the unsaturated aquifer thickness, a Heaviside function of the pressure (value from 0 to 1) has been considered along the topographic surface generating a hydraulic head dependent boundary (Dirichlet type-condition), according to Chui and Freyberg [54]. The latter condition cancels possible undesired fluxes towards the unsaturated zone, very common during the numerical analysis. Boundary conditions are represented in Figure 4.

Figure 4

Hydrogeological section for numerical purposes from Centamore et al. [37] modified on the basis of additional data [23–25, 44] and imposed boundary condition. Red lines indicate the faults. An inner gradient of 0.5% has been set from Cairo Mt. to Gari Spring (G1, G2, and G3) and from Camino Mt. to Peccia Springs (P), no flow condition (Neumann condition at the bottom and at the lateral edges), and a head dependent boundary (Cauchy condition) along the topographic surface. Dolostone and limestone units (DU and CU) represent the aquifer and terrigenous unit (TU) and Quaternary unit (QU) the aquiclude and the aquitard of the system.

3. Results and Discussion

3.1. Hydraulic Conductivity Dataset from the Hydrostructural Approach

The geomechanical analysis of the rock masses pointed out a wide dataset of fractures, collected in correspondence with the dolostone and limestone outcrops. The dolostones are located at the northeast boundary of the Cassino plain between Cervaro and San Vittore (Figure 2(a)). In this area, the dolostone layers are tectonically disturbed and weathered. During the surveys, a systematic geomechanical analysis of the fractures was not possible. The fracture density is very high and the rock-mass has a low compactness, assuming characteristic of a loose, cemented sand. Surveys pointed out fine silt in the main fractures. The hydraulic conductivity measure of the dolostone ( m/s) coming from the Cinquina water-well (Figure 2(a)) is in accordance with the geomechanical observation. Zalaffi [55] described the same hydrostructural setting of the Matese Mts. dolostones (Figure 1).

Conversely, the hydrostructural approach has been fully applied to the limestone layers, between Gari and Peccia Springs. The results of 36 geomechanical stations pointed out a well-organized fracturing setting with several joint sets that take origin from the stratigraphic and tectonic features of the limestone rock-mass. The dip-slip measures of the bedding vary from N.90° to N.130° across the studied area. Data are in accordance with the schematization of Boni et al. [4] and Celico [25]. The dip measures are averagely around 30° except at the Montecassino Hill, where a recumbent anticline fold structure induces a local variation of the dip between 10° and 80°, in accordance with the geological data of Saroli et al. [44]. The spacing (i.e., thickness of every strata) is averagely around 15–40 cm but locally can reach values of 1 m. The opening of these joints is very low, around 0.1 mm, but it is locally higher when the fractures are affected by flexural slip. The bedding can reach an opening of 1 mm. The tectonic joints instead have not a homogeneous dip-slip, and a schematization in different joint sets is required. Collected data shows a diffuse high-angle dip, from 70° to 90°. According to the geomechanical standards [45], all the fractures, with high-angle dip-slip, have been gathered in the same group isolating four tectonic joint sets: the fractures with dip angle in the interval between N.345-N.15 and N.165-N.195 (tectonic joint-set 1); N.65-N.950 and N.255-N.285 (tectonic joint-set 2); N.30-N.60 and N.210-N.240 (tectonic joint-set 3); N.120-N150 and N.300-N330 (tectonic joint-set 4). The spacing and the opening of these tectonic joint sets are variable and depend on the local heterogeneities that affect the rock-mass. After the individuation of the joint sets, the spacing and opening dataset have been divided into classes, to reproduce different degrees of the fracturing condition: in total three spacing classes and four opening classes. A summary of the collected data is in Table 1.

Assuming the simultaneous presence of all the joint sets, cubes with a discrete fracture network have been constructed. The orientation of the fracture has been considered constant, assuming an average value. Then, the three spacing sets have been crossed with the four opening-sets, for a total of twelve numerical analyses. According to the boundary condition, simulated fluxes occur only through the planes that represent the discrete fractures. The numerical simulation pointed out a homogeneous dissipation of the imposed hydraulic head across the opposite face of cubes along the -axis (Figure 5(a)), -axis (Figure 5(b)), and -axis (Figure 5(c)).

(a)

(b)

(c)

Black arrows (Figure 5) outline how the simulated flux is parallel to the fracture orientation, respecting the boundary condition. The ratio between the discharge along the selected direction and the imposed hydraulic head is proportional to the hydraulic conductivity value, if an equivalent porous medium is considered in accordance with the Darcy Law [49]. The simulated hydraulic conductivity range is very wide, from 2 to [m/s] (Figure 6). According to the model, the hydraulic conductivity strongly varies in function of the chosen opening set and subordinately in function of the spacing set. Indeed, inside a single fracture each increment of the opening corresponds to a square increment of the velocity. A wide dataset of hydraulic conductivity up to nine orders of magnitude is in accordance with the huge heterogeneity of the fractured karst media.

Figure 6

Hydraulic conductivity histogram achieved from the hydrostructural approach. Colours indicate the different spacing sets. Dataset is very high and hydraulic conductivity varies from m/s (spacing set 1, opening set 1) to m/s (spacing set 3, opening 3). Calculated value increase with the increase of the opening and the spacing. Hydrostructural approach points out a quasi-isotropic hydraulic conductivity tensor with and slightly lower than . An average anisotropic rate can be valuated around 0.5.

With a wide dataset of results (Figure 6), it would be senseless deriving a reliable hydraulic conductivity value for the limestone unit (CU), based just on the hydrostructural results. However, in every cube and are almost equal and is higher (Figure 6) for all the spacing-opening combinations. Adding more, the ratio remains approximately constant around 0.51. Thus, the results point out a value of the anisotropic rate at regional scale that pumping tests cannot provide. The 36 pumping tests show a hydraulic conductivity magnitude between and m/s, as well as the hydrostructural approach of a ratio around 0.5. By the combination of the hydrostructural approach and the pumping tests, we can assume as acceptable at regional scale the following hydraulic conductivity dataset for the limestone unit (CU): = m/s; = m/s; = m/s.

3.2. Limits of the Hydrostructural Approach

The hydrostructural approach is based on the reconstruction of discrete fracture network and further numerical simulation. Geomechanical analysis reconstructs the heterogeneities of the limestone unit (CU) by punctual measures that cannot be homogeneously distributed. Indeed, vegetation or detritus mostly covers the limestone unit, especially along the piedmont band. Adding more, in absence of explorable karst conduits and caves, all the measures have been achieved just along the topographic surface. Thus, the described setting could not be representative of the fracturing condition in depth. Adding more, the cubic law [46] can be assumed valid under laminar conditions of the flux. Though in the studied area karst processes are young and not well-developed, the flux through wider fractures (karst conduits) and related turbulent conditions has been not taken into account, due to the regional character of the work and the hydrogeological setting of the Central Apennine [29, 30]. However, conduits can influence the results at local scale. In this case specific laboratory tests and field-tests [48–50] are necessary to investigate the relationship between fractures and conduits and the transition from laminar to turbulent flow.

3.3. Validation of the Methodology via Regional Numerical Analysis

The estimated hydraulic conductivity dataset of the Cassino limestone unit (CU) comes from the combination between the hydrostructural approach and pumping tests results. A regional numerical simulation based on the estimated hydraulic conductivity values is a powerful instrument to validate the adopted methodology and the most recent conceptual model [4, 37]. The hydraulic conductivity values of the karst layers derive from the previous results: dolostone unit [m/s], limestone unit [m/s], , and . The remaining hydraulic conductivity values have been set via a literature analysis and pumping tests analyses: terrigenous unit [m/s] and Quaternary unit [m/s]. The constructed numerical section of Figure 4, correlated with the proposed boundary condition, perfectly represents the conceptual model. Indeed, the regional scale simulations outline a groundwater flux from the Cairo Mt. to the Gari Spring and from the Camino Mt. (Venafro Mts.) to the Peccia Springs, entirely developed through the limestone unit (CU). In the numerical model, the simulated velocity field varies from to m/s (Figure 7(a)) with maximum values at the contact between the limestone unit (CU) and the Quaternary unit (QU), near the topographic surface. Adding more, the simulated water springs occur in correspondence with the main Gari Springs (G.1 Spring = 40 m a.s.l., G.2 Spring = 33 and 32 m a.s.l., and G3 Spring = 31 and 30 m a.s.l.), accurately reproducing the hydrogeological setting of this area (Figure 7(a)). A similar setting is observed on the opposite side of the Cassino plain, whereas the velocity field simulates the Peccia Springs (Figure 7(a)).

Figure 7

Numerical analysis of the Boni model [4, 37] (section in Figure 4). Key to the legend: Quaternary unit (QU); terrigenous unit (TU); limestone unit (CU); dolostone unit (DU). Cassino urban Springs G1 (elevation of 41 m a.s.l.); Terme Varroniane Springs G2 (elevation of 32-33 m a.s.l.): Borgo Mastronardi Springs G3 (elevation of 31-30 m a.s.l.); Peccia Spring P (elevation of 25 m a.s.l.). Seepage verifies the Boni model from Cairo Mt. to Gari Spring (G1, G2, and G3) and from Camino Mt. (Venafro Mts.) to Peccia spring (P); (b) in the central sector, seepage from Gari to Peccia Springs is simulated; (c) the upstream and downstream groundwater velocity profile from Gari to Peccia Spring. The downstream/upstream discharge ratio is around 9%, calibrated under the best assumptions, and is not sufficient to justify the Boni conceptual model [4, 37].

Nevertheless, the central part of the section shows very low velocities and a groundwater flux from Gari to Peccia Springs, through the Trocchio and Porchio Mts. (Figure 7(b)). Because the groundwater velocity field and the spring distribution perfectly matches, the next step is the validation of the conceptual model from a quantitatively point of view. In a 2D model, to verify the existence of a 4 m³/s through the Trocchio and Porchio Mts., we based our analysis on the simulated ratio between the downstream and the upstream Gari Spring. The latter is around 9%.

According to the Boni conceptual model [4], there is a flux of 18 m³/s from the Cairo Mt. to the Gari Spring and of 4 m³/s from Gari to Peccia Spring via the Trocchio-Porchio Mts. alignment, for a total discharge of 22 m³/s. Thus, the ratio between downstream and upstream discharge should be higher, around 18% at least. This ratio increases up to 30%, if we consider the minimum values of the Gari Spring discharge (13 m³/s). Consequently, the Boni et al. [4] conceptual model is not numerically verified, along the considered section. The numerical model takes into account favourable conditions to maximize the flux from Gari to Peccia Springs, since there are many possible solutions, compatible with hydrogeological constraints. According to Boni et al. [4], the alignment between the Gari and Peccia spring would represent a continuous groundwater flow direction similar to a conduit, laterally bounded by a thick blanket of terrigenous deposits (TU). Since no lateral groundwater exchanges occur, the 2D section model neglects possible hydraulic head loss related to the real 3D geometry. In order to reduce the hydraulic head loss from Cairo to Peccia Springs, the maximum thickness of the limestone unit (CU) has been assumed in the model. Hydraulic conductivity has been considered constant with depth up to −1500 m a.s.l., neglecting the natural permeability reduction related to the fracture closure; neither lateral variation of the permeability, related to mylonitic bands in depth, nor tectonic-tightening of the limestone unit thickness (CU) has been considered. From a hydraulic point of view, the maximum hydraulic gradient from the highest elevation of Gari Spring to the lowest of Peccia spring has been calculated. No spatial variation of the hydraulic conductivity generates a loss of the hydraulic head, disfavouring the communication between the Gari and Peccia Springs. Thus, as a precaution, a homogeneous value of the limestone unit (CU) has been set. However, small variations of the hydraulic conductivity do not provide any significant variation on the discharge ratio. The only way to increase the upstream/downstream discharge ratio close to the theoretical value of 18% is to switch off the imposed gradient (Dirichlet condition) from Camino Mt., which have a function of backpressure to the flow from Gari to Peccia Spring (Figure 7(c)). This last condition is unacceptable and not in accordance with Celico [23–25] and Boni et al. [4].

4. Conclusions

The proposed methodology, applied to fractured karst media of Cassino plain (Figure 2), allowed defining the average hydraulic conductivity tensor at regional scale. The proposed methodology starts from the geomechanical analysis, defining the average fracturing parameters. A 1 m scale FEM analysis allowed estimating a wide hydraulic conductivity dataset, successively combined with pumping tests. Based on the literature review, application of the proposed methodology, and validation via a regional flow model, the following conclusions can be delineated:(i)1 m scale analysis and the subsequent FEM analysis are appropriate to investigate the hydraulic conductivity anisotropy but simulated magnitude range is too wide, reflecting the heterogeneities of the fractured karst mass at local scale. On the other hand, pumping tests investigate a higher volume but assuming the media homogeneous and isotropic.(ii)Thus, a valid approach to determine the hydraulic conductivity tensor in fractured media, at regional scale, could be obtained by an accurate combination between a detailed scale analysis, based on hydrostructural scale measurements, and average values from pumping tests.(iii)The estimated hydraulic conductivity tensor is valid at regional scale, under Darcian conditions. The numerical model perfectly simulates the distribution of the karst springs along the topographic surface, according to Boni et al. [4] and Celico [23–25]. Nevertheless, under the best assumptions, the model is not quantitatively verified and further data are necessary to clarify the debated conceptual model of an important resource of Central Italy.

The proposed methodology represents an improvement of the continuum approach and can be applied in fractured media, especially when an average value of the hydraulic conductivity tensor is required. For fractured karst media, the Cassino plain and, in general, the Central Apennine region represents a favourable case, since the groundwater preferentially flows in the fracture networks and not in the karst conduits [30]. As future direction, methodology can be implemented with the superimposition of conduits, investigating the transition from laminar to turbulent flow. With this modification, the latter can be applied on every kind of karst media, coupling a discrete and continuum approach.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

References

P. W. Huntoon, “Gradient Controlled Caves, Trapper‐Medicine Lodge Area, Bighorn Basin, Wyoming,” Groundwater, vol. 23, no. 4, pp. 443–448, 1985.
View at: Publisher Site | Google Scholar
P. Galvão, T. Halihan, and R. Hirata, “The karst permeability scale effect of Sete Lagoas, MG, Brazil,” Journal of Hydrology, vol. 532, pp. 149–162, 2016.
View at: Publisher Site | Google Scholar
D. Ford and P. Williams, Karst Hydrogeology and Geomorphology, Wiley, Chichester, UK, 2007.
View at: Publisher Site
C. F. Boni, P. Bono, and G. Capelli, “Schema idrogeologico dellItalia Centrale, (Hydrogeological schemes of Central Italy),” Memorie della Società Geologica Italiana, vol. 35, pp. 991–1012, 1986.
View at: Google Scholar
N. Goldscheider, J. Meiman, M. Pronk, and C. Smart, “Tracer tests in karst hydrogeology and speleology,” International Journal of Speleology, vol. 37, no. 1, pp. 27–40, 2008.
View at: Publisher Site | Google Scholar
F. Fiorillo, “Spring hydrographs as indicators of droughts in a karst environment,” Journal of Hydrology, vol. 373, no. 3-4, pp. 290–301, 2009.
View at: Publisher Site | Google Scholar
P. J. Moore, J. B. Martin, and E. J. Screaton, “Geochemical and statistical evidence of recharge, mixing, and controls on spring discharge in an eogenetic karst aquifer,” Journal of Hydrology, vol. 376, no. 3-4, pp. 443–455, 2009.
View at: Publisher Site | Google Scholar
T. McCormack, Y. O'Connell, E. Daly, L. W. Gill, T. Henry, and M. Perriquet, “Characterisation of karst hydrogeology in Western Ireland using geophysical and hydraulic modelling techniques,” Journal of Hydrology: Regional Studies, vol. 10, pp. 1–17, 2017.
View at: Publisher Site | Google Scholar
T. C. Atkinson, “Diffuse flow and conduit flow in limestone terrain in the Mendip Hills, Somerset (Great Britain),” Journal of Hydrology, vol. 35, no. 1-2, pp. 93–110, 1977.
View at: Publisher Site | Google Scholar
J. Andersson and B. Dverstorp, “Conditional simulations of fluid flow in three‐dimensional networks of discrete fractures,” Water Resources Research, vol. 10, pp. 1876–1886, 1987.
View at: Publisher Site | Google Scholar
B. Berkowitz, “Characterizing flow and transport in fractured geological media: a review,” Advances in Water Resources, vol. 25, no. 8–12, pp. 861–884, 2002.
View at: Publisher Site | Google Scholar
A. J. Maramathas and A. G. Boudouvis, “Manifestation and measurement of the fractal characteristics of karst hydrogeological formations,” Advances in Water Resources, vol. 29, no. 1, pp. 112–116, 2006.
View at: Publisher Site | Google Scholar
E. Pardo-Igúzquiza, P. A. Dowd, C. Xu, and J. J. Durán-Valsero, “Stochastic simulation of karst conduit networks,” Advances in Water Resources, vol. 35, pp. 141–150, 2012.
View at: Publisher Site | Google Scholar
A. Borghi, P. Renard, L. Fournier, and F. Negro, “Stochastic fracture generation accounting for the stratification orientation in a folded environment based on an implicit geological model,” Engineering Geology, vol. 187, pp. 135–142, 2012.
View at: Publisher Site | Google Scholar
S. Bauer, R. Liedl, and M. Sauter, “Modeling of karst aquifer genesis: Influence of exchange flow,” Water Resources Research, vol. 39, p. 1285, 2003.
View at: Google Scholar
R. Liedl, M. Sauter, D. Hückinghaus, T. Clemens, and G. Teutsch, “Simulation of the development of karst aquifers using a coupled continuum pipe flow model,” Water Resources Research, vol. 39, no. 3, p. 1057, 2003.
View at: Google Scholar
C. Langevin, “Sthocastic groundwater flow simulation with a fracture zone continuum mode,” Groundwater, vol. 41, pp. 587–601, 2003.
View at: Publisher Site | Google Scholar
C. Masciopinto, A. Volpe, D. Palmiotta, and C. Cherubini, “A combined PHREEQC-2/parallel fracture model for the simulation of laminar/non-laminar flow and contaminant transport with reactions,” Journal of Contaminant Hydrology, vol. 117, no. 1-4, pp. 94–108, 2010.
View at: Publisher Site | Google Scholar
J. Kordilla, M. Sauter, T. Reimann, and T. Geyer, “Simulation of saturated and unsaturated flow in karst systems at catchment scale using a double continuum approach,” Hydrology and Earth System Sciences, vol. 16, no. 10, pp. 3909–3923, 2012.
View at: Publisher Site | Google Scholar
C. Cherubini, C. Giasi, and N. Pastore, “Fluid flow modeling of a coastal fractured karstic aquifer by means of a lumped parameter approach,” Environmental Earth Sciences, vol. 70, no. 5, pp. 2055–2060, 2013.
View at: Publisher Site | Google Scholar
C. Masciopinto and D. Palmiotta, Flow and Transport in Fractured Aquifers: New Conceptual Models Based on Field Measurments, vol. 96, pp. 117-163, 2003.
M. Sauter, T. Geyer, A. Kovács, and G. Teutsch, “Modeling of the hydraulics of karst aquifers - an overview,” Grundwasser, vol. 11, no. 3, pp. 143–156, 2006.
View at: Publisher Site | Google Scholar
P. Celico, “Indagini Idrogeologiche per la progettazione delle opere di parziale captazione dei deflussi sorgivi del Fiume Gari (Cassino), Costruzione dell’Acquedotto della Campania Occidentale,” in Hydrogeological Investigations for Gari Spring Water Captation, p. 200, Cassa per il Mezzogiorno-Servizi di Intervento Straordinari, 1976.
View at: Google Scholar
P. Celico, “Schema idrogeologico dell’Appennino carbonatico centro-meridionale (Hydrogeological scheme of the central-southern carbonate Apennine),” Memorie e Note dell’Istituto di Geologia Applicata, vol. 14, pp. 1–97, 1978.
View at: Google Scholar
P. Celico, “Idrogeologia dei massicci carbonatici, delle piane quaternarie e delle aree vulcaniche dell’Italia centro meridionale (Marche, Lazio meridionale, Abruzzo, Molise e Campania),” in Hydrogeology of Carbonate Basins, Quaternary Plains and Volcanic Areas of Central-Southern Italy-Marche, Southern Lazio, Molise and Campania Regions, vol. 4 (2), pp. 1–203, Quaderni della Cassa per il Mezzogiorno, 1983.
View at: Google Scholar
CASMEZ, “Schemi Idrici Intersettoriali del Lazio Meridionale, Tronto, Abruzzo, Molise e Campania (Hydrogeological schemes for Tronto River and Southern Lazio, Abruzzo, Molise and Campania),” Progetto Speciale 29, Catalogo, Bari, 1978.
View at: Google Scholar
M. Bakalowicz, “Karst groundwater: A challenge for new resources,” Hydrogeology Journal, vol. 13, no. 1, pp. 148–160, 2005.
View at: Publisher Site | Google Scholar
D. Cosentino, P. Cipollari, M. Marsili, and D. Scrocca, “Geology of the Central Italy: a regional review,” in The geology of Italy; tectonics and life along plate margins, M. Meltrando, A. Peccerillo, M. Mattei, S. Conticelli, and C. Doglioni, Eds., vol. 36, pp. 1–37, Journal of Virtual Explorer, 2010.
View at: Google Scholar
M. Petitta, “Hydrogeology of the Middle Valley of the Velino River and of the S. Vittorino plain (Rieti, Central Italy),” Italian Journal of Engineering Geology and Environment, vol. 1, pp. 157–181, 2009.
View at: Google Scholar
M. Petitta, P. Primavera, P. Tuccimei, and R. Aravena, “Interaction between deep and shallow groundwater systems in areas affected by Quaternary tectonics (Central Italy): A geochemical and isotope approach,” Environmental Earth Sciences, vol. 63, no. 1, pp. 11–30, 2011.
View at: Publisher Site | Google Scholar
M. Surette and M. D. Allen, “Quantifying heterogeneities in variably fractured rock using a hydrostructural domain,” Geological Society of America Bullettin, vol. 120, pp. 225–237, 2008.
View at: Google Scholar
M. Surrette, D. M. Allen, and M. Journeay, “Regional evaluation of hydraulic properties in variably fractured rock using a hydrostructural domain approach,” Hydrogeology Journal, vol. 16, no. 1, pp. 11–30, 2008.
View at: Publisher Site | Google Scholar
N.-Y. Ko, S.-H. Ji, Y.-K. Koh, and J.-W. Choi, “Evaluation of two conceptual approaches for groundwater flow simulation for a rock domain at the block-scale for the Olkiluoto site, Finland,” Engineering Geology, vol. 193, pp. 297–304, 2015.
View at: Publisher Site | Google Scholar
G. Accordi and F. Carbone, “Sequenze carbonatiche Meso-Cenozoiche, (Meso-Cenozoic carbonatic sequences),” in in Carta delle Litofacies del Lazio-Abruzzo aree limitrofe, G. Accordi and F. Carbone, Eds., vol. 111 of Quaderni della Ricerca Scientifica, pp. 11–92, Roma, CNR, 1988.
View at: Google Scholar
B. Accordi, A. Angelucci, and G. Sirna, “Note illustrative della carta geologica d’Italia FF. 159-160,” in Illustrative notes of the geological map of Italy, FF.159-160, p. 77, Nuova Tecnica Grafica, Roma, Italy, 1967.
View at: Google Scholar
G. Devoto, “Lacustrine Pleistocene in the Lower Liri Valley,” Geologica Romana, vol. 4, pp. 291–336, 1965.
View at: Google Scholar
E. Centamore, D. Rossi, and P. Di Manna, “New data on the kinematic evolution of the Volsci Range,” Italian Journal of Geosciences, vol. 126, pp. 159–172, 2007.
View at: Google Scholar
M. Parotto and A. Praturlon, “Geological summary of the central Apennines,” in in Structural Model of Italy, vol. 90 of the Quaderni della Ricerca Scientifica, pp. 257–300, CNR, Roma, Italy, 1976.
View at: Google Scholar
G. P. Cavinato and M. Sirna, “Elementi di tettonica transpressiva lungo la linea di Atina (Lazio Meridionale), (Transpressive elements along the Atina tectonic line),” Memorie della Società Geologica Italiana, vol. 41, pp. 1179–1190, 1988.
View at: Google Scholar
M. Zalaffi, “Su alcune piccole strutture della piana di Cassino, (The small carbonatic structures of the Cassino plain),” Memorie della Società Geological Italiana, vol. 4, pp. 1–14, 1964.
View at: Google Scholar
B. Accordi, “La componente traslativa nella tettonica dellAppennino laziale-abruzzese (The contractional tectonics in the Lazio-Abruzzi Apennine),” Memorie della Società Geologica Italiana, vol. 5, pp. 355–406, 1966.
View at: Google Scholar
C. F. Boni and P. Bono, “Segnalazione di un gruppo di grandi sorgenti nel bacino del Fiume Peccia, aggluente del Garigliano, (Highlighting of huge springs in the Peccia-Garigliano River basin),” Geologica Romana, vol. 12, pp. 227–242, 1973.
View at: Google Scholar
P. Celico and V. Stanganelli, “Sulla struttura idrogeologica dei Monti di Venafro (Italia Meridionale), (The hydrogeological structure of the Venafro Mts.-Southern Italy),” Bollettino della Società dei Naturalisti, vol. 85, pp. 153–178, 1976.
View at: Google Scholar
M. Saroli, M. Lancia, M. Albano, G. Modoni, M. Moro, and G. Scarascia Mugnozza, “New geological data on the Cassino intermontane basin, central Apennines, Italy,” Rendiconti Lincei, vol. 25, no. 2, pp. 189–196, 2014.
View at: Publisher Site | Google Scholar
ISRM, “Suggested methods for the quantitative description of dicontinuities in rock masses,” International Journal of Rock Mechanics and Mining Science, vol. 15, pp. 319–368, 1978.
View at: Google Scholar
P. A. Witherspoon, J. S. Y. Wang, K. Iwai, and J. E. Gale, “Validity of cubic law for fluid flow in a deformable rock fracture,” Water Resources Research, vol. 16, no. 6, pp. 1016–1024, 1980.
View at: Publisher Site | Google Scholar
C. Louis, “A study of groundwater flow in jointed rock and its influence on the stability of rock masses,” Rock Mechanics Research Report, vol. 10, Imperial College, London, p. 90, 1964.
View at: Google Scholar
R. W. Zimmerman, A. Al-Yaarubi, C. C. Pain, and C. A. Grattoni, “Non-linear regimes of fluid flow in rock fractures,” International Journal of Rock Mechanics and Mining Sciences, vol. 41, no. 1, pp. 1–7, 2004.
View at: Publisher Site | Google Scholar
P. G. Ranjith and W. Darlington, “Nonlinear single-phase flow in real rock joints,” Water Resources Research, vol. 43, no. 9, Article ID W09502, pp. 1–9, 2007.
View at: Publisher Site | Google Scholar
C. Cherubini, C. I. Giasi, and N. Pastore, “Evidence of non-Darcy flow and non-Fickian transport in fractured media at laboratory scale,” Hydrology and Earth System Sciences, vol. 17, pp. 2599–2611, 2013.
View at: Google Scholar
J. Rutqvist and O. Stephansson, “The role of the hydromechanical coupling in fractured rock engineering,” Hydrogeology Journal, vol. 11, pp. 11–40, 2003.
View at: Google Scholar
C. Klimczak, R. A. Schultz, R. Parashar, and D. M. Reeves, “Cubic law with aperture-length correlation: Implications for network scale fluid flow,” Hydrogeology Journal, vol. 18, no. 4, pp. 851–862, 2010.
View at: Publisher Site | Google Scholar
H. Darcy, Les Fountaine publiques de la Ville de Dijon , (The public springs of the Dijion villa, pp. 647, Paris, Dalmont, 1956.
T. F. M. Chui and D. Freyberg, “The use of COMSOL for integrated hydrogeological models,” in Proceedings of the COMSOL conference, pp. 217–223, COMSOL INC., Boston, MA, USA, 2007.
View at: Google Scholar
M. Zalaffi, “Osservazioni su alcuni affioramenti di farina di dolomia al bordo meridionale del Matese, (Observations of tectonized dolostone along the Matese Mts. southern border),” Bollettino della Società Geologica Italiana, vol. 88, pp. 161–170, 1968.
View at: Google Scholar

Copyright

Copyright © 2018 M. Lancia et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

1074

Downloads

1350

Citations