System dynamics simulation of crop yield under different irrigation water quality and quantity

According to the increasing population, the requirements for crop production grow bigger every day. And the availability of freshwater is a major limiting factor in sustainable agriculture (Ladeiro 2012).

The main reasons for crop loss during the growing season are drought and salinity stress. Therefore, the problem could be largely solved by applying appropriate management techniques (Houshmand et al. 2005; Sharma & Minhas 2005). In order to study and analyze various types of irrigation management, experiments can be conducted in the field under natural settings. But the main limitations in this method are that they are time consuming and have high costs. In addition, maintaining constant test conditions during the test period is very difficult. So, many attempts have been made to develop computer models. As an example, the DRAINMOD software was developed to simulate the performance of drainage and related water management systems in regions with shallow water tables (Skaggs 1980). Although this model is capable of simulating all aspects, it still faces some limitations such as the high numbers of inputs (Kroes & Van Dam 2008).

Also, Honar et al. (2012) used the CRPSM model as a dynamic and mechanistic model for simulating the yield of canola plants. The results showed that CRPSM is a reliable model in predicting the yield of canola plants. Azamathulla et al. (2008) used Genetic Algorithm and Linear Programming for known total irrigation demand and the optimal allocation of water available to crops at the farm scale. The results showed that the GA model gives better performance than the LP model.

According to studies, there are many numerical models that can be used to simulate solute transport in porous media and evaluate product performance. But due to the sensitivity of the models to boundary conditions, their usage in field crop conditions is controversial (Chops & Hopmans 2002; Jorenush & Sepaskhah 2003). The contribution of irrigation water quality and quantity to crop yield and over time cannot be determined in most of the existing analytical models. Therefore, there is a need to explore new tools to represent the complex relationships found in systems.

One of the most effective methods to evaluate complex systems is system dynamics, which is one of the most powerful and visual methods for simulation (Nozari et al. 2014). This method was first developed by Forrester (1961) to understand the structure and dynamics of complex systems. It applies to dynamic problems arising in complex social, managerial, economic or ecological systems - literally, any dynamic systems characterized by interdependence, mutual interaction, information feedback, and circular causality (Richardson 1991). Fletcher (1998) used this method as a decision analysis method in managing dehydration.

Pourfallah Koushali et al. (2015) used the system dynamics approach for water resources modeling of the Sefīd-Rūd Watershed in the north of Iran. In this study, the river catchment basin, dam reservoir and the irrigation system were studied under three different scenarios.

Matinzadeh et al. (2017) used the system dynamics approach to develop a simulation model of nitrogen dynamics in a shallow water table. The results showed that this approach is appropriate for the nonlinear behavior understanding of complex systems over time and could be used for modeling fertilizer management to prevent environmental pollution. Saysel & Barlas (2001) examined salinity effect on irrigated lands in southeastern Turkey by the system dynamics method. These researchers showed that the system dynamics approach provides a good way to manage environmental problems.

Khan (2004) developed a model to simulate utilization of water resources in irrigation networks in the Qiantang river basin of China. In this study, the influences of effective factors on the management of water resources in irrigation systems, including changes in irrigation efficiency, crop pattern change and land level change for each product were studied by using the system dynamics approach. He found qualitative and quantitative agreements.

Nozari & Liaghat (2014) used the system dynamics approach to simulate the performance of a drainage system in the unsteady state condition. The model test was conducted using data collected at an experimental site located in the Khuzestan province of Iran. The agreement between the observed and simulated variables was very good. They found the model can be used to manage and control saline drainage water and also to prevent environmental damage.

Gonçalves & Giorgetti (2013) simulated water quality of river in steady and unsteady state conditions using the VENSIM software. The researchers simulated several qualitative parameters such as DO, BOD, pH, nitrate, nitrite, temperature, and so on. They stated that using VENSIM as a fundamental tool is a very good solution in developing environmental models. Sadeghi Khalegh Abadi et al. (2014) used the system dynamics approach to develop a comprehensive dynamic model for the Karkheh Basin, which included the population, economic activities, water demand and water resources. Furthermore, Ahmadi et al. (2014) showed that the system dynamics approach is a good alternative to conventional simulation systems, and has some advantages such as increasing model development speed, modifying its structure, the ability for sensitivity analysis and the effective relationship between components.

According to studies, it is observed that many researches have been done on development of computer models and a number of analytical and numerical models have been proposed to calculate agricultural products' performance. The use of these models for different field conditions is controversial, due to the large number of input parameters, the long model run time and the sensitivity of the numerical models to boundary conditions. Most of the existing analytical models are not able to evaluate the influence of quantity and quality of irrigation water on root zone soil salinity over time and on agricultural production. So, considering the importance of using the system dynamics approach in problem solving (such as high speed, user's access to program lines, easily removing or adding a parameter and investigation of its effect on a whole system) we decided to use the system dynamics approach for simulation of water and soil systems in the presence of a subsurface drainage system. For this purpose, water and soil system and crop yield were simulated by the system dynamics model, and the effect of quality and quantity of irrigated water on crop yield, soil moisture and salinity was evaluated over time.

In order to estimate the amount of storage of water and salt in the soil, mass balance and volume balance equations were used, respectively. Therefore, we divided the soil into three layers from the surface down to the underlying bedrock. The first layer is the active root zone. The second layer consists from the bottom of the root zone to subsurface drain depth. The third layer continued from subsurface drain depth to bedrock. Each layer was defined as a state variable. Rainfall, irrigation, lateral flow and upward flow add water to the root zone, and evapotranspiration, drainage lateral leakage and percolation losses remove water from the root zone.

System dynamics is a computer method to understand nonlinear trends in complex systems over time. It is based on feedback processes that are influenced by their past behavior. The design of the causal loop diagram is one of the basic processes of this method. It shows the cause and effect relationships between the variables of a system. Each causal link is assigned a polarity, either positive (+) or negative (−) to indicate the direction of causal relationship. A positive causal link means that the two nodes change in the same direction: if the node in which the link starts decreasing, the other node also decreases. Similarly, if the node in which the link starts increasing, the other node increases as well. A negative causal link means the two nodes change in opposite directions, if the node in which the link starts increasing, the other node decreases and vice versa.

In addition, a systems dynamics model consists of stocks and flows. A stock is an accumulation of material that has built up in a system over time. It indicates system status and decisions so the activities of the system are based on it. A flow is a material that enters or leaves a stock over a period of time.

In this study, the system dynamics tool VENSIM is used to conceptualize, simulate and analyze the system. It allows stock-and-flow diagrams to be drawn before entering equations.

The causal loop structure of crop irrigation system is shown in Figure 1. In loop 1, soil water potential is depleted by actual evapotranspiration. As soil water potential increases, the soil moisture and then the soil moisture factor increase, which in turn increases actual evapotranspiration. Thus soil water potential and soil moisture is decreased. Loop one is negative and represents actual evapotranspiration. Loop two is a positive feedback loop. The larger the actual evapotranspiration, the larger the upward flow of groundwater level, then the larger the soil water potential, the soil moisture and the soil moisture factor, which in turn increases actual evapotranspiration, completing the positive loop. The third feedback loop represents the interaction between soil water potential and deep percolation. The larger the soil water potential, the larger the soil moisture, then lack of moisture and percolation will increase, which in turn decreases soil water potential, completing the negative loop. On the other hand, with increasing deep percolation, the water table rises; as it does, the upward flow rate increases. As a result, soil water in the unsaturated zone increases (loop 4). In loop 5, which also has a rising water table, the upward flow rate increases and the water table reduces, forming a negative loop. In loop 6, which has increased soil salinity, the leaching requirement increases; and for root zone leaching, the irrigation water required increases. Because of the increase in irrigation water, soil salinity in the root zone decreases, so this loop is defined as negative. Finally loop 7, increased groundwater salinization increases the upward flow of salinity, increasing the soil salinity.

The state and flow diagram of soil salinity is shown in Figure 2. Soil salinity as a state variable and rate of salt changes with two flow variables. Soil salinity is increased by evapotranspiration, and is depleted by leaching and infiltration.

Figure 3 also shows the state and flow diagram of the crop yield presented as an example.

Where T is the time, i and j are the time steps, which may not be the same for λ and β periods.

The model was validated using data collected from the Right Abshar Irrigation Network located at the hydrological Zayandehrud River, Iran with a longitude of 51°50′ to 52°09′ and latitude of 32°42′ to 32°67′ and an average elevation of 1,500 m above mean sea level. The monthly temperature ranges from 3 to 30 °C and the area has an arid and semiarid climate with average annual rainfall of 120 mm. Nine treatments of sugar beet cultivated as a randomized complete block statistical design split into two plots (split-split plot) were used to investigate the crop yield. Sugar beet crop was cultivated on 30/04/2007 and harvested on 20/11/2007.

Different treatments have different conditions in terms of irrigation water quantity, irrigation water salinity, soil saturation extract salinity and soil moisture. Tables 1 and 2 show the quantity and quality of irrigation water for each of treatments, respectively (Rezaverdinejad 2010).

Soil characteristics, soil moisture and some crop parameters that were used in the modeling are shown in Table 3.

Where g is the parameter indicating soil capacity to upward flow, P is the discharge coefficient of soil water for evapotranspiration 5 mm per day, θ_Fc is the water content of field capacity, θ_wp is the water content of the wilting point, θ_s is the saturated water content, D_r is the developed deep root, and b is the percent yield loss for each dS·m⁻¹ when salinity is increased.

As mentioned, the model was simulated for the sugar beet for the different treatments of irrigation given in Tables 1 and 2. The crop yield values in different treatments and under different irrigation conditions and soil salinity are presented in Table 4. Table 4 shows that crop yield decreases with increasing soil salinity. Also, in high salinity irrigation water, leaching increased the performance of sugar beet and with the increasing soil salinity, the effect of leaching was so much in increasing sugar beet performance. Yields simulated by the model and measured values are presented in Figure 4.

The value of RMSE and SE indices for performance of this product are 2,777 kg·ha⁻¹ and 0.07, respectively. So, there is a good agreement between the measured and simulated yields. Rezaverdinejad (2010) has assessed the SWAP performance on simulated sugar beet yield. The RMSE to 5,370 kg·ha⁻¹ for sugar beet plants was evaluated to be good.

The F-test was used to determine whether the two variances are equal. The results showed that for α = 0.05, the P-value was 0.470. So, the variance values of the simulated data are not significantly different from those of the measured data. Therefore, we used two-tail t-test with the same variances to demonstrate the difference between the two series. The P-value was 0.474, which indicated that there is no significant difference between simulated and measured data series at 5% significance level. So, it can be concluded that the present model has good accuracy in simulating sugar beet yield. Tables 5 and 6 show F-test for variances and t-test for two-sample assuming equal variances, respectively.

One of the important factors in plant yield is the available moisture in the soil. The comparisons of the simulated and measured soil moistures during the growing season for treatments is shown in Figure 5. Olivera-Guerra (2019) has stated that in the simulation of irrigation based on the modeling of soil moisture dynamics from the water balance, significant uncertainties can be obtained.

The values of RMSE and SE of soil moisture for all treatments are shown in Table 7. According to Table 7 the mean values of RMSE and SE for all treatments are equal to 0.026 and 0.09 respectively. Zand-parsa et al. (2016) examined the soil moisture by the WSM model, and obtained the SE value of soil moisture in different layers of root zone was 0.14, and concluded that the result of this model can be considered acceptable. So, in this study there is a good agreement between the simulated and measured soil moisture.

Also F-test was used for variance analysis of simulated and measured data series, and the P-value was equal to 0.068. According to the P-value, the variance values of the simulated data are not significantly different from those of the measured data at the 5% significance level. So, we used the two-tail t-test with the same variances to understand the difference between two series. The P-value was 0.794, which indicated that there are no significant differences between simulated and measured data series of soil moisture in the root zone at the 5% significance level. It can be concluded that the model has acceptable accuracy in simulating soil moisture. The results of F-test and t-test are shown in Tables 8 and 9 respectively.

The other reference parameter in this study is the amount of soil saturation extract salinity in root zone. The simulated soil salinity versus the measured soil salinity for different treatments is shown in Figure 6 for a large spectrum of salinity levels from 3.5 to 14 dS·m⁻¹. The results show that all simulated data have a standard deviation of less than 8 percent, without an apparent trend for over-prediction or under-prediction. The values of RMSE and SE for all treatments are presented in Table 10 which shows an acceptable agreement between simulated and measured results. According to Table 10, the mean values of RMSE and SE for all treatments are equal to 0.54 dS·m⁻¹ and 0.08, respectively. Thus, according to the research results with the AquaCrop model, if the SE is no more than 0.12 for soil salinity, the results can be considered acceptable (Heydari-nia et al. 2017).

F-test was used for variance analysis of two data series of soil salinity, and the P-value was obtained as 0.470. The results showed that the two standard deviation are equal at significance level of 5% (Table 11). Therefore, the two-tail t-test was used with the same variances to delineate the difference between simulated and measured data series. The P-value was 0.946, which indicated that there are no significant differences between two data series at the 5% significance level. So, the present model has good accuracy in simulating soil salinity in root zone. The results of the F-test for variances and t-test for two-sample assuming equal variances are shown in Tables 11 and 12 respectively.

In this paper, a model of crop yield simulation for saline water irrigation has been developed by using a system dynamics approach. Simulated results of the variations in amount of irrigation water and its salinity on the performance of sugar beet crop yield were consistent with field experimental results, indicating that the system dynamics approach could be applied to simulate crop response to water stress and saline stress. The mean values of RMSE and SE for all crop yield treatments were 2,777 kg·ha⁻¹ and 0.07, respectively, which showed there is good agreement between simulated and actual crop yield. Also simulation of crops' performance under different managements systems of irrigation using the system dynamics approach was examined by Nozari et al. (2014) and the results showed system dynamics was found to be an efficient approach for simulation of the water and soil system and crop yield.

The mean values of RMSE and SE for soil moisture were 0.026 Dmnl and 0.09 and those for soil salinity were 0.54 dS·m⁻¹ and 0.08, respectively. The results show that the model could be used for saline irrigation water management.

However, important merits of the model are the high accuracy, the increased speed of model development, the ability to simply modify model structure in response to changes in the system, the ability to simulate the interactions between the soil and water components and the ability to undertake sensitivity analysis.

None.

All relevant data are included in the paper or its Supplementary Information.