The study area is the experimental station (CDA 236) located in the Mnasra area, about 30 km north of the city of Kenitra, in northwestern Morocco (Fig. 1).

Mnasra region is a strip of land parallel to the Atlantic coast, of about 50,000 ha, located on the right bank of the Sebou River, 7 to 14 km wide and over about 50 km in length extending from the North of the Kenitra City. It is characterized by a diversity of crops, such as cereals, sugar beet, fodder and vegetables. The experimental site is equipped with a weather station that provided various data, such as rainfall and minimum and maximum temperatures (Fig. 1).

The experiment was carried out on a 200-m² field with flat topography, subdivided into ten basic 4 × 5m² plots. This study is focused on uncultivated plots on which are applied a Br and carbofuran flow tracer, and on an uncultivated and untreated control plot.

The site sandy soil was classified as an Alfisol, according to U S Soil Taxonomy (Soil Survey and Staff, 1999) and a sandy loam according to the USDA (USDA-SCS, 1975). Its physicochemical characteristics are presented in Table 1.

The tracer used is an anionic tracer (Br^–) form (KBr). It was applied with a dose of 61.77 g/m². The carbofuran C₁₂H₁₅NO₃ (2,3-dihydro-2,2-dimethylbenzofuran-7-yl-N-methyl carbamate) of the commercial product Furadan was applied with a dose of 4 g/m² based on active ingredient. Carbofuran has a solubility of 700 mg/L and a low vapor pressure, i.e., 3.41 × 10⁻⁶ mmHg at 25 °C. The half-life of carbofuran reaches values ranging from 25.7 to 107 days in soils with pH values and temperatures close to those of our experiments (El Mrabet, 2002). An injection mode was used in a slot corresponding to the injection of a constant concentration for a certain time. Spatial variability was considered by many measurement reiterations on average three times. The solutes were injected using a pulverizer, advancing at a constant speed while keeping the same nozzle flow rates to ensure a homogeneous application. The application of solutes was followed by a first irrigation during the monitoring period. The cumulative amounts of precipitation and irrigation were 62 mm and 160 mm, respectively. Irrigation was carried out with a sprinkler. Evapotranspiration was calculated on a daily basis from observed in situ data.

The spatio-temporal monitoring of physical and chemical parameters, such as humidity, temperature, bromide and pesticide concentrations in the soil was conducted for 90 days.

Samples were collected using an auger from the following depths: 0–20 cm, 20–40 cm and 40–60 cm. No trace of carbofuran or bromide was detected at all levels before treatment.

Determination of bromide in soil samples was carried out by ion chromatography using a Metrohm 881 Compact IC Pro, with a detection limit of 0.1 mg/L (Tirumalesh, 2008).

It should be noted that another method of bromide analysis was carried out by XRF following the procedure of Margui et al. (2009). The comparison between the two methods was statistically evaluated, leading to similar results (Hmimou et al., 2011).

The carbofuran was analyzed by liquid chromatography using a HPLC device: LCQ Advantage Max, Brand Thermo Electron equipped with a type BDS Hypersil C18 (150 × 4, 6 mm × 5 μm) column.

The hydrodynamic characterization of the plot soil is performed in situ using a disc infiltrometer with the negative charge method. Several infiltration measurements are conducted. For the retention curve h(θ), we use the expression of Van Genuchten (1980):

Water movement in the unsaturated soil can be described by Richard's equation as follows:

The transport of an interactive and non-conservative solute in the vadose zone is described on the basis of the advection–dispersion equation assumed in the model. It also incorporates a homogeneous reaction, i.e., pesticide degradation described by first-order kinetics:

R=1+ρKdθ and μ=μw+ρKdθμs, (5)

where R is the retardation factor, a dimensionless quantity, C is the solute concentration (M L⁻³), D_a is the apparent hydrodynamic dispersion coefficient (L² T⁻¹), q_z is the Darcian water flux (LT⁻¹), μ is the degradation coefficient (T⁻¹), ρ is the bulk density of the porous medium (M L⁻³), K_d is the distribution coefficient (L³ M⁻¹). μ_w and μ_s are the degradation constants in the solid and the liquid phase, respectively.

Before the beginning of the experiment, we measured the water contents for each soil profile, since we are dealing with a homogeneous soil. Consequently, we can define an initial uniform water content profile with an average value of θ₀ = 0.11 cm³/cm³. The solute is not present initially in the field. The initial conditions are expressed as:

At the upper boundary of the flow model, we applied a Neuman type flux boundary condition:

At the upper boundary of the transport model, the total advective and dispersive mass fluxes of solutes were prescribed as functions of time:

For the lower limit, we consider a gravity flow condition:

The solute convective flux is, in this limit, expressed by:

Richard equations and advection–dispersion equation (ADE) (Eqs. (3) and (4)) accompanied by their initial and boundary conditions were solved by a numerical finite difference scheme (Crank–Nicolson). Once all discretisations were established, a tridiagonal, and symmetric linear system of equations is obtained, which was resolved by the (SOR) algorithm. The non-linearity of equation (3) was resolved by initializing, at each time step t_i, the c(h) term at time t_i, and then doing as many iterations as necessary to reach an invariable solution.

At each time step, we resolved first equation (3), and then we used its solution in equation (4) to resolve it.

Parameters H_g, n, k_s and η were calibrated by comparing the simulated and measured water contents. Parameter λ was calibrated by comparing the simulated and measured bromide concentrations, while k_d, μ_s and μ_w were calibrated by comparing the simulated and measured concentrations of carbofuran.

The calibration was done by minimizing the RMSE error between measurement and simulation:

After several infiltrometer measurements at different pressures, we were able to reconstruct a part of the retention curve θ (h) (Fig. 2), the value θ_r was taken equal to zero, since the soil is very sandy and has a very low capacity to retain water in dry-state conditions. In Fig. 2, the experimental points are aligned with the fitted curve, hence justifying our choice of Van Genuchten's model (Van Genuchten, 1980). The results of this adjustment are presented in Table 2.

The measured K (h) during infiltration (Fig. 3) has a less dispersed curved shape. Table 2 shows the values of K_s and η obtained by fitting the Brooks and Corey model (Brooks and Corey, 1964) onto the measurement points; these values are close to the average values measured at Mnasra (Saâdi and Maslouhi, 2003).

It should be noted that the adjustment of Figs. 2 and 3 was carried out following the Levenberg–Marquardt Algorithm (LMA). The shape parameter η is often estimated in the literature as (2 + 3L)/L, where 0.2 < L < 3 (Løvoll et al., 2011); in our case, the L value was 0.7. The η value estimated in this study, comprised between 5.8 and 7.4, is very close to the results of Saâdi et al. (1999) and Tamoh and Maslouhi (2004) in the coastal part of Mnasra.

Solving Richard's equation yields solutions in h, which are converted into water content via relation θ (h). Simulations show that the water content of the soil could never mathematically go below its residual value θ_r, even if the internal evaporation phenomenon is taken into consideration. This confirms our consideration of an almost null residual water content (Bagarello and Iovino, 2012; Coppola et al., 2009; Mubarak et al., 2009; Saâdi and Maslouhi, 2003; Tamoh and Maslouhi, 2004). Moreover, considering the residual water content as null when it is very small allows us to cut down the number of parameters to be calibrated, and therefore to avoid over-parameterization. We performed a sensitivity analysis of the model results with respect to different values of the residual water content assigned to θ_r = 0.02, 0.04, 0.06. We found that the variation of the residual water content has only a very limited impact on the mass balance of carbofuran effect.

We measured and simulated water content profiles during three months at three soil depths, i.e., 0–20, 20–40, and 40–60 cm. The volumetric water content is characterized by small variations, especially in the lower depths. It is due to the soil sandy texture, characterized by low retention capacity, a potential drainage, and low initial volumetric water content (see Fig. 4). Fig. 4A shows that the model simulations are in good agreement with the experimental observations. In Fig. 4 B and C, the observed values are simulated with only a slight overestimation for the last sampling dates, probably due to the heterogeneity of the soil or error-related measurements.

The calibration error of h_g, m and η parameters (RMSE) gives a satisfactory value of 0.09 (Table 3). The k_s variation, within the limits of possible values for a soil like ours, does not affect the simulation results and, therefore, the value of RMSE. This allows us to consider that the hydrodynamic parameters are well represented. Furthermore, we find out that the model gives good predictions of the water content values, especially on the 0–20-cm profile.

Fig. 5 shows the temporal evolution of bromide concentrations for the three soil horizons studied. Therefore, the dispersivity value λ is 0.115 m (Table 3). The fact that the value of λ seems high may be due to physical phenomena neglected during modelling, including the three-dimensional nature of the flow across the field, the soil heterogeneity along the vertical profile, and the presence of preferential flows.

However, such a λ value is not surprising (Saâdi and Maslouhi, 2003), since according to the work of Nielsen et al. (1986), λ can take values ranging from 0.05 cm in a laboratory test, to 10 cm and higher for field conditions.

Table 3 shows an irrelevant λ = 0.48 value of the calibration error (RMSE), which is due to a slight fluctuation of the calculated values in the soil layers, especially in the lower ones. Nevertheless, the model yields results consistent with the observations. The calibration method seems therefore qualitatively acceptable, even for a heterogeneous soil.

Simulated and observed concentrations of carbofuran, for levels 0–20 cm, 20–40 cm, and 40–60 cm are shown in Fig. 6. A few days after application, and after 25 days of leaching, carbofuran is practically distributed throughout the soil matrix, and the maximum rate of recovery is in the 0–20 cm layer. This observation suggests the existence of important phenomena of dispersion on the sandy substrate, which appear more pronounced on the topsoil horizon. The adsorption process is probably related to a slight increase in the organic matter content in this level (Table 1).

In addition, the model provides satisfactory results in terms of pesticide dynamics, except for the increase of the RMSE calibration error to 0.54 (Table 3). This increase is probably due to changes in distribution coefficients K_d and to μ degradation in each layer. Bearing in mind that we considered a homogeneous soil as the modelling approach, the physicochemical characteristics of the pesticide were the same in all soil profiles, whereas the soil is actually heterogeneous. The values of the rate of degradation (μ_s and μ_w) inferred by calibration, and those calculated by the equation (μ_s = μ_w = ln2/t_1/2 = 0.693/t_1/2) (Sparks, 1989) show some similarities, as the half-life (t_1/2) of carbofuran derived from a previous study for a sandy soil similar to our soil gave a similar carbofuran half-life (El Mrabet, 2002).

Differences observed between simulated and measured concentrations are probably due to the following reasons:

• • errors in experimental measurements;
• • influence of other physical phenomena related to transport not considered in modeling;
• • we considered 1D model approach, whereas the real phenomena are in 3D.

It is appropriate to specify that the volatilization is not considered, since carbofuran has low volatility and is characterized by a low vapor pressure (3.41 × 10⁻⁶ mmHg at 25 °C). According to Lardier (Lardier, 1987), the volatilization of carbofuran is limited to 0.5% for sandy soils after two months of treatment.

A sensibility analysis shows little influence of θ_r on the water movement, and to a lesser degree, on carbofuran mass balance, especially on the carbofuran leaching flux.

Therefore, we used the model developed to simulate the transport of carbofuran for two different irrigation regimes, i.e., 222 mm and 350 mm. Fig. 7A and B present the simulation results of carbofuran transport in the unsaturated soil down to a 1-m depth in both irrigation schemes (222 mm and 350 mm) area. Fig. 7 shows the distribution of the pesticide 3, 26, 31, and 63 days after application.

According to both irrigation regimes, the maximum carbofuran concentration in soil reaches 300 mg/L and 84 mg/L, respectively. These peaks appear three days after carbofuran surface application. We can observe a decrease in the maximum concentration and a slight advance of the concentration front according to the increase of irrigation amount. The same phenomenon is observed for these concentration profiles 26 days after treatment.

The two profile graphs of carbofuran concentrations appearing 31 and 63 days after application are flat and wide, including the largest irrigation system (Fig. 7B). The curves show that the dispersive transport is enhanced by irrigation. The analysis of the results shows some migration of carbofuran down into in soil. The downward movement of carbofuran is more pronounced in Fig. 7B, which coincides with the increase in irrigation, i.e., 350 mm in total. The carbofuran mobility is also facilitated by its high solubility in water, i.e., 700 mg/L, and the large amounts of water provided by the irrigation system.