Water quality changes and organic matter removal using natural bank infiltration at a boreal lake in Finland

changes and organic matter removal using natural bank infiltration at a boreal lake


Introduction
Managed aquifer recharge (MAR) is a method used to store or treat water by infiltrating it into an aquifer.In MAR, water is purified in the subsurface by natural processes, and it is extracted from production wells after a sufficient retention time ensuring purification.Two main methods of MAR are (1) induced bank infiltration and (2) enhanced infiltration above an aquifer from a surface water reservoir.In induced bank infiltration, groundwater production wells are located close to a river or a lake bank.Surface-water infiltrates from the surface-water body to groundwater when water levels are lowered in the aquifer, and the hydraulic gradient is directed toward the groundwater extraction wells.This is the most widespread method in Europe, used in 57% of the currently active MAR sites (Sprenger et al. 2017).In enhanced infiltration, surface water from a lake or from a river is pumped and infiltrated in an aquifer by ponds, sprinklers, or injection wells; it is used at 39% of the active European MAR sites (Sprenger et al. 2017).MAR has been in use for over two centuries in the USA and in Europe (Ray et al. 2003;Schubert 2002;Sprenger et al. 2017), and from 1965 to 2015, the MAR water production capacity worldwide has grown from 1 to 10 km 3 /year, currently the largest producer being India (Dillon et al. 2019).
Managed aquifer recharge is primarily used for drinking water production, and it is often applied in combination with engineered pre-and post-treatments (e.g., Dillon et al. 2019;Kuehn and Mueller 2000;Kurki et al. 2013).Other uses are production of process water for industry, agricultural irrigation, water storage for dry seasons, mine water treatment, wastewater Electronic supplementary material The online version of this article (https://doi.org/10.1007/s10040-020-02127-9)contains supplementary material, which is available to authorized users.treatment, and urban storm water reuse (Dillon et al. 2019;Kurki et al. 2013).In drinking water production, MAR is used to remove suspended solids, particles, biodegradable compounds, bacteria, viruses and parasites (Hiscock and Grischenk 2002).In Finland and Sweden, the main purpose of MAR is removal of natural organic matter (NOM) from surface waters for drinking water production (Kurki et al. 2013).NOM originates from decomposition products of plants, and it is a mixture of thousands of organic compounds with different molecular sizes, structures, and functional groups (Matilainen et al. 2011;Thurman 1985).The different features of NOM can be characterized by several analytical methods including spectroscopic methods (fluorescence, UV-vis, FTIR, 1 H NMR, 13 C NMR 15 N NMR, 2-D NMR) chromatographic methods (HP-SEC, FIFFF), mass spectrometric methods (LC-MS, FTICR MS, GC-MS) and bulk parameters-total organic carbon (TOC), dissolved organic carbon (DOC), specific ultraviolet adsorbance (SUVA) (Matilainen et al. 2011).NOM in surface waters consist mainly of hydrophobic and hydrophilic components.More than half of the total amount of organic matter are hydrophobic acids (Thurman 1985).
At surface-water treatment plants with NOM removal by coagulation and flocculation, the hydrophobic fraction and high molar mass compounds are removed more efficiently than the hydrophilic fraction and low molar mass compounds (Matilainen et al. 2010;Parsons et al. 2004;Vuorio et al. 1998).Seasonal variations in NOM composition affect the NOM removal efficiency, with better removal percentages in winter, as the share of hydrophobic organic matter is higher during the cold season (Sharp et al. 2006).
In soil water, organic matter is removed by mechanical filtering, sorption and degradation.In MAR, the main NOM removal mechanism is primarily initial adsorption on soil particles followed by biodegradation (Gross-Wittke et al. 2010;Diem et al. 2013;Grünheid et al. 2005).In biodegradation NOM is degraded to dissolved inorganic carbon, which is not accumulated in the soil.On the other hand, in a column test by Kolehmainen et al. (2009) only 50% of NOM removal was estimated to be due to biodegradation.According to this value, NOM accumulation by sorption and mechanical filtering in the aquifer is possible.As biodegradation of NOM to dissolved inorganic carbon consumes dissolved oxygen (DO) (McMahon and Chapelle 2008), it can significantly affect the chemical composition of groundwater.If the conditions in an aquifer become anoxic, nitrate (NO 3

−
) is reduced to nitrite (NO 2 − ), and manganese (Mn) and iron (Fe) in soil are reduced from their insoluble forms, i.e., Mn(III) and Fe(II), to the soluble forms Mn(II) and Fe(II).Managing and monitoring redox conditions, especially DO levels in MAR, is essential to secure good groundwater quality.Fe and Mn concentrations cause aesthetic and technical problems in groundwater use, and high Mn concentrations also are a health risk (World Health Organization 2011).
Studies on existing MAR sites and column tests have shown that the NOM content in infiltrated surface waters decreases rapidly, early on along the flowpath as the more easily degraded or adsorbable fraction is rapidly removed (Frycklund 1995;Grünheid et al. 2005;Hiscock and Grischenk 2002;Kolehmainen et al. 2009Kolehmainen et al. , 2006;;Lindroos et al. 2002;Niinikoski et al. 2016).Hydrogeology (e.g., sediment porosity, permeability and groundwater flow), geochemical and nutrient conditions, temperature, and redox conditions are known to have an impact on the NOM removal rate in MAR (Maeng et al. 2011).Similar to surface-water treatment plants, the average molecular size of organic matter in water decreases in MAR (Lindroos et al. 2002;Nissinen et al. 2001).This is due to a better mobility of small molecules (e.g., McInnis et al. 2014), but also degradation of NOM, which modifies the organic matter towards lower molecular weight (Kim et al. 2006).MAR removes hydrophilic acids and hydrophobic acids efficiently, whereas the decrease in hydrophobic neutral compounds is smaller (Lindroos et al. 2002).
Increased TOC concentrations have been widely reported in the surface waters of the Northern hemisphere (Monteith et al. 2007;de Wit et al. 2016), which is a concern for domestic water producers at both surface-water treatment plants and at MAR sites.Because many of the processes connected to removal of organic matter are temperature-dependent, climate warming also affects the water quality in MAR in the future (Diem et al. 2013).
Recent research on MAR has concentrated on existing plants or on column tests simulating the MAR process (Dillon et al. 2019); thus, the sustainability of MAR over a long time span is not known.In the long term, the essential question is how the capacity to remove organic matter evolves in the aquifer (Kurki et al. 2013).To fill the knowledge gap in the long-term impacts of MAR, water quality was explored in this study by using a natural lake-aquifer system as a surrogate to MAR.A natural bank infiltration site was selected to assess aquifer recharge over a greatly longer time scale than can be achieved by studying existing MAR sites.The main objective was to quantify NOM removal and possible NOM accumulation, and to assess the development of redox conditions in an existing aquifer where natural bank infiltration from a lake has been going on for thousands of years.The portion of lake infiltrate in the aquifer was estimated using the stable isotopes of hydrogen and oxygen in water.Groundwater retention time in the aquifer was estimated by groundwater flow and heat transport modeling.Water samples were collected from the lake and from observation wells within the aquifer for the analysis of total organic carbon (TOC), dissolved organic carbon (DOC), chemical oxygen demand (COD Mn ), pH, DO, Fe, Mn and electrical conductivity (EC).

Study site
To work as a surrogate to MAR, the natural lake-aquifer system had to meet the following criteria: (1) a high percentage of groundwater in the aquifer originates from a surface-water body, (2) surface water is the principal source of NOM in the aquifer, (3) water quality in the surface water and in the aquifer represent the natural, long-term water quality, i.e., human activities in the catchment have not significantly impacted the water quality, (4) the age of the lake-aquifer system can be determined and is much longer than at constructed MAR sites, (5) the infiltration zone from the lake to the aquifer could be located, the main flow direction determined and the aquifer dimensions measured to make a reliable estimate of retention times in the aquifer.Lake Iso Tiilijärvi and the adjacent aquifer in southern Finland (Lat.61°0′ 28.064″ Long.25°30′ 10.982″, 144.5 m amsl, Fig. 1) met the criteria as described in this study.The area of lake Iso Tiilijärvi is 0.534 km 2 (Fig. 1), and the lake receives inflow from surrounding slopes (0.5 km 2 ), and from Lake Keski-Tiilijärvi (0.089 km 2 ) and Lake Vähä-Tiilijärvi (0.1 km 2 ) (Finnish Environment Institute 2017).The mean annual temperature at the site is 4.5 °C; the mean annual precipitation is 636 mm (Pirinen et al. 2012); and the mean annual evapotranspiration from class A pan is 420 mm-observation from Mikkeli, Karila (Korhonen and Haavanlammi 2012).Iso Tiilijärvi is covered with ice from November to April.It has no stream outlet, and outflow only occurs through the lake bank into a groundwater aquifer.The lake depth is mainly 5-10 m, with a local maximum of 15 m close to observation site L2 in Fig. 1.
Groundwater level in the downstream aquifer rises above the soil surface at several locations, forming springs at the Rätäksuo bog (138 m amsl) northwest of the lake Iso Tiilijärvi (Fig. 1a, b).Another bog (Soisalmensuo) is located at the south-eastern end of the lake, and it discharges to the lake.Land-use in the Iso Tiilijärvi catchment consists mostly of forested outdoor recreation areas.About 20 cottages and permanent dwellings are located by the lake.The urban constructed area of Hollola (Fig. 1b) is mainly outside the catchment.As there is no intensive land-use within the catchment, the effect of human activities on the lake-water quality is expected to be small.Lake Iso Tiilijärvi is located on a large formation of glaciofluvial deposits within the First Salpausselkä endmoraine (Saarnisto and Saarinen 2001;Lunkka et al. 2004).The formation mainly consists of assorted sediment material including sand, gravel and sandy till (Benn and Evans 2010;Sauramo 1929).In the surroundings of the lake, the glaciofluvial deposits are tens of meters thick, and they are underlain by a crystalline bedrock.From the northern to northeastern side, the lake is bordered by a bedrock hill (Tiirismaa in Fig. 1b) mainly covered by a few meters thick till deposits (Geological Survey of Finland 2019).The crystalline bedrock in the area is of Precambrian origin, consisting of granodiorite and porphyritic granodiorite.The Tiirismaa hill is quartzite (Geological Survey of Finland 2019).Iron oxides that increase the ability of soil to adsorb organic matter are present in the aquifer, as the mineral material mainly inherits the mineralogical composition of the surrounding bedrock.Magnetite is present as a primary accessory mineral in granodiorite, and biotite in granodiorite weathers to iron oxides.Lake Iso Tiilijärvi was formed approximately 12,500 years ago (cal.BP), at the same time as the final moraine of the First Salpausselkä over the deglaciation of the Scandinavian Ice Sheet (Rinterknecht et al. 2004;Saarnisto and Saarinen 2001).

Hydrological mapping
To estimate the total water amount flowing through the Tiilijärvi aquifer, discharge at the Rätäksuo bog outlet was measured using a V-notch (90 o ) weir (location shown in Fig. 1b) between 5 May 2017 and 11 June 2018.The pressure head at the weir was measured manually once every 1-2 weeks and automatically every hour using a pressure sensor Solinst Levelogger (precision 0.05% FS).The water level at the weir was obtained by compensating for barometric pressure (Solinst Barologger, precision ±0.05 kPa).The aquifer thickness was measured by drilling through soil to bedrock at sampling site 4 (Fig. 1a).The soil type was estimated by visual inspection from eight soil samples taken from sites 40, 50 and 60 during well installation.

Isotopic separation of infiltration fractions
Isotopic tracing was conducted to quantify the two groundwater sources in the aquifer, namely, precipitation and lake water.Stable isotopes of hydrogen (H) and oxygen (O) were used to estimate the fraction of infiltrated lake water at each sampling site, as the isotopic composition of surface-water bodies differs from the composition of groundwater originating from direct infiltration of precipitation (Clark and Fritz 1997;Kendall and McDonnell 1998).
Samples were collected from site L1 at Lake Iso Tiilijärvi, from groundwater wells 2, 3, 4, 5, 6, 10, 30, 40, 50 and 60, and from springs 1 and 2 (Fig. 1a).Additionally, a well located 6 km southeast from the research area (coordinates Lat.60°59′ 10.115″ Long.25°3 6′ 23.277″) was sampled to assess the background quality of groundwater originating from precipitation.Samples from wells 2, 3, 4, 5 and 30 were taken with a suction pump (Honda) and from wells 40, 50 and 60 with a submersible pump (Mini Monsoon).Before sampling, a minimum of 100 L of water was pumped from each groundwater well to clear the water in order to get representative samples.Samples from wells 6, 10 and the background well were taken with a bailer, as the diameters of the wells were too small for a submersible pump, and the water table was too deep for a suction pump.Surface-water samples from the lake were taken with a Limnos water sampler.Water discharging from springs 1 and 2 was taken directly in the sample bottle.Samples were collected in 100-ml PEH bottles and stored cold and dark before analysis.
Sampling was conducted between January 2016 and February 2018, the main part being conducted every second month between April 2017 and February 2018.Stable isotopes were analyzed at the Department of Geosciences and Geography at the University of Helsinki by Picarro cavity ring-down spectrometer (CRDS; precision 0.1‰ for O and 0.5‰ for H).Two laboratory standards calibrated against international VSMOW (Vienna Standard Mean Ocean Water) and SLAP (Standard Light Antarctic Precipitation) standards, as well as quality control samples were analyzed together with the unknown samples.All samples were filtered prior to analysis.
The results of the isotopic analyses were given as δ values, which were calculated from the isotope ratios where R sample is the sample isotope ratio 18 O/ 16 O or 2 H/ 1 H and R std is the isotope ratio in a standard.The values here were given against the Vienna Standard Mean Ocean Water (VSMOW).
The portion of the lake infiltrate in each groundwater sample was estimated from where a is the lake infiltrate fraction in the sample, δ 18 O sample is the isotopic composition of the groundwater sample, δ 18 O Lake is the average isotopic composition of the lake water, and δ 18 O BG is the average isotopic composition of the background samples, representing the isotopic composition of groundwater originating from precipitation.The portion of the lake-water infiltrate can be estimated from Eq. (2) when the system wherein mixing occurs is not subject to evaporation (Clark and Fritz 1997), which is the case in most groundwater systems.

Retention
The retention time and the flow velocity of the lake infiltrate to groundwater well 30 were estimated from a temperature breakthrough curve.The three-dimensional heat transport equation, including convection and conduction in the aquifer, can be written as follows (modified from Stallman 1963, Vandenbohede andVan Houtte 2012): where T is the temperature (°C), θ is the effective porosity (dimensionless), a is the thermal dispersivity tensor (m), q is the Darcy velocity (m/s), t is the time (days), ρ s is the density of soil particles (kg/m 3 ), c s is the specific heat capacity of soil particles (kJ/kg/°C), ρ w is the density of water (kg/m 3 ), c w is the specific heat capacity of water (J/(kg °C), and k z is defined as follows: where k z is the total thermal conductivity of the aquifer (kW/ m/°C), k s is the thermal conductivity of soil particles, and k w is the thermal conductivity of water.
The retention time can be calculated from: where t is the retention time (day), d is the distance from the lake bank to well 30 (m), and q/θ is the seepage flow velocity (m/day).Equation ( 3) was solved with the explicit finite difference method as one-dimensional (1D) with a time step of 1 h and a spatial step of 0.2 m.In the 1D solution, the dispersivity tensor was reduced to only longitudinal thermal dispersivity.The transport and thermal parameters used in the model are given in Table 1.
The water temperature was measured daily in groundwater well 30 and at the observation site L1 (Fig. 1a) at 1 m depth with Solinst Levelogger (precision ±0.05 FS).Lake-water temperatures from spring 2017 and 2018 were used as the upstream boundary condition.The model was calibrated against temperatures measured in well 30.Darcy velocity and longitudinal dispersivity were used as calibration parameters.The modelling results were compared against retention times predicted by a solute transport model (see the following).
Heat convection is mathematically analogous to advection, and conduction is analogous to the combined effect of diffusion and dispersion on solute transport (Vandenbohede and Van Houtte 2012).Due to equilibration of temperature between water and the solid soil material, the temperature front is retarded, or travels slower than the averaged seepage velocity of water.
At observation wells other than well 30, the temporal variation of groundwater temperature was too small to trace the groundwater flow velocity; thus, the retention time was calculated at these locations using 1D parametrization of flow and transport models MODFLOW 2005 (Harbaugh 2005), MODPATH (Pollock 1998) and MT3DMS (Zheng and Wang 1999).The groundwater flow was simplified in the model to 1D along the main flowpath because of the dominance of the lake infiltration, the small model area and the homogenous land cover (Fig. 1b), precipitation, and aquifer material.First, a steady-state groundwater flow model from the lakeshore to the springs at Rätäksuo bog was constructed.Thegrid area of 160 m × 300 m corresponded to the area of the aquifer between Lake Iso Tiilijärvi and Rätäksuo bog.The finite difference grid consisted of one horizontal soil profile divided into 80 vertical columns representing the model grid cells (Fig. 2).
In the model area, the water balance was as follows: where Q in is the total inflow in the catchment from all water sources, Q out is the outflow from springs,Q lake is the infiltration from lake, Q prec is the recharge from precipitation within the grid area, andQ catch is the recharge from precipitation in the catchment outside the grid area.Infiltration from the lake was implemented in the model as an elevated recharge rate on the first cell next to the lake (Fig. 2).Recharge induced by precipitation was spread on all model grid cells.As the model had no-flow boundaries, unflow from the catchment outside the grid area was included in the model  The groundwater originating from precipitation and the water infiltrated from the lake were separated by marking the infiltrated lake water in the model with a constant concentration of 100 mg/L of a conservative solute in the first cell.The concentration of the recharge from precipitation was set to 0 mg/L.Mixing of the infiltrate and the water originating from precipitation was observed as changes in the concentration of each cell.Based on literature values (Airaksinen 1978), the soil saturated hydraulic conductivity was set to 10 −4 m/s, longitudinal dispersivity to 1 m and effective porosity to 0.25.
The recharge from precipitation in the model area and in the catchment, together with the lake-water inflow were used as calibration parameters.The model was calibrated by trial and error against the groundwater discharge from the aquifer and against the measured infiltrate fractions at observation wells 30, 40, 60, 5, 3 and spring 1.In the model calibration, the lake infiltration in the first cell was first adjusted to correspond to the measured share of lake water at the aquifer outlet.Then, recharge on the whole modeled area was calibrated to correspond to the share of groundwater formed from precipitation at the calibration points and the measured amount of groundwater originating from recharge at the outflow.
The model was validated against the infiltration fractions in observation wells 10, 50, 6, 4 and 2. As surface runoff temporarily increased the discharge at the weir, and as the infiltration rate from the lake to the groundwater body varied with the lake-water table level, two estimates of the groundwater discharge were used: the 25th percentile at Rätäksuo weir of measurements from 5 May 2017 to 31 October 2017 and the 25th percentile in the winter before snowmelt, i.e., from 1 November 2017 to 30 April 2018.The infiltrate fraction estimates were based on the composition of stable O isotopes in the samples (see the preceding details).

Quality of lake surface-water and groundwater
To characterize the infiltrating surface-water quality, Lake Iso Tiilijärvi was sampled with a Limnos bailer irregularly between August 2016 and November 2018.From the samples, major cations and anions (Ca, Na, K, Mg, HCO, SO 4 , Cl), Fe, Mn, turbidity, DO, EC, pH, nitrogen compounds (NH 4 , NO 3 , NO 2 ) and organic matter (TOC, DOC, COD Mn ) were analyzed.The samples were taken from observation site L1 at the depth of 1 m.Site L2 was sampled once at the depth of 7 m to assess the horizontal variability in the lake-water quality.The number of analyses related to each observation site are given in Table S1 of the electronic supplementary material (ESM).For frequently sampled observation sites (N > 5), samples were taken evenly in all seasons.
The water quality in the aquifer was investigated by sampling groundwater wells 2, 3, 4, 5, 30, 40, 50, and 60 bimonthly between December 2017 and November 2018, six times in total.From the samples, DO, EC, pH, Fe, Mn, TOC, DOC and COD Mn were analyzed.The sampling equipment and the methods were the same as in the sampling for stable isotopes.Correlation between DO levels with Fe, Mn and OM removal (ΔTOC = TOC sample -TOC lake ) was tested with Spearman's rank correlation coefficient.The significance of changes between the TOC, DO, Fe and Mn concentrations in the lake and groundwater samples was tested with the nonparametric Wilcoxon rank sum test.Tests were performed in MATLAB R2018b and Python 3.6.3.

Analysis methods
Ca, Na, K, Mg, Fe and Mn were analyzed following the standard SFS-EN-ISO 17294-2.Before the analyses, samples were filtered using a 0. Both DO and EC were measured in the field with a YSI ProDSS multiparameter water-quality meter.The measurements were conducted in a flow cell to avoid gas exchange with the atmosphere and to guarantee reliable DO measurements in low concentrations.pH was measured following the standard SFS-EN ISO 10523 and turbidity following the standard SFS-EN ISO 7027:2000 and TOC and DOC were analyzed following the standard SFS-EN 1484, with 0.2 mg/L detection limit.Before the DOC analyses, the samples were filtered in the laboratory using a 0.45-μm filter (see the preceding details).COD Mn was analyzed following the standard SFS 3036.

Hydrological properties of the aquifer
Average discharge at the Rätäksuo bog outlet was 1,760 m 3 / day (Fig. S1 of the ESM).The discharge in the summer season (5 May 2017-31 October 2017) was 1,160 m 3 /day (lower quartile), 1,210 m 3 /day (median) and 1,250 m 3 /day (upper quartile), and in the winter season (1 November 2017-30 April 2018) 1,790 m 3 /day (lower quartile), 2,090 m 3 /day (median) and 2,300 m 3 /day (upper quartile).Using an estimate of 1.70 km 2 for the drainage basin area upstream of the Rätäksuo measurement weir, the runoff was 378 mm/year.The soil type in the aquifer northeast of Lake Iso Tiilijärvi was sand and gravelly sand, and the aquifer thickness was 38.5 m, measured at groundwater well 4.

Isotopic separation
The stable isotope composition of H and O in the water samples formed a mixing line between two end members (Fig. 3a): (1) the background composition of the local groundwater representing the isotopic composition of precipitation, and (2) the composition of Lake Iso Tiilijärvi, influenced by lake evaporation.The background samples of groundwater formed from precipitation corresponded to the isotopic composition of the local meteoric water line (LMWL, Kortelainen 2007), while the surface-water samples strongly deviated from the LMWL.The difference in the average isotope values between the local groundwater originating from precipitation (background in Fig. 3a) and the infiltrated lake water (L1) was 6.0‰ for δ 18 O and 31.9‰ for δH.
The average fraction of lake water was estimated to be 0.95-0.96 in groundwater wells 30 and 40 near the shore, and 0.75-0.76 at the springs 1 and 2 at the end of the groundwater flow path (Fig. 3b).At the northeastern measurement locations 10 and 50, the estimates of the lake-water fractions were 0.28 and 0.65, respectively.The estimates in the northeastern locations reflect the increasing share of groundwater inflow formed from recharge in the Tiirismaa hill (Fig. 1).

Retention times
Figure 4 shows the calibration results of the heat transport model.The seasonal pattern of lake-water temperature was reflected with a delay in the temperature pattern in groundwater well 30.The model calibration yielded seepage velocities of 0.2 m/day in spring 2017 and 0.4 m/day in spring 2018.According to the model, the retention time from the lake to observation well 30 was 15 days in spring 2017 and 7 days in spring 2018.Thermal dispersivity was calibrated to the value of 0.1 m.The heat transport model was compared against the MODFLOW-based estimate of the retention time (Table 2).The average of the retention times in the heat transport model was 29% larger than the average of the retention times in the MODFLOW-based model.As the differences between the models were small, the models produce a consistent picture of typical flow velocities in the aquifer.The sharp spring-time temperature changes in the lake and groundwater, observed in the measurements, were linked to the ice breakup timing of Lake Iso Tiilijärvi (Fig. 4).
Table 2 summarizes  inflow.As two values for groundwater discharge were used for the calibration, i.e., 1,160 m 3 /day for the summer season 2017 and 1,793 m 3 /day for the winter 2017-2018, two sets of modeling results are given.It was observed from the isotopic data that, along the main flow path, the share of groundwater originating from precipitation increased nearly linearly with the distance from the lake (Fig. 3b).Therefore, the best fit during calibration of recharge from precipitation was achieved with a uniform distribution of recharge on the cells.Wells 10 and 50, which did not follow the linear relationship, were outside the main bank infiltration zone and the main flow path.
For the calibration wells, the error between infiltrate fractions estimated by oxygen isotopes and modelled infiltrate fractions was 0-9%-units.For the validation wells, the error in the wells located along the main flow path was 0-2%-units.Outside the intensive infiltration zone, the error was higher: 31%-units in well 50 and 71%-units in well 10.This part of the aquifer was strongly influenced by groundwater inflow from the Tiirismaa hill (Fig. 1).For the summer, the model yielded an inflow of 870 m 3 /day for the lake infiltrate and 320 m 3 /day for the precipitation recharge.The corresponding winter values were 1,265 and 475 m 3 /day.Retention times for the groundwater sampling locations, as calculated with the model, are listed in Table 2.The retention times varied from 7 to 10 days (well 30) to ca. 1.5 years (spring 1).

Water quality
Water quality in Lake Iso Tiilijärvi remained stable during the study period (Table S1 of the ESM).Only small changes (coefficient of variability <0.1 in Table S1 of the ESM) were observed in 11 variables out of 17 (N > 1).The stable isotope composition of δ 18 O and δ 2 H showed the lowest coefficient of variability (0.01).There were no clear differences in the water quality of samples taken from the observation sites L1 (1 m depth) and L2 (7 m depth), and therefore their results were joined.
Figure 5 shows the water quality development in the studied lake-aquifer system.The NOM in terms of TOC, DOC and COD Mn concentrations decreased rapidly when the lake water infiltrated in the aquifer (Fig. 5a-c).In the lake, the TOC concentration was 3.0 ± 0.20 mg/L (average and standard deviation), and in well 30 it was a Observation sites 30, 40, 60, 5, 3 and Spring 1 were used for model calibration; the other sites were used for model validation 1.6 ± 0.4 mg/L.This corresponded to a 46% TOC removal.For DOC, the lake-water concentration was 2.6 ± 0.2 mg/L, and in well 30 the concentration was 1.5 ± 0.3 mg/L, corresponding to a 44% DOC removal.In the most distant groundwater well 2, the TOC concentration was 0.3 ± 0.05 mg/L and the DOC concentration was 0.3 ± 0.1 mg/L, corresponding to a 90% removal of both TOC and DOC.In all samples, TOC and DOC concentrations were above the detection limit (0.2 mg/L).
The lake DO level (Fig. 5d) was 11.0 ± 1.5 mg/L (average and standard deviation), which is close to saturation.As oxygen saturation is temperature-dependent, higher DO concentrations were observed during the winter season and lower concentrations in the summer.In the groundwater, DO concentrations showed spatial and seasonal variability.DO levels at the most intensive infiltration location, i.e., in well 30, varied from anoxic (Min.0.2 mg/L) to almost saturated (Max.10.8 mg/L).The DO concentration was below 0.3 mg/L from June to October and above 9.5 mg/L from December to February, the average value being 5.4 ± 5.0 mg/L.Similar changes with smaller variability were observed in groundwater well 40 with a DO concentration of 5.5 ± 3.8 mg/L.In other wells, DO levels varied from almost saturated (well 50) to almost anoxic (wells 2, 3, 4 and 5), but the amplitude of DO variation in all these wells was smaller than the variation in wells 30 and 40.Considering the entire data set, the highest DO concentrations were observed in well 50, where a high share (35%) of groundwater was formed from precipitation.Observed ΔTOC in all wells showed weak correlation to DO (r = 0.46, p = 0.002; Fig. 6a).
The Wilcoxon rank sum test (Fig. 7) demonstrated that TOC concentrations in all groundwater observation wells, and as DO concentration in all wells excluding well 50, were different from the observed concentrations in the lake.The significance level was 5%.The test showed no statistically significant changes in Fe concentrations from the lake to wells 40, 50 and 60.At the end of the flow path, no significant changes in the TOC or DO concentrations were observed between wells 2, 3 and 4. The similarity of the DO level in well 30 with almost all other wells was best explained by the large variation of DO in well 30.

Flow conditions and scaling to MAR sites
This study was designed to fill knowledge gaps in the long-term development of water quality in MAR, with a focus on NOM removal.The infiltrated fraction in groundwater was calculated from measured δ 18 O values in the lake and in the groundwater.As the seasonal variation in the isotopic composition of both endmembers was small compared to the difference between them, the calculated fractions of lake water and groundwater at each observation were considered to be precise.
Retention times estimated with the MODFLOW-based flow model and from the temperature breakthrough curve were consistent with each other near the lake shore.The assumption behind the 1D model for estimating the mixing ratios was justified on the main flow path, but the results were uncertain at the aquifer boundaries (wells 10 and 50).Due to the low share of lake infiltrate, well 10 was not sampled for later research.As the total discharge from the model area and the mixing ratio at several cells from the two inflow sources were known, detailed data of precipitation or catchment properties were not needed for the model setup.Mustafa et al. (2018) demonstrated that the input uncertainty in the catchment properties affects the model predictions and parameter distributions, and this is the dominant source of uncertainty in the groundwater flow prediction.To limit the uncertainties, the model was based on more accurate spring Fig. 7 Comparison of median TOC, DO, Fe, and Mn concentration between the observation sites.Results of the Wilcoxon rank sum test are shown between the data series from the lake and groundwater observation wells (numbered 2-60).The colors refer to similarity between the data series, i.e., the null hypothesis of the test could not be rejected at the 5% significance level discharge data and isotopic data, instead of inaccurate inflow estimations from the catchment or the lake bank.
The retention time from the lake to the first observation well 30 was 1-2 weeks, and to the most distant observation sites, well 2 and spring 1, around 1-1.5 years.The total retention time exceeded the retention time of a typical MAR system, which may vary from 5 to 100 days (Jokela et al. 2017); Kuehn and Mueller 2000.The estimated surface-water inflow rate of 870-1,265 m 3 /day corresponds to a small MAR site.In Finland, there are 26 enhanced infiltration sites with the production capacity of 350-100,000 m 3 /day (Jokela and Kallio 2015).Additionally 212 groundwater intake plants are estimated to receive, or to potentially receive, water from bank infiltration (Vienonen et al. 2012).This is around 10% of all Finnish groundwater intakes.However, many of the observed processes are universal, suggesting that the results can be generalized to wider scales of MAR designs, as long as the retention times and the quality of recharged water are comparable.

Surface-water quality
The infiltrating surface water (Fig. 5a-c) was characterized by moderate NOM levels as compared to reference studies.Meybeck (1993) reported a global estimate of the average DOC concentration in river waters to be 5.3 mg/L, but the variation is large.In Scandinavia, TOC concentrations between <0.2 to 30 mg/L are measured in lakes (Skjelkvåle et al. 2001).At the largest MAR plant in Finland, water is pre-treated to a TOC concentration of approximately 3 mg/L before infiltration (Niinikoski et al. 2016), but at several sites water with a TOC concentration of 6.5-11 mg/L is infiltrated in MAR (Jokela et al. 2017).In the Lake Iso Tiilijärvi water, 87% of organic carbon was in the dissolved form, which is a typical share in lake waters (Thurman 1985).The turbidity of lake water was low (Table S1 of the ESM), whereby mechanical clogging with suspended solids was not expected to occur in the infiltration.NO 3 − levels were below the detection limit (Table S1 of the ESM), indicating that NO 3 was not to buffer Fe or Mn dissolution when the DO concentration was low.

Removal of natural organic matter
The TOC concentration decreased 46% and the DOC concentration 44% during the first 7-15 days of transport after infiltration, and a 90% reduction of both TOC and DOC was reached by the end of the main flow path (Fig. 5a-c).These results are calculations from the actual measured values.The assumption is that the local groundwater has low TOC concentrations close to 0 mg/L, which reflects available operational groundwater quality data (TOC <1 mg/L, Finnish Environment Institute 2017).If the uncertainty of a 25% dilution at the end of the flow path with local groundwater is taken into account, the total TOC reduction would still be 80-90% in the infiltrated water.
The observed reduction values are on the same level or higher than values reported from existing MAR sites, summarized in Table 3.In these reference studies, DOC removal during the first 15 days after infiltration varied between 33 and 62% and the total removal between 37 and 100%.The large variation in removal efficiency is probably due to variable environmental conditions, aquifer properties and NOM quality (Maeng et al. 2011;Gross-Wittke et al. 2010;Abdelrady et al. 2018).As the reduction percentage at the studied Iso Tiilijärvi infiltration site was high, the results suggest that an aquifer can effectively remove NOM without losing its capacity over the long term, as long as the surface-water quality and the aquifer properties are in a range comparable to the studied conditions.The observed rapid decrease of TOC and DOC levels near the shore was similar to earlier constructed MAR sites (Bourg and Bertin 1993;Brugger et al. 2001;Diem et al. 2013;Frycklund 1995;Grünheid et al. 2005;Hiscock and Grischenk 2002;Kortelainen and Karhu 2006;Niinikoski et al. 2016).This decrease was due to a rapid removal of easily adsorbable and biodegradable NOM fractions, which likely were the hydrophobic, high molecular weight fractions (e.g., Lindroos et al. 2002).In all, 10-20% of NOM remained in groundwater even though the retention time at the Iso Tiilijärvi aquifer was estimated to be longer than 1 year.Statistically significant changes were not observed in the TOC concentrations between the last observation wells 2-5 (Fig. 7), which indicates that processes in soil were not able to remove all NOM, and that a residual fraction travelled through the system.According to earlier studies, the residual organic matter consists of hydrophilic, neutral, low molecular weight organics (e.g., Dietrich et al. 2013;Lindroos et al. 2002;Parsons et al. 2004).At Iso Tiilijärvi, NOM removal efficiency was not correlated to temperature (Jylhä-Ollila et al. 2019), which is in line with earlier results (e.g., Brugger et al. 2001;Diem et al. 2013;Massmann et al. 2006), whereas on the other hand, Abdelrady et al. (2018) reported a correlation between the NOM removal and temperature and Diem et al. (2013) suggested that the easily degradable NOM is decayed in soil at all temperatures within the first hours of infiltration.Rapid decay could explain the absence of seasonality in NOM removal even though bacterial activity is known to be temperature-dependent.
Adsorption of NOM has been the topic of many studies for several decades and the adsorptive nature of NOM on sandy materials is widely known (e.g., Lippold and Lippmann-Pipke 2014).
The sorptive behavior of NOM is dependent on the NOM structure, especially on its molecular weight (e.g., Dietrich et al. 2013), and other structural properties like aromaticity, carboxyl group content, and amino acid residues, that partly correlate with the molecular weight (e.g., Aiken et al. 1985;McKnight et al. 1992).The sorption of NOM is influenced by mineral surface properties, flow rate, and solution conditions such as NOM concentration, pH, ionic strength, and concentrations of multi-valent cations (Dietrich et al. 2013 and references therein).The main sorption mechanism is ligand exchange, with additionally Van-der-Waals attraction, hydrogen bonding and metal bridges (Lippold and Lippmann-Pipke 2014 and references therein).
If sorption or filtering without biodegradation were the dominating processes of NOM removal in MAR, signs of NOM accumulation in soil would be expected at the Iso Tiilijärvi site.These include elevated TOC and DOC concentrations in response to pumping, NOM break through to outflow springs as the sorption sites were fully occupied, lowered hydraulic conductivity due to clogging of soil pores, or NOM accumulation in the visual inspection of soil samples.As none of these accumulation signs were observed, the results support the assumption that biodegradation is the dominating process of NOM removal at our site, and the role of sorption is small.This is promising for the long-term sustainability of MAR.

Dissolved oxygen
Measured DO concentrations revealed that not only biodegradation of NOM in infiltrating surface water but also several other processes affected the redox conditions in the aquifer.Aerobic biodegradation leading to the mineralization of NOM to inorganic carbon can be described as follows: where CH 2 O represents organic matter and CO 2 carbon dioxide and other species of inorganic carbon.According to Eq. ( 7), one mole of O 2 is consumed to oxidize one mole of CH 2 O.This means that biodegradation of organic carbon at the concentration of 3 mg/L from the lake to the average level of ca.0.5 mg/L (ΔTOC -2.5 mg/L) observed in the aquifer would consume 6.5 mg of DO from the infiltrating water; thus, a DO decrease from 11.2 mg/L in the lake to 4.5 mg/L in the groundwater would be expected.The theoretical decrease in the DO level according to Eq. ( 7) is illustrated in Fig. 6a as a solid line.As observed in Fig. 6a, the relationship between DO and ΔTOC in the data did not follow the predicted line, and anoxic conditions in the aquifer were common.Overall, the correlation between NOM removal and DO was not strong (r = 0.46 p = 0.002), which indicates that other significant or even dominant processes affected the DO levels more than the degradation of organic carbon in the infiltrating water.
Earlier column tests and studies at bank infiltration sites (e.g., Diem et al. 2013;Greskowiak et al. 2006;Grischek and Paufler 2017;Henzler et al. 2016;Massmann et al. 2006;Sharma et al. 2012) have demonstrated an elevated DO consumption early along the infiltration flow path.The amount of DO consumption is temperature dependent, typically leading to anoxic conditions in the summer and aerobic conditions in the winter in the aquifer.The elevated oxygen consumption is explained by the microbiological activity in the lake bottom sediment (Diem et al. 2013).As could be expected, DO levels at Lake Iso Tiilijärvi were clearly correlated to temperature at the closest observation site from the lake bank, i.e., groundwater well 30 (r = −0.69,p = 0.001; Jylhä-Ollila et al. 2019).The seasonally observed DO depletion in the infiltrating water revealed that the degradation process of organic matter in the lake bottom sediment had a significant, even dominant impact on the infiltrating water quality.Anyhow, the decrease in DO levels at long retention times (174-440 days) could not be explained by the processes at the lake bank (Fig. 5d).In addition to the previous observations, a high DO content was observed in well 50 due to a high share (35%) of local groundwater.
In the studied area, Mn and Fe demonstrated a typical behavior of dissolution from mineral soil material in reducing conditions.No nitrogen compounds were present to buffer the reduction of Mn and Fe compounds.Mn and Fe correlated negatively to DO levels, and they were mainly dissolved in the aquifer when DO concentrations were below 3-4 mg/L (Fig. 6b-c).
As NOM levels in the surface water in this study were below the estimated average of rivers (Meybeck 1993), further work is required to validate the results in a research with high NOM concentrations in the infiltrating water.Additionally, more research is needed to quantify the processes affecting DO levels in the aquifer, as DO concentration is the key for predicting Mn and Fe dissolution.New analytical methods, e.g., FT-ICR-MS, could be applied to characterize the molecular structure of NOM in the molecular structure of NOM in different phases of infiltration, as it could help to predict the total NOM removal in MAR in variable conditions.

Conclusions
A natural lake-aquifer system was observed to continue removing NOM from infiltrating water for thousands of years, without losing its efficiency.In the studied aquifer, 44-46% of organic matter was removed within the first meters of infiltration, over a retention time of 7-15 days.The total removal was from 3.0 ± 0.20 mg/L (average and standard deviation) TOC to 0.3 ± 0.05 mg/L after a retention for 1 year, corresponding to a removal of 80-90%.The efficiency of organic matter removal in the natural system was on the same level as reported in earlier column tests and studies of MAR sites.A 10-20% residual concentration of nonbiodegradable, nonadsorbable fraction of organic matter was left in groundwater, which suggests that despite the long-term filtering capacity, processes in soil are not able to remove all organic matter from the infiltrating water.
No signs of accumulation of organic matter was observed in the aquifer.This strengthened the assumption that biodegradation was the main process in the removal of NOM, and the role of sorption and mechanical filtering were small or negligible, which is positive in terms of long-term sustainability of MAR.
At least four different processes in both the aquifer and the lake sediments affecting DO levels in the aquifer could be identified: (1) DO consumption due to biodegradation of NOM in the groundwater, (2) DO consumption of organic matter in the lake bed, (3) unexplained DO consumption observed in the wells located far from the lake bank, after a retention time of about 1 year and (4) DO addition due to mixing of groundwater formed from precipitation with lake infiltrate.These processes led to spatial and seasonal changes in the redox conditions.It is important to take DO into account when planning or modeling MAR, as Fe and Mn dissolution was clearly dependent on DO levels.In particular, oxygenconsuming biological processes in the lake sediment under summer temperatures had a negative impact on the final groundwater quality.At MAR sites, where DO is depleted and Fe and Mn are dissolved in water, changing the infiltration method from bank infiltration to enhanced infiltration, could improve the final water quality.

Fig. 1 a
Fig.1a Groundwater and surface-water sampling sites in the study area in Finland, b location of Lake Iso Tiilijärvi and the V-notch weir and land use around the lake.The map contains data from the National Land Survey ofFinland (2019)

Fig. 2
Fig. 2 Structure of the MODFLOW-based 1D flow model Figure4shows the calibration results of the heat transport model.The seasonal pattern of lake-water temperature was reflected with a delay in the temperature pattern in groundwater well 30.The model calibration yielded seepage velocities of 0.2 m/day in spring 2017 and 0.4 m/day in spring 2018.According to the model, the retention time from the lake to observation well 30 was 15 days in spring 2017 and 7 days in spring 2018.Thermal dispersivity was calibrated to the value of 0.1 m.The heat transport model was compared against the MODFLOW-based estimate of the retention time (Table2).The average of the retention times in the heat transport model was 29% larger than the average of the retention times in the MODFLOW-based model.As the differences between the models were small, the models produce a consistent picture of typical flow velocities in the aquifer.The sharp spring-time temperature changes in the lake and groundwater, observed in the measurements, were linked to the ice breakup timing of Lake Iso Tiilijärvi (Fig.4).Table2summarizes the calibration results of the MODFLOW-based flow and solute transport model in terms of the measured and modelled percentage of lake-water

Fig. 5
Fig. 5 Water-quality development in groundwater during natural bank infiltration.Boxplots for a total organic carbon (TOC), b dissolved organic carbon (DOC), c chemical oxygen demand (COD M ), d dissolved oxygen (DO), e Fe, f Mn, g electrical conductivity (EC) and

Fig. 6
Fig. 6 Scatter plots between DO and removal of organic matter as a ΔTOC (ΔTOC = TOC sample -TOC lakewater ), b iron and c manganese in groundwater wells 30, 40, 50, 60, 2, 3, 4 and 5.The line (see a) shows the path where the samples would be located if only oxygen consumption according to Eq. (6) is assumed to affect DO concentrations in groundwater.Observations from well 50 are colored in blue because of a low bank infiltrate fraction (65%)

Table 1
Parameter values of the thermal transport model

Table 2
Measured water temperature at well 30 and Lake Iso Tiilijärvi, and modelled temperature at well 30 produced by the calibrated heat transport model Measurement-based and modelled shares of lake-water recharge to the aquifer, and modelled retention times at different distances (sampling sites) from the bank.Modelled values were produced with the MODFLOW-based flow model

Table 3
Retention time and removal of organic matter at existing MAR sites with different raw water quality d Niinikoski et al. 2016 e Lindroos et al. 2002 f Jokela et al. 2017 g Maeng et al. 2011