# A mixed integer programming reformulation of the mixed fruit-vegetable crop allocation problem

Sara Maqrot, Simon de Givry, Gauthier Quesnel, Marc Tchamitchian

► To cite this version:

Sara Maqrot, Simon de Givry, Gauthier Quesnel, Marc Tchamitchian. A mixed integer programming reformulation of the mixed fruit-vegetable crop allocation problem. 30th International Conference on Industrial, Engineering and Other Applications of Applied Intelligent Systems, IEA/AIE 2017, Jun 2017, Arras, France. 10.1007/978-3-319-60045-1_26 . hal-01605661

HAL Id: hal-01605661

https://hal.science/hal-01605661

Submitted on 2 Jun 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# A Mixed Integer Programming Reformulation of the Mixed Fruit-Vegetable Crop Allocation Problem

Sara Maqrot<fnref n="1" />, Simon de Givry<fnref n="1" />, Gauthier Quesnel<fnref n="1" />, and Marc Tchamitchian<fnref n="2" />

**Abstract.** Mixed fruit-vegetable cropping systems are a promising way of ensuring environmentally sustainable agricultural production systems in response to the challenge of being able to fulfill local market requirements. Indeed, they combine productions and they make a better use of biodiversity. These agroforestry systems are based on a complex set of interactions modifying the utilization of light, water and nutrients. Thus, designing such a system must optimize the use of these resources: by maximizing positive interactions (facilitation) and minimizing negative ones (competition). To attain these objectives, the system's design has to include the spatial and temporal dimensions, taking into account the evolution of above- and belowground interactions over a time horizon. For that, we define the Mixed Fruit-Vegetable Crop Allocation Problem (MFVCAP) using a discrete representation of the land and the interactions between vegetable crops and fruit trees. First, we give a direct formulation as a binary quadratic program (BQP). Then we reformulate the problem using a Benders decomposition approach. The master problem has 0/1 binary variables and deals with tree positioning. The sub-problem deals with crop quantities. The BQP objective function becomes linear in the continuous subproblem by exploiting the fact that it depends only on the quantity of crops assigned to land units having shade, root, or nothing. This problem decomposition allows us to reformulate the MFVCAP into a Mixed Integer linear Program (MIP). The detailed spatial-temporal crop allocation plan is easy to obtain after solving the MIP. Experimental results show the efficiency of our approach compared to a direct solving of the original BQP formulation.

**Keywords:** agroecology, spatial and temporal crop allocation problem, binary quadratic programming, mixed integer programming, Benders decomposition

## 1 Introduction

Agroforestry systems are one of the sustainable approaches that have received considerable research attention over the past with a view to ensure high productions, ecosystem services and environmental benefits [11]. These systems combine two principal land-use sciences: agriculture and forestry. However, to our knowledge, none of the several studies conducted on agroforestry systems has specifically examined the combination of vegetable crops with fruit trees using a modeling approach. Accordingly, the aim of the present study is to design mixed fruit-vegetable cropping systems, represented as a spatial-temporal crop allocation problem. Nonetheless, unlike existing studies in which

<footnote n="1">¹UR 875 MIAT, Université de Toulouse, INRA, Castanet-Tolosan, France</footnote>

<footnote n="2">²UR 767 Ecodéveloppement, INRA, Avignon, France</footnote>

<footnote>`{firstname.lastname}@inra.fr`</footnote>
allocation concerns only annual crops (see Section 2), our approach allocates both annual vegetable crops and perennial fruit trees on the same land while optimizing above- and belowground interactions resulting from this combination [4].

To assess the validity of different modeling choices to design mixed fruit-vegetable cropping systems, we built a first prototype using a Binary Quadratic Programming formulation of the Mixed Fruit-Vegetable Crop Allocation Problem [15]. This model allocates fruit trees and vegetable crops while optimizing interactions between them and minimizing the dispersion of vegetable crops. Such crops must be allocated so that to change their positions between two successive periods (i.e., crop rotation). The aim of this preliminary work was to examine the ability of state-of-the-art exact solver IBM ILOG CPLEX in solving MFVCAP in order to support farmers in their crop allocation strategies. Depending on various modeling simplifications (without crop rotation nor crop dispersion), CPLEX was able to solve in about one hour CPU time a small piece of land divided into 10 × 10 land units, each one contains either a tree or a crop (among 6 possible ones including bare soil), over a time horizon divided into three periods of four seasons (except the first period with only one season starting in autumn).

We further improved the results by exploiting a Benders decomposition of the problem. The master problem deals with tree positioning and is a 0/1 linear program. The BQP objective function becomes linear in the continuous subproblem by exploiting the fact that it depends on crop quantities assigned to land units having shade, root or nothing. The final spatial-temporal crop allocation is obtained in post-processing. This decomposition allows us to reformulate the MFVCAP into a Mixed Integer linear Program.

Section 2 recaps related work. Section 3 gives an informal definition of MFVCAP followed by its BQP formulation. Section 4 describes our main contribution by reformulating MFVCAP into MIP. Section 5 presents computational experiments for various scenarios and land sizes. We conclude in Section 6.

## 2 Related work

Every year, farmers have to allocate their lands to different crops with respect to physical (water use) and chemical (fertilizers, pesticides, etc.) soil properties. To support farmers in these complex decisions, numerous studies [7] have been conducted on cropping plan and crop rotation decisions since they have a major impact on crop productions. Mathematical programming is widely used in this area [9]. For example, Linear Programming (LP) flow models for agricultural planning are proposed in [8, 10, 6], to find an optimal crop rotation that maximizes profits for a given selection of crops on a given piece of land. Compared to these LP flow models, [2] proposes a MIP flow model that minimizes the surface area needed to cover crop demands (aspect of sustainability) that vary over time rather than maximizing incomes. This problem was proven to be strongly NP-hard in [3] and was reformulated as a 0/1 LP compact formulation based on crop-sequence graphs, then as an extended formulation with a polynomial-time pricing problem and a Branch-and-Price-and-Cut algorithm with adapted branching rules and cutting planes. Many other papers propose a column generation approach [13] to model cropping plan and cropping rotations. In this sense, column generation with
Dantzig-Wolfe decomposition was applied in [18] to maximize the land use subject to neighborhood and succession restrictions for crops of the same botanic family. Column generation was also used in [17], but for a more complex MIP: a crop rotation problem with lands divided into plots where the continuous variables represent the surface areas assigned to a given rotation.

Artificial intelligence methods were also employed for solving cropping plan and
crop rotation problems. In this direction, [19] used the classical Constraint Satisfaction
Problem (CSP) formalism to develop whole-farm crop rotation plans on a specific farm,
whereas [1] proposed a Weighted Constraint Satisfaction Problem (WCSP) formalism
to solve the spatial-temporal crop allocation planning problem for a medium-size virtual
farm. Agent-based simulations are also used in this area. In this sense, [20] presents an
application of a multi-agent architecture based on the BDI (Belief-Desire-Intention)
paradigm and the belief theory to simulate cropping plan decision-making.

Indeed, to the best of our knowledge, none of the several studies conducted on
cropping plan has specifically examined the allocation of annual crops (vegetable) with
perennial ones (trees) using a modeling approach. Hence the originality of our model.

3 Mixed Fruit-Vegetable Crop Allocation Problem as BQP

MFVCAP is a spatial-temporal crop allocation problem in which crops are assigned to a piece of land over a time horizon. We discretize the land as a square of $l \times l$ cells; every cell represents a unit land area allocated to a fruit tree, a vegetable crop or a bare soil. In this study, we consider only one type of fruit tree, *apple tree*, which we plant at the beginning of the time horizon. For vegetable crops, we chose a selection that would allow for variable planting dates along the year and would need different cropping durations: *lettuce*, *tomato*, *onion*, *melon* and *carrot*. We added a *green manure* as it is a required practice for the restoration of soil fertility. Indeed, to represent these crops as time passes, we consider a time horizon divided into three periods of four seasons according to tree growth stages: period P1 corresponds to young trees (sprouts), P2 to intermediate growing trees not yet producing fruits (saplings), and P3 to mature trees. Each tree growing period represents a one-year rotation, except P1 which is represented by only one season for the fact that there are no vegetable crops at this period, only green manure to fertilize the soil. Therefore, we have nine time intervals, and thus nine grids of $l \times l$ land units.

We consider classical interactions observed generally in agroforestry systems, due to root extension dynamics [23, 22] as belowground interactions, and microclimate modifications and crop sensitivity to shade (based on solar radiation interception simulation of an apple tree) as above interactions. We model the evolution of the tree root system in three periods as shown in Fig. 1. The shade of a tree takes place only in the spring and summer seasons and we model it as a fixed area orientated to the North of the tree (see Fig. 1).

Moreover, we define operation rules related to the organization of crops in the field: a minimal space between trees allowing them a non-conflict growth and a better light distribution (see Fig. 2), and a minimum (and maximum) number of land units to be allocated to each crop, in order to diversify the food at a given season (see Fig. 3).
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Surface view of root (a) and shade (b) extensions. Blue cells host trees, and checked cells host roots/shade in periods P1, P2, and P3.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Minimum space between trees. The blue cell hosts a tree, and checked cells are forbidden for other trees.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Minimum and maximum crop quantities as a percentage of the total grid area.</figcaption></figure>

We consider the following binary variables in order to formulate the MFVCAP as a Binary Quadratic Program (BQP):

* $crop_{x,y}^{t,c} = 1$ iff land unit of coordinates $x, y$ has crop $c$ at time $t$,

* $tree_{x,y} = 1$ iff land unit of coordinates $x, y$ has a tree,

* $shade_{x,y} = 1$ iff land unit of coordinates $x, y$ has shade,

* $root_{x,y}^p = 1$ iff land unit $x, y$ has roots at period $p$.

Let $L = \{1, \dots, l\}$ be a set of horizontal/vertical land unit positions of a square piece of land discretized as a $l \times l$ grid with top-left (i.e., North-West) corner having coordinates $(x = 1, y = 1)$, $\zeta$ a set of crops, $T$ a set of time steps (corresponding to seasons, starting in autumn), and $P$ a set of periods of tree growth stages, with $T_p$ the set of time steps at each period $p \in P$. In our experiments, we vary $l$ and use a fixed $\zeta$, $\zeta = \{\text{onion}, \text{melon}, \text{lettuce}, \text{carrot}, \text{tomato}, \text{green manure}, \text{bare soil}\}$, $T = \{1, \dots, 9\}$, $P = \{1, 2, 3\}$, $T_1 = \{1\} \sim \{\text{autumn}\}$, $T_2 = \{2, 3, 4, 5\}$, $T_3 = \{6, 7, 8, 9\} \sim \{\text{winter}, \text{spring}, \text{summer}, \text{autumn}\}$.

Concerning the objective function, we minimize the following quadratic function (minimization instead of maximization as in [15] and without crop dispersion):

$$ (1) \quad \min \sum_{p \in P, t \in T_p, c \in \zeta, x \in L, y \in L} [A^{t,c} \times shade_{x,y} + B^{t,c} \times root_{x,y}^p + C^{t,c}] \times crop_{x,y}^{t,c} $$

with coefficients $A^{t,c}$ related to the degree of crop sensitivity to shade at a particular season (without shade at the first period, i.e., $A^{1,c} = 0$), $B^{t,c}$ related to interactions generated by root systems of trees and vegetable crops at a particular season (competition or sharing for water), and $C^{t,c}$ related to the crop selection whatever presence or
absence of shade and tree roots at a given land unit. Specific values of these coefficients have been previously defined by agronomists in [21, 15].

To complete the mathematical model, we define the following linear constraints:

A crop or a tree in each land unit

$$ (2) \qquad \text{tree}_{x,y} + \sum_{c \in C} \text{crop}_{x,y}^{t,c} = 1 \quad (\forall t \in T, \forall x, y \in L^2) $$

Minimal space between trees (see Fig.2)

$$ (3) \qquad \text{tree}_{x,y} + \text{tree}_{x+1,y} \le 1 \quad (\forall x \in L - \{l\}, \forall y \in L) $$

$$ (4) \qquad \text{tree}_{x,y} + \text{tree}_{x,y+1} \le 1 \quad (\forall x \in L, \forall y \in L - \{l\}) $$

No trees on East and West borders

$$ (5) \qquad \text{tree}_{x,l} = 0 \quad (\forall x \in L) $$

$$ (6) \qquad \text{tree}_{1,y} = 0 \quad (\forall y \in L) $$

Definition of shade in the North, East, and West of a tree (Fig.1)

$$ (7) \qquad \sum_{\substack{i \in \{\max(1-x,-1), \dots, \min(l-x,1)\} \\ j \in \{0, \dots, \min(l-y,1)\}}} \text{tree}_{x+i,y+j} - 6\text{shade}_{x,y} \le 0 \quad (\forall x, y \in L^2) $$

$$ (8) \qquad \text{shade}_{x,y} - \sum_{\substack{i \in \{\max(1-x,-1), \dots, \min(l-x,1)\} \\ j \in \{0, \dots, \min(l-y,1)\}}} \text{tree}_{x+i,y+j} \le 0 \quad (\forall x, y \in L^2) $$

Evolution of tree roots (Fig.1)

$$ (9) \qquad \text{root}_{x,y}^1 - \text{tree}_{x,y} = 0 \quad (\forall x, y \in L^2) $$

$$ (10) \qquad \sum_{\substack{i \in \{\max(1-x,-1), \dots, \min(l-x,1)\} \\ j \in \{\max(1-y,-1), \dots, \min(l-y,1)\}}} \text{root}_{x+i,y+j}^{p-1} - 9\text{root}_{x,y}^{p} \le 0 \quad (\forall x, y \in L^2, \forall p \in \{2, 3\}) $$

$$ (11) \qquad \text{root}_{x,y}^{p} - \sum_{\substack{i \in \{\max(1-x,-1), \dots, \min(l-x,1)\} \\ j \in \{\max(1-y,-1), \dots, \min(l-y,1)\}}} \text{root}_{x+i,y+j}^{p-1} \le 0 \quad (\forall x, y \in L^2, \forall p \in \{2, 3\}) $$

Vertical symmetry breaking constraint

$$ (12) \qquad \sum_{x \in \{(l-l/2)+1), \dots, l\}, y \in L} \text{tree}_{x,y} - \sum_{x \in \{1, \dots, \lfloor l/2 \rfloor\}, y \in L} \text{tree}_{x,y} \le 0 $$

Minimum and maximum crop quantities (Fig.3)

$$ (13) \qquad \minBalance^{t,c} \le \sum_{x,y} \text{crop}_{x,y}^{t,c} \le \maxBalance^{t,c} \quad (\forall t \in T, \forall c \in \zeta) $$

Crops growing during two consecutive seasons

$$ (14) \qquad \text{crop}_{x,y}^{t,c} - \text{crop}_{x,y}^{t+1,c} = 0 \quad (\forall p \in P, t=f_c(T_p), \forall x,y \in L^2, \forall c \in \{\text{onion}, \text{tomato}, \text{carrot}, \text{manure}\}) $$

where $f_c(T_p)$ returns if available the (unique) time step corresponding to the crop c planting season in $T_p$ at period p for crops growing during two consecutive seasons (manure is planted in autumn, e.g., $f_{manure}(T_1) = 1$, onion and tomato in spring, $f_{onion \lor tomato}(T_2) = 3$, and carrot in summer, $f_{carrot}(T_3) = 8$). Compared to the original formulation in [15], we simplify the problem by removing crop rotation constraints and crop dispersion, and we add a symmetry breaking constraint (Eq. (12)).
## 4 Reformulation as a MIP by using a Benders Decomposition with Crop Quantity Variables

We apply a Benders decomposition approach [5, 16] in order to separate the MFV-CAP into two parts: (i) tree positioning and (ii) crop production optimizing above- and belowground interactions over time. Problem (i) will correspond to the reduced master problem (MP) in the Benders approach and problem (ii) to the subproblem (SP). The master problem keeps the binary variables $tree_{x,y}$, $shade_{x,y}$, and $root_{x,y}^p$, and the associated constraints (3)-(12). The remaining constraints (13)-(14) belong to the subproblem. The master problem communicates with the subproblem by introducing new integer variables (and constraints) corresponding to the number of trees, i.e., $trees = \sum_{x,y \in L^2} tree_{x,y}$, the number of shade land units for crops, $shades = \sum_{x,y \in L^2} shade_{x,y} - trees$, and the number of root land units for crops at each period, $roots_p = \sum_{x,y \in L^2} root_{x,y}^p - trees$. Recall that there are no shade nor roots land units available for crops at the first period (i.e., $roots_1 = 0$).

These extra variables $trees$, $shades$, $roots_p$ will be taken as constants in the subproblem. Instead of having $crop_{x,y}^{t,c}$ binary variables here, we replace them by crop quantities for each type of land unit configuration. Because at any land unit, shade at periods 2 and 3 implies roots (see Eq. (7)-(11) and Fig. 1), we distinguish three possible configurations: land units having shade and roots ($sr$), land units having roots only ($\bar{sr}$), and land units having no shade and no roots ($\bar{sr}$). The corresponding crop quantity non-negative real variables are respectively: $q_{crop}^{t,c}_{sr}$, $q_{crop}^{t,c}_{sr}$, and $q_{crop}^{t,c}_{sr}$, which defines a partition of all the land units excluding trees (the sum of which equals $q_{crop}^{t,c}$, see Eq. (16)-(17)). Adding Eq. (18)-(19) ensures a feasible allocation w.r.t. land units having shade or roots and no trees. The subproblem has now a linear objective function,

$$ (15) \quad \min \sum_{t \in T, c \in \zeta} (A^{t,c} + B^{t,c}) q_{crop}^{t,c}_{sr} + B^{t,c} q_{crop}^{t,c}_{sr} + C^{t,c} q_{crop}^{t,c} $$

such that,

$$ (16) \quad q_{crop}^{t,c}_{sr} + q_{crop}^{t,c}_{sr} + q_{crop}^{t,c}_{sr} - q_{crop}^{t,c} = 0 \quad (\forall t \in T, c \in \zeta) $$

$$ (17) \quad \sum_{c \in \zeta} q_{crop}^{t,c} = l^2 - \text{trees} \quad (\forall t \in T) $$

$$ (18) \quad \sum_{c \in \zeta} q_{crop}^{t,c}_{sr} = \text{shades} \quad (\forall t \in T \setminus \{1\}) $$

$$ (19) \quad \sum_{c \in \zeta} q_{crop}^{t,c}_{sr} + q_{crop}^{t,c}_{sr} = \text{roots}_p \quad (\forall p \in P, t \in T_p) $$

$$ (20) \minBalance^{t,c} \le q_{crop}^{t,c} \le maxBalance^{t,c} \quad (\forall t \in T, c \in \zeta) $$

$$ (21) \quad q_{crop}^{t,c} - q_{crop}^{t+1,c} = 0 \quad (\forall p \in P, t=f_c(T_p), \forall c \in \{\text{onion, tomato, carrot, manure}\}) $$

$$ (22) \quad q_{crop}^{t,c}_{sr} - q_{crop}^{t+1,c}_{sr} = 0 \quad (\forall p \in P \setminus \{1\}, t=f_c(T_p), \forall c \in \{\text{onion, tomato, carrot, manure}\}) $$

$$ (23) \quad q_{crop}^{t,c}_{sr} - q_{crop}^{t+1,c}_{sr} = 0 \quad (\forall p \in P \setminus \{1\}, t=f_c(T_p), \forall c \in \{\text{onion, tomato, carrot}\}) $$

$$ (24) \quad q_{crop}^{t, \text{manure}}_{sr} - q_{crop}^{t+1, \text{manure}}_{sr} \le 0 \quad (\forall p \in P, t=f_{\text{manure}}(T_p)) $$
Because variables $q_{crop\bar{s}\bar{r}}^{t,c}$ appear only in Eq. 16, such as slack variables ($q_{crop\bar{s}\bar{r}}^{t,c} = q_{crop}^{t,c} - q_{crop\bar{s}\bar{r}}^{t,c} - q_{crop\bar{s}\bar{r}}^{t,c}$), we can remove them, replacing Eq. 16 by

$$ (25) \qquad q_{crop\bar{s}\bar{r}}^{t,c} + q_{crop\bar{s}\bar{r}}^{t,c} - q_{crop}^{t,c} \le 0 \quad (\forall t \in T, c \in \zeta) $$

Notice that the number of land units with shade or roots does not change inside a given period, thus the quantities of crops growing on two consecutive seasons within the same period (except the first) that are allocated to shade or roots must remain the same (Eq. (22)-(23)). Conversely, two-season duration crops like manure that overlap two periods can be allocated to more land units with roots (and no shade) in their second period than in their first period (Eq. (24)). It makes a supplementary partition between land units having no roots at period 2 and roots at period 3. This will be taken into account when reconstructing a complete allocation plan (see Alg. 1, line 5 and below).

In conclusion, the subproblem deals with crop production over time but not with crop allocation on every land units, resulting in a linear objective function instead of the original quadratic one. It is easy to verify that Eq. (15) is equivalent to Eq. (1).

Without constraints (21)-(24) dealing with two-season duration crops, the subproblem can easily be cast as a set of independent minimum-cost flow problems for every season with maximum flow equal to $l^2 - \text{trees}$, and with demands (minBalance) and capacities (maxBalance, $l^2 - \text{trees}$, shades, roots$_p$) given as integer constants. See an example in Fig. 4. It is therefore totally unimodular, with integer optimal solutions. Two-season crops add a chain of equality/inequality binary constraints between crop quantities of consecutive seasons that preserves the integrality property.

To ensure that the subproblem is always feasible, the master problem has to assign trees so as to keep enough land units to be allocated to crops (in the subproblem) satisfying the crop production demands (Eq. 20). This is expressed by the following equations:

$$ (26) \qquad \text{trees} \le l^2 - \sum_{c \in \zeta} \minBalance^{t,c} \quad (\forall t \in T) $$

$$ (27) \qquad \sum_{c \in \zeta} \minBalance^{t,c} \le l^2 \le \sum_{c \in \zeta} \maxBalance^{t,c} \quad (\forall t \in T) $$

To summarize, we have decomposed the original problem in a master and subproblem by adding extra channeling variables (trees, shades, roots$_p$) also resulting in a linear objective function. Because these variables occur as linear terms in MP and SP (Eq.(17)-(19)), we can reformulate the initial problem by the conjunction of the master (MP) and subproblem (SP), i.e., Eq.(3)-(12);(17)-(26) with linear objective Eq.(15), resulting in a Mixed Integer linear Program. This MIP can be solved directly or using a Benders decomposition approach with integer binary variables in the master problem and continuous variables in the subproblem as shown above and implemented in the last version of IBM ILOG cplex solver (v12.7).

## Building a complete allocation plan

It remains to build from the MIP solution a complete crop allocation plan with spatial-temporal crop variables $crop_{x,y}^{t,c}$ correctly assigned. The MIP solution gives the positions of trees, shades, and roots at every period. It also gives the crop quantities depending on land unit configurations. In principle, we could allocate randomly $q_{crop\bar{s}\bar{r}}^{t,c}$
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Crop production in summer at period 2, regardless of two-season crop constraints, for a 10x10 grid as a minimum-cost flow problem with non-zero edge costs and in parentheses, demands and capacities. We define $rt_2 = \text{roots}_2 - \text{shades}$ and $nt_2 = 100 - \text{trees} - \text{roots}_2$.</figcaption></figure>

(respectively $q_{\text{crop}\bar{s}_r}^{t,c}$, $q_{\text{crop}\bar{s}_r}^{t,c}$) land units of crop c among free land units having shade and roots (resp. roots only, no shade nor roots) at each time step t and for every crop c. This is not as simple for two-season duration crops. For these crops c, every assignment of a particular land unit x, y at planting time t must be assigned to the same crop at time t + 1, i.e., $\text{crop}_{x,y}^{t,c} = \text{crop}_{x,y}^{t+1,c}$. Moreover, for manure, which spans over two consecutive periods at t and t + 1, a specific assignment procedure must be done. In case we have $q_{\text{crop}\bar{s}_r}^{t,\text{manure}} < q_{\text{crop}\bar{s}_r}^{t+1,\text{manure}}$, then we need to find $q_{\text{crop}\bar{s}_r}^{t+1,\text{manure}} - q_{\text{crop}\bar{s}_r}^{t,\text{manure}}$ land units without roots at time t and with roots at time t + 1 when allocating $q_{\text{crop}\bar{s}_r}^{t,\text{manure}}$ land units to manure at time t in order to preserve a constant overall quantity of manure at times t and t + 1. The remaining land units $q_{\text{crop}\bar{s}_r}^{t,\text{manure}} - (q_{\text{crop}\bar{s}_r}^{t+1,\text{manure}} - q_{\text{crop}\bar{s}_r}^{t,\text{manure}})$ are chosen arbitrary among land units without roots. Because the two-season crops impose a chain of binary equality/inequality constraints (Eq. (21)-(24)) spanning from spring to winter, we need to build the allocation plan following a temporal order. We use a decreasing order starting in autumn of period 3 in order to satisfy the more complex manure constraint (Eq. (24)) before the other two-season crop equality constraints (see Fig. 5 and Alg. 1 line 3). By doing so, our greedy algorithm, Algorithm 1, building the complete allocation plan, is correct and does not require sophisticated constraint propagation techniques.

Instead of allocating crops to land units randomly, we look for the nearest<fnref n="1" /> free land unit to the top-left corner of the piece of land when doing the assignment of a given crop quantity to land units. It has the effect to reduce horizontal crop dispersion. Furthermore, we allocate crops following a reverse order between periods 3 and 2 (see Fig. 5 and Alg. 1 line 4). It tends to favor crop rotation between two successive periods.

# 5 Experimental Results

To assess the impact of above and belowground interactions, we performed simulations on three scenarios: scenario Above gives a significant importance to above interactions by

<footnote n="1">¹ Minimizing the y-axis distance first, then the x-axis. See Alg. 1 lines 1-2.</footnote>
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Rules to build a complete allocation plan. Reverse order of time steps starting from autumn of period 3 to satisfy the equality/inequality constraints of two-season crops. Reverse order of crops between period 2 and period 3 to favor crop rotation.</figcaption></figure>

multiplying the effect values related to shade by 10, similarly scenario Below attaches a great importance to roots by multiplying their effect values by 10 and the last scenario Equilibrate ascribes equal importance to the impact of both above and belowground interactions.

We solved these scenarios using IBM ILOG cplex v12.7 with its new fully-automatic Benders decomposition strategy for the MIP formulation and with default options for the BQP formulation, except saving the branch-and-bound tree to disk and using 4 cores of an Intel Xeon CPU E5-2680 at 2.5 GHz with 256GB running Linux Debian 8.7. By comparing the real time taken by BQP formulation and MIP formulation to find an optimal solution and prove its optimality (see Fig. 6), we observed that the MIP is several orders of magnitude faster than the BQP in all scenarios.

<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Solving time in seconds taken by BQP and MIP formulations to solve the MFVCAP.</figcaption></figure>

To analyze the obtained results, Figure 7 shows as an example the crop allocation plan of a piece of land of 15x15 land units for scenarios Equilibrate, Below, and Above. We note that these solutions satisfy the constraints on minimal space between trees
<caption>**Algorithm 1:** Build a complete crop allocation plan.</caption>

Input: $tree_{x,y}, shade_{x,y}, root_{x,y}^p, q_{crop_{sr}}^{t,c}, q_{crop_{sr}}^{t,c}, q_{crop_{sr}}^{t,c}, q_{crop_{sr}}^{t,c}$
Output: $crop_{x,y}^{t,c}$
Procedure Allocate ($q_{crop_{type}}^{t,c}, t, c, p, type, p', type'$)
    if ($t < \max(T) \land twoseason(c) \land f_c(T_p) = t$) then return;
    for $y \leftarrow 1$ to $l$ do
        for $x \leftarrow 1$ to $l$ do
            if ($q_{crop_{type}}^{t,c} = 0$) then return;
            if ($free_{x,y}^t \land Type(p, x, y) = type \land Type(p', x, y) = type'$) then
                if ($free_{x,y}^t \leftarrow false$;
                    $q_{crop_{type}}^{t,c} \leftarrow q_{crop_{type}}^{t,c} - 1$;
                    $crop_{x,y}^{t,c} \leftarrow 1, \forall c' \in \zeta \setminus \{c\}, crop_{x,y}^{t,c'} \leftarrow 0$;
                    if ($twoseason(c) \land f_c(T_p) < t$) then
                        if ($free_{x,y}^{t-1} \leftarrow false$;
                            $crop_{x,y}^{t-1,c} \leftarrow 1, \forall c' \in \zeta \setminus \{c\}, crop_{x,y}^{t-1,c'} \leftarrow 0$;
                        Special Boolean array $free_{x,y}^t$ to identify which land units remain available for crops;
                        $\forall x \in L, y \in L, t \in T, \text{if } (tree_{x,y} = 1) \text{ then } free_{x,y}^t \leftarrow \text{false else } free_{x,y}^t \leftarrow \text{true}$;
                        for $t \leftarrow \max(T)$ to 2 do
                            Reverse order of crops between periods P2 and P3;
                            if ($t \ge \min(T_3)$) then $\zeta \leftarrow \zeta_3$; $p \leftarrow 3$ else $\zeta \leftarrow \zeta_2$; $p \leftarrow 2$;
                            foreach ($c \in \zeta$) do
                                Allocate ($q_{crop_{sr}}^{t,c}, t, c, p, sr, p, sr$);
                                if ($c = \text{manure} \land t < \max(T)$) then
                                    Allocate ($q_{crop_{sr}}^{t-1,c}, t, c, p, \bar{s}r, p - 1, \bar{s}r$);
                                    Allocate ($q_{crop_{sr}}^{t,c} - q_{crop_{sr}}^{t-1,c}, t, c, p, \bar{s}r, p - 1, \bar{s}r$);
                                else Allocate ($q_{crop_{sr}}^{t,c}, t, c, p, \bar{s}r, p, \bar{s}r$);
                                Allocate ($q_{crop_{sr}}^{t,c}, t, c, p, \bar{s}r, p, \bar{s}r$);

and consecutive crops. The main differences between these scenarios are the number and position of trees and crops. For scenario Equilibrate, we observe a high density of trees, caused by the overall negative effect of interactions. Therefore, to reduce these negative effects while ensuring the minimum allowed production of crops, the solver assigns many land units to trees (not penalized in the objective function) at the expense of vegetable crops. For scenario Below, trees are grouped together as in a forest to have less land units related to belowground interactions and thus minimum production of vegetable crops with negative effects. Consequently, this leads to a clear separation between the vegetable garden and the orchard, with insertion of some vegetable crops between trees. This clear separation results in shorter solving times for the MIP formulation, with a piece of land of 20x20 land units solved in less than 100 seconds. Concerning scenario Above, we notice a spaced disposition of apple trees in sparse groves. One reason is to ensure a sufficient number of land units of shade for the benefit of vegetable crops (e.g., tomato is at maximum allowed production in summer).
<figure><img src="image_8.png" /><figcaption>Fig. 7. Allocation of trees and crops on 15 × 15 land units over two periods for all scenarios.</figcaption></figure>

# 6 Conclusion and Perspectives

We have proposed a MIP reformulation of the Mixed Fruit-Vegetable Crop Allocation Problem, finding optimal 1-fruit-tree/6-vegetable crop allocation plans for 15 × 15 land units and 9 time steps in less than 1 hour for scenario Equilibrate.

Further research remains to be done in order to extend our model, e.g., by taking into account crop rotations in order to avoid soil depletion and to increase pest and disease natural regulation [7]. To be more realistic w.r.t. real land sizes, we have to solve larger instances (going from 20-by-20 to roughly 100-by-100 land units). A possible direction is to try to generalize solutions found at small size to larger sizes by identifying patterns of tree configurations to be repeated. Although it seems easy to generalize for scenario Below, it looks more complex for scenario Equilibrate and Above. Such a generalization may be further improved by a matheuristic procedure combining local search and linear programming as it is done on block layout design problems [12, 14].
References

1. Akplogan, M., De Givry, S., Metivier, J.P., Quesnel, G., Joannon, A., Garcia, F.: Solving the crop allocation problem using hard and soft constraints. RAIRO-Operations Research 47(2), 151–172 (2013)

2. Alfandari, L., Lemalade, J., Nagih, A., Plateau, G.: A MIP flow model for crop-rotation planning in a context of forest sustainable dev. Ann. Operat. R. 190(1), 149–164 (2011)

3. Alfandari, L., Plateau, A., Schepler, X.: A branch-and-price-and-cut approach for sustainable crop rotation planning. European Journal of Operational Research 241(3), 872–879 (2015)

4. Batish, D., Kohli, R., Jose, S., Singh, H.: Ecological basis of agroforestry. CRC Press (2007)

5. Benders, J.F.: Partitioning procedures for solving mixed-variables programming problems. Numerische mathematik 4(1), 238–252 (1962)

6. Detlefsen, N.K., Jensen, A.L.: Modelling optimal crop sequences using network flows. Agri-cultural Systems 94(2), 566–572 (2007)

7. Dury, J., Schaller, N., Garcia, F., Reynaud, A., Bergez, J.E.: Models to support cropping plan and crop rotation decisions. A review. Agronomy for sustainable dev. 32(2), 567–580 (2012)

8. El-Nazer, T., McCarl, B.A.: The choice of crop rotation: A modeling approach and case study. American Journal of Agricultural Economics 68(1), 127–136 (1986)

9. Glen, J.: Mathematical models in farm planning: a survey. Operat. R. 35(5), 641–666 (1987)

10. Haneveld, W.K., Stegeman, A.W.: Crop succession requirements in agricultural production planning. European Journal of Operational Research 166(2), 406–429 (2005)

11. Jose, S., Gillespie, A.R., Pallardy, S.G.: Interspecific interactions in temperate agroforestry. In: New Vistas in Agroforestry, pp. 237–255. Springer (2004)

12. Kulturel-Konak, S., Konak, A.: Linear programming based genetic algorithm for the unequal area facility layout problem. Int. J. of Production Research 51(14), 4302–4324 (2013)

13. Lbbecke, M.E., Desrosiers, J.: Selected topics in column generation. Operations Research 53(6), 1007–1023 (2005)

14. Maniezzo, V., Sttzle, T. (eds.): Matheuristics. IRIDIA, Brussels, Belgium (2016)

15. Maqrot, S., de Givry, S., Quesnel, G., Tchamitchian, M.: Designing mixed fruit-vegetable cropping systems by integer quadratic programming. In: Proc. of iEMSs. Toulouse, France (2016)

16. Rahmaniani, R., Crainic, T.G., Gendreau, M., Rei, W.: The Benders Decomposition Algorithm: A Literature Review. European Journal of Operational Research (2016)

17. dos Santos, L.M.R., Costa, A.M., Arenales, M.N., Santos, R.H.S.: Sustainable vegetable crop supply problem. European Journal of Operational Research 204(3), 639–647 (2010)

18. dos Santos, L.M.R., Michelon, P., Arenales, M.N., Santos, R.H.S.: Crop rotation scheduling with adjacency constraints. Annals of Operations Research 190(1), 165–180 (2011)

19. Stone, N., Buick, R., Roach, J., Scheckler, R., Rupani, R.: The planning problem in agriculture: Farm-level crop rotation planning as an example. AI applications in natural resource management (USA) (1992)

20. Taillandier, P., Therond, O., Gaudou, B.: A new BDI agent architecture based on the belief theory. Application to the modelling of cropping plan decision-making. In: International environmental modelling and software society (iEMSs). Leipzig, Germany (2012)

21. Tchamitchian, M., Godin, E.: Designing mixed horticultural systems. Building Organic Bridges 1, 179–182 (2014)

22. Vercambre, G., Pags, L., Doussan, C., Habib, R.: Architectural analysis and synthesis of the plum tree root system in an orchard using a quantitative modelling approach. Plant and Soil 251(1), 1–11 (2003)

23. Weaver, J.E., Bruner, W.E.: Root Development Of Vegetable Crops Ist Edn. Mcgraw-Ill Book Co.; London (1927)