Next Article in Journal
Correction: Mayer, A., et al. A Comprehensive Modelling Approach to Assess Water Use Efficiencies of Different Irrigation Management Options in Rice Irrigation Districts of Northern Italy. Water 2019, 11, 1833
Next Article in Special Issue
Do Land Use Changes Balance out Sediment Yields under Climate Change Predictions on the Sub-Basin Scale? The Carpathian Basin as an Example
Previous Article in Journal
HEM Impoundment—A Numerical Prediction Tool for the Water Framework Directive Assessment of Impounded River Reaches
Previous Article in Special Issue
Climate Change Impact on Surface Water and Groundwater Recharge in Northern Thailand
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Evaluation of the Influence of Farming Practices and Land Use on Groundwater Resources in a Coastal Multi-Aquifer System in Puck Region (Northern Poland)

1
Faculty of Civil and Environmental Engineering, Gdańsk University of Technology, ul. Narutowicza 11/12, 80-233 Gdańsk, Poland
2
Institute of Oceanology of the Polish Academy of Sciences, ul. Powstańców Warszawy 55, 81-712 Sopot, Poland
*
Author to whom correspondence should be addressed.
Water 2020, 12(4), 1042; https://doi.org/10.3390/w12041042
Submission received: 21 January 2020 / Revised: 9 March 2020 / Accepted: 2 April 2020 / Published: 7 April 2020

Abstract

:
This study focuses on the modeling of groundwater flow and nitrate transport in a multi-aquifer hydrosystem in northern Poland, adjacent to Puck Bay (Baltic sea). The main goal was to investigate how changes in land use and farming practices may affect groundwater recharge and submarine groundwater discharge (SGD) to the sea and the associated N-NO3 fluxes. An integrated modelling approach has been developed, which couples the SWAT hydrologic model, MODFLOW-NWT groundwater flow model, and MT3DMS transport model. Transient simulations were performed for a 10 y period, assuming 10 different scenarios of land use (farming, grassland, forest) and crop types. Both recharge and SGD showed a distinct pattern of seasonal time variability. In terms of the average flow rates, the effect of varying crop type was more significant than that of land use change, with the minimum recharge and SGD corresponding to winter wheat and the maximum for peas and potatoes. Nitrate loads were strongly affected by both land use and crop type, with minimum values obtained for grassland and maximum values for canola.

1. Introduction

The development of integrated surface–subsurface flow and transport models is often based on the coupling of existing computer codes, capable of simulating processes in various hydrological compartments. The combination of SWAT hydrologic model [1] with MODFLOW/MT3D groundwater flow and transport models [2,3] has been used by a number of authors over the last 20 years [4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24]. Such an approach enhances the capacities of both models because the groundwater component of SWAT is simplified and does not offer the accuracy required in many applications, especially involving complex hydrological conditions. On the other hand, SWAT is well-suited for computations of recharge rates and the associated nitrate loads varying in space and time over the watershed. Coupling between computer codes can be implemented as one-way sequential coupling, where spatially distributed recharge rates and nitrate loads from SWAT (and possibly also river stages) are passed to MODFLOW and MT3DMS, which calculate groundwater flow and transport independently from the SWAT groundwater subroutine. Alternatively, the groundwater procedure in SWAT can be replaced by the MODFLOW code, in principle allowing for a more consistent representation of surface–subsurface water interaction, at the cost of increased complexity of the computer code. Recently, a fully coupled SWAT-MODFLOW-RT3D software was developed in the form of a single executable program [15].
In previous works SWAT was coupled with MODFLOW for watersheds varying in size from a few km2 [24] to several thousand km2 [19]. However, there are only a few studies that used the combination of SWAT, MODFLOW, and MT3D to evaluate the quantity and quality of the submarine groundwater discharge (SGD) from land aquifers to coastal marine waters. SGD has received increasing attention of scientists due to its role played in water and nutrient cycles [25,26,27,28,29]. Galbiati et al. [7] developed an integrated model consisting of the SWAT, MODFLOW, MT3D, and stream water quality model QUAL2E in order to estimate nutrient load discharged to an Adriatic coastal lagoon. In their case, however, the nutrient influx was caused mainly by pumping from channels, not SGD, because the land surface was mostly below the sea level. More recently, [24] applied recharge data from SWAT as an input to the MODFLOW/MT3D model describing local groundwater flow in a small coastal watershed on Samoa. On the other hand, there is only limited experience in applying SWAT to complex multi-aquifer groundwater models [14,19] and, to the best of our knowledge, there were, until now, no attempts to use SWAT coupled with a multi-aquifer groundwater model to estimate SGD. The research reported here has been motivated by the need to predict nutrient loads in submarine groundwater discharge to the Puck Bay (northern Poland), which is a shallow part of the Baltic Sea, very vulnerable to eutrophication (Figure 1) [28,30,31,32]. Our study is a part of a larger project, which combines SWAT, MODFLOW, and MT3D with a hydrodynamic and ecological model of the Puck Bay (EcoPuckBay [33,34]). The integrated set of models can be used as a decision supporting tool and constitutes a part of the web-based service WaterPUCK (waterpuck.pl), which gives local stakeholders (administration and farmers) tools for predicting changes in the quality and quantity of local water resources. Below, we list the objectives of our work and show how they address the existing gaps of knowledge:
  • Development of a numerical model describing transient groundwater flow and nitrate transport in the multi-aquifer system adjacent to the Puck Bay. Earlier groundwater modeling studies in this area focused only on the rate of SGD, without contaminant transport [35,36,37].
  • Coupling of groundwater model with SWAT hydrological model to obtain detailed information on groundwater recharge rate and nitrate loads leached from soil. SWAT has already been used in a number of studies on the hydrological balance and nutrient transport in the Baltic region [38,39,40,41,42,43], including the Puck Bay region [44,45]. However, there has been no attempt to apply SWAT-MODFLOW coupling on the Polish Baltic coast.
  • Evaluation of the time variability of the groundwater recharge rate and SGD rate and the associated nitrate loads under different land use scenarios. While the impact of crop type and land use on water and nutrient fluxes has been extensively studied for different watersheds using SWAT, there have been no such investigations in the context of submarine groundwater discharge. In this work, we evaluate 10 scenarios, corresponding to six types of crops and four types of land use.

2. Materials and Methods

2.1. Study Area

The area under investigation is located in northern Poland (Figure 1). The land represents a typical young glacial landscape with relatively high relief. It consists of isolated fragments of a moraine plateau separated from each other with deeply cut ice marginal valleys. The entire area is covered with Quaternary deposits. The thickness of the sediments ranges from 40 m in the east to 95–100 m in the west. Such a significant difference is probably a result of the varied land surface or large range in altitude of the top of older layers [46]. The Quaternary deposits consist of moraine glacial till with sand and gravel layers, and glacifluvial or river sand and silty sand in the valleys. The dominant soil types are sandy loam (glacial till), sandy loam covered by loamy sand (weathered glacial till), sand (of glacifluvial origin), and peat (in larger river valleys).
Groundwater forms a complex system with several aquifer units (Figure 2). Shallow groundwater occurs in form of small perched aquifers and in sand lenses enclosed in moraine deposits (Qp in Figure 1). These aquifers are situated at depths up to 5–10 m below the ground level. They do not play an important role in water supply, but they are exploited locally by farms and summer houses. Two deeper Quaternary aquifers span most of the area: The upper (sub-moraine) aquifer (Q1) and the lower inter-moraine aquifer (Q2). They were formed in glacifluvial deposits (sand and gravel) separated with a layer of moraine till. These two large aquifers are hydraulically connected. On the major part of the area, they are confined (piezometric head reaches 45 m above mean sea level), except areas in the valleys where the till cover was eroded and the upper aquifer is exposed. In those places, the upper aquifer is unconfined and indirectly connected with surface waters. Shallow groundwater is recharged by the infiltration of rainwater. Deeper aquifers are influenced by inflow from the west and take the seepage water from the layers above. The entire hydrosystem is drained mainly by the Baltic Sea (Puck Bay), either directly via SGD or indirectly via streams and rivers. Aquifers Q1 and Q2 are the main source of water supply in the region.
Proximity to the Baltic Sea results in a specific marine climate of this region typically with moderate winters and mild summers. The annual temperature is about 7.4 °C and the average precipitation in the period of 2000–2009 was 620 mm/y.
The population of the area is about 25,000 people, of which 11,000 live in the largest town, Puck. Land use is divided between agricultural (51%), forests (28%), permanent meadows and grasslands (12%, mostly on peat soils in river valleys), and others, mostly low-density urban (9%). The main crops are (in order of decreasing importance): Winter wheat, winter triticale (cross-over between wheat and rye, typical for Poland), silage corn, canola, mixed spring cereals, potatoes, and peas (Pisum). Most of the farmers use crop rotation [47].

2.2. SWAT Model

A single SWAT model has been built for the investigated area, comprising the watersheds of 3 streams flowing to Puck Bay (Gizdepka, Płutnica, Błądzikowski Potok), as well as a small part of the watershed of Reda river, adjacent to its outflow to the bay (Figure 1). The model has been set-up using the QGIS SWAT interface. The area was divided into 17 subbasins and 353 hydrological response units (HRU—basic units for SWAT calculations). HRUs were defined based on the land use map and soil map. Screening tests showed that including ground slope as a criterion significantly increases the number of HRUs, without affecting the results of simulations. Weather data for the period 2000–2009 were obtained from the station in Żelistrzewo, in the SE part of the area (daily precipitation) and from Climate Forecast System Reanalysis [globalweather.tamu.edu] (daily values of humidity, wind, temperature, and radiation). Potential evapotranspiration was estimated with the Penman–Monteith model. Soil parameters (available water content, vertical hydraulic conductivity, etc.) were taken as average values in each textural class, as reported in the Polish literature [48]. The agricultural management practices for crops and grassland were based on the results of a poll carried out among local farmers [47]. In all agricultural areas, the same pattern of crop rotation was specified, including silage corn, canola, and winter wheat. Similarly, all grassland was assigned the same set of practices (permanent plant cover, Fescue species, hay cut twice a year). Forests were modeled as permanent plant cover. Pine was selected as the representative species.
Due to the lack of sufficient data on stream discharges, the calibration of SWAT has been performed by a manual trial-and-error procedure, using the available data on the components of the hydrological budget and crop yield (soft calibration). The SWAT model was calibrated in parallel with the steady-state groundwater flow model, which provided an additional constraint in the calibration process. As the first objective for SWAT calibration, we set the limits of the calculated groundwater recharge in any HRU from 3% to 30% of precipitation, which is commonly accepted for the climate and soil conditions in Poland [49,50]. The lower limit corresponded to HRU’s on peat and heavy clay soils, while the upper limit - to clean glaci-fluvial sands. As the second objective, we set the value of the average yearly evapotranspiration, which was reported to be in the range from 450 to 495 mm for a location in Poland with very similar climate and land use conditions (629 mm precipitation, canola with different fertilization schemes) [51]. The third objective was the annual biomass production in the forests, which should be close to 7 metric tons of dry mass per hectare per year [54,55].
In the calibration process, the following parameters were adjusted: SCN curve numbers (increased), maximum canopy storage (CANMX, increased in forests), and REVAP (increased in forests), as well as the parameters describing pine growth in the “plant.dat” file. A satisfactory agreement was obtained between the reference values (calibration goals) and the results obtained from the calibrated model (Table 1). In the same table, we report additional validation checks, performed on model outputs that were not used previously in the calibration procedure. The yield of basic crop types (winter wheat, canola, and silage corn) and the biomass harvested as hay on grasslands were compared to the data obtained from local farmers, from the statistics for the Pomeranian region in Poland and reported in the literature (all values given in Table 1 represent dry mass). The model predictions are quite accurate for canola, while there is some underestimation of yield for silage corn and some overestimation for winter wheat. Moreover, both the total runoff and the surface runoff/total runoff ratio are in agreement with previous studies on the hydrology of the region. Overall, the performance of the model versus the available data can be considered satisfactory.
In groundwater simulations described below, we used, as input data, the values of groundwater recharge and N-NO3 load reaching the groundwater table, obtained from SWAT simulations for the period 2000–2009, with a 2 y warm-up period (1998–1999). Note that SWAT uses a simplified procedure to describe the movement of water and nitrate through the vadose zone between the bottom of the soil profile and the groundwater table [1]. The resulting estimations of the vadose zone time lag may be inaccurate, but at the current stage of model development, this simplification was accepted.

2.3. Groundwater Flow Model

Groundwater flow simulations were carried out with MODFLOW-NWT numerical code [57]. The model was set-up with GMS 10.4 [58] and ModelMuse [59] graphical interfaces, while the FloPy [60] library was used to prepare data for transient simulations.
In line with the finite difference discretization of MODFLOW, we used a rectangular grid with 6 layers, 410 rows, and 296 columns. The total number of cells was 728,160, of which 340,496 were active. The spacing in horizontal directions was uniform and equal to 50 m in order to capture the topographic features and minimize discretization errors. The 6 layers were defined as follows: (I) Sandy, perched aquifers, occurring locally (Qp), (II) low-permeability till deposits, (III) sub-moraine aquifer in glacifluvial sand (Q1), (IV) low permeable till layer, (V) inter-moraine aquifer in glacifluvial sand (Q2), (VI) low-permeability layer representing Quarternary glacial till and Neogene clay and silt. The range of hydraulic conductivity values for each layer is presented in Table 2.
In all steady and transient simulations, we used the same set of boundary conditions and source terms, assumed to be constant in time, except for the recharge rate. The sea boundary was represented by a constant head boundary in the Q1 aquifer (direct hydraulic contact with water in Puck Bay) and by a general head boundary (GHB, third type) in the Q2 aquifer, which represented the outflow through a distant outcrop at the sea bottom (the discharge zones are situated in the bottom of the Puck Bay, at a distance of 2–3 km from the coast). The GHB boundary was also used to simulate lateral flow across the land part of the boundary in Q1, Q2, and those perched aquifers (Qp), which were intersected by the model boundary. Wells were assumed to be constantly in operation and the pumping rates were identified based on the available reports from groundwater administration. Rivers were represented by the third-type boundary conditions (RIVER package). The river stages were assumed to be constant in time and represented average values from measurements and earlier reports.
The results obtained from SWAT were processed by scripts written in the Python language using the FloPy library, in order to set the recharge values in the MODFLOW model. In steady-state simulations, we used the average values in each HRU from the period of 2000–2009. For transient simulations, monthly averages were calculated for each HRU from SWAT daily results. We followed the approach described by [15] to transfer HRU-based SWAT results to the grid-based MODFLOW input.
Calibration of the groundwater flow model was based on the steady-state solution, representing average conditions in the period of 2000–2009. We used groundwater head measurements from 95 points in the area. A good agreement was obtained (mean error 0.29 m, mean absolute error 1.10 m, and standard error 1.79 m). The comparison of calculated and measured groundwater head values for all three aquifers (Qp, Q1, and Q2) is shown in Figure 3. In order to validate the model, a number of simulations based on historical records of pumping tests were performed. They confirmed that the model was calibrated adequately.

2.4. Nitrate Transport Model

Based on the groundwater flow model developed in MODFLOW, a simulation of N-NO3 transport was carried out. MT3DMS numerical code was applied to solve the advection–reaction–dispersion equation. Transient calculations were performed with the third-order total-variational diminishing (TVD) numerical method for the advective term, while for steady-state solution, the standard finite difference discretization was applied.
Porosity and longitudinal dispersivity were assigned depending on the soil material. The values were selected based on the literature. Longitudinal dispersivity was assumed to be 1.5 m in permeable layers [61] and 0.75 m in the case of low-permeability layers [62]. Values of porosity were assigned equal to 0.2 in sand and gravel and 0.1 in glacial till. Denitrification of nitrates was described with the first-order kinetic reaction [63,64,65,66]. We used a constant reaction rate of biodegradation equal to 1 × 10−5 1/h, which was established by model calibration.
In the model, we considered only diffuse sources of nitrates. Concentrations of N-NO3 in water reaching the groundwater table were calculated from the results of SWAT simulation for the years of 2000–2009 for each HRU. As in the case of groundwater recharge, the transfer of HRU-based SWAT data to grid-based MT3DMS data was carried out with the aid of Python scripts, making use of the FloPy library. In the areas of the groundwater model not covered by SWAT, we used a concentration value representing an average taken over the whole SWAT model area in the considered time period.
The model of nitrate transport was calibrated manually by adjusting the reaction rate in order to obtain a satisfactory agreement between the calculated N-NO3 concentrations and the results of chemical analysis of the samples taken from the aquifers. Calibration was based on a steady-state flow and transport simulation, using average values of recharge and nitrate loads from the period 2000–2009 and the current land use and agricultural practices. The range of N-NO3 concentrations in layer Q1 obtained from the numerical simulation was 3.8 × 10−6–62.88 mg/dm3, while the range of measured concentrations was <0.02–20.78 mg/dm3 (0.02 is the detection limit). The values obtained from numerical modeling for layer Q2 were from 0.0012 to 1 mg/dm3, while the measured values ranged between <0.02 and 0.23 mg/dm3.

2.5. Scenarios for Transient Flow Simulations

In order to investigate the effects of different land use on the time variability of groundwater recharge and groundwater discharge to the Puck Bay, and the associated N-NO3 loads, we defined 10 scenarios (see Table 3). Scenario S1 (baseline) corresponds to the current land use and agricultural practices. In agricultural areas, crop rotation is used with alternating growing of winter wheat, silage corn, and canola, which are the main crops in the area.
In scenarios S2–S7, we assumed that the land use structure corresponds to the current state, but there is only one type of crop grown on agricultural lands, without any rotation. The crop types are as follows: S2: Winter wheat (also representing winter triticale), S3: Silage corn, S4: Canola, S5: Mixture of spring cereals (represented by barley), S6: Potatoes, S7: Peas (Pisum).
Scenarios S8–S10 represent major shifts in land use structure: S8: Whole area covered by agriculture, S9: Whole area covered by meadows, S10: Whole area covered by forest (in each case, we kept the current location of urban areas). While no such large changes are predicted in the near future, these scenarios can be considered as limit cases, defining the maximum range of expected variations in groundwater fluxes and nitrate loads.

3. Results and Discussion

3.1. Steady-State Flow Simulation

The components of the steady-state groundwater budget are shown in Figure 2. According to the model results, the main components of the budget are recharge from precipitation, external lateral inflow to Q1 and Q2, exchange with streams, and submarine groundwater discharge.
The considered groundwater system is typical for young glacial areas, with flow taking part in a multiple stacked aquifer, partly connected with each other. A typical feature of such systems is significant vertical seepage through glacial till layers, separating the aquifers. Even though glacial till is often considered a weakly permeable deposit, it plays an important role in the circulation of groundwater on glacial areas. These results are consistent with previous research, e.g., [50].
Submarine groundwater discharge constitutes a major component in groundwater budget. The amount of water discharged from the lower aquifer Q2 is approximately twice as large as in the case of the upper aquifer Q1. This can be explained by significant lateral inflow received by Q2 aquifer from the western boundary. In our model, the total groundwater discharge to the sea from Q1 and Q2 is 1286 m3/h, which corresponds to 92 m3/h per kilometer of the coastline. The flow rate varies from 58 m3/h/km in the upland areas to 152 m3/h/km in the river valleys. These values are in a good agreement with the results of previous studies related to groundwater outflow to the Puck Bay [35,36,37,67]. Ref. [35] reported the SGD rate in the same area ranging from 21 to 165 m3/h/km, with an average of 98 m3/h/km, while [37] estimated the SGD rate in the Reda river valley as 165 m3/h/km.

3.2. Transient Simulations of Flow and Nitrogen Transport

A summary of the main results obtained from transient flow simulations is presented in Table 4. We report the values of groundwater discharge and the corresponding load of N-NO3 leached from soil to groundwater, as well as the rate of groundwater discharge to Puck Bay and the associated loads of N-NO3, separately from Q1 and Q2 aquifers. All values represent averages for the whole simulation period (120 months).
Comparison of the average values of recharge and SGD for different scenarios indicates that the effect of changing the crop type within the current agricultural area (S2–S7) can be more pronounced than large-scale changes in land use (S8–S10). SWAT simulations predicted the minimum recharge for winter wheat (S2: 62 mm/y) and the maximum for peas (S7: 92 mm/y), which corresponds to a relative difference of about 50%. In contrast, the difference between land use types is smaller, with the maximum value of 85 mm/y corresponding to grassland (S9) and the minimum value of 71 mm/y to the forest (S10). The average SGD flow rates for various scenarios show a similar variability pattern as the values of recharge. The relative difference between maximum value (S7) and minimum value (S2) is almost 14% in the Q1 aquifer and 11% in the Q2 aquifer.
The yearly nitrate loads obtained from SWAT for various crops (Table 3) are more varied than recharge values, ranging from 7.4 kg/ha/y for winter wheat to 30.5 kg/ha/y for canola. While in the case of winter wheat, the minimum load corresponds to minimum recharge, the opposite is not true for canola. Therefore, there is no simple relationship between the amount of recharge and nitrate load to groundwater, because the latter one strongly depends on crop type and agricultural practices. Overall, the amount of N leached to groundwater according to SWAT simulations is similar to the values reported in the literature for Poland [68] (5–15 kgN/ha/y for grasslands, 5–10 for forests, 20–45 for crops).
Considering nitrate loads entering Puck Bay from SGD, one can see that the relative magnitude of the load for various scenarios corresponds to the relative magnitude of nitrate leaching from soil. Variation in the N-NO3 load carried by SGD between scenarios is much larger that the variation in the SGD flow rate. The major part of the load (about 90%) is discharged from the upper aquifer Q1, despite the fact that its flow rate is about two times smaller than in Q2. This can be easily explained, as Q1 is closer to the ground surface and directly receives NO3 leached from soil.
The results of transient simulations show significant time variation of the components of the groundwater budget and the associated nitrate loads. In Figure 4, Figure 5, and Figure 6, we show average monthly water and nitrogen fluxes in recharge and SGD, for scenarios S1 (baseline), S2 (winter wheat: Low recharge and nitrate load), and S3 (silage corn: High recharge and nitrate load). In most cases, a consistent pattern of seasonal changes can be distinguished, with maximum values in late winter/early spring and minimum values in early autumn. Such a variability is typical for Polish climatic conditions—the maximum corresponds to the period of thawing snow and the minimum to the late vegetation season. An exception can be observed in the very dry years of 2005 and 2006, for which there is no early spring maximum in 2006 in S1 and S2, due to higher water demands of plants in these scenarios.
Seasonal changes in SGD (Figure 4, Figure 5 and Figure 6) correspond approximately to seasonal changes in groundwater recharge; however, the relative differences in monthly values are significantly smaller. This can be seen particularly in the discharge from Q2, which is less affected by variations in recharge and is driven to a larger extent by lateral inflow from the western boundary. The relative fluctuations of SGD from Q1 reach 120% of the lowest value, while in Q2, they do not exceed 31% of the lowest value.
Concentrations of N-NO3 in SGD flux obtained from numerical simulations were compared to the results of chemical analysis of groundwater samples (Table 4). In both cases, the values are similar considering that the SGD concentrations from numerical modeling represent the entire area of investigations, and field research was conducted at points distributed across the site.
As shown in scenarios S8–S10, nitrate loads leaching from soil and transported by SGD are much more sensitive to land use changes than recharge and SGD flow rates. Our simulations show that converting forests and meadows to agriculture (S8) results in nitrate loads increased by about 75% with respect to the baseline scenario. This is consistent with the increase in agricultural area from 51% to 91% of the watershed. The large nitrate load is due to growing corn and canola, as part of crop rotation, assumed in the scenario.
Conversely, if crops and forests are replaced by grassland (S9), the nitrate loads are reduced by about six times compared to the baseline scenario. Grasslands are characterized by a high uptake of nitrogen, promoted by cutting hay twice a year, and the fertilized doses are also smaller.
The introduction of forests to the whole area (S10) also leads to the reduction of nitrate loads, albeit much less significant than in scenario S9. This result is caused by a smaller uptake of nitrogen in the grown-up forest as compared to growing grass. Even though there is no fertilization in the forest, a significant amount of nitrogen is released from peat soils. In other scenarios, this effect is offset with higher yields of crops and hay computed for peat areas. Further research is needed to quantify the actual rate of NO3 leaching from peat soil in the considered area.
The values of recharge and N-NO3 loads obtained from SWAT show significant spatial variability. Figure 7 presents the spatial distribution of values for a month with relatively high recharge and nitrate load (March 2007). The amount of recharge depends on the soil type and land use. Larger values are observed in sandy soils, mostly in the south-west part of the area. Lowest values occur predominantly in peat-covered river valleys. There is no simple relationship between recharge and nitrate load. For example, in the forests in the SW part of the area, nitrate leaching is limited, despite the moderate-to-high recharge. In contrast, in the southernmost part of the model domain, peat soil with a high content of organic matter releases large amounts of nitrate from mineralization of organic nitrogen, even though the recharge rate is small.

4. Conclusions

SWAT, MODFLOW-NWT, and MT3DMS computer codes were sequentially coupled, in order to develop a model of transient groundwater flow and nitrate transport in a coastal multi-aquifer hydrosystem. The SWAT model allowed the capture of effects of hypothetical changes in land use and crop type on groundwater fluxes (including SGD) and the associated nitrate loads. The main findings can be summarized as follows:
  • Groundwater recharge, SGD, and the corresponding nitrate loads show a distinct time variable pattern, with maximum recharge rates and NO3 leaching in late winter/early spring.
  • The average values of recharge and SGD fluxes are influenced more significantly by crop type grown on farmlands than by the changes in land use. The maximum relative difference between the 10 y average of SGD flux between different scenarios did not exceed 12%. In contrast, nitrate leaching from soil and nitrate transport via SGD shows a larger variability, strongly depending on crop type and land use.
  • The lowest N-NO3 load in SGD occurred for the hypothetical scenario with all land converted to grassland, and it was three times smaller than the largest load, corresponding to converting all land to growing crops.
While we believe that the general trends observed in our simulations are realistic, we are also aware of the limitations of the presented approach. Sequential coupling applied in this work does not allow for a consistent representation of aquifer–stream interactions in SWAT and MODFLOW models. In addition, the simplified representation of water and nitrate movement in the vadose zone between the root zone and groundwater table used by SWAT may lead to inaccuracies in the calculated time variability of contaminant transport, especially in complex geological conditions of the investigated area [69,70,71]. Further research will focus on an improved description of these processes in our model.

Author Contributions

Conceptualization, A.S., B.J.-S. and L.D.-G.; Methodology, A.Sz., B.J.-S. and M.P.-C.; Investigation, A.S., B.J.-S., D.P. and A.G.-K.; Writing—Original Draft Preparation, A.Sz., A.G.-K. and D.P.; Funding acquisition, L.D.-G. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the National Centre for Research and Development, Poland, in the framework of the project BIOSTRATEG3/343927/3/NCBR/2017 "Modelling of the impact of the agricultural holdings and land-use structure on the quality of inland and coastal waters of the Baltic Sea set up on the example of the Municipality of Puck region - Integrated info-prediction Web Service WaterPUCK", BIOSTRATEG Programme.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Neitsch, S.L.; Arnold, J.G.; Kiniry, J.R.; Williams, J.R. Soil and Water Assessment Tool Theoretical Documentation Version 2009; Texas Water Resources Institute: College Station, TX, USA, 2011. [Google Scholar]
  2. Harbaugh, A.W. MODFLOW-2005, the US Geological Survey Modular Ground-Water Model: The Ground-Water Flow Process; US Department of the Interior, US Geological Survey: Reston, VA, USA, 2005; pp. A6–A16.
  3. Zheng, C.; Wang, P.P. MT3DMS: A Modular Three-Dimensional Multispecies Transport Model for Simulation of Advection, Dispersion, and Chemical Reactions of Contaminants in Groundwater Systems; Documentation and User’s Guide; Alabama University: Tuscaloosa, AL, USA, 1999. [Google Scholar]
  4. Perkins, S.P.; Sophocleous, M. Development of a comprehensive watershed model applied to study stream yield under drought conditions. Groundwater 1999, 37, 418–426. [Google Scholar] [CrossRef]
  5. Sophocleous, M.A.; Koelliker, J.K.; Govindaraju, R.S.; Birdie, T.; Ramireddygari, S.R.; Perkins, S.P. Integrated numerical modeling for basin-wide water management: The case of the Rattlesnake Creek basin in south-central Kansas. J. Hydrol. 1999, 214, 179–196. [Google Scholar] [CrossRef]
  6. Conan, C.; Bouraoui, F.; Turpin, N.; de Marsily, G.; Bidoglio, G. Modeling flow and nitrate fate at catchment scale in Brittany (France). J. Environ. Qual. 2003, 32, 2026–2032. [Google Scholar] [CrossRef] [PubMed]
  7. Galbiati, L.; Bouraoui, F.; Elorza, F.J.; Bidoglio, G. Modeling diffuse pollution loading into a Mediterranean lagoon: Development and application of an integrated surface–subsurface model tool. Ecol. Model. 2006, 193, 4–18. [Google Scholar] [CrossRef]
  8. Kim, N.W.; Chung, I.M.; Won, Y.S.; Arnold, J.G. Development and application of the integrated SWAT–MODFLOW model. J. Hydrol. 2008, 356, 1–16. [Google Scholar] [CrossRef]
  9. Chung, I.M.; Kim, N.W.; Lee, J.; Sophocleous, M. Assessing distributed groundwater recharge rate using integrated surface water-groundwater modelling: Application to Mihocheon watershed, South Korea. Hydrogeol. J. 2010, 18, 1253–1264. [Google Scholar] [CrossRef]
  10. Lin, H.T.; Ke, K.Y.; Tan, Y.C.; Wu, S.C.; Hsu, G.; Chen, P.C.; Fang, S.T. Estimating pumping rates and identifying potential recharge zones for groundwater management in multi-aquifers system. Water Resour. Manag. 2013, 27, 3293–3306. [Google Scholar] [CrossRef]
  11. Pisinaras, V.; Petalas, C.; Tsihrintzis, V.A.; Karatzas, G.P. Integrated modeling as a decision-aiding tool for groundwater management in a Mediterranean agricultural watershed. Hydrol. Process. 2013, 27, 1973–1987. [Google Scholar] [CrossRef]
  12. Ke, K.Y. Application of an integrated surface water-groundwater model to multi-aquifers modeling in Choushui River alluvial fan, Taiwan. Hydrol. Process. 2014, 28, 1409–1421. [Google Scholar] [CrossRef]
  13. Guzman, J.A.; Moriasi, D.N.; Gowda, P.H.; Steiner, J.L.; Starks, P.J.; Arnold, J.G.; Srinivasan, R. A model integration framework for linking SWAT and MODFLOW. Environ. Model. Softw. 2015, 73, 103–116. [Google Scholar] [CrossRef]
  14. Pulido-Velazquez, M.; Peña-Haro, S.; García-Prats, A.; Mocholi-Almudever, A.F.; Henríquez-Dole, L.; Macian-Sorribes, H.; Lopez-Nicolas, A. Integrated assessment of the impact of climate and land use changes on groundwater quantity and quality in the Mancha Oriental system (Spain). Hydrol. Earth Syst. Sci. 2015, 19, 1677–1693. [Google Scholar] [CrossRef] [Green Version]
  15. Bailey, R.T.; Wible, T.C.; Arabi, M.; Records, R.M.; Ditty, J. Assessing regional-scale spatio-temporal patterns of groundwater–surface water interactions using a coupled SWAT-MODFLOW model. Hydrol. Process. 2016, 30, 4420–4433. [Google Scholar] [CrossRef]
  16. Ehtiat, M.; Mousavi, S.J.; Vaghefi, S.A.; Ghaheri, A. Analysis of recharge conceptualization in inverse groundwater modelling. Hydrol. Sci. J. 2016, 61, 2789–2801. [Google Scholar] [CrossRef] [Green Version]
  17. Ehtiat, M.; Mousavi, S.J.; Srinivasan, R. Groundwater modeling under variable operating conditions using SWAT, MODFLOW and MT3DMS: A catchment scale approach to water resources management. Water Resour. Manag. 2018, 32, 1631–1649. [Google Scholar] [CrossRef]
  18. Aliyari, F.; Bailey, R.T.; Tasdighi, A.; Dozier, A.; Arabi, M.; Zeiler, K. Coupled SWAT-MODFLOW model for large-scale mixed agro-urban river basins. Environ. Model. Softw. 2019, 115, 200–210. [Google Scholar] [CrossRef]
  19. Chunn, D.; Faramarzi, M.; Smerdon, B.; Alessi, D.S. Application of an Integrated SWAT–MODFLOW Model to Evaluate Potential Impacts of Climate Change and Water Withdrawals on Groundwater–Surface Water Interactions in West-Central Alberta. Water 2019, 11, 110. [Google Scholar] [CrossRef] [Green Version]
  20. Gao, F.; Feng, G.; Han, M.; Dash, P.; Jenkins, J.; Liu, C. Assessment of Surface Water Resources in the Big Sunflower River Watershed Using Coupled SWAT–MODFLOW Model. Water 2019, 11, 528. [Google Scholar] [CrossRef] [Green Version]
  21. Semiromi, M.T.; Koch, M. Analysis of spatio-temporal variability of surface–groundwater interactions in the Gharehsoo river basin, Iran, using a coupled SWAT-MODFLOW model. Environ. Earth Sci. 2019, 78, 201. [Google Scholar] [CrossRef]
  22. Wei, X.; Bailey, R.T. Assessment of System Responses in Intensively Irrigated Stream–Aquifer Systems Using SWAT-MODFLOW. Water 2019, 11, 1576. [Google Scholar] [CrossRef] [Green Version]
  23. Wei, X.; Bailey, R.T.; Records, R.M.; Wible, T.C.; Arabi, M. Comprehensive simulation of nitrate transport in coupled surface-subsurface hydrologic systems using the linked SWAT-MODFLOW-RT3D model. Environ. Model. Softw. 2019, 122, 104242. [Google Scholar] [CrossRef]
  24. Welch, E.M.; Dulai, H.; El-Kadi, A.; Shuler, C. Submarine groundwater discharge and stream baseflow sustain pesticide and nutrient fluxes in Faga’alu Bay, American Samoa. Front. Environ. Sci. 2019, 7, 162. [Google Scholar] [CrossRef] [Green Version]
  25. Kim, G.; Swarzenski, P.W. Submarine groundwater discharge (SGD) and associated nutrient fluxes to the coastal ocean. In Carbon and Nutrient Fluxes in Continental Margins; Springer: Berlin/Heidelberg, Germany, 2010; pp. 529–538. [Google Scholar]
  26. Kotwicki, L.; Grzelak, K.; Czub, M.; Dellwig, O.; Gentz, T.; Szymczycha, B.; Böttcher, M.E. Submarine groundwater discharge to the Baltic coastal zone: Impacts on the meiofaunal community. J. Mar. Syst. 2014, 129, 118–126. [Google Scholar] [CrossRef] [Green Version]
  27. Lipka, M.; Böttcher, M.E.; Wu, Z.; Sültenfuß, J.; Jenner, A.K.; Westphal, J.; Dellwig, O.; Escher, P.; Schmiedinger, I.; Winde, V.; et al. Ferruginous groundwaters as a source of P, Fe, and DIC for coastal waters of the southern Baltic Sea:(Isotope) hydrobiogeochemistry and the role of an iron curtain. In E3S Web of Conferences; EDP Sciences: Les Ulis, France, 2018; Volume 54, p. 00019. [Google Scholar]
  28. Kłostowska, Ż.; Szymczycha, B.; Lengier, M.; Zarzeczańska, D.; Dzierzbicka-Głowacka, L. Hydrogeochemistry and magnitude of SGD in the Bay of Puck, southern Baltic Sea. Oceanologia 2019. [Google Scholar] [CrossRef]
  29. Jiao, J.; Post, V. Coastal Hydrogeology; Cambridge University Press: Cambridge, UK, 2019. [Google Scholar]
  30. Wojciechowska, E.; Nawrot, N.; Matej-Łukowicz, K.; Gajewska, M.; Obarska-Pempkowiak, H. Seasonal changes of the concentrations of mineral forms of nitrogen and phosphorus in watercourses in the agricultural catchment area (Bay of Puck, Baltic Sea, Poland). Water Supply 2019, 19, 986–994. [Google Scholar] [CrossRef]
  31. Wojciechowska, E.; Pietrzak, S.; Matej-Łukowicz, K.; Nawrot, N.; Zima, P.; Kalinowska, D.; Wielgat, P.; Obarska-Pempkowiak, H.; Gajewska, M.; Dembska, G.; et al. Nutrient loss from three small-size watersheds in the southern Baltic Sea in relation to agricultural practices and policy. J. Environ. Manag. 2019, 252, 109637. [Google Scholar] [CrossRef]
  32. Szymczycha, B.; Kłostowska, Ż.; Lengier, M.; Dzierzbicka-Głowacka, L. Significance of nutrient fluxes via submarine groundwater discharge in the Bay of Puck, southern Baltic Sea. Oceanologia 2020. [Google Scholar] [CrossRef]
  33. Dzierzbicka-Głowacka, L.; Janecki, M.; Dybowski, D.; Szymczycha, B.; Obarska-Pempkowiak, H.; Wojciechowska, E.; Zima, P.; Pietrzak, S.; Pazikowska-Sapota, G.; Jaworska-Szulc, B.; et al. A new approach for investigating the impact of pesticides and nutrient flux from agricultural holdings and land-use structures on the coastal waters of the Baltic Sea. Pol. J. Environ. Stud. 2019, 28, 2531–2539. [Google Scholar] [CrossRef]
  34. Dybowski, D.; Jakacki, J.; Janecki, M.; Nowicki, A.; Rak, D.; Dzierzbicka-Glowacka, L. High-Resolution Ecosystem Model of the Puck Bay (Southern Baltic Sea)—Hydrodynamic Component Evaluation. Water 2019, 11, 2057. [Google Scholar] [CrossRef] [Green Version]
  35. Piekarek-Jankowska, H. Zatoka Pucka jako obszar drenażu wód podziemnych (The Bay of Puck as A Groundwater Drainage Area); University of Gdańsk: Gdańsk, Poland, 1994; p. 104. (In Polish) [Google Scholar]
  36. Kryza, J.; Kryza, H. Analityczna i modelowa ocena bezpośredniego dopływu podziemnego do Bałtyku na terytorium Polski. Geologos 2006, 10, 153–165. [Google Scholar]
  37. Jaworska-Szulc, B. Groundwater flow modelling of multi-aquifer systems for regional resources evaluation: The Gdansk hydogeological system, Poland. Hydrogeol. J. 2009, 17, 1521–1542. [Google Scholar] [CrossRef]
  38. Piniewski, M.; Szcześniak, M.; Kardel, I.; Berezowski, T.; Okruszko, T.; Srinivasan, R.; Vikhamar Schuler, D.; Kundzewicz, Z.W. Hydrological modelling of the Vistula and Odra river basins using SWAT. Hydrol. Sci. J. 2017, 62, 1266–1289. [Google Scholar] [CrossRef]
  39. Thodsen, H.; Farkas, C.; Chormanski, J.; Trolle, D.; Blicher-Mathiesen, G.; Grant, R.; Engebretsen, A.; Kardel, I.; Andersen, H.E. Modelling nutrient load changes from fertilizer application scenarios in six catchments around the Baltic Sea. Agriculture 2017, 7, 41. [Google Scholar] [CrossRef] [Green Version]
  40. Tamm, O.; Maasikamäe, S.; Padari, A.; Tamm, T. Modelling the effects of land use and climate change on the water resources in the eastern Baltic Sea region using the SWAT model. Catena 2018, 167, 78–89. [Google Scholar] [CrossRef]
  41. Trolle, D.; Nielsen, A.; Andersen, H.E.; Thodsen, H.; Olesen, J.E.; Børgesen, C.D.; Refsgaard, J.C.; Sonnenborg, T.O.; Karlsson, I.B.; Christensen, J.P.; et al. Effects of changes in land use and climate on aquatic ecosystems: Coupling of models and decomposition of uncertainties. Sci. Total Environ. 2019, 657, 627–633. [Google Scholar] [CrossRef]
  42. Olesen, J.E.; Børgesen, C.D.; Hashemi, F.; Jabloun, M.; Bar-Michalczyk, D.; Wachniew, P.; Zurek, A.J.; Bartosova, A.; Bosshard, T.; Hansen, A.L.; et al. Nitrate leaching losses from two Baltic Sea catchments under scenarios of changes in land use, land management and climate. Ambio 2019, 48, 1252–1263. [Google Scholar] [CrossRef]
  43. Terskii, P.; Kuleshov, A.; Chalov, S.; Terskaia, A.; Pluntke, T.; Karthe, D. Assessment of Water Balance for Russian Subcatchment of Western Dvina River Using SWAT Model. Front. Earth Sci. 2019, 7, 241. [Google Scholar] [CrossRef]
  44. Marcinkowski, P.; Piniewski, M.; Kardel, I.; Gielczewski, M.; Okruszko, T. Modelling of discharge, nitrate and phosphate loads from the Reda catchment to the Puck Lagoon using SWAT. Annals of Warsaw University of Life Sciences-SGGW. Land Reclam. 2013, 45, 2. [Google Scholar]
  45. Piniewski, M.; Kardel, I.; Giełczewski, M.; Marcinkowski, P.; Okruszko, T. Climate change and agricultural development: Adapting Polish agriculture to reduce future nutrient loads in a coastal watershed. Ambio 2014, 43, 644–660. [Google Scholar] [CrossRef] [Green Version]
  46. Jereczek-Korzeniowska, K.; Jegliński, W. The Late Glacial and holocene development of valley network in the puck morainic plateau. Geologija 2011, 53, 10–14. [Google Scholar] [CrossRef] [Green Version]
  47. Dzierzbicka-Głowacka, L.; Pietrzak, S.; Dybowski, D.; Białoskórski, M.; Marcinkowski, T.; Rossa, L.; Urbaniak, M.; Majewska, Z.; Juszkowska, D.; Nawalany, P.; et al. Impact of agricultural farms on the environment of the Puck Commune: Integrated agriculture calculator—CalcGosPuck. PeerJ 2019, 7, e6478. [Google Scholar]
  48. Zawadzki, S. Gleboznawstwo [Soil Science]; Państwowe Wydawnictwa Rolnicze i Leśne: Warszawa, Poland, 1999. [Google Scholar]
  49. Duda, R.; Winid, B.; Zdechlik, R.; Stępień, M. Metodyka Wyboru Optymalnej Metody Wyznaczania Zasięgu Stref Ochronnych Ujęć Zwykłych Wód Podziemnych z Uwzględnieniem Warunków Hydrogeologicznych Obszaru RZGW w Krakowie [Methodology of Selecting the Optimal Method of the Wellhead Protection Area Delineation Taking into Account the Hydrogeological Conditions in Areas Administered by the Regional Water Management Board in Cracow]; Akademia Górniczo-Hutnicza: Kraków, Poland, 2013; ISBN 9788388927294. (In Polish) [Google Scholar]
  50. Jaworska-Szulc, B. Formowanie się zasobów wód podziemnych w młodoglacjalnym, wielopoziomowym systemie wodonośnym na przykładzie Pojezierza Kaszubskiego. [Formation of Groundwater Resources in the Young Glacial Multiaquifer System, of the Kashubian Lake District]; GUT (Gdańsk University of Technology) Publishing House Monograph: Gdańsk, Poland, 2015. [Google Scholar]
  51. Czyżyk, F.; Steinhoff-Wrześniewska, A. Zróżnicowanie ewapotranspiracji niektórych gatunków roślin uprawnych w warunkach różnego nawożenia [Variability of evapotranspiration of some crop species under different fertilizing schemes]. Woda-Środowisko-Obszary Wiejskie 2017, 17, s25–s36. [Google Scholar]
  52. Bogdanowicz, R.; Cysewski, A. Przestrzenna i czasowa zmienność transportu zanieczyszczeń w wybranych ciekach Nadmorskiego Parku Krajobrazowego [Spatial and temporal variability of pollutants transport in selected streams of Nadmorski Park Krajobrazowy]. In Wody na obszarach chronionych; Pociask-Karteczka, J., Partyka, J., Eds.; Instytut Geografii i Gospodarki Przestrzennej UJ, Ojcowski Park Narodowy, Komisja Hydrologiczna PTG: Kraków, Poland, 2008; pp. 91–100. [Google Scholar]
  53. Stachy, J. Atlas hydrologiczny Polski [Hydrological atlas of Poland] Volume 1; IMGW, Wydawnictwa Geologiczne: Warszawa, Poland, 1987; T 1. [Google Scholar]
  54. Orzeł, S.; Forgiel, M.; Ochał, W.; Socha, J. Nadziemna biomasa i roczna produkcja drzewostanów sosnowych Puszczy Niepołomickiej [Above-ground biomass and yearly production of pine stands in Puszcza Niepołomicka]. Sylwan 2006, 150, 16–32. [Google Scholar]
  55. Wejherowo Forest Administration; Wejherowo, Poland. Personal communication.
  56. Jadczyszyn, T.; Kowalczyk, J.; Lipiński, W. Zalecenia nawozowe dla roślin uprawy polowej i trwałych użytków zielonych [Recommendations for fertilizing crops and perennial grasslands]. In Materiały szkoleniowe; Nr 95; IUNG-PIB: Puławy, Poland, 2010; p. 24. [Google Scholar]
  57. Niswonger, R.G.; Panday, S.; Ibaraki, M. MODFLOW-NWT, a Newton formulation for MODFLOW-2005. Us Geol. Surv. Tech. Methods 2011, 6, 44. [Google Scholar]
  58. Aquaveo.com. Available online: https://www.aquaveo.com/software/gms-groundwater-modeling-system-introduction (accessed on 6 April 2020).
  59. Winston, R.B. ModelMuse: A Graphical User Interface for MODFLOW-2005 and PHAST; US Geological Survey: Reston, VA, USA, 2009; p. 52.
  60. Bakker, M.; Post, V.; Langevin, C.D.; Hughes, J.D.; White, J.T.; Starn, J.J.; Fienen, M.N. Scripting MODFLOW model development using Python and FloPy. Groundwater 2016, 54, 733–739. [Google Scholar] [CrossRef]
  61. Papadopulos, S.S.; Larson, S.P. Aquifer storage of heated water: II. Numerical simulation of field results. Groundwater 1978, 16, 242–248. [Google Scholar] [CrossRef]
  62. Sykes, J.F.; Pahawa, S.B.; Ward, D.S.; Lantz, D.S. The Validation of SWENT, a Geosphere Transport Model; W: Scientific Computing, red. R. Stapleman i in; IMAES/North Holland: Amsterdam, The Netherlands, 1983; pp. 351–361. [Google Scholar]
  63. Molenat, J.; Gascuel-Odoux, C.H. Modelling flow and nitrate transport in groundwater for the prediction of water travel times and of consequences of land use evolution on water quality. Hydrol. Process. 2002, 16, 479–492. [Google Scholar] [CrossRef]
  64. Almasri, M.N.; Kaluarachchi, J.J. Modeling nitrate contamination of groundwater in agricultural watersheds. J. Hydrol. 2007, 343, 211–229. [Google Scholar] [CrossRef]
  65. Florida Onsite Sewage Nitrogen Reduction Strategies Study. Literature Review of Nitrogen Fate and Transport Modeling. TaskD.2. Final Report; Florida Department of Health, Division of Environmental Health: Tallahassee, FL, USA, 2010. [Google Scholar]
  66. Psarropoulou, E.T.; Karatzas, G.P. Pollution of nitrates–contaminant transport in heterogeneous porous media: A case study of the coastal aquifer of Corinth, Greece. Glob. Nest J. 2014, 16, 9–23. [Google Scholar]
  67. Lidzbarski, M. Identyfikacja systemu krążenia wód podziemnych w procesie ustalania zasobów odnawialnych na przykładzie zlewni Redy i Zagórskiej Strugi. [Identification of groundwater circulation system during assessment renewable resources for example of the Reda and Zagórska Struga catchment.]. Prz. Geol. 2015, 63, 893–900. [Google Scholar]
  68. Ilnicki, P. Polskie Rolnictwo a Ochrona Środowiska [Polish agriculture and environmental protection]; Wydawnictwo Akademia Rolnicza Poznań: Poznań, Poland, 2004. [Google Scholar]
  69. Potrykus, D.; Gumuła-Kawęcka, A.; Jaworska-Szulc, B.; Pruszkowska-Caceres, M.; Szymkiewicz, A. Assessing groundwater vulnerability to pollution in the Puck region (denudation moraine upland) using vertical seepage method. E3S Web. Conf. 2018, 44, 00147. [Google Scholar] [CrossRef]
  70. Szymkiewicz, A.; Gumuła-Kawęcka, A.; Potrykus, D.; Jaworska-Szulc, B.; Pruszkowska-Caceres, M.; Gorczewska-Langner, W. Estimation of conservative contaminant travel time through vadose zone based on transient and steady flow approaches. Water 2018, 10, 1417. [Google Scholar] [CrossRef] [Green Version]
  71. Szymkiewicz, A.; Savard, J.; Jaworska-Szulc, B. Numerical Analysis of Recharge Rates and Contaminant Travel Time in Layered Unsaturated Soils. Water 2019, 11, 545. [Google Scholar] [CrossRef] [Green Version]
Figure 1. Map of land use in the study area with boundaries of SWAT and MODFLOW/MT3D models.
Figure 1. Map of land use in the study area with boundaries of SWAT and MODFLOW/MT3D models.
Water 12 01042 g001
Figure 2. Simplified conceptual hydrogeological model with components of groundwater budget in steady state.
Figure 2. Simplified conceptual hydrogeological model with components of groundwater budget in steady state.
Water 12 01042 g002
Figure 3. Plots of observed versus computed groundwater heads in steady-state model for aquifers Qp, Q1, and Q2.
Figure 3. Plots of observed versus computed groundwater heads in steady-state model for aquifers Qp, Q1, and Q2.
Water 12 01042 g003
Figure 4. Average monthly values of selected groundwater budget components and N-NO3 loads for scenario S1 (baseline).
Figure 4. Average monthly values of selected groundwater budget components and N-NO3 loads for scenario S1 (baseline).
Water 12 01042 g004
Figure 5. Average monthly values of selected groundwater budget components and N-NO3 loads for scenario S2 (winter wheat).
Figure 5. Average monthly values of selected groundwater budget components and N-NO3 loads for scenario S2 (winter wheat).
Water 12 01042 g005
Figure 6. Average monthly values of selected groundwater budget components and N-NO3 loads for scenario S3 (silage corn).
Figure 6. Average monthly values of selected groundwater budget components and N-NO3 loads for scenario S3 (silage corn).
Water 12 01042 g006
Figure 7. Spatial distribution of monthly recharge values (A) and the associated N-NO3 loads (B) in the baseline scenario, March 2007.
Figure 7. Spatial distribution of monthly recharge values (A) and the associated N-NO3 loads (B) in the baseline scenario, March 2007.
Water 12 01042 g007
Table 1. Comparison of calibrated model outputs with reference values (a: Calibration goals, b: Validation checks, c: Average for the whole model area, d: Variability range depending on 353 hydrological response units (HRU) type).
Table 1. Comparison of calibrated model outputs with reference values (a: Calibration goals, b: Validation checks, c: Average for the whole model area, d: Variability range depending on 353 hydrological response units (HRU) type).
Output TypeUnitModel ValuesReference Values
Groundwater recharge (a)mm/y36 to 146 (d)19 to 186 [49,50]
Evapotranspiration, including REVAP (a)mm/y459 (c)450 to 495 [51]
Total runoff (b)mm/y162 (c)47 to 268 [52]
Surface runoff / total runoff ratio (b)-0.55 (c)0.5 [53]
Forest biomass production (a)t/ha/y6.1 to 8.5 (d)6.5 to 7.5 [54,55]
Yield: Winter wheat (b)t/ha/y5.9 to 7.3 (d)5.5 [47]
Yield: Canola (b)t/ha/y2.6 to 3.4 (d)3.4 [47]
Yield: Silage corn (b)t/ha/y9.2 to 10.3 (d)13.5 [47]
Yield: Hay (b)t/ha/y4.7 to 7.0 (d)4.0 to 10.0 [56]
Table 2. Values of hydraulic conductivity and specific yield assigned for each layer of the model.
Table 2. Values of hydraulic conductivity and specific yield assigned for each layer of the model.
Layer No.TypeHydraulic Conductivity [m/s]
MINMAX
1aquifer (Qp)5.56 × 10−51.39 × 10−4
2aquitard1.25 × 10−92.50 × 10−8
3aquifer (Q1)1.94 × 10−54.72 × 10−4
4aquitard7.78 × 10−91.94 × 10−8
5aquifer (Q2)8.33 × 10−55.56 × 10−4
6aquitard7.78 × 10−97.78 × 10−9
Table 3. Average flow rates and nitrate loads obtained from transient simulations.
Table 3. Average flow rates and nitrate loads obtained from transient simulations.
ScenarioRecharge [mm/y]N-NO3 Leaching from Soil [kg/ha/y]Discharge to Puck Bay [m3/h]N-NO3 Load to Puck Bay [kg/h]
Q1Q2Q1Q2
S1 (baseline)7319.73868750.950.10
S2 (winter wheat)627.43678390.620.03
S3 (silage corn8423.94039061.120.14
S4 (canola)7530.53898811.220.15
S5 (summer cereals)7312.73848720.790.07
S6 (potatoes)9026.04149271.210.17
S7 (peas) 9219.54179331.020.12
S8 (only farmland)7933.93918811.750.14
S9 (only grassland)853.43988970.590.03
S10 (only forest)7114.73798581.190.04
Table 4. Calculated and measured N-NO3 concentration
Table 4. Calculated and measured N-NO3 concentration
LayerN-NO3 in SGD Flux in Numerical Simulations
[mg/dm3]
N-NO3 Measured in Groundwater in Coastal Area
[mg/dm3]
Upper aquifer (Q1)1.48–5.03<0.02–20.78
Lower aquifer (Q2)0.03–0.38<0.02–0.23

Share and Cite

MDPI and ACS Style

Szymkiewicz, A.; Potrykus, D.; Jaworska-Szulc, B.; Gumuła-Kawęcka, A.; Pruszkowska-Caceres, M.; Dzierzbicka-Głowacka, L. Evaluation of the Influence of Farming Practices and Land Use on Groundwater Resources in a Coastal Multi-Aquifer System in Puck Region (Northern Poland). Water 2020, 12, 1042. https://doi.org/10.3390/w12041042

AMA Style

Szymkiewicz A, Potrykus D, Jaworska-Szulc B, Gumuła-Kawęcka A, Pruszkowska-Caceres M, Dzierzbicka-Głowacka L. Evaluation of the Influence of Farming Practices and Land Use on Groundwater Resources in a Coastal Multi-Aquifer System in Puck Region (Northern Poland). Water. 2020; 12(4):1042. https://doi.org/10.3390/w12041042

Chicago/Turabian Style

Szymkiewicz, Adam, Dawid Potrykus, Beata Jaworska-Szulc, Anna Gumuła-Kawęcka, Małgorzata Pruszkowska-Caceres, and Lidia Dzierzbicka-Głowacka. 2020. "Evaluation of the Influence of Farming Practices and Land Use on Groundwater Resources in a Coastal Multi-Aquifer System in Puck Region (Northern Poland)" Water 12, no. 4: 1042. https://doi.org/10.3390/w12041042

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop