Comparison of Approaches for Irrigation Scheduling Using AquaCrop and NSGA-III Models under Climate Uncertainty

In the face of increased competition for water resources, optimal irrigation scheduling is necessary for sustainable development of irrigated agriculture. However, optimal irrigation scheduling is a nonlinear problem with many competing and conflicting objectives and constraints, and deals with an environment in which conditions are uncertain. In this study, a multi-objective optimization problem for irrigation scheduling was presented in which maximization of net benefits and water use efficiency and minimization of risk were the objectives. The presented optimization problem was solved using four different approaches, all of which used the AquaCrop model and nondominated sorting genetic algorithm III. Approach 1 used dynamic climate data without adaption; Approach 2 used dynamic climate data with adaption; Approach 3 used static climate data without adaption; and Approach 4 used static climate data with adaption. The dynamic climate data were generated using the bootstrap resampling of original climate data. A case study of maize production in north Jiangsu Province of China was used to evaluate the proposed approaches. Under the multi-objective scenario presented and other conditions of the study, Approach 4 gave the best results, and showed that irrigation depths of 400, 325, and 200 mm were required to produce a maize crop in a dry, normal, and wet year, respectively.


Introduction
In the face of increased competition for water resources, optimal irrigation scheduling is necessary for the sustainable development of irrigated agriculture. Irrigated agriculture faces challenges of reduced water and land resources because of increased urbanization, industrial development, the need to protect the environment, requirements for a better quality of life, and many others. Therefore, agriculture must utilize water resources more efficiently while meeting other necessary objectives. Irrigation scheduling is a key element in on-farm irrigation water management and is therefore an important consideration in this endeavor. However, optimal irrigation scheduling is not a trivial problem because it is a nonlinear problem with many competing and conflicting objectives and constraints. It is further influenced by many factors in the soil-plant-atmosphere system, and the best solution for any particular case has to be found from a very large number of other possible solutions, thus requiring a systematic means for locating it [1].
Field experiments conducted to establish crop yield response to irrigation amount and timing require much time and resources. Simple empirical models such as that of Doorenbos and Kassam [2] Sustainability 2020, 12, 7694 2 of 20 have been applied to estimate crop yield response to irrigation water but these models are often only applicable to sites where they were derived. Dynamic crop simulation models which are based on current understanding of crop physiological responses to environmental variables can simulate crop development and yield more realistically than simple empirical models. After calibration, these models are a more cost-effective means for evaluating crop yield response to irrigation under various field conditions than field experiments. While field experiments are important and cannot be done away with, combining field experiments with dynamic crop modeling cannot only reduce the cost of field experiments by narrowing down viable strategies to test, but can also enrich the field experiment studies by introducing alternatives not previously considered for testing. However, many crop simulation models require a large number of input parameters, many of which are difficult to measure in the field, thereby limiting application of the models to researchers who are more familiar with the required input parameters and leaving out farmers and other stakeholders who are usually not familiar with them. The AquaCrop model [3][4][5][6] was developed with a view to balance simplicity, robustness, and accuracy and is a more suitable model for farmers to use than the more complex models available.
Evolutionary optimization algorithms offer many advantages over traditional optimization techniques. Their capability to link with simulation models makes them particularly attractive for many applications. Genetic algorithms, a type of evolutionary optimization algorithm, has been applied in many optimal irrigation scheduling studies. Recent studies include Jamal et al. [7], who presented a real-time decision support modeling framework for irrigation scheduling using probabilistic seasonal weather forecasts, a crop simulation model, and genetic algorithms. Li et al. [8] presented a simulation-optimization model for irrigation scheduling of spring wheat under uncertainty using the AquaCrop model, genetic algorithms, and a bootstrap technique for resampling climate data. Wen et al. [9] used a crop water production function, a soil water balance model, and genetic algorithms to optimize irrigation scheduling of spring wheat in northwest China. Lalehzari et al. [10] combined the nondominated sorting genetic algorithm II (NSGA-II), particle swarm optimization algorithm, and the Doorenbos and Kassam [2] crop production function to evaluate optimal cropping patterns, allocation of surface and groundwater resources for irrigation, and irrigation scheduling. They expressed the problem they were solving as a constrained biobjective optimization problem.
The value in combining an optimization algorithm with a crop simulation model to solve irrigation scheduling optimization problems was demonstrated by Ioslovich et al. [11] who showed that an optimization model consisting of a heuristic optimization algorithm and AquaCrop model generated an irrigation schedule that gave a higher crop yield than one which was generated by the AquaCrop model alone even though the irrigation amounts were roughly the same. Other recent studies on optimal irrigation scheduling include Rafiee and Shourian [12], who presented a simulation-optimization approach for combining the Soil and Water Assessment Tool (SWAT) with the harmony search algorithm to find optimal irrigation and cropping patterns. Linker [13] also presented a model-based optimization procedure for optimal allocation of cropping areas and water by utilizing the AquaCrop model and a constrained particle swarm optimization solver. Their procedure first determined optimal land and water allocations, and then optimal crop yield and irrigation scheduling. In that study, the AquaCrop model was not used within the main optimization procedure but only to generate water productivity functions which were the ones used in the main optimization procedure. Linker and Sylaios [14] presented an efficient model-based procedure for generating near-optimal irrigation schedules for real-time applications utilizing the AquaCrop model and the epsilon-moga evolutionary algorithm. Linker et al. [15] presented an optimization scheme for generating sub-optimal irrigation schedules for three crops using the AquaCrop model and the OptQuest Nonlinear Program (OQNLP) solver in combination with the qlcAssign procedure for global nonlinear search. Xu et al. [16] presented a procedure for optimal allocation of irrigation water resources to multiple crops and multiple agricultural subareas using a fuzzy random simulation-based particle swarm optimization and an empirical crop production function. In their study, the decision variables were the area and irrigation water volume allocated to each crop in the season, and the state variable was the available water Sustainability 2020, 12, 7694 3 of 20 stored in the reservoir at any period of the season. Seasonal reservoir inflow and rainfall were treated as fuzzy random variables. Shan et al. [17] used a water-yield-quality-benefit relationship to establish an irrigation scheduling optimization framework. Their model was designed to create a balance between yield quantity and quality, and to find optimal irrigation schedules for various market conditions. Pereira et al. [18] developed a multi-objective model to solve an optimal irrigation water resources allocation problem, considering economic, social, and environmental benefits as objectives. Nguyen et al. [1] presented a simulation-optimization framework for optimal irrigation and fertilizer scheduling. García-Vila and Fereres [19] presented a model to assist farmers' plan cropping patterns and irrigation strategies using the AquaCrop model and an economic model under various scenarios of irrigation water constraints, agricultural policies, changes in crop product and water prices, and communication delays to farmers of seasonal irrigation water allocations. Sadati et al. [20] presented a nonlinear programming optimization model with an integrated soil water balance to determine optimal reservoir release policies and cropping patterns. They solved the model using a genetic algorithm. Li et al. [21] presented an inexact fuzzy stochastic simulation-optimization programming model for the optimal irrigation water allocation of two water sources. Their model combined a crop water model and a field water cycle model with an uncertainty optimization model. Zhai et al. [22] used the AquaCrop crop model and entropy-cloud model to develop an irrigation scheduling optimization model for rice.
Most of the studies cited above treated irrigation scheduling optimization as a single objective optimization problem, with maximization of net benefits as the most common objective. However, in situations when rainfall is insufficient for crop production, achieving a high crop production per unit of water used might be a policy that is more likely to be followed than one that seeks to merely maximize crop yield or benefit [23]. In this case, optimizing crop production per unit of water used will involve both the production and water use simultaneously, requiring optimization of both the amount of water used and the benefit derived. Therefore, some studies (e.g., [10,14,24]) have considered maximization of net benefits or yield, and maximization of relative water use efficiency or minimization of total irrigation water used or minimization of the total number of irrigations, simultaneously, to achieve this aim. However, some risk-averse farmers may not prefer an irrigation schedule that would give them the highest expected net benefits and water use efficiency if the risk involved is very high. Such farmers may be happy to choose an irrigation schedule with low expected net benefits if the risk involved is also lower than one with a high expected net benefit but at high risk. Therefore, in such cases it would be useful to have risk minimization as an additional objective in the optimization problem. Therefore, in this study, irrigation scheduling optimization was treated as a multi-optimization problem with maximization of net benefits and water use efficiency and minimization of risk as the objectives.
The aims of the current study were to (1) solve a multi-objective optimization problem for irrigation scheduling of maize in north Jiangsu Province of China using four different approaches, (2) determine the best performing approach, and (3) use the best approach to determine optimal irrigation schedules for a dry, normal, and wet year. All the approaches used the AquaCrop model and nondominated sorting genetic algorithm III (NSGA-III).

The AquaCrop Model
The AquaCrop model [4][5][6] simulates water-driven plant growth and yield. It is suitable for evaluating the effect of irrigation schedules on crop yield. The model calculates the soil water balance, considering rainfall, irrigation, capillary rise, runoff, evaporation, transpiration, and deep percolation. To simulate plant growth and yield requires climate, soil, crop, and field management characteristics to be specified in the model. In the model, crop yield is determined from biomass produced via a harvest index. Biomass produced is related to plant transpiration via water productivity normalized for climate and atmospheric CO 2 concentration as: where bio k is biomass produced on day k; w is the normalized water productivity; T r k is plant transpiration on day k; and ET o k is reference crop evapotranspiration on day k. The model calculates canopy cover instead of leaf area index. The canopy cover is defined as the fraction of ground covered by the crop canopy when observed from above the crop. During the first part of the growth season, canopy cover of the crop is modeled as: where c is the canopy cover at period k; c o is canopy cover at 90% seedling emergence; and GC is the canopy growth coefficient (expressed as fraction per growing degree day (GDD) or per day under optimal conditions but adjusted for stresses). During the second half of the plant growth, c is expressed as: where c x is the maximum canopy cover under optimal conditions. Deterioration in green canopy cover as a result of leaf senescence is expressed as: where DC is the canopy decline coefficient (in fraction per GDD or per day). Transpiration is determined from: where c is green canopy cover adjusted for microadvective effects and Kcb x is the crop coefficient when c = 1. Crop yield is expressed as: where Y is crop yield; HI is harvest index; and bio is biomass. More details on the AquaCrop model calculation procedures can be found in the literature (e.g., Hsiao et al. [4], Raes et al. [5], Steduto et al. [6], Vanuytrecht et al. [25]). The AquaCrop plug-in model version 6.0 was used in the current study to evaluate crop yield. It has identical calculation procedures to the standard AquaCrop version program but does not have a user interface and is more ideal for applications requiring iterative simulations as in the current study.

The Nondominated Sorting Genetic Algorithm III
The NSGA-III [20,21] is a multi-objective optimization algorithm based on the NSGA-II [24] framework. Multi-objective optimization problems are a feature of many real-life problems and optimal irrigation scheduling is not an exception. Where the NSGA-II used crowding distance to maintain diversity among population members, the NSGA-III uses reference points.
Reference points are placed on a normalized hyper-plane which is equally inclined to all objective axes, with an intercept of one on each axis. To move from one generation to the next, nondominated solutions are emphasized as well as solutions associated with the reference points. A wide distribution of the reference points on the hyper-plane creates diversity in the population members.
The total number of reference points is calculated as: Sustainability 2020, 12, 7694

of 20
where r is the number of reference points; G is the number of objectives; D iv is the number of divisions along each objective; and () represents the binomial coefficient. The minimum value for each objective function ever found in the population is identified to constitute the ideal point. The ideal point is then subtracted from the objective values in the population set. The extreme point for each objective axis is then identified to constitute the G-dimensional hyper-plane. The objective values are normalized as: where where f w is the objective value in P t ; v min w is the minimum value for each objective function ever found; P t is the parent population at the tth generation; S t is the new population; f n w is the normalized wth objective; and c w is the intercept of the wth objective axis and the linear hyper-plane.
The normalization procedure and creation of the hyper-plane is done at each generation using all extreme points ever found from the start of the simulation, thereby adaptively maintaining a diversity in the members of S t at every generation. After normalizing each objective based on the extent of members in S t , each population member is associated with a reference point. A reference line corresponding to each reference point on the hyper-plane is defined by joining the reference point with the origin. The perpendicular distance of each population member of S t from each of the reference lines is calculated and the reference point whose reference line is closest to a population member in the normalized objective space is taken to be associated with the population member. A new population P t+1 is formed by elitist selection of solutions and by emphasizing solutions closest to the reference line of each reference point. An offspring population is then created by the usual genetic operators of crossover and mutation, and randomly selecting parents from P t+1 . The NSGA-III algorithm also requires the usual genetic algorithm parameters such as population size, termination criteria, and crossover and mutation probabilities. More details on the NSGA-III can be found in Deb and Jain [26] and Jain and Deb [27].

Ranking Method for Selecting Solutions from Pareto Optimal Front
A drawback of the Pareto optimality concept used in multi-objective optimization algorithms is that Pareto optimal alternatives usually make making a selection of one solution difficult. Gogodze [25] proposed an approach for selecting one solution from the Pareto optimal front. The approach was based on ranking-theory methods used to rank participants in competitive sports tournaments. This approach was adopted in the current study. The approach involved making a count of the number of times a particular solution was better than the other solutions in all objectives and selecting the solution with the largest count. The set of all solutions was a union of all Pareto fronts ever found, i.e.,: where A is the set of all solutions (A = [a 1 , . . . , a M ]); F 1 is a set of Pareto front solutions in the tth generation; and t m is the last generation. The total count for solution a i is given by: where c i is the total count for solution a i ; c ij is the count of times a i is superior to other solutions regarding objective j; G is the number of objectives; and M is the total number of solutions.

Bootstrap Resampling
Theoretical probability distribution functions were found to not fit well to reference crop transpiration and rainfall data in this study. This result was also observed by Prats and Picó [28]. Therefore, a bootstrapping technique was used to generate climate sequences for the model approaches 1 and 2 in order to simulate climatic uncertainty. A 21-year (1999-2019) record of daily historical climate data of the study area was used in the current study. The climate data consisted of maximum and minimum temperature ( • C), rainfall (mm), average air pressure (kPa), average wind speed (m/s), average relative humidity (%), and sunshine hours (h). Reference evapotranspiration (ETo) data were calculated according to Allen et al. [29].
The climatic data were divided into wet, normal, and dry years, according to the annual rainfall distribution. A wet year was considered to be one whose annual rainfall was above the third quartile (i.e., wet year > third quartile); a normal year was one whose annual rainfall was between the first quartile and the third quartile (i.e., first quartile < normal year < third quartile); and a dry year was one whose annual rainfall was below the first quartile (i.e., dry year < first quartile). The probabilities of the occurrence of a dry year (P 1 ), normal year (P 2 ), and wet year (P 3 ) were estimated using the historical climate data. The crop season was divided into 10-day periods. Climate data corresponding to each 10-day period from all years in with a particular climate (i.e., dry, normal, or wet year) were compiled in a set so that each 10-day period in a season with a particular climate had a set of possible climatic values for that climate. Within each set, the values were randomly shuffled. To generate dynamic climate data for Approaches 1 and 2, the following steps were followed: Step 1: A random value z was generated from the standard uniform distribution in the open interval [0,1].
, wet year climate was selected.
Step 3: 10 random whole numbers between 1 and total number of values in each 10-day set for the climate selected in Step 2 were generated, inclusively, and values whose positions in the sets corresponded with these numbers were selected.
Step 4: The 10-day weather sets generated in Step 3 were joined end to end to create seasonal weather.
Step 5: Steps 1 to 4 were repeated 20 times to create 20 seasons.

Optimal Irrigation Scheduling Approaches
A flow chart of the proposed approaches for solving the optimal irrigation scheduling problem are shown in Figure 1. Following are the steps followed in each of the four different approaches, which are illustrated in Figure 1.

Approach 1
Step 1: Maximum number of iterations was selected. This was the number of generations through which the solutions would evolve to become optimal solutions by the genetic operators of selection, crossover, and mutation.

Step 2:
The population size was selected. This was the number of the solutions that would be evolving to become optimal solutions. Step 3: The reference points were created. These were used for selecting solutions when moving from one iteration to the next.

Step 4:
Percentage of crossover was selected. This was the percentage of the population that would undergo crossover to produce new population members or offspring.

Step 5:
Percentage of mutation was selected. This was the percentage of the population that would undergo mutation to produce new population members.

Step 6:
Complete climate files for 20 crop seasons were created by the bootstrap resampling of the original climate data and input into the AquaCrop model.

Step 7:
A population of solutions of the size from Step 2 was randomly generated. Each population member was a vector of real numbers representing maximum allowable soil water depletion levels for each 10-day period in the crop season.

Step 8:
A population member from Step 7 was input into the AquaCrop model and the AquaCrop model was run to evaluate crop yield and total water use for each of the 20 crop seasons from Step 6.
Step 9: Results from Step 8 were retrieved from the AquaCrop model and used to evaluate the objective function values corresponding to the population member (or solution) selected in Step 8.
Step 10: Steps 8 and 9 were repeated for all the population members from Step 7.
Step 11: Some population members from Step 7 were randomly selected for crossover.
Step 12: Steps 8 and 9 were repeated for offspring from Step 11.
Step 13: Some population members from Step 7 were selected for mutation.

Step 14:
Steps 8 and 9 were repeated for new members from Step 13.

Step 15:
Population members were sorted according to nondomination and selected for the next iteration.

Step 16:
In the next iteration, Steps 11 to 14 were repeated but with the population from Step 15.

Step 18:
Steps 11 to 17 were repeated until the total number of iterations from Step 1 was completed.
Step 19: The ranking method was used to find the best solution.
Step 6: Current 10-day period of the season for which irrigation was being determined was selected. For example, the first period considered at the beginning was the first 10-day period of the season. Afterwards, the second 10-day period was considered. The periods were considered one after the other and the next period was not considered unless the current period's solution was determined.
Step 21: The maximum allowable depletion level for the 10-day period from Step 6 was saved to a file.

Step 22:
The actual season climate for the 10-day period from Step 6 was used to update climate files in Step 7.
Step 23: Step 8 was repeated but with the allowable depletion level for the period chosen in Step 6 as the given value.
Step 24: Step 6 was repeated by selecting the next 10-day period.
Step 25: Steps 9 to 24 were repeated until the total number of 10-day periods in the season was covered.

Steps 6:
Climate files were created using original climate data and were input into the AquaCrop model. Steps 7-16: Same as Steps 7-16 for Approach 1, respectively.

Steps 17:
Steps 11 to 16 were repeated until the total number of iterations from Step 1 was completed.

Steps 18:
The ranking method was used to find the best solution.

Step 6:
Current 10-day period of the season for which irrigation was being determined was selected.

Step 7:
Climate files were created using original climate data and were input into the AquaCrop model. Steps 8-17: Same as Steps 7-16 in Approach 1, respectively.
Step 18: Steps 12 to 17 were repeated until the total number of iterations from Step 1 was completed. Flow chart of the four proposed approaches for optimal irrigation scheduling. TMP is temperature file, ETo is reference evapotranspiration file, PLU is rainfall file, CLI is climate file, NB is net benefit, R is risk, WUE is water use efficiency, D is maximum allowable soil water depletion level, NSGA-III is nondominated sorting genetic algorithm III, N is no, and Y is yes.

Approach 1
Step 1: Maximum number of iterations was selected. This was the number of generations through which the solutions would evolve to become optimal solutions by the genetic operators of selection, crossover, and mutation.

Figure 1.
Flow chart of the four proposed approaches for optimal irrigation scheduling. TMP is temperature file, ETo is reference evapotranspiration file, PLU is rainfall file, CLI is climate file, NB is net benefit, R is risk, WUE is water use efficiency, D is maximum allowable soil water depletion level, NSGA-III is nondominated sorting genetic algorithm III, N is no, and Y is yes.

Irrigation Scheduling Optimization Problem
The multi-objective irrigation scheduling optimization problem considered in the paper was expressed as: Find: such that: subject to: where [D 1 , D 2 , . . . , D n ] are the maximum allowable depletion levels (% AW) in the root zone; n is the number of 10-day periods in the crop season; AW is available water, defined as the difference between field capacity and wilting point in the root zone (mm); NB is net benefit (RMB/ha), where RMB is the currency of the People's Republic of China; WUE is the water use efficiency (kg/m 3 ); R is the risk (t/ha); T is the total number of irrigation days in a season; I k is irrigation depth on day k in the season (mm); I m is the allocated seasonal irrigation quota (mm); θ k is the soil water content in the root zone (m 3 /m 3 ); θ D k k is the soil water content in the root zone corresponding to allowable depletion level D k (m 3 /m 3 ); I c is irrigation system capacity (mm/day); d l , d m are days during the irrigation season on which irrigation is applied; l, m are indices denoting different days in the season; δ is the allowable irrigation interval (days); D L and D U refer to lower and upper bounds set for maximum allowable soil water depletion (%); and D i is the maximum allowable depletion level for any 10-day period of the season.
Approaches 1 and 3 only use the forecasted climate in the optimization model. Therefore, the following constraint on climate was added for Approaches 1 and 3: where ρ k is the climate used in the optimization model for day k and ρ f k is the forecasted climate. For Approaches 2 and 4, the additional constraints added were: where k * is current day in the season; I * k is irrigation already applied in the season before day k * ; and ρ * k is observed climate. For Approaches 1 and 2, NB and WUE were expressed as: where Y is the average crop yield (ton/ha); S P is the selling price of the crop produce (RMB/ton); m 0 is the fixed cost of production (RMB/ha); m 1 is the cost of irrigation water (RMB/m 3 ); m 2 is the fixed cost of pumping (RMB/mm·ha); E f is the efficiency of the irrigation system (%); and V k is the rainfall on day k (mm). For Approaches 3 and 4, NB and WUE were expressed as: where subscripts i, j, and k are for climate (dry year, normal year, and wet year), season, and day; N c is the number of climate years (i.e., wet, normal, or dry year); N si is the number of crop seasons in the climate i; P i is the probability of season j to be of climate i according to the climate data; Y ij is crop yield from crop season j in climate i; and V k is precipitation (mm). For all the approaches, risk was expressed as [19]: where R is risk (RMB) and NB U , NB L are net benefits corresponding to years with satisfactory and unsatisfactory weather, respectively. Irrigation schedules for dry, normal, and wet years were determined according to the following steps: Step 1: Optimal D was determined for each year having a particular climate.
Step 2: The optimal D obtained was used to evaluate the corresponding optimal irrigation schedule in the AquaCrop model for that year.
Step 3: The determined irrigation schedules were applied to other years with the same climate.
Step 4: The irrigation schedule which resulted in the best average objective values was selected as the irrigation schedule for the climate.

Case Study
The case study was of north Jiangsu Province of China (34 • 17 N, 117 • 09 E, elevation 41 m). The soil was sandy loam, with mean soil dry bulk density of 1.5 gcm −3 , mean saturated water content of 0.39 cm 3 cm −3 , mean field capacity (FC) of 0.29 cm 3 cm −3 , mean permanent wilting point of 0.12 cm 3 cm −3 , and mean saturated hydraulic conductivity of 560 mmd −1 for the 0-120 cm soil layer. The applicability of the AquaCrop model to simulate water consumption of summer maize in China has been evaluated by Liu and Shen [30] who conducted field experiments at Luancheng Agricultural Ecosystem Experimental Station of the Chinese Academy of Sciences and concluded that the model could be used to evaluate the water consumption and water use efficiency of summer maize in the North China Plain. Ran et al. [31] and Ran et al. [32] also showed that with minor adjustments to the model parameters, the AquaCrop model could simulate the growth and yield of seed maize in northwestern China with acceptable accuracy. Parameter values established by Liu and Shen [30] for summer maize were used to calibrate the AquaCrop model used in the current study. The crop season considered was from May to mid-August and simulations were performed for this period only. The planting date was 1 May. The length of the irrigation season was 100 days. The decision variable was a vector of 10 real numbers representing maximum allowable soil water depletion levels as percentages of AW. Irrigation was to be applied when soil water depletion levels reached these values. The irrigation application depth was 25 mm per irrigation during the first 30 days of the season, and 50 mm per irrigation afterwards until the end of the season. A constant irrigation application depth instead of a variable one presents some management advantages for the irrigator. In genetic algorithms, a population size of~3 times the number of variables is recommended, therefore, in this study, the population size was 30. The number of iterations was 200. No attempt was made to optimize the probabilities of crossover and mutation, and they were both given a value of 0.3. A minimum irrigation interval (δ) of 5 days was adopted. The maximum allowable seasonal irrigation depth (I m ) was taken as 450 mm. The lower (D L ) and upper limits (D U ) of maximum allowable soil water depletion levels were taken as 10% and 40% of AW, respectively. The initial soil water content was assumed to be at 80% of field capacity. The price of maize, fixed costs of production, fixed cost of irrigation water, variable cost of irrigation water, and efficiency of the irrigation system were assumed to be RMB33000/ton, RMB3750/ha, RMB30/ha, RMB0.157/m 3 , and 75%, respectively, and were assumed constant. In the climate data, 5 years were dry years, 11 years were normal years, and 5 years were wet years. The simulation-optimization model was developed and run on the MATLAB platform (MATLAB 2016a).

Results
The four approaches were applied to generate an irrigation schedule for the 2019 crop season. Therefore, the 2019 crop season climate data were not used in Approaches 1 and 3 completely, except in Approaches 2 and 4 for updating the climate files as the season progressed. Figure 2 shows the minimum and maximum air temperature and rainfall for the study area during the 2019 crop season.
summer maize were used to calibrate the AquaCrop model used in the current study. The crop season considered was from May to mid-August and simulations were performed for this period only. The planting date was 1 May. The length of the irrigation season was 100 days. The decision variable was a vector of 10 real numbers representing maximum allowable soil water depletion levels as percentages of AW. Irrigation was to be applied when soil water depletion levels reached these values. The irrigation application depth was 25 mm per irrigation during the first 30 days of the season, and 50 mm per irrigation afterwards until the end of the season. A constant irrigation application depth instead of a variable one presents some management advantages for the irrigator. In genetic algorithms, a population size of ~3 times the number of variables is recommended, therefore, in this study, the population size was 30. The number of iterations was 200. No attempt was made to optimize the probabilities of crossover and mutation, and they were both given a value of 0.3. A minimum irrigation interval ( ) of 5 days was adopted. The maximum allowable seasonal irrigation depth ( ) was taken as 450 mm. The lower ( ) and upper limits ( ) of maximum allowable soil water depletion levels were taken as 10% and 40% of AW, respectively. The initial soil water content was assumed to be at 80% of field capacity. The price of maize, fixed costs of production, fixed cost of irrigation water, variable cost of irrigation water, and efficiency of the irrigation system were assumed to be RMB33000/ton, RMB3750/ha, RMB30/ha, RMB0.157/m 3 , and 75%, respectively, and were assumed constant. In the climate data, 5 years were dry years, 11 years were normal years, and 5 years were wet years. The simulation-optimization model was developed and run on the MATLAB platform (MATLAB 2016a).

Results
The four approaches were applied to generate an irrigation schedule for the 2019 crop season. Therefore, the 2019 crop season climate data were not used in Approaches 1 and 3 completely, except in Approaches 2 and 4 for updating the climate files as the season progressed. Figure 2 shows the minimum and maximum air temperature and rainfall for the study area during the 2019 crop season. Improvement in objective value with iteration number for Approach 3, as an example, is shown in Figure 3. The shapes of the plots were typical for all the approaches. There were many Pareto optimal solutions at each iteration number but only one solution was selected at each iteration and Improvement in objective value with iteration number for Approach 3, as an example, is shown in Figure 3. The shapes of the plots were typical for all the approaches. There were many Pareto optimal solutions at each iteration number but only one solution was selected at each iteration and these are the ones shown in Figure 3. Generally, NB and WUE values increased with iteration number while R values decreased. Pareto front optimal solutions are mathematically considered equivalent to each other, and one objective cannot improve without causing a corresponding worsening of the other objectives. Table 1 shows the allowable depletion levels for the 2019 crop season as determined using the four different approaches.
these are the ones shown in Figure 3. Generally, NB and WUE values increased with iteration number while R values decreased. Pareto front optimal solutions are mathematically considered equivalent to each other, and one objective cannot improve without causing a corresponding worsening of the other objectives. Table 1 shows the allowable depletion levels for the 2019 crop season as determined using the four different approaches.
Approaches 2 and 4 had the lowest value for R and Approach 4 had the highest value for WUE. Approach 4 was the best in two objectives (WUE and R). Approach 2 and 3 were the best in one objective each, namely, R and NB, respectively. Therefore, under the conditions of this study, Approach 4 produced the best results, followed by Approaches 2 and 3. Approach 4 was used to determine irrigation schedules for a dry, normal, and wet year in the study area. Approaches 2 and 4 had the lowest value for R and Approach 4 had the highest value for WUE. Approach 4 was the best in two objectives (WUE and R). Approach 2 and 3 were the best in one objective each, namely, R and NB, respectively. Therefore, under the conditions of this study, Approach 4 produced the best results, followed by Approaches 2 and 3. Approach 4 was used to determine irrigation schedules for a dry, normal, and wet year in the study area.
Since weather varied between years even among those with the same climate, each year's weather resulted in a unique irrigation schedule for the same maximum allowable soil water depletion levels. Moreover, Approach 4 resulted in each year having a unique set of optimal allowable soil water depletion levels. Each year's crop seasonal irrigation schedule was tested on the other years sharing the same climate. The irrigation schedule resulting in the best average objective values was selected as the irrigation schedule for that climate. The selection of the best irrigation schedule was done using the ranking method. As an example, Figure 6 shows the objective values resulting from five unique irrigation schedules determined for the five dry years. For this particular example, the irrigation schedule which resulted in the best average objective values was the one in the 5th position and it was selected and is shown in Figure 7a. Shown also in Figure 7 are irrigation schedules for normal and wet year climates, determined in a similar manner. Since weather varied between years even among those with the same climate, each year's weather resulted in a unique irrigation schedule for the same maximum allowable soil water depletion levels. Moreover, Approach 4 resulted in each year having a unique set of optimal allowable soil water depletion levels. Each year's crop seasonal irrigation schedule was tested on the other years sharing the same climate. The irrigation schedule resulting in the best average objective values was selected as the irrigation schedule for that climate. The selection of the best irrigation schedule was done using the ranking method. As an example, Figure 6 shows the objective values resulting from five unique irrigation schedules determined for the five dry years. For this particular example, the irrigation schedule which resulted in the best average objective values was the one in the 5th position and it was selected and is shown in Figure 7a. Shown also in Figure 7 are irrigation schedules for normal and wet year climates, determined in a similar manner.
A dry year required 11 irrigations with a total seasonal irrigation depth of 400 mm. A normal year required nine irrigations with a total seasonal irrigation depth of 325 mm. A wet year required six irrigations with a total seasonal irrigation depth of 200 mm. A dry year required 11 irrigations with a total seasonal irrigation depth of 400 mm. A n year required nine irrigations with a total seasonal irrigation depth of 325 mm. A wet year re six irrigations with a total seasonal irrigation depth of 200 mm.

Discussion
In this study, a multi-objective irrigation scheduling optimization problem was solved four approaches. The crop which was investigated was maize. Figure 3 confirmed the ability optimization model to search for more optimal solutions as there was a clear trend of objective improvement with an increase in iteration number. The large dip in WUE between iterations 10 was due to sub-optimal solutions which applied more water than was necessary. As a resu yields were adversely impacted because there was poor soil aeration, causing a high risk of tot failure. Figure 2b showed that most of the rainfall for the 2019 crop season occurred towards July and in August. As a result, in Figure 4, none of the approaches recommended irrigatio about July 20. This was to maximize the effectiveness of rainfall since maximization of WUE w of the objectives of the optimization problem. The optimal depletion levels, such as the ones in Table 1, varied with crop season for Approaches 2 and 4 because of adaption, but were co for all seasons under Approaches 1 and 3. These optimal depletion levels were used in the Aqu model to determine irrigation dates. These depletion levels, as well as preferred irrigation appl depths, were specified in AquaCrop and, with this information and that of the climate management, the soil, and the crop, the model was able to generate irrigation schedule AquaCrop model has a daily time step and simulates the soil water balance to determin development and final yield.
Considering Figure 6, Approach 4, which was employed, involved adaption. As a result, fo dry year, there was a unique set of optimal depletion levels. This fact, and because of w variation among the dry years, resulted in a different irrigation schedule for each year. This w case for normal and wet years as well.
The authors hope that this study has demonstrated the usefulness of the presented fram for solving a multi-objective optimization problem of irrigation scheduling. Only three obj were considered, namely, maximizing NB and WUE and minimizing R. However, other obj could also have been considered or added, such as minimizing end-of-season soil salinity in t zone, minimizing the total number of irrigations in a season, minimizing water losses through and deep percolation, etc. The AquaCrop model, which was utilized to represent the co interactions of the soil-water-plant-atmosphere system, could also have been replaced by a suitable crop model.

Discussion
In this study, a multi-objective irrigation scheduling optimization problem was solved using four approaches. The crop which was investigated was maize. Figure 3 confirmed the ability of the optimization model to search for more optimal solutions as there was a clear trend of objective value improvement with an increase in iteration number. The large dip in WUE between iterations 0 and 10 was due to sub-optimal solutions which applied more water than was necessary. As a result, crop yields were adversely impacted because there was poor soil aeration, causing a high risk of total crop failure. Figure 2b showed that most of the rainfall for the 2019 crop season occurred towards end of July and in August. As a result, in Figure 4, none of the approaches recommended irrigation after about July 20. This was to maximize the effectiveness of rainfall since maximization of WUE was one of the objectives of the optimization problem. The optimal depletion levels, such as the ones shown in Table 1, varied with crop season for Approaches 2 and 4 because of adaption, but were constant for all seasons under Approaches 1 and 3. These optimal depletion levels were used in the AquaCrop model to determine irrigation dates. These depletion levels, as well as preferred irrigation application depths, were specified in AquaCrop and, with this information and that of the climate, field management, the soil, and the crop, the model was able to generate irrigation schedules. The AquaCrop model has a daily time step and simulates the soil water balance to determine crop development and final yield.
Considering Figure 6, Approach 4, which was employed, involved adaption. As a result, for each dry year, there was a unique set of optimal depletion levels. This fact, and because of weather variation among the dry years, resulted in a different irrigation schedule for each year. This was the case for normal and wet years as well.
The authors hope that this study has demonstrated the usefulness of the presented framework for solving a multi-objective optimization problem of irrigation scheduling. Only three objectives were considered, namely, maximizing NB and WUE and minimizing R. However, other objectives could also have been considered or added, such as minimizing end-of-season soil salinity in the root zone, minimizing the total number of irrigations in a season, minimizing water losses through runoff and deep percolation, etc. The AquaCrop model, which was utilized to represent the complex interactions of the soil-water-plant-atmosphere system, could also have been replaced by another suitable crop model. This study borrowed from other studies and attempted to develop a unique approach for solving a multi-objective irrigation scheduling optimization problem. Garg and Dadhich [33] also developed a model to optimize deficit levels for a number of crops. However, their study incorporated a multiplicative crop water production function of the Doorenbos and Kassam [2] type to simulate crop response to water, and only one value of the allowable maximum depletion level was determined for the whole season for each crop they considered. Li et al. [8] bootstrapped climate data to fit theoretical probability distributions to the data. They used the probability distribution functions to find parameter intervals. In the current study, Approaches 1 and 2 utilized the bootstrap resampling technique to generate climate data, but theoretical probability distributions were not fitted to the generated data as they were not found to fit well to the climate data used, a result also observed by Prats and Picó [28]. Interval mathematical programming was also not used to handle uncertainties in the climate data. Instead, the approach of Prats and Picó [28] was adopted in the current study. Prats and Picó [28] generated climate data by bootstrapping original climate data in a Monte Carlo simulation and did not fit a theoretical probability distribution function to the generated data. Linker and Kisekka [24] presented an approach for optimal irrigation scheduling by optimizing the maximum allowable soil water depletion levels. They used CERES-Maize as a surrogate crop and the AquaCrop model in the optimization procedure to determine the optimal soil water depletion levels, but only used static climate data. Jamal et al. [7] also used static climate data in their optimization models. Xu et al. [16] used a dynamic simulation model approach to optimize water resource allocation for irrigation but they utilized an empirical crop water production function to express the relationship between irrigation water amount and crop yield.
In the current study, Approaches 3 and 4 are similar to the explicit single-stage and explicit two-stage approaches presented by Jamal et al. [7] though that study used genetic algorithms, the Soil Water Atmosphere Plant (SWAP) model, the concept of aggregation periods, and only considered the single objective of maximizing net benefits.
Anvari et al. [34] presented an approach for optimal irrigation scheduling using stochastic dynamic programming. They also divided the crop season into 10-day periods and their study considered uncertainty in crop evapotranspiration during the season, but they applied a crop water production function of the type presented by Doorenbos and Kassam [2] to simulate crop yield response to irrigation.
In the current study, Approaches 2 and 4 used adaption and historical weather as forecasted weather. If reasonably accurate weather forecasts for the next few days were available, these could have been incorporated into the model. Approaches 1 and 3 required much less time to produce results for the whole season than Approaches 2 and 4. However, since the irrigator would only be concerned with the next few days and not the whole season at any given time, Approaches 2 and 4 would be comparable with Approaches 1 and 3 in terms of time required to produce results. Moreover, for Approaches 2 and 4, historical climate data, which were used as forecast climate data for the next 10 days, could be replaced with actual weather forecasts.
The results obtained in this study for Approach 4, when compared to Approach 3, are consistent with what [7] found in that the approach involving adaption produced better results than the one without adaption. Approaches 1 and 2 in the current study used climate data obtained by resampling original climate data and only 200 iterations were conducted. In future studies, the effect of the number of iterations on the performance of Approaches 1 and 2 should be established.
The results of the study are limited by the accuracy of the crop model in representing the soil-plant-atmosphere system. However, the developed framework can be used with any crop model. Furthermore, results of the study were only tested by simulation. Future studies can test the results in the field. Future studies can also explore the usefulness of utilizing more than one crop model in the optimization model in comparison to using only one crop model.

Conclusions
A multi-objective optimization problem of irrigation scheduling was presented and solved using four different approaches. All the approaches utilized the AquaCrop model and nondominated sorting genetic algorithm III in a simulation-optimization framework. Approach 1 used dynamic climate data without adaption. Approach 2 used dynamic climate data with adaption. Approach 3 used static climate data without adaption. Approach 4 used static climate data with adaption. The bootstrap technique was used for resampling original climate data to create the dynamic climate data. The results showed that Approach 4 produced the best results and that 400, 325, and 200 mm of irrigation water were required in a dry, normal, and wet year to produce a maize crop. The study demonstrated the usefulness of the presented simulation-optimization framework in solving a multi-objective optimization problem of irrigation scheduling. While the results of the study were limited by the accuracy of the crop model in representing the soil-water-plant-atmosphere system, the developed framework can be utilized with any other crop model. Furthermore, the results of the study were only tested by simulation. Future studies can test the results in the field. Future studies can also explore the usefulness of utilizing more than one crop model in the optimization model in comparison to using only one crop model.