Glyphosate dynamics prediction in a soil under conventional and no-tillage systems during the crop cycle

Simulation models are efficient tools to predict the fate of different solutes in agricultural soils. This work aimed to compare measured and predicted glyphosate and AMPA (aminomethyl phosphonic acid; its main metabolite) contents in a soil under no-tillage (NT), and conventional tillage (CT); and to compare the predictions considering constant and time-variable hydraulic properties. Additionally, we evaluated the ability of the model to predict glyphosate and AMPA accumulation during the crop cycle. Hydrus 1-D code was used to predict the glyphosate and AMPA dynamics, considering constant and time-variable hydraulic properties during the studied crop cycle. In general, the prediction of glyphosate and AMPA distribution along the soil profile using HYDRUS 1-D was satisfactory; however, an overestimation of both compounds was observed in the first 0.20 m of the soil probably because of the preferential flow. Additionally, the accumulation process of glyphosate and AMPA in the soil during the crop cycle was underestimated by HYDRUS 1-D, as compared with the observed field data. Simulated data show that higher values of K0 increase the risk of glyphosate and AMPA vertical transport. The inclusion of temporal variation of hydraulic properties in glyphosate and AMPA simulation did not improve the simulation performance, showing that the model is more sensitive to the parameters related to the solutes. From the obtained results, HYDRUS 1-D code allowed to predict glyphosate and AMPA dynamics reasonably well in agricultural soils of the Argentinean Pampas region and is a potential model to give support in the analysis of the environmental risk of leaching and soil contamination.


INTRODUCTION
The prediction of the impact of different tillage systems on the environment and the fate of different pesticides is a complex and necessary task to improve water and soil management. However, there is a knowledge gap in the prediction and modeling of the fate of different pesticides in agricultural soils under different tillage systems (Elias et al., 2018). Glyphosate (N-[phosphonomethyl] glycine) is a broad-spectrum herbicide, used non-selectively in agriculture to control weeds and herbaceous plants, especially under no-tillage (NT) management systems. In Argentina, the use of this herbicide has increased drastically since its introduction, and it is the most commonly used herbicide in the country (Primost et al., 2017). In the last years, several authors reported the presence and the potential mobility of glyphosate and its main metabolite (aminomethyl phosphonic acid, AMPA) in Argentinian agricultural soils (Peruzzo et al., 2008;Lupi et al., 2015;Okada et al., 2016;Soracco et al., 2018a); however, these results are site-specific. Furthermore, the procedure for the quantification of these molecules is time-consuming and expensive (Baer and Calvet, 1999). In this sense, the use of different simulation models provides an efficient alternative, which can be adapted to different soil types and management (van Genuchten et al., 1974;Walker, 1987;Armstrong et al., 2000;Jarvis et al., 2000;Pang et al., 2000;Worrall and Kolpin, 2004;Šimůnek et al., 2005;Köhne et al., 2006;Šimůnek and van Genuchten, 2008;Gupta et al., 2012).
In recent years, different models have been used to predict water flux and solute transport in soil. In general, these models solve the continuity of Richards' equation and the convection-dispersion equation. A large number of simulation studies have been developed to predict the dynamics of a wide range of agrochemicals such as pyrethroids, Cyanazine, Metribuzin, Atrazine, isoproturon, among others (Sadeghi et al., 1995;Ahuja et al., 2000;Bayless et al., 2008;Luo and Zhang, 2011;Filipović et al., 2016). Additionally, agrochemicals registration and leaching risk are supported by different simulation models, such as PRMZ and PEARL (Dubus et al., 2003). In this sense, simulation models are a useful tool to support authorities in decisions concerning the approval of pesticide registration at the European Union level and in the United States (Scorza Júnior and Boesten, 2005).
For glyphosate and AMPA modeling, Mamy et al. (2008) studied their lixiviation under field conditions, comparing measured and predicted data obtained using the PRZM model, obtaining a good correlation. Klier et al. (2008) reported similar results in sandy soils from data obtained using LEACHP code. These authors mentioned that the risk of glyphosate lixiviation is very low due to its high adsorption coefficient and degradation rate. Candela et al. (2007) also suggested strong bounds of glyphosate and soil in their study based on laboratory test modeled with a non-commercial software; however, later in 2010 they detected, in a field condition experiment, the presence of glyphosate and AMPA in unexpectedly high depths suggesting preferential flow and/or colloidal mediated transport as the potential causes (Candela et al., 2010). Laitinen et al. (2007) reported that the observed and simulated glyphosate residues by PEARL model in soil did not correlate. These authors attributed their results to the fact that glyphosate residues must originate from exudation from plant roots, and that the translocation process could not be described by the model. In the Pampas Region, Peruzzo et al. (2008) used the SoilFug model to predict the glyphosate distribution in agricultural soils. These authors found a very good agreement between the model prediction and the measured field data.
In recent years, HYDRUS code (Šimůnek et al., 2013) has been successfully used to model and predict soil water movement (Garg and Ali, 2000) and several solutes transport (Coquet et al., 2005;Dousset et al., 2007;Kodešová et al., 2008;Valdés-Abellán et al., 2014;Alletto et al., 2015) in different scenarios. The popularity of the HYDRUS code within soil physics researches is reflected by its increasing use in a variety of applications and Rev Bras Cienc Solo 2020;44:e0190130 publications (Šimůnek et al., 2016). Kodešová et al. (2008) studied the transport of the herbicide chlorotoluron in three different soil types. These authors found a good correlation between predicted and measured data, only when a dual-porosity model was included. Köhne et al. (2006) studied the transport of isoproturon and terbuthylazine (systemic herbicides) and observed a good prediction of their behavior using HYDRUS 1-D. On the other hand, Filipović et al. (2014), in a 6 years' field trial, reported that HYDRUS 1-D predicted the isoproturon vertical transport, except when events of preferential flux were observed.
The use of this model requires a good description of the soil water retention curve (WRC) and the hydraulic conductivity function [K (h)]. In general, these functions are taken as constant in time during the crop cycle in most simulation studies (Schwen et al., 2011a). However, several authors emphasized that soil physical properties show temporal variation (Alletto and Coquet, 2009;Hu et al., 2009;Schwen et al., 2011a;Jirků et al., 2013;Lozano et al., 2016;Villarreal et al., 2017;Soracco et al., 2018b) and that this temporal variation has often overshadowed any measured differences between management treatments and situations (Strudley et al., 2008). These temporal changes, especially in the near saturation range, could have a great impact on soil water and solute dynamics. Or et al. (2000) introduced a model to describe the temporal changes of the water retention properties after the tillage operation, based on the soil pore distribution. Xu and Mermoud (2003) used an empirical function to describe the decay of saturated hydraulic conductivity (K 0 ) in the soil water movement simulation. Beyond the studies mentioned above, the inclusion of hydraulic properties temporal variability is, in general, not common when solute transport is predicted. In addition, despite the wide use of HYDRUS code to model and predict different solutes transport, there are few studies about the performance of this code on glyphosate and AMPA modeling under real field conditions. In this sense, the introduction of the temporal variation of soil hydraulic properties during the crop cycle could improve the performance of HYDRUS code to predict glyphosate and AMPA dynamics in agricultural soils. We hypothesized that it is possible to predict glyphosate and AMPA contents using HYDRUS 1-D, and that the prediction is improved considering the time-variable hydraulic properties. The objectives of this study were to compare measured and predicted glyphosate and AMPA content values in a soil under NT and CT; to compare the predictions considering constant and time-variable hydraulic properties. Additionally, we evaluated the ability of the model to predict glyphosate and AMPA accumulation through the crop cycle.

Site and treatments
The experiment was carried out near the city of Chascomús, Argentina (located at 35° 44' 37" south and 58° 03' 10" west). The soil was classified as a fine, illitic, thermic abruptic Argiudoll (Soil Survey Staff, 2006), which corresponds to a Luvic Phaeozem (IUSS Working Group WRB, 2007). The climate in the region is temperate without frost. The mean annual precipitation is 946 mm and the mean annual potential reference evapotranspiration is 929 mm. Daily precipitation and air temperature were recorded during the experiment period (June 2015-August 2016) (Figure 1). The total rainfall during the experiment period was 1,178 mm.
Before the treatments were applied, the plots were under CT and with the same crop rotation for more than 20 years. In the year 2000, an experimental design with two treatments (plots of 30 m wide and 50 m long for each treatment) was applied: (a) no tillage (NT), in which only a narrow (0.05 m) strip of the soil was drilled to deposit crop seeds; (b) conventional tillage (CT) in which the soil was ploughed (disc plough+tooth Rev Bras Cienc Solo 2020;44:e0190130 harrow) at 0.20 m depth, and later smoothed using the tooth harrow each year in October. The experiment was developed during a glyphosate-resistant soybean period.
During the studied period, glyphosate Roundup ® UltraMax was applied three times (1.6 kg ha -1 active ingredient) on September 5th (before the CT plough), November

Soil hydraulic parameters
The K (h) was measured in the field in both treatments using a disk tension infiltrometer with a 12.5 cm-diameter base. Infiltration tests were carried out at three soil water tensions: -6, -3, and 0 cm. The saturated hydraulic conductivity (K 0 ) was calculated from the steady-state infiltration rates in accordance with the multiple-head method (Ankeny et al., 1991). Samples with a volume of 98 cm 3 (5 cm height and 5 cm diameter) were taken from the A horizon, avoiding rows and visible wheel tracks to determine the soil WRC. Ten replicates from each treatment and moments were collected. Water retention data at tensions, h, values of 0, -10, -30, -50, -70, -100, -300, and -15000 cm were determined on the undisturbed soil cores using a sand box apparatus for h values between 0 and -100, and a pressure chamber for h values < -300 cm. The θ s , α, and n parameters (called van Genuchten parameters, VG) of the soil water retention curve were optimized by using the RETC software (van Genuchten et al., 1991) by fitting the measured retention and hydraulic conductivity data. To reduce the amount of unknown variables, tortuosity parameter, l, and residual soil water content, θ r , were set constant at 0.5 and 0, respectively. The average bulk density of the A horizon was measured by using five 500 cm 3 soil cores (9 cm by 8.4 cm diameter) taken from each treatment and sampling moment. Water content in the field was measured at four layers: 0.00-0.10, 0.10-0.20, 0.20-0.30, and 0.30-0.40 m in each treatment and sampling date.

Glyphosate and AMPA determinations
To determine glyphosate and AMPA contents, the first 0.40 m of the A horizon was sampled and divided in four layers: 0.00-0.10, 0.10-0.20, 0.20-0.30, and 0.30-0.40 m. Glyphosate and AMPA quantifications were carried out according to Aparicio et al. (2013). The results were expressed as μg of glyphosate or AMPA per kg of dry soil. A full description of the methodology can be seen in Soracco et al. (2018a).

Numerical simulation
Glyphosate and AMPA transport simulations were carried out using HYDRUS 1-D code (Šimůnek et al., 2008) that predicts one-dimensional water flow and solute transport in unsaturated porous media. Non-steady water flux is described by the one dimensional Richards' equation (Equation 1): in which θ represents volumetric water content (L 3 L -3 ); h is pressure head (L); z is the vertical coordinate [L]; t is time [T], K is the unsaturated hydraulic conductivity (L T -1 ), and A (L T -1 ) represents a sink term (root water uptake).
The one dimensional solute transport is described by the advection-dispersion equation (Equation 2): in which: in which: c is the solute concentration in the liquid phase (M L -3 ); s is the solute concentration in the solid phase (M M −1 ), ρ is the soil bulk density (M L -3 ), q is the volumetric water flux density (L T -1 ), D is the effective dispersion coefficient (L 2 T -1 ), μw is the first-order degradation rate in the liquid phase (T -1 ), μ s is the first-order degradation rate in the solid phase (T -1 ), and K d is the adsorption coefficient (L 3 M -1 ). A Galerkin-type linear finite element scheme is used for the spatial distribution and an implicit finite difference scheme is used for the temporal variation of the variables in equations 1 and 2 (Gupta et al., 2012).
For each simulation, the considered soil depth was 0.50 m, dividing the soil profile in 49 elements with 50 nodes. To guarantee numerical stability, smaller size elements were used close to the topsoil, where water conditions vary more rapidly. The initial time step was taken as 0.001 day. The time period for simulations was from 12th June 2016 to 4th August 2017, with a daily temporal discretization.
The initial condition was the measured initial soil water content in the soil profile at time t = 0 at different depths. An atmospheric boundary condition was selected at the top of the soil profile and free drainage was considered at the bottom of the soil profile. Initial solute (Glyphosate and AMPA) concentration was set equal to the field measured values at time t = 0. Solute flux was applied as the upper boundary condition and zero concentration gradient was selected as the lower boundary condition (Gupta et al., 2012). Glyphosate uptake by roots was neglected, and decay was considered independent of soil temperature.

Water flow and solute transport parameters
Soil hydraulic parameters required as input data for HYDRUS 1-D code were determined from the soil WRC and field infiltration data. For each treatment, two different water dynamics and solute transport simulations were carried-out: with constant and time-variable hydraulic properties during the studied crop cycle. The constant hydraulic properties simulations were performed with the parameters measured in June 2015; the variable hydraulic properties simulations consisted of four different simulations corresponding to periods between sampling moments (Table 1). Obtained hydraulic output data from each simulation were connected with a subsequent simulation according to Schwen et al., (2011b), while the initial glyphosate and AMPA soil content were the final values obtained from the previous simulation to allow a continuous description. Adsorption coefficients for glyphosate and AMPA were calculated using the pedotransfer function proposed by De Gerónimo et al. (2018). Degradation rates for glyphosate and AMPA were obtained from Okada (2014) and Mamy et al. (2008), respectively (Table 1). Column miscible displacement experiments were carried out to obtain the dispersivity parameter (λ, L). Stainless-steel columns with a length of 5 cm and a diameter of 2 cm were used. The internal wall of columns was covered with a non-reactive material to provide good contact between the soil and the column wall. At both ends, 25 μm filters were used to avoid soil loss during the experiments.
The columns were uniformly packed under vibration with air-dried and 2 mm sieved soil (OECD, 2000) and slowly saturated upward with the electrolyte solution (KNO 3 0.05 mol L -1 ). A peristaltic pump was connected with polytetrafluoroethylene tubing to the columns, and a stainless-steel valve was used to switch between different inflow solutions. The experiments were conducted at a flux density of 3 cm h -1 , similar to hydraulic conductivity values obtained in the field. A volume equivalent to the column porosity of LiBr [0.01 mol L -1 (C 0 )] as a conservative tracer was injected in each experiment. Effluent solutions from the column were collected in 2 mL aliquots using an automatic fraction collector every 16 min. The electric conductivities of each effluent sample were measure using an EC electrode (Eijkelkamp ® ). From obtained electrical conductivity values of each sample and the LiBr solution, Brconcentration was determined following Larsbo et al. (2014). Three replicates were carried out. Finally, breakthrough curves (BTCs) were constructed from measured effluent concentrations. The Brconcentration values for different times were adjusted using HYDRUS 1-D code to obtain the soil λ.
All obtained BTCs were symmetric with a single peak, showing that bromide did not suffer any mechanism of adsorption, precipitation, nor anion exclusion during transport through the soil column (Candela et al., 2007). The single peak also certified the absence of preferential flow in the experimental setup.
Mass recovery was over 95 %, showing that bromide behaved as an ideal conservative tracer. The obtained values of λ were in the same order than previous reports for similar soils (Ersahin et al., 2002;Montoya et al., 2006;Candela et al., 2007;Bedmar et al., 2008;Okada et al., 2014) (Table 1).

Model performance analysis
From field measured values of glyphosate and AMPA contents during the studied crop cycle, the simulation performance of HYDRUS 1-D code was evaluated using different statistical indicators, according to Mamy et al. (2008): Correlation Coefficient (r) which is a measure of the degree of association between simulation and measurement: Eq. 4 in which O i and S i are the observed and simulated values respectively, Ô and Ŝ are the mean observed and predicted values respectively, and n is the number of sampling dates.
Root mean square error (RMSE), calculated as: Coefficient of residual mass (CRM), which indicates the consistent errors in the distribution of all simulated values across all measurements with no consideration of the order of the measurements:

Soil hydraulic properties during the crop cycle
The temporal variation of soil hydraulic properties during the crop cycle and its relationship with glyphosate and AMPA content can be seen in detail in Soracco et al. (2018a). In general, higher values of K 0 under CT during the crop cycle were observed, as compared with NT (Table 1). The K 0 showed no significant differences between treatments during the fallow periods (June 2015 and August 2016), while this variable was significantly higher under CT during the growing season from October 2015 to March 2016 (Table 1), especially after the tillage (October 2015). Under NT, K 0 values increased during the vegetative period. Regarding the VG parameters, CT treatment showed higher values of α as compared with NT only in October 2015, after tillage (Table 1). On the other hand, the n parameter showed constant values during the crop cycle, without differences between tillage systems (Table 1).

Measured glyphosate and AMPA contents
Overall, >60 % of glyphosate and 80 % of AMPA content along the studied soil profile for all sampling dates was found between 0.00-0.20 m (Figures 2 and 3). However, comparing two consecutive measuring dates, vertical transport was detected, especially when high precipitations occurred near the application. Glyphosate accumulation during the crop cycle under both tillage systems was observed due to high rates of glyphosate applications. A more detailed description can be found in Soracco et al. (2018a).

Glyphosate and AMPA simulation
In general, glyphosate and AMPA vertical distributions were well described by the model (Figures 2 and 3). For both tillage systems, was observed a good description of the soil profile distributions of glyphosate and AMPA (r >0.7), except in March 2016 under CT and in August 2016 under CT and NT (Figures 2 and 3). Under both tillage systems, the lower values of r were observed towards the end of the crop cycle, for both studied molecules, showing a low prediction capacity in long simulation periods. Furthermore, the values of the CRM indicator were lower than zero, showing an overestimation of glyphosate and AMPA contents (Table 2). This overestimation was higher for the glyphosate content under CT, due to the absence of measured glyphosate in August 2016. For both molecules, the overestimation was observed in the first 0.20 m of the soil profile, while between 0.20 and 0.40 m, the simulated concentrations were lower than the measured data. Lixiviation events, where the measured glyphosate contents were higher in the deeper depths as compared with the top soil, were not correctly predicted by the model (Figures 2d,  2e, 2f, 2j, and 2k). The RMSE values observed for glyphosate and AMPA were relatively high, attributed to the overestimation observed in October 2015 and August 2016 for glyphosate and AMPA, respectively.
Simulations with time-variable soil hydraulic properties did not improve the simulation performance, as compared with time-constant parameters, showing similar RMSE, CRM, and r values (Table 2 and Figures 2 and 3). However, a little higher performance in the prediction of glyphosate and AMPA contents when lixiviation occurred was observed  results are in disagreement with field measured data, which indicated an increase of 54 and 82 % during the crop cycle for NT and CT, respectively.

Soil hydraulic properties during the crop cycle
Higher values of K 0 under CT as compared with NT after the tillage can be attributed to the loosening effect, also reported by several authors (Angulo-Jaramillo et al., 1997;Strudley et al., 2008;Villarreal et al., 2017). These results are in agreement with several reports from the Pampas region, which showed that tillage effects on soil pore system configuration do not persist after harvest (Sasal et al., 2006;Soracco et al., 2010;Villarreal et al., 2017). The decrease of K 0 values in both treatments after harvest can be attributed to high traffic intensity associated with the harvest operations, which have been shown to damage the soil structure (Soracco et al., 2012). Under NT, increasing K 0 values during the vegetative period, which is in agreement with Schwen et al. (2011a) who mentioned that K 0 and water-conducting macroporosity increase in spring and summer due to higher biological activity and root growth. In addition, the decrease of K 0 values observed under both treatments after harvest can be attributed to high traffic intensity associated with the harvest operations, which have been shown to damage the soil structure (Soracco et al., 2012).
The behavior observed for VG parameters are in agreement with previous reports. Several authors reported that the α parameter is related to soil structure and it is affected by soil loosening (Schwen et al., 2011a;Jirků et al., 2013;Peña-Sancho et al., 2016). On the other hand, the n parameter showed constant values during the crop cycle, without differences between tillage systems (Table 1). This parameter is related to soil texture, showing small temporal variation (Jirků et al., 2013;Peña-Sancho et al., 2016). A more detailed discussion of the temporal variation of soil hydraulic properties can be found in Soracco et al. (2018a).

Measured glyphosate and AMPA contents
Glyphosate and AMPA vertical transport observed during the studied period, especially when high precipitations occurred near the application, is in agreement with several previous reports (Borggaard and Gimsing, 2008;Peruzzo et al., 2008). The results showed that the temporal dynamics of glyphosate and AMPA were related to the temporal variation of soil hydraulic properties. High values of K 0 may lead to lower retention of glyphosate and AMPA in the topsoil, favoring vertical transport to deeper soil layers. Glyphosate accumulation during the crop cycle under both tillage systems was observed due to high rates of glyphosate applications. A more detailed discussion can be found in Soracco et al. (2018a).

Glyphosate and AMPA simulation
In general, glyphosate and AMPA vertical distributions were well predicted by the model (Figures 2 and 3). The RMSE values observed for glyphosate and AMPA were relatively high, attributed to the overestimation observed in October 2015 and August 2016 for glyphosate and AMPA, respectively, as mentioned before.
The overestimation in the topsoil and the underestimation in the subsoil was also reported by Klier et al. (2008). Lixiviation events, when the measured glyphosate contents were higher in the greater depths as compared with the topsoil, were not correctly predicted by the model (Figures 2d, 2e, 2f, 2j, and 2k), probably because the model did not consider the preferential fluxes. Several authors have reported similar results and mentioned that isoproturon (Filipović et al., 2014) and glyphosate and AMPA (Mamy et al., 2008) distribution along the soil profile were well described by different models, except when preferential flux was observed. However, preferential flow may also keep the concentrations in the topsoil layer higher, since the preferential flow of water bypasses the soil matrix where glyphosate and AMPA are found. Furthermore, in this study, input parameters related to solute transport (i.e., K d , degradation rates and partition coefficients) were obtained from pedotransfer functions and the literature, and may lead to the observed underestimation.
Another possible reason is the underestimation of the λ values, due to the small length of the soil columns as compared with the studied soil profile. Vanderborght and Vereecken (2007) mentioned that soil dispersivities increase with increasing transport distance and scale of the experiment.
From the total simulated masses in each sampling date (Table 3), it was observed that during the applications, glyphosate increased, followed by a decrease until the end of the crop cycle. Coupled with the glyphosate decay, simulate AMPA values increased due to the degradation. However, this prediction by HYDRUS did not correlate with the measured values. Glyphosate and AMPA increments observed in March 2016, despite the absence of application, could be attributed to the release of retained herbicide in soybean and weeds residues (Mamy et al., 2016) and the roots exudation. Laitinen et al. (2007), using the PEARL model, observed that the simulated and observed values of glyphosate did not correlated, mentioning the importance of the translocation process in glyphosate fate. Also, these authors suggested that plant translocation of glyphosate to roots should be included in simulation models.
The slight improvement in the prediction with time-variable hydraulic properties when lixiviation occurred can be related to the higher values of K 0 and macroporosity. Higher values of K 0 and macroporosity may lead to lower retention of glyphosate and AMPA in the topsoil, favoring vertical transport to deeper soil layers during the crop cycle (de Jonge et al., 2000;Kjaer et al., 2005;Stone and Wilson, 2006;Soracco et al., 2018a).
Regarding glyphosate and AMPA accumulation, the HYDRUS 1-D predictions are in disagreement with field measured data. These results show that the model underestimates the glyphosate and AMPA accumulation along with the soil profile, as reported by several authors in different soil types under field conditions (Bento et al., 2016;Primost et al., 2017;Soracco et al., 2018a).
Although HYDRUS 1-D predictions overestimated glyphosate and AMPA values in the topsoil, and underestimated those values in the subsoil, in general, the predicted values were satisfactory, showing similar values of r and RMSE as compared with other authors (Mamy et al., 2008). On the other hand, the performance of the simulation did not improve with time-variable hydraulic parameters. Several authors found that the performance of the simulation of water content dynamics improved using time-variable hydraulic parameters (Schwen et al., 2011b;Feki et al., 2018). However, Klier et al. (2008) mentioned that solute transport models are more sensitive to the solute transport parameters than to hydraulic parameters. In addition, this predictive capacity decreased towards the end of the studied period, together with an underestimation of the accumulation process of glyphosate and AMPA along with the soil profile. However, despite not taking into account important processes such as preferential flow or glyphosate degradation dependence with soil temperature, the HYDRUS 1-D code, in general, allowed to predict reasonably well glyphosate and AMPA dynamics in agricultural soils of the Argentinean Pampas region.

CONCLUSIONS
HYDRUS 1-D code represents a simple and useful tool to study and predict the glyphosate and AMPA dynamics in agricultural soils under different tillage systems. Simulated data show that higher values of K 0 increase the risk of glyphosate and AMPA vertical transport. However, the prediction did not improve considering the time-variable hydraulic properties, indicating that the model is more sensitive to the parameters related to the solutes.
HYDRUS 1-D, the most popular simulation code within soil physics research, is a potential model to give support in the analysis of the environmental risk of leaching and soil contamination.