in the 525 km² Geropotamos basin located in central-southern part of the island. Nerantzaki et al. [161] also used Karst-SWAT to forecast the hydrological response of the Crete region under climate change scenarios up to the year 2098, considering different irrigation sources in SWAT. Moreover, Tapoglou et al. [159] applied Karst-SWAT to predict the impact of climate change on the hydrological cycle and the frequency of extreme hydrological and meteorological events in Crete. Nerantzaki et al. [163] further expanded the research work conducted in the Koiliaris River basin with Karst-SWAT by assessing: (1) the uncertainty of the watershed runoff and karstic flow simulations due to the parameter uncertainty in SWAT and Karst-SWAT, and (2) the impact of internal variability (or stochastic uncertainty) of the meteorological input data on the flow simulations for the reference period and under the climate change scenarios. The uncertainty of the flow models was estimated by combining the Sequential Uncertainty Fitting Version 2 (SUFI-2) in the SWAT Calibration and Uncertainty Program (SWAT-CUP) interface and the @RISK by PALISADE software, while the effect of input internal variability on the flow output was evaluated using Monte Carlo simulations.

Within the framework of studying the hydrological and geochemical processes in the Koiliaris European Critical Zone Observatory, Lilli et al. [164] used Karst-SWAT to simulate the hydrological budget of the Koiliaris River basin and gain insight on the hydrological pathways and response of the karst during extreme events. Additionally, Karst-SWAT was applied to simulate surface flow and spring flow in the Koiliaris River basin, which were required to determine the design flows and flood frequency within the framework of developing nature-based solutions for the riparian forest restoration and flood protection project at the Koiliaris Critical Zone Observatory [162].

Next, Malagò et al. [165] developed a two-reservoir modeling approach by linking SWAT-B&B [153] and Karst-SWAT [157]. The resultant hybrid model, called KSWAT, was used in conjunction with SWAT to simulate the water balance, spring flow, and total discharge in the Island of Crete in Greece. The study area extends over 8336 km² of which 2730 km² is karst.

A SWAT model of the Crete Island was set up and the modified model KSWAT was applied only in the karst subbasins of the region. The daily aquifer recharge from the karst subbasins was simulated using SWAT-B&B. The area of the subbasins contributing to the recharge of a particular spring or group of springs was identified based on the local geological maps and dominant karst soils. Recharge from the soil profile bottom, stream losses, and seepage from other water bodies to the deep aquifer were maximized by: (1) setting the deep aquifer percolation fraction and minimum groundwater delay to 1, (2) adjusting the groundwater coefficient of capillary rise to 0.1 to prevent the upward movement of water to the unsaturated zone, and (3) minimizing the return flow from the shallow and deep aquifers in SWAT. The deep aquifer recharge time series generated by SWAT-B&B were then input to Karst-SWAT in order to simulate and calibrate the discharge of the springs. The parameters of the Karst-SWAT model were adjusted based on daily spring discharge data from 47 gauging stations in Crete.

The hydrologic model in SWAT was adjusted using a step-wise calibration with monthly streamflow data at 15 stream gauging stations. Snow, surface runoff, lateral flow, and baseflow parameters were first calibrated separately in order to adjust the timing of the runoff signal and the discharge values (peak flow and baseflow). Then, the model was recalibrated based on the streamflow-related parameters combined with the adjusted variables of the other water budget components. The final near optimal parameter set of the calibrated subbasins was transferred to the ungauged subbasins using the hydrological similarity approach with the Partial Least Squares Regression, in order to identify similar subbasins based on the correlation between the watershed and discharge characteristics.

Subsequently, the calibrated spring discharge time series from Karst-SWAT were added to the Crete SWAT model as point sources in order to predict the total monthly runoff across the island, and a final calibration was performed to adjust discharge. The results of the performance indicators showed that only 40% of the calibrated gauging
stations scored NSE values greater than 0.5, while 50% had R² values higher than 0.5 and 64% reached PBIAS lower than 25%.

Moreover, Nguyen et al. [166] added a two-reservoir karstic flow model to the original groundwater module of SWAT. The improved SWAT code, termed SWAT_IFG, consists of two conceptual groundwater models compiled in a single executable file: (1) the standard SWAT one-reservoir model applied to non-karst terrains, and (2) the modified two-reservoir model used in karst areas. The two-reservoir groundwater model of SWAT_IFG represents a variant of the Karst-flow model by Nikolaidis et al. [157]. In SWAT_IFG, the matrix reservoir receives diffuse recharge as a linear function of daily infiltration from the soil bottom, considering the time delay of flow in the unsaturated zone (Equations (A15) and (A16) of Appendix H). The conduit reservoir receives another fraction of the soil water seepage as concentrated recharge with infiltration losses from sinking streams. It is also fed by diffuse discharge from the matrix reservoir, which represents the flow exchange mechanism between the two karst domains (Equations (A18) and (A19) in Appendix H). Groundwater outflow from the matrix storage reservoir to the conduit is modeled using a linear storage-discharge equation with a matrix recession coefficient (Equation (A17) in Appendix H). Outflow from the conduit reservoir to the spring is also modeled via a linear storage-discharge relationship adjusted for the total recharge volume to the conduits and a conduit recession coefficient (Equation (A20) in Appendix H). The total discharge of the basin where the spring is located is then simulated as the sum of the direct runoff (surface runoff and lateral flow) and the outflow from the conduit reservoir to the spring (Equation (A21) in Appendix H). Consequently, SWAT_IFG can simulate surface runoff and subsurface flow in non-karst and karst areas. It can also account for the dual recharge and storage functions of the matrix and conduits. and model spring discharge in regions where the karst aquifer boundaries do not coincide with the surface subbasin areas.

SWAT_IFG was applied to simulate discharge in the drainage basin of the karst dominated southern Harz rim and non-karst southwest Harz Mountains in Northern Germany. The watershed covers an area of 384 km² and has one river outlet and a main spring outlet (the Rhume spring). The spring is mainly fed by allogenic recharge and river transmission losses from upstream subbasins via a connected network of losing streams in the area, with only 4% of the spring discharge originating from autogenic recharge and nearly 96% from IGF.

When applying SWAT_IFG, an aquifer classification map with information about the aquifer types and spring recharge area in the study site was incorporated into the model to delineate karst and non-karst HRUs. Subsequently, the suitable conceptual reservoir model could be assigned (the two-reservoir model for the karst HRUs and the one-reservoir model for the non-karst HRUs), and recharge from the extended karst area could be routed to the spring outlet.

SWAT_IFG was run for 14 years at the daily time step, with 3 years of warm-up, 6 years of calibration, and 5 years of validation. The model parameters were fitted based on multi-site daily streamflow data and satellite-derived actual evapotranspiration records (the Moderate Resolution Imaging Spectroradiometer MOD16 ETa). A multi-criteria NSE objective function was used to assess the overall performance of the model simulation, with equal weights allocated to the multi-gauge streamflow observations and evapotranspiration data. Results showed that the use of MOD16 ETa data in the calibration did not affect the model performance. The flow simulation at the spring outlet improved with multi-gauge calibration, as the NSE values varied from 0.75–0.48 (calibration-validation) with the single gauge calibration to 0.69–0.62 under the multi-site calibration. The model performance for all remaining streamflow gauging stations also improved with multi-site calibration, and NSE values of 0.54–0.91 and 0.6–0.91 were reached for the calibration and validation periods, respectively. Additionally, the model prediction uncertainty was reduced. The PBIAS values calculated at the different gages fell below 10%, while KGE values ranged between 0.68 and 0.91. Yet, the observed and simulated streamflow hydrographs showed that SWAT_IFG underestimated the high and low flows, which is a property inherited
stations scored NSE values greater than 0.5, while 50% had R² values higher than 0.5 and 64% reached PBIAS lower than 25%.

Moreover, Nguyen et al. [166] added a two-reservoir karstic flow model to the original groundwater module of SWAT. The improved SWAT code, termed SWAT_IFG, consists of two conceptual groundwater models compiled in a single executable file: (1) the standard SWAT one-reservoir model applied to non-karst terrains, and (2) the modified two-reservoir model used in karst areas. The two-reservoir groundwater model of SWAT_IFG represents a variant of the Karst-flow model by Nikolaidis et al. [157]. In SWAT_IFG, the matrix reservoir receives diffuse recharge as a linear function of daily infiltration from the soil bottom, considering the time delay of flow in the unsaturated zone (Equations (A15) and (A16) of Appendix H). The conduit reservoir receives another fraction of the soil water seepage as concentrated recharge with infiltration losses from sinking streams. It is also fed by diffuse discharge from the matrix reservoir, which represents the flow exchange mechanism between the two karst domains (Equations (A18) and (A19) in Appendix H). Groundwater outflow from the matrix storage reservoir to the conduit is modeled using a linear storage-discharge equation with a matrix recession coefficient (Equation (A17) in Appendix H). Outflow from the conduit reservoir to the spring is also modeled via a linear storage-discharge relationship adjusted for the total recharge volume to the conduits and a conduit recession coefficient (Equation (A20) in Appendix H). The total discharge of the basin where the spring is located is then simulated as the sum of the direct runoff (surface runoff and lateral flow) and the outflow from the conduit reservoir to the spring (Equation (A21) in Appendix H). Consequently, SWAT_IFG can simulate surface runoff and subsurface flow in non-karst and karst areas. It can also account for the dual recharge and storage functions of the matrix and conduits. and model spring discharge in regions where the karst aquifer boundaries do not coincide with the surface subbasin areas.

SWAT_IFG was applied to simulate discharge in the drainage basin of the karst dominated southern Harz rim and non-karst southwest Harz Mountains in Northern Germany. The watershed covers an area of 384 km² and has one river outlet and a main spring outlet (the Rhume spring). The spring is mainly fed by allogenic recharge and river transmission losses from upstream subbasins via a connected network of losing streams in the area, with only 4% of the spring discharge originating from autogenic recharge and nearly 96% from IGF.

When applying SWAT_IFG, an aquifer classification map with information about the aquifer types and spring recharge area in the study site was incorporated into the model to delineate karst and non-karst HRUs. Subsequently, the suitable conceptual reservoir model could be assigned (the two-reservoir model for the karst HRUs and the one-reservoir model for the non-karst HRUs), and recharge from the extended karst area could be routed to the spring outlet.

SWAT_IFG was run for 14 years at the daily time step, with 3 years of warm-up, 6 years of calibration, and 5 years of validation. The model parameters were fitted based on multi-site daily streamflow data and satellite-derived actual evapotranspiration records (the Moderate Resolution Imaging Spectroradiometer MOD16 ETa). A multi-criteria NSE objective function was used to assess the overall performance of the model simulation, with equal weights allocated to the multi-gauge streamflow observations and evapotranspiration data. Results showed that the use of MOD16 ETa data in the calibration did not affect the model performance. The flow simulation at the spring outlet improved with multi-gauge calibration, as the NSE values varied from 0.75–0.48 (calibration-validation) with the single gauge calibration to 0.69–0.62 under the multi-site calibration. The model performance for all remaining streamflow gauging stations also improved with multi-site calibration, and NSE values of 0.54–0.91 and 0.6–0.91 were reached for the calibration and validation periods, respectively. Additionally, the model prediction uncertainty was reduced. The PBIAS values calculated at the different gages fell below 10%, while KGE values ranged between 0.68 and 0.91. Yet, the observed and simulated streamflow hydrographs showed that SWAT_IFG underestimated the high and low flows, which is a property inherited
from the original SWAT model. Nonetheless, the model successfully simulated IGF and transmission losses from the rivers contributing to the spring discharge.

### 4.2.3. Conceptual Linear Three-Reservoir Model

This subsection describes the modified SWAT models with a linear three-reservoir groundwater module to simulate flows of the main karst aquifer components (i.e., epikarst, matrix and conduits). These models include the modified SWAT by [167,168].

Wang et al. [167] coupled SWAT (v2012) with a linear three-reservoir model. The modified model consists of: (1) an upper reservoir that reproduces the regulation and storage function of the epikarst and is recharged by percolation from the soil bottom, (2) a middle reservoir that represents the conduits system fed by infiltration from the epikarst, depressions, and avens, and (3) a lower reservoir corresponding to the matrix system recharged by the epikarst and conduit reservoirs. Daily infiltration to the upper reservoir is simulated as a function of the saturation moisture content in the epikarst system, its water-holding content, and saturated hydraulic conductivity (Equations (A22) and (A23) of Appendix I). To simulate the intercompartment fluxes, a proportionality coefficient ($a_1$: 0.5–1) is introduced to separate the recharge from the epikarst reservoir between the quick flow and slow flow reservoirs, based on the degree of karstification of the watershed. Another coefficient ($a_2$: 0.1–0.5) is used to split the discharge from the conduit reservoir between the slow flow reservoir and the basin outlet (Equations (A24) and (A25) of Appendix I). The outflows from the conduit and matrix reservoirs are modeled using the standard attenuation functions of SWAT (See Equations (A26) and (A27) of Appendix I). Finally, the total karst flow is calculated as the sum of the discharge components of the matrix and conduit reservoirs (Equation (A28) of Appendix I), which is then added to the surface runoff to simulate the total discharge at the watershed outlet.

The original and modified SWAT models were applied to predict daily runoff for a fully karstified watershed of 26.76 km² located in Hunan Province, China, with a calibration period of 180 days and a validation period of 100 days. The study area is primarily covered by Devonian and Carboniferous limestone and exhibits karst depressions, caves, and underground rivers. The SWAT three-reservoir model yielded a streamflow simulation that was significantly better than that obtained by the standard SWAT model in both calibration and validation, with NSE values increasing from 0.57–0.63 to 0.81–0.83 and R² values increasing from 0.58–0.62 to 0.82–0.84.

Geng et al. [168] later modified the SWAT model proposed by Wang et al. [167] in order to improve the simulation of rapid recharge to the epikarst reservoir through direct water percolation from the soil bottom without attenuation. The modelers also added a discharge component from the epikarst reservoir to the river channels. Three coefficients of proportionality were thus introduced to the three-reservoir groundwater model in order to separate the flow from the epikarst reservoir between a surface runoff component and two recharge components to the matrix and conduit reservoirs (See Equations (A29) to (A31) of Appendix J). Rapid infiltration through sinkholes, ponors, and fractures was also replaced by pond leakage with concentrated (fast) recharge similar to the computational method proposed by Baffaut and Benson [153]. The remaining intercompartment fluxes were modeled similarly to the model by Wang et al. [167].

The modified SWAT model was applied to simulate the hydrological cycle processes in the Daotian River basin, including the contribution of the streamflow and baseflow components to the runoff at the watershed outlet. The study site is situated in the Guizhou Province, China, and has a temperate monsoon climate. It covers a total area of 99.21 km² of which ~53% is dolomite, ~38% limestones, and ~9% clastic rocks. The karst landscape is characterized by karst depressions, sinkholes, and well-connected networks of conduits of high hydraulic conductivity, particularly in the limestone area. Due to karst effects, the watershed recharge boundaries extend by 24.75 km² beyond its surface drainage area, and the additional water is discharged into the watershed through underground conduits. The areas outside the topographic drainage divides and the flow paths of the karst subterranean
from the original SWAT model. Nonetheless, the model successfully simulated IGF and transmission losses from the rivers contributing to the spring discharge.

### 4.2.3. Conceptual Linear Three-Reservoir Model

This subsection describes the modified SWAT models with a linear three-reservoir groundwater module to simulate flows of the main karst aquifer components (i.e., epikarst, matrix and conduits). These models include the modified SWAT by [167,168].

Wang et al. [167] coupled SWAT (v2012) with a linear three-reservoir model. The modified model consists of: (1) an upper reservoir that reproduces the regulation and storage function of the epikarst and is recharged by percolation from the soil bottom, (2) a middle reservoir that represents the conduits system fed by infiltration from the epikarst, depressions, and avens, and (3) a lower reservoir corresponding to the matrix system recharged by the epikarst and conduit reservoirs. Daily infiltration to the upper reservoir is simulated as a function of the saturation moisture content in the epikarst system, its water-holding content, and saturated hydraulic conductivity (Equations (A22) and (A23) of Appendix I). To simulate the intercompartment fluxes, a proportionality coefficient ($a_1$: 0.5–1) is introduced to separate the recharge from the epikarst reservoir between the quick flow and slow flow reservoirs, based on the degree of karstification of the watershed. Another coefficient ($a_2$: 0.1–0.5) is used to split the discharge from the conduit reservoir between the slow flow reservoir and the basin outlet (Equations (A24) and (A25) of Appendix I). The outflows from the conduit and matrix reservoirs are modeled using the standard attenuation functions of SWAT (See Equations (A26) and (A27) of Appendix I). Finally, the total karst flow is calculated as the sum of the discharge components of the matrix and conduit reservoirs (Equation (A28) of Appendix I), which is then added to the surface runoff to simulate the total discharge at the watershed outlet.

The original and modified SWAT models were applied to predict daily runoff for a fully karstified watershed of 26.76 km² located in Hunan Province, China, with a calibration period of 180 days and a validation period of 100 days. The study area is primarily covered by Devonian and Carboniferous limestone and exhibits karst depressions, caves, and underground rivers. The SWAT three-reservoir model yielded a streamflow simulation that was significantly better than that obtained by the standard SWAT model in both calibration and validation, with NSE values increasing from 0.57–0.63 to 0.81–0.83 and R² values increasing from 0.58–0.62 to 0.82–0.84.

Geng et al. [168] later modified the SWAT model proposed by Wang et al. [167] in order to improve the simulation of rapid recharge to the epikarst reservoir through direct water percolation from the soil bottom without attenuation. The modelers also added a discharge component from the epikarst reservoir to the river channels. Three coefficients of proportionality were thus introduced to the three-reservoir groundwater model in order to separate the flow from the epikarst reservoir between a surface runoff component and two recharge components to the matrix and conduit reservoirs (See Equations (A29) to (A31) of Appendix J). Rapid infiltration through sinkholes, ponors, and fractures was also replaced by pond leakage with concentrated (fast) recharge similar to the computational method proposed by Baffaut and Benson [153]. The remaining intercompartment fluxes were modeled similarly to the model by Wang et al. [167].

The modified SWAT model was applied to simulate the hydrological cycle processes in the Daotian River basin, including the contribution of the streamflow and baseflow components to the runoff at the watershed outlet. The study site is situated in the Guizhou Province, China, and has a temperate monsoon climate. It covers a total area of 99.21 km² of which ~53% is dolomite, ~38% limestones, and ~9% clastic rocks. The karst landscape is characterized by karst depressions, sinkholes, and well-connected networks of conduits of high hydraulic conductivity, particularly in the limestone area. Due to karst effects, the watershed recharge boundaries extend by 24.75 km² beyond its surface drainage area, and the additional water is discharged into the watershed through underground conduits. The areas outside the topographic drainage divides and the flow paths of the karst subterranean
rivers that exist in the watershed were determined by conducting a karst survey and an artificial tracer test prior to hydrological modeling. After determining the flow paths of the subterranean river based on the spatial distribution of the subterranean river inlet, ponors, and sinkholes, the DEM data were modified to convert the subterranean river into a surface river. By adopting this approach, the actual catchment boundaries of the watershed were correctly identified by SWAT for the subbasins delineation.

The modified SWAT model was calibrated using daily streamflow measurements at the watershed outlet for 6 years (1993–1998). Two validation periods were considered under various annual precipitation patterns: the first from 1999 to 2002 (normal, dry, and wet years) and the second from 2003 to 2006 (normal and dry years). The performance of the modified SWAT model was compared to that of a previous model run in SWAT at monthly time scale. Both models had satisfactory simulations of monthly discharge. Yet, the modified SWAT model had a better prediction efficiency than the original SWAT, scoring NSE values of 0.87–0.83/0.85, R² of 0.88–0.84/0.86, and PBIAS of 2.5%–(–1.9%/–15%) for the calibration and two validation periods, respectively. The three-reservoir model also improved the simulation of the karst water cycle by increasing the groundwater recharge and return flow components. As a result, the NSE for the baseflow simulation of the modified SWAT model was 0.09 higher than that of the original SWAT model, which underestimated flows below 0.7 m³.s⁻¹ in the dry periods and overestimated runoff during wet periods.

### 4.2.4. Conceptual Non-Linear One-Reservoir Model

To the best of our knowledge, a non-linear groundwater reservoir has been previously implemented few times in SWAT to: (1) model groundwater flow in a karstic aquifer [169], (2) estimate baseflow in rivers dominated by snow and glacier melt [179], and (3) improve the prediction of streamflow in watersheds with complex groundwater processes under heavy irrigation pumping [180].

Wang and Brubaker [169] replaced the linear reservoir model in SWAT with a single non-linear reservoir based on the algorithm of Wittenberg [181] (see Equations (A32) and (A33) of Appendix K), providing a modified SWAT version called ISWAT. The ISWAT model was tested in the Shenandoah River watershed of the Potomac River basin (USA), which drains a large karstic area of 7607 km². It was calibrated using 13 years of daily discharge records across 14 gauging stations, with 2 years of warm-up, and validated for 4 years. To account for the spatial variability of the geology and soils in the watershed during calibration, the parameters of the non-linear reservoir (i.e., groundwater recession coefficients and exponents) were grouped by soil type.

The ISWAT model performance was assessed against that of the linear SWAT model by comparing the simulated and observed streamflow hydrographs at the different gauging stations, and the recession curves in the low flow periods with the baseflow taken as the lowest 30% of daily discharge. The non-linear ISWAT model reproduced low flow discharge and recession curves better than SWAT, but simulated peak flows with a comparable accuracy to SWAT. The NSE (modified) and R² indices improved with the use of the non-linear model at the level of eight and ten observational river reaches, respectively, with values of 0.5 and 0.6. The non-linear model also lowered the overall relative bias of the simulations by 3%, with the majority of the observational river reaches scoring a bias less or equal to 10%.

### 4.2.5. Modified Crack Flow with Conceptual Linear One-Reservoir Model

This subsection summarizes the application of the modified SWAT models by [36], which were adapted to simulate fast recharge in karstified watersheds based on the crack flow approach in SWAT.

Eini et al. [36] modified SWAT (v2012) to increase groundwater recharge in karst areas by: (1) adjusting the groundwater recharge function in SWAT to increase infiltration in karst HRUs (the SWAT-ML model), and (2) expanding the crack flow module in SWAT to
rivers that exist in the watershed were determined by conducting a karst survey and an artificial tracer test prior to hydrological modeling. After determining the flow paths of the subterranean river based on the spatial distribution of the subterranean river inlet, ponors, and sinkholes, the DEM data were modified to convert the subterranean river into a surface river. By adopting this approach, the actual catchment boundaries of the watershed were correctly identified by SWAT for the subbasins delineation.

The modified SWAT model was calibrated using daily streamflow measurements at the watershed outlet for 6 years (1993–1998). Two validation periods were considered under various annual precipitation patterns: the first from 1999 to 2002 (normal, dry, and wet years) and the second from 2003 to 2006 (normal and dry years). The performance of the modified SWAT model was compared to that of a previous model run in SWAT at monthly time scale. Both models had satisfactory simulations of monthly discharge. Yet, the modified SWAT model had a better prediction efficiency than the original SWAT, scoring NSE values of 0.87–0.83/0.85, R² of 0.88–0.84/0.86, and PBIAS of 2.5%–(–1.9%/–15%) for the calibration and two validation periods, respectively. The three-reservoir model also improved the simulation of the karst water cycle by increasing the groundwater recharge and return flow components. As a result, the NSE for the baseflow simulation of the modified SWAT model was 0.09 higher than that of the original SWAT model, which underestimated flows below 0.7 m³.s⁻¹ in the dry periods and overestimated runoff during wet periods.

### 4.2.4. Conceptual Non-Linear One-Reservoir Model

To the best of our knowledge, a non-linear groundwater reservoir has been previously implemented few times in SWAT to: (1) model groundwater flow in a karstic aquifer [169], (2) estimate baseflow in rivers dominated by snow and glacier melt [179], and (3) improve the prediction of streamflow in watersheds with complex groundwater processes under heavy irrigation pumping [180].

Wang and Brubaker [169] replaced the linear reservoir model in SWAT with a single non-linear reservoir based on the algorithm of Wittenberg [181] (see Equations (A32) and (A33) of Appendix K), providing a modified SWAT version called ISWAT. The ISWAT model was tested in the Shenandoah River watershed of the Potomac River basin (USA), which drains a large karstic area of 7607 km². It was calibrated using 13 years of daily discharge records across 14 gauging stations, with 2 years of warm-up, and validated for 4 years. To account for the spatial variability of the geology and soils in the watershed during calibration, the parameters of the non-linear reservoir (i.e., groundwater recession coefficients and exponents) were grouped by soil type.

The ISWAT model performance was assessed against that of the linear SWAT model by comparing the simulated and observed streamflow hydrographs at the different gauging stations, and the recession curves in the low flow periods with the baseflow taken as the lowest 30% of daily discharge. The non-linear ISWAT model reproduced low flow discharge and recession curves better than SWAT, but simulated peak flows with a comparable accuracy to SWAT. The NSE (modified) and R² indices improved with the use of the non-linear model at the level of eight and ten observational river reaches, respectively, with values of 0.5 and 0.6. The non-linear model also lowered the overall relative bias of the simulations by 3%, with the majority of the observational river reaches scoring a bias less or equal to 10%.

### 4.2.5. Modified Crack Flow with Conceptual Linear One-Reservoir Model

This subsection summarizes the application of the modified SWAT models by [36], which were adapted to simulate fast recharge in karstified watersheds based on the crack flow approach in SWAT.

Eini et al. [36] modified SWAT (v2012) to increase groundwater recharge in karst areas by: (1) adjusting the groundwater recharge function in SWAT to increase infiltration in karst HRUs (the SWAT-ML model), and (2) expanding the crack flow module in SWAT to
retain the formation of cracks independently of soil moisture conditions (the SWAT-CF model). These modifications were based on the premise that preferential flow through the soil sinks and cracks can be representative of rapid recharge in karst landforms.

Both SWAT-ML and SWAT-CF were applied in the Maharlu Lake, a large watershed of 4270 km² in Southwest Iran, of which 37% is covered by extensive karst areas (limestone and dolomite). Several karst-fissured aquifers are well developed in these areas due to lithology, climate, and tectonic activity.

In SWAT-ML, infiltration from non-karst areas was calculated using the standard SWAT recharge equation (Equation (A34) of Appendix L), while fast infiltration from karst areas was modeled by dividing the delay variable in the original groundwater recharge equation by a new non-dimensional calibration parameter (X) (Equation (A35) of Appendix L). This parameter can be increased depending on the volume and numbers of cracks in a karst HRU.

In SWAT-CF, the standard SWAT function that calculates crack volume during the crack flow process was modified by adding a new parameter of the crack volume based on new moisture conditions (Equation (A36) of Appendix L) to achieve a year-round crack formation in the soil matrix. As a result, SWAT-CF can simulate crack flow in karst HRUs both in dry and wet soil conditions during surface flow events, as opposed to standard SWAT, which only models crack volume on drier days as a function of crop capacity and soil moisture.

The hydrological models developed with SWAT, SWAT-ML, and SWAT-CF were run using streamflow data recorded at 3 gauging stations, with 3 years of warm-up, 26 years of calibration, and 4 years of validation. Both modified SWAT models outperformed SWAT in simulating monthly streamflow at different stations, with SWAT-ML having the best overall accuracy. The average NSE value increased from 0.64 with SWAT to 0.67 using SWAT-CF and 0.69 using SWAT-ML, while the average R² value varied from 0.70 to 0.69 and 0.72 under SWAT-CF and SWAT-ML, respectively. The modified models also increased the prediction accuracy of the baseflow and water budget components.

### 4.2.6. Variable Source Area Hydrology with Conceptual Linear One-Reservoir Model

Topo-SWAT (initially termed SWAT-VSA) is a modified version of SWAT that was applied to simulate flow, sediments yield, and nutrients loads in a watershed with variable source area (VSA) hydrology [170]. Compared to SWAT, Topo-SWAT incorporates a topographic wetness index (TI) that indicates the saturation potential of a landscape unit and the subsequent likelihood of runoff generation. Ten equal-area wetness classes ranging from 1 to 10 (1 being 10% of the watershed with the lowest runoff potential, and 10 being 10% of the watershed with the highest runoff potential) were used to overlay a wetness class layer with the soil map of the study site and generate a single GIS layer and associated lookup tables for the SWAT slope class and soil layers.

Topo-SWAT was tested in the Spring Creek watershed in northeastern USA. Spring Creek has a surface water watershed area of 370 km² but is defined by a groundwater recharge boundary of 450 km², which is characterized by saturation excess surface runoff from VSAs (e.g., perched and losing streams in the headwater regions of the watershed, low surface runoff in the forested uplands due to quick infiltration through shallow soils, overland flow generated at the base of hillslopes). Some of the adjustments made to the model parameters to accurately represent karst hydrology in the study watershed included reducing the initial curve number, restricting the groundwater delay factor to 1 day, setting the baseflow recession factor to 0.011 day based on observed daily streamflow records, and introducing the contribution of the springs that recharge outside Spring Creek but discharge inside the watershed as point sources.

Topo-SWAT outperformed SWAT in modeling daily streamflow for a 12-year simulation period, with NSE of 0.73–0.79, PBIAS of –2.8 to –3.7%, and R² of 0.71–0.77 in the calibration and validation periods. Moreover, Topo-SWAT successfully reproduced the VSA
retain the formation of cracks independently of soil moisture conditions (the SWAT-CF model). These modifications were based on the premise that preferential flow through the soil sinks and cracks can be representative of rapid recharge in karst landforms.

Both SWAT-ML and SWAT-CF were applied in the Maharlu Lake, a large watershed of 4270 km² in Southwest Iran, of which 37% is covered by extensive karst areas (limestone and dolomite). Several karst-fissured aquifers are well developed in these areas due to lithology, climate, and tectonic activity.

In SWAT-ML, infiltration from non-karst areas was calculated using the standard SWAT recharge equation (Equation (A34) of Appendix L), while fast infiltration from karst areas was modeled by dividing the delay variable in the original groundwater recharge equation by a new non-dimensional calibration parameter (X) (Equation (A35) of Appendix L). This parameter can be increased depending on the volume and numbers of cracks in a karst HRU.

In SWAT-CF, the standard SWAT function that calculates crack volume during the crack flow process was modified by adding a new parameter of the crack volume based on new moisture conditions (Equation (A36) of Appendix L) to achieve a year-round crack formation in the soil matrix. As a result, SWAT-CF can simulate crack flow in karst HRUs both in dry and wet soil conditions during surface flow events, as opposed to standard SWAT, which only models crack volume on drier days as a function of crop capacity and soil moisture.

The hydrological models developed with SWAT, SWAT-ML, and SWAT-CF were run using streamflow data recorded at 3 gauging stations, with 3 years of warm-up, 26 years of calibration, and 4 years of validation. Both modified SWAT models outperformed SWAT in simulating monthly streamflow at different stations, with SWAT-ML having the best overall accuracy. The average NSE value increased from 0.64 with SWAT to 0.67 using SWAT-CF and 0.69 using SWAT-ML, while the average R² value varied from 0.70 to 0.69 and 0.72 under SWAT-CF and SWAT-ML, respectively. The modified models also increased the prediction accuracy of the baseflow and water budget components.

### 4.2.6. Variable Source Area Hydrology with Conceptual Linear One-Reservoir Model

Topo-SWAT (initially termed SWAT-VSA) is a modified version of SWAT that was applied to simulate flow, sediments yield, and nutrients loads in a watershed with variable source area (VSA) hydrology [170]. Compared to SWAT, Topo-SWAT incorporates a topographic wetness index (TI) that indicates the saturation potential of a landscape unit and the subsequent likelihood of runoff generation. Ten equal-area wetness classes ranging from 1 to 10 (1 being 10% of the watershed with the lowest runoff potential, and 10 being 10% of the watershed with the highest runoff potential) were used to overlay a wetness class layer with the soil map of the study site and generate a single GIS layer and associated lookup tables for the SWAT slope class and soil layers.

Topo-SWAT was tested in the Spring Creek watershed in northeastern USA. Spring Creek has a surface water watershed area of 370 km² but is defined by a groundwater recharge boundary of 450 km², which is characterized by saturation excess surface runoff from VSAs (e.g., perched and losing streams in the headwater regions of the watershed, low surface runoff in the forested uplands due to quick infiltration through shallow soils, overland flow generated at the base of hillslopes). Some of the adjustments made to the model parameters to accurately represent karst hydrology in the study watershed included reducing the initial curve number, restricting the groundwater delay factor to 1 day, setting the baseflow recession factor to 0.011 day based on observed daily streamflow records, and introducing the contribution of the springs that recharge outside Spring Creek but discharge inside the watershed as point sources.

Topo-SWAT outperformed SWAT in modeling daily streamflow for a 12-year simulation period, with NSE of 0.73–0.79, PBIAS of –2.8 to –3.7%, and R² of 0.71–0.77 in the calibration and validation periods. Moreover, Topo-SWAT successfully reproduced the VSA
hydrology of the watershed using the wetness class distribution approach and predicted water quality adequately.

The calibrated Topo-SWAT model of the Spring Creek watershed was later used to simulate nutrient and sediment loadings under four dairy cropping scenarios of different land areas, feed production, and nutrient input strategies [171]. The model was also applied to evaluate the impacts of agricultural best management practices on nutrient water quality in the watershed [172]. Gunn et al. [173] further modified SWAT-VSA by integrating a daily dynamic time series of CO₂ into the model and implementing changes to the plant subroutines to additionally include flexible stomatal conductance and LAI parameters. The generated models, namely SWAT-VSA_CO₂ and SWAT-VSA_CO₂+Plant, were used to predict the impacts of climate change with increasing atmospheric CO₂ on the water balance of the Spring Creek watershed.

### 4.2.7. SWAT + Water Accounting Plus (WA+) Framework

Delavar et al. [174] coupled SWAT with the Water Accounting Plus (WA+) framework, providing a hybrid tool to analyze water resources and support macro and micro water planning in watersheds, based on past, current and future trends in water demand and supply. WA+ uses four sheets to assess water resources in a basin, including resource base, evapotranspiration, productivity, and withdrawal components. In order to populate these sheets with input data generated by hydrological simulation in SWAT, the SWAT source code was modified to: (1) simulate and report the daily groundwater level changes, and (2) model the interactions and exchanges between the aquifers located in different subbasins by overlaying the subbasin HRUs layer and aquifer boundaries during HRUs definition. The SWAT-WA+ tool was used to evaluate the trends in water supply and consumption in the Tashk-Bakhtegan karst watershed (27,520 km²) in Iran, where 60% of the irrigation demand is met by groundwater. Delavar et al. [175] remodified SWAT (v2012) and linked it with WA+, adding the crack flow module proposed by Eini et al. [36] with the other modifications previously implemented by Delavar et al. [174] to link SWAT and WA+. The model was applied to assess different conditions of water supply and demand under wet and dry periods in the karstified Karkheh River basin (42,267 km²) in Iran.

The customized SWAT-WA+ frameworks, namely SWAT-FARS for the Tashk-Bakhtegan basin and SWAT-Karkheh for Karkheh River basin, were calibrated for streamflow, groundwater levels, evapotranspiration, and crop yields using a multi-stage calibration process. Both models scored NSE and R² values higher than 0.5 for the streamflow and groundwater levels simulations in the calibration and validation periods. The outputs of the modified SWAT models were then used to run the WA+ framework. The results of the WA+ assessments revealed that the Tashk-Bakhtegan basin has suffered a decline in the volume of manageable water by more than 20% while irrigation demand increased by more than 50%, and that the volume of manageable water in the Karkheh River basin has dropped while groundwater abstraction increased by 17% due to climate change.

### 4.2.8. Overall Performance of Modified SWAT Models

Based on the reviewed literature, a total of 18 modified SWAT models have been developed and applied across 25 studies in watersheds characterized by karst geology (Table 2). Models that were run at daily and monthly time intervals reported a higher prediction efficiency of the flow at the monthly scale than the daily scale, both in calibration and validation periods [152,157,163,170]. Around 80% of the studies that used NSE to evaluate the hydrological model performance at the daily time step reported NSE values greater than 0.5 for the calibration and validation periods. In comparison, more than 80% of the studies that used the monthly time step scored NSE values higher than 0.5 for calibration, and more than 90% reported NSE values greater than 0.5 for validation (Figure 3).
hydrology of the watershed using the wetness class distribution approach and predicted water quality adequately.

The calibrated Topo-SWAT model of the Spring Creek watershed was later used to simulate nutrient and sediment loadings under four dairy cropping scenarios of different land areas, feed production, and nutrient input strategies [171]. The model was also applied to evaluate the impacts of agricultural best management practices on nutrient water quality in the watershed [172]. Gunn et al. [173] further modified SWAT-VSA by integrating a daily dynamic time series of CO₂ into the model and implementing changes to the plant subroutines to additionally include flexible stomatal conductance and LAI parameters. The generated models, namely SWAT-VSA_CO₂ and SWAT-VSA_CO₂+Plant, were used to predict the impacts of climate change with increasing atmospheric CO₂ on the water balance of the Spring Creek watershed.

### 4.2.7. SWAT + Water Accounting Plus (WA+) Framework

Delavar et al. [174] coupled SWAT with the Water Accounting Plus (WA+) framework, providing a hybrid tool to analyze water resources and support macro and micro water planning in watersheds, based on past, current and future trends in water demand and supply. WA+ uses four sheets to assess water resources in a basin, including resource base, evapotranspiration, productivity, and withdrawal components. In order to populate these sheets with input data generated by hydrological simulation in SWAT, the SWAT source code was modified to: (1) simulate and report the daily groundwater level changes, and (2) model the interactions and exchanges between the aquifers located in different subbasins by overlaying the subbasin HRUs layer and aquifer boundaries during HRUs definition. The SWAT-WA+ tool was used to evaluate the trends in water supply and consumption in the Tashk-Bakhtegan karst watershed (27,520 km²) in Iran, where 60% of the irrigation demand is met by groundwater. Delavar et al. [175] remodified SWAT (v2012) and linked it with WA+, adding the crack flow module proposed by Eini et al. [36] with the other modifications previously implemented by Delavar et al. [174] to link SWAT and WA+. The model was applied to assess different conditions of water supply and demand under wet and dry periods in the karstified Karkheh River basin (42,267 km²) in Iran.

The customized SWAT-WA+ frameworks, namely SWAT-FARS for the Tashk-Bakhtegan basin and SWAT-Karkheh for Karkheh River basin, were calibrated for streamflow, groundwater levels, evapotranspiration, and crop yields using a multi-stage calibration process. Both models scored NSE and R² values higher than 0.5 for the streamflow and groundwater levels simulations in the calibration and validation periods. The outputs of the modified SWAT models were then used to run the WA+ framework. The results of the WA+ assessments revealed that the Tashk-Bakhtegan basin has suffered a decline in the volume of manageable water by more than 20% while irrigation demand increased by more than 50%, and that the volume of manageable water in the Karkheh River basin has dropped while groundwater abstraction increased by 17% due to climate change.

### 4.2.8. Overall Performance of Modified SWAT Models

Based on the reviewed literature, a total of 18 modified SWAT models have been developed and applied across 25 studies in watersheds characterized by karst geology (Table 2). Models that were run at daily and monthly time intervals reported a higher prediction efficiency of the flow at the monthly scale than the daily scale, both in calibration and validation periods [152,157,163,170]. Around 80% of the studies that used NSE to evaluate the hydrological model performance at the daily time step reported NSE values greater than 0.5 for the calibration and validation periods. In comparison, more than 80% of the studies that used the monthly time step scored NSE values higher than 0.5 for calibration, and more than 90% reported NSE values greater than 0.5 for validation (Figure 3).
<figure><img src="image_7.png" /><figcaption>**Figure 3.** Number of modified SWAT-based studies in karst watersheds under the NSE performance ratings recommended by Moriasi et al. [102] for daily and monthly discharge simulation.</figcaption></figure>

Several applications showed that accounting for external flows from sinkholes and IGF, in conjunction with the implementation of reservoir models in SWAT, is needed to achieve an adequate representation of karstic flows [155,157,158,163,166,167,170,174]. Overall, coupling SWAT with karst flow models and adding external functions that reproduce karst features and processes to the SWAT source code have resulted in semi-distributed karstic flow models that simulated discharge with a comparable or better prediction efficiency than standard SWAT [36,155,156,167–169] while accounting for the dynamics of the different components in a karst system (when applicable).

However, poor daily and monthly performance statistics were still reported by Afinowicz et al. [152], Baffaut and Benson [153], and Malagò et al. [165] following the modification of the SWAT code. These results which suggest that modified approaches applied to the groundwater recharge and reservoir functions in SWAT may not always guarantee a successful simulation of the flow in complex karstic environments. For instance, Afinowicz et al. [152] indicated that additional modifications to the SWAT flow modules could be required to adequately simulate the large volumes of surface runoff and return flow during flood events as opposed to baseflow during low flow periods. Additionally, we could not identify published studies in which the modified SWAT models were applied across different watersheds, with the exception of the Karst-SWAT model [157] which has been used in the Koiliaris River Basin and the Island of Crete for a multitude of applications (i.e., hydrological and geochemical analyses, climate change impacts, management practices). Thus, one could infer that these models have not been widely tested in other karst basins or they only worked for watersheds with comparable geomorphological and hydrogeological characteristics to the basins in which they were initially applied. Consequently, the modified SWAT models listed in Table 2 can be further improved to have a better representation of karstic heterogeneity and non-linearity in their structure.

# 5. Recommendations and Perspectives for Future Research

## 5.1. *Constraints of Modified SWAT Models in Representing Karst Flow Characteristics*

Highlighting the constraints of the modified SWAT models (Table 2) would be the initial step to enhance their adaptability to other karst watersheds with complex surface water-groundwater hydrodynamics.

First, the Karst-SWAT [157] and KSWAT [165] two-reservoir models do not consider the function of the epikarst and do not explicitly include the diffuse and concentrated recharge components of infiltration from the karstic soils to the deep aquifer reservoir in
<figure><img src="image_7.png" /><figcaption>**Figure 3.** Number of modified SWAT-based studies in karst watersheds under the NSE performance ratings recommended by Moriasi et al. [102] for daily and monthly discharge simulation.</figcaption></figure>

Several applications showed that accounting for external flows from sinkholes and IGF, in conjunction with the implementation of reservoir models in SWAT, is needed to achieve an adequate representation of karstic flows [155,157,158,163,166,167,170,174]. Overall, coupling SWAT with karst flow models and adding external functions that reproduce karst features and processes to the SWAT source code have resulted in semi-distributed karstic flow models that simulated discharge with a comparable or better prediction efficiency than standard SWAT [36,155,156,167–169] while accounting for the dynamics of the different components in a karst system (when applicable).

However, poor daily and monthly performance statistics were still reported by Afinowicz et al. [152], Baffaut and Benson [153], and Malagò et al. [165] following the modification of the SWAT code. These results which suggest that modified approaches applied to the groundwater recharge and reservoir functions in SWAT may not always guarantee a successful simulation of the flow in complex karstic environments. For instance, Afinowicz et al. [152] indicated that additional modifications to the SWAT flow modules could be required to adequately simulate the large volumes of surface runoff and return flow during flood events as opposed to baseflow during low flow periods. Additionally, we could not identify published studies in which the modified SWAT models were applied across different watersheds, with the exception of the Karst-SWAT model [157] which has been used in the Koiliaris River Basin and the Island of Crete for a multitude of applications (i.e., hydrological and geochemical analyses, climate change impacts, management practices). Thus, one could infer that these models have not been widely tested in other karst basins or they only worked for watersheds with comparable geomorphological and hydrogeological characteristics to the basins in which they were initially applied. Consequently, the modified SWAT models listed in Table 2 can be further improved to have a better representation of karstic heterogeneity and non-linearity in their structure.

# 5. Recommendations and Perspectives for Future Research

## 5.1. *Constraints of Modified SWAT Models in Representing Karst Flow Characteristics*

Highlighting the constraints of the modified SWAT models (Table 2) would be the initial step to enhance their adaptability to other karst watersheds with complex surface water-groundwater hydrodynamics.

First, the Karst-SWAT [157] and KSWAT [165] two-reservoir models do not consider the function of the epikarst and do not explicitly include the diffuse and concentrated recharge components of infiltration from the karstic soils to the deep aquifer reservoir in
SWAT. Both models use the exponential decay weighting function to simulate recharge and outflows of the matrix and conduit reservoirs. However, the non-linear models are generally more suitable than the linear ones in representing the hydraulic behavior of karst systems, particularly during low flow periods and flood events [182–184]. Moreover, the two models follow the watershed surface delineation in SWAT to determine the recharge area of the spring, which may not always coincide with the groundwater recharge boundaries. This requires a tedious assessment of the karst areas recharging the springs outside the watershed in tandem with the introduction of the spring flow time series contributing to the watershed discharge as point sources in SWAT to simulate IGF.

Some of the Karst-SWAT and KSWAT modeling constraints were improved in the two-reservoir SWAT_IGF model [166], which simulates the hydrological processes in non-karst and karst regions as well as IGF in a single executable file. Although SWAT_IGF considers the dual recharge and storage functions in karst systems, it uses the linear storage-discharge relationship to model the outflows of the matrix and conduit reservoirs, and does not account for the function of the epikarst either. In SWAT_IGF, the exchange flow rate between the matrix and conduits is simulated as a diffuse net unidirectional flow from the matrix to the conduit. Yet, flow can be transferred from the conduit to the matrix, and vice versa, based on the change in the water level gradient between the two mediums. Additionally, SWAT_IGF models the spring flow contribution to discharge in karst areas as a single outflow from the conduit reservoir, whereas both the matrix and conduits can contribute to the karst discharge with different flow regimes (slow matrix discharge during low flow periods and fast conduit discharge during heavy rainfall seasons).

Next, the three-reservoir models developed by Wang et al. [167] and Geng et al. [168] incorporate the epikarst, matrix, and conduit functions, and thus represent a complete underground karst system. Yet, their main constraint is that fluxes between the reservoirs are simulated using a linear storage-discharge relationship.

SWAT-B&B [153], SWAT-karst [154], KarstSWAT [155], and the karst model developed by Zhou et al. [156] can simulate fast infiltration in karst watersheds dominated by spring flow fed by sinkholes. However, they all apply the linear reservoir of SWAT to model groundwater flow without considering the different storage and recharge functions of the main karst components (epikarst, matrix, and conduits). In comparison with SWAT-B&B, SWAT-karst, and the modified SWAT model by Zhou et al. [156], which all rely on the pond and wetland modules of SWAT for the simulation of flow through the sinkholes, only KarstSWAT can simulate IGF. With KarstSWAT, groundwater basins that drain the aquifer water to the spring can be identified during watershed delineation and included in a user defined input file to route the total recharge to the spring. Nonetheless, the model may only be applied in karst watersheds where sinkholes are solely located along the river streambed. Moreover, KarstSWAT assumed that flow from the losing streams directly recharges the deep aquifer through the sinkholes and emerges at the spring within a day, as the length of the flow path from the aquifer to the spring in the study watershed was unknown. Although this assumption was suitable for a specific study basin, additional data on storage, time of travel, and flow diversion pathways would have been required to simulate discharge in other watersheds. The study by Yactayo [154] corroborates this finding, as field investigations were needed to improve the model performance by determining whether the sinkholes route the flow within the study watershed or divert it outside of the watershed.

The SWAT-ML and SWAT-CF dedicated to watersheds with crack/preferential flow [36], the Topo-SWAT specific to watersheds with variable surface area hydrology [170], and the SWAT-WA+ models [174,175] may be directly applied to basins affected by karst hydrology or other rapid infiltration phenomena. Nevertheless, they do not represent the underground flow dynamics (epikarst-matrix-conduits) of karst aquifers either.

Finally, the non-linear ISWAT model [169] does not account for the diffuse/slow recharge and concentrated/rapid recharge into karst aquifer systems. In addition, it does not explicitly represent the storage and discharge functions of three main subsystems in
SWAT. Both models use the exponential decay weighting function to simulate recharge and outflows of the matrix and conduit reservoirs. However, the non-linear models are generally more suitable than the linear ones in representing the hydraulic behavior of karst systems, particularly during low flow periods and flood events [182–184]. Moreover, the two models follow the watershed surface delineation in SWAT to determine the recharge area of the spring, which may not always coincide with the groundwater recharge boundaries. This requires a tedious assessment of the karst areas recharging the springs outside the watershed in tandem with the introduction of the spring flow time series contributing to the watershed discharge as point sources in SWAT to simulate IGF.

Some of the Karst-SWAT and KSWAT modeling constraints were improved in the two-reservoir SWAT_IGF model [166], which simulates the hydrological processes in non-karst and karst regions as well as IGF in a single executable file. Although SWAT_IGF considers the dual recharge and storage functions in karst systems, it uses the linear storage-discharge relationship to model the outflows of the matrix and conduit reservoirs, and does not account for the function of the epikarst either. In SWAT_IGF, the exchange flow rate between the matrix and conduits is simulated as a diffuse net unidirectional flow from the matrix to the conduit. Yet, flow can be transferred from the conduit to the matrix, and vice versa, based on the change in the water level gradient between the two mediums. Additionally, SWAT_IGF models the spring flow contribution to discharge in karst areas as a single outflow from the conduit reservoir, whereas both the matrix and conduits can contribute to the karst discharge with different flow regimes (slow matrix discharge during low flow periods and fast conduit discharge during heavy rainfall seasons).

Next, the three-reservoir models developed by Wang et al. [167] and Geng et al. [168] incorporate the epikarst, matrix, and conduit functions, and thus represent a complete underground karst system. Yet, their main constraint is that fluxes between the reservoirs are simulated using a linear storage–discharge relationship.

SWAT-B&B [153], SWAT-karst [154], KarstSWAT [155], and the karst model developed by Zhou et al. [156] can simulate fast infiltration in karst watersheds dominated by spring flow fed by sinkholes. However, they all apply the linear reservoir of SWAT to model groundwater flow without considering the different storage and recharge functions of the main karst components (epikarst, matrix, and conduits). In comparison with SWAT-B&B, SWAT-karst, and the modified SWAT model by Zhou et al. [156], which all rely on the pond and wetland modules of SWAT for the simulation of flow through the sinkholes, only KarstSWAT can simulate IGF. With KarstSWAT, groundwater basins that drain the aquifer water to the spring can be identified during watershed delineation and included in a user defined input file to route the total recharge to the spring. Nonetheless, the model may only be applied in karst watersheds where sinkholes are solely located along the river streambed. Moreover, KarstSWAT assumed that flow from the losing streams directly recharges the deep aquifer through the sinkholes and emerges at the spring within a day, as the length of the flow path from the aquifer to the spring in the study watershed was unknown. Although this assumption was suitable for a specific study basin, additional data on storage, time of travel, and flow diversion pathways would have been required to simulate discharge in other watersheds. The study by Yactayo [154] corroborates this finding, as field investigations were needed to improve the model performance by determining whether the sinkholes route the flow within the study watershed or divert it outside of the watershed.

The SWAT-ML and SWAT-CF dedicated to watersheds with crack/preferential flow [36], the Topo-SWAT specific to watersheds with variable surface area hydrology [170], and the SWAT-WA+ models [174,175] may be directly applied to basins affected by karst hydrology or other rapid infiltration phenomena. Nevertheless, they do not represent the underground flow dynamics (epikarst-matrix-conduits) of karst aquifers either.

Finally, the non-linear ISWAT model [169] does not account for the diffuse/slow recharge and concentrated/rapid recharge into karst aquifer systems. In addition, it does not explicitly represent the storage and discharge functions of three main subsystems in
karst (epikarst, matrix, and conduits) due to the lumped feature of the reservoir model in SWAT.

## 5.2. *Future Research Areas for the Application of Modified SWAT Models*

Climate change effects on karst hydrology and water quality were investigated with the modified SWAT models proposed by Nikolaidis et al. [157], Nerantzaki et al. [158], Tapoglou et al. [159], Nerantzaki et al. [161], and Gunn et al. [173]. Additionally, some studies evaluated water resources in a karst watershed for different trends in water supply and consumption [174] as well as the joint impacts of climate change and groundwater withdrawal on water resources availability [175]. Other studies simulated in-stream water quality under different agricultural management practices [170–172].

However, the use of modified SWAT models for the integrated understanding of the critical zone processes and the quantification of the impacts of evapotranspiration and vegetation cover change on karst water resources are still lacking. Among the studies that applied a modified SWAT code in a karst basin (Table 2), Lilli et al. [164] confirmed the conceptualization of the two reservoir well-mixed karstic system with Karst-SWAT, using geomorphologic and tidal analyses. Then, the authors applied Karst-SWAT to simulate the hydrological budget and pathways in the critical zone and investigate the response of karst during extreme events.

On the other hand, Nguyen et al. [166] investigated the impact of evapotranspiration on the hydrological performance of SWAT_IGF by using satellite-derived ETa (Moderate Resolution Imaging Spectroradiometer MOD16 ETa at 8-day time step and 1 km² spatial resolution) in tandem with multi-site streamflow records to calibrate the model. ET was simulated in SWAT_IGF using the Penman-Monteith approach. Results of the NSE index showed that the use of ETa as an additional variable for the calibration of discharge had little to no effect on the model performance in the study watershed, compared with the case in which multi-gage streamflow calibration was carried out separately. In addition, a strong positive correlation was established between the MOD16 ETa data and SWAT-simulated ET with only streamflow used in calibration. As the model performance for ETa was shown to improve with the improvement of the model performance for streamflow, the use of MOD16 ETa for model calibration was disregarded. The authors maintained, nonetheless, that these findings should not be generalized to other remote sensing products and to studies in other karst areas, considering more research on the use of ETa for the calibration of karst hydrological models and streamflow estimation is needed.

Moreover, only Afinowicz et al. [152] used a modified SWAT code to predict the impact of land-use change on the hydrological cycle of a karst area. In particular, the authors developed scenarios for brush control strategies in function of slope, rangeland cover density, and soil depth, to determine the most favorable areas for increasing the watershed water yield. All scenarios of the brush reduction cover resulted in a decrease in evapotranspiration and increases in surface runoff, baseflow, and deep aquifer recharge, with the greatest effect observed on recharge.

Therefore, future research in the realm of karst hydrological modeling should integrate spatially distributed ET data from remote sensing models that account for the dynamics of the land use [33] with multi-source precipitation data derived from ground-based observations or satellite products [147,148]. This approach could improve the spatial distribution of aquifer recharge and the overall rainfall-discharge relationship in karstic watersheds.

Other areas of future research could include testing the capabilities of the newly released SWAT+ version in simulating discharge in karst watersheds, particularly extreme flows (peak and low flows), and comparing the performance efficiency of SWAT+ to previous SWAT versions. The performance of SWAT should also be compared to other modeling approaches used in karst hydrological applications [28] in order to improve its representation of the high and low flows sustained by karst features. Additionally, it is recommended to model the rainfall-discharge relationship in highly dynamic karst aquifers using subdaily time intervals (e.g., hourly time step) in order to reach a better
karst (epikarst, matrix, and conduits) due to the lumped feature of the reservoir model in SWAT.

## 5.2. *Future Research Areas for the Application of Modified SWAT Models*

Climate change effects on karst hydrology and water quality were investigated with the modified SWAT models proposed by Nikolaidis et al. [157], Nerantzaki et al. [158], Tapoglou et al. [159], Nerantzaki et al. [161], and Gunn et al. [173]. Additionally, some studies evaluated water resources in a karst watershed for different trends in water supply and consumption [174] as well as the joint impacts of climate change and groundwater withdrawal on water resources availability [175]. Other studies simulated in-stream water quality under different agricultural management practices [170–172].

However, the use of modified SWAT models for the integrated understanding of the critical zone processes and the quantification of the impacts of evapotranspiration and vegetation cover change on karst water resources are still lacking. Among the studies that applied a modified SWAT code in a karst basin (Table 2), Lilli et al. [164] confirmed the conceptualization of the two reservoir well-mixed karstic system with Karst-SWAT, using geomorphologic and tidal analyses. Then, the authors applied Karst-SWAT to simulate the hydrological budget and pathways in the critical zone and investigate the response of karst during extreme events.

On the other hand, Nguyen et al. [166] investigated the impact of evapotranspiration on the hydrological performance of SWAT_IGF by using satellite-derived ETa (Moderate Resolution Imaging Spectroradiometer MOD16 ETa at 8-day time step and 1 km² spatial resolution) in tandem with multi-site streamflow records to calibrate the model. ET was simulated in SWAT_IGF using the Penman-Monteith approach. Results of the NSE index showed that the use of ETa as an additional variable for the calibration of discharge had little to no effect on the model performance in the study watershed, compared with the case in which multi-gage streamflow calibration was carried out separately. In addition, a strong positive correlation was established between the MOD16 ETa data and SWAT-simulated ET with only streamflow used in calibration. As the model performance for ETa was shown to improve with the improvement of the model performance for streamflow, the use of MOD16 ETa for model calibration was disregarded. The authors maintained, nonetheless, that these findings should not be generalized to other remote sensing products and to studies in other karst areas, considering more research on the use of ETa for the calibration of karst hydrological models and streamflow estimation is needed.

Moreover, only Afinowicz et al. [152] used a modified SWAT code to predict the impact of land-use change on the hydrological cycle of a karst area. In particular, the authors developed scenarios for brush control strategies in function of slope, rangeland cover density, and soil depth, to determine the most favorable areas for increasing the watershed water yield. All scenarios of the brush reduction cover resulted in a decrease in evapotranspiration and increases in surface runoff, baseflow, and deep aquifer recharge, with the greatest effect observed on recharge.

Therefore, future research in the realm of karst hydrological modeling should integrate spatially distributed ET data from remote sensing models that account for the dynamics of the land use [33] with multi-source precipitation data derived from ground-based observations or satellite products [147,148]. This approach could improve the spatial distribution of aquifer recharge and the overall rainfall-discharge relationship in karstic watersheds.

Other areas of future research could include testing the capabilities of the newly released SWAT+ version in simulating discharge in karst watersheds, particularly extreme flows (peak and low flows), and comparing the performance efficiency of SWAT+ to previous SWAT versions. The performance of SWAT should also be compared to other modeling approaches used in karst hydrological applications [28] in order to improve its representation of the high and low flows sustained by karst features. Additionally, it is recommended to model the rainfall-discharge relationship in highly dynamic karst aquifers using subdaily time intervals (e.g., hourly time step) in order to reach a better
prediction of the flood peak discharge during high rainfall events. Assessing the discharge at lower time series can improve the mitigation of karst flash floods at the spring outlet and the management of groundwater storage for future water supply [185]. Finally, future studies could focus on developing solute transport models that incorporate the different components and flow dynamics in karst hydrosystems.

**Author Contributions:** Conceptualization, I.A.K., D.L. and L.B.; methodology, I.A.K., D.L. and L.B.; data curation, I.A.K.; writing—original draft preparation, I.A.K.; writing—review and editing, I.A.K., L.B. and D.L.; supervision, D.L. and L.B.; funding acquisition, D.L. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research has been supported through the grant EUR TESS N°ANR-18-EURE-0018 in the framework of the Programme des Investissements d'Avenir and by the Ministry of Higher Education, Research and Innovation (Ministère de l'Enseignement Supérieur, de la Recherche et de l'Innovation, MESRI) of France.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare no conflict of interest.

## Appendix A

<table><caption>**Table A1.** List of SWAT Variables in the Main Text.</caption><thead><tr><td>Variable</td><td>Definition</td></tr></thead><tbody><tr><td>SW</td><td>Soil water content</td></tr><tr><td>P</td><td>Precipitation</td></tr><tr><td>Q<sub>surf</sub></td><td>Surface runoff</td></tr><tr><td>PET</td><td>Potential evapotranspiration</td></tr><tr><td>ET<sub>a</sub></td><td>Actual evapotranspiration</td></tr><tr><td>W<sub>seep</sub></td><td>Percolation and bypass flow drained from the soil profile</td></tr><tr><td>Q<sub>gw</sub></td><td>Return flow</td></tr><tr><td>Q<sub>lat</sub></td><td>Soil lateral flow</td></tr><tr><td>Q<sub>tile</sub></td><td>Tile flow</td></tr><tr><td>T<sub>loss</sub></td><td>Channel transmission losses</td></tr><tr><td>WYLD</td><td>Water yield</td></tr><tr><td>Irr</td><td>Irrigation</td></tr><tr><td>LAI</td><td>Lead area index</td></tr><tr><td>λ</td><td>Latent heat of vaporization</td></tr><tr><td>E</td><td>Depth rate evaporation</td></tr><tr><td>Δ</td><td>Slope of the saturation vapor pressure-temperature curve</td></tr><tr><td>H<sub>net</sub></td><td>Net radiation</td></tr><tr><td>G</td><td>Heat flux density to the ground</td></tr><tr><td>ρ<sub>air</sub></td><td>Air density</td></tr><tr><td>c<sub>p</sub></td><td>Specific heat at constant pressure</td></tr><tr><td>e<sub>0</sub></td><td>Saturation vapor pressure of air at height z</td></tr><tr><td>e<sub>z</sub></td><td>Water vapor pressure of air at height z</td></tr><tr><td>γ</td><td>Psychometric constant</td></tr><tr><td>r<sub>c</sub></td><td>Plant canopy resistance</td></tr><tr><td>r<sub>a</sub></td><td>Aerodynamic resistance</td></tr><tr><td>S</td><td>Soil moisture retention parameter</td></tr></tbody></table>
prediction of the flood peak discharge during high rainfall events. Assessing the discharge at lower time series can improve the mitigation of karst flash floods at the spring outlet and the management of groundwater storage for future water supply [185]. Finally, future studies could focus on developing solute transport models that incorporate the different components and flow dynamics in karst hydrosystems.

**Author Contributions:** Conceptualization, I.A.K., D.L. and L.B.; methodology, I.A.K., D.L. and L.B.; data curation, I.A.K.; writing—original draft preparation, I.A.K.; writing—review and editing, I.A.K., L.B. and D.L.; supervision, D.L. and L.B.; funding acquisition, D.L. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research has been supported through the grant EUR TESS N°ANR-18-EURE-0018 in the framework of the Programme des Investissements d'Avenir and by the Ministry of Higher Education, Research and Innovation (Ministère de l'Enseignement Supérieur, de la Recherche et de l'Innovation, MESRI) of France.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare no conflict of interest.

## Appendix A

<table><caption>**Table A1.** List of SWAT Variables in the Main Text.</caption><thead><tr><td>Variable</td><td>Definition</td></tr></thead><tbody><tr><td>SW</td><td>Soil water content</td></tr><tr><td>P</td><td>Precipitation</td></tr><tr><td>Q<sub>surf</sub></td><td>Surface runoff</td></tr><tr><td>PET</td><td>Potential evapotranspiration</td></tr><tr><td>ET<sub>a</sub></td><td>Actual evapotranspiration</td></tr><tr><td>W<sub>seep</sub></td><td>Percolation and bypass flow drained from the soil profile</td></tr><tr><td>Q<sub>gw</sub></td><td>Return flow</td></tr><tr><td>Q<sub>lat</sub></td><td>Soil lateral flow</td></tr><tr><td>Q<sub>tile</sub></td><td>Tile flow</td></tr><tr><td>T<sub>loss</sub></td><td>Channel transmission losses</td></tr><tr><td>WYLD</td><td>Water yield</td></tr><tr><td>Irr</td><td>Irrigation</td></tr><tr><td>LAI</td><td>Lead area index</td></tr><tr><td>λ</td><td>Latent heat of vaporization</td></tr><tr><td>E</td><td>Depth rate evaporation</td></tr><tr><td>Δ</td><td>Slope of the saturation vapor pressure-temperature curve</td></tr><tr><td>H<sub>net</sub></td><td>Net radiation</td></tr><tr><td>G</td><td>Heat flux density to the ground</td></tr><tr><td>ρ<sub>air</sub></td><td>Air density</td></tr><tr><td>c<sub>p</sub></td><td>Specific heat at constant pressure</td></tr><tr><td>e<sub>0</sub></td><td>Saturation vapor pressure of air at height z</td></tr><tr><td>e<sub>z</sub></td><td>Water vapor pressure of air at height z</td></tr><tr><td>γ</td><td>Psychometric constant</td></tr><tr><td>r<sub>c</sub></td><td>Plant canopy resistance</td></tr><tr><td>r<sub>a</sub></td><td>Aerodynamic resistance</td></tr><tr><td>S</td><td>Soil moisture retention parameter</td></tr></tbody></table>
<table><caption>Table A1. Cont.</caption><thead><tr><th>Variable</th><th>Definition</th></tr></thead><tbody><tr><td>CN</td><td>Curve number</td></tr><tr><td>I<sub>a</sub></td><td>Initial water abstraction prior to runoff</td></tr><tr><td>I(t)</td><td>Infiltration</td></tr><tr><td>I<sub>acc</sub>(t)</td><td>Cumulative infiltration after ponding</td></tr><tr><td>ψ</td><td>Wetting front matric potential</td></tr><tr><td>Δθ</td><td>Change in soil moisture content</td></tr><tr><td>K<sub>e</sub></td><td>Effective hydraulic conductivity</td></tr><tr><td>q<sub>peak</sub></td><td>Peak runoff rate</td></tr><tr><td>α<sub>tc</sub></td><td>Fraction of daily rainfall during the time of concentration</td></tr><tr><td>A</td><td>Subbasin area</td></tr><tr><td>t<sub>conc</sub></td><td>Time of concentration for the subbasin</td></tr><tr><td>K<sub>ch</sub></td><td>Effective hydraulic conductivity of the channel alluvium</td></tr><tr><td>L<sub>ch</sub></td><td>Channel length</td></tr><tr><td>P<sub>ch</sub></td><td>Wetted perimeter in the channel</td></tr><tr><td>TT</td><td>Flow travel time</td></tr><tr><td>K<sub>sat</sub></td><td>Saturated hydraulic conductivity</td></tr><tr><td>SAT</td><td>Saturation water content</td></tr><tr><td>FC</td><td>Field capacity water content</td></tr><tr><td>W<sub>perc,ly</sub></td><td>Water percolation to an underlying soil layer</td></tr><tr><td>W<sub>perc,ly=n</sub></td><td>Water percolating out of the lowest soil layer</td></tr><tr><td>SW<sub>ly,excess</sub></td><td>Drainable volume of water in the soil layer</td></tr><tr><td>TT<sub>perc</sub></td><td>Travel time through the soil layer</td></tr><tr><td>W<sub>crk,btm</sub></td><td>Crack flow past the lower boundary of the soil profile</td></tr><tr><td>crk</td><td>Total crack volume for the soil profile</td></tr><tr><td>crk<sub>ly=nn</sub></td><td>Crack volume for the deepest soil layer</td></tr><tr><td>depth<sub>ly=nn</sub></td><td>Depth of the deepest soil layer</td></tr><tr><td>slp</td><td>Increase in elevation per unit distance</td></tr><tr><td>Ø<sub>d</sub></td><td>Drainable (residual) porosity of the soil layer</td></tr><tr><td>L<sub>hill</sub></td><td>Hillslope length</td></tr><tr><td>E<sub>e,ly</sub></td><td>Evaporation drawn from the soil layer</td></tr><tr><td>E<sub>t,ly</sub></td><td>Transpiration drawn from the soil layer</td></tr><tr><td>δ<sub>gw</sub></td><td>Delay time for recharge to reach the aquifers</td></tr><tr><td>W<sub>rchrg</sub></td><td>Recharge to the aquifers</td></tr><tr><td>W<sub>rchrg,sh</sub></td><td>Recharge to the shallow aquifer</td></tr><tr><td>W<sub>rchrg,dp</sub></td><td>Recharge to the deep aquifer</td></tr><tr><td>Q<sub>gw,sh</sub></td><td>Return flow from the shallow aquifer</td></tr><tr><td>Q<sub>gw,dp</sub></td><td>Groundwater flow from the deep aquifer</td></tr><tr><td>α<sub>gw,sh</sub></td><td>Groundwater recession constant of the shallow aquifer</td></tr><tr><td>α<sub>gw,dp</sub></td><td>Groundwater recession constant of the deep aquifer</td></tr><tr><td>β<sub>dp</sub></td><td>Coefficient of percolation to the deep aquifer</td></tr><tr><td>aq<sub>sh</sub></td><td>Water level in the shallow aquifer</td></tr><tr><td>aq<sub>shthr,q</sub></td><td>Threshold water level in the shallow aquifer for return flow</td></tr><tr><td>ET<sub>rvp</sub></td><td>Actual groundwater evapotranspiration by plant uptake/capillary rise</td></tr><tr><td>ET<sub>rvp,max</sub></td><td>Maximum water removed from the shallow aquifer</td></tr><tr><td>β<sub>rvp</sub></td><td>Groundwater evaporation coefficient</td></tr><tr><td>aq<sub>shthr,rvp</sub></td><td>Threshold water level in the shallow aquifer for groundwater evaporation to occur</td></tr><tr><td>Q<sub>pump,sh</sub></td><td>Water withdrawn by pumping from the shallow aquifer</td></tr><tr><td>V<sub>seep</sub></td><td>Seepage from the pond or wetland</td></tr><tr><td>A<sub>wet</sub></td><td>Water surface area of the pond or wetland</td></tr><tr><td>V<sub>seep,pot</sub></td><td>Seepage from a pothole</td></tr><tr><td>SA</td><td>Pothole surface area</td></tr></tbody></table>
<table><caption>Table A1. Cont.</caption><thead><tr><th>Variable</th><th>Definition</th></tr></thead><tbody><tr><td>CN</td><td>Curve number</td></tr><tr><td>I<sub>a</sub></td><td>Initial water abstraction prior to runoff</td></tr><tr><td>I(t)</td><td>Infiltration</td></tr><tr><td>I<sub>acc</sub>(t)</td><td>Cumulative infiltration after ponding</td></tr><tr><td>ψ</td><td>Wetting front matric potential</td></tr><tr><td>Δθ</td><td>Change in soil moisture content</td></tr><tr><td>K<sub>e</sub></td><td>Effective hydraulic conductivity</td></tr><tr><td>q<sub>peak</sub></td><td>Peak runoff rate</td></tr><tr><td>α<sub>tc</sub></td><td>Fraction of daily rainfall during the time of concentration</td></tr><tr><td>A</td><td>Subbasin area</td></tr><tr><td>t<sub>conc</sub></td><td>Time of concentration for the subbasin</td></tr><tr><td>K<sub>ch</sub></td><td>Effective hydraulic conductivity of the channel alluvium</td></tr><tr><td>L<sub>ch</sub></td><td>Channel length</td></tr><tr><td>P<sub>ch</sub></td><td>Wetted perimeter in the channel</td></tr><tr><td>TT</td><td>Flow travel time</td></tr><tr><td>K<sub>sat</sub></td><td>Saturated hydraulic conductivity</td></tr><tr><td>SAT</td><td>Saturation water content</td></tr><tr><td>FC</td><td>Field capacity water content</td></tr><tr><td>W<sub>perc,ly</sub></td><td>Water percolation to an underlying soil layer</td></tr><tr><td>W<sub>perc,ly=n</sub></td><td>Water percolating out of the lowest soil layer</td></tr><tr><td>SW<sub>ly,excess</sub></td><td>Drainable volume of water in the soil layer</td></tr><tr><td>TT<sub>perc</sub></td><td>Travel time through the soil layer</td></tr><tr><td>W<sub>crk,btm</sub></td><td>Crack flow past the lower boundary of the soil profile</td></tr><tr><td>crk</td><td>Total crack volume for the soil profile</td></tr><tr><td>crk<sub>ly=nn</sub></td><td>Crack volume for the deepest soil layer</td></tr><tr><td>depth<sub>ly=nn</sub></td><td>Depth of the deepest soil layer</td></tr><tr><td>slp</td><td>Increase in elevation per unit distance</td></tr><tr><td>Ø<sub>d</sub></td><td>Drainable (residual) porosity of the soil layer</td></tr><tr><td>L<sub>hill</sub></td><td>Hillslope length</td></tr><tr><td>E<sub>e,ly</sub></td><td>Evaporation drawn from the soil layer</td></tr><tr><td>E<sub>t,ly</sub></td><td>Transpiration drawn from the soil layer</td></tr><tr><td>δ<sub>gw</sub></td><td>Delay time for recharge to reach the aquifers</td></tr><tr><td>W<sub>rchrg</sub></td><td>Recharge to the aquifers</td></tr><tr><td>W<sub>rchrg,sh</sub></td><td>Recharge to the shallow aquifer</td></tr><tr><td>W<sub>rchrg,dp</sub></td><td>Recharge to the deep aquifer</td></tr><tr><td>Q<sub>gw,sh</sub></td><td>Return flow from the shallow aquifer</td></tr><tr><td>Q<sub>gw,dp</sub></td><td>Groundwater flow from the deep aquifer</td></tr><tr><td>α<sub>gw,sh</sub></td><td>Groundwater recession constant of the shallow aquifer</td></tr><tr><td>α<sub>gw,dp</sub></td><td>Groundwater recession constant of the deep aquifer</td></tr><tr><td>β<sub>dp</sub></td><td>Coefficient of percolation to the deep aquifer</td></tr><tr><td>aq<sub>sh</sub></td><td>Water level in the shallow aquifer</td></tr><tr><td>aq<sub>shthr,q</sub></td><td>Threshold water level in the shallow aquifer for return flow</td></tr><tr><td>ET<sub>rvp</sub></td><td>Actual groundwater evapotranspiration by plant uptake/capillary rise</td></tr><tr><td>ET<sub>rvp,max</sub></td><td>Maximum water removed from the shallow aquifer</td></tr><tr><td>β<sub>rvp</sub></td><td>Groundwater evaporation coefficient</td></tr><tr><td>aq<sub>shthr,rvp</sub></td><td>Threshold water level in the shallow aquifer for groundwater evaporation to occur</td></tr><tr><td>Q<sub>pump,sh</sub></td><td>Water withdrawn by pumping from the shallow aquifer</td></tr><tr><td>V<sub>seep</sub></td><td>Seepage from the pond or wetland</td></tr><tr><td>A<sub>wet</sub></td><td>Water surface area of the pond or wetland</td></tr><tr><td>V<sub>seep,pot</sub></td><td>Seepage from a pothole</td></tr><tr><td>SA</td><td>Pothole surface area</td></tr></tbody></table>
## Appendix B. Afinowicz et al. [152]

Baseflow is computed using Equation (A1):

$$Q_{gw,i} = Q_{gw, i-1} \times e^{-\alpha_{gw} \times \Delta t} + (w_{rchrg} - w_{deep}) \times (1 - e^{-\alpha_{gw} \times \Delta t}) \quad (A1)$$

where $Q_{gw,i}$ and $Q_{gw, i-1}$ are the baseflow values for the current and previous day (mm H₂O.day⁻¹), $\alpha_{gw}$ is baseflow recession constant, $\Delta t$ is the time interval (days), $w_{rchrg}$ is the water percolated past the root zone (mm H₂O.day⁻¹), and $w_{deep}$ is the water percolated to the deep aquifer (mm H₂O.day⁻¹).

## Appendix C. Baffaut and Benson [153]

Slow recharge from the bottom of the soil profile and fast recharge in karst areas are computed by Equations (A2) and (A3), respectively:

$$rchr g_{seep(t)} = \left(1 - e^{\frac{-1}{g_{delay}}}\right) soil_{seep(t)} + \left(e^{\frac{-1}{g_{delay}}}\right) rchr g_{seep(t-1)} \quad (A2)$$

$$rchr g_{karst(t)} = \left(1 - e^{\frac{-1}{karst_{delay}}}\right) karst_{seep(t)} + \left(e^{\frac{-1}{karst_{delay}}}\right) rchr g_{karst(t-1)} \quad (A3)$$

where $rchr g_{seep(t)}$ and $rchr g_{karst(t)}$ are the slow recharge from the soil layers and fast recharge from the sinkholes and losing streams to the aquifer on a given day (t) (mm H₂O.day⁻¹), respectively, $soil_{seep(t)}$ and $karst_{seep(t)}$ are the percolation from the soil bottom and losses from the sinkholes and losing streams (mm H₂O.day⁻¹), respectively, $gw_{delay}$ and $karst_{delay}$ represent the time delay for the water percolating from the soil bottom and water infiltrating from the sinkholes and riverbeds to reach the aquifer (days), respectively, and $rchr g_{seep(t-1)}$ and $rchr g_{karst(t-1)}$ are the recharge values of the previous day (mm H₂O.day⁻¹).

Daily return flow is calculated using the standard equation of SWAT as a function of groundwater flow of the previous day and aquifer recharge of that day, as follows:

$$\left\{ \begin{aligned} Q_{gw,t} &= rchr g_{(t)} \times (1 - e^{-\alpha_{gw} \times \Delta t}) + Q_{gw,t-1} \times e^{-\alpha_{gw} \times \Delta t} \\ rchr g_{(t)} &= rchr g_{seep(t)} + rchr g_{karst(t)} \end{aligned} \right. \quad (A4)$$

where $Q_{gw,t}$ is the daily return flow (mm H₂O.day⁻¹), $rchr g_{(t)}$ is the total aquifer recharge calculated as the sum of the slow recharge from the soil layers ($rchr g_{seep(t)}$) and fast recharge from sinkholes and losing streams to the aquifer ($rchr g_{karst(t)}$) (mm H₂O.day⁻¹), and $\alpha_{gw}$ is the baseflow recession coefficient (days⁻¹).

## Appendix D. Yactayo [154]

In non-karst regions, the recharge to the unconfined aquifer is calculated using Equation (A5):

$$\left\{ \begin{aligned} rchr g_{seep(t)} &= \left(1 - e^{\frac{-1}{gw_{delay}}}\right) seep + \left(e^{\frac{-1}{gw_{delay}}}\right) rchr g_{seep(t-1)} \\ seep &= sepbm_{(t)} + twlpnd + twlwet \end{aligned} \right. \quad (A5)$$

where $rchr g_{seep(t)}$ and $rchr g_{seep(t-1)}$ represent the recharge from the water percolating from the soil bottom to the aquifer on a given day (t) and the day before (mm H₂O.day⁻¹), respectively, $seep$ is percolating water and seepage from impoundments on a given day (mm H₂O.day⁻¹), $sepbm_{(t)}$ is the water drained from the bottom of the soil profile (mm H₂O.day⁻¹), $twlpnd$ is the seepage from the ponds (mm H₂O.day⁻¹), $twlwet$ is the seepage
## Appendix B. Afinowicz et al. [152]

Baseflow is computed using Equation (A1):

$$Q_{gw,i} = Q_{gw, i-1} \times e^{-\alpha_{gw} \times \Delta t} + (w_{rchrg} - w_{deep}) \times (1 - e^{-\alpha_{gw} \times \Delta t}) \quad (A1)$$

where $Q_{gw,i}$ and $Q_{gw, i-1}$ are the baseflow values for the current and previous day (mm H₂O.day⁻¹), $\alpha_{gw}$ is baseflow recession constant, $\Delta t$ is the time interval (days), $w_{rchrg}$ is the water percolated past the root zone (mm H₂O.day⁻¹), and $w_{deep}$ is the water percolated to the deep aquifer (mm H₂O.day⁻¹).

## Appendix C. Baffaut and Benson [153]

Slow recharge from the bottom of the soil profile and fast recharge in karst areas are computed by Equations (A2) and (A3), respectively:

$$rchr_{g_{seep}(t)} = \left(1 - e^{-\frac{1}{g_{w_{delay}}}t}\right) soil_{seep(t)} + \left(e^{-\frac{1}{g_{w_{delay}}}t}\right) rchr_{g_{seep}(t-1)} \quad (A2)$$

$$rchr_{g_{karst}(t)} = \left(1 - e^{-\frac{1}{karst_{delay}}t}\right) karst_{seep(t)} + \left(e^{-\frac{1}{karst_{delay}}t}\right) rchr_{g_{karst}(t-1)} \quad (A3)$$

where $rchr_{g_{seep}(t)}$ and $rchr_{g_{karst}(t)}$ are the slow recharge from the soil layers and fast recharge from the sinkholes and losing streams to the aquifer on a given day (t) (mm H₂O.day⁻¹), respectively, $soil_{seep(t)}$ and $karst_{seep(t)}$ are the percolation from the soil bottom and losses from the sinkholes and losing streams (mm H₂O.day⁻¹), respectively, $gw_{delay}$ and $karst_{delay}$ represent the time delay for the water percolating from the soil bottom and water infiltrating from the sinkholes and riverbeds to reach the aquifer (days), respectively, and $rchr_{g_{seep}(t-1)}$ and $rchr_{g_{karst}(t-1)}$ are the recharge values of the previous day (mm H₂O.day⁻¹).

Daily return flow is calculated using the standard equation of SWAT as a function of groundwater flow of the previous day and aquifer recharge of that day, as follows:

$$\left\{ \begin{aligned} Q_{gw,t} &= rchr_{g_{(t)}} \times (1 - e^{-\alpha_{gw} \times \Delta t}) + Q_{gw,t-1} \times e^{-\alpha_{gw} \times \Delta t} \\ rchr_{g_{(t)}} &= rchr_{g_{seep(t)}} + rchr_{g_{karst(t)}} \end{aligned} \right. \quad (A4)$$

where $Q_{gw,t}$ is the daily return flow (mm H₂O.day⁻¹), $rchr_{g_{(t)}}$ is the total aquifer recharge calculated as the sum of the slow recharge from the soil layers ($rchr_{g_{seep(t)}}$) and fast recharge from sinkholes and losing streams to the aquifer ($rchr_{g_{karst(t)}}$) (mm H₂O.day⁻¹), and $\alpha_{gw}$ is the baseflow recession coefficient (days⁻¹).

## Appendix D. Yactayo [154]

In non-karst regions, the recharge to the unconfined aquifer is calculated using Equation (A5):

$$\left\{ \begin{aligned} rchr_{g_{seep}(t)} &= \left(1 - e^{-\frac{1}{gw_{delay}}t}\right) seep + \left(e^{-\frac{1}{gw_{delay}}t}\right) rchr_{g_{seep}(t-1)} \\ seep &= sepbm_{(t)} + twlpnd + twlwet \end{aligned} \right. \quad (A5)$$

where $rchr_{g_{seep(t)}}$ and $rchr_{g_{seep(t-1)}}$ represent the recharge from the water percolating from the soil bottom to the aquifer on a given day (t) and the day before (mm H₂O.day⁻¹), respectively, $seep$ is percolating water and seepage from impoundments on a given day (mm H₂O.day⁻¹), $sepbm_{(t)}$ is the water drained from the bottom of the soil profile (mm H₂O.day⁻¹), $twlpnd$ is the seepage from the ponds (mm H₂O.day⁻¹), $twlwet$ is the seepage
from the wetlands (mm H₂O.day⁻¹), and gwdelay is the time delay for the water percolating from the soil bottom to reach the aquifer (days).

If a sinkhole is simulated as an HRU in karst regions (sink > 0), surface runoff and lateral flow are added to karst seepage, and karst recharge via direct conduits is calculated as follows:

$$ \left\{ \begin{aligned} rchrg_{karst(t)} &= \left(1 - e^{-\frac{1}{gw_{delay}/10}}\right)seep_{direct} + \left(e^{-\frac{1}{gw_{delay}/10}}\right)rchrg_{karst(t-1)} \\ seep_{direct} &= sink \times (surfq + latq) \end{aligned} \right. \quad (A6) $$

where rchrgkarst(t) and rchrgkarst(t-1) represent the recharge from sinkholes and losing streams via direct conduits to the aquifer on a given day and the day before (mm H₂O.day⁻¹), respectively, seepdirect is the seepage from the sinking streams, ponds, and sinkholes on a given day from the HRU (mm H₂O.day⁻¹), sink is the sinkhole partitioning coefficient (0–1), and surfq + latq is the sum of surface runoff and lateral flow from the HRU (mm H₂O.day⁻¹).

Water routed to the deep (confined) aquifer is assumed to be lost from the system and is calculated as follows:

$$ deepst = gw_{seep} + (1 - sink) \times (surfq + latq) \quad (A7) $$

where deepst is the depth of water in the deep aquifer for the day from the HRU (mm H₂O.day⁻¹), gwseep is the water recharging the deep aquifer from the HRU (mm H₂O.day⁻¹), sink is the sinkhole partitioning coefficient (0–1), and surfq + latq is the sum of surface runoff and lateral flow from the HRU (mm H₂O.day⁻¹).

## Appendix E. Palanisamy and Workman [155]

The capacity of a sinkhole is computed using the head-discharge relationship shown in Equation (A8):

$$ Q = C_d A \sqrt{2gH} \quad (A8) $$

where Q is the capacity of the sinkhole (L³/T), Cd is a coefficient of discharge, A is the area of the orifice (L²), g is the acceleration due to gravity (L/T²), and H is the head (water level) at the orifice mouth, set as the depth of water leaving the reach segment upstream of the sinkhole (L).

To compute the diversion of flow from a stream channel located upstream of a sinkhole, the daily flow of water in the reach is input into the sinkhole, and water leaving the reach is set to zero if the incoming simulated streamflow is lower than the capacity of the sinkhole. Conversely, the sinkhole inflow is set to its full capacity, and the remaining flow is transferred to the next reach if the streamflow is higher than the capacity of the sinkhole.

## Appendix F. Zhou et al. [156]

The daily recharge to the aquifer is computed using a modified version of the recharge function of SWAT, shown in Equation (A9):

$$ rchrg_{(j)} = \left(1 - e^{-\frac{1}{delay}}\right) (sepbtm_{(j)} + gwq_{ru(j)}) + \left(1 - e^{-\frac{1}{delay}}\right) \times rchrg_{karst} + \left(e^{-\frac{1}{delay}}\right) rchrg_{(j-1)} \quad (A9) $$

where rchrg(j) is the recharge to the shallow aquifer on a given day (j) (mm H₂O.day⁻¹), delay is the groundwater delay time required for the water to infiltrate from the soil bottom to the aquifer (days), sepbtm(j) is the daily percolation from the bottom of the soil profile (mm H₂O.day⁻¹), gwqru(j) is the daily groundwater contribution to streamflow (mm H₂O.day⁻¹), rchrgkarst is the amount of water seeping through the ponds (mm H₂O.day⁻¹), and rchrg(j-1) is the recharge from the previous day (mm H₂O.day⁻¹).
from the wetlands (mm H₂O.day⁻¹), and gwdelay is the time delay for the water percolating from the soil bottom to reach the aquifer (days).

If a sinkhole is simulated as an HRU in karst regions (sink > 0), surface runoff and lateral flow are added to karst seepage, and karst recharge via direct conduits is calculated as follows:

$$ \left\{ \begin{aligned} rchrg_{karst(t)} &= \left(1 - e^{-\frac{1}{gw_{delay}/10}}\right)seep_{direct} + \left(e^{-\frac{1}{gw_{delay}/10}}\right)rchrg_{karst(t-1)} \\ seep_{direct} &= sink \times (surfq + latq) \end{aligned} \right. \quad (A6) $$

where rchrgkarst(t) and rchrgkarst(t-1) represent the recharge from sinkholes and losing streams via direct conduits to the aquifer on a given day and the day before (mm H₂O.day⁻¹), respectively, seepdirect is the seepage from the sinking streams, ponds, and sinkholes on a given day from the HRU (mm H₂O.day⁻¹), sink is the sinkhole partitioning coefficient (0–1), and surfq + latq is the sum of surface runoff and lateral flow from the HRU (mm H₂O.day⁻¹).

Water routed to the deep (confined) aquifer is assumed to be lost from the system and is calculated as follows:

$$ deepst = gw_{seep} + (1 - sink) \times (surfq + latq) \quad (A7) $$

where deepst is the depth of water in the deep aquifer for the day from the HRU (mm H₂O.day⁻¹), gwseep is the water recharging the deep aquifer from the HRU (mm H₂O.day⁻¹), sink is the sinkhole partitioning coefficient (0–1), and surfq + latq is the sum of surface runoff and lateral flow from the HRU (mm H₂O.day⁻¹).

## Appendix E. Palanisamy and Workman [155]

The capacity of a sinkhole is computed using the head-discharge relationship shown in Equation (A8):

$$ Q = C_d A \sqrt{2gH} \quad (A8) $$

where Q is the capacity of the sinkhole (L³/T), Cd is a coefficient of discharge, A is the area of the orifice (L²), g is the acceleration due to gravity (L/T²), and H is the head (water level) at the orifice mouth, set as the depth of water leaving the reach segment upstream of the sinkhole (L).

To compute the diversion of flow from a stream channel located upstream of a sinkhole, the daily flow of water in the reach is input into the sinkhole, and water leaving the reach is set to zero if the incoming simulated streamflow is lower than the capacity of the sinkhole. Conversely, the sinkhole inflow is set to its full capacity, and the remaining flow is transferred to the next reach if the streamflow is higher than the capacity of the sinkhole.

## Appendix F. Zhou et al. [156]

The daily recharge to the aquifer is computed using a modified version of the recharge function of SWAT, shown in Equation (A9):

$$ rchrg_{(j)} = \left(1 - e^{-\frac{1}{delay}}\right) (sepbtm_{(j)} + gwq_{ru(j)}) + \left(1 - e^{-\frac{1}{delay}}\right) \times rchrg_{karst} + \left(e^{-\frac{1}{delay}}\right) rchrg_{(j-1)} \quad (A9) $$

where rchrg(j) is the recharge to the shallow aquifer on a given day (j) (mm H₂O.day⁻¹), delay is the groundwater delay time required for the water to infiltrate from the soil bottom to the aquifer (days), sepbtm(j) is the daily percolation from the bottom of the soil profile (mm H₂O.day⁻¹), gwqru(j) is the daily groundwater contribution to streamflow (mm H₂O.day⁻¹), rchrgkarst is the amount of water seeping through the ponds (mm H₂O.day⁻¹), and rchrg(j-1) is the recharge from the previous day (mm H₂O.day⁻¹).
## Appendix G. Nikolaidis et al. [157]

The daily inflows of water to the upper and lower reservoirs are calculated using the following equations:

$$Q_{in,up} = \alpha_1 \times Q_{in\_deepGW} \quad (A10)$$

$$Q_{in,low} = (1 - \alpha_1) \times Q_{in\_deepGW} + \alpha_2 Q_{up} \quad (A11)$$

where $Q_{in,up}$ and $Q_{in,low}$ represent the water inflows to the upper and lower reservoirs (mm H₂O.day⁻¹), respectively, $Q_{in\_deepGW}$ is the deep groundwater flow from SWAT (mm H₂O.day⁻¹), $\alpha_1$ is the fraction of deep groundwater flow entering the upper reservoir, $\alpha_2$ is the fraction of flow from the upper reservoir to the lower reservoir, and $Q_{up}$ is the outflow of the upper reservoir (mm H₂O.day⁻¹).

The outflows from the reservoirs and the total spring discharge are calculated as follows:

$$Q_{up} = Q_{up-1}e^{-k_u t} + (\alpha_1 Q_{in\_deepGW}) (1 - e^{-k_u t}) \quad (A12)$$

$$Q_{low} = Q_{low-1}e^{-k_l t} + [(1 - \alpha_1)Q_{in\_deepGW} + \alpha_2 Q_{up}] (1 - e^{-k_l t}) \quad (A13)$$

$$Q_T = (1 - \alpha_2) Q_{up} + Q_{low} \quad (A14)$$

where $Q_T$ is the total spring discharge (mm H₂O.day⁻¹), $Q_{up}$ and $Q_{low}$ represent the outflows of the upper and lower reservoirs (mm H₂O.day⁻¹), respectively, $Q_{up-1}$ and $Q_{low-1}$ are the values of $Q_{up}$ and $Q_{low}$ at the previous time step (mm H₂O.day⁻¹), respectively, $k_u$ and $k_l$ are the recession coefficients of the upper and lower reservoirs (day⁻¹), respectively, $Q_{in\_deepGW}$ is the deep groundwater flow from SWAT (mm H₂O.day⁻¹), $\alpha_1$ is the fraction of deep groundwater flow entering the upper reservoir, $\alpha_2$ is the fraction of flow from the upper reservoir to the lower reservoir, and $t$ is the time step (1 day).

## Appendix H. Nguyen et al. [166]

The diffuse recharge (mm H₂O.day⁻¹) from the bottom of the soil profile to the matrix storage reservoir on day i is calculated as follows:

$$W_{rd,i} = \left[ \left( 1 - e^{-\frac{1}{\delta_{gw}}} \right) W_{seep,i} \right] \times \beta + \left( e^{-\frac{1}{\delta_{gw}}} \right) W_{rd,i-1} \quad (A15)$$

where $W_{rd,i}$ and $W_{rd,i-1}$ (mm H₂O.day⁻¹) represent the amount of diffuse recharge to the matrix reservoir on day i and i − 1, respectively, $\delta_{gw}$ is the delay time for the infiltrated water to reach the matrix storage reservoir (days), $\beta$ is a recharge separation factor (0–1), $W_{seep,i}$ is the total amount of water exiting the bottom of the soil profile on day i (mm H₂O.day⁻¹).

The total volume of diffuse recharge $Q_{in\_matrix,i}$ (m³ H₂O.day⁻¹) to the matrix reservoir on day i is computed as follows:

$$Q_{in\_matrix,i} = \sum_{j=1}^{nhrus} W_{rd,i,j} \times a_j \times 10^{-3} \quad (A16)$$

where $W_{rd,i,j}$ (mm H₂O.day⁻¹) and $a_j$ (m²) are the diffuse recharge and area of HRU number j, respectively, $10^{-3}$ is a unit conversion factor from mm H₂O to m H₂O, and $nhrus$ is the number of HRUs in the recharge area.

The outflow from the matrix storage reservoir is simulated using the linear storage-discharge relationship, as follows:

$$Q_{matrix,i} = Q_{matrix,i-1} (e^{-\alpha_{matrix}\Delta t}) + (1 - e^{-\alpha_{matrix}\Delta t}) \sum_{j=1}^{nhrus} W_{rd,i,j} \times a_j \times 10^{-3} \quad (A17)$$
## Appendix G. Nikolaidis et al. [157]

The daily inflows of water to the upper and lower reservoirs are calculated using the following equations:

$$Q_{in,up} = \alpha_1 \times Q_{in\_deepGW} \quad (A10)$$

$$Q_{in,low} = (1 - \alpha_1) \times Q_{in\_deepGW} + \alpha_2 Q_{up} \quad (A11)$$

where $Q_{in,up}$ and $Q_{in,low}$ represent the water inflows to the upper and lower reservoirs (mm H₂O.day⁻¹), respectively, $Q_{in\_deepGW}$ is the deep groundwater flow from SWAT (mm H₂O.day⁻¹), $\alpha_1$ is the fraction of deep groundwater flow entering the upper reservoir, $\alpha_2$ is the fraction of flow from the upper reservoir to the lower reservoir, and $Q_{up}$ is the outflow of the upper reservoir (mm H₂O.day⁻¹).

The outflows from the reservoirs and the total spring discharge are calculated as follows:

$$Q_{up} = Q_{up-1}e^{-k_u t} + (\alpha_1 Q_{in\_deepGW}) (1 - e^{-k_u t}) \quad (A12)$$

$$Q_{low} = Q_{low-1}e^{-k_l t} + [(1 - \alpha_1)Q_{in\_deepGW} + \alpha_2 Q_{up}] (1 - e^{-k_l t}) \quad (A13)$$

$$Q_T = (1 - \alpha_2) Q_{up} + Q_{low} \quad (A14)$$

where $Q_T$ is the total spring discharge (mm H₂O.day⁻¹), $Q_{up}$ and $Q_{low}$ represent the outflows of the upper and lower reservoirs (mm H₂O.day⁻¹), respectively, $Q_{up-1}$ and $Q_{low-1}$ are the values of $Q_{up}$ and $Q_{low}$ at the previous time step (mm H₂O.day⁻¹), respectively, $k_u$ and $k_l$ are the recession coefficients of the upper and lower reservoirs (day⁻¹), respectively, $Q_{in\_deepGW}$ is the deep groundwater flow from SWAT (mm H₂O.day⁻¹), $\alpha_1$ is the fraction of deep groundwater flow entering the upper reservoir, $\alpha_2$ is the fraction of flow from the upper reservoir to the lower reservoir, and $t$ is the time step (day).

## Appendix H. Nguyen et al. [166]

The diffuse recharge (mm H₂O.day⁻¹) from the bottom of the soil profile to the matrix storage reservoir on day $i$ is calculated as follows:

$$W_{rd,i} = \left[ \left( 1 - e^{-\frac{1}{\delta_{gw}}} \right) W_{seep,i} \right] \times \beta + \left( e^{-\frac{1}{\delta_{gw}}} \right) W_{rd,i-1} \quad (A15)$$

where $W_{rd,i}$ and $W_{rd,i-1}$ (mm H₂O.day⁻¹) represent the amount of diffuse recharge to the matrix reservoir on day $i$ and $i-1$, respectively, $\delta_{gw}$ is the delay time for the infiltrated water to reach the matrix storage reservoir (days), $\beta$ is a recharge separation factor (0–1), $W_{seep,i}$ is the total amount of water exiting the bottom of the soil profile on day $i$ (mm H₂O.day⁻¹).

The total volume of diffuse recharge $Q_{in\_matrix,i}$ (m³ H₂O.day⁻¹) to the matrix reservoir on day $i$ is computed as follows:

$$Q_{in\_matrix,i} = \sum_{j=1}^{nhrus} W_{rd,i,j} \times a_j \times 10^{-3} \quad (A16)$$

where $W_{rd,i,j}$ (mm H₂O.day⁻¹) and $a_j$ (m²) are the diffuse recharge and area of HRU number $j$, respectively, $10^{-3}$ is a unit conversion factor from mm H₂O to m H₂O, and $nhrus$ is the number of HRUs in the recharge area.

The outflow from the matrix storage reservoir is simulated using the linear storage-discharge relationship, as follows:

$$Q_{matrix,i} = Q_{matrix,i-1} (e^{-\alpha_{matrix}\Delta t}) + (1 - e^{-\alpha_{matrix}\Delta t}) \sum_{j=1}^{nhrus} W_{rd,i,j} \times a_j \times 10^{-3} \quad (A17)$$
where $Q_{matrix,i}$ and $Q_{matrix,i-1}$ are the outflows from the matrix storage reservoir on day $i$ and $i-1$ (m$^3$ H$_2$O.day$^{-1}$), respectively, $\alpha_{matrix}$ is the recession constant of the matrix storage reservoir (day$^{-1}$), $\Delta t$ is the time step (1 day), $W_{rd,i,j}$ (mm H$_2$O.day$^{-1}$) and $a_j$ (m$^2$) are the diffuse recharge and area of HRU number $j$, respectively, $10^{-3}$ is a unit conversion factor from mm H$_2$O to m H$_2$O, and $nhrus$ is the number of HRUs in the recharge area.

The concentrated recharge $W_{rc,i}$ (mm H$_2$O.day$^{-1}$) from closed depressions, fractures, and sinkholes to the conduit storage reservoir on day $i$ is calculated as follows:

$$W_{rc,i} = (1 - \beta)W_{seep,i} \qquad (A18)$$

where $\beta$ is a recharge separation factor (0-1), and $W_{seep,i}$ is the total amount of water exiting the bottom of the soil profile on day $i$ (mm H$_2$O.day$^{-1}$).

The total volume of concentrated recharge $Q_{in,conduit}$ (m$^3$ H$_2$O.day$^{-1}$) to the conduit reservoir on day $i$ is calculated as follows:

$$Q_{in,conduit} = \sum_{j=1}^{nhrus} W_{rc,i,j} \times a_j \times 10^{-3} + rttlc_i + Q_{matrix,i} \qquad (A19)$$

where $Q_{matrix,i}$ is the outflow from the matrix storage reservoir on day $i$, $rttlc_i$ is the recharge from losing streams on day $i$ (m$^3$ H$_2$O.day$^{-1}$), $\Delta t$ is the time step (1 day), $W_{rc,i,j}$ (mm H$_2$O.day$^{-1}$) and $a_j$ (m$^2$) are the concentrated recharge and area of HRU number $j$, respectively, $10^{-3}$ is a unit conversion factor from mm H$_2$O to m H$_2$O, and $nhrus$ is the number of HRUs in the recharge area.

The outflow from the conduit storage reservoir is simulated using the linear storage-discharge relationship, as follows:

$$Q_{conduit,i} = Q_{conduit,i-1}(e^{-\alpha_{conduit}\Delta t}) + \left[ \sum_{j=1}^{nhrus} W_{rc,i,j} \times a_j \times 10^{-3} + rttlc_i + Q_{matrix,i} \right] (1 - e^{-\alpha_{conduit}\Delta t}) \qquad (A20)$$

where $Q_{conduit,i}$ and $Q_{conduit,i-1}$ are the outflows from the conduit storage reservoir on day $i$ and $i-1$ (m$^3$ H$_2$O.day$^{-1}$), respectively, $\alpha_{conduit}$ is the recession constant of the conduit storage reservoir (day$^{-1}$), $rttlc_i$ is the amount of recharge from losing streams on day $i$ (m$^3$ H$_2$O.day$^{-1}$), $\Delta t$ is the time step (1 day), $W_{rc,i,j}$ (mm H$_2$O.day$^{-1}$) and $a_j$ (m$^2$) are the concentrated recharge and area of HRU number $j$, respectively, $10^{-3}$ is a unit conversion factor from mm H$_2$O to m H$_2$O, and $nhrus$ is the number of HRUs in the recharge area.

The total runoff $Q_{river,i}$ (m$^3$ H$_2$O.day$^{-1}$) of the basin at the location of the spring is calculated as follows:

$$Q_{river,i} = Q_{conduit,i} + Q_{direct,i} \qquad (A21)$$

where $Q_{direct,i}$ is the daily direct runoff calculated as the sum of the surface runoff and lateral flow from the basin where the spring is located (m$^3$ H$_2$O.day$^{-1}$), and $Q_{conduit,i}$ is the outflow from the conduit storage reservoir on day $i$ (m$^3$ H$_2$O.day$^{-1}$).

## Appendix I. Wang et al. [167]

According to the model proposed by Wang et al. [167], water added to the epikarst is first computed as a function of the water percolated from the bottom of the soil profile minus the water recharge from the depression and aven.

The attenuation function of the epikarst is represented as follows:

$$TT = \frac{sat - f_c}{k} \qquad (A22)$$

where $TT$ is the attenuation coefficient, $sat$ is the saturation moisture content in the epikarst system (mm H$_2$O), $f_c$ is the water holding capacity in the epikarst system (mm H$_2$O), and $k$ is the saturated hydraulic conductivity (mm.h$^{-1}$).
where $Q_{matrix,i}$ and $Q_{matrix,i-1}$ are the outflows from the matrix storage reservoir on day $i$ and $i-1$ (m$^3$ H$_2$O.day$^{-1}$), respectively, $\alpha_{matrix}$ is the recession constant of the matrix storage reservoir (day$^{-1}$), $\Delta t$ is the time step (1 day), $W_{rd,i,j}$ (mm H$_2$O.day$^{-1}$) and $a_j$ (m$^2$) are the diffuse recharge and area of HRU number $j$, respectively, $10^{-3}$ is a unit conversion factor from mm H$_2$O to m H$_2$O, and $nhrus$ is the number of HRUs in the recharge area.

The concentrated recharge $W_{rc,i}$ (mm H$_2$O.day$^{-1}$) from closed depressions, fractures, and sinkholes to the conduit storage reservoir on day $i$ is calculated as follows:

$$W_{rc,i} = (1 - \beta)W_{seep,i} \qquad (A18)$$

where $\beta$ is a recharge separation factor (0-1), and $W_{seep,i}$ is the total amount of water exiting the bottom of the soil profile on day $i$ (mm H$_2$O.day$^{-1}$).

The total volume of concentrated recharge $Q_{in,conduit}$ (m$^3$ H$_2$O.day$^{-1}$) to the conduit reservoir on day $i$ is calculated as follows:

$$Q_{in,conduit} = \sum_{j=1}^{nhrus} W_{rc,i,j} \times a_j \times 10^{-3} + rttlc_i + Q_{matrix,i} \qquad (A19)$$

where $Q_{matrix,i}$ is the outflow from the matrix storage reservoir on day $i$, $rttlc_i$ is the recharge from losing streams on day $i$ (m$^3$ H$_2$O.day$^{-1}$), $\Delta t$ is the time step (1 day), $W_{rc,i,j}$ (mm H$_2$O.day$^{-1}$) and $a_j$ (m$^2$) are the concentrated recharge and area of HRU number $j$, respectively, $10^{-3}$ is a unit conversion factor from mm H$_2$O to m H$_2$O, and $nhrus$ is the number of HRUs in the recharge area.

The outflow from the conduit storage reservoir is simulated using the linear storage-discharge relationship, as follows:

$$Q_{conduit,i} = Q_{conduit,i-1} (e^{-\alpha_{conduit}\Delta t}) + \left[ \sum_{j=1}^{nhrus} W_{rc,i,j} \times a_j \times 10^{-3} + rttlc_i + Q_{matrix,i} \right] (1 - e^{-\alpha_{conduit}\Delta t}) \qquad (A20)$$

where $Q_{conduit,i}$ and $Q_{conduit,i-1}$ are the outflows from the conduit storage reservoir on day $i$ and $i-1$ (m$^3$ H$_2$O.day$^{-1}$), respectively, $\alpha_{conduit}$ is the recession constant of the conduit storage reservoir (day$^{-1}$), $rttlc_i$ is the amount of recharge from losing streams on day $i$ (m$^3$ H$_2$O.day$^{-1}$), $\Delta t$ is the time step (1 day), $W_{rc,i,j}$ (mm H$_2$O.day$^{-1}$) and $a_j$ (m$^2$) are the concentrated recharge and area of HRU number $j$, respectively, $10^{-3}$ is a unit conversion factor from mm H$_2$O to m H$_2$O, and $nhrus$ is the number of HRUs in the recharge area.

The total runoff $Q_{river,i}$ (m$^3$ H$_2$O.day$^{-1}$) of the basin at the location of the spring is calculated as follows:

$$Q_{river,i} = Q_{conduit,i} + Q_{direct,i} \qquad (A21)$$

where $Q_{direct,i}$ is the daily direct runoff calculated as the sum of the surface runoff and lateral flow from the basin where the spring is located (m$^3$ H$_2$O.day$^{-1}$), and $Q_{conduit,i}$ is the outflow from the conduit storage reservoir on day $i$ (m$^3$ H$_2$O.day$^{-1}$).

## Appendix I. Wang et al. [167]

According to the model proposed by Wang et al. [167], water added to the epikarst is first computed as a function of the water percolated from the bottom of the soil profile minus the water recharge from the depression and aven.

The attenuation function of the epikarst is represented as follows:

$$TT = \frac{sat - f_c}{k} \qquad (A22)$$

where $TT$ is the attenuation coefficient, $sat$ is the saturation moisture content in the epikarst system (mm H$_2$O), $f_c$ is the water holding capacity in the epikarst system (mm H$_2$O), and $k$ is the saturated hydraulic conductivity (mm.h$^{-1}$).
Infiltration water through the epikarst is calculated as follows:

$$Q_i = W_{st} \left(1 - e^{-t/T}\right) \qquad (A23)$$

where $Q_i$ is the water infiltrated through the epikarst on a given day (mm H₂O.day⁻¹), $W_{st}$ is the water content in the epikarst which varies in function of the daily percolation from the soil bottom (mm H₂O.day⁻¹), and $t$ is the simulation time step (1 day), and $TT$ is the attenuation coefficient.

The daily recharge values to the upper and lower reservoirs are computed as follows:

$$Q_{in,up} = \alpha_1 Q_i + Q_C \qquad (A24)$$

$$Q_{in,low} = (1 - \alpha_1) Q_i + \alpha_2 Q_{up} \qquad (A25)$$

where $Q_i$ is the water infiltrated through the epikarst on a given day (mm H₂O.day⁻¹), $Q_C$ is the injection volume from the depression and aven (mm H₂O.day⁻¹).

The discharge values of the upper and lower reservoirs are calculated using Equations (A26) and (A27), respectively, as follows:

$$Q_{up} = Q_{up-1} \left(e^{-k_u t}\right) + \left(\alpha_1 Q_i + Q_C\right) \left(1 - e^{-k_u t}\right) \qquad (A26)$$

$$Q_{low} = Q_{low-1} \left(e^{-k_l t}\right) + \left[(1 - \alpha_1) Q_i + \alpha_2 Q_{up}\right] \left(1 - e^{-k_l t}\right) \qquad (A27)$$

where $Q_{up}$ and $Q_{low}$ are the daily discharge values of the upper and lower reservoirs (mm H₂O.day⁻¹), $Q_{up-1}$ and $Q_{low-1}$ are the discharge values of the upper and lower reservoirs on the previous day (mm H₂O.day⁻¹), $Q_i$ is the water infiltrated through the epikarst on a given day (mm H₂O.day⁻¹), $Q_C$ is the recharge from the depression and aven (mm H₂O.day⁻¹), $\alpha_1$ is the coefficient of proportionality for infiltration from the epikarst to the upper reservoir, and $\alpha_2$ is the coefficient of proportionality for infiltration from the upper to the lower reservoir, and $k_u$ and $k_l$ are the recession coefficients of the upper and lower reservoirs (day⁻¹), respectively.

The total contribution of the reservoirs to the streamflow $Q_T$ (mm H₂O.day⁻¹) is then calculated as the sum of outflows from both reservoirs:

$$Q_T = (1 - \alpha_2) Q_{up} + Q_{low} \qquad (A28)$$

where $Q_{up}$ and $Q_{low}$ represent the discharge components from the upper reservoir and lower reservoirs (mm H₂O.day⁻¹), respectively, and $\alpha_2$ is the coefficient of proportionality for infiltration from the upper to the lower reservoir.

## Appendix J. Geng et al. [168]

In the model proposed by Geng et al. [168], the upper reservoir is recharged by water infiltration from the bottom of the soil profile. The fast recharge through the sinkholes and cracks is modeled using pond leakage, as follows:

$$Q_{kr} = (1 - K_d) twl pnd + K_d Q_{kr0} \qquad (A29)$$

where $Q_{kr}$ and $Q_{kr0}$ represent the recharge of karst groundwater on a given day and the day before (mm H₂O.day⁻¹), respectively, $twl pnd$ represents pond leakage (mm H₂O.day⁻¹), and $K_d$ is the flow delay coefficient in the karst groundwater recharge from sinkholes (days).
Infiltration water through the epikarst is calculated as follows:

$$Q_i = W_{st} (1 - e^{-t/T}) \qquad (A23)$$

where $Q_i$ is the water infiltrated through the epikarst on a given day (mm H₂O.day⁻¹), $W_{st}$ is the water content in the epikarst which varies in function of the daily percolation from the soil bottom (mm H₂O.day⁻¹), and $t$ is the simulation time step (1 day), and $TT$ is the attenuation coefficient.

The daily recharge values to the upper and lower reservoirs are computed as follows:

$$Q_{in,up} = \alpha_1 Q_i + Q_C \qquad (A24)$$

$$Q_{in,low} = (1 - \alpha_1) Q_i + \alpha_2 Q_{up} \qquad (A25)$$

where $Q_i$ is the water infiltrated through the epikarst on a given day (mm H₂O.day⁻¹), $Q_C$ is the injection volume from the depression and aven (mm H₂O.day⁻¹).

The discharge values of the upper and lower reservoirs are calculated using Equations (A26) and (A27), respectively, as follows:

$$Q_{up} = Q_{up-1} (e^{-k_u t}) + (\alpha_1 Q_i + Q_C) \times (1 - e^{-k_u t}) \qquad (A26)$$

$$Q_{low} = Q_{low-1} (e^{-k_l t}) + [(1 - \alpha_1) Q_i + \alpha_2 Q_{up}] \times (1 - e^{-k_l t}) \qquad (A27)$$

where $Q_{up}$ and $Q_{low}$ are the daily discharge values of the upper and lower reservoirs (mm H₂O.day⁻¹), $Q_{up-1}$ and $Q_{low-1}$ are the discharge values of the upper and lower reservoirs on the previous day (mm H₂O.day⁻¹), $Q_i$ is the water infiltrated through the epikarst on a given day (mm H₂O.day⁻¹), $Q_C$ is the recharge from the depression and aven (mm H₂O.day⁻¹), $\alpha_1$ is the coefficient of proportionality for infiltration from the epikarst to the upper reservoir, and $\alpha_2$ is the coefficient of proportionality for infiltration from the upper to the lower reservoir, and $k_u$ and $k_l$ are the recession coefficients of the upper and lower reservoirs (day⁻¹), respectively.

The total contribution of the reservoirs to the streamflow $Q_T$ (mm H₂O.day⁻¹) is then calculated as the sum of outflows from both reservoirs:

$$Q_T = (1 - \alpha_2) Q_{up} + Q_{low} \qquad (A28)$$

where $Q_{up}$ and $Q_{low}$ represent the discharge components from the upper reservoir and lower reservoirs (mm H₂O.day⁻¹), respectively, and $\alpha_2$ is the coefficient of proportionality for infiltration from the upper to the lower reservoir.

## Appendix J. Geng et al. [168]

In the model proposed by Geng et al. [168], the upper reservoir is recharged by water infiltration from the bottom of the soil profile. The fast recharge through the sinkholes and cracks is modeled using pond leakage, as follows:

$$Q_{kr} = (1 - K_d) twl pnd + K_d Q_{kr0} \qquad (A29)$$

where $Q_{kr}$ and $Q_{kr0}$ represent the recharge of karst groundwater on a given day and the day before (mm H₂O.day⁻¹), respectively, $twl pnd$ represents pond leakage (mm H₂O.day⁻¹), and $K_d$ is the flow delay coefficient in the karst groundwater recharge from sinkholes (days).
The discharge for each reservoir is calculated as follows:

$$
\left\{
\begin{array}{l}
Q_i = Q_{i,0} (e^{-\alpha_i \Delta t}) + Q_{in,i} (1 - e^{-\alpha_i \Delta t}) \\
Q_{in,up} = Q_{seep} \\
Q_{in, mid} = Q_{kr} + \beta_1 Q_{up} \\
Q_{in,low} = \beta_2 Q_{up} + \beta_3 Q_{mid} \\
i = up, mid, low
\end{array}
\right.
\qquad (A30)
$$

where $Q_i$ is the discharge of reservoir $i$ on a given day (mm H₂O.day⁻¹), $i$ denotes the reservoir (one of the upper, middle, and lower reservoirs), $Q_{in,i}$ is the daily recharge to reservoir $i$ (mm H₂O.day⁻¹), $\alpha_i$ is the recession constant of reservoir $i$ (day⁻¹), $\Delta t$ is the time step (1 day), $Q_{i,0}$ is the discharge of reservoir $i$ on the previous day (mm H₂O.day⁻¹), $Q_{kr}$ is the karst groundwater recharge on a given day (mm H₂O.day⁻¹), $Q_{seep}$ the daily infiltration recharge from the soil bottom (mm H₂O.day⁻¹), and $\beta_1$, $\beta_2$, and $\beta_3$ are coefficients of proportionality.

The total discharge of karst groundwater $Q_T$ (mm H₂O.day⁻¹) is then calculated as the sum of outflows from the three reservoirs, as follows:

$$Q_T = (1 - \beta_1 - \beta_2)Q_{up} + (1 - \beta_3)Q_{mid} + Q_{low} \quad (A31)$$

where $Q_{up}$, $Q_{mid}$ and $Q_{low}$ represent the outflows from the upper, middle, and lower reservoirs (mm H₂O.day⁻¹), respectively, and $\beta_1$, $\beta_2$, and $\beta_3$ are coefficients of proportionality.

## Appendix K. Wang and Brubaker [169]

The return flow is estimated as an explicit function of the shallow aquifer storage based on the non-linear algorithm of Wittenberg [180], as follows:

$$S_{sh} - S_{sh,min} = \frac{1}{\alpha_{gw}} Q_{gw}^{\beta_{gw}} \qquad (A32)$$

$$Q_{gw} = (\alpha_{gw} (S_{sh} - S_{sh,min}))^{1/\beta_{gw}} \qquad (A33)$$

where $S_{sh}$ is the shallow aquifer storage (L), $S_{sh,min}$ is the minimum storage for groundwater flow to occur (L), $\alpha_{gw}$ is a scale parameter ($T^\beta L^{3(1-\beta)}$), $\beta_{gw}$ is a dimensionless coefficient. If $\beta_{gw} = 1$, the non-linear model becomes linear.

## Appendix L. Eini et al. [36]

### Appendix L.1. SWAT-ML

In SWAT-ML, aquifer recharge from the non-karst areas and recharge from the karst areas are computed separately using Equations (A34) and (A35), respectively, each with a different delay time parameter.

$$rchr g_{(j)} = (1 - g w_{delay(j)}) \times (sepbtm_{(j)} + g w q_{ru(j)}) + g w_{delay(j)} \times rchr g_{1} \quad (A34)$$

$$rchr g_{karst(j)} = \left(1 - e^{\frac{-1}{\frac{8 w_{delay(j)}}{X}}}\right) (sepbtm_{(j)} + g w q_{ru(j)}) + \left(e^{\frac{-1}{\frac{8 w_{delay(j)}}{X}}}\right) rchr g_{karst1} \quad (A35)$$

where $rchr g_{(j)}$ and $rchr g_{karst(j)}$ represent the daily recharge values from non-karst and karst HRUs (mm H₂O.day⁻¹), respectively, $j$ denotes the HRU number, $gw_{delay(j)}$ is the recharge delay time (days), $sepbtm_{(j)}$ is the daily percolation from the bottom of the soil profile (mm H₂O.day⁻¹), $gwq_{ru(j)}$ represents the seepage from the lakes, wetlands, and riverside branches (mm H₂O.day⁻¹), $rchr g_1$ and $rchr g_{karst1}$ represent the recharge values from the previous day (mm H₂O.day⁻¹), and $X$ is a non-dimensional calibration parameter (1, +∞) used to adjust infiltration rates in karstic HRUs.
The discharge for each reservoir is calculated as follows:

$$
\left\{
\begin{array}{l}
Q_i = Q_{i,0} (e^{-\alpha_i \Delta t}) + Q_{in,i} (1 - e^{-\alpha_i \Delta t}) \\
Q_{in,up} = Q_{seep} \\
Q_{in, mid} = Q_{kr} + \beta_1 Q_{up} \\
Q_{in,low} = \beta_2 Q_{up} + \beta_3 Q_{mid} \\
i = up, mid, low
\end{array}
\right.
\qquad (A30)
$$

where $Q_i$ is the discharge of reservoir $i$ on a given day (mm H₂O.day⁻¹), $i$ denotes the reservoir (one of the upper, middle, and lower reservoirs), $Q_{in,i}$ is the daily recharge to reservoir $i$ (mm H₂O.day⁻¹), $\alpha_i$ is the recession constant of reservoir $i$ (day⁻¹), $\Delta t$ is the time step (1 day), $Q_{i,0}$ is the discharge of reservoir $i$ on the previous day (mm H₂O.day⁻¹), $Q_{kr}$ is the karst groundwater recharge on a given day (mm H₂O.day⁻¹), $Q_{seep}$ the daily infiltration recharge from the soil bottom (mm H₂O.day⁻¹), and $\beta_1, \beta_2, \beta_3$ are coefficients of proportionality.

The total discharge of karst groundwater $Q_T$ (mm H₂O.day⁻¹) is then calculated as the sum of outflows from the three reservoirs, as follows:

$$Q_T = (1 - \beta_1 - \beta_2)Q_{up} + (1 - \beta_3)Q_{mid} + Q_{low} \quad (A31)$$

where $Q_{up}$, $Q_{mid}$ and $Q_{low}$ represent the outflows from the upper, middle, and lower reservoirs (mm H₂O.day⁻¹), respectively, and $\beta_1, \beta_2$, and $\beta_3$ are coefficients of proportionality.

## Appendix K. Wang and Brubaker [169]

The return flow is estimated as an explicit function of the shallow aquifer storage based on the non-linear algorithm of Wittenberg [180], as follows:

$$S_{sh} - S_{sh,min} = \frac{1}{\alpha_{gw}} Q_{gw}^{\beta_{gw}} \qquad (A32)$$

$$Q_{gw} = (\alpha_{gw} (S_{sh} - S_{sh,min}))^{1/\beta_{gw}} \qquad (A33)$$

where $S_{sh}$ is the shallow aquifer storage (L), $S_{sh,min}$ is the minimum storage for groundwater flow to occur (L), $\alpha_{gw}$ is a scale parameter ($T^\beta L^{3(1-\beta)}$), $\beta_{gw}$ is a dimensionless coefficient. If $\beta_{gw} = 1$, the non-linear model becomes linear.

## Appendix L. Eini et al. [36]

### Appendix L.1. SWAT-ML

In SWAT-ML, aquifer recharge from the non-karst areas and recharge from the karst areas are computed separately using Equations (A34) and (A35), respectively, each with a different delay time parameter.

$$rchr g_{(j)} = (1 - g w_{delay(j)}) \times (sepbtm_{(j)} + g w q_{ru(j)}) + g w_{delay(j)} \times rchr g_{1} \quad (A34)$$

$$rchr g_{karst(j)} = \left(1 - e^{\frac{-1}{\frac{8 w_{delay(j)}}{X}}}\right) (sepbtm_{(j)} + g w q_{ru(j)}) + \left(e^{\frac{-1}{\frac{8 w_{delay(j)}}{X}}}\right) rchr g_{karst1} \quad (A35)$$

where $rchr g_{(j)}$ and $rchr g_{karst(j)}$ represent the daily recharge values from non-karst and karst HRUs (mm H₂O.day⁻¹), respectively, $j$ denotes the HRU number, $gw_{delay(j)}$ is the recharge delay time (days), $sepbtm_{(j)}$ is the daily percolation from the bottom of the soil profile (mm H₂O.day⁻¹), $gwq_{ru(j)}$ represents the seepage from the lakes, wetlands, and riverside branches (mm H₂O.day⁻¹), $rchr g_1$ and $rchr g_{karst1}$ represent the recharge values from the previous day (mm H₂O.day⁻¹), and $X$ is a non-dimensional calibration parameter (1, +∞) used to adjust infiltration rates in karstic HRUs.
## Appendix L.2. SWAT-CF

In the standard SWAT, crack volume is set to zero for wet conditions. Conversely, the crack volume in SWAT-CF is calculated using Equation (A36), which allows cracks to form in wet soils, thus removing the distinction of dry and wet soil conditions during the crack flow process.

$$ volcrk_{(1,j)} = crlag \times volcrk_{(1,j)} + (1 - crlag) \times volcrknew \quad (A36) $$

where volcrk is the crack volume for soil layer (mm), crlag is a daily lag factor for crack development, volcrknew is the crack volume for soil layer based on new moisture conditions (mm), j is the HRU number, and l is the counter.

## References

1. Chen, Z.; Auler, A.S.; Bakalowicz, M.; Drew, D.; Griger, F.; Hartmann, J.; Jiang, G.; Moosdorf, N.; Richts, A.; Stevanovic, Z.; et al. The World Karst Aquifer Mapping project: Concept, mapping procedure and map of Europe. *Hydrogeol. J.* **2017**, *25*, 771–785. [CrossRef]

2. Biondić, R.; Meaški, H.; Biondić, B.; Loborec, J. Karst Aquifer Vulnerability Assessment (KAVA) Method—A Novel GIS-Based Method for Deep Karst Aquifers. *Sustainability* **2021**, *13*, 3325. [CrossRef]

3. Auler, A.S.; Stevanović, Z. Preface: Five decades of advances in karst hydrogeology. *Hydrogeol. J.* **2021**, *29*, 1–6. [CrossRef]

4. Goldscheider, N.; Chen, Z.; Auler, A.S.; Bakalowicz, M.; Broda, S.; Drew, D.; Hartmann, J.; Jiang, G.; Moosdorf, N.; Stevanovic, Z.; et al. Global distribution of carbonate rocks and karst water resources. *Hydrogeol. J.* **2020**, *28*, 1661–1677. [CrossRef]

5. Zeiger, S.J.; Owen, M.R.; Pavlowsky, R.T. Simulating nonpoint source pollutant loading in a karst basin: A SWAT modeling application. *Sci. Total Environ.* **2021**, *785*, 147295. [CrossRef]

6. Dal Soglio, L.; Danquigny, C.; Mazzilli, N.; Emblanch, C.; Massonnat, G. Taking into Account both Explicit Conduits and the Unsaturated Zone in Karst Reservoir Hybrid Models: Impact on the Outlet Hydrograph. *Water* **2020**, *12*, 3221. [CrossRef]

7. Wang, Z.; Yin, J.-J.; Pu, J.; Wang, P.; Liang, X.; Yang, P.; He, Q.; Gou, P.; Yuan, D. Integrated understanding of the Critical Zone processes in a subtropical karst watershed (Qingmuguan, Southwestern China): Hydrochemical and isotopic constraints. *Sci. Total Environ.* **2020**, *749*, 141257. [CrossRef] [PubMed]

8. Chen, X.; Zhang, Z.; Soulsby, C.; Cheng, Q.; Binley, A.; Jiang, R.; Tao, M. Characterizing the heterogeneity of karst critical zone and its hydrological function: An integrated approach. *Hydrol. Process.* **2018**, *32*, 2932–2946. [CrossRef]

9. Kovačič, G.; Petrič, M.; Ravbar, N. Evaluation and Quantification of the Effects of Climate and Vegetation Cover Change on Karst Water Sources: Case Studies of Two Springs in South-Western Slovenia. *Water* **2020**, *12*, 3087. [CrossRef]

10. Stroj, A.; Briški, M.; Oštrić, M. Study of Groundwater Flow Properties in a Karst System by Coupled Analysis of Diverse Environmental Tracers and Discharge Dynamics. *Water* **2020**, *12*, 2442. [CrossRef]

11. Bauer, S.; Liedl, R.; Sauter, M. Modeling the influence of epikarst evolution on karst aquifer genesis: A time-variant recharge boundary condition for joint karst-epikarst development. *Water Resour. Res.* **2005**, *41*, W09416. [CrossRef]

12. Paiva, I.; Cunha, L. Characterization of the hydrodynamic functioning of the Degracias-Sicó Karst Aquifer, Portugal. *Hydrogeol. J.* **2020**, *28*, 2613–2629. [CrossRef]

13. Hartmann, A.; Liu, Y.; Olarinoye, T.; Marx, V. Integrating field work and large-scale modeling to improve assessment of karst water resources. *Hydrogeol. J.* **2021**, *29*, 315–329. [CrossRef]

14. Giese, M.; Reimann, T.; Bailly-Comte, V.; Maréchal, J.-C.; Sauter, M.; Geyer, T. Turbulent and Laminar Flow in Karst Conduits Under Unsteady Flow Conditions: Interpretation of Pumping Tests by Discrete Conduit-Continuum Modeling. *Water Resour. Res.* **2018**, *54*, 1918–1933. [CrossRef]

15. de Rooij, R.; Graham, W. Generation of complex karstic conduit networks with a hydrochemical model. *Water Resour. Res.* **2017**, *53*, 6993–7011. [CrossRef]

16. Zhao, L.-J.; Yang, Y.; Cao, J.-W.; Wang, Z.; Luan, S.; Xia, R.-Y. Applying a modified conduit flow process to understand conduit-matrix exchange of a karst aquifer. *China Geol.* **2022**, *5*, 26–33. [CrossRef]

17. Le Mesnil, M.; Charlier, J.-B.; Moussa, R.; Caballero, Y.; Dörfliger, N. Interbasin groundwater flow: Characterization, role of karst areas, impact on annual water balance and flood processes. *J. Hydrol.* **2020**, *585*, 124583. [CrossRef]

18. Gutiérrez, F.; Parise, M.; De Waele, J.; Jourde, H. A review on natural and human-induced geohazards and impacts in karst. *Earth-Sci. Rev.* **2014**, *138*, 61–88. [CrossRef]

19. Goldscheider, N. A holistic approach to groundwater protection and ecosystem services in karst terrains. *Carbonates Evaporites* **2019**, *34*, 1241–1249. [CrossRef]

20. Iván, V.; Mádl-Szőnyi, J. State of the art of karst vulnerability assessment: Overview, evaluation and outlook. *Environ. Earth. Sci.* **2017**, *76*, 112. [CrossRef]

21. Entezari, M.; Yamani, M.; Jafari Aghdam, M. Evaluation of intrinsic vulnerability, hazard and risk mapping for karst aquifers, Khorein aquifer, Kermanshah province: A case study. *Environ. Earth Sci.* **2016**, *75*, 435. [CrossRef]
## Appendix L.2. SWAT-CF

In the standard SWAT, crack volume is set to zero for wet conditions. Conversely, the crack volume in SWAT-CF is calculated using Equation (A36), which allows cracks to form in wet soils, thus removing the distinction of dry and wet soil conditions during the crack flow process.

$$ volcrk_{(l,j)} = crlag \times volcrk_{(l,j)} + (1 - crlag) \times volcrknew \quad (A36) $$

where volcrk is the crack volume for soil layer (mm), crlag is a daily lag factor for crack development, volcrknew is the crack volume for soil layer based on new moisture conditions (mm), j is the HRU number, and l is the counter.

## References

1. Chen, Z.; Auler, A.S.; Bakalowicz, M.; Drew, D.; Griger, F.; Hartmann, J.; Jiang, G.; Moosdorf, N.; Richts, A.; Stevanovic, Z.; et al. The World Karst Aquifer Mapping project: Concept, mapping procedure and map of Europe. *Hydrogeol. J.* **2017**, 25, 771–785. [CrossRef]

2. Biondić, R.; Meaški, H.; Biondić, B.; Loborec, J. Karst Aquifer Vulnerability Assessment (KAVA) Method—A Novel GIS-Based Method for Deep Karst Aquifers. *Sustainability* **2021**, 13, 3325. [CrossRef]

3. Auler, A.S.; Stevanović, Z. Preface: Five decades of advances in karst hydrogeology. *Hydrogeol. J.* **2021**, 29, 1–6. [CrossRef]

4. Goldscheider, N.; Chen, Z.; Auler, A.S.; Bakalowicz, M.; Broda, S.; Drew, D.; Hartmann, J.; Jiang, G.; Moosdorf, N.; Stevanovic, Z.; et al. Global distribution of carbonate rocks and karst water resources. *Hydrogeol. J.* **2020**, 28, 1661–1677. [CrossRef]

5. Zeiger, S.J.; Owen, M.R.; Pavlowsky, R.T. Simulating nonpoint source pollutant loading in a karst basin: A SWAT modeling application. *Sci. Total Environ.* **2021**, 785, 147295. [CrossRef]

6. Dal Soglio, L.; Danquigny, C.; Mazzilli, N.; Emblanch, C.; Massonnat, G. Taking into Account both Explicit Conduits and the Unsaturated Zone in Karst Reservoir Hybrid Models: Impact on the Outlet Hydrograph. *Water* **2020**, 12, 3221. [CrossRef]

7. Wang, Z.; Yin, J.-J.; Pu, J.; Wang, P.; Liang, X.; Yang, P.; He, Q.; Gou, P.; Yuan, D. Integrated understanding of the Critical Zone processes in a subtropical karst watershed (Qingmuguan, Southwestern China): Hydrochemical and isotopic constraints. *Sci. Total Environ.* **2020**, 749, 141257. [CrossRef] [PubMed]

8. Chen, X.; Zhang, Z.; Soulsby, C.; Cheng, Q.; Binley, A.; Jiang, R.; Tao, M. Characterizing the heterogeneity of karst critical zone and its hydrological function: An integrated approach. *Hydrol. Process.* **2018**, 32, 2932–2946. [CrossRef]

9. Kovačič, G.; Petrič, M.; Ravbar, N. Evaluation and Quantification of the Effects of Climate and Vegetation Cover Change on Karst Water Sources: Case Studies of Two Springs in South-Western Slovenia. *Water* **2020**, 12, 3087. [CrossRef]

10. Stroj, A.; Briški, M.; Oštrić, M. Study of Groundwater Flow Properties in a Karst System by Coupled Analysis of Diverse Environmental Tracers and Discharge Dynamics. *Water* **2020**, 12, 2442. [CrossRef]

11. Bauer, S.; Liedl, R.; Sauter, M. Modeling the influence of epikarst evolution on karst aquifer genesis: A time-variant recharge boundary condition for joint karst-epikarst development. *Water Resour. Res.* **2005**, 41, W09416. [CrossRef]

12. Paiva, I.; Cunha, L. Characterization of the hydrodynamic functioning of the Degracias-Sicó Karst Aquifer, Portugal. *Hydrogeol. J.* **2020**, 28, 2613–2629. [CrossRef]

13. Hartmann, A.; Liu, Y.; Olarinoye, T.; Marx, V. Integrating field work and large-scale modeling to improve assessment of karst water resources. *Hydrogeol. J.* **2021**, 29, 315–329. [CrossRef]

14. Giese, M.; Reimann, T.; Bailly-Comte, V.; Maréchal, J.-C.; Sauter, M.; Geyer, T. Turbulent and Laminar Flow in Karst Conduits Under Unsteady Flow Conditions: Interpretation of Pumping Tests by Discrete Conduit-Continuum Modeling. *Water Resour. Res.* **2018**, 54, 1918–1933. [CrossRef]

15. de Rooij, R.; Graham, W. Generation of complex karstic conduit networks with a hydrochemical model. *Water Resour. Res.* **2017**, 53, 6993–7011. [CrossRef]

16. Zhao, L.-J.; Yang, Y.; Cao, J.-W.; Wang, Z.; Luan, S.; Xia, R.-Y. Applying a modified conduit flow process to understand conduit-matrix exchange of a karst aquifer. *China Geol.* **2022**, 5, 26–33. [CrossRef]

17. Le Mesnil, M.; Charlier, J.-B.; Moussa, R.; Caballero, Y.; Dörfliger, N. Interbasin groundwater flow: Characterization, role of karst areas, impact on annual water balance and flood processes. *J. Hydrol.* **2020**, 585, 124583. [CrossRef]

18. Gutiérrez, F.; Parise, M.; De Waele, J.; Jourde, H. A review on natural and human-induced geohazards and impacts in karst. *Earth-Sci. Rev.* **2014**, 138, 61–88. [CrossRef]

19. Goldscheider, N. A holistic approach to groundwater protection and ecosystem services in karst terrains. *Carbonates Evaporites* **2019**, 34, 1241–1249. [CrossRef]

20. Iván, V.; Mádl-Szőnyi, J. State of the art of karst vulnerability assessment: Overview, evaluation and outlook. *Environ. Earth. Sci.* **2017**, 76, 112. [CrossRef]

21. Entezari, M.; Yamani, M.; Jafari Aghdam, M. Evaluation of intrinsic vulnerability, hazard and risk mapping for karst aquifers, Khorein aquifer, Kermanshah province: A case study. *Environ. Earth Sci.* **2016**, 75, 435. [CrossRef]
22. Lukač Reberski, J.; Rubinić, J.; Terzić, J.; Radišić, M. Climate Change Impacts on Groundwater Resources in the Coastal Karstic Adriatic Area: A Case Study from the Dinaric Karst. *Nat. Resour. Res.* **2020**, *29*, 1975–1988. [CrossRef]

23. Hartmann, A.; Goldscheider, N.; Wagener, T.; Lange, J.; Weiler, M. Karst water resources in a changing world: Review of hydrological modeling approaches. *Rev. Geophys.* **2014**, *52*, 218–242. [CrossRef]

24. Sivelle, V.; Jourde, H.; Bittner, D.; Mazzilli, N.; Tramblay, Y. Assessment of the relative impacts of climate changes and anthropogenic forcing on spring discharge of a Mediterranean karst system. *J. Hydrol.* **2021**, *598*, 126396. [CrossRef]

25. Mudarra, M.; Hartmann, A.; Andreo, B. Combining experimental methods and modeling to quantify the complex recharge behavior of karst aquifers. *Water Resour. Res.* **2019**, *55*, 1384–1404. [CrossRef]

26. Gill, L.W.; Schuler, P.; Duran, L.; Morrissey, P.; Johnston, P.M. An evaluation of semidistributed-pipe-network and distributed-finite-difference models to simulate karst systems. *Hydrogeol. J.* **2021**, *29*, 259–279. [CrossRef] [PubMed]

27. Dwarakish, G.S.; Ganasri, B.P. Impact of land-use change on hydrological systems: A review of current modeling approaches. *Cogent Geosci.* **2015**, *1*, 1115691. [CrossRef]

28. Jeannin, P.-Y.; Artigue, G.; Butscher, C.; Chang, Y.; Charlier, J.-B.; Duran, L.; Gill, L.; Hartmann, A.; Johannet, A.; Jourde, H.; et al. Karst modelling challenge 1: Results of hydrological modelling. *J. Hydrol.* **2021**, *600*, 126508. [CrossRef]

29. Doummar, J.; Sauter, M.; Geyer, T. 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.* **2012**, *426–427*, 112–123. [CrossRef]

30. Sarrazin, F.; Hartmann, A.; Pianosi, F.; Rosolem, R.; Wagener, T. V2Karst V1.1: A parsimonious large-scale integrated vegetation–recharge model to simulate the impact of climate and land cover change in karst regions. *Geosci. Model Dev.* **2018**, *11*, 4933–4964. [CrossRef]

31. Bittner, D.; Narany, T.S.; Kohl, B.; Disse, M.; Chiogna, G. Modeling the hydrological impact of land-use change in a dolomite-dominated karst system. *J. Hydrol.* **2018**, *567*, 267–279. [CrossRef]

32. Ruggieri, G.; Allocca, V.; Borfecchia, F.; Cusano, D.; Marsiglia, P.; De Vita, P. Testing Evapotranspiration Estimates Based on MODIS Satellite Data in the Assessment of the Groundwater Recharge of Karst Aquifers in Southern Italy. *Water* **2021**, *13*, 118. [CrossRef]

33. Ollivier, C.; Olioso, A.; Carrière, S.D.; Boulet, G.; Chalikakis, K.; Chanzy, A.; Charlier, J.-B.; Combemale, D.; Davi, H.; Emblanch, C.; et al. An evapotranspiration model driven by remote sensing data for assessing groundwater resource in karst watershed. *Sci. Total Environ.* **2021**, *781*, 146706. [CrossRef]

34. Yang, W.; Chen, L.; Chen, X.; Chen, H. Subdaily precipitation-streamflow modelling of the karst dominated basin using an improved grid-based distributed Xinanjiang hydrological model. *J. Hydrol. Reg. Stud.* **2022**, *42*, 101125. [CrossRef]

35. Kibii, J.K.; Kipkorir, E.C.; Kosgei, J.R. Application of Soil and Water Assessment Tool (SWAT) to Evaluate the Impact of Land-use and Climate Variability on the Kaptagat Catchment River Discharge. *Sustainability* **2021**, *13*, 1802. [CrossRef]

36. Eini, R.M.; Javadi, S.; Delavar, M.; Gassman, P.W.; Jarihani, B. Development of alternative SWAT-based models for simulating water budget components and streamflow for a karstic-influenced watershed. *CATENA* **2020**, *195*, 104801. [CrossRef]

37. Gassman, P.W.; Reyes, M.R.; Green, C.H.; Arnold, J.G. The soil and water assessment tool: Historical development, applications, and future research directions. *Trans. ASABE* **2007**, *50*, 1211–1250. [CrossRef]

38. CARD. SWAT Literature Database for Peer-Reviewed Journal Articles; Center for Agricultural and Rural Development, Iowa State University: Ames, IA, USA, 2022; Available online: https://www.card.iastate.edu/swatжиd articles/ (accessed on 24 November 2022).

39. Douglas-Mankin, K.R.; Srinivasan, R.; Arnold, J.G. Soil and Water Assessment Tool (SWAT) Model: Current Developments and Applications. *Trans. ASABE* **2010**, *53*, 1423–1431. [CrossRef]

40. Tuppad, P.; Douglas-Mankin, K.R.; Lee, T.; Srinivasan, R.; Arnold, J.G. Soil and Water Assessment Tool (SWAT) Hydrologic/Water Quality Model: Extended Capability and Wider Adoption. *Trans. ASABE* **2011**, *5*, 1677–1684. [CrossRef]

41. Francesconi, W.; Srinivasan, R.; Pérez-Miñana, E.; Willcock, S.P.; Quintero, M. Using the Soil and Water Assessment Tool (SWAT) to model ecosystem services: A systematic review. *J. Hydrol.* **2016**, *535*, 625–636. [CrossRef]

42. Wang, Y.; Jiang, R.; Xie, J.; Zhao, Y.; Yan, D.; Yang, S. Soil and water assessment tool (SWAT) model: A systemic review. *J. Coast. Res.* **2019**, *93*, 22–30. [CrossRef]

43. Brighenti, T.M.; Bonumá, N.B.; Srinivasan, R.; Chaffe, P.L.B. Simulating sub-daily hydrological process with SWAT: A review. *Hydrol. Sci. J.* **2019**, *64*, 1415–1423. [CrossRef]

44. Karki, R.; Srivastava, P.; Veith, T.L. Application of the Soil and Water Assessment Tool (SWAT) at Field Scale: Categorizing Methods and Review of Applications. *Trans. ASABE* **2020**, *63*, 513–522. [CrossRef]

45. Tan, M.L.; Gassman, P.W.; Yang, X.; Haywood, J. A review of SWAT applications, performance and future needs for simulation of hydro-climatic extremes. *Adv. Water Resour.* **2020**, *143*, 103662. [CrossRef]

46. Aloui, S.; Mazzoni, A.; Elomri, A.; Auissi, J.; Boufekane, A.; Zghibi, A. A review of Soil and Water Assessment Tool (SWAT) studies of Mediterranean catchments: Applications, feasibility, and future directions. *J. Environ. Manag.* **2023**, *326*, 116799. [CrossRef] [PubMed]

47. Arnold, J.G.; Srinivasan, R.; Muttiah, R.S.; Williams, J.R. Large area hydrologic modeling and assessment Part I: Model development. *J. Am. Water Resour. Assoc.* **1998**, *34*, 73–89. [CrossRef]
22. Lukač Reberski, J.; Rubinić, J.; Terzić, J.; Radišić, M. Climate Change Impacts on Groundwater Resources in the Coastal Karstic Adriatic Area: A Case Study from the Dinaric Karst. *Nat. Resour. Res.* **2020**, *29*, 1975–1988. [CrossRef]

23. Hartmann, A.; Goldscheider, N.; Wagener, T.; Lange, J.; Weiler, M. Karst water resources in a changing world: Review of hydrological modeling approaches. *Rev. Geophys.* **2014**, *52*, 218–242. [CrossRef]

24. Sivelle, V.; Jourde, H.; Bittner, D.; Mazzilli, N.; Tramblay, Y. Assessment of the relative impacts of climate changes and anthropogenic forcing on spring discharge of a Mediterranean karst system. *J. Hydrol.* **2021**, *598*, 126396. [CrossRef]

25. Mudarra, M.; Hartmann, A.; Andreo, B. Combining experimental methods and modeling to quantify the complex recharge behavior of karst aquifers. *Water Resour. Res.* **2019**, *55*, 1384–1404. [CrossRef]

26. Gill, L.W.; Schuler, P.; Duran, L.; Morrissey, P.; Johnston, P.M. An evaluation of semidistributed-pipe-network and distributed-finite-difference models to simulate karst systems. *Hydrogeol. J.* **2021**, *29*, 259–279. [CrossRef] [PubMed]

27. Dwarakish, G.S.; Ganasri, B.P. Impact of land-use change on hydrological systems: A review of current modeling approaches. *Cogent Geosci.* **2015**, *1*, 1115691. [CrossRef]

28. Jeannin, P.-Y.; Artigue, G.; Butscher, C.; Chang, Y.; Charlier, J.-B.; Duran, L.; Gill, L.; Hartmann, A.; Johannet, A.; Jourde, H.; et al. Karst modelling challenge 1: Results of hydrological modelling. *J. Hydrol.* **2021**, *600*, 126508. [CrossRef]

29. Doummar, J.; Sauter, M.; Geyer, T. 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.* **2012**, *426–427*, 112–123. [CrossRef]

30. Sarrazin, F.; Hartmann, A.; Pianosi, F.; Rosolem, R.; Wagener, T. V2Karst V1.1: A parsimonious large-scale integrated vegetation–recharge model to simulate the impact of climate and land cover change in karst regions. *Geosci. Model Dev.* **2018**, *11*, 4933–4964. [CrossRef]

31. Bittner, D.; Narany, T.S.; Kohl, B.; Disse, M.; Chiogna, G. Modeling the hydrological impact of land-use change in a dolomite-dominated karst system. *J. Hydrol.* **2018**, *567*, 267–279. [CrossRef]

32. Ruggieri, G.; Allocca, V.; Borfecchia, F.; Cusano, D.; Marsiglia, P.; De Vita, P. Testing Evapotranspiration Estimates Based on MODIS Satellite Data in the Assessment of the Groundwater Recharge of Karst Aquifers in Southern Italy. *Water* **2021**, *13*, 118. [CrossRef]

33. Ollivier, C.; Olioso, A.; Carrière, S.D.; Boulet, G.; Chalikakis, K.; Chanzy, A.; Charlier, J.-B.; Combemale, D.; Davi, H.; Emblanch, C.; et al. An evapotranspiration model driven by remote sensing data for assessing groundwater resource in karst watershed. *Sci. Total Environ.* **2021**, *781*, 146706. [CrossRef]

34. Yang, W.; Chen, L.; Chen, X.; Chen, H. Subdaily precipitation-streamflow modelling of the karst dominated basin using an improved grid-based distributed Xinanjiang hydrological model. *J. Hydrol. Reg. Stud.* **2022**, *42*, 101125. [CrossRef]

35. Kibii, J.K.; Kipkorir, E.C.; Kosgei, J.R. Application of Soil and Water Assessment Tool (SWAT) to Evaluate the Impact of Land-use and Climate Variability on the Kaptagat Catchment River Discharge. *Sustainability* **2021**, *13*, 1802. [CrossRef]

36. Eini, R.M.; Javadi, S.; Delavar, M.; Gassman, P.W.; Jarihani, B. Development of alternative SWAT-based models for simulating water budget components and streamflow for a karstic-influenced watershed. *CATENA* **2020**, *195*, 104801. [CrossRef]

37. Gassman, P.W.; Reyes, M.R.; Green, C.H.; Arnold, J.G. The soil and water assessment tool: Historical development, applications, and future research directions. *Trans. ASABE* **2007**, *50*, 1211–1250. [CrossRef]

38. CARD. SWAT Literature Database for Peer-Reviewed Journal Articles; Center for Agricultural and Rural Development, Iowa State University: Ames, IA, USA, 2022; Available online: https://www.card.iastate.edu/swatжиd articles/ (accessed on 24 November 2022).

39. Douglas-Mankin, K.R.; Srinivasan, R.; Arnold, J.G. Soil and Water Assessment Tool (SWAT) Model: Current Developments and Applications. *Trans. ASABE* **2010**, *53*, 1423–1431. [CrossRef]

40. Tuppad, P.; Douglas-Mankin, K.R.; Lee, T.; Srinivasan, R.; Arnold, J.G. Soil and Water Assessment Tool (SWAT) Hydrologic/Water Quality Model: Extended Capability and Wider Adoption. *Trans. ASABE* **2011**, *5*, 1677–1684. [CrossRef]

41. Francesconi, W.; Srinivasan, R.; Pérez-Miñana, E.; Willcock, S.P.; Quintero, M. Using the Soil and Water Assessment Tool (SWAT) to model ecosystem services: A systematic review. *J. Hydrol.* **2016**, *535*, 625–636. [CrossRef]

42. Wang, Y.; Jiang, R.; Xie, J.; Zhao, Y.; Yan, D.; Yang, S. Soil and water assessment tool (SWAT) model: A systemic review. *J. Coast. Res.* **2019**, *93*, 22–30. [CrossRef]

43. Brighenti, T.M.; Bonumá, N.B.; Srinivasan, R.; Chaffe, P.L.B. Simulating sub-daily hydrological process with SWAT: A review. *Hydrol. Sci. J.* **2019**, *64*, 1415–1423. [CrossRef]

44. Karki, R.; Srivastava, P.; Veith, T.L. Application of the Soil and Water Assessment Tool (SWAT) at Field Scale: Categorizing Methods and Review of Applications. *Trans. ASABE* **2020**, *63*, 513–522. [CrossRef]

45. Tan, M.L.; Gassman, P.W.; Yang, X.; Haywood, J. A review of SWAT applications, performance and future needs for simulation of hydro-climatic extremes. *Adv. Water Resour.* **2020**, *143*, 103662. [CrossRef]

46. Aloui, S.; Mazzoni, A.; Elomri, A.; Auissi, J.; Boufekane, A.; Zghibi, A. A review of Soil and Water Assessment Tool (SWAT) studies of Mediterranean catchments: Applications, feasibility, and future directions. *J. Environ. Manag.* **2023**, *326*, 116799. [CrossRef] [PubMed]

47. Arnold, J.G.; Srinivasan, R.; Muttiah, R.S.; Williams, J.R. Large area hydrologic modeling and assessment Part I: Model development. *J. Am. Water Resour. Assoc.* **1998**, *34*, 73–89. [CrossRef]
48. Kim, J.; Choi, J.; Choi, C.; Park, S. Impacts of changes in climate and land-use/land cover under IPCC RCP scenarios on streamflow in the Hoeya River Basin, Korea. *Sci. Total Environ.* **2013**, *452–453*, 181–195. [CrossRef]

49. Wang, R.; Kalin, L.; Kuang, W.; Hanqin, T. Individual and combined effects of land-use/cover and climate change on Wolf Bay watershed streamflow in southern Alabama. *Hydrol. Process.* **2014**, *28*, 5530–5546. [CrossRef]

50. Tan, M.L.; Ibrahim, A.L.; Yusop, Z.; Duan, Z.; Ling, L. Impacts of land-use and climate variability on hydrological components in the Johor River Basin, Malaysia. *Hydrol. Sci. J.* **2015**, *60*, 873–889. [CrossRef]

51. Mittal, N.; Bhave, A.G.; Mishra, A.; Singh, R. Impact of human intervention and climate change on natural flow regime. *Water Resour. Manag.* **2016**, *30*, 685–699. [CrossRef]

52. Setegn, S.G.; Darfahi, B.; Srinivasan, R.; Melesse, A.M. Modeling of Sediment Yield From Anjeni-Gauged Watershed, Ethiopia using SWAT Model. *J. Am. Water Resour. Assoc.* **2010**, *46*, 514–526. [CrossRef]

53. De Girolamo, A.M.; Bouraoui, F.; Buffagni, A.; Pappagallo, G.; Lo Porto, A. Hydrology under climate change in a temporary river system: Potential impact on water balance and flow regime. *River Res. Appl.* **2017**, *33*, 1219–1232. [CrossRef]

54. Marhaento, H.; Martijn, J.B.; Hoekstra, A.Y. Hydrological response to future land-use change and climate change in a tropical catchment. *Hydrol. Sci. J.* **2018**, *63*, 1368–1385. [CrossRef]

55. Boufala, M.; El Hmaidı, A.; Essahlaoui, A.; Chadli, K.; El Ouali, A.; Lahjouj, A. Assessment of the best management practices under a semi-arid basin using SWAT model (case of M’dez Watershed, Morocco). *Model. Earth Syst. Environ.* **2022**, *8*, 713–731. [CrossRef]

56. Zango, B.-S.; Seidou, O.; Sartaj, M.; Nakhael, N.; Stiles, K. Impacts of urbanization and climate change on water quantity and quality in the Carp River Watershed. *J. Water Clim. Change* **2022**, *13*, 786–816. [CrossRef]

57. Schilling, K.E.; Mount, J.; Suttles, K.M.; McLellan, E.L.; Gassman, P.W.; White, M.J.; Arnold, J.G. An Approach for Prioritizing Natural Infrastructure Practices to Mitigate Flood and Nitrate Risks in the Mississippi-Atchafalaya River Basin. *Land* **2023**, *12*, 276. [CrossRef]

58. Bennour, A.; Jia, L.; Menenti, M.; Zheng, C.; Zeng, Y.; Asenso Barnieh, B.; Jiang, M. Calibration and Validation of SWAT Model by Using Hydrological Remote Sensing Observables in the Lake Chad Basin. *Remote Sens.* **2022**, *14*, 1511. [CrossRef]

59. Yesuf, H.M.; Assen, M.; Alamirew, T.; Melesse, A.M. Modeling of sediment yield in Maybar gauged watershed using SWAT, northeast Ethiopia. *CATENA* **2015**, *127*, 191–205. [CrossRef]

60. Dile, Y.T.; Daggupati, P.; George, C.; Srinivasan, R.; Arnold, J. Introducing a new open source GIS user interface for the SWAT model. *Environ. Model. Softw.* **2016**, *85*, 129–138. [CrossRef]

61. Golmohammadi, G.; Prasher, S.; Madani, A.; Rudra, R. Evaluating Three Hydrological Distributed Watershed Models: MIKE-SHE, APEX, SWAT. *Hydrology* **2014**, *1*, 20–39. [CrossRef]

62. Mehdi, B.; Ludwig, R.; Lehner, B. Evaluating the impacts of climate change and crop land-use change on streamflow, nitrates and phosphorus: A modeling study in Bavaria. *J. Hydrol. Reg. Stud.* **2015**, *4*, 60–90. [CrossRef]

63. Bieger, K.; Georg Hörmann, G.; Fohrer, N. Detailed spatial analysis of SWAT-simulated surface runoff and sediment yield in a mountainous watershed in China. *Hydrol. Sci. J.* **2015**, *60*, 784–800. [CrossRef]

64. Neupane, R.P.; Kumar, S. Estimating the effects of potential climate and land-use changes on hydrologic processes of a large agriculture dominated watershed. *J. Hydrol.* **2015**, *529*, 418–429. [CrossRef]

65. Uniyal, B.; Dietrich, J. Modifying automatic irrigation in swat for plant water stress scheduling. *Agric. Water Manag.* **2019**, *223*, 105714. [CrossRef]

66. Thomas, T.; Ghosh, N.C.; Sudheer, K.P. Optimal reservoir operation—A climate change adaptation strategy for Narmada basin in central India. *J. Hydrol.* **2021**, *598*, 126238. [CrossRef]

67. Shao, G.; Zhang, D.; Guan, Y.; Xie, Y.; Huang, F. Application of SWAT Model with a Modified Groundwater Module to the Semi-Arid Hailiutu River Catchment, Northwest China. *Sustainability* **2019**, *11*, 2031. [CrossRef]

68. Ayivi, F.; Jha, M.K. Estimation of water balance and water yield in the Reedy Fork-Buffalo Creek Watershed in North Carolina using SWAT. *Int. Soil Water Conserv. Res.* **2018**, *6*, 203–213. [CrossRef]

69. Neitsch, S.L.; Arnold, J.G.; Kiniry, J.R.; Williams, J.R. *Soil and Water Assessment Tool Theoretical Documentation Version 2009*; Texas Water Resources Institute: College Station, TX, USA, 2011.

70. Monteith, J.L. Evaporation and environment. *Symp. Soc. Exp. Biol.* **1965**, *19*, 205–234. [PubMed]

71. Priestley, C.H.B.; Taylor, R.J. On the assessment of surface heat flux and evaporation using large-scale parameters. *Mon. Weather Rev.* **1972**, *100*, 81–92. [CrossRef]

72. Hargreaves, G.H.; Hargreaves, G.L.; Riley, J.P. Agricultural benefits for Senegal River Basin. *J. Irrig. Drain. Eng.* **1985**, *111*, 113–124. [CrossRef]

73. Strauch, M.; Volk, M. SWAT plant growth modification for improved modeling of perennial vegetation in the tropics. *Ecol. Model.* **2013**, *269*, 98–112. [CrossRef]

74. Ferreira, A.D.N.; de Almeida, A.; Koide, S.; Minoti, R.T.; Siqueira, M.B.B.D. Evaluation of Evapotranspiration in Brazilian Cerrado Biome Simulated with the SWAT Model. *Water* **2021**, *13*, 2037. [CrossRef]

75. Alemayehu, T.; van Griensven, A.; Woldegiorgis, B.T.; Bauwens, W. An improved SWAT vegetation growth module and its evaluation for four tropical ecosystems. *Hydrol. Earth Syst. Sci.* **2017**, *21*, 4449–4467. [CrossRef]
48. Kim, J.; Choi, J.; Choi, C.; Park, S. Impacts of changes in climate and land-use/land cover under IPCC RCP scenarios on streamflow in the Hoeya River Basin, Korea. *Sci. Total Environ.* **2013**, *452–453*, 181–195. [CrossRef]

49. Wang, R.; Kalin, L.; Kuang, W.; Hanqin, T. Individual and combined effects of land-use/cover and climate change on Wolf Bay watershed streamflow in southern Alabama. *Hydrol. Process.* **2014**, *28*, 5530–5546. [CrossRef]

50. Tan, M.L.; Ibrahim, A.L.; Yusop, Z.; Duan, Z.; Ling, L. Impacts of land-use and climate variability on hydrological components in the Johor River Basin, Malaysia. *Hydrol. Sci. J.* **2015**, *60*, 873–889. [CrossRef]

51. Mittal, N.; Bhave, A.G.; Mishra, A.; Singh, R. Impact of human intervention and climate change on natural flow regime. *Water Resour. Manag.* **2016**, *30*, 685–699. [CrossRef]

52. Setegn, S.G.; Darfahi, B.; Srinivasan, R.; Melesse, A.M. Modeling of Sediment Yield From Anjeni-Gauged Watershed, Ethiopia using SWAT Model. *J. Am. Water Resour. Assoc.* **2010**, *46*, 514–526. [CrossRef]

53. De Girolamo, A.M.; Bouraoui, F.; Buffagni, A.; Pappagallo, G.; Lo Porto, A. Hydrology under climate change in a temporary river system: Potential impact on water balance and flow regime. *River Res. Appl.* **2017**, *33*, 1219–1232. [CrossRef]

54. Marhaento, H.; Martijn, J.B.; Hoekstra, A.Y. Hydrological response to future land-use change and climate change in a tropical catchment. *Hydrol. Sci. J.* **2018**, *63*, 1368–1385. [CrossRef]

55. Boufala, M.; El Hmaidı, A.; Essahlaoui, A.; Chadli, K.; El Ouali, A.; Lahjouj, A. Assessment of the best management practices under a semi-arid basin using SWAT model (case of M’dez Watershed, Morocco). *Model. Earth Syst. Environ.* **2022**, *8*, 713–731. [CrossRef]

56. Zango, B.-S.; Seidou, O.; Sartaj, M.; Nakhael, N.; Stiles, K. Impacts of urbanization and climate change on water quantity and quality in the Carp River Watershed. *J. Water Clim. Change* **2022**, *13*, 786–816. [CrossRef]

57. Schilling, K.E.; Mount, J.; Suttles, K.M.; McLellan, E.L.; Gassman, P.W.; White, M.J.; Arnold, J.G. An Approach for Prioritizing Natural Infrastructure Practices to Mitigate Flood and Nitrate Risks in the Mississippi-Atchafalaya River Basin. *Land* **2023**, *12*, 276. [CrossRef]

58. Bennour, A.; Jia, L.; Menenti, M.; Zheng, C.; Zeng, Y.; Asenso Barnieh, B.; Jiang, M. Calibration and Validation of SWAT Model by Using Hydrological Remote Sensing Observables in the Lake Chad Basin. *Remote Sens.* **2022**, *14*, 1511. [CrossRef]

59. Yesuf, H.M.; Assen, M.; Alamirew, T.; Melesse, A.M. Modeling of sediment yield in Maybar gauged watershed using SWAT, northeast Ethiopia. *CATENA* **2015**, *127*, 191–205. [CrossRef]

60. Dile, Y.T.; Daggupati, P.; George, C.; Srinivasan, R.; Arnold, J. Introducing a new open source GIS user interface for the SWAT model. *Environ. Model. Softw.* **2016**, *85*, 129–138. [CrossRef]

61. Golmohammadi, G.; Prasher, S.; Madani, A.; Rudra, R. Evaluating Three Hydrological Distributed Watershed Models: MIKE-SHE, APEX, SWAT. *Hydrology* **2014**, *1*, 20–39. [CrossRef]

62. Mehdi, B.; Ludwig, R.; Lehner, B. Evaluating the impacts of climate change and crop land-use change on streamflow, nitrates and phosphorus: A modeling study in Bavaria. *J. Hydrol. Reg. Stud.* **2015**, *4*, 60–90. [CrossRef]

63. Bieger, K.; Georg Hörmann, G.; Fohrer, N. Detailed spatial analysis of SWAT-simulated surface runoff and sediment yield in a mountainous watershed in China. *Hydrol. Sci. J.* **2015**, *60*, 784–800. [CrossRef]

64. Neupane, R.P.; Kumar, S. Estimating the effects of potential climate and land-use changes on hydrologic processes of a large agriculture dominated watershed. *J. Hydrol.* **2015**, *529*, 418–429. [CrossRef]

65. Uniyal, B.; Dietrich, J. Modifying automatic irrigation in swat for plant water stress scheduling. *Agric. Water Manag.* **2019**, *223*, 105714. [CrossRef]

66. Thomas, T.; Ghosh, N.C.; Sudheer, K.P. Optimal reservoir operation—A climate change adaptation strategy for Narmada basin in central India. *J. Hydrol.* **2021**, *598*, 126238. [CrossRef]

67. Shao, G.; Zhang, D.; Guan, Y.; Xie, Y.; Huang, F. Application of SWAT Model with a Modified Groundwater Module to the Semi-Arid Hailiutu River Catchment, Northwest China. *Sustainability* **2019**, *11*, 2031. [CrossRef]

68. Ayivi, F.; Jha, M.K. Estimation of water balance and water yield in the Reedy Fork-Buffalo Creek Watershed in North Carolina using SWAT. *Int. Soil Water Conserv. Res.* **2018**, *6*, 203–213. [CrossRef]

69. Neitsch, S.L.; Arnold, J.G.; Kiniry, J.R.; Williams, J.R. *Soil and Water Assessment Tool Theoretical Documentation Version 2009*; Texas Water Resources Institute: College Station, TX, USA, 2011.

70. Monteith, J.L. Evaporation and environment. *Symp. Soc. Exp. Biol.* **1965**, *19*, 205–234. [PubMed]

71. Priestley, C.H.B.; Taylor, R.J. On the assessment of surface heat flux and evaporation using large-scale parameters. *Mon. Weather Rev.* **1972**, *100*, 81–92. [CrossRef]

72. Hargreaves, G.H.; Hargreaves, G.L.; Riley, J.P. Agricultural benefits for Senegal River Basin. *J. Irrig. Drain. Eng.* **1985**, *111*, 113–124. [CrossRef]

73. Strauch, M.; Volk, M. SWAT plant growth modification for improved modeling of perennial vegetation in the tropics. *Ecol. Model.* **2013**, *269*, 98–112. [CrossRef]

74. Ferreira, A.D.N.; de Almeida, A.; Koide, S.; Minoti, R.T.; Siqueira, M.B.B.D. Evaluation of Evapotranspiration in Brazilian Cerrado Biome Simulated with the SWAT Model. *Water* **2021**, *13*, 2037. [CrossRef]

75. Alemayehu, T.; van Griensven, A.; Woldegiorgis, B.T.; Bauwens, W. An improved SWAT vegetation growth module and its evaluation for four tropical ecosystems. *Hydrol. Earth Syst. Sci.* **2017**, *21*, 4449–4467. [CrossRef]
76. Ma, T.; Duan, Z.; Runkui Li, R.; Song, X. Enhancing SWAT with remotely sensed LAI for improved modelling of ecohydrological process in subtropics. *J. Hydrol.* **2019**, *570*, 802–815. [CrossRef]

77. Sinnathamby, S.; Douglas-Mankin, K.R.; Craige, C. Field-scale calibration of crop-yield parameters in the Soil and Water Assessment Tool (SWAT). *Agric. Water Manag.* **2017**, *180*, 61–69. [CrossRef]

78. Luo, Y.; He, C.; Sophocleous, M.; Yin, Z.; Hongrui, R.; Ouyang, Z. Assessment of crop growth and soil water modules in SWAT 2000 using extensive field experiment data in an irrigation district of the Yellow River Basin. *J. Hydrol.* **2008**, *352*, 139–156. [CrossRef]

79. Abiodun, O.O.; Guan, H.; Post, A.V.E.; Batelaan, O. Comparison of MODIS and SWAT evapotranspiration over a complex terrain at different spatial scales. *Hydrol. Earth Syst. Sci.* **2018**, *22*, 2775–2794. [CrossRef]

80. Aouissi, J.; Benabdallah, S.; Chabaâne, Z.L.; Cudennec, C. Evaluation of potential evapotranspiration assessment methods for hydrological modelling with SWAT—Application in data-scarce rural Tunisia. *Agric. Water Manag.* **2016**, *174*, 39–51. [CrossRef]

81. Campbell, A.; Pradhanang, S.M.; Anbaran, S.K.; Sargent, J.; Palmer, Z.; Audette, M. Assessing the impact of urbanization on flood risk and severity for the Pawtuxet watershed, Rhode Island. *Lake Reserv. Manag.* **2018**, *34*, 74–87. [CrossRef]

82. Bacopoulos, P.; Tang, Y.; Wang, D.; Hagen, S.C. Integrated Hydrologic-Hydrodynamic Modeling of Estuarine-Riverine Flooding: 2008 Tropical Storm Fay. *J. Hydrol. Eng.* **2017**, *22*, 04017022. [CrossRef]

83. Her, Y.; Jeong, J.; Arnold, J.; Gosselink, L.; Glick, R.; Jaber, F. A new framework for modeling decentralized low impact developments using Soil and Water Assessment Tool. *Environ. Model. Softw.* **2017**, *96*, 305–322. [CrossRef]

84. Her, Y.; Jeong, J.; Bieger, K.; Rathjens, H.; Arnold, J.; Srinivasan, R. Implications of Conceptual Channel Representation on SWAT Streamflow and Sediment Modeling. *J. Am. Water Resour. Assoc.* **2017**, *53*, 725–747. [CrossRef]

85. Malagò, A.; Vigiak, O.; Bouraoui, F.; Pagliero, L.; Franchini, M. The Hillslope Length Impact on SWAT Streamflow Prediction in Large Basins. *J. Environ. Inform.* **2018**, *32*, 82–97. [CrossRef]

86. Cunge, J.A. On the subject of a flood propagation computation method (Muskingum method). *J. Hydraul. Res.* **1969**, *7*, 205–230. [CrossRef]

87. Williams, J.R. Flood routing with variable travel time or variable storage coefficients. *Trans. ASABE* **1969**, *12*, 0100–0103. [CrossRef]

88. Nguyen, V.T.; Dietrich, J.; Uniyal, B.; Tran, D.A. Verification and Correction of the Hydrologic Routing in the Soil and Water Assessment Tool. *Water* **2018**, *10*, 1419. [CrossRef]

89. Holvoet, K.; van Griensven, A.; Gevaert, V.; Seuntjens, P.; Vanrolleghem, P.A. Modifications to the SWAT code for modelling direct pesticide losses. *Environ. Model. Softw.* **2008**, *23*, 72–81. [CrossRef]

90. Baffaut, C.; John Sadler, E.; Ghidey, F.; Anderson, S.H. Long-Term Agroecosystem Research in the Central Mississippi River Basin: SWAT Simulation of Flow and Water Quality in the Goodwater Creek Experimental Watershed. *J. Environ. Qual.* **2015**, *44*, 84–96. [CrossRef]

91. Rahbeh, M.; Srinivasan, R.; Mohtar, R. Numerical and conceptual evaluation of preferential flow in Zarqa River Basin, Jordan. *Ecohydrol. Hydrobiol.* **2019**, *19*, 224–237. [CrossRef]

92. Mapes, K.L.; Pricope, N.G. Evaluating SWAT Model Performance for Runoff, Percolation, and Sediment Loss Estimation in Low-Gradient Watersheds of the Atlantic Coastal Plain. *Hydrology* **2020**, *7*, 21. [CrossRef]

93. Fu, C.; James, A.L.; Yao, H. SWAT-CS: Revision and testing of SWAT for Canadian Shield catchments. *J. Hydrol.* **2014**, *511*, 719–735. [CrossRef]

94. Qi, J.; Zhang, X.; McCartyc, G.W.; Sadeghi, A.M.; Cosh, M.H.; Zeng, X.; Gao, F.; Daughtry, C.S.T.; Huang, C.; Lang, M.W.; et al. Assessing the performance of a physically-based soil moisture module integrated within the Soil and Water Assessment Tool. *Environ. Model. Softw.* **2018**, *109*, 329–341. [CrossRef]

95. Luo, Y.; Arnold, J.; Allen, P.; Chen, X. Baseflow simulation using SWAT model in an inland river basin in Tianshan Mountains, Northwest China. *Hydrol. Earth Syst. Sci.* **2012**, *16*, 1259–1267. [CrossRef]

96. Uniyal, B.; Dietrich, J.; Vu, N.Q.; Jha, M.K.; Arumí, J.L. Simulation of regional irrigation requirement with SWAT in different agro-climatic zones driven by observed climate and two reanalysis datasets. *Sci. Total Environ.* **2019**, *649*, 846–865. [CrossRef] [PubMed]

97. Xie, H.; Longuevergne, L.; Ringler, C.; Scanlon, B.R. Integrating groundwater irrigation into hydrological simulation of India: Case of improving model representation of anthropogenic water use impact using GRACE. *J. Hydrol. Reg. Stud.* **2020**, *29*, 100681. [CrossRef]

98. Nguyen, V.T.; Dietrich, J. Modification of the SWAT model to simulate regional groundwater flow using a multicell aquifer. *Hydrol. Process.* **2018**, *32*, 939–953. [CrossRef]

99. Phiri, W.K.; Vanzo, D.; Banda, K.; Nyirenda, E.; Nyambe, I.A. A pseudo-reservoir concept in SWAT model for the simulation of an alluvial floodplain in a complex tropical river system. *J. Hydrol. Reg. Stud.* **2021**, *33*, 100770. [CrossRef]

100. Rahman, M.M.; Thompson, J.R.; Flower, R.J. An Enhanced SWAT Wetland Module to Quantify Hydraulic Interactions between Riparian Depressional Wetlands, Rivers and Aquifers. *Environ. Model. Softw.* **2016**, *84*, 263–289. [CrossRef]

101. Dash, S.D.; Sahoo, B.; Raghuwanshi, N.S. A novel embedded pothole module for Soil and Water Assessment Tool (SWAT) improving streamflow estimation in paddy-dominated catchments. *J. Hydrol.* **2020**, *588*, 125103. [CrossRef]

102. Moriasi, D.N.; Gitau, M.W.; Pai, N.; Daggupati, P. Hydrologic and Water Quality Models: Performance Measures and Evaluation Criteria. *Trans. ASABE* **2015**, *58*, 1763–1785. [CrossRef]
76. Ma, T.; Duan, Z.; Runkui Li, R.; Song, X. Enhancing SWAT with remotely sensed LAI for improved modelling of ecohydrological process in subtropics. *J. Hydrol.* **2019**, *570*, 802–815. [CrossRef]

77. Sinnathamby, S.; Douglas-Mankin, K.R.; Craige, C. Field-scale calibration of crop-yield parameters in the Soil and Water Assessment Tool (SWAT). *Agric. Water Manag.* **2017**, *180*, 61–69. [CrossRef]

78. Luo, Y.; He, C.; Sophocleous, M.; Yin, Z.; Hongrui, R.; Ouyang, Z. Assessment of crop growth and soil water modules in SWAT 2000 using extensive field experiment data in an irrigation district of the Yellow River Basin. *J. Hydrol.* **2008**, *352*, 139–156. [CrossRef]

79. Abiodun, O.O.; Guan, H.; Post, A.V.E.; Batelaan, O. Comparison of MODIS and SWAT evapotranspiration over a complex terrain at different spatial scales. *Hydrol. Earth Syst. Sci.* **2018**, *22*, 2775–2794. [CrossRef]

80. Aouissi, J.; Benabdallah, S.; Chabaâne, Z.L.; Cudennec, C. Evaluation of potential evapotranspiration assessment methods for hydrological modelling with SWAT—Application in data-scarce rural Tunisia. *Agric. Water Manag.* **2016**, *174*, 39–51. [CrossRef]

81. Campbell, A.; Pradhanang, S.M.; Anbaran, S.K.; Sargent, J.; Palmer, Z.; Audette, M. Assessing the impact of urbanization on flood risk and severity for the Pawtuxet watershed, Rhode Island. *Lake Reserv. Manag.* **2018**, *34*, 74–87. [CrossRef]

82. Bacopoulos, P.; Tang, Y.; Wang, D.; Hagen, S.C. Integrated Hydrologic-Hydrodynamic Modeling of Estuarine-Riverine Flooding: 2008 Tropical Storm Fay. *J. Hydrol. Eng.* **2017**, *22*, 04017022. [CrossRef]

83. Her, Y.; Jeong, J.; Arnold, J.; Gosselink, L.; Glick, R.; Jaber, F. A new framework for modeling decentralized low impact developments using Soil and Water Assessment Tool. *Environ. Model. Softw.* **2017**, *96*, 305–322. [CrossRef]

84. Her, Y.; Jeong, J.; Bieger, K.; Rathjens, H.; Arnold, J.; Srinivasan, R. Implications of Conceptual Channel Representation on SWAT Streamflow and Sediment Modeling. *J. Am. Water Resour. Assoc.* **2017**, *53*, 725–747. [CrossRef]

85. Malagò, A.; Vigiak, O.; Bouraoui, F.; Pagliero, L.; Franchini, M. The Hillslope Length Impact on SWAT Streamflow Prediction in Large Basins. *J. Environ. Inform.* **2018**, *32*, 82–97. [CrossRef]

86. Cunge, J.A. On the subject of a flood propagation computation method (Muskingum method). *J. Hydraul. Res.* **1969**, *7*, 205–230. [CrossRef]

87. Williams, J.R. Flood routing with variable travel time or variable storage coefficients. *Trans. ASABE* **1969**, *12*, 0100–0103. [CrossRef]

88. Nguyen, V.T.; Dietrich, J.; Uniyal, B.; Tran, D.A. Verification and Correction of the Hydrologic Routing in the Soil and Water Assessment Tool. *Water* **2018**, *10*, 1419. [CrossRef]

89. Holvoet, K.; van Griensven, A.; Gevaert, V.; Seuntjens, P.; Vanrolleghem, P.A. Modifications to the SWAT code for modelling direct pesticide losses. *Environ. Model. Softw.* **2008**, *23*, 72–81. [CrossRef]

90. Baffaut, C.; John Sadler, E.; Ghidey, F.; Anderson, S.H. Long-Term Agroecosystem Research in the Central Mississippi River Basin: SWAT Simulation of Flow and Water Quality in the Goodwater Creek Experimental Watershed. *J. Environ. Qual.* **2015**, *44*, 84–96. [CrossRef]

91. Rahbeh, M.; Srinivasan, R.; Mohtar, R. Numerical and conceptual evaluation of preferential flow in Zarqa River Basin, Jordan. *Ecohydrol. Hydrobiol.* **2019**, *19*, 224–237. [CrossRef]

92. Mapes, K.L.; Pricope, N.G. Evaluating SWAT Model Performance for Runoff, Percolation, and Sediment Loss Estimation in Low-Gradient Watersheds of the Atlantic Coastal Plain. *Hydrology* **2020**, *7*, 21. [CrossRef]

93. Fu, C.; James, A.L.; Yao, H. SWAT-CS: Revision and testing of SWAT for Canadian Shield catchments. *J. Hydrol.* **2014**, *511*, 719–735. [CrossRef]

94. Qi, J.; Zhang, X.; McCartyc, G.W.; Sadeghi, A.M.; Cosh, M.H.; Zeng, X.; Gao, F.; Daughtry, C.S.T.; Huang, C.; Lang, M.W.; et al. Assessing the performance of a physically-based soil moisture module integrated within the Soil and Water Assessment Tool. *Environ. Model. Softw.* **2018**, *109*, 329–341. [CrossRef]

95. Luo, Y.; Arnold, J.; Allen, P.; Chen, X. Baseflow simulation using SWAT model in an inland river basin in Tianshan Mountains, Northwest China. *Hydrol. Earth Syst. Sci.* **2012**, *16*, 1259–1267. [CrossRef]

96. Uniyal, B.; Dietrich, J.; Vu, N.Q.; Jha, M.K.; Arumí, J.L. Simulation of regional irrigation requirement with SWAT in different agro-climatic zones driven by observed climate and two reanalysis datasets. *Sci. Total Environ.* **2019**, *649*, 846–865. [CrossRef] [PubMed]

97. Xie, H.; Longuevergne, L.; Ringler, C.; Scanlon, B.R. Integrating groundwater irrigation into hydrological simulation of India: Case of improving model representation of anthropogenic water use impact using GRACE. *J. Hydrol. Reg. Stud.* **2020**, *29*, 100681. [CrossRef]

98. Nguyen, V.T.; Dietrich, J. Modification of the SWAT model to simulate regional groundwater flow using a multicell aquifer. *Hydrol. Process.* **2018**, *32*, 939–953. [CrossRef]

99. Phiri, W.K.; Vanzo, D.; Banda, K.; Nyirenda, E.; Nyambe, I.A. A pseudo-reservoir concept in SWAT model for the simulation of an alluvial floodplain in a complex tropical river system. *J. Hydrol. Reg. Stud.* **2021**, *33*, 100770. [CrossRef]

100. Rahman, M.M.; Thompson, J.R.; Flower, R.J. An Enhanced SWAT Wetland Module to Quantify Hydraulic Interactions between Riparian Depressional Wetlands, Rivers and Aquifers. *Environ. Model. Softw.* **2016**, *84*, 263–289. [CrossRef]

101. Dash, S.D.; Sahoo, B.; Raghuwanshi, N.S. A novel embedded pothole module for Soil and Water Assessment Tool (SWAT) improving streamflow estimation in paddy-dominated catchments. *J. Hydrol.* **2020**, *588*, 125103. [CrossRef]

102. Moriasi, D.N.; Gitau, M.W.; Pai, N.; Daggupati, P. Hydrologic and Water Quality Models: Performance Measures and Evaluation Criteria. *Trans. ASABE* **2015**, *58*, 1763–1785. [CrossRef]
103. Guse, B.; Pfannerstill, M.; Gafurov, A.; Kiesel, J.; Lehr, C.; Fohrer, N. Identifying the connective strength between model parameters and performance criteria. *Hydrol. Earth Syst. Sci.* **2017**, *21*, 5663–5679. [CrossRef]

104. Spruill, C.; Workman, S.; Taraba, J. Simulation of daily and monthly stream discharge from small watersheds using the SWAT model. *Trans. ASABE 2000*, *43*, 1431. [CrossRef]

105. Coffey, M.E.; Workman, S.R.; Taraba, J.L.; Fogle, A.W. Statistical procedures for evaluating daily and monthly hydrologic model predictions. *Trans. ASABE 2004*, *47*, 59–68. [CrossRef]

106. Benham, B.L.; Baffaut, C.; Zeckoski, R.W.; Mankin, K.R.; Pachepsky, Y.A.; Sadeghi, A.M.; Brannan, K.M.; Soupir, M.L.; Habersack, M.J. Modeling bacteria fate and transport in watersheds to support TMDLs. *Trans. ASABE 2006*, *49*, 987–1002. [CrossRef]

107. Amatya, D.M.; Jha, M.; Edwards, A.E.; Williams, T.M.; Hitchcock, D.R. SWAT-Based Streamflow and Embayment Modeling of Karst-Affected Chapel Branch Watershed, South Carolina. *Trans. ASABE 2011*, *4*, 1311–1323. [CrossRef]

108. Amatya, D.M.; Jha, M.K.; Williams, T.M.; Edwards, A.E.; Hitchcock, D.R. SWAT Model Prediction of Phosphorus Loading in a South Carolina Karst Watershed with a Downstream Embayment. *J. Environ. Prot.* **2013**, *4*, 75–90. [CrossRef]

109. Williams, T.M.; Amatya, D.M.; Hitchcock, D.R.; Edwards, A.E. Streamflow and Nutrients from a Karst Watershed with a Downstream Embayment: Chapel Branch Creek. *J. Hydrol. Eng.* **2014**, *19*, 428–438. [CrossRef]

110. Wilson, G.L.; Dalzell, B.J.; Mulla, D.J.; Dogwiler, T.; Porter, P.M. Estimating water quality effects of conservation practices and grazing land use scenarios. *J. Soil Water Conserv.* **2014**, *69*, 330–342. [CrossRef]

111. Jain, S.; Ale, S.; Munster, C.L.; Ansley, R.J.; Kiniry, J.R. Simulating the Hydrologic Impact of *Arundo donax* Invasion on the Headwaters of the Nueces River in Texas. *Hydrology* **2015**, *2*, 134–147. [CrossRef]

112. Sunde, M.; He, H.S.; Hubbart, J.A.; Scroggins, G. Forecasting streamflow response to increased imperviousness in an urbanizing Midwestern watershed using a coupled modeling approach. *Appl. Geogr.* **2016**, *72*, 14–25. [CrossRef]

113. Sunde, M.G.; He, H.S.; Hubbart, J.A.; Urban, M.A. Integrating downscaled CMIP5 data with a physically based hydrologic model to estimate potential climate change impacts on streamflow processes in a mixed-use watershed. *Hydrol. Process.* **2017**, *31*, 1790–1803. [CrossRef]

114. Sunde, M.G.; He, H.S.; Hubbart, J.A.; Urban, M.A. An integrated modeling approach for estimating hydrologic responses to future urbanization and climate changes in a mixed-use midwestern watershed. *J. Environ. Manag.* **2018**, *220*, 149–162. [CrossRef]

115. Sarkar, S.; Yonce, H.N.; Keeley, A.; Canfield, T.J.; Butcher, J.B.; Paul, M.J. Integration of SWAT and HSPF for Simulation of Sediment Sources in Legacy Sediment-Impacted Agricultural Watersheds. *J. Am. Water Resour. Assoc.* **2019**, *55*, 497–510. [CrossRef]

116. Merriman, K.R.; Daggupati, P.; Srinivasan, R.; Hayhurst, B. Assessment of site-specific agricultural Best Management Practices in the Upper East River watershed, Wisconsin, using a field-scale SWAT model. *J. Great Lakes Res.* **2019**, *45*, 619–641. [CrossRef]

117. Sullivan, T.P.; Gao, Y.; Reimann, T. Nitrate transport in a karst aquifer: Numerical model development and source evaluation. *J. Hydrol.* **2019**, *573*, 432–448. [CrossRef]

118. Chen, H.; Wang, S.; Wang, Y.; Zhu, J. Probabilistic projections of hydrological droughts through convection-permitting climate simulations and multimodel hydrological predictions. *J. Geophys. Res. Atmos.* **2020**, *125*, e2020JD032914. [CrossRef]

119. Al Aamery, N.; Adams, E.; Fox, J.; Husic, A.; Zhu, J.; Gerlitz, M.; Agouridis, C.; Bettel, L. Numerical model development for investigating hydrologic pathways in shallow fluviokarst. *J. Hydrol.* **2021**, *593*, 125844. [CrossRef]

120. Karki, R.; Srivastava, P.; Kalin, L.; Mitra, S.; Singh, S. Assessment of impact in groundwater levels and stream-aquifer interaction due to increased groundwater withdrawal in the lower Apalachicola-Chattahoochee-Flint (ACF) River Basin using MODFLOW. *J. Hydrol. Reg. Stud.* **2021**, *34*, 100802. [CrossRef]

121. Salerno, F.; Tartari, G. A coupled approach of surface hydrological modelling and Wavelet Analysis for understanding the baseflow components of river discharge in karst environments. *J. Hydrol.* **2009**, *376*, 295–306. [CrossRef]

122. Vale, M.; Holman, I.P. Understanding the hydrological functioning of a shallow lake system within a coastal karstic aquifer in Wales, UK. *J. Hydrol.* **2009**, *376*, 285–294. [CrossRef]

123. Tzoraki, O.; Cooper, D.; Kjeldsen, T.; Nikolaidis, N.P.; Gamvroudis, C.; Froebrich, J.; Querner, E.; Gallart, F.; Karalemas, N. Flood generation and classification of a semi-arid intermittent flow watershed: Evrotas river. *Int. J. River Basin Manag.* **2013**, *11*, 77–92. [CrossRef]

124. Palazón, L.; Navas, A. Sediment production of an alpine catchment with SWAT. *Z. Geomorphol* **2013**, *57*, 69–85. [CrossRef]

125. Palazón, L.; Navas, A. Modeling sediment sources and yields in a Pyrenean catchment draining to a large reservoir (Ésera River, Ebro Basin). *J. Soils Sediments* **2014**, *14*, 1612–1625. [CrossRef]

126. Sellami, H.; La Jeunesse, I.; Benabdallah, S.; Baghdadi, N.; Vanclooster, M. Uncertainty analysis in model parameters regionalization: A case study involving the SWAT model in Mediterranean catchments (Southern France). *Hydrol. Earth Syst. Sci.* **2014**, *18*, 2393–2413. [CrossRef]

127. Gamvroudis, C.; Nikolaidis, N.; Tzoraki, O.; Papadoulakis, V.; Karalemas, N. Water and sediment transport modeling of a large temporary river basin in Greece. *Sci. Total Environ.* **2015**, *508*, 354–365. [CrossRef] [PubMed]

128. Malagò, A.; Pagliero, L.; Bouraoui, F.; Franchini, M. Comparing calibrated parameter sets of the SWAT model for the Scandinavian and Iberian peninsulas. *Hydrol. Sci. J.* **2015**, *60*, 949–967. [CrossRef]

129. Sellami, H.; Benabdallah, S.; La Jeunesse, I.; Vanclooster, M. Quantifying hydrological responses of small Mediterranean catchments under climate change projections. *Sci. Total Environ.* **2016**, *543*, 924–936. [CrossRef]
103. Guse, B.; Pfannerstill, M.; Gafurov, A.; Kiesel, J.; Lehr, C.; Fohrer, N. Identifying the connective strength between model parameters and performance criteria. *Hydrol. Earth Syst. Sci.* **2017**, *21*, 5663–5679. [CrossRef]

104. Spruill, C.; Workman, S.; Taraba, J. Simulation of daily and monthly stream discharge from small watersheds using the SWAT model. *Trans. ASABE 2000*, *43*, 1431. [CrossRef]

105. Coffey, M.E.; Workman, S.R.; Taraba, J.L.; Fogle, A.W. Statistical procedures for evaluating daily and monthly hydrologic model predictions. *Trans. ASABE 2004*, *47*, 59–68. [CrossRef]

106. Benham, B.L.; Baffaut, C.; Zeckoski, R.W.; Mankin, K.R.; Pachepsky, Y.A.; Sadeghi, A.M.; Brannan, K.M.; Soupir, M.L.; Habersack, M.J. Modeling bacteria fate and transport in watersheds to support TMDLs. *Trans. ASABE 2006*, *49*, 987–1002. [CrossRef]

107. Amatya, D.M.; Jha, M.; Edwards, A.E.; Williams, T.M.; Hitchcock, D.R. SWAT-Based Streamflow and Embayment Modeling of Karst-Affected Chapel Branch Watershed, South Carolina. *Trans. ASABE 2011*, *4*, 1311–1323. [CrossRef]

108. Amatya, D.M.; Jha, M.K.; Williams, T.M.; Edwards, A.E.; Hitchcock, D.R. SWAT Model Prediction of Phosphorus Loading in a South Carolina Karst Watershed with a Downstream Embayment. *J. Environ. Prot.* **2013**, *4*, 75–90. [CrossRef]

109. Williams, T.M.; Amatya, D.M.; Hitchcock, D.R.; Edwards, A.E. Streamflow and Nutrients from a Karst Watershed with a Downstream Embayment: Chapel Branch Creek. *J. Hydrol. Eng.* **2014**, *19*, 428–438. [CrossRef]

110. Wilson, G.L.; Dalzell, B.J.; Mulla, D.J.; Dogwiler, T.; Porter, P.M. Estimating water quality effects of conservation practices and grazing land use scenarios. *J. Soil Water Conserv.* **2014**, *69*, 330–342. [CrossRef]

111. Jain, S.; Ale, S.; Munster, C.L.; Ansley, R.J.; Kiniry, J.R. Simulating the Hydrologic Impact of *Arundo donax* Invasion on the Headwaters of the Nueces River in Texas. *Hydrology* **2015**, *2*, 134–147. [CrossRef]

112. Sunde, M.; He, H.S.; Hubbart, J.A.; Scroggins, G. Forecasting streamflow response to increased imperviousness in an urbanizing Midwestern watershed using a coupled modeling approach. *Appl. Geogr.* **2016**, *72*, 14–25. [CrossRef]

113. Sunde, M.G.; He, H.S.; Hubbart, J.A.; Urban, M.A. Integrating downscaled CMIP5 data with a physically based hydrologic model to estimate potential climate change impacts on streamflow processes in a mixed-use watershed. *Hydrol. Process.* **2017**, *31*, 1790–1803. [CrossRef]

114. Sunde, M.G.; He, H.S.; Hubbart, J.A.; Urban, M.A. An integrated modeling approach for estimating hydrologic responses to future urbanization and climate changes in a mixed-use midwestern watershed. *J. Environ. Manag.* **2018**, *220*, 149–162. [CrossRef]

115. Sarkar, S.; Yonce, H.N.; Keeley, A.; Canfield, T.J.; Butcher, J.B.; Paul, M.J. Integration of SWAT and HSPF for Simulation of Sediment Sources in Legacy Sediment-Impacted Agricultural Watersheds. *J. Am. Water Resour. Assoc.* **2019**, *55*, 497–510. [CrossRef]

116. Merriman, K.R.; Daggupati, P.; Srinivasan, R.; Hayhurst, B. Assessment of site-specific agricultural Best Management Practices in the Upper East River watershed, Wisconsin, using a field-scale SWAT model. *J. Great Lakes Res.* **2019**, *45*, 619–641. [CrossRef]

117. Sullivan, T.P.; Gao, Y.; Reimann, T. Nitrate transport in a karst aquifer: Numerical model development and source evaluation. *J. Hydrol.* **2019**, *573*, 432–448. [CrossRef]

118. Chen, H.; Wang, S.; Wang, Y.; Zhu, J. Probabilistic projections of hydrological droughts through convection-permitting climate simulations and multimodel hydrological predictions. *J. Geophys. Res. Atmos.* **2020**, *125*, e2020JD032914. [CrossRef]

119. Al Aamery, N.; Adams, E.; Fox, J.; Husic, A.; Zhu, J.; Gerlitz, M.; Agouridis, C.; Bettel, L. Numerical model development for investigating hydrologic pathways in shallow fluviokarst. *J. Hydrol.* **2021**, *593*, 125844. [CrossRef]

120. Karki, R.; Srivastava, P.; Kalin, L.; Mitra, S.; Singh, S. Assessment of impact in groundwater levels and stream-aquifer interaction due to increased groundwater withdrawal in the lower Apalachicola-Chattahoochee-Flint (ACF) River Basin using MODFLOW. *J. Hydrol. Reg. Stud.* **2021**, *34*, 100802. [CrossRef]

121. Salerno, F.; Tartari, G. A coupled approach of surface hydrological modelling and Wavelet Analysis for understanding the baseflow components of river discharge in karst environments. *J. Hydrol.* **2009**, *376*, 295–306. [CrossRef]

122. Vale, M.; Holman, I.P. Understanding the hydrological functioning of a shallow lake system within a coastal karstic aquifer in Wales, UK. *J. Hydrol.* **2009**, *376*, 285–294. [CrossRef]

123. Tzoraki, O.; Cooper, D.; Kjeldsen, T.; Nikolaidis, N.P.; Gamvroudis, C.; Froebrich, J.; Querner, E.; Gallart, F.; Karalemas, N. Flood generation and classification of a semi-arid intermittent flow watershed: Evrotas river. *Int. J. River Basin Manag.* **2013**, *11*, 77–92. [CrossRef]

124. Palazón, L.; Navas, A. Sediment production of an alpine catchment with SWAT. *Z. Geomorphol* **2013**, *57*, 69–85. [CrossRef]

125. Palazón, L.; Navas, A. Modeling sediment sources and yields in a Pyrenean catchment draining to a large reservoir (Ésera River, Ebro Basin). *J. Soils Sediments* **2014**, *14*, 1612–1625. [CrossRef]

126. Sellami, H.; La Jeunesse, I.; Benabdallah, S.; Baghdadi, N.; Vanclooster, M. Uncertainty analysis in model parameters regionalization: A case study involving the SWAT model in Mediterranean catchments (Southern France). *Hydrol. Earth Syst. Sci.* **2014**, *18*, 2393–2413. [CrossRef]

127. Gamvroudis, C.; Nikolaidis, N.; Tzoraki, O.; Papadoulakis, V.; Karalemas, N. Water and sediment transport modeling of a large temporary river basin in Greece. *Sci. Total Environ.* **2015**, *508*, 354–365. [CrossRef] [PubMed]

128. Malagò, A.; Pagliero, L.; Bouraoui, F.; Franchini, M. Comparing calibrated parameter sets of the SWAT model for the Scandinavian and Iberian peninsulas. *Hydrol. Sci. J.* **2015**, *60*, 949–967. [CrossRef]

129. Sellami, H.; Benabdallah, S.; La Jeunesse, I.; Vanclooster, M. Quantifying hydrological responses of small Mediterranean catchments under climate change projections. *Sci. Total Environ.* **2016**, *543*, 924–936. [CrossRef]
130. Palazón, L.; Navas, A. Case Study: Effect of Climatic Characterization on River Discharge in an Alpine-Prealpine Catchment of the Spanish Pyrenees Using the SWAT Model. *Water* **2016**, *8*, 471. [CrossRef]

131. Vigiak, O.; Malagó, A.; Bouraoui, F.; Vanmaercke, M.; Obreja, F.; Poesen, P.; Habersack, H.; Fehér, J.; Grošelj, S. Modelling sediment fluxes in the Danube River Basin with SWAT. *Sci. Total Environ.* **2017**, *599–600*, 992–1012. [CrossRef]

132. Efthimiou, N. Hydrological simulation using the SWAT model: The case of Kalamas River catchment. *J. Appl. Water Eng. Res* **2018**, *6*, 210–227. [CrossRef]

133. Martínez-Salvador, A.; Conesa-García, C. Suitability of the SWAT Model for Simulating Water Discharge and Sediment Load in a Karst Watershed of the Semiarid Mediterranean Basin. *Water Resour. Manag.* **2020**, *34*, 785–802. [CrossRef]

134. Senent-Aparicio, J.; Alcalá, F.J.; Liu, S.; Jimeno-Sáez, P. Coupling SWAT Model and CMB Method for Modeling of High-Permeability Bedrock Basins Receiving Interbasin Groundwater Flow. *Water* **2020**, *12*, 657. [CrossRef]

135. Busico, G.; Ntona, M.M.; Carvalho, S.C.P.; Patrikaki, O.; Voudouris, K.; Kazakis, N. Simulating Future Groundwater Recharge in Coastal and Inland Catchments. *Water Resour. Manag.* **2021**, *35*, 3617–3632. [CrossRef]

136. Sánchez-Gómez, A.; Martínez-Pérez, S.; Pérez-Chavero, F.M.; Molina-Navarro, E. Optimization of a SWAT model by incorporating geological information through calibration strategies. *Optim. Eng.* **2022**, *23*, 2203–2233. [CrossRef]

137. Jiang, R.; Li, Y.; Wang, Q.; Kuramochi, K.; Hayakawa, A.; Woli, K.P.; Hatano, R. Modeling the water balance processes for understanding the components of river discharge in a non-conservative watershed. *Trans. ASABE* **2011**, *54*, 2171–2180. [CrossRef]

138. Tian, Y.; Wang, S.; Bai, X.; Luo, G.; Xu, Y. Trade-offs among ecosystem services in a typical Karst watershed, SW China. *Sci. Total Environ.* **2016**, *566*, 1297–1308. [CrossRef] [PubMed]

139. Bucak, T.; Trolle, D.; Andersen, H.E.; Thodsen, H.; Erdoğan, S.; Levi, E.E.; Filiz, N.; Jeppesen, E.; Beklioğlu, M. Future water availability in the largest freshwater Mediterranean lake is at great risk as evidenced from simulations with the SWAT model. *Sci. Total Environ.* **2017**, *581–582*, 413–425. [CrossRef]

140. Hou, W.; Gao, J. Simulating runoff generation and its spatial correlation with environmental factors in Sancha River Basin: The southern source of the Wujiang River. *J. Geogr. Sci.* **2019**, *29*, 432–448. [CrossRef]

141. Jakada, H.; Chen, Z. An approach to runoff modelling in small karst watersheds using the SWAT model. *Arab J. Geosci.* **2020**, *13*, 318. [CrossRef]

142. Mo, C.; Zhang, M.; Ruan, Y.; Qin, J.; Wang, Y.; Sun, G.; Xing, Z. Accuracy Analysis of IMERG Satellite Rainfall Data and Its Application in Long-term Runoff Simulation. *Water* **2020**, *12*, 2177. [CrossRef]

143. Hou, W.; Gao, J.; Wu, S. Quantitative Analysis of the Influencing Factors and Their Interactions in Runoff Generation in a Karst Basin of Southwestern China. *Water* **2020**, *12*, 2898. [CrossRef]

144. Gao, J.; Jiang, Y.; Anker, Y. Contribution analysis on spatial tradeoff/synergy of Karst soil conservation and water retention for various geomorphological types: Geographical detector application. *Ecol. Indic.* **2021**, *125*, 107470. [CrossRef]

145. Jiang, Y.; Gao, J.; Yang, L.; Wu, S.; Dai, E. The interactive effects of elevation, precipitation and lithology on karst rainfall and runoff erosivity. *CATENA* **2021**, *207*, 105588. [CrossRef]

146. Chang, W.; Li, W.; Ma, H.; Wang, D.; Bandala, E.R.; Yu, Y.; Rodrigo-Comino, J. An integrated approach for shaping drought characteristics at the watershed scale. *J. Hydro.* **2022**, *604*, 127248. [CrossRef]

147. Zhang, J.; Zhang, P.; Song, Y. Comparative Water Environment Simulation Study of Two Typical Models with BMPs in a Karst Basin. *Agriculture* **2022**, *12*, 69. [CrossRef]

148. Mo, C.; Chen, X.; Lei, X.; Wang, Y.; Ruan, Y.; Lai, S.; Xing, Z. Evaluation of Hydrological Simulation in a Karst Basin with Different Calibration Methods and Rainfall Inputs. *Atmosphere* **2022**, *13*, 844. [CrossRef]

149. Yuan, J.; Li, R.; Huang, K. Driving factors of the variation of ecosystem service and the trade-off and synergistic relationships in typical karst basin. *Ecol. Indic.* **2022**, *142*, 109253. [CrossRef]

150. Zettam, A.; Taleb, A.; Sauvage, S.; Boithias, L.; Belaidi, N.; Sánchez-Pérez, J.M. Modelling Hydrology and Sediment Transport in a Semi-Arid and Anthropized Catchment Using the SWAT Model: The Case of the Tafna River (Northwest Algeria). *Water* **2017**, *9*, 216. [CrossRef]

151. Zaibak, I.; Meddi, M. Simulating streamflow in the Cheliff basin of west northern Algeria using the SWAT model. *J. Earth Syst. Sci.* **2022**, *131*, 25. [CrossRef]

152. Afinowicz, J.D.; Munster, C.L.; Wilcox, B.P. Modeling effects of brush management on the rangeland water budget: Edwards Plateau, Texas. *J. Am. Water Resour. Assoc.* **2005**, *41*, 181–193. [CrossRef]

153. Baffaut, C.; Benson, V. Modeling flow and pollutant transport in a karst watershed with SWAT. *Trans. ASABE* **2009**, *52*, 469–479. [CrossRef]

154. Yactayo, G.A. Modification of the SWAT Model to Simulate Hydrologic Processes in a Karst-Influenced Watershed. Master's Thesis, Virginia Tech, Blacksburg, VA, USA, 2009.

155. Palanisamy, B.; Workman, S.R. Hydrologic modeling of flow through sinkholes located in streambeds of cane run stream, Kentucky. *J. Hydrol. Eng.* **2014**, *20*, 04014066. [CrossRef]

156. Zhou, Y.; Zhao, L.; Cao, J.; Wang, Y. Using an Improved SWAT Model to Simulate Karst Sinkholes: A Case Study in Southwest China. *Front. Environ. Sci.* **2022**, *10*, 950098. [CrossRef]
130. Palazón, L.; Navas, A. Case Study: Effect of Climatic Characterization on River Discharge in an Alpine-Prealpine Catchment of the Spanish Pyrenees Using the SWAT Model. *Water* **2016**, *8*, 471. [CrossRef]

131. Vigiak, O.; Malagó, A.; Bouraoui, F.; Vanmaercke, M.; Obreja, F.; Poesen, P.; Habersack, H.; Fehér, J.; Grošelj, S. Modelling sediment fluxes in the Danube River Basin with SWAT. *Sci. Total Environ.* **2017**, *599–600*, 992–1012. [CrossRef]

132. Efthimiou, N. Hydrological simulation using the SWAT model: The case of Kalamas River catchment. *J. Appl. Water Eng. Res* **2018**, *6*, 210–227. [CrossRef]

133. Martínez-Salvador, A.; Conesa-García, C. Suitability of the SWAT Model for Simulating Water Discharge and Sediment Load in a Karst Watershed of the Semiarid Mediterranean Basin. *Water Resour. Manag.* **2020**, *34*, 785–802. [CrossRef]

134. Senent-Aparicio, J.; Alcalá, F.J.; Liu, S.; Jimeno-Sáez, P. Coupling SWAT Model and CMB Method for Modeling of High-Permeability Bedrock Basins Receiving Interbasin Groundwater Flow. *Water* **2020**, *12*, 657. [CrossRef]

135. Busico, G.; Ntona, M.M.; Carvalho, S.C.P.; Patrikaki, O.; Voudouris, K.; Kazakis, N. Simulating Future Groundwater Recharge in Coastal and Inland Catchments. *Water Resour. Manag.* **2021**, *35*, 3617–3632. [CrossRef]

136. Sánchez-Gómez, A.; Martínez-Pérez, S.; Pérez-Chavero, F.M.; Molina-Navarro, E. Optimization of a SWAT model by incorporating geological information through calibration strategies. *Optim. Eng.* **2022**, *23*, 2203–2233. [CrossRef]

137. Jiang, R.; Li, Y.; Wang, Q.; Kuramochi, K.; Hayakawa, A.; Woli, K.P.; Hatano, R. Modeling the water balance processes for understanding the components of river discharge in a non-conservative watershed. *Trans. ASABE* **2011**, *54*, 2171–2180. [CrossRef]

138. Tian, Y.; Wang, S.; Bai, X.; Luo, G.; Xu, Y. Trade-offs among ecosystem services in a typical Karst watershed, SW China. *Sci. Total Environ.* **2016**, *566*, 1297–1308. [CrossRef] [PubMed]

139. Bucak, T.; Trolle, D.; Andersen, H.E.; Thodsen, H.; Erdoğan, S.; Levi, E.E.; Filiz, N.; Jeppesen, E.; Beklioğlu, M. Future water availability in the largest freshwater Mediterranean lake is at great risk as evidenced from simulations with the SWAT model. *Sci. Total Environ.* **2017**, *581–582*, 413–425. [CrossRef]

140. Hou, W.; Gao, J. Simulating runoff generation and its spatial correlation with environmental factors in Sancha River Basin: The southern source of the Wujiang River. *J. Geogr. Sci.* **2019**, *29*, 432–448. [CrossRef]

141. Jakada, H.; Chen, Z. An approach to runoff modelling in small karst watersheds using the SWAT model. *Arab J. Geosci.* **2020**, *13*, 318. [CrossRef]

142. Mo, C.; Zhang, M.; Ruan, Y.; Qin, J.; Wang, Y.; Sun, G.; Xing, Z. Accuracy Analysis of IMERG Satellite Rainfall Data and Its Application in Long-term Runoff Simulation. *Water* **2020**, *12*, 2177. [CrossRef]

143. Hou, W.; Gao, J.; Wu, S. Quantitative Analysis of the Influencing Factors and Their Interactions in Runoff Generation in a Karst Basin of Southwestern China. *Water* **2020**, *12*, 2898. [CrossRef]

144. Gao, J.; Jiang, Y.; Anker, Y. Contribution analysis on spatial tradeoff/synergy of Karst soil conservation and water retention for various geomorphological types: Geographical detector application. *Ecol. Indic.* **2021**, *125*, 107470. [CrossRef]

145. Jiang, Y.; Gao, J.; Yang, L.; Wu, S.; Dai, E. The interactive effects of elevation, precipitation and lithology on karst rainfall and runoff erosivity. *CATENA* **2021**, *207*, 105588. [CrossRef]

146. Chang, W.; Li, W.; Ma, H.; Wang, D.; Bandala, E.R.; Yu, Y.; Rodrigo-Comino, J. An integrated approach for shaping drought characteristics at the watershed scale. *J. Hydro.* **2022**, *604*, 127248. [CrossRef]

147. Zhang, J.; Zhang, P.; Song, Y. Comparative Water Environment Simulation Study of Two Typical Models with BMPs in a Karst Basin. *Agriculture* **2022**, *12*, 69. [CrossRef]

148. Mo, C.; Chen, X.; Lei, X.; Wang, Y.; Ruan, Y.; Lai, S.; Xing, Z. Evaluation of Hydrological Simulation in a Karst Basin with Different Calibration Methods and Rainfall Inputs. *Atmosphere* **2022**, *13*, 844. [CrossRef]

149. Yuan, J.; Li, R.; Huang, K. Driving factors of the variation of ecosystem service and the trade-off and synergistic relationships in typical karst basin. *Ecol. Indic.* **2022**, *142*, 109253. [CrossRef]

150. Zettam, A.; Taleb, A.; Sauvage, S.; Boithias, L.; Belaidi, N.; Sánchez-Pérez, J.M. Modelling Hydrology and Sediment Transport in a Semi-Arid and Anthropized Catchment Using the SWAT Model: The Case of the Tafna River (Northwest Algeria). *Water* **2017**, *9*, 216. [CrossRef]

151. Zaibak, I.; Meddi, M. Simulating streamflow in the Cheliff basin of west northern Algeria using the SWAT model. *J. Earth Syst. Sci.* **2022**, *131*, 25. [CrossRef]

152. Afinowicz, J.D.; Munster, C.L.; Wilcox, B.P. Modeling effects of brush management on the rangeland water budget: Edwards Plateau, Texas. *J. Am. Water Resour. Assoc.* **2005**, *41*, 181–193. [CrossRef]

153. Baffaut, C.; Benson, V. Modeling flow and pollutant transport in a karst watershed with SWAT. *Trans. ASABE* **2009**, *52*, 469–479. [CrossRef]

154. Yactayo, G.A. Modification of the SWAT Model to Simulate Hydrologic Processes in a Karst-Influenced Watershed. Master's Thesis, Virginia Tech, Blacksburg, VA, USA, 2009.

155. Palanisamy, B.; Workman, S.R. Hydrologic modeling of flow through sinkholes located in streambeds of cane run stream, Kentucky. *J. Hydrol. Eng.* **2014**, *20*, 04014066. [CrossRef]

156. Zhou, Y.; Zhao, L.; Cao, J.; Wang, Y. Using an Improved SWAT Model to Simulate Karst Sinkholes: A Case Study in Southwest China. *Front. Environ. Sci.* **2022**, *10*, 950098. [CrossRef]