Closed-Loop Scheduling and Control for Precision Irrigation

(5) McCarthy et al. used crop production models in MPC to optimize irrigation ...... T. C. A Receding Horizon Control algorithm for adaptive managemen...
0 downloads 0 Views 1MB Size
Subscriber access provided by WEBSTER UNIV

Process Systems Engineering

Closed-loop scheduling and control for precision irrigation Jannatun Nahar, Su Liu, Yawen Mao, Jinfeng Liu, and Sirish Lalji Shah Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b06184 • Publication Date (Web): 01 Mar 2019 Downloaded from http://pubs.acs.org on March 2, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Closed-loop scheduling and control for precision irrigation† Jannatun Nahar,





Su Liu,

Yawen Mao,

‡,¶

∗,‡

Jinfeng Liu,

and Sirish L. Shah



‡Department of Chemical & Materials Engineering, University of Alberta, Edmonton, AB T6G 2V4, Canada

¶Key Laboratory of Advanced Process Control for Light Industry, School of Internet of Things Engineering, Jiangnan University, Wuxi 214122, China E-mail: [email protected]

Abstract In agriculture irrigation management, irrigation scheduling is typically performed in an open-loop fashion and is done only once at the beginning of a growing season. In this work, we study whether closed-loop scheduling with closed-loop control can lead to improved performance in terms of crop yield and water conservation in agriculture irrigation. The interaction between soil water, crop (maize in this work) and atmosphere is described by an agro-hydrological model, which is a partial dierential equation. In the proposed scheduling and control scheme, both the scheduler and the controller are designed using model predictive control (MPC). The scheduler uses a long horizon (with a sampling period of one day) that covers the entire crop growth season and the horizon shrinks as time moves. The primary objective of the scheduler is to maximize the crop yield. The controller uses a much shorter prediction horizon and a much ner sampling †

J. Nahar and J. Liu would like to thank Prof.

Shah for his supervision and mentorship over the past

few years.

1

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

period. The primary objective of the controller is to track the soil moisture reference calculated by the scheduler. To alleviate the computational complexity of the scheduler and the controller, a linear parameter varying (LPV) model is identied for the scheduler and controller, respectively. The performance of the closed-loop scheduling scheme is evaluated against the traditional open-loop scheduling scheme under dierent scenarios. Keywords:

Irrigation scheduling, closed-loop irrigation, model identication, water conser-

vation, yield deciency.

1

Introduction

Water scarcity is an increasingly pressing issue globally. It is critical to increase the eciency of water usage especially in agriculture since it consumes the largest portion of fresh water (70% approximately). 1,2 Traditionally, agriculture irrigation decisions are primarily made in open-loop based on empirical or heuristic knowledge and no real-time feedback information (e.g., soil moisture) from the eld is used. These open-loop irrigation methods are not precise and typically lead to signicant over irrigation. In the past decade, there has been an increasing interest in lower-level closed-loop agriculture irrigation for precision irrigation. Dierent control algorithms based on real-time eld feedback like soil moisture measurement have been reported in the literature. For example, Goodchild et al. proposed a modied proportional-integral-derivative (PID) control to improve the irrigation precision. 3 Kim et al. developed a closed-loop irrigation control system for a self-propelled lateral-move sprinkler irrigation system. 4 Model predictive control (MPC) has also been investigated in the optimal control of agriculture irrigation. Park et al. used MPC to maintain the eld at desired soil moisture and salt levels. 5 McCarthy et al. used crop production models in MPC to optimize irrigation time and amount. 6 Delgoda et al. implemented MPC to minimize root zone soil moisture decit and irrigation amount with a limit on water supply. 7 Recently, Mao et al. developed a zone MPC algorithm to 2

ACS Paragon Plus Environment

Page 2 of 35

Page 3 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

maintain soil moisture within a target zone instead of a set-point. This approach was shown to lead to further water conservation. 8 The above mentioned irrigation control studies all focus on short time (hourly or daily) irrigation management. A typical objective of these studies is to maintain the soil moisture or root zone soil moisture decit at a pre-determined set-point or zone. Irrigation scheduling is another important subject in agriculture irrigation management, which considers water management over a much longer time period (e.g., the entire growing season of a crop). In general, a water balance accounting for plant water consumption, evaporation and drainage loss, historical weather data and avaiable soil moisture is used in irrigation scheduling. 911 Dierent objectives were considered in irrigation scheduling including, for example, prot maximization, 12 crop yield maximization 13 and irrigation uniformity optimization. 14 One common feature of the above scheduling studies is that the scheduling optimization problem is only solved once at the beginning of the crop growing season and the solution is used for the entire season. In this work, we introduce feedback into irrigation scheduling and integrate it with irrigation control to form a hierarchical closed-loop scheduling and control scheme. It is well recognized in chemical process industries that the closed-loop scheduling and control can bring improved economic performance. 1519 The primary objective of this work is to study the closed-loop scheduling and control in agriculture irrigation management and investigate whether closed-loop scheduling brings any benets to irrigation management. In the proposed closed-loop scheduling and control scheme, a scheduler operates in the top layer with the objective to maximize the crop yield. The scheduler calculates the target soil moisture value and passes it to a lower layer irrigation controller. The lower layer controller's primary objective is to follow the target soil moisture. Both the scheduler and the controller are designed in the framework of MPC. The scheduler uses a prediction horizon that covers the entire growing season of the crop (maize in this work). The scheduler is re-evaluated every week based on updated weather forecast and soil moisture measurement from the eld 3

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

with a shrinking horizon. The controller uses a prediction horizon of one week and calculates irrigation commands to follow the target soil moisture determined by the scheduler. The soil water dynamics of the eld are described by an agro-hydrological model which consists of a partial dierential equation. To alleviate the computational complexity of the scheduler and the controller, a linear parameter varying (LPV) model is identied for the scheduler and the controller, respectively, based on input and output data from the agro-hydrological model. The main contributions of this work include:

• A method to identify an LPV model for the scheduler to describe the interaction between soil water and the crop over the entire growing season that contains dierent growing stages. The LPV model is computationally ecient for online optimization.

• A detailed closed-loop scheduler and controller design to maximize crop yield and water conservation.

• Extensive simulations that investigate the performance of the closed-loop scheduler with respect to the traditional open-loop scheduling.

2

System description

The system considered consists of the soil-water-atmosphere-crop environment. The interaction between the soil water, the atmosphere and the crop is described by an agro-hydrological model. The yield of the crop is characterized using the ratio of the actual evapo-transpiration to the potential evapo-transpiration. The crop considered is maize. A schematic of the considered system is shown in Figure 1.

2.1 Agro-hydrological model We consider an agro-hydrological system that characterizes the interaction between the soil water, the atmosphere and the crop. The system is modeled by a set of empirical and 4

ACS Paragon Plus Environment

Page 4 of 35

Page 5 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Rain/ Irrigation

Transpiration Evaporation

Interception

Surface runoff/ runon

soil layers

Ponding

Infiltration

Rootwater extraction

BABY MAIZE PLANT

INTERMIDIATE MAIZE PLANT

Drainage

MATURED MAIZE PLANT

Figure 1: The system considered consists of the soil-water-atmosphere-crop environment over the entire growing season of the crop maize. physical relations. The dynamics of the water ow inside soil is described using the Richards equation 20 :

∂ ∂h ∂θ = [K(h)( + 1)] − S(h) ∂t ∂z ∂z

(1)

where θ is the soil moisture content (cm/cm), K is the hydraulic conductivity (cm/hr), z is the vertical distance inside soil matrix (cm) taken positive upwards, S includes the source and sink terms (hr−1 ), such as soil water extraction rate by plant roots and h is the matric potential (cm) also known as soil suction. The hydraulic conductivity K can be expressed in terms of soil moisture content as follows 21 : 1

K = Ksat Seλ [1 − (1 − Sem )m ]2

(2)

where Ksat is saturated hydraulic conductivity, λ = 0.5 is the shape parameter 22 and Se is the relative saturation given by:

Se =

θ − θres θsat − θres 5

ACS Paragon Plus Environment

(3)

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 35

where θres is the residual moisture content (cm/cm), θsat is the saturated moisture content (cm/cm). The relation between soil moisture content θ and matric potential h is described by the Van Genuchten model: 23

θ(h) = θres + (θsat − θres )(1 + |αv h|η )−m

(4)

where η and m are empirical shape factors, and αv is the inverse of air entry suction (cm−1 ). In general, m is considered as a dependent variable of η such as m = 1 − η1 . Thus, the Richards equation can be described as follows:

C(h)

∂ ∂h ∂h = [K(h)( + 1)] − S(h) ∂t ∂z ∂z

(5)

where C is the hydraulic capacity, which accounts for the change of soil moisture due to the change of soil matric potential; that is, C = dθ/dh. In the Richards equation (5), S(h) represents the soil water extraction rate by the crop roots which accounts for the transpiration rate T . The top water ux qtop at the soilatmosphere interface is considered in the top boundary calculation of the Richards equation. The top water ux qtop can be expressed as follows:

qtop = I + P − E

(6)

where qtop is the top ux, I is the irrigation rate, P is the precipitation rate and E denotes evaporation. The transpiration T and evaporation E are calculated following the method in SWAP model. 24 The top boundary of this agro-hydrological system is time varying depending on the surface soil moisture condition, and may switch between the ux and the pressure head boundary conditions. Both of these boundary conditions are determined using the qtop value as described in the work of Van Dam and Feddes. 25 In this study the bottom boundary condition is taken as free drainage.

6

ACS Paragon Plus Environment

Page 7 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Following the method used in the work of Van Dam et al. 24 and Nahar et al., 26 equation (5) can be discretized using the implicit nite dierence method with explicit linearization of the hydraulic conductivity K(h). Let us consider that the spatial domain is discretized into ND compartments with the top one as compartment 1 and the bottom one as compartment ND . After discretization, at compartment i (i = 2, . . . , ND − 1), the water content dynamics are described as follows:

di (h(n))hi−1 (n + 1) + bi (h(n))hi (n + 1) + ai (h(n))hi+1 (n + 1) = ei (h(n)) − ui (n)

(7)

where hi is the matric potential of the i-th compartment, n denotes the discretized time with Ki−1/2 (h) Ki+1/2 (h) , di (h) = − 1 , bi (h) = a sampling period of a hour, ai (h) = − 1 (∆zi + ∆zi+1 ) (∆zi−1 + ∆zi ) 2 2 ∆zi ∆zi Ci (hi ) − ai (h) − di (h), ei (h) = Ci (hi ) hi + Ki−1/2 (h) − Ki+1/2 (h), ui (k) = Si (k)∆zi , ∆t ∆t and ∆zi denotes the thickness of compartment i. In equation (7), Ki+1/2 denotes the average hydraulic conductivity between compartment i and i + 1. The boundary conditions of the system dene i = 1 and i = ND . Dening h = [h1 h2 ...hND ]T , combining the equations for the ND compartments together, and splitting the inputs to the system into manipulated inputs and uncontrolled inputs, the overall dynamics can be written in the following compact form:

h(n + 1) = f (h(n), u(n), w(n))

(8)

where u is the manipulated input (irrigation) and w denotes disturbances and uncontrolled inputs such as rain and weather conditions. Note that h is a column vector whose size depends on the number of discretization points ND in the spatial domain. Based on equation (4), the soil moisture is a function of the matric potential h. For the soil moisture at each compartment i, we can use equation (4) to nd the corresponding matric potential hi . For brevity, we write these equations in the following compact form:

θi (n) = g(hi (n)), i = 1, . . . , ND 7

ACS Paragon Plus Environment

(9)

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 35

where the function g represents the relation between θ and h as expressed in equation (4). Note that when i = 1, it indicates the top surface. Equations (8) and (9) dene the system equations for the dynamics of water ow inside soil. For more details of the derivation of equations (8) and (9), the reader may refer to Nahar et al. 26

2.2 Crop yield model The crop yield is modeled as a function of the potential evapo-transpiration and the actual evapo-transpiration as follows: 12 T

X Ya ETa (1 − ) = Ky (k)(1 − ) Yp ET p k=1

(10)

where Ya is the actual yield, Yp is the potential maximum yield, ETa is the actual evapotranspiration and ETp is the potential evapo-transpiration, and Ky (k) is the crop sensitivity factor (also known as yield response factor) for the growing period at time k . Crop sensitivity factor is a measure used to indicate the impact of water stress (at dierent growing stages) on the overall yield of the crop. When the actual yield is equal to the potential yield, Ya equation (10) takes its minimum value zero. In (10), (1 − ) represents the yield reduction Yp ETa or the yield deciency and 1 − represents the relative reduction in evapo-transpiration ETp or ET deciency. The potential and the actual evapo-transpiration values of a crop are dened as follows: 27

ETp = Kc ET0

(11)

ETa = Ks (θ)ETp

(12)

where Kc is an empirical crop factor that describes the relation of ETp to the potential evapo-transpiration of the reference crop (grass) ET0 , and Ks (θ) is the water stress factor which is a function of the soil moisture and relates ETa to ETp .

8

ACS Paragon Plus Environment

Page 9 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

In equation (11), the reference evapo-transpiration ET0 is calculated using PenmanMontieth equation based on hourly weather data. 27 In equation (12), the water stress factor

Ks (θ) is modeled following the method in Feddes's work 28 and will be made clear in Section 3.3.1. It is noted that equation (12) holds when other factors such as nutrient deciency or diseases are not present. Combining equation (10) and (12), the following relation can be obtained to describe the relation between the crop yield and soil moisture: T

X Ya (1 − ) = Ky (k)(1 − Ks (θ(k))) Yp k=1

(13)

In this work, following the work of Van Dam et al., 24 a simple maize model is used to calculate the evapo-transpiration amounts at dierent growing stages. Note that the crop water requirement increases as the crop grows. This implies that more water is required to keep the soil moisture at the same level as crop grows.

3

Problem formulation and proposed design

3.1 Closed-loop scheduling and control problem formulation In this section, we describe the proposed closed-loop scheduling and control framework for precision irrigation. The scheduler and the controller form a hierarchical decision making system. In this work, both the scheduler and the controller will be designed in the framework of MPC. The primary objective of the scheduler is to maximize the crop yield at the end of the entire growing season. It considers a much longer horizon that covers the entire growing season and uses historical weather data, available weekly weather forecast and soil moisture measurement from the eld to calculate the soil moisture reference trajectories for the lower layer controller. The primary objective of the controller is to track the soil moisture reference calculated by the scheduler. It uses a xed prediction horizon of a week and takes into 9

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2: Structure of the closed-loop scheduling and control system. Left: Process block diagram; Right: The design of prediction horizon. account the soil moisture measurement and daily weather forecast in its irrigation command calculation. A schematic of the closed-loop scheduling and control framework is shown in Figure 2.

3.2 Scheduling and control model identication The direct use of the agro-hydrological model, equations (8) and (9), in the scheduler and the controller is computationally challenging due to the stiness of the nonlinlear system when the soil moisture is low, the high dimensions of the discretized system, and the large horizon covered by the scheduler. In order to overcome the computational challenge and to capture the dynamics between the input and the output of the agro-hydrological system, we propose to identify an LPV model for the scheduler and the controller, respectively. The two LPV models can be identied based on the simulated input-output data from the agrohydrological model. LPV models have been widely used to describe nonlinear processes due

10

ACS Paragon Plus Environment

Page 10 of 35

Page 11 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

to its simple structure and its ability to approximate complex nonlinear systems. 2932

3.2.1

Scheduling model

The output of the scheduling model is the average root zone soil moisture content described as follows:

Nr 1 X θi θs = Nr i=1

(14)

where θs is the average root zone soil moisture, θi is determined by equation (9), and Nr < ND is the deepest node that contains the crop roots. Note that Nr changes over time as the length of roots grows. The LPV model consists of a set of linear models that are connected with varying parameters which are functions of scheduling variables. 33 The scheduling variables are the ones that can be used to reect the nonlinearity of the dynamics of the system. For the scheduling model, two scheduling variables are used. They are the soil moisture s1 (s1 = θ) and the root depth s2 . The sampling period used in the scheduling model is one day. For each scheduling variable, a few working points need to be determined and around those working points, linear models are identied. Let us assume that we have H working points for s1 and G working points for s2 . These working points dene a mesh with HG nodes. For each of these nodes, a linear model should be identied. Let us also use (i, j) (i = 1, 2, ..., H , j = 1, 2, ..., G) to denote the node dened by the i-th s1 working point and the j -th s2 working point. For each node (i, j), a linear model of the following form is identied:

xi,j (m) = −

na X

ai,j p xi,j (m

− p) +

p=1

nb X

i,j bi,j q u(m − q) + d

(15)

q=1

where m denotes the discrete time with a sampling period of a day, xi,j is the output of i,j i,j the linear model for (i, j), ai,j are model p (p = 1, 2, ..., na ), bq (q = 1, 2, ..., nb ) and d

parameters, u is the manipulated input. The model orders na and nb are determined prior

11

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 35

to the identication process. This model can also be expressed as follows:

xi,j (m) = ϕTi,j (m)βi,j

(16)

where ϕi,j (m) = [−xi,j (m − 1), ..., −xi,j (m − na ), u(m − 1), ..., u(m − nb ), 1]T , and βi,j = i,j i,j i,j i,j T [ai,j 1 , ..., ana , b1 , ..., bnb , d ] .

With respect to s1 (soil moisture), the linear models are combined together using Gaussian weighting functions as follows:

y∗j (m) =

H X

αi (s1 (m − 1))xi,j (m)

(17)

i=1

where y∗j is an intermediate variable and the weighting function αi is given by:

ωi (s1 (m)) αi (s1 (m)) = PH i=1 ωi (s1 (m))

(18)

where

 (−sign(s1 (m) − s1,i ) + 1) −(s1 (m) − s1,i )2 ωi (s1 (m)) = exp − × 2 2 2σi1   (−sign(s1 (m) − s1,i ) + 1) −(s1 (m) − s1,i )2 × + exp − −1 2 2 2σi2 

with s1,i the soil moisture value of the i-th pre-specied s1 working point, σi1 and σi2 being two parameters that needs to be identied which is related to the validity width of the linear model. 8 With respect to s2 (root length), the outputs obtained from equation (17) are further linearly interpolated to get the nal output of the LPV model:

y(m) =

s2,j+1 − s2 (m) j+1 s2 (m) − s2,j j y∗ (m) + y (m) s2,j+1 − s2,j s2,j+1 − s2,j ∗

12

ACS Paragon Plus Environment

(19)

Page 13 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

given that the root length at time m lies between s2,j and s2,j+1 (i.e., s2,j ≤ s2 (m) ≤ s2,j+1 ). Note that the length of the roots grows according to the time after seeding and is considered as a function of time in this work. In equation (19), y is the average water content within the root zone given by the LPV model. Equation (19) is the model that will be used in the scheduler and may be written equivalently in the following state-space form:

x(m + 1) = Ax(m) + Bus (m) + d

(20)

y(m) = C(s1 (m − 1), s2 (m))x(m)

(21)

where A, B are composed of elements in βi,j (i = 1, 2, . . . , H, j = 1, 2, . . . , G), C depends on the two scheduling variables, d is composed of di,j , and us is composed of the current and previous daily irrigation amounts. The detailed identication of the above LPV model will be discussed in Section 4.

3.2.2

Control model

Since the controller uses a relatively much smaller prediction horizon (a week), we consider an LPV model with only the soil moisture (s1 ) as the scheduling variable. In this model, we consider the root length of the crop is constant within each control prediction horizon. The output of the control model is also the average soil water content within the root zone as dened in equation (14). The sampling period of the LPV model for control purposes is also a day. The LPV takes the following form:

x˜i (m) = −

na X

aip x˜i (m

− p) +

p=1

y˜(m) =

G X

nb X

biq u(m − q) + di

(22a)

q=1

αi (˜ s1 (m − 1))˜ xi (m)

(22b)

i=1

where x ˜i is the output of a sub-model and y˜ is the nal output of the LPV model, aip (p =

1, 2, ..., na ), biq (q = 1, 2, ..., nb ) and di are model parameters, αi is a weighting factor deter13

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 35

1

κ

0 Soil moisture

Figure 3: The relation of soil moisture content and water stress factor 28 . mined following equation (18). At dierent crop growing stages, the above control model should be re-identied to account for the growth of the crop. The above model can also be written equivalently in the following state-space form:

˜x(m) + Bu ˜ c (m) + d˜ x˜(m + 1) = A˜

(23)

˜ 1 (m − 1))˜ y˜(m) = C(s x(m)

(24)

˜ are composed of ai and bi (i = 1, 2, . . . , H ), C˜ depends on the scheduling where A˜, B q p variable, d˜ is composed of di , and uc is composed of the current and previous daily irrigation amounts. The details will be discussed also in Section 4.

3.3 Design of scheduler and controller 3.3.1

Design of the scheduler

As explained earlier, the main objective of the scheduler is to minimize the crop yield deciency while reducing the irrigation water usage. From equation (13), it can be seen that crop yield relates to the water stress factor Ks . The stress factor characterizes the suppression of root water update due to either drought or insucient aeration. Figure 3 shows a typical relation of Ks and soil moisture 28 . When the soil moisture θ is between θL (limiting point) 14

ACS Paragon Plus Environment

Page 15 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

and θa (anaerobic point), there is no water stress (Ks = 1). When the soil moisture θ is below θL but larger than θw (wilting point), water stress takes place due to drought. When the soil moisture θ is above θa but less than θu , there is crop water stress due to insucient aeration; when θ is above θu , crop transpiration stops. For all other soil moisture region, root water extraction is completely suppressed. Therefore, it is also important to ensure that the soil moisture is always between θw and θu . Specically, the value of the water stress factor Ks is determined as follows: 28

 θ − θw   ,   θL − θw      1, Ks = θ − θa   ,   θU − θa      0,

for θw ≥ θ ≤ θL for θL < θ < θa for θa ≥ θ ≤ θU

(25)

for θ ≤ θw or θ ≥ θU

The scheduler optimization problem is evaluated every week. It uses a prediction horizon from the current week to the end of the growing season T . Within each week, the irrigation

15

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 35

is considered to be constant. For week k , the optimization problem is shown below:

( min

u,ε1 ,ε2

Ya 1− Yp

2 Qs1 +

T −1 X

2

u(l) Qs2 +

l=k

T X 

2

2

ε1 (l) Qs3 + ε2 (l) Qs4

) 

k T X X Ya Ky (l)(1 − Ks (θs (l))) + s.t.(1 − ) = Ky (l)(1 − Ks (y(l))) Yp l=1 l=k+1

x(l + 1) = A7 x(l) +

6 X

(26a)

l=k

A6−m Bu(l) +

m=0

6 X

A6−m d, l = k, . . . , T − 1

(26b) (26c)

m=0

y(l) = Cx(l), l = k, . . . , T

(26d)

x(k) = x¯(k)

(26e)

0 ≤ u(l) ≤ umax , l = k, . . . , T − 1

(26f)

y(l) − ε1 < θu , y(l) + ε2 > θw , l = k, . . . , T

(26g)

ε1 (l) ≥ 0, ε2 (l) ≥ 0, l = k, . . . , T

(26h)

In the optimization problem (26), equation (26a) denes the cost function to be minimized. In equation (26a), the rst term indicates the crop yield deciency at the end of the growing season, the second term measures the irrigation amount in the remaining growing season, the last summation term includes two slack variables that are used to handle constraints on the soil moisture in equation (26g), and Qs1 , Qs2 , Qs3 are positive weighting factors. Equation (26b) is the model used to calculate the crop yield deciency. Note that in equation (26b) the rst part on the right-hand-side is calculated based on actual soil water content measurements and the second part is calculated based on soil water content predicted in the scheduler. Equations (26c)-(26d) are the LPV model for scheduling evaluated for a week with constant input. Note that the sampling period of the model identied in the previous subsection is a day. The scheduler optimization problem, however, evaluates once every week. Within a week, the input is held constant. Equation (26c) reects this type of consideration. Equation (26e) denes the initial condition of the scheduling optimization problem and x ¯(k) denotes the current state measurement. Note that we assume that the state is measured. 16

ACS Paragon Plus Environment

Page 17 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Equation (26f) requires that all the irrigation should be non-negative and is less than the maximum allowed limit umax . Equation (26g) imposes the constraint such that θw < y < θu as soft constraints through the introduction of two slack variables ε1 and ε2 . Equation (26h) implies that the two slack variables are non-negative. The weighting factors Qs3 and Qs4 should be tuned to be relatively higher to avoid constraint violation. The solution to the above scheduling optimization problem is a sequence of optimal irrigation amount from week k to T . Let us denote the optimal weekly irrigation as u∗ (l|k), l =

k, . . . , T −1. Let us also denote the corresponding weekly optimal soil moisture as y ∗ (l|k), l = k, . . . , T − 1. The rst two values of the optimal irrigation values u∗ (k|k), u∗ (k + 1|k) and the corresponding soil moisture values y ∗ (k|k), y ∗ (k + 1|k) are sent to the controller as references. At the next week, the scheduling optimization problem is re-evaluated with new soil moisture measurement.

3.3.2

Design of the controller

In the controller design, the control LPV model is rst augmented with an output disturbance to account for model-plant mismatch. Then, an observer is designed to estimate the modelplant mismatch. The augmented model is as follows:

˜x(m) + Bu ˜ c (m) + d˜ x˜(m + 1) = A˜

(27)

p(m + 1) = p(m)

(28)

y˜(m) = C˜ x˜(m) + p(m)

(29)

17

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 35

Based on the above augmented model, the following observer is designed to estimate the output disturbance:

xˆ(m|m) = xˆ(m|m − 1) + Lx (m)(¯ y (m) − yˆ(m|m − 1))

(30)

pˆ(m|m) = pˆ(m|m − 1) + Lp (m)(¯ y (m) − yˆ(m|m − 1))

(31)

˜x(m|m) + Bu ˜ c (m) + d˜ xˆ(m + 1|m) = Aˆ

(32)

pˆ(m + 1|m) = pˆ(m)

(33)

yˆ(m|m) = C˜ xˆ(m + 1|k) + pˆ(m + 1|m)

(34)

where x ˆ, pˆ, yˆ and y¯ denote the estimated state, estimated disturbance, estimated output, and output measurement, respectively. The observer gains Lx and Lp are determined as discussed in Mao et al. 8 The primary objective of the controller is to track the soil moisture and total irrigation water consumption references calculated by the scheduler. The controller is also designed in the framework of MPC. The controller is evaluated every day with updated soil moisture measurements. The controller uses a xed prediction horizon of a week (i.e., 7). At day m, the controller is formulated as follows:

18

ACS Paragon Plus Environment

Page 19 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

min

u,ξ1 ,ξ2

m+7 X

[ˆ y (l) − y ∗ (l)]2 Qc +

l=m

"m+6 X

u(l) −

l=m

m+6 X

#2 u∗ (l)

Rc +

l=m

m+7 X

  ξ1 (l)2 Qc1 + ξ2 (l)2 Qc2 (35a)

l=m

˜ l = m, . . . , m + 6 ˜x(l) + Bu(l) ˜ s.t.x ˜(l + 1) = A˜ + d,

(35b)

p(m + 1) = p(m)

(35c)

y˜(l) = C˜ x˜(l) + p(l), l = m, . . . , m + 7

(35d)

x˜(m) = xˆ(m|m),

(35e)

p(m) = pˆ(m|m)

0 ≤ u(l) ≤ umax /7, l = m, . . . , m + 6

(35f)

y˜(l) − ξ2 (l) ≤ θu , y˜(l) + ξ1 (l) ≥ θw , l = m, . . . , m + 7

(35g)

ξ1 (l) ≥ 0, ξ2 (l) ≥ 0, l = m, . . . , m + 7

(35h)

In the optimization problem (35), equation (35a) denes the cost function to be minimized. In equation (35a), y ∗ and u∗ denote the references from the scheduler. The rst term requires that the predicted average soil moisture content within the root zone should be close to the reference from the scheduler, the second term requires that the total water consumption over a week is close to the predicted consumption of the scheduler, and the last terms contain two slack variables. In equation (35a), Qc , Qc1 , Qc2 and R are positive weighting factors. Equation (35b) is the LPV model for control. Equation (35c) implies that a constant output disturbance is considered in the optimization of the controller and the disturbance value is given by the observer as specied in equation (35e). This constant output disturbance is added to the output equation as specied in equation (35d). Equation (35e) denes the initial condition of the control optimization problem and x ˆ is from the observer. Equation (35f) requires that all the irrigation should be non-negative and is smaller than or equal to one seventh of the weekly upper limit. Equation (35g) imposes the constraint such that θw < y˜ < θu as soft constraints. Equation (35h) implies that the two slack variables are non-negative.

19

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4: Root length of maize according to time after seeding. The solution to the optimization problem (35) is a sequence of optimal irrigation amount for time l = m, . . . , m + 6 denoted by u∗ (l|m), l = m, . . . , m + 6. At time m, only the rst element from the optimal daily irrigation input is applied to the eld. The rest are discarded. At the next day m + 1, the optimization problem is solved again

4

Results and discussion

4.1 LPV model identication and validation In the identication of the scheduling LPV model, three working points are considered for each of the two scheduling variables (H = 3, G = 3). In the input-output data generation, the three s2 working points are s2 = [s2,1 , s2,2 , s2,3 ]T = [10.3571, 55.3571, 99.881]T cm. These working points are picked according to the root length of maize at dierent growing stages (see Figure 4). Specically, the rst working point of s2 corresponds to the root length about ten days after seeding and is picked to represent the baby plant. The second working point is the mean value of the root length of the rapid growing stage and the last working point is the root length at the fully matured stage. For each working point of s2 , 20

ACS Paragon Plus Environment

Page 20 of 35

Page 21 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 5: Input-output data for model identication with two scheduling variables. three dierent s1 working points are considered. For s2,1 , the three s1 working points are {1}

{1}

{1}

{2}

{2}

{2}

{3}

{3}

{3}

{1}

= [s1,1 , s1,2 , s1,3 ]T = [0.0755, 0.0936, 0.1792]T . For s2,2 , the three s1 working points are

{2}

= [s1,1 , s1,2 , s1,3 ]T = [0.0789, 0.1211, 0.1988]T . For s2,3 , the three s1 working points are

{3}

= [s1,1 , s1,2 , s1,3 ]T = [0.0835, 0.1176, 0.1564]T . For each s2 working point, the three s1

s1

s1

s1

working points are determined according to three (small, medium, large) steady-state irrigation amounts, which are picked based on extensive simulations and trail-and-error. In determining these irrigation amounts, one consideration is the possible soil moisture range of the system (about 0.05 - 0.35) and another consideration is the nonlinearity feature of the agro-hydrological system. The orders of the local models na and nb are both picked to be one based on the results of Mao et al. 8 The input and output data for the LPV identication generated from the agro-hydrological model is shown in Figure 5. In generating the data, the daily weather condition is kept constant with no rain. Binary input signals are used around the working points. During the transition from one working point to another of s1 (soil moisture), constant inputs are used. 21

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 35

In the top plot of Figure 5, CS1, CS2, CS3 correspond to the three working points of s2 (root depth). During these regions, binary inputs are used. Region CV 1 and CV 2 correspond to transition periods from s2,1 to s2,2 and from s2,2 to s2,3 . In these two regions, constant inputs are used to show the eect of crop root zone on soil moisture. As crop root length increases (from CS1 to CS3), more irrigation is required to maintain the soil moisture as shown in the top plot of Figure 5. Following the approach in Mao et al., 8 the parameters of the LPV model are identied and shown in Table 1. The bottom plot of Figure 5 shows the predicted Table 1: Estimated local model parameters and weights for the scheduling LPV model. Soil moisture \Root length 1 2 3

Gaussian weights

a1 b1 d1 a1 b1 d1 a1 b1 d1 σ12 2 σ21 2 σ22 σ32

1 -0.7749 0.0152 0.0144 -0.8745 0.0542 -0.0045 -0.5701 0.0264 0.0575 0.0201 0.0094 0.0200 0.0205

2 -0.6337 0.0074 0.0235 -0.7334 0.0083 0.0106 -0.8873 0.0065 -0.0072 0.0193 0.0129 0.0210 0.0203

3 -0.4202 0.0022 0.0433 -0.1006 0.0054 0.0643 -0.3665 0.0036 0.0556 0.0194 0.0157 0.0209 0.0186

trajectory and the actual trajectory generated by the agro-hydrological model. The LPV model overall has a good t over the entire period with the percentage dierence between the predicted value of the LPV model and the actual value of the agro-hydrological model about 5% or less except in the region CV 1 which corresponds to the rapid growing stage (whose percentage dierence is about 18%). To further validate the developed LPV model over the crop growing season, an open-loop simulation of the agro-hydrological model is performed. Fixed weather condition is considered with no rain, and an input irrigation of 2.5 cm/day is considered. The result is shown in Figure 6. The overall percentage dierence between the predicted value given by the LPV model and the value generated by the agro-hydrological model is about 27% with most of it contributed from the period between day 23 and day 22

ACS Paragon Plus Environment

Page 23 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 6: Validation of the scheduling LPV model. 45 which corresponds to the rapid growing stage. In the rest of the period, the percentage dierence is around 15%. Given the complexity of the agro-hydrological system and the use of feedback, the mismatch between the LPV model and the agro-hydrological model is considered to be acceptable. Note that the performance of the LPV model may be improved by including more working points. In the controller, the three local linear models of each root length working point are weighted using the corresponding Gaussian weights in Table 1 and then is used in the controller design. For example, in CS1 and up to the middle of CV 1, the local models of root length working point 1 are weighted and used in the MPC controller. The switching between root length working points is based on the root length.

4.2 Scheduling and control simulation results The proposed closed-loop scheduling and control is applied to the agro-hydrological system under dierent scenarios. The proposed closed-loop scheduling approach is compared with an open-loop scheduling scheme in which the scheduling optimization is only solved once 23

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Table 2: The values of crop sensitivity factor. Crop sensitivity factor (Ky ) 0.4000 0.4000 0.9000 0.9000 1.5000 1.5000 2.3000 2.3000 0.5000 0.5000 0.2000 0.2000

Development stage (DVS) 0 0.2375 0.2500 0.5625 0.5750 0.7500 0.7625 1.3750 1.3875 1.6250 1.6375 2.0000

at the beginning and applied throughout the growing season. In the open-loop scheduling scheme, the same closed-loop controller is used. The crop growing season is considered to be 23 weeks from seeding. The values of the crop sensitivity factor Ky at dierent crop growth stages (from 0 to 2 with 0 indicating initial seeding stage and 2 implying mature stage) for maize are shown in Table 2. 34 The weekly irrigation upper limit is umax = 70 cm/week , which implies that the daily irrigation limit is 10 cm.

4.2.1

Scenario 1: Favorable weather

In the rst scenario, we consider that the amount of rain in the current crop season is similar to that of the historical weather data. In this scenario, we consider that the actual total rain in the crop growing season is 164.05 cm but the historical weather data suggests the total rain over the season to be 154.50 cm. In the closed-loop scheduling, the scheduler is updated at the end of each week with an update weather forecast for the next two weeks and the current soil moisture measurement. While in the open-loop scheduling, only the historical weather data is used at the beginning of the season. As for the controller, in both schemes, the eld soil moisture is updated every day with eld measurements and the weather data is also updated for the next seven days.

24

ACS Paragon Plus Environment

Page 24 of 35

Page 25 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 7: Simulation results of the closed-loop scheduling scheme under favorable weather scenario with more weighting on water conservation. Uncertainties are considered in weather forecast. First, we consider a water conservation case. That is, in the scheduler design, more emphasis is put on water conservation. In this case, the weighting factors are chosen as

Qs1 = 10002 , Qs2 = 100, Qs3 = 109 and Qs4 = 1. The value of Qs3 is chosen very high to ensure that the constraint on the lower soil moisture limit is not violated. The weight on the upper soil moisture limit is relatively low since it is less important and can be controlled to some extent by the irrigation amount. Compared with the next case to be studied, Qs2 = 100 puts more weighting on water consumption. The values for θu and θw are taken as 0.41cm/cm and 0.15 cm/cm. The wilting value θw for the system is considered to be 0.09 cm/cm. Since there is model-plant mismatch, a lower limit that is higher than the wilting point is used. For the MPC controller, the weighting factors are Qc = 100, Rc = 0.001, Qc1 = 1000 and Qc2 = 0.01. In the controller, a dierent lower limit 0.13 cm/cm on the soil moisture is used since the controller gets daily feedback of soil moisture from the farm and it can

25

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 8: Simulation results of the open-loop scheduling scheme under favorable weather scenario with more weighting on water conservation. handle uncertainty better. The simulation results are shown in Figure 7 for the closed-loop scheduling scheme and in Figure 8 for the open-loop scheduling scheme. For simplicity, only one day weather forecast is shown in the gures. In both schemes, the MPC controller overall is able to track the soil moisture targets well. The over prediction by the LPV model subplot (c) of Figures 7 and 8 is relatively signicant in the crop rapid growing stage which is due the relatively big mismatch of the LPV model from the actual system in the rapid growing stage. The peak at around day 15 is due to the sudden increase in the actual rain that is not captured in the weather forecast. Table 3 summarizes the results of this case. From the table, it can be seen that the open-loop scheme has a slightly (1.5%) better crop yield deciency while the closed-loop scheme has a slightly reduced (1.9%) irrigation amount. Second, we consider a case that we put less emphasis on water conservation (i.e., more on crop yield reduction). In the scheduler, the weights are Qs1 = 10004 , Qs2 = 1, Qs3 = 109 and

Qs4 = 1. The lower soil moisture constraint in the scheduler is decreased to θw = 0.13cm/cm

26

ACS Paragon Plus Environment

Page 26 of 35

Page 27 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

whereas for the controller the lower limit is increased and set to be the set-point from the scheduler. The weights in the controller cost function are kept same as in the previous case. Thus the controller in this case tries to track the soil moisture set point from some higher values thereby ensuring less yield deciency. A summary of the results is given in Table 3. In this case, we can see that the closed-loop and the open-loop schemes give very similar crop yield while the closed-loop scheme uses slightly (0.8%) less water. Table 3: Simulation results under favorable weather scenario.

Yield deciency Irrigation implemented by controller (cm3 /cm2 ) Irrigation prescribed by scheduler(cm3 /cm2 )

4.2.2

More water conservation Closed-loop Open-loop 12.47 12.29

Less water conservation Closed-loop Open-loop 5.71 5.71

990.10

1009.30

1222.40

1211.90

1020.50

993.44

1266.80

1242.00

Scenario 2: Dry weather

In this scenario, we consider that the amount of rain in the current crop season is signicantly less than the total historical rain data. In this set of simulations, we consider that there is

79.44 cm of rain in the crop growing season. As in the previous scenario, we rst consider a water saving mode in the scheduler. The scheduler and controller settings are the same as in the previous scenario. A summary of the results is given in Table 4. It can be observed that the closed-loop scheme saves more (2.3%) water relative to the open-loop scheme at the cost of (8.2%) less yield. We also consider a similar case like in the previous scenario to put less weight on water conservation. The results are summarized in Table 4. From the table, a similar conclusion can be drawn. That is, the closed-loop scheme is able to use (2.5%) less water but at a cost of (5.9%) less yield. To understand the eects of weights in the scheduling optimization, let us look at the results of the same scheduling scheme (closed-loop or open-loop) in dierent operating modes 27

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 35

Table 4: Simulation results under dry weather scenario.

Yield deciency Irrigation implemented by controller (cm3 /cm2 ) Irrigation prescribed by scheduler(cm3 /cm2 )

More water conservation Closed-loop Open-loop 17.94 16.57

Less water conservation Closed-loop Open-loop 9.73 9.12

1055.30

1079.60

1234.50

1266.80

1092.50

993.56

1273.70

1245.10

(more water conservation and less water conservation). From Table 3 (columns 2 & 4), it can be observed that when the objective in the closed-loop scheduler is changed to consider less water conservation, an increase of 54.24% in yield can be achieved with the cost of a 23.46% increase of water usage under the favorable weather condition. From Table 3 (columns 3 & 5), it can be seen that an increase of 53.52% in yield can be achieved with the cost of a 20.07% increase of water usage for the open-loop scheduler under the favorable weather condition. This may imply that the yield is more sensitive to the weights and a relatively smaller change in the consumption of irrigation water may result a much larger change in the yield. A similar pattern can be observed from the results of Table 4 for the dry weather condition. In addition to the above cases, we consider another case in which the water price varies. The water price is increased in the relatively dry season between week 10 and 15 corresponding to day 70 and 105. The corresponding Qs2 value changes according to the price as follows:

Qs2 (k) =

  

1, if k < 10 or k > 15

(36a)

if 10 ≤ k ≤ 15

  100,

The results are shown in Figure 9 for the closed-loop integrated scheme and in Figure 10 for the open-loop irrigation scheme. A summary of the result is given in Table (5). The result shows that less water is used relatively to open-loop scheme. In the gures, it can be observed that a tendency to use more water when water price is low presents. The soil 28

ACS Paragon Plus Environment

Page 29 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 9: Simulation results of the closed-loop scheduling scheme under dry weather scenario with varying water price. Table 5: Simulation results under dry weather scenario with varying water price. Yield deciency Irrigation implemented by controller (cm3 /cm2 ) Irrigation prescribed by scheduler (cm3 /cm2 )

Closed-loop 14.9129 1163.6 1236.8

Open-loop 14.1672 1228.9 1167.7

moisture set point before day 60 is around 0.3 cm/cm. This creates a water storage which later is used when the water price increases. It can be observed that the actual irrigation during day 59 drops to zero since the soil moisture set point gets reduced. In this case, the closed-loop scheme leads to a 5.3% water saving but at the cost of 5.2% crop yield reduction.

5

Conclusions and future work

In this work, we investigated closed-loop scheduling together with closed-loop control for agriculture irrigation. In the proposed design, LPV models were identied and used in the design of the scheduler and the controller. The designed hierarchical scheduling and control 29

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 10: Simulation results of the open-loop scheduling scheme under dry weather scenario with varying water price. system is applied to a eld with maize. In the simulations, a couple of dierent weather scenarios were considered. Within each scenario, a few dierent cases were studied including dierent weights on water conservation and varying water prices. From the simulation results, we can observe that the gain (either in water conservation or crop yield) in closed-loop scheduling in most of the cases comes with a cost (either lower crop yield or increased water consumption). The results in the simulations considered do not show a very clear benet of using closed-loop scheduling. One possible reason for the results is that the crop maize considered in this work requires signicant water in the growing season. It may have a much stronger impact on the irrigation pattern than the impacts of uncertainty in weather and model. This may have caused the irrigation to be high in most of the growing season as can be seen in the simulation results and thus makes the dierence between closed-loop scheduling and open-loop scheduling not obvious. Moreover, in this work, the spatial heterogeneity is not considered and the eld only includes maize. When spatial heterogeneity and more

30

ACS Paragon Plus Environment

Page 30 of 35

Page 31 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

types of crops are considered, the benets of a closed-loop scheduling may be more obvious. Further, in this work, we have only considered the soil moisture as a controlled variable, when other manipulated inputs (e.g., salinity) are considered, closed-loop scheduling may bring more benets. In our future work, we will consider dierent crops other than maize. Spatial heterogeneity and other manipulated inputs will also be considered.

References (1) Guan, D.; Hubacek, K. Assessment of regional trade and virtual water ows in China.

Ecological Economics

2007,

61, 159170.

(2) Mubako, S.; Lahiri, S.; Lant, C. Input-output analysis of virtual water transfers: Case study of California and Illinois. Ecological Economics

2013,

93, 230238.

(3) Goodchild, M.; Kühn, K.; Jenkins, M.; Burek, K.; Button, A. A method for precision closed-loop irrigation using a modied PID control algorithm. Sensors & Transducers 2015,

188, 6168.

(4) Kim, Y.; Evans, R.; Iversen, W. Evaluation of closed-loop site-specic irrigation with wireless sensor network. Journal of Irrigation and Drainage Engineering

2009,

135,

2531. (5) Park, Y.; Shamma, J. S.; Harmon, T. C. A Receding Horizon Control algorithm for adaptive management of soil moisture and chemical levels during irrigation. Environ-

mental Modelling & Software

2009,

24, 11121121.

(6) McCarthy, A. C.; Hancock, N. H.; Raine, S. R. Simulation of irrigation control strategies for cotton using Model Predictive Control within the VARIwise simulation framework.

Computers and Electronics in Agriculture

2014,

101, 135147.

(7) Delgoda, D.; Malano, H.; Saleem, S. K.; Halgamuge, M. N. Irrigation control based on 31

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 35

model predictive control (MPC): Formulation of theory and validation using weather forecast data and AQUACROP model. Environmental Modelling & Software

2016,

78,

4053. (8) Mao, Y.; Liu, S.; Nahar, J.; Liu, J.; Ding, F. Soil moisture regulation of agrohydrological systems using zone model predictive control. Computers and Electronics

in Agriculture

2018,

154, 239247.

(9) Thorp, K. R.; Hunsaker, D. J.; Bronson, K. F.; Andrade-Sanchez, P.; Barnes, E. M. Cotton irrigation scheduling using a crop growth model and FAO-56 methods: Field and simulation studies. Transactions of the ASABE

2017,

60, 20232039.

(10) Cahoon, J.; Ferguson, J.; Edwards, D.; Tacker, P. A microcomputer-based irrigation scheduler for the humid mid-south region. Applied Engineering in Agriculture

1990,

6,

289295. (11) Jensen, M. E.; Robb, D. C. N.; Franzoy, C. E. Scheduling irrigations using climatecrop-soil data. Journal of the Irrigation and Drainage Division

1970,

96, 2538.

(12) Bras, R. L.; Cordova, J. R. Intraseasonal water allocation in decit irrigation. Water

Resources Research

1981,

17, 866874.

(13) Wardlaw, R.; Barnes, J. Optimal allocation of irrigation water supplies in real time.

Journal of Irrigation and Drainage Engineering

1999,

125, 345354.

(14) Hassan-Esfahani, L.; Torres-Rua, A.; McKee, M. Assessment of optimal irrigation water allocation for pressurized irrigation system using water balance approach, learning machines, and remotely sensed data. Agricultural Water Management

2015,

153, 42

50. (15) Engell, S.; Harjunkoski, I. Optimal operation: Scheduling, advanced control and their integration. Computers & Chemical Engineering 32

2012,

ACS Paragon Plus Environment

47, 121133.

Page 33 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(16) Touretzky, C. R.; Baldea, M. Integrating scheduling and control for economic MPC of buildings with energy storage. Journal of Process Control

2014,

24, 12921300.

(17) Charitopoulos, V. M.; Papageorgiou, L. G.; Dua, V. Closed-loop integration of planning, scheduling and multi-parametric nonlinear control. Computers & Chemical Engineering 2018,

https://doi.org/10.1016/j.compchemeng.2018.06.021.

(18) Zhuge, J.; Ierapetritou, M. G. Integration of scheduling and control for batch processes using multi-parametric model predictive control. AIChE Journal

2014,

60, 31693183.

(19) Pistikopoulos, E. N.; Diangelakis, N. A. Towards the integration of process design, control and scheduling: Are we getting closer? Computers & Chemical Engineering 2016,

91, 8592.

(20) Richards, L. A. Capillary conduction of liquids through porous mediums. Journal of

Applied Physics

1931,

1(5), 318333.

(21) Mualem, Y. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resources Research

1976,

12, 513522.

(22) Schaap, M. G.; Van Genuchten, M. T. A modied Mualemvan Genuchten formulation for improved description of the hydraulic conductivity near saturation. Vadose Zone

Journal

2006,

5, 2734.

(23) Genuchten, M. T. V. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Science Society of America Journal

1980,

44, 892898.

(24) Van Dam, J. C.; Groenendijk, P.; Hendriks, R. F. A.; Jacobs, C. M. J. SWAP version 3.2: Theory description and user manual. Wageningen UR, Alterra: Postbus 47, 6700 AA Wageningen, The Netherlands, 2008. (25) Van Dam, J. C.; Feddes, R. A. Numerical simulation of inltration, evaporation and

33

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 35

shallow groundwater levels with the Richards equation. Journal of Hydrology

2000,

233, 7285. (26) Nahar, J.; Liu, J.; Shah, S. L. Parameter and state estimation of an agro-hydrological system based on system observability analysis. Computers and Chemical Engineering 2019,

121, 450464.

(27) Allen, R. G.; Pereira, L. S.; Raes, D.; Smith, M. FAO Irrigation and drainage paper No. 56. Rome: Food and Agriculture Organization of the United Nations

1998,

56, e156.

(28) Feddes, R. A. Simulation of eld water use and crop yield ; Pudoc, 1982. (29) Jin, X.; Huang, B.; Shook, D. S. Multiple model LPV approach to nonlinear process identication with EM algorithm. Journal of Process Control

2011,

21, 182193.

(30) Bolea, Y.; Puig, V.; Blesa, J. Linear parameter varying modeling and identication for real-time control of open-ow irrigation canals. Environmental Modelling & Software 2014,

53, 8797.

(31) Yang, X.; Huang, B.; Gao, H. A direct maximum likelihood optimization approach to identication of LPV time-delay systems. Journal of the Franklin Institute

2016,

353,

18621881. (32) Xu, Z.; Zhao, J.; Qian, J.; Zhu, Y. Nonlinear MPC using an identied LPV model.

Industrial & Engineering Chemistry Research

2009,

48, 30433051.

(33) Huang, J.; Ji, G.; Zhu, Y.; van den Bosch, P. Identication of multi-model LPV models with two scheduling variables. Journal of Process Control

2012,

22, 11981208.

(34) Steduto, P.; Hsiao, T. C.; Fereres, E.; Raes, D. FAO irrigation and drainage paper No. 66 “Crop yield response to water”. FAOFood and Agriculture Organization of

the United Nations, Rome

2012,

34

ACS Paragon Plus Environment

Page 35 of 35 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Table of Contents Graphic

35

ACS Paragon Plus Environment