Groundwater recharge influenced by ephemeral river flow and land use in the semiarid Limpopo Province of South Africa

Determining the sustainability of groundwater use in drylands with high climate variability is complex. Central to this determination is an understanding of groundwater recharge and associated processes and controls. Groundwater recharge in drylands can occur by diffuse and focused recharge (focused recharge being associated with intense episodic rainfall events and ephemeral river flow, predicted to increase and intensify with climate change). This study evaluated the relative significance and dominant controls on these two recharge processes. Ten groundwater hydrographs with multidecadal observations were collated from the Limpopo Province, South Africa, based on their proximity to river channels and rain gauges, representing diversity in local climate, landscape, vegetation, and hydrogeological conditions. The hydrographs showed that groundwater-level rises are sensitive to rainfall intensity during the rainy season, with generally larger increases after years with large episodic rainfall events, which disproportionately contribute to groundwater replenishment. Recharge processes and annual recharge volumes were quantified using the water-table fluctuation method and the numerical model HYDRUS-1D. This allowed for the inference of additional recharge contributions from focused recharge in proximity to ephemeral rivers, up to a factor of five relative to diffuse recharge. The analysis revealed synchronicity and linear correlation between annual river discharge and recharge close to the river, substantiating the importance of focused recharge close to the river network. The study showed that recharge in drylands is subject to large spatial and temporal variation and that consideration of focused and episodic recharge is critically important for managing groundwater resources at various scales in these regions.


Introduction
Drylands, defined as areas with an aridity index less than 0.65 (UNEP-WCMC 2007), are characterized by high spatial and temporal variabilities in rainfall, which is often concentrated in episodic, intensive storm events interspersed with longer dry periods.These storm events are projected to increase in intensity and frequency under future climate change, whereas dry periods are expected to become more protracted, especially in the tropics (Allan et al. 2010;Taylor et al. 2013;Fisher and Knutti 2016;Kendon et al. 2019).River systems in drylands are often ephemeral, with intermittent flow occurring in response to storm events (Taylor et al. 2013;Wekesa et al. 2020;Acworth et al. 2021).The periodic nature of river flow in drylands challenges its use as a reliable source of water supply.Consequently, communities in many of these regions are highly dependent on groundwater resources to meet freshwater demands for household and agricultural purposes (Ebrahim et al. 2019;Issoufou Ousmane et al. 2023).In contrast to surface waters, groundwater is widely available throughout the landscape, offers natural interannual and interseasonal storage, and is often of higher quality than surface-water sources.
However, the long-term renewability of groundwater withdrawals in drylands is dependent upon recharge under both current and future climate conditions (Scanlon et al. 2006;Gleeson et al. 2020;MacDonald et al. 2021).Concurrently, rising freshwater demand under climate change is intensifying dependence upon groundwater resources to shield communities and local economies against the threat of water and food insecurity (MacDonald et al. 2012;Gaye and Tindimugaya 2018).Enhancing the sustainable use of groundwater could reduce rural poverty and increase water and food security in sub-Saharan Africa (Braune and Xu 2010;Pavelic et al. 2013;Villholth 2013;Cobbing and Hiller 2019;Bellwood-Howard et al. 2022).
Past research has highlighted substantial spatial and temporal (seasonal to interannual) variations in recharge in drylands mainly due to variations in rainfall (Vries and Simmers 2002;Scanlon et al. 2006;Healy 2010).Groundwater replenishment in drylands can occur principally by diffuse and focused recharge (Goni et al. 2021), potentially occurring simultaneously (Scanlon et al. 2006).Diffuse recharge takes place throughout the landscape in a distributed manner as rainfall infiltrates the soil surface, leaves the root zone as net rainfall, and subsequently percolates to the water table.This recharge pathway is influenced by land use, topography, soil and geological conditions, and spatio-temporal variations in rainfall.
Rainfall varies on timescales that are hourly to seasonal and interannual variability to long-term climate-driven patterns.Interannual variability is often an important factor as several years of low recharge due to a scarcity of large, intensive rainfall events may be followed by a year with one or more extreme rainfall events leading to potential disproportionate increases in recharge and associated abrupt rises in groundwater levels.
Focused recharge occurs locally beneath ephemeral rivers and ponds in response to temporary accumulation of water in river courses and topographical depressions from intense rainfall (Favreau et al. 2009;Taylor et al. 2013;Seddon et al. 2021).This recharge component depends on rainfall intensity and duration, topography, geomorphology, catchment-wide surface and subsurface characteristics, including infiltration capacity of the soil, as well as the hydraulic properties of the river system and the river bed (Scanlon et al. 2006).It can be significant at the local scale (Scanlon et al. 2006) but the overall net effect on groundwater storage depends on several factors including the depth of the unsaturated zone and links between any shallow aquifer systems underlying ephemeral water bodies and deeper, larger aquifers (Walker et al. 2018;Kotchoni et al. 2019;Acworth et al. 2021).Due to the complexity of the processes involved, the focused recharge component, though predominantly close to ephemeral water bodies, is often more difficult to quantify (Wheater et al. 2010;Cuthbert et al. 2016;Acworth et al. 2021;Sorensen et al. 2021).
Under the ongoing intensification of rainfall, focused recharge is expected to increase in its magnitude relative to diffusive recharge.This is further accentuated by anthropogenic effects of urbanization (reducing diffuse recharge and increasing storm flow), catchment degradation (increasing overland flow) and increasing depletion of surface-water systems resulting in an increasing number of ephemeral systems (Blackman et al. 2021).This transition calls for a better understanding of the processes related to dryland ephemeral systems including the relative significance of focused recharge in terms of aquifer replenishment and storage as well as hydrological and ecosystem support functions (Messager et al. 2021;Blackman et al. 2021).
A widely used method to estimate groundwater recharge is the water-table fluctuation (WTF) method (Healy and Cook 2002;Healy, 2010).The approach was used by Cuthbert et al. (2019) to estimate recharge based on multidecadal groundwater level records from 14 wells from hyper-arid to humid environments across sub-Saharan Africa.Relations between annual values of recharge and rainfall were assessed alongside characterization of aquifer response mechanisms.Though specific local controls on recharge beyond rainfall were not thoroughly addressed, it was concluded that focused recharge plays a significant role in drylands.Kotchoni et al. (2019) similarly used the WTF method and explored the relationships between annual rainfall and diffuse recharge in different geological environments in humid Benin.They found linear relationships between groundwater recharge and rainfall exceeding a threshold as well as significant control of hydrogeology and depth of unsaturated zone on the recharge volumes relative to rainfall.
Physically based numerical unsaturated flow models such as HYDRUS-1D (Simunek et al. 2005) characterize diffuse recharge in a one-dimensional (1D) framework.Using observed data for precipitation and reference evapotranspiration to drive the model and parameterizing the model based on vegetation and soil characteristics, diffusive recharge can be simulated as the water flux at the groundwater table.It is not possible to directly consider focused recharge by this method.
This study combines numerical simulations of diffuse recharge by the HYDRUS-1D model and estimations of total recharge by the WTF method to estimate the relative magnitude of focused recharge and to evaluate dominant environmental controls on diffuse and focused recharge in the Limpopo Province of South Africa.Existing comprehensive national databases are exploited to extract time series of observed rainfall, pan evaporation, groundwater level, and river discharge from a range of locations in the region.

Study area
The study area is located in the Limpopo Province, which covers much of the northeastern part of South Africa (Fig. 1) and is part of the transboundary Limpopo River Basin that extends over parts of Botswana, Mozambique, and Zimbabwe.The climate in the study area, located around the city of Polokwane (Fig. 1), is semiarid with a mean annual rainfall of about 600 mm and potential evapotranspiration of ~1,600 mm (aridity index of 0.38; Zhu and Ringler 2012).The rainy season extends from October to April with substantial interannual variability in rainfall and weather extremes that give rise to both flooding and severe droughts, often correlated to the ENSO Index (Weldon and Reason 2014;Kolusu et al. 2019;Sorensen et al. 2021).River discharge in this region is mostly ephemeral and flood events are episodic.Such a river regime imposes severe challenges in sustainably supplying the growing population and economic sector with water for domestic, agricultural, and industrial uses while maintaining ecosystem function (Petrie et al. 2015).Although groundwater resources in the area are naturally limited in terms of their replenishment and storage, they have traditionally been relied upon to meet most needs including irrigation for agriculture.Evidence of groundwater depletion, especially around intensively irrigated areas and in periods of less intense rainfall has, however, been observed (Ebrahim et al. 2019;Sorensen et al. 2021).
The geology of the study area is dominated by weathered and fractured hard rocks, primarily in the form of crystalline gneiss (Muchingami et al. 2021).The subsurface is intersected by subvertical dolorite dykes, which act as barriers to perpendicular groundwater flow while facilitating parallel flows and storage accumulation up-gradient of dykes.Due to fracturing and higher permeability zones adjacent to the dykes, these are often target areas for groundwater exploration and abstraction (Ebrahim et al. 2019).Over large parts of the study area, the hydrogeology comprises an upper aquifer of unconsolidated weathered gneiss (regolith) to depths ranging from 12 to 50 m below the surface that, in turn, overlies a semiconfined aquifer of fractured gneiss bedrock extending to depths of ~120 m.Under this, unfractured bedrock is present.The fractured bedrock is more transmissive than the regolith and is the predominant target for most groundwater development.An additional significant source of groundwater exists in local alluvial aquifers associated with the ephemeral river network (Holland and Witthüser 2011;Mpala et al. 2020).The hydraulic connectivity between the hard rock and the alluvial aquifers is, however, not well understood.
Land use within the region is overall divided into natural vegetation and agriculture, the latter including commercial and small-scale cropping as well as extensive livestock rearing.The Limpopo Province is well known as an important agricultural region within South Africa, growing a wide range of crops (Cai et al. 2017).Apart from smallholder irrigation in smaller areas, intensive irrigation takes place throughout the year with the use of center pivot schemes (Ebrahim et al. 2019).The global database ISRIC-World Soil Information (2023) indicates that the soil types of the region are sandy to sandy clay loam.

Data
To investigate recharge dynamics and their relation to climatic, hydrological, and landscape conditions, a range of datasets was employed and are outlined in the following.Multidecadal observations from groundwater monitoring wells, rainfall, and river gauging stations in the Limpopo Province were collected from the Department for Water and Sanitation (DWS).Daily monitoring data were available although gaps were present in the records.Hydrogeological information included the name of the well, coordinates, length of time series, time of measurements, groundwater levels measured in meters below ground level, and a quality control code.The geological setting of each well was determined from a combination of information from geological maps and specific borehole logs from the online National Groundwater Archive hosted by DWS.Daily rainfall data were also retrieved from DWS including the name of the station, length of the time series, rainfall amount, and a quality control code.Additional daily rainfall data were collected from the South African Weather Services (SAWS) to expand the number of rainfall stations, thereby increasing the likelihood of having site locations for both rain and well gauges in proximity.In total, data from 209 monitoring wells, 98 rainfall gauges, and 7 river gauging stations were collected.Bulk figures of groundwater abstraction were retrieved from the literature (Ebrahim et al. 2019) since abstraction data for individual wells were not available.Daily evaporation measurements from class A evaporation pans were retrieved from DWS and SAWS.To convert the measurements to values for reference evapotranspiration, a standard pan coefficient of 0.7 was used (Dingman 2002).

Selection of observation stations
To analyze the dynamics of and relationships between rainfall and groundwater levels, monitoring well and rainfall gauge should ideally be co-located.Given the availability of observation locations, separation distances were minimized.The target was initially a separation distance of 10 km as suggested by Kotchoni et al. (2019) and Cuthbert et al. (2019), though there is no consistent evidence-based methodology in or close to the tropics.Here, this value was increased to 20 km to enable the inclusion of wells with longer and more reliable time series.Further, a criterion for time series selection was that the interpretability of the groundwater hydrographs should not be extensively affected by abstraction.
To capture focused recharge, the observation wells need to be in close proximity to the river.The distance between the wells and the river sections across the sites, except site 4, varied between 80 and 950 m (Table 1), which may be appropriate as the impact of recharge from a losing river was estimated to extend perpendicularly in the order of 1-2 km in similar semiarid environments (Cuthbert et al. 2016).
Ultimately, the analysis identified six observation sites where the integrity and interpretability of records were deemed reasonable and representative of the diversity in climate, landscape, and hydrogeological conditions.Each site included one rainfall station, at least one well, and preferably a river gauging station all within 20 km distances between measuring points.Due to the sparse availability of river gauging stations in the area, two of the sites (sites 2 and 4) did not include river flow measurements (Table 1).
Figure 1 shows the location of the six sites, mostly located within or close to the Hout/Sand River catchments, forming the basis for the analysis here.Although the procedure involved the removal of measuring stations with obvious anomalies, such as erratic measurements, large gaps, and impacts of other factors such as groundwater abstraction, use of observational records often involves compromises in the selection of data that are coherent while fulfilling minimum interpretability criteria (Taylor et al. 2019).Particularly it is noted that the records from gauging station A7H007 (site 3) had many quality flags for missing data.
Due to a lack of consistent evaporation data across the time series, in which some years had very intermittent data, an average value for each day of the year was estimated and the applied throughout the entire time series.Table 1 lists the relevant characteristics of the monitoring stations of the six observation sites.

Recharge estimation methods
It is assumed that recharge estimates based on the WTF method are composed of both diffuse and focused components and by comparing to the simulations using HYDRUS-1D, which only represents the diffuse part, it is possible to separately estimate the two components and infer the relative significance of focused recharge from the ephemeral rivers.

The water-table fluctuation method
The WTF method estimates recharge based on the groundwater level responses embodied in the groundwater hydrographs (Healy and Cook 2002;Cuthbert et al. 2019).The methodology assumes that recharge can be determined from a water balance calculation where recharge is estimated by the incremental change in groundwater level, lateral drainage of groundwater, and drawdown due to groundwater abstraction (Eq.1): where R t  2019) is used for the analysis.The analysis was carried out to estimate annual recharge values based on annual pickpoints from October 1st representing the beginning of the hydrological year and the rainy season to September 30th the following year, representing the end of the dry season.Estimation of a representative groundwater drainage constant (D) for each well applicable for all years was derived from observed dry-season groundwater recessions by calculating the rate of approximated linear decline of the groundwater level during this time of the year and taking the mean for all years.This value was assumed to be applicable also in the rainy season.Groundwater abstraction values are unknown but assumed less important for sites selected, and are thus not considered in the calculation (i.e., s t = 0).Applied values for specific yields for the various aquifers were derived from literature based on the geological characteristics (Table 2).

The HYDRUS-1D method
Estimation of diffuse recharge was carried out using the HYDRUS-1D model (Šimunek et al. 2005).The model simulates the vertical flow of water through the unsaturated zone by solving Richards' equation using a finite element numerical method: where (positive upwards) and S is the sink term representing the water uptake by the root system [T -1 ].
The soil hydraulic functions are represented by the van Genuchten (1980) expressions given by The parameter values of the hydraulic functions applied in the model simulations were retrieved from a database contained in the HYDRUS-1D software package where the entry point was the soil class.
The sink term S representing the water uptake by roots is described by the expression developed by Feddes et al. (1978) where the root uptake water stress response function α 1 (h) is a dimensionless function of soil water pressure head h (0 ≤ α 1 ≤ 1) and S p [T -1 ] is potential water uptake.
The potential root water uptake S p is assumed to be linearly distributed over the root zone where T p [L/T] is the potential transpiration rate and L R [L] is the depth of the root zone, which can vary over the season.
Crop-specific values for the parameters in this model were retrieved from the database contained in the software package.
HYDRUS-1D was applied to the individual sites based on the following model setup and assumptions/simplifications: 1.Given the conceptualization of the Hydrus-1D model, single-domain Darcy type of flow is assumed to occur only in a vertical direction from the soil surface to the groundwater table, i.e. lateral flow components are not considered.2. The soil profile is assumed homogenous from the soil surface to the groundwater table.3. Daily observed values of rainfall and reference evapotranspiration are drivers of the model.The reference evapotranspiration is specified as transpiration and thus soil evaporation and interception are not considered explicitly. (4)

Groundwater hydrographs and rainfall
In Fig. 2, observed groundwater levels (hydrographs) are shown together with daily rainfall (hyetographs) from the six sites; two well records are available for sites 1-4, whereas only one well record is available for each of sites 5 and 6.Significant variability exists in the frequency of the groundwater level measurements, with intervals varying on hourly, daily, or monthly scales.Sharp, sudden declines followed by rapid rises in groundwater levels at wells A6N0059 (2015) (site 1), A7N0549 (2011) (site 3), and A7N0561 (2015) (site 3) are assumed to be the result of pumping (or pumping tests) or other anomaly and were excluded from consideration in the analysis.
At site 1 (Fig. 1), two wells (A6N0059 and A6N0079) are close to the ephemeral Mogalakwena River (370 and 250 m, respectively) and installed within a sandy alluvial aquifer (screened at 25 and 21 m depths, respectively).The hydrographs for both wells show clear evidence of episodic recharge events by rapid and high increases in groundwater levels in response to high rainfall events, occurring with recurrence intervals of several years.The two wells have similarities in the timing and magnitude of episodic recharge with significant rises in for example in 1988 and 1995-1996.Significant recharge events also occur in 2005-2006, 2008-2009, and 2009-2010 in well A6N0059, but are less pronounced in A6N0079.There are years with an apparent lack of recharge despite high rainfall, for example well A6N0059 in 1985 and well A6N0079 in 2018, and in 1993 for both wells.These cases may be explained by factors such as rainfall not arriving in intense, concentrated storms; the aquifer being fully replenished, and recharge being rejected; and an antecedent, very dry year, respectively.
At site 2 (Fig. 1), hydrographs for two wells (B5N0011 and B5N0013) are presented that are both proximate to the Frisgewaag River (260 and 950 m, respectively) and screened within fractured gneiss at depths of 70 and 85 m, respectively.These two hydrographs also exhibit episodic recharge patterns but B5N0011 has a more dynamic response with larger events recorded and often followed by larger declines in groundwater levels.Both wells show a dramatic increase in levels in response to the 1995-1996 season.Whereas B5N0011 responds extensively in 2000, this is muted in the hydrograph of B5N0013 as explained in the following.Conversely, sharp responses in B5N0013 occurring during 1974-1975 and 1988-1989 are more subdued in B5N0011.During 1999-2000, a major event occurs, particularly in B5N0011, in this case with an extraordinary groundwater level rise of over 50 m followed by a long period of fairly stable levels.Incidentally, both wells have a much less dynamic behavior after 2000.In B5N0013, the aquifer responds very differently if the groundwater level is above and below approximately 36 m depth.The different modes of behavior above and below this depth are likely due to shifts in specific yield and hydraulic conductivity in the aquifer at this level.Once below, levels continue to drop until they are restored by large episodic events, whilst above, the recharge responses are much smaller with groundwater levels remaining nearly constant at 36 m once recession has returned levels to this depth.Similar dynamics are also visible in other wells but not as pronounced as in B5N0013.Another noticeable observation is the decline in groundwater levels over several years during the 1980s and 1990s particularly for B5N0011 signifying that groundwater discharge is intermittently in excess of recharge for more than 10 years.
At site 3, two wells (A7N0549 and A7N0561) are located in the Sand River subcatchment (Fig. 1) close to this ephemeral river (100 and 530 m, respectively) within a weathered gneiss aquifer (regolith); wells are screened at depths of 26 and 20 m, respectively.Both wells show relatively regular annual recharge responses yet larger episodic responses are visible in A7N0549 during the 1995-1996, 1996-1997, and 1999-2000 rainy seasons following an almost 15 years overall decline in groundwater levels.Episodic responses are much less pronounced in A7N0561, and the time series shows a regular seasonal/annual pattern with no longer-term declines.Large groundwater level increases in A70549 during episodic events are more pronounced when the aquifer has been depleted, which indicates that areas without long-term depletion could be affected by rejected recharge (i.e. the phenomenon that recharge is reduced due to lack of available subsurface storage).
At site 4, two wells (A7N0524 and A7N0529) are located in the Hout River subcatchment (Fig. 1), remotely from rivers (17 and 11 km, respectively).Well A7N0524 is located in an area with natural vegetation and thus with a permanent land cover over the year, while well A7N0529 is located in agricultural land with seasonality in land cover.Wells are screened at depths 50 and 60 m, respectively.Groundwater levels in both wells occur at depths ranging from 20 to 50 m but present contrasting groundwater-level responses to rainfall over time.The time series for A7N0524 is characterized by episodic recharge events in 1981-1982 and 1999-2000 that maintain relatively higher groundwater levels for a few years followed by long-term steady declines with little seasonal fluctuations.In A7N0529, there is a strong seasonal/ annual dynamic groundwater level response that similarly peaks in 1981-1982.A sharp decline occurs after 1992, more or less coinciding with the very steady drop in well A7N0524, but with continuous annual fluctuations.
Site 5 is represented by well A7N0019, which is screened within the shallow alluvial aquifer at a depth of 30 m and is in close proximity (170 m) to the Sand River (Fig. 1).For mid-80s to mid-90s, there is a long-term decline in groundwater levels with short-term intermittent increases in response to recharge events.These are not, however, sufficient to compensate for the longer-term trend.More substantial recharge events occur in the 1995-1996, 1999-2000, and 2013-2014 rainy seasons.
At site 6, well A8N0508 is located close to the upper reaches of the Nzhelele River in an aquifer intersected by intruding dolorite dikes and screened at 110 m depth.This hydrograph shows a consistent annual response to rainfall in groundwater level rise.The largest rises coincide with rainfall events of over 100 mm/day in 2010 with levels returning to ~5 m below ground level every year.The record from this site is significantly shorter (~14 years) than the other sites.
Interestingly, the sites do not show excessive net drawdown over a long time horizon of almost 50 years, from the 1970s to late 2010s; indeed, for some like site 2, a net increase is observed despite concerns of overexploitation and groundwater depletion regionally.

Groundwater and river hydrographs
Figure 3 shows observed groundwater and river discharge hydrographs for four of the sites (two sites are missing river discharge observations).The ephemeral nature of the rivers can be observed from the highly variable nature of flows recorded that involve large individual peaks reported for some days, which are followed by long stretches across years of little to no-flow.For all sites in Fig. 3, which are all close to the river, there is a consistent concurrence between the timing of high river flows and episodic recharge events indicated by the sharp rises in the groundwater hydrographs.Note that the discharge time series for the various sites are of different lengths.

Estimates of diffusive groundwater recharge
At site 4, groundwater recharge reflected by the hydrographs can be assumed to be solely attributed to localized diffuse processes as the two wells are located more than 10 km away from the river and therefore responses to focused recharge are presumed to be negligible although infiltration of water temporarily accumulated in topographical depressions could potentially occur.The two hydrographs at site 4 show two distinct types of responses to rainfall; one fluctuates on a seasonal/annual basis (well A7N0529), whereas the other is dominated by episodic responses (well A7N0524; Fig. 2).These differences imply sensitivities to local controlling factors including climate, vegetation, soil, and geology, though pumping cannot be excluded as a contributing factor.Given that the two wells are located not too far from each other, it is assumed that climate, soil, and geology are similar and that land use is a dominant factor explaining the difference in responses.
A7N0529 is also located in agricultural land where a seasonal variation in vegetation cover and root depth is assumed with a maximum depth of 100 cm during the midst of the growing season.The soil data suggest that its textural class is sandy and based on a sensitivity analysis the sandy soil type is selected from the soil catalog in HYDRUS-1D as the most appropriate representation of the soil characteristics.The hydraulic parameter values for this soil type are listed in Table 3.
Based on these model assumptions, the simulation results of groundwater recharge shown in Fig. 4a are obtained.As shown in the figure, the dynamics of the simulated recharge fluxes at the water table overall are in reasonable agreement with the observed groundwater-level fluctuations, suggesting that the model gives a credible description of the flow dynamics at the site; simulated mean annual diffusive recharge for this land cover is 150 mm/year.Well A7N0524 at site 4 is located in an area with natural vegetation, implying that the land cover is permanent over the year and with other vegetation and root zone characteristics distinct from those of well For the simulation of the flow conditions at this site, a sandy soil type produces the best fit.Consistent with the natural vegetation of trees and shrubs, the root depth is kept constant at 200 cm throughout the year.This parameterization introduces a higher water retention capacity of the root zone and the simulations produce more pronounced peaks only in response to heavy rainfall events reflected in the episodic nature of groundwater level rises (Fig. 4b).Significant variation exists in the relative size of the peaks that do not always correlate to the magnitude of groundwater level rises in the hydrograph such as during the 2011-2012 season.Also, the long recession in groundwater levels after the apparent recharge input during 1991-1992 is not captured by the simulations.It may be explained by abstraction taking place in the vicinity during the 1990s not accounted for in the model.Nevertheless, a reasonable overall agreement is found between the episodic behavior of measured groundwater levels and simulated recharge.The simulated mean annual diffuse recharge rate for A7N0524 is 57 mm/year.Lower recharge at this site is consistent with higher evapotranspiration from the natural perennial and deep-rooted vegetation.For comparison, recharge for the two wells is computed using the WTF method.Since the focused recharge contribution is assumed to be insignificant in this case, the two methods are expected to give similar results.To display the results, annual recharge is plotted against rainfall within a hydrological year to see whether both estimates provide consistency in the overall levels and trends given that the WTF method is calculated independently of rainfall measurements (Fig. 5).Rainfall is given as annual values within the hydrological year, which is justified given that the vast majority of events fall within the rainy season.
For well A7N0529, the HYDRUS-1D recharge simulations for (hydrological) years of relatively low rainfall appear lower than the WTF estimations (Fig. 5a); this relationship reverses for two events of very high annual rainfall.From a visual inspection, the HYDRUS-1D estimates seem to increase relatively linearly with increasing rainfall.WTF estimates increase more steeply with rainfall for rainfall less than 300 mm and then taper off to remain between 150-250 mm/year for a wide range of annual rainfall values between 300 and 800 mm.
At A7N0524, the HYDRUS-1D model appears to estimate greater recharge than the WTF model, particularly for rainfall less than 450 mm (Fig. 5b).Two high recharge events are seen in the WTF estimates with values of more than 100 mm of annual recharge.The corresponding HYDRUS-1D estimates for these years diverge from these values with one much larger and another much lower.Such differences are also observed in Fig. 4b, as the peaks of the simulated diffuse recharge fluxes do not always correspond to an observed rise in groundwater level.

Estimates of groundwater recharge for wells in proximity to rivers
Figure 6 compares recharge estimates based on the HYDRUS-1D model and the WTF method for the wells located near river stretches.The land use in proximity of rivers is largely controlled by natural vegetation and the soil is predominantly sandy and thus the same parameterization as for well A7N0524 is applied.
For wells A6N0059 and A6N0079 at site 1 (Fig. 6a) values for the diffuse recharge simulated by HYDRUS-1D gradually increase as annual rainfall increases.Similarly, there is an upward trend in the WTF estimates with increasing rainfall.The recharge estimates for higher rainfall are generally representative of years with increased groundwater level rises and episodic recharge as seen in the hydrographic records.There is not a perfect trend because there are still years of low WTF recharge estimates and high rainfall totals, a reflection also exhibited in the combined groundwater hydrographs and hyetographs (Fig. 2).
Site 2 (Fig. 6b) shows a somewhat similar behavior as site 1.The WTF recharge estimates show an increasing tendency with increasing rainfall although this is less evident for well B5N0011 as reflected by a lower coefficient of determination.For increasing annual rainfall, there is indication of an increasing divergence for both wells between the two estimation methods, which is more evident for site 2 than for site 1.
At site 3 (Fig. 6c), there is less divergence between the two WTF estimates and WTF estimates are overall smaller than for sites 1 and 2. The results show a linear and shallow gradient in recharge estimations in both the diffuse component from HYDRUS-1D and the WTF estimations for both wells.For well A7N0549, four large recharge events are estimated for rainfall above 600 mm.
For the well at site 5 (Fig. 6d), large differences exist in recharge estimated by HYDRUS-1D and the WTF method, suggesting groundwater level rises at this site are overwhelmingly dominated by focused recharge mechanisms.
At site 6 (Fig. 6e), a mix of natural vegetation and agricultural land is observed and is then tested in HYDRUS-1D using the characteristics of both land surfaces, respectively.The results show an increase in recharge modelled by HYDRUS-1D with rainfall, whereas recharge computed by WTF data remains at a low level.

Relation between groundwater recharge and river flow
In Fig. 7, the relationship between annual recharge estimated by the WTF method and observed river discharge is shown for sites 1, 3, 5, and 6.The graphs show an apparent correlation between the two variables as reflected by the R 2 values, which generally are higher than the correlations between recharge and rainfall (Fig. 6).

Discussion
Most of the groundwater hydrographs for the investigated wells show clear evidence of episodic recharge with disproportionate increases in groundwater table after episodically high rainfall events or seasons followed by at times several years of relatively stable or gradually declining levels.The shapes of the hydrographs are dependent on the distance of the wells from the river segments, the lithological conditions including any stratification, the depth of the well screen and the land cover.Dominant episodic events are particularly seen in wells located in alluvial sand and close to the river segments (A6N0059, A6N0079 and A7N0019, Fig. 2), which is indicative of the groundwater levels impacted by focused recharge from the river.This is substantiated by comparing the groundwater and river hydrographs (Fig. 3), which clearly show a concurrence between high river flow events and rapid increases in groundwater levels in wells close to the rivers that are indicative of the dissipation of focused recharge from the river.
To differentiate the components of diffusive and focused recharge the WTF method and the numerical model HYDRUS-1D were used in combination.The 1 3 HYDRUS-1D model was driven by daily values for rainfall and reference evapotranspiration.These were defined from gauge observations, whereas for reference evapotranspiration, spatial and temporal averaging were performed.Uncertainties in model parameterization, primarily vegetation characteristics (variation in root depth over the season) and soil properties (retention and hydraulic conductivity functions) constrained the analysis.Predefined sets of parameter values for vegetation and soil characteristics as provided in the modelling software were adopted and the parameter values providing the best simulations of the dynamics of groundwater recharge in comparison with observed groundwater hydrographs at the two locations at site 4 were applied.Both are sufficiently remote from river discharge for it to be assumed that recharge results entirely from diffuse recharge.The data availability, lacking measurements of soil moisture and actual evapotranspiration, did not allow for a stringent calibration; nevertheless, the adopted parameterization in the HYDRUS-1D model was used as a baseline for simulating diffusive recharge for agricultural crop and natural vegetation, respectively, and compared these to recharge simulated by the WTF method.To estimate the diffuse fluxes next to the river, the same soil type was assumed, rendering it possible to interpret a difference in contribution of focused recharge in settings close to and away from the river, notwithstanding expected variations soil type in the floodplains closer to and further away from the river.
The WTF method is also based on assumptions and estimations subject to uncertainties.The greatest uncertainty lies in the applied value of the specific yield S y for the aquifer (Eq.1).However, even if this value was estimated too high, realistic lower values for the parameter would still lead to recharge estimates by the WTF method that would be indicative of focused recharge (data not shown).Another unaccounted factor in the analysis is the vertical heterogeneity of specific yield values resulting from geological stratification The results by the two methods are compared at the two wells at site 4 (Fig. 4).Both methods are in agreement, indicating that recharge is smaller for well A7N0524, the natural vegetation site, and explained by higher evapotranspiration due to perennial deep-rooted vegetation.Accepting that both methods are subject to uncertainties as outlined in the preceding, they nevertheless are in reasonable agreement in terms of absolute magnitude of recharge for the two cases representing two different land uses.It is thus contended that the results simulated by the HYDRUS-1D are indicative of the magnitude of diffuse recharge.Furthermore, the WTF model provides comparable results to HYDRUS-1D (Fig. 5), suggesting that there is no substantial contribution of focused recharge, as expected, at this site.The unexpected higher diffuse recharge simulated by HYDRUS-1D for well A7N0524 compared to the WTF method could be an artifact if significant overland flow occurs in this upland location (contributing ephemeral flow to the river downstream).This is not accounted for in the model, which assumes no overland flow, and hence possibly overestimates recharge from a larger pool of potential infiltration.
A comparison of estimated recharge for the wells at site 1 (Fig. 6a) shows that, with increasing rainfall, there is a tendency for increasing divergence between the HYDRUS-1D and WTF estimates.These differences are attributed to an additional source of recharge, focused recharge from the adjacent river system that is accounted for in the WTF method but not in the HYDRUS-1D model.This divergence, and hence the component of recharge, increases with rainfall, which is theoretically explainable from an increase in intermittent river flow as a function of rainfall and hence more riverine leakage to groundwater.Interestingly, the divergent recharge estimates exceed values of annual rainfall, indicating that this additional water derives from focused recharge.The WTF recharge values exceeding 2,000 mm are approximately double that of the corresponding annual rainfall.In this case, the focused recharge is simulated to be about a factor of 4-5 greater than the diffuse recharge in the same area.Similar inference of focused recharge can be made at site 2 (Fig. 6b).
For site 3 (Fig. 6c), the difference between the two estimates is less pronounced, indicating that the additional input from focused recharge is less significant than that seen at sites 1 and 2. Both wells are located near to the river, with well A7N0549 lying closer at 100 m and displays generally larger recharge estimates compared to well A7N0561, which is located 500 m away from the river.These observations suggest that the degree of impact of focused recharge depends on the distance from the river, although for wells close to a leaking river, mounding effects and temporary river bank storage may lead to overestimation of recharge (Favreau et al. 2009).There is also a slight difference between the response in the two wells in the sense that the well closest has some high events, while the furthest away is more linear, indicative of the dissipation of focused recharge with distance from the river.
For site 5 (Fig. 6d), the estimated recharge far exceeds the rainfall amount.Similar to site 1, this site is located in an alluvial aquifer and the well is screened at shallow depth (30 m), which promotes high recharge volumes with estimates well above the annual rainfall values in most years, and a maximum over 3,500 mm/year.
The results for site 6 (Fig. 6e) differ from the other sites as the recharge predicted by HYDRUS-1D for both land use types exceeds the WTF estimate.However, the well screen at this site is placed at 110 m depth in intrusive rock.Thus, the impact from focused recharge may be small due to the presence of less permeable and potentially confining sediments lying above the screen intake.Such geological conditions are not considered in the HYDRUS-1D model leading to an overestimation of diffuse recharge at these depths.
Comparing the relationship between annual recharge estimated by the WTF method and observed river discharge (Fig. 7) with those of WTF-estimated recharge and observed rainfall (Fig. 6), it was found that although the number of data points behind the two analyses are different, the obtained correlation values indicate that river discharge is a stronger predictor of recharge than rainfall close to the river, and that focused recharge (deriving from river loss) is a critical recharge process around the rivers.Seddon et al. (2021) reported similar stronger correlations between recharge estimates and river discharge observations than between recharge estimates and rainfall recorded close to piezometers for drylands in Tanzania.Sorensen et al. (2021) also found that groundwater storage is more correlated to river flow dynamics than rainfall dynamics in the South African part of the Limpopo River Basin.
In this study, recharge events over an annual timescale (hydrological years) were analyzed with more or less annual, or less frequent, recurrent monocyclic rises in groundwater levels coinciding with the onset of the rainy season and with an assumed linear groundwater drainage and associated groundwater level recession until the next rise.As a result, the resolution of the WTF method was set on an annual basis starting on October 1st.This allows the groundwater level to recede to background levels and a recharge value can be estimated (Cuthbert et al. 2019).A higher resolution on a monthly and daily basis was also attempted, but this led to inaccurate estimations due to limitations in the regularity of groundwater measurements, which are highly variable within and across the well dataset available.
The total number of wells used in the study was much smaller than the number initially collected.Many were excluded in terms of interpretability due to long gaps in the data series.These gaps have limited the comprehensiveness and the spatial scale at which interpretations can be made when considering the vastness of the Limpopo Province.In some hydrographs, large gaps in the order of years may have resulted in significant and smaller recharge events being missed.Greater inclusion of wells, river and rainfall gauges across the study area would have been desirable to enhance the representativeness of the study.Importantly, many monitoring wells are impacted by groundwater abstraction from neighboring wells which makes their validity for recharge estimation doubtful or at best results uncertain (Sorensen et al. 2021).Unfortunately, groundwater abstraction data are very limited.
The distance between rainfall and river discharge gauges and wells is an important parameter and hence, an initial effort was made to keep this distance below 10 km.However, this criterion was not possible to fulfill in all cases.The majority of the wells included were in close proximity to the river.This was useful in terms of estimating focused recharge close to the river but led to a shortage of representation in terms of the gradient from focused to diffuse recharge dominance.
All wells are assumed to be placed in unconfined aquifers; lithologies were determined through a combination of geological maps and borehole data.Screen depths are assumed to be the same as the depth of the well.Uncertainty exists due to unknown hydrogeological conditions throughout the profiles, for example, at the well at site 6 which extends to a depth of 110 m.

Conclusions
Multidecadal records for groundwater level, rainfall and river discharge from the semiarid Limpopo Province of South Africa, up to 50 years from the early 1970s, were used to analyze the spatial and temporal variation of groundwater replenishment in the area and how it relates to environmental factors such as rainfall, river discharge, land use, soil and geology.Many of the groundwater hydrographs showed clear evidence of episodic recharge events; recharge dynamics varied greatly among the different sites.
No significant and consistent long-term groundwater storage decline or accumulation trends were observed over the period, though intermittent decadal persistent declines were not uncommon, followed by abrupt large recharge events influenced by extreme rainfall and river flow events, which fully or partially contributed to recovering aquifer storage volumes.
Two methods were used to estimate groundwater recharge, the water-table fluctuation method (WTF) method and the numerical model (HYDRUS-1D), which enabled differentiating the contributions from the two principal recharge mechanisms occurring in drylands, diffuse and focused recharge.
The results showed that with increasing rainfall, and in particular with the prevalence of episodic recharge events, the WTF method generally produced higher recharge estimates at localities next to the river than HYDRUS-1D.Further, focused recharge from ephemeral river systems manifested in synchronous rises of groundwater levels with river discharge.Estimated total annual river-proximal recharge displayed a closer correlation with accumulated annual river discharge than with rainfall.
These results imply that increased contributions derive from an additional recharge pathway in the form of focused recharge from the river systems that is not considered/represented in HYDRUS-1D.Inferred focused recharge fluxes disproportionately replenish aquifers in close proximity to ephemeral rivers, up to a factor of five relative to diffuse recharge and, hence, augment the availability of groundwater resources in these areas, which is often disregarded in conventional water resources assessments but of critical importance for groundwater management and use.
Focused recharge from ephemeral river system manifests in near-synchronous rises of groundwater levels with river discharge; estimated total annual recharge displayed a closer correlation with accumulated annual river discharge than with rainfall.The HYDRUS-1D simulations of diffusive recharge for natural vegetation showed more episodic recharge events than for agricultural crop cover as soil moisture deficits were not overcome in many years in these natural settings.
The analysis has only documented the dynamics and significance of focused recharge at the local scale.To estimate the mutual significance and relative magnitude of diffuse and focused recharge at a larger scale including the catchment scale, development and calibration of an integrated hydrological model that considers all hydrological processes are required.Such a model will also allow for an integrated analysis of the influence of environmental factors such as climate, land use, soil, geology and distance from the river.

Fig. 2
Fig.2Multidecadal groundwater hydrographs (red graph) and hyetographs (blue graph) for observation sites 1-6 in the Limpopo Province in South Africa.Gaps in rainfall data are indicated by grey-shaded areas.Groundwater level is shown as depth, in meters below ground level(mbgl)

Fig. 3
Fig. 3 Groundwater and river hydrographs for a site 1, b site 3, c site 5, and d site 6 in the Limpopo Province in South Africa

Fig. 4
Fig. 4 Simulated recharge and observed groundwater hydrographs for site 4 for a well A7N0529 (agricultural crop) and b well A7N0524 (natural vegetation)

Fig. 5
Fig. 5 Comparison of annual recharge estimates using the WTF method and HYDRUS-1D model for a well A7N0529 (agricultural crop) and b well A7N0524 (natural vegetation)

Fig. 6
Fig. 6 Comparison of recharge estimates using the WTF method and the HYDRUS-1D model for a site 1, b site 2 c site 3, d site 5, and e site 6.Coefficients of determination R 2 for linear regression between

Table 1
Characteristics of monitoring stations (Observation period is the time period of overlapping rainfall and groundwater level measurements; Mix indicates a mix of natural and agricultural land cover)

Table 2
Ѳ r is residual water content [L/L], Ѳ s is water content at full saturation [L/L], K s is saturated hydraulic conductivity [L/T], S e is effective saturation, and α [L -1 ] and n [-] are empirical parameters.
All rainfall is assumed to infiltrate, i.e. no overland flow is considered, as the soil type used in the simulations is of the sandy type. 5.The most prevalent land use cover of each site is defined from Google Earth photos.6. Root depths are assumed to be variable over the season for agricultural crops and constant for natural vegetation (no specific data are available).7. The lower boundary condition is zero pressure head at the location of the water table, which is specified at the observed mean depth of the groundwater table for the individual wells and kept constant in time.The seasonal variation of the groundwater level is not of importance as the variable of interest is the flow at this approximate level.8.The model simulates the flux at the water table, which is assumed to correspond to groundwater recharge.9.The model simulates actual evapotranspiration based on specified reference evapotranspiration and root water uptake.Unfortunately, to the knowledge of the authors, no data are available for verifying the modelling results of evapotranspiration.10.The model is not subject to a dedicated calibration, but is fine-tuned manually by testing soil and vegetation characteristics for the categories of the built-in databases in the model to best reflect the dynamics of the observed hydrographs in a way that the simulated fluxes are synchronous with the observed increases in groundwater levels.Neither measurements of soil moisture nor evapotranspiration are available for calibration purposes.

Table 3
Soil characteristics (van Genuchten parameter values) used in the HYDRUS-1D simulations (K s saturated hydraulic conductivity, Θ s water content at saturation, Θ r residual water content, α and n empirical parameters)