# Advances, challenges and perspective in modelling the functioning of karst systems: a review

Hervé Jourde, Xiaoguang Wang

► To cite this version:

Hervé Jourde, Xiaoguang Wang. Advances, challenges and perspective in modelling the functioning of karst systems: a review. Environmental Earth Sciences, 2023, 82 (17), pp.396. 10.1007/s12665-023-11034-7. hal-04875957

HAL Id: hal-04875957

https://hal.science/hal-04875957v1

Submitted on 9 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Advances, challenges and perspective in modelling the functioning of karst systems: a review

Hervé Jourde<fnref n="1" /> · Xiaoguang Wang<fnref n="1" />,<fnref n="2" />,<fnref n="3" />

Received: 3 January 2023 / Accepted: 24 June 2023 / Published online: 11 August 2023

© The Author(s) 2023

## Abstract

We present a discussion of the state-of-the art on modelling geometrical characteristics, hydrogeological behavior and geochemical evolution of karst aquifers of meteoric origin. The considered key karst processes include: (1) the recharge processes, (2) the various hydrodynamic processes in the vadose and phreatic zones, (3) the related transport processes and (4) the speleogenesis processes. Different types of approaches for modelling geometrical characteristics of karst conduit networks are summarized. Integrated numerical studies on hydrogeological behavior of karst aquifers with functional and physically based models are then reviewed. Challenging issues in karst modelling are further discussed based on previous modelling progresses. The paper finally offers recommendations for advancing the modelling of hydrogeological behavior of karst systems and lists several open questions for future research.

**Keywords** Karst · Numerical models · Hydrodynamics · Solute transport

## Introduction

Karst aquifers are highly heterogeneous media that are composed of numerous flow paths of high permeability formed by the dissolution of carbonate rocks. The primary porosity of the carbonate rocks can play a major role in the karst genesis processes, as well as the fractures and bedding planes within the rock mass, which are often enlarged by dissolution, forming hierarchically organized flow networks. Modelling key karst processes is thus a quite challenging issue. This difficulty ultimately hinders our ability in water management of karst groundwater resources and prediction, prevention and mitigation of aquifer contamination.

The first challenge in modelling karst is to solve the complex hydrodynamics problem, i.e., the physical processes

of fluid flow and transport in karstic rocks. Karst aquifers consist of multiple media that conduct flow (e.g., conduit, fracture and matrix) with distinct hydraulic properties. Accordingly, the hydraulic conductivity may vary by several orders of magnitude over a short distance where it exhibits a sharp contrast and strong anisotropy (Kuniansky et al. 2012). In addition, the various flow features are hierarchically organized such that multi-scale flow behaviors are very common. As a result, the behavior of karst aquifers is characterized by the striking characteristic of duality (Goldscheider and Drew 2014). Although numerical techniques exist for numerical modelling of flow processes in individual compartments, the coupling of different flow processes remains challenging. The second fundamental challenge in modelling karst is the geometrical representation of complex three-dimensional karstic conduit systems. In most cases, karst forms by dissolving limestone fractures, which occur at all scales. Though some large cavities may be detected by geophysical tools, direct observation of the detailed 3D structure/geometry of karst conduit networks in the subsurface is generally not possible. The traditional approach of direct data sampling relies on cave exploration, outcrop cross-section mapping and borehole logging. Data collections by caving and outcrop mapping are time consuming and remain limited. Considering the technological advances in well-logging and in-situ testing, a detailed, continuous,

<footnote>✉ Xiaoguang.Wang
xiaoguang.wang@umontpellier.fr</footnote>

<footnote>Hervé Jourde
herve.jourde@umontpellier.fr</footnote>

<footnote n="1">¹ HydroSciences Montpellier (HSM), Univ. Montpellier, CNRS, IRD, Montpellier, France</footnote>

<footnote n="2">² College of Energy, Chengdu University of Technology, Sichuan, China</footnote>

<footnote n="3">³ Tianfu Yongxing Laboratory, Chengdu, China</footnote>
characterization of the distribution of discontinuous features (e.g., fractures, bedding planes, vugs, karstic conduits, etc.) is possible. However, the resulting data are largely limited to one-dimensional and only within the regions adjacent to the observation boreholes. Seismic surveys may be able to provide spatial dimensions about the 3D large-scale features such as faults and cavities, but their application in karst media is often limited due to a resolution threshold. Consequently, the 3D description of karst geometries requires extrapolating 1D/2D observations in small samples of the karst reservoir to the whole study domain in 3D. Therefore, the question of how to create realistic karst networks remains an unresolved issue.

Another critical issue of karst modelling concerns the hydrological and hydrochemical functioning of karst media under complex artificial or natural hydraulic conditions. Groundwater movement and contaminant migration through karst systems is predominately governed by conduits and caves. Recharge is often highly focused at sinkholes or stream networks and therefore is generally less dependent on evapotranspiration and soil cover than in other hydrogeologic settings. Discharge from karst aquifers primarily occurs at conduit-flow karst springs, and can be very large and highly variable both spatially and temporarily. The focused discharge from springs, however, represents an integrated response of water movement through different flow paths and mixing processes that occurred in various locations within the karst system. The understanding of hydrochemical changes in karst aquifers is also essential for groundwater protection and aquatic ecosystem management. Given these complexities and challenges, many numerical approaches have been proposed in the last decades to simulate the main karst processes. All methods have both advantages and disadvantages. There is no a systematic and complete numerical model that can accurately describe the flow and transport processes in karst aquifers. Relevant reviews of numerical approaches for modelling karst flow processes exist (Ghasemizadeh et al. 2012; Kovács and Sauter 2007). Their main focus was to classify karst flow models and the various numerical approaches that, in most reviews, were divided into two general classes based on the type of parameterization: the lumped and distributed parameter models. Even though distributed lumped parameters approaches were proposed to account for the spatial heterogeneity of either epikarst (Hartmann et al. 2012, 2015) or land cover land uses (Bittner et al. 2018) at the scale of the catchment’s scale, the lumped parameter models generally do not consider spatial variability of the karst system hydrodynamic properties and rather focus on the temporal variation of global precipitation–recharge–discharge behavior of the system. In contrast, the distributed parameter models aim at mimicking the distributed physical properties that controls both

the spatial and temporal variations of hydraulic response of the system. This classification was useful for a general category but may become irrelevant for the new generation of models. Indeed, some models now couple a functional approach describing the flow process in one system component and a physically based approach that captures the spatial variability of another component. In addition, novel modelling methods that introduce simple spatial discretization (often one-dimensional variation) into the functional approach also emerges in recent years. For these reasons, in this review, we do not refer to “lumped” and “distributed parameter models” but rather adopt the terminology of ‘functional approach’ and ‘physically based approach’ when classifying models.

In previous reviews, the various models were discussed separately, which did not allow for a thorough presentation of the links between them, let alone the possibility for model and data integration. Here, we attempt to provide an integrated point of view on karst modelling and to better describe the links between the various models. Accordingly, we put emphasis on the connections and combinations of existing models for an improved characterization of karst systems, as well as on the use of data from various sources through the state-of-the-art inversion techniques. The karst modelling approaches that will be discussed in this review paper are summarized in Fig. 1. Two major group of models are considered: (1) models for karst hydrology and hydraulics simulation and (2) models for karst generation and geomorphology simulation. Both forward and inverse modelling techniques are reviewed. We highlight the models and techniques that emerge in recent years (those indicated by red font color in Fig. 1). We also discuss the specific issues than each model can address, and their potential impact on the management decision making.

This paper focuses on discussing numerical modelling approaches for karst processes simulation. It is structured in the following way. “Key karst processes and related conceptual models” details the key processes occurring in epigene karst systems and highlights the importance of accurately describing them in conceptual models before capturing them in numerical models. In “Modelling of key karst processes” and “Modelling of karst conduits network”, we present recent advances in numerical modelling of key karst processes and karst network geometry, respectively. In “Modelling of karst hydrogeological functioning”, we show what can be achieved in terms of simulation of karst hydrogeological functioning with the current modelling tools and discuss the limitations of these tools. We finish this review by providing several perspectives for future modeling approaches, with the aim to better understand, characterize and reproduce key flow and transport processes that occur in natural epigenic karst systems.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /><figcaption>**Fig. 1** Summary of karst modelling approaches that will be discussed in this review paper. Black font color indicates classical approaches; red font color indicates the approaches, variables or processes more recently considered in the modelling</figcaption></figure>

## Key karst processes and related conceptual models

The first step in karst modelling is the building of an appropriate conceptual model (Fig. 2). The conceptual model provides a collective visualization of the hydrodynamic cycle: recharge, flow/storage and discharge. A well-structured conceptual model should be readily translated into a mathematical model formulated with a set of coupled differential equations. There exist several conceptual models for karst aquifers in the literature. A classical conceptual model originates from the dual-medium model developed for fractured petroleum reservoirs (Drogue and Bidaux 1992). In this model, the main flow paths correspond to fracture and karst conduits while the fissured matrix contributes to the main storage. The model was quite successful for representing flow processes in the phreatic zones of karst aquifers but it is hardly applicable for karst systems where the unsaturated zone processes are significant. By examining the spatial and temporal variation of water table level and tracer concentration, it was found that the shallow portion of the karst aquifer, i.e., the soil zone and epikarst, may play an important role in regulating flow and transport processes (Perrin et al. 2003). The concept of epikarst as a shallow, highly karstified zone below the soil zone was introduced by Mangin (1975) and Williams (1983). The epikarst was also

described as playing a major role in the storage of karst aquifer (Klimchouk 2000; Sauter 1992; Williams 1983) and in the sequence of vertical infiltration from the soil and epikarst zone, towards vadose zone and phreatic zone (Ghasemizadeh et al. 2012). Accordingly, all subsequent models adopted a similar divisional structure formed by four main components: the soil zone, epikarst zone, vadose/unsaturated and phreatic/saturated zones. Such models adequately describe the structure of karst systems formed by epigene processes.

In this paper, we focus primarily on karst systems formed by epigene and/or hypergene processes; we do not discuss hypogene karst whose genesis and development may be driven by very different mechanisms. In epigene karst systems, the hydrological interactions between meteoric water, surface water, soil water and karst groundwater are complex. During rainfall events, the meteoric water first enters the soil zone or/and the epikarst where soils are absent. Within this zone, rapid flow may occur in vertical fractures and conduits (sinkhole), while slower flow takes place in matrix. The epikarst zone can also play a major role in concentrating water flow and thus generating intense hydraulic response that may propagate through the system and be monitored at karst springs. In shallow karst system, the epikarst processes thus can dominate the spring discharge and chemical signal (Perrin et al. 2003). Although the intensity of fissures decreases rapidly with depth, concentrated water may continue to
<figure><img src="image_4.png" /><figcaption>**Fig. 2** a Conceptual model of karst systems (modified from Donsimoni et al. 2008), b example of karstic conduits developed along an inception feature, c a karst underground river</figcaption></figure>

percolate quickly to the unsaturated zone through major tectonic fractures and shafts. Within this zone, a horizontal redistribution of the infiltrated water may become important. The unsaturated zone is mainly seen as a transmissive zone connecting the epikarst and phreatic zones, although it may be highly capacitive, especially when thick vadose zone exists and control the organization of the flows with the depth in the unsaturated zone of the karst system (Emblanch et al. 2003). Flow in the unsaturated zone can also dominate the transport processes and thus modify karst spring chemographs (Perrin et al. 2007). Once the water flows in the phreatic zone, the sub-horizontal flow component predominates and, although diffusive flow may still exist, groundwater is mainly transmitted by the network of karst conduits towards springs and outlets.

Another important characteristic of karst aquifer at larger time scale is that their hydrological functioning varies constantly due to the dissolution of carbonate rocks that continuously modify the properties of the flow path networks. The evolution process due to dissolution starts when meteoric water mixes with CO₂ in the soil zone. The majority of limestone dissolution occurs in the epikarst zone and to a lesser extent in the unsaturated zone and phreatic/saturated zone. Indeed, the water often reaches carbonate saturation when transmitted to the phreatic zone. The temporal variability of karst systems hydrological functioning is a key feature that distinguish them from classical aquifers where changes can

occur over time, but not as much and as fast as in carbonate aquifers.

## Modelling of key karst processes

The terminology of ‘functional approach’ and ‘physically based approach’ is considered when referring to the models of karst processes. In this section, we focus on the key karst processes that play important roles in the system’s functioning: recharge, hydrodynamics, transport and speleogenesis. The commonly used methods for modelling these processes, either according to the ‘functional approach’ or ‘physically based approach’, are then reviewed.

### Recharge process

Estimating groundwater recharge in karst is mostly relying on the modelling of the dynamics of the soil, epikarst and the unsaturated zone (Hartmann et al. 2012). The soil and epikarst have an important hydrological role, due to their high storage capacity and thus highly regulate aquifer recharge. After infiltrating through the epikarst, the recharge is generally divided into two major portions: the vadose percolation, which is the part of the diffusive flow through fissures and matrix, and the vadose fast flow, which is associated to piston movement of water flux through shafts, karstic

<figure><img src="image_5.png" /></figure>
conduits and fractures with large openings. Groundwater recharge can be modeled by means of diverse methods. The conventional karst recharge models are based on water balance or budget calculations (Blavoux et al. 1992; Scanlon et al. 2002). The methods rely on the fact that water entry should be equal to the amount discharged plus or minus the variation in the volume of water that is stored. By combining with geographical information system (GIS), an estimate of the spatial distribution of recharge can be obtained (Andreo et al. 2008; Malard et al. 2016; Radulovic et al. 2012). The water balance approaches have several advantages: they are easy to implement, but they suffer from several important drawbacks. First, they are based on empirical interpolations and estimates which may be erroneous. Also, they rely essentially on site-specific information, so that their transferability to other regions or towards the future is limited (Hartmann et al. 2014).

An alternative approach to karst recharge modelling relies on physically based distributed methods (Ireson and Butler 2011; Kiraly et al. 1995; Perrin et al. 2003). In this method, groundwater flow equations are solved under saturated or unsaturated conditions. Along this line, the pioneering work by Kiraly et al. (1995) considers a nested 2D model to simulate the epikarst, which redistribute flows toward a underlying 3D karst network model. Outflow from the epikarst is considered as recharge inflow at sinkholes or dolines and computed as a proportion of concentrated and diffuse recharge, respectively. The assumption of partial saturation is often adopted to model the flow in the saturated zone, which is under a steady-state regime. To deal with the condition of the total saturation, a time variant model may be used, but it requires much more data (Andreo et al. 2008). In this situation, recharge estimation with the functional approaches exhibit promising results and show good potential for estimating recharge and capture realistic flow behaviors in karst aquifers.

Karst recharge models can also be distinguished as a function of the considered approach, i.e. functional or physically based, accounting for unsaturated zone flow processes. In the case of limited data availability and/or aquifer scale modelling purposes, functional recharge models are generally preferred (Fleury et al. 2007; Jukić and Denić-Jukić 2009; Tritz et al. 2011). However, contrary to physically based approaches (Andreu et al. 2011; Dvory et al. 2016; Hughes et al. 2008; Ireson and Butler 2011; Kiraly et al. 1995; Perrin et al. 2003), such functional approaches do not allow considering the spatial variability of recharge and can only provide a mean recharge, spatially averaged over the karst catchment. In the case of physically based approaches, spatial variability of soil and epikarst hydrodynamic properties can be measured in the field and considered either for the explicit modelling of water movement in the vadose zone (Kordilla et al. 2012; Rimmer and Salingar 2006) or

the estimation of groundwater recharge, distributed over the karst catchment.

## Hydrodynamic processes

### Dual-medium flow

Karst aquifers functioning is often characterized by the dual flow behavior, which results from high-conductive flow paths made of karst conduits and fractures, within a low conductive matrix. Such conceptualization of flow is also well-suited to describe flow in either the epikarst or vadose zone. However, due to multiple phase flow processes (mainly air and water) in these compartments, it was mainly considered for the modelling of single-phase flow in the saturated zone. The first attempt to simulate groundwater flow and solute transport in a karst aquifer with a hybrid discrete-continuum numerical model was proposed by Kiraly (1998), who considered conduit flow in discrete pipes embedded in a continuum porous medium with Darcian flow. This hybrid approach was implemented for the development of the CAVE (Carbonate Aquifer Void Evolution) model (Clemens et al. 1996; Liedl et al. 2003).

### Filled conduit flow/laminar and turbulent flows

Early flow models only accounted for laminar flow, although turbulent or free surface flow may occur in conduits or other enlarged flow paths. Turbulent processes were considered in MODFLOW-CFP (Hill et al. 2010; Reimann and Hill 2009) while considering the hydraulic gradient in conduits and adjusting the hydraulic conductivity so that it accounts for turbulence (Shoemaker et al. 2008) when the hydraulic gradient exceeds a critical value (Fig. 3a). Further developments were proposed by Reimann et al. (2011, 2014) to simulate conduit-associated drainable storage and unsaturated conduit flow, and account for time variable boundary conditions (Fig. 3b). MODFLOW-2005 considers laminar flow in the continuum (calculated by Darcy's equation), while CFPM1 (Conduit Flow Process Mode 1) considers laminar flow in the continuum (calculated by Darcy's equation) but also laminar or turbulent flow in a discrete pipe network (calculated by the Hagen-Poiseuille and the Colebrook-White equations). CFPM2 (Conduit Flow Process Mode 2) considers both laminar flow (calculated by Darcy's equation) and turbulent flow in the continuum (estimated by adjusting the hydraulic conductivity as described by Shoemaker et al. 2008). Finally, CTFc (Conduit-Type Flow in Continuum Cells) also considers laminar flow (calculated by Darcy's equation) and turbulent flow in the continuum (estimated a power law).

The simulations proposed by Reimann et al. (2011) consider nonlinear flow at the outlet of a karst conduit
<figure><img src="image_6.png" /><figcaption>**Fig. 3** Simulation of turbulent flow in karst aquifer. a Specific discharge as a function of hydraulic gradient for laminar flow and turbulent flow, with corresponding linear and nonlinear equations. Re<sub>c1</sub> denotes the transition from laminar to turbulent flow conditions, Re<sub>c2</sub> denotes the transition from turbulent to laminar flow conditions. b Normalized discharge simulated with different codes (MODFLOW-2005, CFPM1 CFPM2 CTFc) that consider nonlinear flow at the outlet of a karst conduit embedded in a porous matrix, for a direct recharge pulse injected in the conduit and a constant diffuse recharge in the matrix (after Reimann et al. 2011)</figcaption></figure>

embedded in a porous matrix, for a direct recharge pulse
injected in the conduit and a constant diffuse recharge in the
matrix. These simulations show that the non-consideration
of turbulent flow processes may lead to an overestimation
of the discharge at the outlet of the karst conduit (Fig. 3b).

Partially filled conduit flow/free-surface flow

The simulation of strongly variable transient flow processes,
such as flood, requires considering free-surface flow, wave
propagation in karst conduits, and changes between pressur-
ized and non-pressurized conduit flow. The Saint–Venant
equations for free-surface flow simulation in karst conduits
were implemented in the coupled overland and groundwater
flow model MODBRANCH, to simulate unsteady and non-
uniform flow processes in karst conduits. The Saint-Venant
equations allow describing unsteady flow in partially filled
conduits (de Saint-Venant 1871) and are therefore often

used to simulate discharge in conduit-flow-dominated karst
systems. In the unsaturated zone, free-surface conditions
control the conduit flow, so that Manning’s equations are
generally considered, although they may underestimate head
losses due to turbulent flows (Jeannin and Maréchal 1995).

As the cross-sectional area of karst conduits is finite, both
pressurized and non-pressurized flow may occur. Equations
were therefore adapted for pressurized conditions through
considering the Preissmann slot model in MODBRANCH
package (Reimann et al. 2011), or the Darcy–Weisbach
equation in pipe-flow models (Jeannin 2001), such as the
Storm Water Management Model (SWMM) (Campbell and
Sullivan 2002; Chen and Goldscheider 2014; Peterson and
Wicks 2006). However, the use of such models requires
detailed distributions of hydraulic parameters which are
often unavailable due to the unknown location and geometry
of karst conduits.

Transport processes

The highly heterogeneous flow fields in karst media originat-
ing from the varieties of flow paths can result in complex
transport phenomena in the system (Fig. 4a, b). Indeed, the
range of residence time of karst groundwater can be consid-
erably contrasted, spanning from hours to decades. Moreo-
ver, the mean residence time typically exhibit a strong tem-
poral variability and is highly dependent on the hydraulic
condition, which greatly controls transport and speleogen-
esis processes. Accordingly, various computational models
have been developed to simulate transport dynamics in karst
formations (Fig. 4c, d).

Most interpretation methods of artificial tracer data are
based on the classical one-dimensional advection–disper-
sion equation (ADE). Field observations have shown that
the breakthrough curves from karst aquifer are most likely
asymmetric and exhibit a sudden rising limb, a long tail,
and multiple peaks (Goldscheider et al. 2008; Moreno and
Tsang 1991). The traditional ADE performs poorly in char-
acterizing these BTCs. To handle the tailing effect, several
transport models have been proposed, which may generate
more adjustable tails (Morales et al. 2010). The Mobile-
Immobile Model (MIM) originally proposed by Coats and
Smith (1964) partitions the aquifer into mobile and immo-
bile regions. The method was later successfully applied to
laboratory column experiment data by van Genuchten and
Wierenga (1977) and field data (Atkinson et al. 1973; Brown
and Ford 1971; Field and Pinsky 2000; Göppert and Gold-
scheider 2008; Massei et al. 2007, 2006; Mull et al. 1988;
Toride and Leij 1996), and incorporated into standard sol-
ute transport codes such as the popular MT3D (Zheng and
Wang 1999). Skopp et al. (1981) proposed dual-permeability
model which modifies the MIM by considering both regions
are mobile and ADE is valid in both regions; this model
<figure><img src="image_7.png" /><figcaption>**Fig. 4** Numerical models for solute and thermal transport processes simulation. **a** simulated versus observed solute transport break-through curves (after Morales et al. 2010); **b** simulated hydro-chemo-thermograph in response to precipitation in a karst conduit (modified from Raeisi et al. 2007); **c** conceptual model for solute transport interpretation; **d** conceptual model for mass and thermal transport calculation</figcaption></figure>

is referred to as Multiple Region Advection-Dispersion (MRAD) by Majdalani et al. (2018). Another possibility for simulating the long tails of karst BTCs is the Continuous Time Random Walk (CTRW) method (Berkowitz et al. 2006). The CTRW approach visualizes solute transport as a series of particle jumps or transitions with spatially changing velocities. Recently, it has been applied to tracer test interpretation in the context of karst aquifers by Goeppert et al. (2020).

To handle the multi-peak effect, some progresses have also been made. The dual-peaked BTCs were successfully modeled by the application of the MIM (Field and Leij 2012; Maloszewski et al. 1994, 1992; Werner et al. 1997). Field and Leij (2012) compared the ability of dual-advection dispersion equation (DADE) and the weighted sum advection-dispersion equation (WSADE) to produce dual-peaked BTCs. Recently, Wang et al. (2020) applied the DADE model and successfully fitted the BTCs generated by controlled laboratory conduit network experiments.

Theoretical studies have been conducted to investigate the problem of solute transport in a system of a single conduit embedded in the background porous media (Li et al. 2008).

Analytical solutions have been derived for these problems. Although the hydrodynamic interaction between the conduit and matrix is considered, the mass exchange at the conduit-matrix interface is assumed constant and does not vary according to the hydrodynamic condition (Li and Field 2016; Li and Loper 2011). Zoghbi and Basha (2020) extended this model to incorporate a fully coupled conduit-matrix exchange behavior, accounting for the change in hydraulic gradient/flow and transport between the conduit and matrix. The application of this model to field data obtained on two karst aquifers demonstrated its effectiveness in simulating karst BTCs and estimating key transport parameters.

When one attempts to model spatially variable output phenomena, physically based models may provide more accurate representation of some of the transport processes (Birk et al. 2006, 2005; Faulkner et al. 2009; Oehlmann et al. 2015). In these models, spatially variable velocities were estimated first from flow simulations and then used in the transport simulation. Thus, depending on the parameterization of flow models, the type of distributed transport model includes equivalent porous medium (Zheng and Wang 2012).
1999), dual-porosity (Gerke and Van Genuchten 1993), discrete conduit network (Peterson and Wicks 2006) and hybrid (Faulkner et al. 2009; Oehlmann et al. 2015) approaches.

Many characteristics/processes of karst transport dynamics have not been explicitly considered in the existing karst transport models. For instance, inertial flow and dynamic transition between laminar and turbulent flow regimes that frequently occur in karst conduits may have important effects on the transport process but are difficult to handle with traditional numerical models. To simulate these processes, computational fluid dynamics (CFD) simulations need to be performed (Hauns et al. 2001). However, the use of CFD models is limited for the understanding of general transport behavior based on scenario models, as they require significant computational resources. Therefore, their application for regional-scale karst problems is prohibitive.

The lattice Boltzmann methods (LBMs) offers an alternative numerical approach for simulating the effect of these complex hydrodynamics on conduit transport properties while consuming acceptable computational resources (Sukop et al. 2008; Anwar and Sukop 2009).

## Speleogenesis

Earliest speleogenesis models were proposed by Dreybrodt (1990, 1996) and Palmer (1991). These one-dimensional (1D) models allowed demonstrating the importance of non-linear dissolution processes and the existence of positive feedback loops while considering the evolution of a single fracture (Fig. 5a). They also allowed calculating dissolution length scales and gave rise to the concept of breakthrough time; the breakthrough time corresponds to the time when

<figure><img src="image_8.png" /><figcaption>**Fig. 5** Speleogenesis simulations: a aperture enlargement due to dissolution in a single fracture (after Dreybrodt and Gabrovšek 2018); b karst conduit's aperture and pattern at a mature stage of karstification (after Hubinger and Birk 2011); c wormhole formation in a single fracture highlighted by the normalized undersaturation $(C_{sat} - C)/C_{sat}$ contour plots (after Szymczak and Ladd 2011) plotted at successive times (from left to right); $t = 25\tau$, $t = 50\tau$, $t = 75\tau$, and $t = 100\tau$ ($\tau$ is the</figcaption></figure>

time it takes for the initial fracture aperture at the inlet to double); d concentration and flow discharge evolution during karst generation in a 20 m long fracture for linear and non-linear dissolution kinetics, in which $q_0$ is the initial flux in the fracture, $l_p$ is the penetration length, and $C_{sat}$ is the concentration at saturation (after Szymczak and Ladd 2011)
a steep increase in the flow rate is observed, due to the increasing widening of the fracture aperture at the outlet Dreybrodt (1990). These models were the basis for more complex 2D models (Dreybrodt 2012; Groves and Howard 1994; Kaufmann et al 2012; Siemens and Dreybrodt 1998), with more complex dynamics related to 2D networks of fractures (Fig. 5b). In such cases, the region with dissolution of first-order kinetics penetrates downstream, as the discharge through the fracture network increases (Hubinger and Birk 2011), which results in a main branch of the karst conduits network fanning out at the discharge side (Fig. 5b).

At the scale of a single fracture, 2D models allowed demonstrating that aperture heterogeneity can be at the origin of preferential flow paths and thus generate earlier breakthrough than in the 1D case (Hanna and Rajaram 1998). Flow exchange between matrix and fracture (Bauer et al. 2003; Liedl et al. 2003; Kaufmann and Braun 2000), as well as flow exchange between prominent fractures and fissures (Gabrovšek et al. 2004), can also generate earlier breakthrough. Fingering of the dissolution front within a single fracture (Fig. 5c) was further demonstrated as the consequence of the high instability of the dissolution front propagation (Szymczak and Ladd 2011). In such modelling, breakthrough time and concentration profiles, as well as flow discharge evolution during speleogenesis highly depend on the considered dissolution kinetics (Fig. 5d). Further modelling was performed to study the influence of mixing corrosion and the competition between different flow paths and the effect of CO₂ sources (Gabrovšek et al. 2000), flank margin caves formation (Dreybrodt et al. 2009; Dreybrodt and Romanov 2007), karstification around dam sites (Dreybrodt et al. 2002), as well as buoyant convection (Chaudhuri et al. 2013). 2D fracture networks dissolution models allowed analyzing the horizontal and vertical evolution of caves (Gabrovšek and Dreybrodt 2001; Groves and Howard 1994; Siemens and Dreybrodt 1998). Dual-porosity models were further developed (Kaufmann and Braun 2000; Liedl et al. 2003), while considering flow through discrete conduits and through the porous rock matrix as well as exchange between both media. 3D model accounting for both speleogenesis and landscape evolution was proposed by Kaufmann (2009) to assess karst evolution over time. A complete review of fracture network speleogenesis models has been proposed by Dreybrodt et al. (2005) and Dreybrodt and Gabrovšek (2018).

The topology of the initial flow path network, which often consists of discrete fracture networks, may also influence significantly the dissolution pattern (Aliouache et al. 2019). When the flow is parallel to the persisting fracture orientation and becomes channelized (Fig. 6a–c), the initial hydrogeological setting highly controls the dissolution pattern. On the contrary, intermittent interruptions of the positive feedback loop between flow and dissolution occur when the

flow direction is normal to the persisting fracture orientation (Fig. 6d–f), and the breakthrough is delayed.

## Summary with a view towards hydrogeological modelling

From the above discussion, it can be seen that the modelling of key karst processes can be achieved by either functional or physically based approaches, which have important differences in the conceptualization of karst systems. Table 1 presents a comparison of the functional and physically based models. The functional modelling approach allows bypassing the geometrical complexity using equivalent behavioral parameters that reflect the overall hydrologic and hydrodynamic functioning of the karst system. On the other hand, physically based modelling scheme treats the system as an assemblage of interacting individual components and allows for the explicit integration of complex physical laws for hydraulic processes and interactions. However, some of the model parameters (e.g., conduit friction coefficient, exchange coefficient) may need to be obtained by numerical calibration rather than from direct measurements from the field. Furthermore, the computational time for solving physically based models can be substantially larger than that for functional models.

Since the important characteristics of hydraulic processes are preserved in the conventional functional and physically based modelling, some hydrogeological characteristics of karst aquifers can be reproduced. However, simplifications are commonly adopted in order for the problem to be tractable. Many aspects may be improved for a better modelling of karst hydrological functioning. To improve recharge simulation methods, a promising direction for future research is the development of hybrid models that assimilate the advantages of both functional and physically based models. Some of the current recharge estimation models have already exhibited some of such characteristics. For example, Bradley et al. (2010) combined several flow paths and storages within a functional model. In the subsequent stochastic recharge model proposed by Hartmann et al. (2012, 2015), statistical distribution functions were used to represent the soil and epikarst heterogeneities controlling storage and flow paths. Most physically based models that were developed in the past focused on the behavior of the karst aquifer in the phreatic zone. Despite recent progress in the identification of karst vadose zone dynamics, the appropriate characterization of process parameters and related conceptual models remain to be initiated for an efficient implementation of the typical processes of karst vadose zone into respective physically based models. Complete models should therefore account for unsaturated flow in the vadose zone, which would require considering film flow in wide fractures and shafts, and saturated storm flow in smaller fractures. These models should
<figure><img src="image_9.png" /><figcaption>**Fig. 6** Incipient karst formation in the same initial natural fracture network (modified after Aliouache et al. 2019) considering directional flow along X and Y directions. **a, d** flow rate distribution at breakthrough for x- and y-direction flow, respectively; **b, e** aperture distribution at breakthrough for x- and y-direction flow, respectively.</figcaption></figure>

<table><caption>**Table 1** Comparison of functional and physically based approaches for karst process modelling</caption><thead><tr><td>Numerical models</td><td>Key inputs</td><td>Strengths</td><td>Limitations</td></tr></thead><tbody><tr><td>Functional approaches</td><td>Effective mathematical relations; conceptual model structure; equivalent hydraulic properties</td><td>(i) Simple model structure<br>(ii) Efficient computation<br>(iii) Suitability for basin-scale applications</td><td>(i) Difficulty in deriving hydrodynamic properties<br>(ii) Oversimplification of conduit hydrodynamics<br>(iii) No consideration of complex flow processes related to conduit-matrix interaction</td></tr><tr><td>Physically based approaches</td><td>Geometric and hydraulic properties for conduits</td><td>(i) Explicit integration of karst networks<br>(ii) Explicit representation of hydraulic in conduits<br>(iii) Integration of localized boundary conditions (e.g. infiltration, pumping...)</td><td>(i) Large computation time<br>(ii) Limited data on conduit network geometry</td></tr></tbody></table>

also consider the dampening effect of the epikarst related to temporary storage and retardation of precipitation. In general, the simulation of coupled saturated and unsaturated flow is still a challenge in hydrogeology, especially for karst systems (Reimann et al. 2011; de Rooij et al. 2013). The dual continuum approach may be adopted for such a purpose (Sauter et al. 2006; Teutsch and Sauter 1991). However, a

sacrifice for reality of aquifer geometry must be accepted (Kordilla et al. 2012). A promising way of addressing the difficulties in modelling complex karst hydraulics is the employment of Lattice Boltzmann Methods (LBMs). The LBM-based karst models are particularly attractive as they are capable of simulating inertial flows, incorporating complex conduit geometries, and solving coupled Darcian and
Navier-Stokes flow problems. They also allow solving sol-
ute transport problems characterized by laminar and eddy
flows in conduits or Darcian flow in porous matrix, and tack-
ling multiphase flow and transport problems (Anwar and
Sukop 2009; Sukop et al. 2008).

Modelling of karst conduits network

In karst aquifers, the majority of flow is conducted through
the conduit network. The distribution of conduits within the
aquifer influences the hydrodynamic functioning of karst
aquifer. Thus, they should be ideally incorporated into
numerical models, particularly when the physical mecha-
nism of flow and transport in karst systems is concerned.
The rigorous analysis of parameter uncertainty requires a
large number of model-runs with an ensemble of plausible
karst conduit networks. Therefore, methods for generating
realistic karst conduit networks are needed. In practice, karst
networks may be generated (1) based on field-mapped infor-
mation of conduit geometry; (2) by stochastic simulation or
(3) by sophisticated transport models accounting for com-
plex geochemical reactions and dissolution processes. In the
literature of karst hydrology, the stochastic and geochemi-
cal approaches are also referred to as structure-imitating
and process-imitating approaches, following the classifica-
tion of numerical methods proposed originally by Kolter-
mann and Gorelick (1996) for simulating heterogeneity in
sedimentary deposits (Borghi et al. 2012; Pardo-Iguzquiza
et al. 2011). The structure-imitating approaches relies on
a statistical approach aiming to reproduce the structure of
conduit network without considering physical and chemical
processes, while the process-imitating approaches rely on

dissolution processes that control conduit network evolution,
i.e., speleogenesis.

Field-mapped (and interpolated) karst conduits
networks

Conduit profile and pattern may be determined from detailed
speleological investigations (e.g., Jouves et al. 2017) which
allow reconstructing the three-dimensional geometry of
explored cave systems. During such a field survey, three-
dimensional information may be acquired at a sequence of
topographic stations before being processed numerically to
obtain discrete conduits with elliptical or rectangular sec-
tion shapes. A 3D conduit network may then be constructed/
simulated by interpolating the two-dimensional sections.
This task can be solved with common speleological pro-
grams such as Survex, Visual Topo or GHTopo. Rendering
with more realistic section shapes is also available in these
programs. LiDAR technologies also allow to precisely map
karst conduits and caves (Jaillet et al. 2011). This method is
adapted to explore individual conduits.

As only a small portion of conduits can be surveyed,
modelling techniques are required to reconstruct 3D karst
conduit networks based on one- and two-dimensional meas-
urements. Henrion et al. (2010) proposed a new method
(object-distance simulation method—ODSIM) to gener-
ate cylindrical and lens shaped envelops along an irregu-
lar curve considered as the skeleton of the karst conduit, in
order to reproduce its 3D geometry. Rongier et al. (2014)
further improved this method by considering a fast-marching
method to account for the different geological properties
influencing conduit shapes (Fig. 7).

The method produced conduit network with more realis-
tic section shapes. The use of these generated complex 3D

<figure><img src="image_10.png" /><figcaption>**Fig. 7** Simulation of 3D conduit networks (after Rongier et al. 2014)</figcaption></figure>
conduit network for flow and transport simulation at aquifer scale remains to be demonstrated, considering the high computational cost of 3D simulations to couple geometry modelling with modern simulators of flow and transport.

## Stochastically-generated karst conduits networks

The complete geometrical, geomorphological and topological properties of karst systems cannot be fully assessed due to the difficulty to identify karst conduit networks at all scales. This favored the development of stochastic approaches that rely on the hierarchy and scaling of natural karst systems, which can be adequately characterized by fractal geometry or power law models (Hendrick and Renard 2016a; Pardo-Iguzquiza et al. 2011).

The spatial organization of natural karst networks, which often exhibit a hierarchical structure as highlighted by extensive field observations, can be characterized by a fractal dimension, *D*, which quantifies the manner of how fractals spread and cluster in the Euclidean space. The reported fractal dimension for natural karst networks range from 1.043 to 1.67 (Hendrick and Renard 2016b; Pardo-Igúzquiza et al. 2012).

The fractal characteristics may be incorporated in stochastically-generated karst networks. One example is the diffusion-limited aggregation approach of Pardo-Iguzquiza et al. (2012). In this approach, the geometry of conduit segments (e.g., orientation and length) is first simulated by sampling with replacement from known conduits (conduit templates). Then, the topology (connectivity between conduit segments) is generated stochastically using the diffusion-limited aggregation method. Although the method is able to incorporate the fractal character of karst conduit geometry, the generated connectivity is not explicitly based on statistical input parameters. Thus, it is not guaranteed that the resulting network will have a realistic connectivity pattern, even though this approach was successfully applied for either conditional or non-conditional stochastic simulation of karst conduit networks.

Apart from the fractal approach, other numerical methods may be used to simulate karst. Borghi et al. (2012) used a fast-marching algorithm to simulate the evolution of complex karst patterns formed under different boundary conditions (Fig. 8a–c). Jaquet et al. (2004) developed a lattice-gas cellular automation approach to simulate the diameter growth of a conduit network and the spatiotemporal evolution of a multiphase speleogenesis system (Fig. 8d). Ronayne (2013) used a non-looping invasion percolation model (Stark 1991) that allows generating dendritic (i.e., non-looping) connectivity patterns. Although dendritic geometries are generally related to stream networks, such connectivity patterns can be generated by a focused recharge, which results in dendritic karst conduit networks (Palmer 1991).

## Geochemically-simulated karst conduits networks

Geochemically based models, which benefited from the increasing knowledge of speleogenesis processes, consider carbonate dissolution kinetics to simulate conduit growth and karst network evolution in response to hydrologic and geomorphic variations. Such model can be used to reproduce natural karst patterns in response to paleo-hydrological condition. One of the approaches is the hydrochemical model of de Rooij and Graham (2017), which has been applied by Henson et al. (2018) to understand the major factors, which contribute to the formation of large springs. However, such methods favor the development of hierarchical networks, which are one particular geometry among many other possible geometries of karst networks (Palmer 1991). The computational cost of these speleogenesis models is also very high and they are very difficult to constrain with field data.

## Summary with a view towards hydrogeological modelling

Each of the three types of conduit network generation models have their strengths and limitations (Table 2). Field-mapped karst networks models allow capturing many realistic features of karst conduits; however these models become limited to characterize phreatic conduit network in deep rocks and caves with small diameter that human cannot access. The stochastic conduit network approaches have the merit to be efficient and also applicable for 3D problems. However, they require a detailed statistical analysis of conduit properties and accurately deriving control parameters (e.g., conduit sizes and fractal dimension). Such a practice requires a large amount of data at all scales, and is therefore barely difficult to apply for realistic problems. The process-imitating karst models, which may capture some genesis characteristics of natural karst, are sensitive to the inferred paleo-flow fields and boundary conditions. These models need to be further improved to take into account geological, thermal and geomechanical effects to reproduce actual karst systems. Recent advances in karst modelling have proposed such evolution but still remain unapplied to realistic karst system.

To generate better karstic conduit networks, a future research direction that is attracting much effort is the development of process-oriented numerical models for catchment scale problems. Several studies adopted such a modelling strategy and have already shown promising results. For example, de Rooij and Graham (2017) used the spring locations to condition a 2D conduit network generated by a hydrochemical model. Aliouache et al. (2019) applied numerical karst genesis modelling in a geologically-mapped fracture network following dissolution kinetic law. Their conduit growth is regulated by the

<figure><img src="image_11.png" /></figure>
<figure><img src="image_12.png" /><figcaption>**Fig. 8** Stochastic conduit network generation with fast marching method (a–c after Borghi et al. 2012)) and with the cellular automation approach (d after Jaquet et al. 2004)</figcaption></figure>

realistic topological properties of natural fractures. Jiang et al. (2022) proposed a hydro-chemical model to simulate incipient karst generation in 3D jointed layered carbonates and highlighted the role of different flow regimes as well as the initial penetration length in the generation of distinct karst morphologies (i.e. pipe-like, stripe-like and sheet-like karst features). Such karst genesis models may thus supply additional information for outlining the geometry of conduit network in realistic modelling practice.

Modelling of karst hydrogeological functioning

Various analytical and numerical approaches have been used to model the hydrogeological functioning of karst aquifers. The preference for a functional or physically based modelling scheme often depends on the objective of modelling application, scale of problem, system complexity and data availability. The functional approach
<table><caption>**Table 2** Comparison of different methods for karst network modelling</caption><thead><tr><td>Numerical models</td><td>Essential input data</td><td>Advantages</td><td>Limitations</td></tr></thead><tbody><tr><td>Geological karst networks</td><td>Speleological characterization, borehole imaging, LIDAR scan or geophysical survey</td><td>(i) Deterministic description of a karst system<br>(ii) Preservation of geological realisms</td><td>(i) Limited feasibility for application to basin scale karst characterization<br>(ii) Constraints from speleological and geophysical investigation extent and scale<br>(iii) Difficult to model 3D karst geometry</td></tr><tr><td>Stochastic karst networks</td><td>Statistical data of conduit lengths, orientation, position, shape and their connection type</td><td>(i) Simplicity and efficiency in network generation<br>(ii) Applicability for both 2D and 3D<br>(iii) Applicability for multiple scales</td><td>(i) Simplification of conduit geometry and topology<br>(ii) Uncertainties in statistical parameters<br>(iii) Requirement of multiple realizations<br>(iv) Negligence of physical processes</td></tr><tr><td>Geochemically-simulated karst networks</td><td>Paleo-hydraulic conditions, rock properties, reaction kinetics</td><td>(i) Linking geometry with physical mechanisms<br>(ii) Correlation between different fracture attributes</td><td>(i) Large computational time<br>(ii) Uncertainties in input properties as well as paleo-hydrological and hydrochemical conditions<br>(iii) Negligence of thermal and mechanical processes</td></tr></tbody></table>

has the advantage of simplicity and its great efficiency to model large-scale system functioning. In contrast, the physically based models can straightforwardly incorporate spatial heterogeneities and model complex hydrodynamic processes. Depending on the discretization type, the physically based models may be further classified as equivalent porous medium (EPM), double continuum (DC), discrete conduit network (DCN) and combined discrete-continuum (CDC) models as proposed in former reviews (Hartmann et al. 2014; Kovács and Sauter 2007). In this section, commonly used models for simulating hydrogeological behaviors of karst systems, including spring responses to natural forcing (i.e., rainfall or recharge) and borehole responses to artificial forcing (i.e., hydraulic tests) will be reviewed.

## Karst spring responses

One peculiarity of karst aquifers is that most of the groundwater discharge occurs through a small number of outlets, and sometimes through a single spring, which can be characterized by hydrographs or chemographs controlled by the whole water moving through the aquifer. Numerical studies therefore focused on the simulation of karst spring responses, with either functional or physically based model, to understand system's dynamics.

## Functional approach

Functional approaches model the karst system's functioning or evolutionary dynamics at the scale of the entire system, from the infiltration of meteoric water into the soil zone

and epikarst to discharge at springs. The dynamic processes, such as infiltration and storage in various aquifer components are modeled with linear or nonlinear empirical relationships. Due to this lumped nature of model structure, the parameters of the functional models only represent effective properties of the whole karst system and may not bear a physical meaning. For this reason, the model parameters cannot be quantified by field measurements and only the dynamic of discharge at springs can give insight about the base flow recession that can be further considered in the model. Instead, such parameters may need to be determined by an indirect process of numerical calibration. Calibrations of functional karst models are important task, particularly when they are used for real applications. The commonly used data for calibrating functional karst models include continuously monitored flow and physicochemical data at karst springs and to a much less extent, in-situ measurements of water levels in monitoring boreholes (Ladouche et al. 2014; Mazzilli et al. 2019, 2011).

Functional models can be applied for single event (or recession) analysis or time series analysis. Recession analysis is widely used to understand the hydrodynamic functioning of karst system (e.g., Amit et al. 2002; Bonacci 1993; Chang et al. 2015; Cinkus et al. 2021; Fiorillo 2011; Florea and Vacher 2006; Forkasiewicz and Paloc 1967; Kovács et al. 2005; Malík and Vojtková, 2012; Mangin 1975). The time series analysis is commonly used to reveal insight about the internal functioning of karst aquifers by analyzing the precipitation-discharge relationship (e.g., Bailly-Comte et al. 2008; Labat et al. 2000; Larocque et al. 1998; Padilla and Pulido-Bosch 1995; Panagopoulos and Lambrakis 2006;
Rahnemaei et al. 2005). Detailed mathematical analysis has been continuously and progressively inserted to increase the efficiency of this method (Gárfias-Soliz et al. 2010; Jukić and Denić-Jukić, 2011; Mathevet et al. 2004).

The early stage of karst functional models focused on reproducing the discharge behavior caused by the duality of flow regime, one of the most important concepts in karst hydrology. The components of fast conduit flow and diffusive matrix flow are simulated using two separate reservoirs or tanks, which are combined to generate the discharge. These models were extended to consider the conduit-matrix coupling by formulating an exchange relationship between the reservoirs. For instance, in the models proposed by Geyer et al. (2008), matrix may recharge the conduit reservoir if a pre-defined hydraulic threshold is reached. Functional models that include ‘two-way’ flow exchanges between matrix and conduits have also been proposed (Mazzilli et al. 2011). Another alternative way of simulating the dual flow characteristics of karst aquifer is to couple a matrix reservoir and a routing function representing the rapid flow component, which is transferred directly to spring discharge (Le Moine et al. 2008).

To account for the influence of soil and epikarst in regulating recharge and providing storage, a separation function that divides recharge into matrix and conduit reservoirs according to a pre-defined fraction ratio has been introduced in classical functional models (Charlier et al. 2012; Fleury et al. 2009). Another approach developed to account for the particular flow dynamics related to soil and epikarst consists in integrating a hysteresis scheme that allows activating and deactivating concentrated recharge (Tritz et al. 2011).

A new class of functional karst models that integrate the hydrochemical modelling into original hydrodynamic analyses have also emerged. Thorough incorporating extra hydrochemical data, such as chemical and heat tracers (Anderson 2005; Covington et al. 2011; Luhmann et al. 2011; Martin and Dean 1999; Saar 2011), water chemistry (Barberá and Andreo 2012; Cholet et al. 2015; Filippini et al. 2018; Mudarra et al. 2014; Raeisi and Karami 2011; Toran et al. 2006), water conductivity (Massei et al. 2006; Wang et al. 2020) in addition to hydrographs, it has been shown that these integrated functional models provides a more reliable characterization of the system functioning of karst aquifers. Hydrochemical information in functional models (Charlier et al. 2012; Hartmann 2013; Hartmann et al. 2016, 2017; Sivelle et al. 2022) were introduced to combine information inferred from both hydrochemical signal and hydrodynamic signal.

The simplicity of lumped model structure may be regarded as an advantageous aspect, because detailed information, such as geometrical and hydraulic properties for dominant processes, is often not available when analyzing complex karst systems, for which physically based modelling

is greatly limited. However, it is still very important to continue improving the realism and accuracy of functional karst models, because the spring responses (i.e., the spatial and temporal variations of discharge and physicochemical parameters) are governed by the recharge pattern and physical processes integrated within the model. Developments are needed towards (1) a more thorough characterization of the flow and transport dynamics respecting more details of real karst system, such as the complex role of epikarst storage (Perrin et al. 2003) in the regulation of recharge, mixing of waters from rainfall and aquifer tributaries, as well as chemical reaction, and (2) a more rational and efficient implementation of model structure respecting the identified dominant processes and accounting for the coupling between some of these processes. For more extensive discussion on the use of functional models for modelling karst spring responses and other improved implementations, the reader is referred to Jeannin and Sauter (1998), Kovács and Sauter (2007) and Hartmann et al. (2014).

## Physically based approach

Extensive studies have been conducted to model karst spring dynamics and understand their relation to spatial heterogeneities and physical processes (Doummar et al. 2012; Koridilla et al. 2012; Liedl et al. 2003; Oehlmann et al. 2015; Xanke et al. 2016). The increased knowledge of the mechanisms of conduit genesis and hydraulics promoted the development of physically based karst models that incorporate the mechanisms governing the flow and transport processes and simulate karst spring signal evolution as a hydrological response to spatial heterogeneity.

Various numerical discretization methods have been proposed and the one based on equivalent porous medium (EPM) concept is frequently used (Ghasemizadeh et al. 2015). In the EPM model, the rock matrix, fissures and conduit networks are represented by a continuum domain with equivalent hydraulic conductivities assigned to individual grid cells. The inclusion of conduit features can be achieved through using lines of drain cells and internal sinks (Quinn et al. 2006; Quinn and Tomasko 2000) or lines of high hydraulic conductivity cells (Lindgren et al. 2005). Since the modelling objective time dependent behaviors such as spring discharge and water table variation, transient simulations are generally required. Storage parameters, such as specific yield, may be assessed by core testing techniques or in some studies were treated as variables to be determined by a calibration process. The transient recharge is difficult to quantify as it depends on various meteoric factors such as rainfall and evapotranspiration as well as reservoir characteristics, such as elevation, runoff and recharge area. Water budget methods are often adopted to identify the recharge type and estimate the recharge parameter values. The initial
and boundary conditions of the karst groundwater system may be identified and approximated by hydrological and hydrogeological characterization of field observations and monitored data (Ghasemizadeh et al. 2012).

In general, the EPM models are more suitable for karst modelling at regional scales where the salient influence of conduit hydraulics is averaged out. Scanlon et al. (2003) developed an EPM model for simulating regional scale groundwater flow in the Edwards aquifer in Austin, Texas and found that the model was able to simulate temporal variations in spring discharge and to reproduce the effect of pumping on discharge with fairly good accuracy. To evaluate the potential impact of conduit flow, Saller et al. (2013) integrated high conductivity zones, which represents karst conduits, into the EPM models developed originally by Long and Putnam (2009) for the Madison aquifer. They observed that the inclusion of high conductivity zones resulted in a better modelling of transient water table but did not improve the capacity of the original EPM model in karst spring discharge simulations. The classic physically based numerical models based on the EPM formulation can capture some aspects of the karst hydrological behavior, but were inadequate to simulate solute transport. Worthington and Ford (2009) created a regional model for the limestone aquifer at Mammoth Cave (Kentucky, USA) and investigated the difference in model response with and without the presence of high-permeability conduits. The authors found that the EPM model, even while considering conduit effects, cannot reproduce the tracer test results.

Apart from the EPM approach, discrete conduit network (DCN) model was used to successfully simulate flow of a karst system in Kentucky, USA (Thrailkill 1974). Jeannin (2001) modelled the hydraulic response of the Holloch cave (Switzerland) with a DCN model, and successfully related it to the local cave geometry. In subsequent studies, the storm water management model (SWMM) that is based on the DCN formulation has been adopted for simulating flow and transport in conduit-dominated karst systems (Campbell and Sullivan 1999; Chen and Goldscheider 2014; Peterson and Wicks 2006; Wu et al. 2008). One of the main focuses of such models was the examination of the influence of conduit hydraulic properties. Peterson and Wicks (2006) performed a parameter sensitivity test using SWMM, and reported that a slight change in the relative conduit roughness resulted in a substantial variation in the simulated spring discharge and solute transport. Because SWMM solve one-dimensional flow and transport equations, the calculation is very efficient even for solving regional problem. It can be coupled to a calibration system to study model parameter sensitivity and fitted to field measurements of long-term time series such as spring discharge and hydrochemical variation (Chen and Goldscheider 2014). Such a calibration capability also offers a possibility of inferring essential geometrical structure of

karst network through inversions of multi-tracer test including transient flow (Sivelle et al. 2020). The DCN methods focus on representing conduit hydraulics but do not address the flow and transport processes in the matrix adequately.

To overcome such a deficiency, a dual continuum model may be adopted. The dual continuum approach is able to simulate the dual flow behavior of karst system without the need to incorporate an explicit representation of complex conduit networks (Teutsch and Sauter 1991). The integration of saturated and unsaturated flow regimes with the dual continuum conceptualization has also been implemented in the HydroGeoSphere model (Therrien et al. 2010).

## Hybrid approach

The strong hydraulic coupling effect of various zones in a karst system can affect significantly the simulated aquifer functioning behavior. For example, the recharge function has a strong impact on spring hydrographs in karst aquifers (Hartmann 2013; Kiraly et al. 1995; Kovács et al. 2005). One promising method is thus to include recharge through the use of a functional reservoir model. In this way, the role of the epikarst zone in regulating and storing rainfalls is considered and the output of the reservoir model gives the recharge input for the physically based model component (Borghi et al. 2016; Chen and Goldscheider 2014). This method can be implemented in all kinds of karst simulation model and has a potential benefit to simulate a real karst aquifer. Oftenderger et al. (2015) extended the simple hybrid models and used a complex hydrological model to compute recharge distribution and the runoff of precipitation. However, this kind of model is very difficult to parameterize given that the number of measurement stations is often limited. On the other hand, the parameter values of recharge models are difficult to identify. One encouraging way to address these difficulties is to incorporate some of the recharge model parameters into an inverse framework and to estimate their values by an optimization process (Mohammadi et al. 2018; Kavousi et al. 2020). The calibrated parameters need to be back-checked with field data to confirm if they are physically meaningful. Thus, their application to real problems is only restricted to well monitored karst systems.

## Numerical representation of hydraulic tests

Karst aquifer characterization is typically based on spring response interpretation. Recently, it has been extended by systematic analysis of pumping tests. Experimental field sites equipped with boreholes and in-situ monitoring sensors may provide useful information about karst aquifer properties at local scales. However, they do not allow determining hydrodynamic properties at catchment scale because of the restricted extent of the investigated area due to limited
pumping rates. In particular cases with exceptional pumping rate of several hundred liters per second, large-scale characterization of karst aquifers has been made possible (Maréchal et al. 2008; Noushabadi et al. 2011; Dausse et al. 2019).

Numerical models with distributed hydrodynamic parameter fields in both matrix and karst conduits can also be used to characterize karst aquifers at large scale. Discrete-continuum models represent a suitable approach to simulate karst aquifers (e.g., Kaufmann 2009; Kiraly 1998, 2003; Sauter et al. 2006) while considering the matrix interaction and highly anisotropic hydraulic parameter fields. Liedl et al. (2003) proposed a model comprising both a discrete-continuum for the simulation of laminar and turbulent pipe flow, and a continuum that represents the matrix. As this model applies steady-state pipe flow equations, it cannot capture the transient water storage within the conduit. To overcome this deficiency in representing realistic storage dynamics, the Saint-Venant equations for variably filled conduit flow can be solved (Reimann et al. 2014). Other important extensions for better reproducing the realistic hydrodynamic process in karstic media include coupling surface flow (de Rooij 2008) and variably saturated matrix flow (Rooij et al. 2013). Giese et al. (2017) used the discrete conduit-continuum model of Reimann et al. (2014) to investigate the effect of different boundary conditions as well as sink/source terms for idealized two-dimensional mixed karst aquifer systems. Giese

et al. (2018) applied the same model and examined the effect of conduit flow regimes (i.e., laminar and non-laminar flow) on the drawdown in conduit and matrix with varied conduit hydraulic properties (i.e., conduit diameter and relative roughness). They found that the turbulent flow may generate important influence on pumping test interpretation for small conduit diameters, while for mature conduit system where the diameter is sufficiently large, the laminar flow equation may be applicable for simulating karstic conduit hydrodynamics even under the condition of fully-developed turbulent flow.

The numerical representation of pumping tests together with the borehole hydraulic head measurements can be used to infer the spatial distribution and hydraulic properties of karst conduit networks. Such an inversion approach is referred to as hydraulic tomography, which can be used to identify hydraulic connectivity of karst aquifers. While considering a quasi-Newton inverse method, Wang et al. (2016) estimated the transmissivity field and thus the major karst conduits and their connectivity at the scale of an experimental field site, which however highly depends on a priori information provided on the inversions (Fig. 9a, b). A stochastic Newton inverse method was further applied by these authors to quantify the uncertainty in the inverted transmissivity fields and address the issue of uniqueness of the inversion (Wang et al. 2017). On the same experimental field site, the cellular automata-based deterministic inversion

<figure><img src="image_13.png" /><figcaption>**Fig. 9** Karst conduit networks constructed by simultaneously inverting multiple cross-hole pumping test data; **a** inversion results with a constant initial hydraulic conductivity distribution; **b** inversion results</figcaption></figure>

with an initial hydraulic conductivity distribution inferred from interference tests; **c** inversion results obtained by the CADI method (after Fischer et al. 2017; Wang et al. 2016)
(CADI) method (Fischer et al. 2017) allowed identifying a high transmissivity discrete conduit network embedded within a lower transmissivity matrix (Fig. 9c). The CADI method minimizes the difference between the observed and modeled hydraulic data by distributing the hydraulic properties along linear structures and iteratively modifying the conduit network geometry.

Further research considered harmonic pumping technique to assess both cross-borehole connectivity and karst conduit networks geometry (Fischer et al. 2018). With this technique, one can assess the degree of connectivity between the pumping well and the observation wells, which depends on the phase offset of the monitored responses, and the conductivity of main flow paths, which is controlled by the amplitude of the monitored response. High conductivity flow paths between observation wells and the pumping well can be identified with high frequency pumping, while dual connectivity between the observation wells and the pumping well (propagation occurring both in the karst conduit networks and in the matrix) can be highlighted by low-frequency pumping. The hierarchical flow path network (i.e., karst conduit, fracture, and matrix) can thus be assessed while considering various harmonic pumping tests frequencies. Tomographic inversions of harmonic pumping tests are also much faster (about 10 times) than the one using constant pumping rates and they require at least two times as few tests, which makes this method a powerful approach to enhance the interpretation of karst systems.

## Comments on the modelling of karst hydrogeological functioning

In the conventional hydrogeological modelling, the karst is regarded as a particular and distinct type of aquifer. Besides, most of these hydrogeological modelling studies focused on the phreatic zone, with flow and transport dynamics in the epikarst and vadose zone either completely ignored or greatly simplified. Such treatments may be acceptable for modelling karst systems with thin and relatively homogeneous vadose zone, for which the effect of flow redistribution and storage of the epikarst and vadose zone may be negligible or approximated by simple recharge function. However, for most karst systems, the vadose zone is quite thick and affects the storage and recharge processes greatly. Under high flow conditions, fluid transfers are likely to occur at the conduit-matrix interface driven by high local pressure gradients. The matrix of the epikarst and vadose zone may thus be filled by conduit flows; while under low flow conditions, water stored in the epikarst and vadose matrix may be released to conduits and thus become the primary source of spring discharge. This regulation process can create retention and releasing of fluid and contaminant and thus affect greatly the hydrological and hydrochemical properties

of karst aquifers and groundwater, as observed in the field (Perrin et al. 2003). Furthermore, the irregular roughness of karst conduits, which are often assumed to be straight line (in 2D) or cylindrical pipe (in 3D), may accommodate enhanced pressure loss under high flow rates, particularly at sites of diameter contractions, which has important effects on flow and transport processes. Hence, towards an improved hydrogeological modelling of karst aquifers, the capability of precisely modelling the coupled saturated-unsaturated hydraulic behavior is prerequisite, for which the important characteristics of flow and transport processes, e.g., conduit-matrix interaction, flow regime transition, need to be appropriately treated.

The development of calibration and validation procedures is a crucial issue for physically based karst models based on field measurements, especially when they are applied for practical purposes. The genetically-generated conduit networks can be examined by comparing the geometrical properties of a simulated conduit network with those of the conduit pattern mapped by speleology exploration. Boreholes and underground research laboratories (Emblanch et al. 2003; Perineau et al. 2011) provide also valuable information on karst conduit sizes that intersect borehole and tunnel walls. These data might be of major importance for the calibration of karst network geometries. Hydrodynamic information, such as water table variations and pressure derivatives in response to pumping (single-well and cross-wells tests), as well as chemical concentration monitored during tracer tests (Ronayne 2013) can also help tuning the hydrological performance of karst conduit network models. Furthermore, the continuous physicochemical data monitored at springs and water table variations recorded in boreholes offer useful references for the calibration of karst models used for hydrological simulation at catchment scale. The newly-developed hydraulic tomography technique can be employed to calculate best-fit hydraulic diameter distributions for karst flow models based on an optimization process using multiple objective functions.

## Outstanding issues and open questions

Modelling the hydrogeological and geochemical characteristics of karst aquifers is a challenging issue. In previous modelling studies, a reasonable understanding of the main components of karst aquifers and their hydrogeological functions has been established and efforts in developing quantitative techniques for modelling key dynamic processes have shown promising results. However, many important issues and scientific questions still exist. These issues/scientific questions, that would be valuable to investigate, offer opportunities for future research in the field of karst modelling.
## Issues of functional karst models

The functional models based on empirical relationships have shown great adaptability and strength in simulating the hydrogeological behavior under different contexts and scales (Hartmann et al. 2015). However, it remains relevant to continue improving their realism and accuracy. Enhancements to the current modelling schemes by assimilating distinctive karst processes are worthy of further exploration. Some functional models have been proposed to condition the model to flow measurements in conduit or matrix (Fleury et al. 2009; Schmidt et al. 2014), to discharge at multiple springs (Rimmer and Salingar 2006), but also to consider the variations of recharge and storage in soil and epikarst in time (Tritz et al. 2011) and space (Hartmann et al. 2012). Understanding the sensitivity and better assess the parametric uncertainty of models (Cousquer and Jourde 2022) or their structural uncertainty is also an important aspect for future research. The predictive capacities of the functional models can be increased if field measurements in addition to spring discharge data are incorporated into the calibration procedure. Some studies have already been attempted to test the feasibility and consequence of integrating hydrochemical (Charlier et al. 2012; Hartmann 2013; Hartmann et al. 2016; Pinault et al. 2001) and gravimetric (Mazzilli et al. 2013) information. The joint use of data from different type of monitoring should help further reducing the ambiguity of functional models. Future studies should go on this direction.

## Issues of physically based karst models

In previous studies, simulation results obtained from physically based karst models have shown consistency with well-controlled laboratory experiments and/or field tests. However, there remains a number of difficulties in applying physically based numerical approaches for modelling karst processes. And a number of important issues still exist in numerical modelling of the hydrogeological behavior of karst aquifers using physically based models.

The first task will be the development of realistic hydrodynamic models that can represent the important transient flow and transport characteristics of natural karst systems (e.g., coupled free surface-Darcian flow systems, transition between saturated and unsaturated flow regimes, dynamic conduit-matrix interactions), and that can be used afterwards to capture the key hydrogeological and hydrochemical characteristics of karst systems (e.g., dual flow regimes, hydraulic thresholds, matrix retention and releasing of solutes).

The geometry of the created karst conduit networks needs to be consistent with the speleogenesis principles and field mapped conduit statistics so that the hierarchically-organized fluid flow can be properly modelled. Developments

are needed towards (1) a more detailed characterization of the underlying statistics, such as multifractals for which more than one scaling exponent is expected (Hendrick and Renard 2016b; Pardo-Iguzquiza et al. 2011); (2) a more precise and efficient model to generate karst conduit network with more details to reproduce the diversity of individual conduit shapes and morphology (Rongier et al. 2014), the topological complexity in karst conduit populations (Pardo-Igúzquiza et al. 2019) and the correlation between conduit geometrical attributes (Jouves et al. 2017); (3) extrapolation of current 2D karst networks generated by the stochastic methods to 3D while preserving geostatistical, topological and hydrological properties.

Another major issue is to propose adequate upscaling approaches to hydrodynamic models to evaluate the large-scale behavior of karst systems. The challenge for such models will be to preserve the hydrological and hydrochemical characteristics observed at the local scale. Some techniques were used to construct large-scale models on the basis of either upscaled permeabilities of local cells (Popov et al. 2009) or multiscale integration methods (Correia et al. 2020). In the future, more efforts should be put in other aspects, such as the development of more advanced parallelization schemes (De Rooij 2019), the simulation of multiphase flow effects on transport processes and, importantly, extension of the models to complex 3D systems.

## General open questions: linking functional and physically based models

As can be seen from above, the hydrogeological modelling of karst aquifers can be achieved by functional or physically based approaches, which differ significantly in the conceptualization of geological media. Because the functional models treat the system as a black box and attempts to simulate the hydrogeological processes using analytical relationships between input signal and output response, the computational demand is much lower comparing to physically based models. They may be easily coupled to an optimization system to calibrate model parameters and predict future aquifer's functioning behavior. However, an obvious drawback is that the calibrated model parameters may not bear any physical meaning. The continuum modelling scheme (i.e., EPM) mainly reflects the dynamic process in a karst from a more overarching view and simplifies the geometrical complexity using equivalent material properties derived from upscaling techniques. Such a treatment allows the continuum models to deal with large-scale problems with great efficiency. However, it cannot consider the effects of conduit flow dynamics, sharp hydraulic contrast in conduit and matrix, complex conduit-matrix interaction, and conduit surface dissolution kinetics.
On the other hand, the physically based models regard the system as a group of interacting individual flow components and permits the integration of complex spatial variability for rock materials and realistic physical processes. Various physically based models have been proposed to handle complexities at different level. A physically based model can be established at a specific investigation scale without presuming the form of recharge-discharge transfer relationships. However, the resolution of discontinuous problems may require considerably larger computational time than the one needed for continuum models, which greatly limits the use of such models in calibration and inversion studies.

A direction for future research is to develop an integrated modelling approach by further exploring the link between functional and physically based models. In practical applications, physically based models that integrate all available data about the modelling site could be established first. It may then be used to drive functional proxies for the physically based models, where the actual flow and transport physics is replaced by a set of reservoir structures characterized by empirical relationships that produce similar system responses. As such, the numerical efficiency will be greatly improved and the functional proxies could be adopted for calibration and sampling purposes at a much lower computational cost. To explore this link, a comparison of a large number of simulations needs to be conducted. For the task to be successful, advanced computing techniques, such as parallelization on GPU and/or CPU, as well as modern data analysis schemes, such as machine-learning and deep learning algorithms, may need to be implemented.

The improved understanding of the link between functional and physically based model parameters may also help to propose more robust karst typology. Efforts to classification of karst system have been attempted previously with only functional approaches. The main focus was on looking for similarities in discharge responses. Although distinct discharge signatures have been recognized, difficulties still exist in extrapolating the results because the system's signature cannot be related to system properties (i.e., geometrical and hydrodynamic properties) which are not considered in functional models; as a consequence, the karst topology is not generalized. The use of physically based models to infer structure and parameters of functional models will infuse physical meanings into the individual or groups of parameters. The development of such a physically based karst typology is important for decision makings at the global scale under constraints of climate change and anthropogenic forcing as well as to improve our general understanding of karst hydrology.

## Concluding remarks

This review paper began by presenting a conceptual model that identifies the key components and processes of karst aquifers. An overview for simulating the key flow and transport processes in epikarst, vadose and phreatic zones were then reviewed. Different techniques for modelling the geometrical and topological characteristics of karst conduit networks were also surveyed. When they properly account for flow and transport dynamics and accurately capture the dynamic regime transition of the studied karst system, modelling approaches can be adopted to reduce uncertainties in the management of karst groundwater. Up to now, many hydraulic behaviors and their activating thresholds that were identified by field observations and measurements have not been incorporated into physically based karst models. Proper parameterization and efficient implementation strategies for these processes are needed to derive advanced karst models for simulating more realistic and complete hydrogeological responses. Other major issues that need to be addressed in the future are those related to development of advanced karst models integrating the complete hydraulic processes, but also the creation of realistic karst networks, the upscaling of physically based models for catchment-scale practical applications, the study of coupled bio-hydro-chemical effects on karst aquifer evolution, linking functional and physically based approaches. Another challenge will be the transfer of models to ungauged karst catchments on the basis of specific models linked to the typology of karst hydrological response.

**Acknowledgements** This work was supported by the European Commission through the Partnership for Research and Innovation in the Mediterranean Area (PRIMA) program under Horizon 2020 (KARMA project, Grant agreement number 01DH19022A). X. W. is supported by the National Key Research and Development Program of China (no. 2020YFC1808300) and National Natural Science Foundation of China (no. 42102300).

**Author contributions** HJ and XW wrote the manuscript text and prepare the figures. All authors reviewed the manuscript.

**Funding** The authors have not disclosed any funding.

## Declarations

**Conflict of interest** The authors have no relevant financial or non-financial interests to disclose.

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will
need to obtain permission directly from the copyright holder. To view a
copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

References

Aliouache M, Wang X, Jourde H, Huang Z, Yao J (2019) Incipient karst formation in carbonate rocks: Influence of fracture network topology. J Hydrol 575:824–837. https://doi.org/10.1016/j.jhydrol.2019.05.082

Amit H, Lyakhovsky V, Katz A, Starinsky A, Burg A (2002) Interpretation of spring recession curves. Groundwater 40:543–551

Anderson MP (2005) Heat as a ground water tracer. Groundwater 43:951–968. https://doi.org/10.1111/j.1745-6584.2005.00052.x

Andreo B, Vías J, Durán JJ, Jiménez P, López-Geta JA, Carrasco F (2008) Methodology for groundwater recharge assessment in carbonate aquifers: application to pilot sites in southern Spain. Hydrogeol J 16:911–925

Andreu JM, Alcalá FJ, Vallejos A, Pulido-Bosch A (2011) Recharge to mountainous carbonated aquifers in SE Spain: different approaches and new challenges. J Arid Environ 75:1262–1270

Anwar S, Sukop MC (2009) Lattice Boltzmann models for flow and transport in saturated karst. Groundwater 47:401–413

Atkinson TC, Smith DI, Lavis JJ, Witaker RJ (1973) Experiments in tracing underground waters in limestones. J Hydrol 19:323–349

Bailly-Comte V, Jourde H, Roesch A, Pistre S, Batiot-Guilhe C (2008) Time series analyses for Karst/River interactions assessment: case of the Coulazou river (southern France). J Hydrol 349:98–114. https://doi.org/10.1016/j.jhydrol.2007.10.028

Barberá JA, Andreo B (2012) Functioning of a karst aquifer from S Spain under highly variable climate conditions, deduced from hydrochemical records. Environ Earth Sci 65:2337–2349

Bauer S, Liedl R, Sauter M (2003) Modelling of karst aquifer genesis: Influence of exchange flow. Water Resour Res. https://doi.org/10.1029/2003WR002218

Berkowitz B, Cortis A, Dentz M, Scher H (2006) Modeling non-Fickian transport in geological formations as a continuous time random walk. Rev Geophys 44:RG2003. https://doi.org/10.1029/2005RG000178

Birk S, Geyer T, Liedl R, Sauter M (2005) Process-based interpretation of tracer tests in carbonate aquifers. Groundwater 43:381–388

Birk S, Liedl R, Sauter M (2006) Karst spring responses examined by process-based modelling. Groundwater 44:832–836

Bittner D, Narany TS, Kohl B, Disse M, Chiogna G (2018) Modeling the hydrological impact of land use change in a dolomite-dominated karst system. J Hydrol 567:267–279. https://doi.org/10.1016/j.jhydrol.2018.10.017

Blavoux B, Mudry J, Puig J-M (1992) The karst system of the Fontaine de Vaucluse (Southeastern France). Environ Geol Water Sci 19:215–225

Bonacci O (1993) Karst springs hydrographs as indicators of karst aquifers. Hydrol Sci J 38:51–62

Borghi A, Renard P, Jenni S (2012) A pseudo-genetic stochastic model to generate karstic networks. J Hydrol 414:516–529

Borghi A, Renard P, Cornaton F (2016) Can one identify karst conduit networks geometry and properties from hydraulic and tracer test data? Adv Water Resour 90:99–115

Bradley C, Baker A, Jex CN, Leng MJ (2010) Hydrological uncertainties in the modelling of cave drip-water δ18O and the implications for stalagmite palaeoclimate reconstructions. Quat Sci Rev 29:2201–2214

Brown MC, Ford DC (1971) Quantitative tracer method for investigation of karst hydrological system. Trans Cave Res Group G.B. 13(1):37–51

Campbell CW, Sullivan SM (1999) Simulating time-varying cave flow and water levels using the storm water management model (SWMM). In: Presented at the hydrogeology and engineering geology of sinkholes and karst, pp 383–388

Campbell CW, Sullivan SM (2002) Simulating time-varying cave flow and water levels using the storm water management model. Eng Geol 65:133–139

Chang Y, Wu J, Liu L (2015) Effects of the conduit network on the spring hydrograph of the karst aquifer. J Hydrol 527:517–530

Charlier J-B, Bertrand C, Mudry J (2012) Conceptual hydrogeological model of flow and transport of dissolved organic carbon in a small Jura karst system. J Hydrol 460:52–64

Charlier J-B, Moussa R, Bailly-Comte V, Danneville L, Desprats J-F, Ladouche B, Marchandise A (2015) Use of a flood-routing model to assess lateral flows in a karstic stream: implications to the hydrogeological functioning of the Grands Causses area (Tarn River, Southern France). Environ Earth Sci 74:7605–7616

Chaudhuri A, Rajaram H, Wiswanathan H (2013) Early-stage hypogene karstification in a mountain hydrologic system: A coupled thermohydrochemical model incorporating buoyant convection. Water Resour Res 49:5880–5899. https://doi.org/10.1002/wrcr.20427

Chen Z, Goldscheider N (2014) Modelling spatially and temporally varied hydraulic behavior of a folded karst system with dominant conduit drainage at catchment scale, Hochifen-Gottesacker, Alps. J Hydrol 514:41–52. https://doi.org/10.1016/j.jhydrol.2014.04.005

Cholet C, Steinmann M, Charlier J-B, Denimal S (2015) Comparative study of the physicochemical response of two karst systems during contrasting flood events in the French Jura Mountains. In: Andreo B, Carrasco F, Durán JJ, Jiménez P, LaMoreaux JW (eds) Hydrogeological and environmental investigations in karst systems, environmental earth sciences. Springer, Berlin, pp 1–9. https://doi.org/10.1007/978-3-642-17435-3_1

Cinkus G, Mazzilli N, Jourde H (2021) Identification of relevant indicators for the assessment of karst systems hydrological functioning: proposal of a new classification. J Hydrol 603:127006

Clemens T, Hückinghaus D, Sauter M, Liedl R, Teutsch G (1996) A combined continuum and discrete network reactive transport model for the simulation of karst development. IAHS Publication, Wallingford, pp 309–320

Coats K, Smith B (1964) Dead-end pore volume and dispersion in porous media. Soc Petrol Eng J 4:73–84

Cornaton F, Perrochet P (2002) Analytical 1D dual-porosity equivalent solutions to 3D discrete single-continuum models. Application to karstic spring hydrograph modelling. J Hydrol 262:165–176

Correia MG, von Hohendorff Filho JC, Schiozer DJ (2020) Multiscale integration for karst-reservoir flow-simulation models. SPE Reserv Eval Eng 23:518–533. https://doi.org/10.2118/195545-PA

Cousquer Y, Jourde H (2022) Reducing uncertainty of karst aquifer modelling with complementary hydrological observations for the sustainable management of groundwater resources. J Hydrol. https://doi.org/10.1016/j.jhydrol.2022.128130

Covington MD, Luhmann AJ, Gabrovšek F, Saar MO, Wicks CM (2011) Mechanisms of heat exchange between water and rock in karst conduits. Water Resour Res. https://doi.org/10.1029/2011WR010683

Dausse A, Leonardi V, Jourde H (2019) Hydraulic characterization and identification of flow-bearing structures based on multi-scale investigations applied to the Lez karst aquifer. J Hydrol Reg Stud. https://doi.org/10.1016/j.ejrh.2019.100627
de Rooij R (2008) Towards improved numerical modelling of karst aquifers: coupling turbulent conduit flow and laminar matrix flow under variably saturated conditions. Université de Neuchâtel, Neuchâtel

de Saint-Venant AJ-C (1871) Théorie du mouvement non-permanent des eaux, avec application aux crues des rivières et à l'introduction des marées dans leur lit. CR Acad Sci Paris 73:5

De Rooij R (2019) Improving accuracy and efficiency in discrete-con tinuum karst models. Environ Earth Sci 78:115. https://doi.org/ 10.1007/s12665-019-8115-5

de Rooij R, Graham W (2017) Generation of complex karstic conduits networks with a hydrochemical model. Water Resour Res 53:6993-7011

de Rooij R, Graham W, Maxwell RM (2013) A particle-tracking scheme for simulating pathlines in coupled surface-subsurface flows. Adv Water Resour 52:7-18

Donsimoni M, Berthier H, Martin J, Nachbaur A (2008) Détermination de la présence de karsts sous les levées domaniales du bassin de la Loire et réduction du risque d'effondrement de laigue lié à la présence de ces conduits souterrains naturels (Tranche 1: analyse géologique et hydrogéologique). Rapport technique, BRGM

Doumarj M, Sauter M, Geyer T (2012) Simulation of flow processes in a large scale karst system with an integrated catchment model (Mike She)—identification of relevant parameters influencing spring discharge. J Hydrol 426:112–123

Dreybrodt W (1990) The role of dissolution kinetics in the develop-
ment of karst aquifers in limestone: a model simulation of karst
evolution. J Geol 98(5):639–655

Dreybrodt W (1996) Principles of early development of karst conduits under natural and man-made conditions revealed by mathematical analysis of numerical models. Water Resour Res 32:2923-2935

Dreybrodt W (2012) Processes in karst systems: physics, chemistry, and geology. Springer Science & Business Media, Berlin

Dreybrodt W, Romanov D, Gabrovsek F (2002) Karstification below dam sites: a model of increasing leakage from reservoirs. Environ Geol 42:518-524

Dreybrodt W, Gabrovšek F, Romanov D (2005) Processes of speleogenesis: a modelling approach, Carsologica 4. Karst Research Institute at ZRC SAZU. ZRC Publishing, Postojna, Ljubljana

Dreybrodt W, Romanov D, Kaufmann G (2009) Evolution of isolated caves in porous limestone by mixing of phreatic water and surface water at the water table of unconfined aquifers: a model approach. J Hydrol 376:200-208

Dreybrodt W, Gabrovšek F (2018) Dynamics of wormhole formation in fractured karst aquifers. Hydrol Earth Syst Sci Discuss. https://doi.org/10.5194/hess-2018-275

Dreybrodt W, Romanov D (2007) Time scales in the evolution of solu-
tion porosity in porous coastal carbonate aquifers by mixing cor-
rosion in the saltwater-freshwater transition zone. Acta Carsol
36(1):25–34

Drogue C, Bidaux P (1992) Structural and hydrogeological origin of tower karst in southern China (Lijiang plain in the Guilin region). Zeitschrift für Geomorphologie 25–36

Dvory NZ, Livshitz Y, Kuznetsov M, Adar E, Yakirevich A (2016) The effect of hydrogeological conditions on variability and dynamic of groundwater recharge in a carbonate aquifer at local scale. J Hydrol 535:480–494

Emblanch C, Zuppi GM, Mudry J, Blavoux B, Batiot C (2003) Carbon 13 of TDIC to quantify the role of the unsaturated zone: the example of the Vaucluse karst systems (Southeastern France). J Hydrol 279:262–274

Faulkner J, Hu BX, Kish S, Hua F (2009) Laboratory analog and numerical study of groundwater flow and solute transport in a

karst aquifer with conduit and matrix domains. J Contam Hydrol
110:34–44

Field MS, Leij FJ (2012) Solute transport in solution conduits exhibit-
ing multi-peaked breakthrough curves. J Hydrol 440:26–35

Field MS, Pinsky PF (2000) A two-region nonequilibrium model for solute transport in solution conduits in karstic aquifers. J Contam Hydrol 44:329-351

Filippini M, Squarzoni G, De Waele J, Fiorucci A, Vigna B, Grillo B, Riva A, Rossetti S, Zini L, Casagrande G, Stumpp C, Gargini A (2018) Differentiated spring behavior under changing hydrological conditions in an alpine karst aquifer. J Hydrol 556:572–584. https://doi.org/10.1016/j.jhydrol.2017.11.040

Filipponi M (2009) Spatial analysis of karst conduit networks and determination of parameters controlling the speleogenesis along preferential lithostratigraphic horizons. Ph.D. Thesis, Verlag nicht ermittelbar. https://doi.org/10.5075/epfl-thesis-4376

Fiorillo F (2011) Tank-reservoir drainage as a simulation of the recession limb of karst spring hydrographs. Hydrogeol J 19:1009-1019

Fischer P, Jardani A, Wang X, Jourde H, Lecoq N (2017) Identifying flow networks in a karstified aquifer by application of the cellular automata-based deterministic inversion method (Lez aquifer, France). Water Resour Res 53:10508-10522

Fischer P, Jardani A, Jourde H, Cardiff M, Wang X, Chedeville S, Lecoq N (2018) Harmonic pumping tomography applied to image the hydraulic properties and interpret the connectivity of a karstic and fractured aquifer (Lez aquifer, France). Adv Water Resour 119:227-244

Fleury P, Plagnes V, Bakalowicz M (2007) Modelling of the func-
tioning of karst aquifers with a reservoir model: Application to
Fontaine de Vaucluse (South of France). J Hydrol 345:38-49

Fleury P, Ladouche B, Conroux Y, Jourde H, Dörfliger N (2009) Mod-
elling the hydrologic functions of a karst aquifer under active
water management—the Lez spring. J Hydrol 365:235–243

Florea LJ, Vacher H (2006) Springflow hydrographs: eogenetic vs. telogenetic karst. Groundwater 44:352–361

Forkasiewicz J, Paloc H (1967) Le régime de tarissement de la Foux de la Vis. Etude préliminaire. In: Hydrologie Des Roches Fis-surées. Actes Du Colloque de Dubrovnik, IASH and UNESCO, Belgium, pp 2123-2228

Gabrovšek F, Dreybrodt W (2001) A model of the early evolution of karst aquifers in limestone in the dimensions of length and depth. J Hydrol 240:206-224

Gabrovšek F, Menne B, Dreybrodt W (2000) A model of early evolu-
tion of karst conduits affected by subterranean CO 2 sources.
Environ Geol 39:531-543

Gabrovšek F, Romanov D, Dreybrodt W (2004) Early karstification in a dual-fracture aquifer: the role of exchange flow between prominent fractures and a dense net of fissures. J Hydrol 299:45-66

Gárfias-Soliz J, Llanos-Acebo H, Martel R (2010) Time series and stochastic analyses to study the hydrodynamic characteristics of karstic aquifers. Hydrol Process 24:300–316. https://doi.org/10.1002/hyp.7487

Gerke HH, Van Genuchten MT (1993) A dual-porosity model for simu-
lating the preferential movement of water and solutes in struc-
tured porous media. Water Resour Res 29:305–319

Geyer T, Birk S, Liedl R, Sauter M (2008) Quantification of temporal distribution of recharge in karst systems from spring hydrographs. J Hydrol 348:452–463

Ghasemizadeh R, Hellweger F, Butscher C, Padilla I, Vesper D, Field M, Alshawabkeh A (2012) Groundwater flow and transport modelling of karst aquifers, with particular reference to the North Coast Limestone aquifer system of Puerto Rico. Hydrogeol J 20:1441–1461

Ghasemizadeh R, Yu X, Butscher C, Hellweger F, Padilla I, Alshawabkeh A (2015) Equivalent porous media (EPM) simulation
of groundwater hydraulics and contaminant transport in karst aquifers. PLoS ONE 10:e0138954. https://doi.org/10.1371/journal.pone.0138954

Giese M, Reimann T, Liedl R, Maréchal J-C, Sauter M (2017) Application of the flow dimension concept for numerical drawdown data analyses in mixed-flow karst systems. Hydrogeol J 25:799–811

Giese M, Reimann T, Bailly-Comte V, Maréchal J-C, Sauter M, Geyer T (2018) Turbulent and laminar flow in karst conduits under unsteady flow conditions: interpretation of pumping tests by discrete conduit-continuum modelling. Water Resour Res 54:1918–1933

Goeppert N, Goldscheider N, Berkowitz B (2020) Experimental and modelling evidence of kilometer-scale anomalous tracer transport in an alpine karst aquifer. Water Res 178:115755

Goldscheider N, Drew D (2014) Methods in karst hydrogeology. CRC Press, Boca Raton

Goldscheider N, Meiman J, Pronk M, Smart C (2008) Tracer tests in karst hydrogeology and speleology. Int J Speleol 37:27–40

Goeppert N, Goldscheider N (2008) Solute and colloid transport in karst conduits under low-and high-flow conditions. Groundwater 46:61–68

Groves CG, Howard AD (1994) Early development of karst systems: 1. Preferential flow path enlargement under laminar flow. Water Resour Res 30:2837–2846. https://doi.org/10.1029/94WR01303

Hanna RB, Rajaram H (1998) Influence of aperture variability on dissolutional growth of fissures in Karst Formations. Water Resour Res 34:2843–2853. https://doi.org/10.1029/98WR01528

Hartmann A, Lange J, Weiler M, Arbel Y, Greenbaum N, Zehe E (2012) A new approach to model the spatial and temporal variability of recharge to karst aquifers. Hydrol Earth Syst Sci 16:2219–2231

Hartmann A, Weiler M, Wagener T, Lange J, Kralik M, Humer F, Mizyed N, Rimmer A, Barberá JA, Andreo B (2013) Process-based karst modelling to relate hydrodynamic and hydrochemical characteristics to system properties. Hydrol Earth Syst Sci 17:3305–3321

Hartmann A, Goldscheider N, Wagener T, Lange J, Weiler M (2014) Karst water resources in a changing world: review of hydrological modelling approaches. Rev Geophys 52:218–242

Hartmann A, Kobler J, Kralik M, Dirnböck T, Humer F, Weiler M (2016) Model-aided quantification of dissolved carbon and nitrogen release after windthrow disturbance in an Austrian karst system. Biogeosciences 13:159–174

Hartmann A, Gleeson T, Wada Y, Wagener T (2017) Enhanced groundwater recharge rates and altered recharge sensitivity to climate variability through subsurface heterogeneity. Proc Natl Acad Sci 114:2842–2847

Hartmann A, Gleeson T, Rosolem R, Pianosi F, Wada Y, Wagener T (2015) A large-scale simulation model to assess karstic groundwater recharge over Europe and the Mediterranean. Geosci Model Dev 8:1729–1746. https://doi.org/10.5194/gmd-8-1729-2015

Hartmann A (2013) Modelling karst hydrology and hydrochemistry at different scales and in different climates considering uncertainty. Ph.D. thesis, Verlag nicht ermittelbar

Hauns M, Jeannin P-Y, Atteia O (2001) Dispersion, retardation and scale effect in tracer breakthrough curves in karst conduits. J Hydrol 241:177–193

Hendrick M, Renard P (2016a) Fractal dimension, walk dimension and conductivity exponent of karst networks around Tulum. Front Phys 4:27

Hendrick M, Renard P (2016b) Subnetworks of percolation backbones to model karst systems around Tulum, Mexico. Front Phys 4:43

Henrion V, Caumon G, Cherpeau N (2010) ODSIM: an object-distance simulation method for conditioning complex natural structures. Math Geosci 42:911–924

Henson WR, de Rooij R, Graham W (2018) What makes a first-magnitude spring?: Global sensitivity analysis of a speleogenesis model to gain insight into karst network and spring genesis. Water Resour Res 54:7417–7434

Hill ME, Stewart MT, Martin A (2010) Evaluation of the MODFLOW-2005 conduit flow process. Groundwater 48:549–559

Hubinger B, Birk, S., (2011) Influence of initial heterogeneities and recharge limitations on the evolution of aperture distributions in carbonate aquifers. Hydrol Earth Syst Sci 15:3715–3729. https://doi.org/10.5194/hess-15-3715-2011

Hughes A, Mansour M, Robins N (2008) Evaluation of distributed recharge in an upland semi-arid karst system: the West Bank Mountain Aquifer, Middle East. Hydrogeol J 16:845–854

Ireson AM, Butler AP (2011) Controls on preferential recharge to Chalk aquifers. J Hydrol 398:109–123

Jaillet S, Sadier B, Arnaud J, Azéma M, Boche E, Cailhol D, Filipponi M, Le Roux P, Varrel E (2011) Topographie, représentation et analyse morphologique 3D de drains, de conduits et de parois du Karst. Collection EDYTEM. Cahiers De Géographie 12:119–130

Jaquet O, Siegel P, Klubertanz G, Benabderrhamane H (2004) Stochastic discrete model of karstic networks. Adv Water Resour 27:751–760

Jeannin P-Y (2001) Modelling flow in phreatic and epiphreatic karst conduits in the Hölloch cave (Muotatal, Switzerland). Water Resour Res 37:191–200

Jeannin P-Y, Sauter M (1998) Analysis of karst hydrodynamic behaviour using global approaches: a review. Bull Hydrogéol (neuchâtel) 16:31–48

Jeannin P-Y, Maréchal JC (1995) Lois de pertes de charge dans les conduits karstiques: base théorique et observations. Bulletin d'hydrogéologie 149–176

Jiang C, Wang X, Pu S, Jourde H (2022) Incipient karst generation in jointed layered carbonates: insights from three-dimensional hydro-chemical simulations. J Hydrol. https://doi.org/10.1016/j.jhydrol.2022.127831

Jouves J, Viseur S, Arfib B, Baudement C, Camus H, Collon P, Guglielmi Y (2017) Speleogenesis, geometry, and topology of caves: a quantitative study of 3D karst conduits. Geomorphology 298:86–106

Jukić D, Denić-Jukić V (2009) Groundwater balance estimation in karst by using a conceptual rainfall-runoff model. J Hydrol 373:302–315

Jukić D, Denić-Jukić V (2011) Partial spectral analysis of hydrological time series. J Hydrol 400:223–233. https://doi.org/10.1016/j.jhydrol.2011.01.044

Kaufmann G (2009) Modelling karst geomorphology on different time scales. Geomorphol Recent Dev Surf Subsurf Karst Geomorphol 106:62–77. https://doi.org/10.1016/j.geomorph.2008.09.016

Kaufmann G, Braun J (2000) Karst Aquifer evolution in fractured, porous rocks. Water Resour Res 36:1381–1391. https://doi.org/10.1029/1999WR900356

Kaufmann G, Romanov D, Dreybrodt W (2012) Modeling of karst aquifers, encyclopedia of caves, 2nd edn. Academic Press, New York, pp 508–512. https://doi.org/10.1016/B978-0-12-383832-2.00073-6 (editor(s): William B. White, David C. Culver, ISBN 9780123838322)

Kavousi A, Reimann T, Liedl R, Raeisi E (2020) Karst aquifer characterization by inverse application of MODFLOW-2005 CFPv2 discrete-continuum flow and transport model. J Hydrol 587:124922

Kiraly L (1998) Modelling karst aquifers by the combined discrete channel and continuum approach. Bulletin Du Centre D'hydrogéologie 16:77–98

Kiraly L (2003) Karstification and groundwater flow. Speleogenesis and Evolution of Karst Aquifers 1:155–192
Kiraly L, Perrochet P, Rossier Y (1995) Effect of the epikarst on the hydrograph of karst springs. Bulletin Du Centre D'hydrogéologie 14:199-220

Klimchouk A (2000) The formation of epikarst and its role in vadose speleogenesis. Speleogenesis: evolution of karst aquifers. National Speleological Society, Huntsville

Koltermann CE, Gorelick SM (1996) Heterogeneity in sedimentary deposits: a review of structure-imitating, process-imitating, and descriptive approaches. Water Resour Res 32:2617-2658

Kordilla J, Sauter M, Reimann T, Geyer T (2012) Simulation of saturated and unsaturated flow in karst systems at catchment scale using a double continuum approach. Hydrol Earth Syst Sci 16:3909-3923

Kovács A, Sauter M (2007) Modelling karst hydrodynamics. Front Karst Res 13:107

Kovács A, Perrochet P, Király L, Jeannin P-Y (2005) A quantitative method for the characterisation of karst aquifers based on spring hydrograph analysis. J Hydrol 303:152-164. https://doi.org/10.1016/j.jhydrol.2004.08.023

Kuniansky EL, Bellino JC, Dixon J (2012) Transmissivity of the upper Floridan aquifer in Florida and parts of Georgia, South Carolina, and Alabama. US Department of the Interior, US Geological Survey

Labat D, Ababou R, Mangin A (2000) Rainfall-runoff relations for karstic springs. Part II: continuous wavelet and discrete orthogonal multiresolution analyses. J Hydrol 238:149-178

Ladouche B, Maréchal J-C, Dorfliger N (2014) Semi-distributed lumped model of a karst system under active management. J Hydrol 509:215-230

Larocque M, Mangin A, Razack M, Banton O (1998) Contribution of correlation and spectral analyses to the regional study of a large karst aquifer (Charente, France). J Hydrol 205:217-231

Le Moine N, Andréassian V, Mathevet T (2008) Confronting surface-and groundwater balances on the La Rochefoucault-Touvre karstic system (Charente, France). Water Resour Res 44:W03403. https://doi.org/10.1029/2007WR005984

Lesparre N, Boyle A, Grychtol B, Cabrera J, Marteau J, Adler A (2016) Electrical resistivity imaging in transmission between surface and underground tunnel for fault characterization. J Appl Geophys 128:163-178. https://doi.org/10.1016/j.jappgeo.2016.03.004

Li G, Loper DE, Kung R (2008) Contaminant sequestration in karstic aquifers: experiments and quantification. Water Resour Res 44:W02429. https://doi.org/10.1029/2006WR005797

Li G, Field MS (2016) Solute migration from the aquifer matrix into a solution conduit and the reverse. Groundwater 54:699-708. https://doi.org/10.1111/gwat.12416

Li G, Loper DE (2011) Transport, dilution, and dispersion of contaminant in a leaky karst conduit. Transp Porous Med 88:31-43. https://doi.org/10.1007/s11242-011-9721-1

Liedl R, Sauter M, Hückinghaus D, Clemens T, Teutsch G (2003) Simulation of the development of karst aquifers using a coupled continuum pipe flow model. Water Resour Res 39:1057. https://doi.org/10.1029/2001WR001206

Lindgren RJ, Dutton AR, Hovorka SD, Worthington SR, Painter S (2005) Conceptualization and simulation of the Edwards aquifer, San Antonio region, Texas. In: Sinkholes and the engineering and environmental impacts of karst, pp 122-130

Long AJ, Putnam LD (2009) Age-distribution estimation for karst groundwater: Issues of parameterization and complexity in inverse modelling by convolution. J Hydrol 376:579-588. https://doi.org/10.1016/j.jhydrol.2009.07.064

Luhmann AJ, Covington MD, Peters AJ, Alexander SC, Anger CT, Green JA, Runkel AC, Alexander EC Jr (2011) Classification of thermal patterns at karst springs and cave streams. Groundwater 49:324-335. https://doi.org/10.1111/j.1745-6584.2010.00737.x

Majdalani S, Guinot V, Delenne C, Gebran H (2018) Modelling solute dispersion in periodic heterogeneous porous media: model benchmarking against intermediate scale experiments. J Hydrol 561:427-443. https://doi.org/10.1016/j.jhydrol.2018.03.024

Malard A, Sinreich M, Jeannin P-Y (2016) A novel approach for estimating karst groundwater recharge in mountainous regions and its application in Switzerland. Hydrol Process 30:2153-2166

Malík P, Vojtková S (2012) Use of recession-curve analysis for estimation of karstification degree and its application in assessing overflow/underflow conditions in closely spaced karstic springs. Environ Earth Sci 65:2245-2257

Maloszewski P, Harum T, Benischke R (1992) Mathematical modelling of tracer experiments in the karst of Lurbach system. Steierische Beiträge Zur Hydrogeologie 43:116-136

Maloszewski P, Benischke R, Harum T, Zojer H (1994) Estimation of solute transport parameters in heterogeneous groundwater system of a karstic aquifer using artificial tracer experiments. Water down under 94: groundwater papers; preprints of papers 105

Mangin A (1975) Contribution à l'étude hydrodynamique des aquifères karstiques. Ph.D. thesis, Université de Dijon

Maréchal J-C, Ladouche B, Dörfliger N, Lachassagne P (2008) Interpretation of pumping tests in a mixed flow karst system. Water Resour Res 44:W05401. https://doi.org/10.1029/2007WR006288

Martin JB, Dean RW (1999) Temperature as a natural tracer of short residence times for groundwater in karst aquifers. Karst Model 5:236-242

Martínez-Santos P, Andreu J (2010) Lumped and distributed approaches to model natural recharge in semiarid karst aquifers. J Hydrol 388:389-398

Massei N, Wang HQ, Field MS, Dupont JP, Bakalowicz M, Rodet J (2006) Interpreting tracer breakthrough tailing in a conduit-dominated karstic aquifer. Hydrogeol J 14:849-858

Massei N, Mahler BJ, Bakalowicz M, Fournier M, Dupont J-P (2007) Quantitative interpretation of specific conductance frequency distributions in karst. Groundwater 45:288-293

Mathevet T, Lepiller M, I, Mangin, A., (2004) Application of time-series analyses to the hydrological functioning of an Alpine karstic system: the case of Bange-L'Eau-Morte. Hydrol Earth Syst Sci 8:1051-1064. https://doi.org/10.5194/hess-8-1051-2004

Mazzilli N, Jourde H, Jacob T, Guinot V, Le Moigne N, Boucher M, Chalikakis K, Guyard H, Legtchenko A (2013) On the inclusion of ground-based gravity measurements to the calibration process of a global rainfall-discharge reservoir model: case of the Durzon karst system (Larzac, southern France). Environ Earth Sci 68:1631-1646

Mazzilli N, Guinot V, Jourde H, Lecoq N, Labat D, Arfib B, Baudement C, Danquigny C, Dal Soglio L, Bertin D (2019) Karst-Mod: a modelling platform for rainfall-discharge analysis and modelling dedicated to karst systems. Environ Model Softw 122:103927. https://doi.org/10.1016/j.envsoft.2017.03.015

Mazzilli N, Jourde H, Guinot V, Bailly-Comte V, Fleury P (2011) Hydrological modelling of a karst aquifer under active groundwater management using a parsimonious conceptual model

Mohammadi Z, Illman WA, Karimi M (2018) Optimization of the hydrodynamic characteristics of a karst conduit with CFPv2 coupled to OSTRICH. J Hydrol 567:564-578

Morales T, Uriarte JA, Olazar M, Antigüedad I, Angulo B (2010) Solute transport modelling in karst conduits with slow zones during different hydrologic conditions. J Hydrol 390:182-189

Moreno L, Tsang CF (1991) Multiple-peak response to tracer injection tests in single fractures: a numerical study. Water Resour Res 27:2143-2150

Moussa R, Bocquillon C (1996) Criteria for the choice of flood-routing methods in natural channels. J Hydrol 186:1-30

Mudarra M, Andrea B, Marín AI, Vadillo I, Barberá JA (2014) Combined use of natural and artificial tracers to determine the
hydrogeological functioning of a karst aquifer: the Villanueva del Rosario system (Andalusia, southern Spain). Hydrogeol J 22:1027–1039

Mull DS, Liebermann TD, Smoot JL, Woosley LH (1988) Application of dye-tracing techniques for determining solute-transport characteristics of ground water in karst terranes. Environmental Protection Agency, Atlanta, GA (United States). Region IV

Noushabadi MJ, Jourde H, Massonnat G (2011) Influence of the observation scale on permeability estimation at local and regional scales through well tests in a fractured and karstic aquifer (Lez aquifer, Southern France). J Hydrol 403:321–336

Oehlmann S, Geyer T, Licha T, Sauter M (2015) Reducing the ambiguity of karst aquifer models by pattern matching of flow and transport on catchment scale. Hydrol Earth Syst Sci 19:893–912

Ofterdinger U, Comte J-C, Cassidy R, Caulfield J, Pilatova K, Nitsche J, Wilson C, Cai Z, Flynn R (2015) Multidisciplinary investigations of poorly productive hard rock aquifers in Ireland: typologies, properties and significance in the Irish water cycle. Extended Abstract, IAH International Conference Hard Rock Aquifers: Up to Date Concepts and Practical Applications, La Roche-Sur-Yon, France, 11/06/2015 – 13/06/2015

Padilla A, Pulido-Bosch A (1995) Study of hydrographs of karstic aquifers by means of correlation and cross-spectral analysis. J Hydrol 168:73–89. https://doi.org/10.1016/0022-1694(94)02648-U

Palmer AN (1991) Origin and morphology of limestone caves. Geol Soc Am Bull 103:1–21

Panagopoulos G, Lambrakis N (2006) The contribution of time series analysis to the study of the hydrodynamic characteristics of the karst systems: application on two typical karst aquifers of Greece (Trifilia, Almyros Crete). J Hydrol 329:368–376. https://doi.org/10.1016/j.jhydrol.2006.02.023

Pardo-Igúzquiza E, Durán-Valsero JJ, Rodríguez-Galiano V (2011) Morphometric analysis of three-dimensional networks of karst conduits. Geomorphology 132:17–28

Pardo-Igúzquiza E, Dowd PA, Xu C, Durán-Valsero JJ (2012) Stochastic simulation of karst conduit networks. Adv Water Resour 35:141–150

Pardo-Igúzquiza E, Collados-Lara AJ, Pulido-Velazquez D (2019) Potential future impact of climate change on recharge in the Sierra de las Nieves (southern Spain) high-relief karst aquifer using regional climate models and statistical corrections. Environ Earth Sci 78:598. https://doi.org/10.1007/s12665-019-8594-4

Perineau A, Danquigny C, Emblanch C, di Borgo EP, Boyer D, Poupeney J (2011) Hydrodynamic organisation of the flows in the unsaturated zone of the Fontaine de Vaucluse karst system. First results. i-DUST 2010 01001

Perrin J, Luetscher M (2008) Inference of the structure of karst conduits using quantitative tracer tests and geological information: example of the Swiss Jura. Hydrogeol J 16:951–967

Perrin J, Jeannin P-Y, Zwahlen F (2003) Epikarst storage in a karst aquifer: a conceptual model based on isotopic data, Milandre test site, Switzerland. J Hydrol 279:106–124

Perrin J, Jeannin P-Y, Cornaton F (2007) The role of tributary mixing in chemical variations at a karst spring, Milandre, Switzerland. J Hydrol 332:158–173

Peterson EW, Wicks CM (2006) Assessing the importance of conduit geometry and physical parameters in karst systems using the storm water management model (SWMM). J Hydrol 329:294–305

Pinault J-L, Plagnes V, Aquilina L, Bakalowicz M (2001) Inverse modelling of the hydrological and the hydrochemical behavior of hydrosystems: characterization of karst system functioning. Water Resour Res 37:2191–2204

Popov P, Qin G, Bi L, Efendiev Y, Ewing RE, Li J (2009) Multiphysics and multiscale methods for modelling fluid flow through

naturally fractured carbonate karst reservoirs. SPE Reserv Eval Eng 12:218–231. https://doi.org/10.2118/105378-PA

Quinn JJ, Tomasko D, Kuiper JA (2006) Modelling complex flow in a karst aquifer. Sed Geol 184:343–351

Quinn JJ, Tomasko D (2000) A numerical approach to simulating mixed flow in karst aquifers. In: Groundwater flow and contaminant transport in carbonate aquifers, pp 147–156

Radulovic M, Stevanovic Z, Radulovic M (2012) A new approach in assessing recharge of highly karstified terrains–Montenegro case studies. Environ Earth Sci 65:2221–2230

Raeisi E, Groves C, Meiman J (2007) Effects of partial and full pipe flow on hydrochemographs of Logsdon river, Mammoth Cave Kentucky USA. J Hydrol 337:1–10. https://doi.org/10.1016/j.jhydrol.2006.11.015

Raeisi E, Karami G (2011) Hydrochemographs of Berghan karst spring as indicators of aquifer characteristics. J Cave Karst Stud Natl Speleol Soc Bull

Rahnemaei M, Zare M, Nematollahi AR, Sedghi H (2005) Application of spectral analysis of daily water level and spring discharge hydrographs data for comparing physical characteristics of karstic aquifers. J Hydrol 311:106–116. https://doi.org/10.1016/j.jhydrol.2005.01.011

Reimann T, Hill ME (2009) Performance review of the MODFLOW-2005 conduit flow process. Ground Water 47:321–325

Reimann T, Giese M, Geyer T, Liedl R, Maréchal J-C, Shoemaker WB (2014) Representation of water abstraction from a karst conduit with numerical discrete-continuum models. Hydrol Earth Syst Sci 18:227–241

Reimann T, Rehrl C, Shoemaker WB, Geyer T, Birk S (2011) The significance of turbulent flow representation in single-continuum models. Water Resour Res 47:W11503. https://doi.org/10.1029/2011WR010446

Rimmer A, Hartmann A (2012) Simplified conceptual structures and analytical solutions for groundwater discharge using reservoir equations. Water Resour Manag Model 2:217–238

Rimmer A, Salingar Y (2006) Modelling precipitation-streamflow processes in karst basin: the case of the Jordan River sources, Israel. J Hydrol 331:524–542

Ronayne MJ (2013) Influence of conduit network geometry on solute transport in karst aquifers with a permeable matrix. Adv Water Resour 56:27–34

Rongier G, Collon-Drouaillet P, Filipponi M (2014) Simulation of 3D karst conduits with an object-distance based method integrating geological knowledge. Geomorphology 217:152–164

Saar MO (2011) Review: geothermal heat as a tracer of large-scale groundwater flow and as a means to determine permeability fields. Hydrogeol J 19:31–52. https://doi.org/10.1007/s10040-010-0657-2

Saller SP, Ronayne MJ, Long AJ (2013) Comparison of a karst groundwater model with and without discrete conduit flow. Hydrogeol J 21:1555–1566

Sauter M, Geyer T, Kovács A, Teutsch G (2006) Modellierung der hydraulik von karstgrundwasserleitern–eine übersicht. Grundwasser 11:143–156

Sauter M (1992) A numerical approach to simulating mixed flow in karst aquifers. In: Sasowsky I, Wicks C (eds) Groundwater Flow and Contaminant Transport in Carbonate Aquifers. A.A. Balkema, Rotterdam, Holland, pp. 147–156

Scanlon BR, Healy RW, Cook PG (2002) Choosing appropriate techniques for quantifying groundwater recharge. Hydrogeol J 10:18–39

Scanlon BR, Mace RE, Barrett ME, Smith B (2003) Can we simulate regional groundwater flow in a karst system using equivalent porous media models? Case study, Barton Springs Edwards aquifer, USA. J Hydrol 276:137–158
Schmidt S, Geyer T, Guttman J, Marei A, Ries F, Sauter M (2014) Characterisation and modelling of conduit restricted karst aquifers—example of the Auja spring, Jordan Valley. J Hydrol 511:750–763

Shoemaker WB, Cunningham KJ, Kuniansky EL, Dixon J (2008) Effects of turbulence on hydraulic heads and parameter sensitivities in preferential groundwater flow layers. Water Resour Res 44:W03501. https://doi.org/10.1029/2007WR006601

Siemers J, Dreybrodt W (1998) Early development of karst aquifers on percolation networks of fractures in limestone. Water Resour Res 34:409–419. https://doi.org/10.1029/97WR03218

Sivelle V, Labat D, Mazzilli N, Massei N, Jourde H (2019) Dynamics of the flow exchanges between matrix and conduits in karstified watersheds at multiple temporal scales. Water 11:569. https://doi.org/10.3390/w11030569

Sivelle V, Renard P, Labat D (2020) Coupling SKS and SWMM to solve the inverse problem based on artificial tracer tests in Karstic Aquifers. Water 12:1139. https://doi.org/10.3390/w12041139

Sivelle V, Pérotin L, Ladouche B, de Montety V, Bailly-Comte V, Champollion C, Jourde H (2022) A lumped parameter model to evaluate the relevance of excess air as a tracer of exchanged flows between transmissive and capacitive compartments of karst systems. Front Water 4:930115. https://doi.org/10.3389/frwa.2022.930115

Skopp J, Gardner WR, Tyler EJ (1981) Solute movement in structured soils: two-region model with small interaction. Soil Sci Soc Am J 45:837–842

Stark CP (1991) An invasion percolation model of drainage network evolution. Nature 352:423–425

Sukop MC, Anwar S, Lee JS, Cunningham KJ, Langevin CD (2008) Modelling ground-water flow and solute transport in karst with lattice Boltzmann methods. In: Proceedings of the US Geological Survey Karst Interest Group workshop, Bowling Green, pp 27–29

Szymczak P, Ladd AJ (2012) Reactive-infiltration instabilities in rocks. Fracture dissolution. arXiv preprint arXiv:1204.0917

Szymczak P, Ladd AJ (2011) The initial stages of cave formation: beyond the one-dimensional paradigm. Earth Planet Sci Lett 301:424–432

Teutsch G, Sauter M (1991) Groundwater modelling in karst terranes: scale effects, data acquisition and field validation. In: Third conference on hydrogeology, ecology, monitoring, and management of ground water in Karst Terranes. National Ground Water Association, Dublin, Ohio, pp 17–35

Therrien R, McLaren RG, Sudicky EA (2010) HydroGeoSphre: a three dimensional numerical model describing fully integrated subsurface and surface flow and solute transport. Groundwater Simulations Group, University of Waterloo, Waterloo, ON

Thrailkill J (1974) Pipe flow models of a Kentucky limestone aquifer. Groundwater 12:202–205

Toran L, Tancredi JH, Herman EK, White WB (2006) Conductivity and sediment variation during storms as evidence of pathways to karst springs. Spec Pap the Geol Soc Am 404:169–176

Toride N, Leij FJ (1996) Convective-dispersive stream tube model for field-scale solute transport: I. moment analysis. Soil Sci Soc Am J 60:342–351. https://doi.org/10.2136/sssaj1996.03615995006000020004x

Tritz S, Guinot V, Jourde H (2011) Modelling the behaviour of a karst system catchment using non-linear hysteretic conceptual model. J Hydrol 397:250–262

Van Genuchten MT, Wierenga PJ (1977) Mass transfer studies in sorbing porous media: II. Experimental evaluation with tritium (3H₂O). Soil Sci Soc Am J 41:272–278

Wang X, Jardani A, Jourde H, Lonergan L, Cosgrove J, Gosselin O, Massonnat G (2016) Characterisation of the transmissivity field of a fractured and karstic aquifer, Southern France. Adv Water Resour 87:106–121. https://doi.org/10.1016/j.advwatres.2015.10.014

Wang X, Jardani A, Jourde H (2017) A hybrid inverse method for hydraulic tomography in fractured and karstic media. J Hydrol Investig Coast Aquif 551:29–46. https://doi.org/10.1016/j.jhydrol.2017.05.051

Wang C, Wang X, Majdalani S, Guinot V, Jourde H (2020) Influence of dual conduit structure on solute transport in karst tracer tests: an experimental laboratory study. J Hydrol 590:125255

Werner A, Hötzl H, Käss W (1997) The interpretation of a high water tracer test in the Danube-Aach-System (Western Swabian Alb, Germany). Sciences & Techniques de l'Environnement, université de Franche-Comté. Mémoire hors-série, pp 187–190

Williams PW (1983) The role of the subcutaneous zone in karst hydrology. J Hydrol 61:45–67

Worthington SRH, Ford DC (2009) Self-organized permeability in carbonate aquifers. Groundwater 47:326–336

Wu Q, Xu H, Pang W (2008) GIS and ANN coupling model: an innovative approach to evaluate vulnerability of karst water inrush in coalmines of north China. Environ Geol 54:937–943. https://doi.org/10.1007/s00254-007-0887-3

Xanke J, Jourde H, Liesch T, Goldscheider N (2016) Numerical long-term assessment of managed aquifer recharge from a reservoir into a karst aquifer in Jordan. J Hydrol 540:603–614

Zheng C, Wang PP (1999) MT3DMS: a modular three-dimensional multispecies transport model for simulation of advection, dispersion, and chemical reactions of contaminants in groundwater systems. Documentation and User's Guide

Zoghbi C, Basha H (2020) Simple transport models for karst systems. J Hydrol 588:125046

**Publisher's Note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.