Spatially distributed simulation of water balance and sediment transport in an agricultural field

Soil & Tillage Research 143 (2014) 26–37

Contents lists available at ScienceDirect

Soil & Tillage Research journal homepage: www.elsevier.com/locate/still

Spatially distributed simulation of water balance and sediment transport in an agricultural field Lassi Warsta a,*, Antti Taskinen b, Maija Paasonen-Kiveka¨s c, Tuomo Karvonen d, Harri Koivusalo e a

Aalto University School of Engineering, Department of Civil and Environmental Engineering, P.O. Box 15500, FI-00076 Aalto, Finland Finnish Environment Institute, P.O. Box 140, FI-00251 Helsinki, Finland Sven Hallin Research Foundation, Simonkatu 12 A 11, FI-00100 Helsinki, Finland d Waterhope, Munkkiniemen Puistotie 20 A 6, FI-00330 Helsinki, Finland e Aalto University School of Engineering, Department of Civil and Environmental Engineering, P.O. Box 15200, FI-00076 Aalto, Finland b c

A R T I C L E I N F O

A B S T R A C T

Article history: Received 12 November 2013 Received in revised form 11 April 2014 Accepted 4 May 2014

Runoff and sediment transport are distinctively three-dimensional (3D) processes and occur through overland, tillage layer and subsurface pathways. The objective was to quantify water balance and sediment concentrations in runoff waters and to assess sediment loads via surface runoff and drainflow in a clayey, subsurface drained field section using the FLUSH model. The model can simulate field scale two-dimensional overland flow and 3D unsaturated and saturated subsurface flow, including preferential flow in macropores. The erosive processes, comprised of hydraulic and raindrop splash erosion, occur in the overland domain while suspended sediment is conveyed from the field surface to subsurface drains via preferential transport in macropores in the subsurface domain. The study site, located in southern Finland, is a clayey, subsurface drained field section with an area of 12 ha and an average slope of 2.8%. The growing seasons and the following autumns of two years with distinctly different rainstorm characteristics were modelled. The simulated sediment loads via tillage layer runoff and drainflow were 85 and 117 kg ha1, respectively in 1988, and 63 and 189 kg ha1 in 1984. Despite the high precipitation in October 1984 (143 mm), erosion in the field area was low due to minimal surface runoff and consequently minimal hydraulic erosion. The suspended sediment at the site was generated by raindrop impacts and the eroded material was transported to the open ditch through tillage layer flow and subsurface drainflow. The model was able to reproduce the measured sediment concentrations in the main open ditch into which both tillage layer runoff and drainflow are discharged from the field. Spatially distributed erosion simulations facilitate the detection of net erosion and deposition locations in the field and could be used to design intensive measurement campaigns and guide erosion control practices in the future. ß 2014 Elsevier B.V. All rights reserved.

Keywords: FLUSH model Clay soil Erosion Sediment concentration Surface runoff Drainflow

1. Introduction Water balance and sediment loads in agricultural fields are tightly coupled in the Nordic countries, where hydrometeorological conditions are characterized by a short growing season with wet soils, and high runoff and sediment load generation potential outside the growing season. Sediment export from cultivated fields into open waterways poses problems due to the environmental effects of sediment itself, such as siltation and increased turbidity,

* Corresponding author. Tel.: +358 505416581. E-mail address: [email protected] (L. Warsta). http://dx.doi.org/10.1016/j.still.2014.05.008 0167-1987/ß 2014 Elsevier B.V. All rights reserved.

and due to the particle bound nutrients, pesticides and heavy metals that are potentially released into receiving rivers and lakes (e.g. Boardman and Poesen, 2006). Suspended sediment is mainly generated by erosive impacts caused by raindrops and by hydraulic erosion on the field surface of cultivated soils in high latitudes (e.g. Boardman and Poesen, 2006). Erosion remains an elusive, spatial and stochastic process, and the processes that lead to large sediment loads in one field might produce only a little erosion in another site (Ba¨rlund et al., 2009). Erosion is affected by climate and weather, cultivation and drainage practices, the topography of the field and the surrounding areas, and the soil type and moisture conditions. The variety of these controlling factors makes it difficult to generalize sediment load estimates computed from

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

measured concentrations and discharge volumes in individual experiments. In empirical studies, sediment loads are estimated from field or small catchment-scale runoff and sediment concentration data using various schemes (e.g. Bechmann and Sta˚lnacke, 2005; Puustinen et al., 2005; Paasonen-Kiveka¨s et al., 2008). The effects of tilling on erosion and nutrient loads in cultivated fields in high latitudes are recognized and investigated in several studies (e.g. Chapman et al., 2005; Puustinen et al., 2005; Turtola et al., 2007; Uusitalo et al., 2007; Muukkonen et al., 2009; Ule´n et al., 2010; Skøien et al., 2012), while research on the effects of seasonal and annual hydrological variations on erosion are more scarce (e.g. Puustinen et al., 2007). In clayey, subsurface drained fields, preferential flow in soil macropores has a large impact on water balance, flow pathways and sediment transport because a notable part of the annual sediment load can be lost via drainflow (e.g. Øygarden et al., 1997; Turtola et al., 2007; Paasonen-Kiveka¨s et al., 2008). Different lines of research have produced evidence for the direct and short routes from the surface to the subsurface drains including Caesium-137 (Laubel et al., 1999; Uusitalo et al., 2001), colloid and particle tracer studies (McKay et al., 1993; Jacobsen et al., 1997; Joel et al., 2012) and matching sediment concentrations in surface runoff and drainflow (Turtola and Paajanen, 1995; Uusitalo et al., 2001; Paasonen-Kiveka¨s et al., 2008). In some cases the soil profile is found to filter suspended sediment from the water moving in the macropores (Turtola and Paajanen, 1995; Jacobsen et al., 1997; Turtola et al., 2007). Mathematical models can be used to explain different factors affecting erosion because models facilitate separation of the water balance and sediment transport into individual components – such as surface and tillage layer runoff, drainflow and groundwater outflow – and thus disclose the underlying governing processes. Sediment loads from cultivated fields have been quantified by applying simple models (e.g. Lundekvam, 2007; Bechmann, 2012) and field scale process-based models (Knisel and Turtola, 2000; Tattari et al., 2001; Taskinen and Bruen, 2007a,b; Ba¨rlund et al., 2009), as well as large-scale distributed models (e.g. Ba¨rlund et al., 2007; Puustinen et al., 2010). Only a few previously published modelling applications quantify sediment load estimates via both surface runoff and drainflow (Knisel and Turtola, 2000; Larsson et al., 2007; Lundekvam, 2007; Warsta et al., 2013b). When the models are calibrated and validated directly against empirically estimated sediment load accumulations, the assessment of the model performance may be biased by the method of comparison between model results and measurements. Firstly, the model can underestimate runoff but overestimate concentrations, which results in a load that is comparable to the load estimated from data for the wrong reasons. Secondly, the number of sediment concentration measurements is usually much lower than the number of runoff measurements, which leads to uncertainty in the derivation of the empirically estimated loads (e.g. Bechmann and Sta˚lnacke, 2005; Puustinen et al., 2005; Paasonen-Kiveka¨s et al., 2008). Thirdly, the comparison between measured and simulated long-term load accumulations may hide periods of load under- and overestimation, which can compensate for each other over longer periods. Finally, erosion and sediment transport are spatially variable processes and both model applications and measurement campaigns should strive for simulating and detecting threedimensional (3D) variables instead of only aggregated measures. Water pollution control measures are implemented within agriculture fields, where the sediment source areas are located. The assessment of efficiency of these measures has been proven to be very difficult (e.g. Granlund et al., 2005). We propose that model simulations based on detailed 3D description of the field are a way forward in the estimation of management impacts on water balance, sediment loads and nutrient export.

27

Ule´n and Bechmann (2012) identified several gaps in soil erosion knowledge in the Nordic countries including (1) the quantification of sediment transport in macropores to subsurface drains and (2) the utilization of spatially distributed models for assessing erosion. In this study we aim to demonstrate through a case study how the newly introduced FLUSH model (Warsta, 2011; Warsta et al., 2013a,b) can be applied to fill in these gaps. The leading objective is to assess the water balance and the sediment loads via different pathways during the growing seasons and following autumns in a subsurface drained field section in southern Finland. While suspended sediment concentrations are routinely measured in empirical measurement campaigns, modelling studies usually present only cumulative sediment loads or total loads during the investigated period (e.g. Knisel and Turtola, 2000; Tattari et al., 2001; Larsson et al., 2007; Lundekvam, 2007). We provide a more detailed picture of water and sediment processes by presenting both hourly simulated suspended sediment concentrations and measured concentrations not included in these earlier studies. Taskinen and Bruen (2007a,b) previously compared simulated sediment concentrations against measurements in surface runoff with an event based model. We will produce an hourly simulation of soil erosion during the growing seasons and the following autumns, and account for sediment discharge via both surface runoff and subsurface drains, focusing on the autumn periods with their high erosion risk after tilling. We also present 3D field scale, spatial distribution of the erosion affected by the topography of the site. Exploring the spatial distribution of field scale erosion processes is an extension of the commonly produced aggregated loads of one-dimensional (1D) or conceptual models lumped into the field scale (e.g. Knisel and Turtola, 2000; Tattari et al., 2001; Larsson et al., 2007; Lundekvam, 2007; Ba¨rlund et al., 2009). 2. Materials and methods The water balance, sediment concentrations in runoff waters and sediment loads of 1984 and 1988 of the Hovi monitoring field, located in southern Finland, are investigated with the FLUSH modelling system (Warsta, 2011; Warsta et al., 2013a,b). These years were chosen because they included recorded rainfall–runoff events with erosive potential in autumns. Temporal data including precipitation, runoff and sediment concentration measurements from the field site are shown in Section 3 in conjunction with the simulation results. 2.1. Site description The Hovi field site is one of the small catchments where water flow and quality, including erosion, are monitored in Finland by the Finnish Environment Institute. The gathered data has been used to estimate sediment and nutrient loads for different land-use types in the whole of Finland (Vuorenmaa et al., 2002). The Hovi field represents clayey, subsurface drained fields, which are abundant in the southern and western parts of Finland. The detailed simulations will produce further evidence of the reliability of the monitoring system and the data. Runoff and water quality data from the Hovi field have been used in the estimation of the water balance, the environmental impacts of drainage, and erosion and nutrient loads in several previous studies (e.g. Seuna and Kauppi, 1981; Bengtsson et al., 1992; Puustinen et al., 2007, 2010). Moreover, Seuna and Kauppi (1981) studied the effects of subsurface drainage installation on the quantity and quality of runoff and Bengtsson et al. (1992) investigated the pathways of melt water movement using oxygen-18 as a tracer concluding that water reached subsurface drains via cracks and vertical sand– gravel drains.

[(Fig._1)TD$IG]

28

[(Fig._2)TD$IG]

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

a)

b)

Fig. 2. The Hovi monitoring field land use areas (LU1 and LU2) in (a) 1984 and (b) 1988.

Fig. 1. The location of Vihti and a map of the Hovi monitoring field. Thin lines inside the field borders are the elevation contours (elevations are presented in metres), the thicker lines are the subsurface drains and the white circles are the measurement points for subsurface drainflow (S) and flow in the main ditch (P).

The Hovi agricultural catchment is located in Vihti (608250 2000 N 248220 700 E, Fig. 1). The climate of the Hovi area is temperate with a mean annual precipitation of 700 mm (uncorrected value) and a mean annual air temperature of +5 8C. The 12 ha field section is L-shaped and embedded within a larger area of arable land. The topography is undulating with a general flow direction spiralling clockwise from south-east to north-east. The average slope is 2.8% (Seuna and Kauppi, 1981). A digital elevation map (DEM) of the field was provided by MTT Agrifood Research Finland, with a pixel resolution of 2.0 m  2.0 m. The DEM was measured with a real time kinematic global positioning system mounted on a tractor. The field topography features a large ovalshaped depression (100  50 m2) in the middle of the field (Fig. 1). The field section is not a natural catchment but defined according to the subsurface drainage system. Embankments were ploughed along the north-east side to direct surface runoff to the open ditch (Fig. 1). A road ditch borders the north-west side of the field. At the south-west side, the field is delimited by a natural water divide and a farmyard. During the simulation years the field was divided into two land use areas (LU1 and LU2) (Fig. 2). The differences in land use in 1984 and 1988 are due to changes in cultivation practices and fallow rotation. The land use information corresponding to these areas is presented in Table 1, i.e. the field section was cultivated normally during the study years. The subsurface drainage system was installed in 1971 with a drain spacing of 20 m and an average drain depth of 1 m (Bengtsson et al., 1992). The drainage area is the whole field section shown in Fig. 1. The lateral plastic drainpipes, with a diameter of 0.055 m, were connected to a single collector pipe, which in turn discharged into the main ditch. Vertical gravel drains (0.1 m2) with a spacing of 10 m were installed simultaneously with the lateral drainpipes (Bengtsson et al., 1992). Subsurface drainflow and total runoff in the main ditch were measured near the north-east corner of the field by using v-notch weirs (Fig. 1).

MTT Agrifood Research Finland recently conducted a campaign in the field to measure the particle size distributions, organic carbon contents and hydraulic properties of the soil from several locations above and between the drain lines in the north-eastern part of the field. The clay content varied between 54 and 77% and silt content between 17 and 27%. The samples were taken from three different depths (0.0–0.2, 0.2–0.35 and 0.35–0.6 m) (Table 2). Clay, silt and organic carbon contents are presented in Table 2. The measured hydraulic properties included saturated hydraulic conductivities (KS), macroporosities (w), and water retention curves (WRCs). Since measurements taken at the drain line locations characterize only a minor fraction of the field, the data gathered between the drain lines were used to describe the general soil properties of the field in the simulations. The measurements were conducted with samples collected from the same depths as the soil texture data (Table 2). The diameter dividing the pore volume between macropores and micropores was assumed to be 300 mm. Because data on the hydraulic properties of the soil were only available from limited locations in the north-east part of the field, arithmetic mean values of the hydraulic parameters were applied to the whole field area. The changes in the hydraulic properties of the soil are more pronounced in the vertical direction, i.e. running vertically down the profile, rather than in horizontal directions within the field area (e.g. Table 4 in Alakukku et al., 2010).

Table 1 Land use information about the Hovi monitoring field. Date

LU1

LU2

3 Oct. 1983–8 Oct. 1983 1 May 1984 17 Aug. 1984 13 Sep. 1984 23 Oct. 1984 29 Sep. 1987 10 May 1988 11 May 1988

Ploughing Seeding of spring wheat Harvesting – Ploughing Ploughing Fallow –

17 Aug. 1988 9 Sep. 1988

– Ploughing, harrowing, seeding of autumn wheat Autumn wheat springing up

Ploughing Fallow Rye seeding Rye springing up – Ploughing Harrowing Seeding of spring wheat Harvesting –

29 Sep. 1988

Ploughing

Table 2 Mean clay, silt and organic carbon contents (g g1) at different soil depths (m). Range is in parentheses (minimum/maximum values). Depth

Clay (0.002 mm)

Silt (0.002–0.02 mm)

Organic carbon

0.0–0.2 0.2–0.35 0.35–0.6

0.54 (0.32/0.75) 0.71 (0.40/0.86) 0.77 (0.64/0.89)

0.27 (0.15/0.38) 0.20 (0.08/0.38) 0.17 (0.06/0.29)

0.02 (0.02/0.03) 0.01 (0.003/0.02) 0.003 (0.002/0.004)

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

The subsurface drainflow and total runoff (drainflow plus tillage layer runoff) data covered the whole year in 1984 and 1988, and they were measured with v-notch weirs. Tillage layer runoff, which is composed of surface runoff and seepage to the open ditches, can be calculated by subtracting the measured drainflow from the total runoff. Suspended sediment concentration measurements were available from the main ditch at irregular intervals throughout April–October in 1988. Only a few measurements from the ditch and subsurface drainflow were available from 1984. Precipitation was recorded with pluviographs with Wild-type collectors (Posch and Rekolainen, 1993; Taskinen, 2002). The hourly precipitation data covered the period from May to October in 1984 and 1988. The measured precipitation was corrected with a constant factor of 1.05, following the recommendation by Førland et al. (1996) for a gage located between a forest and plain at least 10 km from coast. The cumulative daily precipitation values were verified by comparing them to the data from the nearby Maasoja weather station (608250 700 N 248230 5400 E), operated by the Finnish Meteorological Institute and located 2.5 km from the field. Potential evapotranspiration (PET) was calculated from daily Class A evaporation pan results from the Maasoja weather station. The pan results were multiplied with the monthly correction coefficients of Vakkilainen (1982) to derive the daily PET values. Vakkilainen (1982) also estimated how much latent energy was consumed by plant transpiration in each hour in a day during the April–October period, and these data were used to divide the daily PET values into hourly values. 2.2. Description of the model FLUSH is a 3D distributed hydrological model that simulates water flow, soil erosion and sediment transport in clayey, subsurface drained agricultural fields (Warsta, 2011; Warsta et al., 2013a,b). The simulated area is divided into two-dimensional (2D) overland and 3D subsurface domains (including the tillage layer) with curvilinear grids composed of rectangular and hexahedric cells. The governing partial differential equations are solved implicitly with the finite volume method. The field section can have several independent drainage systems composed of subsurface drains and open ditches (Turunen et al., 2013). Overland flow is described with the diffuse wave approximation of the Saint Venant equations (Appendix A, Eqs. (1) and (2)), while subsurface flow is presented with a dual-permeability approach, applying the Richards equation in both macropore and matrix pore systems (Eqs. (3) and (4)). The Manning’s equation is used to describe the flow friction terms in the overland domain (Eq. (2)). A parametric Mualem-van Genuchten (MVG) WRC (van Genuchten, 1980) is applied to describe soil water retention properties and unsaturated hydraulic conductivity in both pore systems. Water exchange between the pore systems is described with the first-order approach, driven by the pressure difference between the pore systems (Warsta et al., 2013a). The model supports soil shrinkage and swelling derived from the SWAP model (Van Dam et al., 2008). Tillage layer runoff includes surface runoff and rapid horizontal flow in the topsoil layers above the tillage pan. The PET value derived from the Class-A pan measurement is divided into the soil profile according to the plant root mass distribution. The root depths are given as an input time series, and the root mass is computed as a linearly decreasing function from the field surface to the prescribed root depth. Evapotranspiration (ET) is decreased with the moisture limiter presented by Feddes et al. (1978). Overland soil erosion is based on the sediment continuity equation combined with hydraulic and raindrop splash erosion, sediment settling due to gravity and transport capacity processes (Appendix A, Eq. (5)). Suspended sediment can infiltrate into the

29

soil macropores in the field surface and be transported with preferential flow in macropores to the subsurface drains (Eq. (6)). Sediment can also be transported with lateral groundwater flow when the groundwater table rises above the drain level. It is assumed that the macropore pathways below the subsurface drains are too small for particle transport. More information about the applied solution methods, and water flow and soil erosion processes are available in Warsta et al. (2013a,b) and Warsta (2011).

3. Simulation results The model parameters were calibrated with data from 1988 and validated with data from 1984 because the former year contained more data than the latter one. The model parametrisation is presented first, and thereafter cumulative runoff results are shown with illustrations of spatially distributed cumulative erosion and deposition in the field area. Finally, hourly results of the simulated runoff and sediment concentrations are presented. 3.1. Model parameterization The simulation periods of 1988 and 1984 started in the beginning of May, when the snow cover had melted in southern Finland, and continued till the end of October. The computational time steps in the simulations varied between 0.5 min and 1.0 h, the shorter values being used during wet periods. The Nash–Sutcliffe model efficiency (NSME) coefficient (Nash and Sutcliffe, 1970), bias and mean absolute error were applied to quantify the quality of the calibration and validation results against measurements. The measures were computed from hourly runoff results. The calibration–validation process was conducted with the trial and error approach. The Hovi field area (Fig. 1) was divided into cells with a 3D computational grid. The horizontal resolution of the grid was 104  144 cells in x- and y-directions and the cell dimensions were 4 m  4 m. The depth of the computational grid was 2.4 m and it was divided into 16 cells. Most of the vertical resolution was concentrated into the volume between the field surface and subsurface drains. The thicknesses of the cell layers were 0.02, 0.05, 0.08, 0.1  9, 0.25  2, 0.35 and 0.5 m. Water levels in the open ditches were not simulated. The depths of the ditches were set to 1.0 m and a constant water depth of 0.3 m was maintained as a boundary condition. No flow from the ditch towards the soil domain was allowed. The subsurface drainage network is presented in Fig. 1. The depth of drains was assumed to be 1.0 m and the radius of the pipes was set to 0.025 m. The flow path length to the subsurface drains was set here to 1.0 m (Turunen et al., 2013). According to the available land use information (Table 1), the field was partially fallow in both 1984 and 1988. On the 23rd of October 1984, LU1 (Fig. 2) was ploughed while LU2 was seeded with rye on the 17th of August. A limited data is available about the operations in LU2 in 1984 but it was assumed that light tillage was conducted on the fallow field before the seeding. In 1988, LU1 and LU2 were ploughed on the 9th and the 29th September, respectively. The values of the raindrop erodibility (kR in Eq. (5) in Appendix A) was increased in LU1 at the time of ploughing and in LU2 at the time of seeding, from 4 to 20 J1 in 1984, to take into account the increase in erodibility caused by the tillage. In 1988, the values of kR in LU1 and LU2 were increased in the same way (from 4 to 20 J1) after ploughing (Table 1). Maximum root depth was set to a 0.6 m depth from the field surface, and the crops reached this depth in August. ET in the fallow and harvested areas were restricted to a 0.05 m depth from the surface.

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

30

Table 3 Statistics of the measured KS (mm h1) and w (m3 m3) of the soil samples at different depths (m). Statistic

w

KS

Depth

0.0–0.20

0.20–0.35

0.35–0.55

0.0–0.20

0.20–0.35

0.35–0.55

Mean Median Min. Max.

163 85 0.1 752

85 23 0.01 545

0.6 0.1 0.01 7

0.019 0.019 0 0.043

0.003 0.002 0 0.013

0.001 0.001 0 0.003

The fitted MVG WRC parameter values are presented in Table 4. The parameter values for the macropore system were derived from Warsta et al. (2013a). Tillage layer WRCs (Table 4) were used in the tillage layer matrix pore system (0–0.25 m), while the bottom soil WRCs were assigned to the soil matrix pore system in the rest of the profile (0.25–2.4 m). The macropore WRC (Table 4) was employed in all the soil horizons. The static macroporosity values (w in Eqs. (3) and (4)) in the upper part of the profile were derived directly from the measurements (Table 3, mean values). The values in the layer 0.35–0.55 m were also applied in the profile below this layer, down to the depth of the subsurface drains (1.0 m). A small w value (1.0  104 m3 m3) was applied to the rest of the profile below the subsurface drains. There were vertical gravel drains within the field but they were too small to be explicitly presented in the computational grid. Their effect was taken into account by increasing w of the soil by 0.001 m3 m3 in the profile above the subsurface drains. The w values were decreased to 40% of the original values in the tillage layer (0–0.25 m) during both years after the tillage in both land use areas (Table 1). The values of w are not affected by soil shrinkage and swelling processes in the model. The Manning’s n (Eq. (2)) for the field surface was calibrated to a high value of 1.0. The value was further increased after the tillage to 3.0. The rest of the model parameters were given the same value as in Warsta et al. (2013a,b) including the particle density (rS = 2650.0 kg m3, see Eq. (3) in Warsta et al., 2013b), the mean particle diameter (DS = 1.5  105 m, see Eq. (3) in Warsta et al., 2013b) the overland flow threshold depth (hW,THR = 0.001 m, Eq. 2), the parameters of the soil shrinkage characteristic curves, the macropore saturated hydraulic conductivity multiplier (KFS,MUL = 80 m h1, see Eq. (6) in Warsta et al., 2013a), the parameters controlling water exchange between the pore systems and the sediment lateral and transverse dispersivities (aL = 0.1 m and aT = 0.01 m, respectively, see Eq. (6) in Warsta et al., 2013b) in the macropore system. In the simulations, the groundwater levels were adjusted to 0.1 m in 1988 and 1.2 m in 1984 by calibration against the measured runoff. No groundwater table level measurements were available from the field to assess the modelled levels. In particular, the 1988 initial groundwater table level had to be higher than in 1984 because, according to measurements, the first precipitation events in the end of spring in 1988 generated runoff, indicating that the profile had a higher moisture content compared to that of the validation year. Thus it was not possible to use the same initial groundwater table depths during 1988 and 1984. Overland water depths and sediment concentrations, as well as subsurface sediment concentrations, were set to zero at the beginning of the simulations. 3.2. Cumulative runoff results and the spatial distribution of erosion in the field area Precipitation was the sole source of water in the simulations, and water was removed from the field area via tillage layer runoff, ET, drainflow and groundwater outflow. The corrected cumulative

Table 4 Parameter values of MVG WRCs for tillage layer and bottom soil in the Hovi monitoring field. Parameter

Unit

Tillage layer

Bottom soil

aMVG

(m1) (–) (m3 m3) (m3 m3)

7.2011 1.0884 0.5294 0.1

1.3807 1.1197 0.5352 0.1

nMVG

uS uR

precipitation amount from May to October was 384 mm in the calibration year and 476 mm in the validation year. In terms of the simulated mass balance, tillage layer runoff was a minor component (6% of the precipitation during both years). The modelled surface runoff fraction of the tillage layer runoff was very small at the site during both years (the cumulative sum was 11 and 17 mm in 1988 and 1984, respectively), even though October 1984 was wet (the corrected precipitation sum was 143 mm). The fraction of modelled tillage layer runoff (24 and 30 mm in 1988 and 1984, respectively) corresponds with the difference between the measured total runoff and measured drainflow (25 and 41 mm in 1988 and 1984, respectively). Modelled ET was the primary water loss component and it comprised 58% and 49% of precipitation in 1988 and 1984, respectively. Simulated groundwater outflow had an important role in the water balance (15% and 16% of precipitation in 1988 and 1984, respectively). Cumulative measured and simulated tillage layer runoff and drainflow in the calibration (1988) and validation (1984) years are presented in Fig. 3. Cumulative drainflow was over predicted in September and October 1988, while in 1984 simulated drainflow was lower in July compared to the measured values (Fig. 3). Cumulative simulated tillage layer runoff was smaller (27%) than the measured tillage layer runoff in October 1984 (Fig. 3b). The NSME values for tillage layer runoff and drainflow were 0.494 and 0.482, respectively in 1988 and 0.715 and 0.705 in 1984. The coefficients were better in the validation year because runoff events were concentrated into a shorter interval at the end of the simulation period (Fig. 3b). The bias values for tillage layer runoff and drainflow were 0.00008 and +0.003 mm h1, respectively in 1988 and 0.003 and +0.002 mm h1 in 1984. The average tillage layer runoff and drainflow rates were 0.006 and 0.019 mm h1, respectively in 1988 and 0.006 and 0.020 mm h1 in 1984. The total entrained sediment amounts in the simulations were 252 and 298 kg ha1 in 1988 and 1984, respectively. The sediment load fractions of the total entrained sediment amount via tillage layer flow, drainflow and groundwater outflow were 34%, 47% and 12%, respectively in 1988 and 21%, 64% and 15% in 1984. It is noteworthy here that groundwater flow was only able to transport sediment when the water level rose above the drain level and water level in the ditches outside of the field. The sediment was transported to the drains with preferential flow in macropores in the simulations. Some suspended sediment was left on the field surface at the end of the simulations (8% and 0.2% of the entrained amounts in 1988 and 1984, respectively). In 1988, 73% of the sediment flux via tillage layer runoff was transported by surface runoff and in 1984, 72% was transported by surface runoff. The suspended sediment was only detached by raindrop impacts in the model. It is emphasized here that although raindrop impacts were responsible for the particle detachment, sediment was still transported on the field surface by overland flow. Spatial distribution of overland net erosion in the field area at the end of the simulation periods is presented in Fig. 4. Negative erosion values denote sediment deposition in the cell (Fig. 4, depicted in blue). The maximum and minimum simulated net erosion results were +0.168 and 3.53 kg m2, respectively in 1988 and +0.243 and 4.24 kg m2 in 1984. The two land use areas (Fig. 2) are visible in Fig. 4. The value of kR was increased after

[(Fig._3)TD$IG]

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

(b) Tillage layer, measured Tillage layer, simulated Drain, measured Drain, simulated

Cumulative runoff (mm)

100 80 60 40 20

Tillage layer, measured Tillage layer, simulated Drain, measured Drain, simulated

100 Cumulative runoff (mm)

(a)

31

80 60 40 20 0

0 5

6

7

8

9

10

11

Time (month)

5

6

7

8

9

10

11

Time (month)

Fig. 3. Measured and simulated cumulative tillage layer runoff and drainflow in (a) 1988 (calibration) and (b) 1984 (validation). Tillage layer runoff is composed of surface runoff and seepage to open ditches through soil layers affected by tillage.

tillage (Table 1) in autumn, which explains the higher erosion rates. Instead of constant erodibility, it is possible that the erodibility of the soil slowly decreases after the tillage operation. The large blue areas visible in Fig. 4a and partly visible in Fig. 4b coincide with the depression in the field surface (Fig. 1), which acts as a sediment deposition area. The subsurface drainage network is faintly visible as lines with less net erosion compared to their surroundings (Fig. 4). This is partly caused by the drier conditions over the drains where sediment is able to settle down on the field surface. 3.3. Hourly calibration results with the 1988 data Hourly measured and simulated tillage layer runoff values in the calibration year 1988 are presented in Fig. 5. Only the end of September and the beginning of October are presented because tillage layer runoff was minimal during the summer and early autumn. The simulated results compared well to the measurements

[(Fig._4)TD$IG]

in terms of NSME (0.501), although a high n value (3.0) had to be applied to reproduce the measured pattern. The measured and simulated hourly drainflow results are presented in Fig. 6. The NSME value (0.499) was reduced by occasional mismatches between the measured and simulated peaks in early May, August and late October 1988 (Fig. 6). In the autumn, soil swelling and tillage can decrease the hydraulic conductivity of the profile by blocking macropore pathways from the field surface to subsurface drains. Although the model simulates soil swelling, the effects of tillage have to be applied manually. In the simulations we decreased w values by 40% after tillage during both years. However, this can lead to problems, such as flattening and widening of drainflow events, visible in the last runoff event in Fig. 6. Because drainflow sediment concentration measurements were not available in 1988, the simulated drainflow concentrations were plotted against concentrations measured in the open ditch (Fig. 1) for comparison (Fig. 7). While Fig. 7 only presents an overview of the events in 1988, it can be seen that the timing and magnitude of

Fig. 4. Simulated net erosion in (a) 1988 and (b) 1984. The colour range from red to blue corresponds to net erosion results of +0.15 to 0.15 kg m2, respectively. Green indicates zero net erosion.

[(Fig._5)TD$IG]

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

32

1.5

Measured Simulated

1

Precipitation

10 20 30

0.5

40

0

50 24 25 26 27 28 29 30 September

[(Fig._6)TD$IG]

Precipitation (mm h-1)

0 Tillage (LU2)

Till. layer runoff (mm h-1)

2

1

2

3

4

5

6 7 8 October

9

10 11 12

Fig. 5. Hourly precipitation and hourly measured and simulated tillage layer runoff in September–October 1988 (calibration).

0.6

Simulated Precipitation

0.4

10

Tillage (LU2)

Measured

Tillage (LU1)

Drainflow (mm h-1)

0.8

20 30 40

0.2

50

Precipitation (mm h-1)

0

1

0 5

6

7

8

9

10

11

Time (month) Fig. 6. Hourly precipitation and hourly measured and simulated drainflow in May–October 1988 (calibration).

the measurements. The mean absolute error between the measured and simulated concentration values was 0.27 g l1 and the bias was 0.111 g l1 (average measured concentration, 0.42 g l1). The simulated peak was higher than the recorded peak at the end of September.

the simulated erosion events are similar to the measurements. The beginning of the simulation period is not shown because only minimal runoff was simulated. The mean absolute error and bias between the measured and simulated concentration values for the period were 0.17 and 0.093 g l1 (average measured concentration, 0.25 g l1), respectively. Some simulated peaks, such as the one in the end of July, were not measured. Measured sediment concentrations in the ditch are presented in Fig. 8, with the simulated concentrations that represent the combined effect of modelled drainflow and tillage layer runoff [(Fig._7)TD$IG]entering the ditch. The timing of the simulated peaks was similar to

3.4. Hourly validation results with the 1984 data The hourly measured and simulated tillage layer runoff results in the validation year 1984 are presented in Fig. 9. Only October and the beginning of November are presented because, as in 1988,

Ditch, measured Drain, simulated

1

Precipitation

0.5

10 20 30 40

0

Precipitation (mm h-1)

1.5

Tillage (LU2)

0 Tillage (LU1)

Sediment conc. (g l-1l)

2

50 16 20 24 28 1 July

5

9 13 17 21 25 29 2 August

6 10 14 18 22 26 30 4 8 12 September October

Fig. 7. Hourly precipitation and simulated sediment concentrations in drainflow in July–October 1988 (calibration). Occasional measured concentrations in the ditch are plotted for comparison because drainflow concentration measurements were not available.

[(Fig._8)TD$IG]

1

10

Measured Simulated

20

Precipitation

30

0.5

40

0

Precipitation (mm h-1)

1.5

0 Tillage (LU2)

Ditch conc. (g l-1l)

2

50 25 26 27 28 29 30 September

1

2

3

4

5

6

7 8 October

9

10 11 12 13

Fig. 8. Hourly precipitation and measured and simulated ditch sediment concentrations in September–October 1988.

[(Fig._9)TD$IG]

33

1.5

Measured Simulated

1

10

Precipitation

20 30

0.5

40

0

Precipitation (mm h-1)

0

2 Tillage (LU1)

Till. layer runoff (mm h-1)

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

50 30

[(Fig._10)TD$IG]

2

4

6

8

10

12

14 16 18 October

20

22

24

26

28

30

1

3

Fig. 9. Hourly precipitation and measured and simulated tillage layer runoff in 1984 (validation).

Measured

0.6

Simulated Precipitation

0.4

10

Tillage (LU1)

0.8

20

30 40

0.2

50

Precipitation (mm h-1)

0 Tillage (LU2)

Drainflow (mm h-1)

1

0 5

6

7

8

9

10

11

Time (month) Fig. 10. Hourly precipitation and measured and simulated drainflow in 1984 (validation).

year 1984 were better than the results from the calibration year 1988. Only four concentration measurements from drainflow were available in 1984 (Fig. 11). The simulated sediment concentrations in drainflow started to rise in the beginning of October. While measured concentrations were lower, they exhibited a similar rising trend. The simulated concentration results in the ditch flow (Fig. 12) were again calculated from the surface layer runoff and drainflow results. The simulated ditch concentration peaks were an order of magnitude higher than the sediment concentrations in the drainflow (Fig. 11 – note the different y-axis scales in Figs. 11 and 12). No sediment concentration peaks were simulated in

runoff was minimal during the summer and early autumn. Simulated tillage layer runoff peaks were lower than the measurements but the timing was similar. The lower simulated peaks could have been caused by too high infiltration rates and high values of n (1.0–3.0). The large drainflow event at the end of June was not reproduced in the simulated drainflow results (Fig. 10), and some of the peaks at the end of October were smaller than the measured peaks. At the end of June the soil was presumably too dry to produce drainflow due to modelled ET. The drop in drainflow intensities at the end of October is caused by soil swelling and the effect of tillage on modelled drainflow. In terms of the NSME values, the tillage layer [(Fig._1)TD$IG]runoff (0.715) and drainflow (0.705) results from the validation

Simulated 0.1

Precipitation

20 30

0.05

40

0

50 15

19 23 27 September

1

5

9

13

17 21 October

25

29

2

Fig. 11. Hourly precipitation and measured and simulated drainflow sediment concentrations in 1984 (validation).

2 1.5

Measured Simulated

1

10

Precipitation

20 30

0.5

40

0

Precipitation (mm h-1)

0 Tillage (LU1)

Ditch conc. (g l-1l)

[(Fig._12)TD$IG]

Precipitation (mm h-1)

Measured

0.15

10

Tillage (LU1)

Drain conc. (g l-1l)

0 0.2

50 15

19

23 27 September

1

5

9

13

17 21 October

25

29

2

Fig. 12. Hourly precipitation and measured and simulated sediment concentrations in the ditch in 1984 (validation).

34

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

drainflow. The difference between the simulated drainflow and ditch concentrations was solely caused by the flow dynamics as sediment sieving in the macropores was not applied in the model. The amount of precipitation in October 1984 (143 mm) was higher than in 1988 (47 mm), although simulations in 1988 already ended on the 23rd October. Apparently, a large share of the suspended sediment detached by raindrops in October 1984 was transported to the open ditch by tillage layer runoff. It also has to be stated that the concentration data available from the ditch and drainflow are not sufficient for model validation and that the concentration measurements do not coincide with the precipitation events. 4. Discussion According to FLUSH simulations for Hovi, the sediment loads via tillage layer runoff and drainflow during the simulated periods (May–October in both years) were 85 and 117 kg ha1, respectively in 1988 and 63 and 189 kg ha1 in 1984. Ba¨rlund et al. (2009) applied the ICECREAM model to differentiate the effects of different soil types on water balance and erosion, and phosphorus transport in cultivated soils. According to their results, annual sediment loads from silty–clay soils were, on average, 500 kg ha1 a1, with minimum and maximum loads of 200 and 600 kg ha1 a1, respectively. Ba¨rlund et al. (2009) noted that their applied model did not include preferential flow or transport processes in macropores and thus only included sediment loads via surface runoff. Vakkilainen et al. (2010) reported recent experimental results of agricultural pollution in southern Finland and found annual erosion of 948–1774 and 320–2663 kg ha1 a1 from clayey fields in Jokioinen and Siuntio, respectively, under normal cultivation practices. Sediment export in these sites was clearly higher in subsurface drain flow than near surface tillage layer runoff. Bechmann (2012) investigated the field and catchment scale effects of soil tillage in a field section and surrounding catchment of silt and silt–clay soils in south-eastern Norway using a measurement campaign and a simple erosion model. The average annual sediment load in the Norwegian site was 940 kg ha1 a1 when the field was ploughed in autumn and 174 kg ha1 a1 when the field was not tilled. Warsta et al. (2013b) studied sediment loads during growing seasons and fallowing autumns and reported large simulated loads (3221 and 4477 kg ha1 a1) via surface runoff and drainflow in a two-year study in a clayey field (Sjo¨kulla) in southern Finland. The Sjo¨kulla area had exceptional rains (129 mm) in October 1998 and record-breaking rains (234 mm) in November 1996 (the average November rain in that particular area of Finland is 90 mm). In our study, the loads were low even though the field area was ploughed in the autumn, except in 1984, when LU2 was presumably lightly cultivated or harrowed after the fallow period before seeding (Table 1). There are several factors that can explain why the sediment loads in the simulations were low compared to other studies. We want to emphasize that our simulation periods included only the growing seasons and the following autumns. Puustinen et al. (2007) reported annual measured sediment loads of approximately 1400 kg ha1 in the Hovi monitoring field, which were mainly generated during the snow melt period in the spring. In the study by Bechmann (2012), soil losses occurred during snowmelt in March–April 1999, when high concentrations and high runoff volumes caused very high sediment losses. However, the amount of rainfall in October 1984 (144 mm) is comparable to the rainfall in October 1998 reported by Warsta et al. (2013b), which produced heavy erosion in their simulations. Our sediment load results are plausible because the model was able to reproduce the measured tillage layer runoff and drainflow amounts (Fig. 3), and sediment concentrations in the open ditch and subsurface drain discharge waters (Figs. 7–11 and 12). Some uncertainties were introduced

into the results due to the sediment concentration measurements conducted at the end of the open ditch (100 m) (Fig. 1). Erosion and deposition occur in agricultural brooks and are affected by channel dimensions and flow and vegetation conditions (e.g. Va¨stila¨ and Ja¨rvela¨, 2011). It is possible that the sediment concentrations changed along the ditch in the Hovi field, which might be worthwhile to investigate in future studies. During our simulation periods, the measured concentrations in the open ditch ranged between 0.007 and 1.3 g l1 (mean, 0.217 g l1) in 1988, and 0.011 and 0.15 g l1 (mean 0.063 g l1) in 1984. In drainflow, measured concentrations varied between 0.11 and 0.0026 g l1 (mean 0.040 g l1). In 1984, simulated sediment concentration peaks in the open ditch in October were an order of magnitude higher (Fig. 12) compared to the concentrations in drainflow (Fig. 11). In the model this difference was caused by the transport of detached particles from the field area to the open ditch by horizontal subsurface flow in the tillage layer. In the Hovi catchment of this study, Puustinen et al. (2007) measured seasonal flow-weighted mean sediment concentrations of 0.430– 0.485, 0.120–0.160, 0.870–1.295 and 0.140–0.145 g l1 in the main ditch during autumn, winter, spring and summer, respectively. Turtola et al. (2007) suggested that low concentrations in drainflow were due to sediment filtering in the soil profile during transport from the field surface to the subsurface drains, which is not described in FLUSH. Puustinen et al. (2005) measured sediment concentration of tillage layer runoff in a clay loam field (Aurajoki experimental field) in western Finland and reported concentrations of 1.298, 0.691 and 0.576 g l1 in autumn, winter and spring, respectively, with normal ploughing and a 8–9% slope. The range of concentrations measured in this study was low compared to the levels listed above, which indicates that the simulation periods represent conditions with low erosion potential. Also, only a few concentration measurements were available and they were partly collected during a low flow period in 1984. The lack of surface runoff was puzzling compared to results from other clayey fields in Finland (e.g. Turunen et al., 2013; Warsta et al., 2013a), because there certainly was enough precipitation in October 1984 to produce overland flow. The Sjo¨kulla field (Warsta et al., 2013a,b) has steeper parts compared to the Hovi field, and is therefore more susceptible for concentrated overland flow on the slopes, and a flat area close to the outlet where water and suspended sediment could accumulate. As a result, hydraulic erosion produced bigger sediment loads in the comparable simulation made in the Sjo¨kulla field (Warsta et al., 2013a,b). Hydraulic erosion due to overland flow was initially included in the simulations but, during the model calibration, hydraulic erosion was excluded and erosion was modelled merely with raindrop splash erosion. The exclusion of hydraulic erosion in the model does not prove that it has no role in erosion processes but the small amount of tillage layer runoff measured and the absence of the direct measurement of surface runoff quality did not support the detailed identification of hydraulic erosion. An inspection of our results revealed that cumulative surface runoff (a component of tillage layer runoff in the model) was only 11 and 17 mm in 1988 and 1984, respectively. Taskinen and Bruen (2007a,b) ended up using a similar assumption in their eventbased sediment transport simulations of the Hovi monitoring field. It is possible that the vertical gravel drains increased the conductivity of the profile and decreased surface runoff in the Hovi field (Bengtsson et al., 1992). The different surface areas of the drainage basins might also cause differences in surface runoff generation. When the total modelled water balance of the Hovi field was assessed, groundwater outflow was found to be a major component of the water balance (15% and 16% of precipitation in 1988 and 1984, respectively), which is somewhat surprising in clayey

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

soils. Warsta et al. (2013a) estimated with FLUSH that 14% and 11% of precipitation were lost via groundwater outflow in 1998 and 1996, respectively (in the Sjo¨kulla field), while Turunen et al. (2013) reported groundwater outflow proportions of 9–15% of the precipitation amounts in another clayey field in southern Finland. Three other experimental, clayey, subdrained field sites – including Sjo¨kulla (Alakukku et al., 2010), Nummela and Ga˚rdskulla Ga˚rd (Vakkilainen et al., 2010) are located in the vicinity of the Hovi field, and it was possible to compare the soil properties of the Hovi field against them. The clay content in the tillage layer (0– 0.2 m) was higher (0.67) in Nummela compared to the other fields, while in the bottom soil (0.45 m) the clay content was lower (0.59) and silt content was higher (0.29) in Ga˚rdskulla. Organic carbon contents were similar between the fields. Saturated hydraulic conductivities (KS) in the tillage layer in the Hovi field were two and three times lower than in Sjo¨kulla and Nummela respectively but, in the bottom soil, the average conductivity in Sjo¨kulla was only half of the Hovi value. However, the macroporosities (w) in Hovi were lower compared to the Sjo¨kulla (Alakukku et al., 2010) and Nummela (Vakkilainen et al., 2010) values in both the tillage layer and bottom soil. The effect of gravel drains was taken into account in this study by increasing w by 0.001 from the field surface to the level of subsurface drains (1.0 m). In our model applications, static w values (macroporosity not affected by soil shrinkage and swelling) were decreased to 40% of the original values during both years after tillage. Warsta et al. (2013a) decreased the static w values in the Sjo¨kulla field to 25% and 67% of the original values in 1998 and 1996, respectively. The comparison of the fields revealed differences between their soil properties but the differences as such did not explain the low summer and autumn sediment loads in Hovi. The use of coarse gravel as trench fill and the construction of gravel traps from drain level to soil surface remain a likely factor affecting the runoff generation mechanism and explaining the low levels of surface runoff and erosion in Hovi. Warsta et al. (2013b) noted that the sediment load was more sensitive to flow model parameters controlling the shares of surface runoff and subsurface drainflow than to erosion model parameters. Jetten et al. (1999) already noted over a decade ago that modellers are usually more concerned with the outgoing fluxes from the source area than with investigating the main processes controlling water movement and erosion within the area, which still seems to be true. Also, Ule´n and Bechmann (2012) found the lack of information on spatially distributed erosion to be a gap in erosion knowledge in high latitudes. 3D distributed models like FLUSH can simulate horizontally occurring erosion on the field surface and also vertical suspended sediment transport in the soil profile to the subsurface drains. We presented spatially distributed erosion results along with total loads, although we did not have spatially distributed data. However, we still found the results useful as the model enabled us to divide the field into different land use areas to describe cultivated and fallow areas during the calibration and validation years. The land use areas had different crop root development phases and were tilled at different times. Also, different features such as subsurface drains and the depression in the field had impact on the spatially distributed erosion results (Fig. 4). The vertical dimension enabled us to simulate preferential transport of suspended sediment in varyingly saturated soil profiles, which is only available in a few 1D models (Larsson et al., 2007). With previously published 1D models (e.g. Knisel and Turtola, 2000; Larsson et al., 2007; Ba¨rlund et al., 2009), it is not possible to describe spatially distributed moisture conditions or erosion, divide the field into separate land-use areas or investigate sediment source areas within the field.

35

Agricultural mitigation measures are implemented in the field scale and FLUSH facilitates future modellers to explicitly parameterize and test the impacts of erosion control measures, such as vegetated buffer strips on sediment sieving, in the field area. The effectiveness of this protection measure can be lower than expected, especially in subsurface drained fields where the notable sediment loads are lost via subsurface drainflow and surface runoff (e.g. Turtola et al., 2007; Vakkilainen et al., 2010; Warsta et al., 2013b). Spatial information provided by the model can aid the experimental design of erosion measurements. For example, aggregated measurements of sediment transport to drainflow could be conducted in subareas of the field. The selection of the subareas can be based on the simulation of erosion hot spots. Topographical changes due to erosion can be directly measured with light detection and ranging (LIDAR) systems (e.g. Lilja, 2010). The results show that considerable erosion can occur during a short interval in response to a single major storm event. LIDAR measurements before and after such an event in erosion prone areas would provide useful information for the validation of the modelled distribution of erosion. Currently, the simulation times of FLUSH are relatively long, especially with the higher resolution grid (104  144  16 cells in x-, y- and z-directions) applied here. Nonetheless, automatic calibration methods are needed for robust parameter sensitivity analysis and the calibration process in FLUSH. Methods such as the ensemble Kalman filter approach (e.g. Erdal et al., 2014) provide promising options for parameter identification in future studies. 5. Conclusions The distributed hydrological model, FLUSH, was applied to a clayey, subsurface drained agricultural field (abundant around the Baltic Sea area) to assess water balance and soil erosion during two growing seasons and the following autumns. The model divides the simulated field into 2D overland and 3D subsurface domains, which facilitated simulation of overland flow and erosion together with subsurface flow and suspended sediment transport in soil macropores. The model was able to decompose the water balance into individual components (tillage layer runoff, ET, drainflow, groundwater outflow and storage changes) and take into account the effects of separate drainage systems in the field area. The dualpermeability subsurface flow model enabled us to use low saturated hydraulic conductivity measurements in the matrix pore system (which is not possible in single pore system models) but still reproduce fast measured drainflow events due to preferential flow in macropores after precipitation events. The simulated sediment loads via tillage layer runoff and drainflow, and consequently erosion, were relatively low during the simulation periods (May–October) compared to the previously published results. According to the simulations, most of the suspended sediment was detached by raindrop impacts as concentrated overland flow, and consequently hydraulic erosion was absent in the field. Vertical gravel drains decreased surface runoff in this field in the simulations. The measurement system, in which sediment concentrations were recorded at the end of the open ditch (100 m), might have modified the measured concentrations and loads, and should be inspected in the future. Simulated spatially distributed erosion results helped us to determine erosion processes occurring in different parts of the field and these results could be used to determine optimal measurement locations within a field for intensive spatially distributed erosion measurements. Our tests with the FLUSH model indicate that the model is well suited for describing fields with complex topography, structured soils (such as clay or till) and several separate drainage systems, which would all be difficult to simulate with 1D models.

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37

36

Acknowledgements

Overland erosion is simulated the with the sediment continuity equation (e.g. Warsta et al., 2013b):

The data sets were provided by several organizations, including the Finnish Environment Institute, MTT Agrifood Research Finland and the Finnish Meteorological Institute. The authors would like to thank Liisa Pesonen, Antti Ristolainen and Jere Kaivosoja for the possibility of using the MTT data set collected from the Hovi monitoring field in this study. We would like to express our gratitude to Prof. Seppo Rekolainen for giving us an opportunity to use the runoff data set from the Hovi field. The authors are grateful to Prof. Pertti Vakkilainen for his support and advice during this study. We would like to thank Prof. Nicholas Jarvis, Prof. Jirˇı´ Sˇimu˚nek and Prof. Eila Turtola for their valuable feedback and advice. This study was partly conducted in an Academy of Finland project ‘‘Flow–vegetation–sediment interaction: physical, chemical, and biological processes in environmental channels’’ led by Dr. Juha Ja¨rvela¨. The study was funded by Maa- ja vesitekniikan tuki ry., the Academy of Finland, Salaojituksen tukisa¨a¨tio¨ (the Finnish Drainage Foundation), the Ministry of Agriculture and Forestry and Sven Hallinin tutkimussa¨a¨tio¨ (the Sven Hallin Research Foundation). We would like to thank the anonymous reviewers for their constructive comments in the preparation of the manuscript. We acknowledge CSC – IT Center for Science Ltd. for the allocation of computational resources.

@ðhW CÞ @ðQ X CÞ @ðQ Y CÞ þ þ ¼ g R þ g H  wS U S C  SS @t @x @y

(5)

where C (kg mS3) is the sediment concentration in the overland water, QX and QY (m2 sS1) are the overland unit flow rates, gR (kg mS2 sS1) is the potential raindrop splash erosion rate, gH (kg mS2 sS1) is the potential hydraulic erosion rate, wS (–) is the coefficient for the particle settling, US (m sS1) is the particle settling velocity and SS (kg mS2 sS1) is the overland sediment source/sink, e.g. open ditches. Sediment transport in the macropore system is depicted with the advection dispersion equation (e.g. Warsta et al., 2013b):       @ðuCÞ @ @C @ @C @ @C @ ¼ þ þ  uDXX uDYY uDZZ @x @y @z @x @t @x @y @z  ðqx cÞ 

@ @ ðq cÞ  ðqz cÞ  sS @y y @z

(6)

where c (kg mS3) is the sediment concentration in the subsurface domain, DXX, DYY and DZZ (m2 hS1) are the dispersion coefficients, q (m hS1) is the unit (Darcian) water flux and sS (kg mS3 hS1) is the sediment source/sink. References

Appendix A. Governing equations Overland water flow is presented with the diffusion wave simplification (e.g. Warsta et al., 2013a):

@hW @ðvX hW Þ @ðvY hW Þ þ þ ¼ sW @t @x @y 1 n

vX ¼ h2=3 W

(1)

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi @zS @hW 1 2=3 @zS @hW   and vY ¼ hW n @x @x @y @y

(2)

where hW (m) is the overland water depth, t (h) is the time, v (m sS1) is the overland flow velocity, sW (m sS1) is the source/sink term in the surface domain, n (–) is the Manning’s coefficient and zS (m) is the elevation of the field surface. hW is calculated as the difference between the total overland water depth HW,TOT (m) and the overland flow threshold depth HW,THR (m). HW,THR depicts overland micro depression storage, i.e. overland flow ensues when HW,TOT S HW,THR > 0.0 m. Subsurface water flow is presented with a dual-permeability approach applying the Richards equation in matrix and macropore systems (e.g. Warsta et al., 2013a): 











@uF @ @HF @ @HF @ @HF GW ¼ K þ K þ K  @t @x FX @x W @y FY @Y @z FZ @z  SWF 

(3) 







@uM @ @HM @ @HM @ @HM ¼ K þ K þ K @t @x MX @x @y MY @Y @z MZ @z GW þ

1W

 SWM



(4)

where u (m3 mS3) is the volumetric water content, K (m hS1) is the unsaturated hydraulic conductivity, H (m) is the hydraulic head, GW (hS1) is the water exchange rate between the pore systems, w (m3 mS3) is the fraction of the macropores of the total porosity (i.e. macroporosity) and SW is the source/sink term in the subsurface domain. Subscripts F and M refer to the macropore and matrix systems, respectively.

Alakukku, L., Nuutinen, V., Ketoja, E., Koivusalo, H., Paasonen-Kiveka¨s, M., 2010. Soil macroporosity in relation to subsurface drain location on a sloping clay field in humid conditions. Soil Till. Res. 106, 275–284. Bechmann, M., 2012. Effect of tillage on sediment and phosphorus losses from a field and a catchment in south eastern Norway. Acta Agric. Scand. B – Soil Plant Sci. 62, 206–216. Bechmann, M., Sta˚lnacke, P., 2005. Effect of policy-induced measures on suspended sediments and total phosphorus concentrations from three Norwegian agricultural catchments. Sci. Total Environ. 344, 129–142. Bengtsson, L., Seuna, P., Lepisto¨, A., Saxena, R.K., 1992. Particle movement of melt water in a subdrained agricultural basin. J. Hydrol. 135, 383–398. Boardman, J., Poesen, J. (Eds.), 2006. Soil Erosion in Europe. Wiley, West Sussex, England, p. 878. Ba¨rlund, I., Kirkkala, T., Malve, O., Ka¨ma¨ri, J., 2007. Assessing SWAT model performance in the evaluation of management actions for the implementation of the Water Framework Directive in a Finnish catchment. Environ. Model. Softw. 22, 719–724. Ba¨rlund, I., Tattari, S., Puustinen, M., Koskiaho, J., Yli-Halla, M., Posch, M., 2009. Soil parameter variability affecting simulated fieldscale water balance, erosion and phosphorus losses. Agric. Food Sci. 18, 402–416. Chapman, A.S., Foster, I.D.L., Lees, J.A., Hodgkinson, R.A., 2005. Sediment delivery from agricultural land to rivers via subsurface drainage. Hydrol. Proc. 19, 2875–2897. Erdal, D., Neuweiler, I., Wollschla¨ger, U., 2014. Using a bias aware EnKF to account for unresolved structure in an unsaturated zone model. Water Resour. Res. 50, 132–147. Feddes, R.A., Kowalik, P.J., Zaradny, H., 1978. Simulation of field water use and crop yield. Wageningen, Netherlands, pp. 189. Førland, E.J., Allerup, P., Dahlstro¨m, B., Elomaa, E., Jo´nsson, T., Madsen, H., Pera¨la¨, J., Rissanen, P., Vedin, H., Vejen, F., 1996. Manual for Operational Correction of Nordic Precipitation Data. Norwegian Meteorol. Inst., Report 24, KLIMA, pp. 66. Granlund, K., Ra¨ike, A., Ekholm, P., Rankinen, K., Rekolainen, S., 2005. Assessment of water protection targets for agricultural nutrient loading in Finland. J. Hydrol. 304, 251–260. Jacobsen, O., Moldrup, P., Larsen, C., Konnerup, L., Petersen, L., 1997. Particle transport in macropores of undisturbed soil columns. J. Hydrol. 196, 185–203. Jetten, V., Roo, A., Davis-Mortlock, D., 1999. Evaluation of field-scale and catchmentscale soil erosion models. Catena 37, 521–541. Joel, A., Wesstro¨m, I., Strock, J., Messing, I., 2012. Method for in situ measurements of water, sediment and phosphorous transport in the upper soil profile. Acta Agric. Scand. B – Soil Plant Sci. 62, 252–259. Knisel, W.G., Turtola, E., 2000. Gleams model application on a heavy clay soil in Finland. Agric. Water Manage. 43, 285–309. Larsson, M.H., Persson, K., Ule´n, B., Lindsjo¨, A., Jarvis, N.J., 2007. A dual porosity model to quantify phosphorus losses from macroporous soils. Ecol. Mod. 205, 123–134. Laubel, A., Jacobsen, O.H., Kronvang, B., Grant, R., Andersen, H.E., 1999. Subsurface drainage loss of particles and phosphorous from field plot experiments and a tile-drained catchment. J. Environ. Qual. 28, 576–584. Lilja, H., 2010. Erosion mapping with Light Detection and Ranging (LIDAR) and RUSLE – method testing at experimental plots and farmers’ fields. In: Turtola, E., Ekholm, P., Chardon, W. (Eds.), Novel methods for reducing agricultural

L. Warsta et al. / Soil & Tillage Research 143 (2014) 26–37 nutrient loading and eutrophication. Meeting of Cost, vol. 869. 14–16 June, Jokioinen, Finland. http://jukuri.mtt.fi/handle/10024/477234. Lundekvam, H.E., 2007. Plot studies and modelling of hydrology and erosion in southeast Norway. Catena 71, 200–209. McKay, L.D., Gillham, R.W., Cherry, J.A., 1993. Field experiments in a fractured clay till 2. Solute and colloid transport. Water Resour. Res. 29, 3879–3890. Muukkonen, P., Hartikainen, H., Alakukku, L., 2009. Effect of soil structure disturbance on erosion and phosphorus losses from Finnish clay soil. Soil Till. Res. 103, 84–91. Nash, J.E., Sutcliffe, J.V., 1970. River flow forecasting through conceptual models. Part I. A discussion of principles. J. Hydrol. 10, 282–290. Øygarden, L., Kværner, J., Jenssen, P.D., 1997. Soil erosion via preferential flow to drainage systems in clay soils. Geoderma 76, 65–86. Paasonen-Kiveka¨s, M., Vakkilainen, P., Karvonen, T., 2008. Nutrient transport through tile drains on a clayey field. In: 10th International Drainage Workshop of ICID Working Group on Drainage, Helsinki–Tallinn, Finland–Estonia, 6–11 July 2008, Conference proceedings, pp. 142–152. Posch, M., Rekolainen, S., 1993. Erosivity factor in the universal soil loss equation estimated from Finnish rainfall data. Agric. Sci. Finland 2, 271–279. Puustinen, M., Koskiaho, J., Peltonen, K., 2005. Influence of cultivation methods on suspended solids and phosphorus concentrations in surface runoff on clayey sloped fields in boreal climate. Agric. Ecosyst. Environ. 105, 565–579. Puustinen, M., Tattari, S., Koskiaho, J., Linjama, J., 2007. Influence of seasonal and annual hydrological variations on erosion and phosphorus transport from arable areas in Finland. Soil Till. Res. 93, 44–55. Puustinen, M., Turtola, E., Kukkonen, M., Koskiahoa, J., Linjama, J., Niinioja, R., Tattari, S., 2010. VIHMA – a tool for allocation of measures to control erosion and nutrient loading from Finnish agricultural catchments. Agric. Ecosyst. Environ. 138, 306–317. Seuna, P., Kauppi, L., 1981. Influence of sub-drainage on water quantity and quality in a cultivated area in Finland. Publications of the Water Research Institute, National Board of Waters, Finland, No. 43, , pp. 32–47. Skøien, S.E., Børresen, T., Bechmann, M., 2012. Effect of tillage methods on soil erosion in Norway. Acta Agric. Scand. B – Soil Plant Sci. 62, 191–198. Taskinen, A., 2002. Mathematical Modelling of Overland Flow, Erosion and Transport of Phosphorus. (Ph.D. thesis), University College Dublin, Dublin, Ireland. Taskinen, A., Bruen, M., 2007a. Incremental distributed modeling investigation in a small agricultural catchment. 1. Overland flow with comparison with the unit hydrograph model. Hydrol. Proc. 21, 80–91. Taskinen, A., Bruen, M., 2007b. Incremental distributed modeling investigation in a small agricultural catchment. 2. Erosion and phosphorus transport. Hydrol. Proc. 21, 92–102. Tattari, S., Ba¨rlund, I., Rekolainen, S., Posch, M., Siimes, K., Tuhkanen, H-R., Yli-Halla, M., 2001. Modeling sediment yield and phosphorus transport in Finnish clayey soils. Trans. ASAE 44, 297–307. Turtola, E., Alakukku, L., Uusitalo, R., 2007. Surface runoff, subsurface drainflow and soil erosion as affected by tillage in a clayey Finnish soil. Agric. Food Sci. 16, 332–351.

37

Turtola, E., Paajanen, A., 1995. Influence of improved subsurface drainage on phosphorus losses and nitrogen leaching from a heavy clay soil. Agric. Water Manage. 28, 295–310. Turunen, M., Warsta, L., Paasonen-Kiveka¨s, M., Nurminen, J., Myllys, M., Alakukku, L., A¨ijo¨, H., Puustinen, M., Koivusalo, H., 2013. Modeling water balance and effects of different subsurface drainage methods on water outflow components in a clayey agricultural field in boreal conditions. Agric. Water Manage. 121, 135–148. Ule´n, B., Bechmann, M., Øygarden, L., Kyllmar, K., 2012. Soil erosion in Nordic countries – future challenges and research needs. Acta Agric. Scand. B – Soil Plant Sci. 62, 176–184. Ule´n, B., Aronsson, H., Bechmann, M., Krogstad, T., Øygarden, L., Stenberg, M., 2010. Soil tillage measures against phosphorus erosion and leaching and potential side-effects – a Scandinavian review. Soil Use Manage. 26, 94–107. Uusitalo, R., Turtola, E., Kauppila, T., Lilja, T., 2001. Particulate phosphorus and sediment in surface runoff and drainflow from clayey soils. J. Environ. Qual. 30, 589–595. Uusitalo, R., Turtola, E., Lemola, R., 2007. Phosphorus losses from a subdrained clayey soil as affected by cultivation practices. Agric. Food Sci. 16, 352–365. Vakkilainen, P., 1982. On the Estimation of Evapotranspiration. (Doctoral dissertation) Series C Technica No. 20, Artes Constructionum No.6, University of Oulu, Oulu, Finland, pp. 146. Vakkilainen, P., Alakukku, L., Koskiaho, J., Myllys, M., Nurminen, J., PaasonenKiveka¨s, M., Peltomaa, R., Puustinen, M., A¨ijo¨, H., 2010. Pellon vesitalouden optimointi, Loppuraportti 2010. Salaojituksen tutkimusyhdistys ry:n tiedote 30, Helsinki 2010. Salaojituksen tutkimusyhdistys ry, pp. 114. Van Dam, J.C., Groenendijk, P., Hendriks, R.F.A., Kroes, J.G., 2008. Advances of modeling water flow in variably saturated soils with SWAP. Vadose Zone J. 7, 640–653. van Genuchten, M.T., 1980. A closed form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J. 44, 892–898. Vuorenmaa, J., Rekolainen, S., Lepisto¨, A., Kentta¨mies, K., Kauppila, P., 2002. Losses of nitrogen and phosphorus from agricultural and forest areas in Finland during the 1980 and 1990. Environ. Monit. Assess. 76, 213–248. Va¨stila¨, K., Ja¨rvela¨, J., 2011. Environmentally preferable two-stage drainage channels: considerations for cohesive sediments and conveyance. JRBM 9, 171–180. Warsta, L., 2011. Modelling Water Flow and Soil Erosion in Clayey, Subsurface Drained Agricultural Fields, vol. 82 (Doctoral dissertations), Aalto University, Espoo, Finland, pp. 212. Warsta, L., Karvonen, T., Koivusalo, H., Paasonen-Kiveka¨s, M., Taskinen, A., 2013a. Simulation of water balance in a clayey, subsurface drained agricultural field with three-dimensional FLUSH model. J. Hydrol. 476, 395–409. Warsta, L., Taskinen, A., Koivusalo, H., Paasonen-Kiveka¨s, M., Karvonen, T., 2013b. Simulation of soil erosion in a clayey, subsurface drained agricultural field with a three-dimensional FLUSH model. J. Hydrol. 498, 132–143.