Unintentional Recharge of Aquifers from Dams and Dykes in Spain. A GIS-Based Approach to Determine a Fractional Volume

1. Introduction

Aquifers can be anthropically recharged intentionally or unintentionally. Unintentional aquifer recharge can occur via diverse processes such as percolation through dams, irrigation returns, leakage from sewer and water supply systems, and wastewater disposal, among others. Intentional aquifer recharge is conducted through a variety of methods that entail infiltrating water through permeable surfaces, the vadose zone or directly into aquifer through injection systems [1]. This sort of intentional recharge is often called managed aquifer recharge (MAR) or artificial recharge.

The amount of anthropogenic recharge can be a considerable component of the hydrological balance and make up a large percentage of the annual recharge of a region, or even a country.

In Spain, a country where water stress is tangible across many surface and groundwater bodies, regional water budgets are of paramount importance to properly manage and allocate water resources. In this context, there have been various attempts to estimate the contribution from dams of various sizes (also referred here as transverse structures and dykes) and MAR systems to groundwater recharge.

Apart from a few large-scale operational schemes located in the regions of Castile and Leon and Catalonia [1], most of the large-scale MAR implementations in Spain have been built with an experimental character. According to the DINAMAR project [1,2], the total volume of water infiltrated in Spain through MAR (considering 32 schemes) and a preliminary inventory of nearly 800 dykes is 380 Mm³/year [1,3,4,5]. During the 2020 CONAMA conference, the Spanish water authorities estimated this figure at about 500 Mm³/year [3]. Other references suggest a lower figure in the order of 50 Mm³/year [6]. Other sources calculated the recharge from MAR and dykes and dams in urban and forested areas, aimed to flood mitigation and groundwater recharge, at 350 Mm³/year [7] and 380 Mm³/year [1,2,7]. The infiltration through the bottom of large dams, which in the Spanish territory exceed 1,400 [8,9,10], is not included in these estimations, and could result in a recharge of about 800 Mm³/year [3].

Given the relatively small recharge volume that MAR represents compared to other sources of unintentional recharge [1,3,7], the mismatch in the figure above is attributed predominantly to the difference in the number of transversal structures considered, and the methodologies employed.

This article has two main objectives: (1) to present the inventory of transverse structures created by the Spanish Ministry for the Ecological Transition and the Demographic Challenge (MITERD); and (2), to provide a calculation of the (un)managed aquifer recharge caused by cross-cutting river structures, such as dykes and small dams over the peninsular Spanish territory.

To the first end, the transversal structure inventory conducted by the MITERD is presented and described [11,12,13], exploring the methodology employed to construct the database, the information provided for each transversal structure and the level of accuracy involved in the determination of such information. To achieve the second objective, a GIS-based approach relying in various information layers is applied.

This analysis does not imply that recharge from dykes and dams is considered an intentional mechanism to increase groundwater storage, and therefore, can be regarded as MAR, but rather such additional recharge constitutes a secondary effect [14].

The results from this research should help disseminating the result from the Spanish Inventory of water course structures; contribute to closing knowledge gaps concerning anthropogenic recharge in Spain; providing a methodology to estimate recharge from a large number of structures (dykes and small dams) at the regional level; and ultimately, provide figures to better allocate and manage water resources.

2. Methodology

The methodology of this article has been subdivided into two sections. In the first one the background and some of the information pertaining to the transverse structure inventory is provided. The second part deals with the estimation of unintentional recharge.

2.1. Inventory of transverse structures

2.1.1. Theory of river channel transverse structures

Rivers channel works are frequently used in hydrological and forestry restoration efforts [11] to reduce the effect of torrential precipitation by preventing floods [12,15]. These works are also aimed at adjusting and controlling the effects of water flows. For instance, they can act on the river bed and banks, preventing the formation of flows with high concentrations of solids and high density, caused by erosion and sedimentation of materials [15,16,17].

River channel works can be either longitudinal or transverse. Longitudinal works are built in torrential watercourses for lateral erosion defence, landslide containment, and flood defence. They have a very limited effect on water retention and infiltration into aquifers [18]. Transverse structures are the most suitable solution to deal with sediment transport, erosion, and flooding processes in torrential watercourses with massive transport of materials, erosion of banks, and destabilisation of adjacent massifs [15,17,18]. These structures are also the best choice when the objective is to retain water.

The technical literature separates "dyke” from "weir", or small dam. The former term is employed for cross-cutting structures, designed predominantly to respond to issues associated with the retention of solid flows and the consolidation of slopes and beds of water channels, such as torrents, ravines and wadis [18,19,20]. The latter term, weir, or small dam, involves impounding water for diverse purposes.

As long as the reservoir basin that they originate stays unlanded, the dam effect causes the dammed water to slow down the rate at which sediment arrives, resulting in the stagnation of coarser sediments and a decrease in the solid proportion of the discharge flow. This has an essential effect on the type of clogging that accumulates at the bottom and constitutes an impediment to water infiltration [21]. The deposits that occur form a siltation that raises the bed of the channel until it reaches the compensation slope, which is always lower than the natural channel´s slope [22]. The compensation slope is formed when the sieve size composition of the bed and the stream flow is equalised, compensating for the solid volumes entering and leaving the siltation [20].

Raising the riverbed leads to a lower slope of the embankment, which reduces its velocity and, consequently, its capacity for dragging and erosion. The embankment wedge attached to the construction site also consolidates [23].

Dykes can be further categorised into consolidation and retention dykes. Consolidation dykes (Figure 1) serve the purpose of decreasing river channel erosion and stabilising the slope. Retention dykes (Figure 2) are built to retain solid or liquid flows. In the first case, i.e., to retain solid flows, they are known as closed or semi-hollow dykes, and, when aimed at retaining water, they are either detention or recharge dykes [24]. The dykes for retention play a special role in lamination and (un)intentional aquifer recharge [20,25,26].

Longitudinal works are used in torrential watercourses as lateral erosion defence works, landslide containment works, and flood defence works. They have a scant or null effect on water infiltration into aquifers [18,25].

The main difficulties of this type of work arise from the analysis and forecasting of phenomena related to sediment transport, its influence on the stability of the modified channel, both in terms of the extent of this influence [17], its evolution over time [18,26], and finally, the "anthropic" seepage variation [1,4,27,28].

In addition to these processes, dams can also contribute to infiltrating water into aquifers by increasing water retention times, working as infiltration basins [15,16,17]. The rate at which recharge proceeds is conditioned by factors such as the wet perimeter, the water stage, and the riverbed grain size distribution. As long as the dam reservoir preserves its volumetric capacity, the dam decreases the water sediment load, slows down the sediment delivery rate, and results in the deposition of fine sediments, a process commonly known as siltation, resulting in a low-permeability layer or clogging layer that decreases infiltration rates [21]. The low-permeability layer raises the channel bed up to the compensation slope [22], which forms when the grain size distribution of the river bed and the stream are equalised and there is a balance between the solid fraction entering and leaving the siltation layer [20,22].

2.1.2. inventory of transverse structures. Background and methology

In 2019, the Spanish Ministry for the Ecological Transition and the Demographic Challenge (MITERD) requested Tragsatec to conduct and digitalise the first inventory of dykes and transverse structures in watercourses in the peninsular Spanish territory (excluding the islands) [11,27,28].

Most of these transversal structures were constructed by former institutions for forest protection and agricultural development (ICONA and IRYDA, respectively) in the 50s and 60s [31], leaving aside the big dams. The documents concerning their construction were scattered across numerous libraries. The construction of the inventory relied on information from some of these documents; the photointerpretation of orthophotos from the National Aerial Orthophotography Plan (PNOA), satellite photography from Google Earth, and panoramic images at street level from Google Street View; and field visits, to corroborate ground truth satellite and remote information.

After three years of intense work, the inventory was consolidated, and the resulting number of transverse structures, i.e., 27,680 units, proved that a preliminary estimation of 8,000 [13] underestimated the reality [27,28].

The inventory and the related GIS information have been publicly available online on different websites (Table 1).

An interactive tool at MITERD's geoportal shows the transverse structures' location (Link provided in Table 1). There are links to download the actual inventory in shape file and KMZ formats and documentation about the transverse structures' attributes

2.2. Estimation of recharge from transverse structures.

2.2.1. General approach

The calculation of groundwater recharge from transversal structures is based on the following information sources: UTM grids to be sampled [31], the national lithological map [33], Spain's surface permeability map [33], design information from the new inventory of transversal structures in Spain [11,23,34,35], and other layers that contain relevant information, such as the DINA-MAR hydro-geoportal [36,37].

Equation 1 is the general formula employed to compute infiltration from transverse structures.

$$\sum _{i=1}^n{WA}_i{K}_{vi}{C}_{cfi}{GWL}_{cfi}$$

(1)

Where WA_i is the wetted area, K_vi is the vertical hydraulic conductivity, C_cfi is the clogging correction factor, and GWL_cfi is the groundwater level correction factor of the i^th transversal structure. The next subsection describes the methodologies used to estimate each variable in Equation 1.

2.2.2. Wetted area (AW_i)

The long-term average wetted area from transverse structures has been estimated through design characteristics obtained from the inventory of transverse structures.

The geodatabase [34] includes a layer of the National Flood Zone Mapping System (in Spanish, SNCZI-IPE) [38], with additional cross-cutting obstacles, such as minor vertical jumps, piped crossings, and crossings over faces.

Two methods to estimate the AW_i variable have been tested.

The first is the formula by Mozzi et al. [39], devised to estimate the recharge from check dams in India. According to these authors, the wetted area can be computed through Equation 2 and Equation 3.

$$A_{surf}=\frac{W_s{H}_s}{\tan \theta }$$

(2)

$$WA=\left[1+\frac{2H_s}{\cos \alpha }\left(1-\sin \alpha \right)\right]A_{surf}$$

(3)

Where A_surf is the water surface, W_s is the structure width, H_s is the structure's height, θ is the stream gradient, and α is the river bank slope. All the necessary input data except the slope were available in the inventory [11]. The stream gradient for each construction was calculated from the 1-meter resolution DTM PNOA-IGN [40], as the average of the upstream and downstream slopes (Figure 3).

Different examples were conducted to check the reliability of the calculation.

An alternative method to approximate WA is proposed. This method considers a simplified equation (4) based on the backwater length and width.

$$WA=\frac{L_{bw}{W}_{bw}}{2}$$

(4)

If a transverse structure misses this data, the average wet area from well-defined structures has been assigned. i.e., 6,914.87 m².

An obstacle with a jump height of more than 10 metres is considered a large dam and included in a specific inventory [9,27]. Dimensions above 3 m height cause sufficiently large backwaters [28]. Both dimensions are not always available to apply the calculation. Therefore, other data had to be taken into account to check the magnitude of the effect of each weir in concrete and the average width, and depth of the channel before the backwater created by the obstacle.

Only 4,171 transverse structures (out of 27,685) count with complete design and backwater dimensions data that were confirmed in the field.

3.3. Estimation of Kv

As groundwater recharge through surface structures is mainly driven by gravity, vertical hydraulic conductivity (K_v) is the required parameter to estimate the recharge from a surface water body into a shallow aquifer.

The methodology to obtain values of K_v at every transverse structure consists of the following steps: 1) assign the predominant lithology to every transverse structure; 2) use this lithology to designate permeability values; and 3) use standard anisotropy coefficients to estimate the vertical hydraulic conductivity.

First, the predominant lithology at every transverse structure has been determined as the lithology from the Spanish lithostratigraphic map with the highest cover percentage within a buffer zone [32] (Figure 4). The buffer zone was defined using the length of the wet perimeter. In case this value was missing, the average was employed.

Subsequently, the assigned lithologies have been translated into horizontal permeability values using Spain's surface permeability map at a scale of 1:200,000 [34]. This map qualitatively displays lithostratigraphic units' permeability, establishing five categories: very high, high, medium, low, and very low permeability. In addition, the map considers seven large groups of lithologies: carbonate, granular, quaternary granular, volcanic, meta-granular, igneous, and evaporitic.

For each permeability category, the average of the range of typical hydraulic conductivity values reported by Terzaghi and Peck [41], and Bendixen et al. [42] were assigned, resulting in the categorisation in Table 2.

Finally, to obtain vertical hydraulic conductivities, the anisotropy coefficient, which relates K_h and K_v, was used [1,43].

The anisotropy coefficients employed (Table 3) were taken from the memory of the Spanish permeability map [33], Bouwer and Rice [43], Chapuis and Ernest [44], and Sridhar et al. [45].

3.5. Groundwater level correction factor (GWL_cfi)

According to Bouwer and Rice [43], the infiltration rate below a structure aimed to infiltrate variates as a function of the depth until the groundwater level underneath. Consequently, an additional correction factor considering the groundwater level depth has been applied. The correction factor selected is 1.05 if the groundwater table depth is > 25 meters deep. This factor is purely bibliographic [43,47].

3.6. Clogging Correction Factor (C_cfi)

Clogging layers often decrease the vertical hydraulic conductivity of river beds and can control infiltration rates completely, sometimes acting as a bottleneck [40,43] that renders the underlying hydrogeological configuration irrelevant [46]. Consequently, the infiltration from transverse structures in riverbeds can depend considerably on the hydraulic impedance of the clogging layer [47], often composed of structure-less loams, sands and clays.

The clogging correction factor varies considerably in the hydrogeological literature. There are average values of about 0.75 that consider different types of clogging (e.g., physical and gas clogging) and are based on fieldwork conducted in sites with different lithologies, including alluvial and mudstone formations [48,49]. There are also more conservative correction factors of 0.5 [49,50] that account for different types of clogging, such as leakage, efficiency reduction of the structures, meteorological factors, etc. [50]. Clogging correction factors of about 0.8 have also been estimated for specific MAR systems (e.g., Aquifer Storage and Recovery (ASR) [51]). A study on check dams in SE Spain found a factor of 0.73-0.78 [52].

To consider the diverse types of lithologies involved in this research, the Dutch protocol 11001 approach [53] has been utilised. This protocol (Equation 5) was initially designed for injection wells and boreholes. Thus, it includes factors such as the number of full operation hours, the specific obstruction rate, and the modified fouling index (MFI).

$$v_{inf}=\mathrm{1,000}{\left(\frac{K}{150}\right)}^{0.6}\sqrt{\frac{v_{verstopping}}{2MFI{u}_{eq}}}$$

(5)

Where v_inf is the infiltration velocity (m/hour), K is the original material hydraulic conductivity (m/day); v_verstopping is the standard value for each lithology (m/j), which is often a constant value of 0.1; MFI is the modified fouling index (s/L2), which habitually is assigned a value of 2 [54]; and u_eq is the number of total operation hours (h/j) [54, appendix 4].

The resulting clogging factors for each lithology are presented in Table 4. These factors depend on the predominant sediments in the clogging layer (gravel, grit, loam, sand and clay), structure, compaction degree, and gas entrainment (v_verstopping). This index refers specifically for infiltration boreholes and wells, not for surface recharge devices. From now on, special attention is paid on extreme values: a conservative correction factors of 0.5, usual for volcanic materials, and another, rather lax, 0.8, more frequently cited in the references [51,52].

4. Results and Discussion

4.1. Analysis of the database.

The public database consists of 27,680 structures represented as points in a shape file (shp) (Figure 5). Each point has 33 attributes [11,29]. An additional database on longitudinal structures along river channels is available on the Internet [15,27]. It has not been discussed, as it is rather irrelevant to the main research question.

The document published by Tragsatec entitled "Data model and storage of information contained in the inventory of transverse barriers and longitudinal defence works" [28], describes the design and content of the database in full detail. The following is a description of its structure, focusing on the general interest attributes and those relevant to the calculations conducted.

The "ObsTrans" layer contains a unique alphanumeric identifier for each structure (autoUID), and a series of identifying codes, which are shown in Table 4.

The type of obstacle is specified in column "TIPO_INFR" through the abbreviations and in the proportions shown in Figure 6. Vertical jumps, crossing with roads, and step over a wall are the most predominant.

Given the uncertainty in the methods employed to classify the structures, namely, satellite imagery, each structure in the database has a "CIERTO" attribute, which determines the degree of certainty in the type of structure assigned (Figure 7).

The "REVISIONES" column specifies whether a structure has been visited in the field and provides the visit date (FECHA REVISION).

The "TRANSGENERALES" table includes two identifiers for each structure, an alphanumeric code (AutoUID), and a code based on the type of structure and water body (ID_CLAVE). This table also includes the code of the surface water body (CODMAS), the code of the hydromorphological section (TRAMO_COD), the Autonomous Community (NOM_CCAA), the province (NOM_PROV), the municipality (NOM_MUNICI), coordinates, and name of the river (NOM_RIO) structure location. The river basin where the structure is located is specified in column "DEMARCACIO" and follows the "EUUMCODE" nomenclature. The number of transverse structures within each major Spanish River basin is presented in Figure 8. The majority of structures are located in central and northern eastern Spain, in basins such as the Duero, Júcar, Ebro, Tajo and Miño-Sil.

The attribute "EXIS_PASO" indicates if the structure has a fish ladder (yes/no), or a lateral diversion (deviation) (Figure 9).

The following columns in the database specify the origin of the information (ORIGEN), whether it was present in previous databases (BBDD_CONF) [27], whether it has been modified in any way (MOD_IMPRESS), and the name of the structure, if any.

The "TRANSTECNICOS" table includes identification data that assigns a main use (USO_P) among 22 options (Figure 10). The main uses are road crossing (10,317 structures), and irrigation-related (3,664 structures).

In certain cases, secondary uses of the infrastructure, conservation quality, construction year, concession date (if applicable), and type of passage are also reported. The construction date is important to estimate the clogging accumulation (in general, the oldest, the more clogged the river bed) and the structural fatigue (the oldest the structure, the greater the amount of stored sediment and the greater the confining pressure).

Another table, "EFECTO BARRERA" (barrier effect), contains information about the design and construction of the transverse structures. The main information is shown in Table 5. Approximately 20% of entries (4,171 structures out of 27,685) contain all the parameters in this table (Figure 11). The "EFECTO BARRERA" table includes most of the constructive data of particular relevance for the target of the study plus environmental aspects.

This information is complemented with attributes concerning the existence of a side channel in the structure (yes, no or unknown), grating, crossing device, the difficulty of ascent at the foot of the structure, call effect, whether there is significant turbulence, the roughness of the face, slopes or changes of gradient, absence of a clear gradient of velocity in the dammed area, passage through turbines, mills or falls of more than 10 m, accessibility of the channel to the downstream crossing device, call effect on the ascent, lift well at the foot of the structure, and accessibility of the crossing device to the upstream channel, among others.

Very relevant to this research is the attribute "COLMATACION", indicating if the structure is clogged (yes or no).

4.3. Estimation of unintentional recharge from transverse structures and validation

Two methods have been used to estimate the wetted area at each transverse structure.

The simplified formula proposed results in a mean wetted area of 6,914.87 m², which is closer to field observations and a total recharge that ranges between 812.5 Mm³/year (applying the conservative clogging correction factor x0.5), and 1,218.8 Mm³/year (using the less restrictive clogging correction factor x0.8).

The formula by Mozzi et al. [39] yields a mean wetted area of 19,747.4 m², and a total unintentional recharge of 2,716.8 and 1,901.8 Mm³/year, applying respectively the clogging correction factors.

These values are in the same order of magnitude as previous estimations. Nonetheless, the calculated volumes considerably exceed previous estimations by the Water General Directorate of the Spanish Ministry (500 Mm³/year, according to own estimations [3]). It is also significantly above the DINAMAR project's figure, 380 Mm³/year [4], considering the combined recharge of operative MAR systems and a limited inventory of about 800 transverse structures [1,2,4]. Therefore, the value obtained using the conservative clogging factor (i.e., 812.5 Mm³/year) is the closest to previous estimations.

These values are also contrasted with the SIMulación Precipitación-Aportación (SIMPA) quasi-distributed hydrological model (SIMulation Precipitation-Input), developed by CEDEX [29,31] for the Ministry of Ecological Transition. Its latest version (the period between hydrological years 1980/81-2017/18) has gridded Spain into approximately 2,000,000 square cells of 500 m on each side, and computes the hydrological balance, resulting in a series of published maps displaying precipitation (PRE), potential evapotranspiration (ETP), actual evapotranspiration (ETR), soil moisture (HUM), infiltration (INF), surface runoff (ASP), groundwater runoff (ASB), and runoff (AES), or total input (APN) [29].

According to the SIMPA model, the river network recharges more than 800 Mm³/year. This model considers direct infiltration from watercourses, and includes infiltration from large transverse structures, which do not include the structures in the inventory presented (see validation chapter).

The figure for the total recharge in Spain from river courses is 1,946 Mm³ [31]. In this respect, the percentage of infiltration from transverse structures ranges between 41.75% and 62.63%, not so far from the authors' 0.5 and 0.8 correction factors. An identical C_cf between 0.50 and 0.80 has also been found in a reference [39] for check dams in India. It is worth noticing that errors of different signs cancel each other out.

Values are also compared to the HEC-RAS modeling of watercourses for the Southeast of the Iberian Peninsula using HEC-RAS modeling on watercourses [52]. The variations in the flow rate through river stretches with an influent behaviour provide information about the infiltration rate from both the bottom and the banks. The model performance was therefore validated against the observed measurements in the semiarid conditions of the SE of Spain (see validation chapter). The HEC-RAS model at the regional scale for SE Spain [52, pages 991-993] found a C_cf between 0.73 and 0.78 for areas around check-dams, which are not far from the estimations through SIMPA, and in the present research.

5.1. Sensitivity analysis

Changes in the parameters considered have relatively little impact on the infiltration rates.

Calculations using the simplified methodology proposed by the authors reflect mean values for the wetted area, which is a factor in close correlation with precipitation, runoff and meteorological variables.

Determining the vertical permeability for each lithology below each transverse structure relies upon methods of general acceptance. Hydrogeological parameters variation is far from being a heavyweight factor for the infiltrated volume, due to the basic statistical rule that defects on one side, counterbalance excesses on the other. Permeability ranges are a classical figure found in the literature since at least 1967 [41], with slight or no modifications since then [42]. Consequently, K_v has a moderate influence on infiltration volumes, and thus, it is not a critical variable.

The anisotropy coefficient is a parameter with wide ranges of variation in the references, having even three orders of magnitude, e.g. for volcanic rocks. It is dependent on the lithologies and relies on the homogeneity and isotropy of the aquifer.

Algebra mapping operations conducted in ARC/GIS are secure to guarantee the certainty of the results.

The accumulation of clogging at the bottom of the wetted area is the factor with the greatest influence on the final weight. Corrections factors met in the references have been obtained, in general, for local or regional conditions.

Mozzi et al. [39] discussed the impact of clogging layers. They mention that as impedance reaches 100 days (30 cm of silt with a hydraulic conductivity of 0.003 m/day), infiltration is reduced between 20%-50% for check dams in India. This range is consistent with the results of the proposed simplified formula.

Regional data may represent an important limitation when extrapolating results to the scale of one country. Future studies on correction factors determined by the composition and development of clogging processes under different environmental conditions would be necessary.

Other factors, such as the catchment area surface, and runoff coefficient, depending on the geometry of the area, have also been reported as low-impact variables [39].

After different tests varying the correction factor in the database, the main conclusion is that C_cf is the parameter more sensitive to the final result. Any reasonable variation of the input parameters hardly modifies the computed value 812.52 Mm³/year, and 1,218.79 Mm³/year, as the final result for the infiltrated volume from dykes and small dams in peninsular Spain.

Limitations and factors to consider in future estimations

• The estimation could be enhanced considering the following factors:
• The meteorological variability can modify the wetted perimeters of the transversal structures.
• Considering the geometrical configuration of the wetted area through using the slopes around the structures to constrain the wetted area.
• The clogging factor could be further improved by considering the year of construction of the structures, assigning a higher value to older structures.
• The mean water table, which can constrain recharge.
• The moisture conditions of the basin
• Use of some cartographic information from previous works by TRAGSA, in particular, the "hydro-geoportal" [37,38], which has detailed information on the vertical permeability (K_v) of the aquifers across Spain.

6. Conclusions

Sites across the Mediterranean are facing similar water stress problems, according to standardised groundwater level indexes for the Mediterranean region [55]. Water resources are becoming increasingly scarcer due to drought and water demand, particularly for irrigation.

MAR is at the forefront of two diverging stresses in the Mediterranean rural context. On the one hand, and following [56], groundwater resources are increasingly scarcer under unsustainable abstractions and drought. According to these authors, groundwater is key to sustaining irrigation where the rural exodus to urban centres poses the threat of declining populations. Within this context, MAR can play an important role in this issue by taking advantage of untapped water resources and buffering variable precipitation and stream flows.

Transverse structures result in considerable and, in most cases, unintentional recharge of aquifers. A new GIS-based calculation has been presented for the Spanish peninsula territory, and the recharged volume ranges between 812.5 - 1,218.8 Mm³/year. The lower range value is more realistic, as it matches the water balances for the Spanish territory according to PHCs, the third stage, CEDEX calculations, which estimated 500 Mm³/year [30], and regional modeling calculus [52].

The volume infiltrated into the aquifers from transverse structures (unintentionally and excluding large dams [9]) exceeds by more than four times the volume from MAR systems in Spain (about 200 Mm³/year [1]).

Additional meteorological and clogging data, and specific corrections, should be included in the future to refine the initial estimation.

Patents. Not applicable.