Novel Route Planning System for Machinery Selection. Case: Slurry Application

: The problem of ﬁnding an optimal solution for the slurry application process is casted as a capacitated vehicle routing problem (CVRP) in which by considering the vehicle’s capacity, it is required to visit all the tracks only once to fully cover the ﬁeld, as well as complying with a speciﬁed targeted application rate. A key objective in this study was to determine an optimized coverage plan in order to minimize the driving distance in the ﬁeld, while at the same time allowing for varying the application rate. The coverage plan includes the optimal sequence of tracks with a speciﬁed application rate for each track. Two algorithms were developed for optimization and simulation of the slurry application cast as capacitated operations. In order to validate the proposed algorithms, a slurry application operation was recorded, and the results of the optimization algorithm were compared with the conventional non-optimized method. The comparison showed that applying the proposed new method reduces the non-working distance by 18.6% and the non-working time by 28.1%.


Introduction
Capacitated field operations are operations that either involve input material flow such as fertilizing, or output material flow such as harvesting. Manure distribution operation requires considerable labor inputs and high-capacity specialized machines and equipment [1]. Furthermore, the odor and nuisance concerns related to the transportation and application of manure as well as the field traffic and soil compaction related to the spreading manure by heavy tankers especially in wet seasons, which impairs crop growth, are among the concerns which require farmers to optimize these operations to reduce the operational time and transportation costs [2,3].
As a solution, coverage planning methods can be used to increase the operational efficiency of the agricultural field operations [4][5][6], cutting the operational costs by improving the path tracking performance of the vehicle and by minimizing the overlapped area in the field [6]. Moreover, vehicle route planning (VRP) can be applied as a tool to reduce the operational time and costs by minimizing the total traveled distance in the field. In agriculture, VRP also has received increased attention and several studies show the application of VRP for planning the in-field operations [7][8][9]. Bochtis et al. (2013) [10] proposed an algorithmic method towards minimizing the total non-working distance traveled by an agricultural machine in the field. This involved expressing the field coverage as the traversal of a weighted graph and showed that the problem of finding an optimal traversal sequence is equivalent to finding the shortest path in the graph. By deploying algorithmically computed optimal sequences, the total non-working distance could be reduced by up to 50%, depending on the operating width and turning radius of the machine, shape of the fields and the type of conventional patterns that compared with the B-patterns. Moreover, it was shown that by exploiting the optimal headland patterns, apart from the reduction of non-working distance (followed by the reduction of the fuel consumption and the non-productive time), a reduction of soil compaction in the headland area will also be achieved [7].
Although the VRP has received considerable attention among computer scientists, operations management specialists, and others researching logistics, solving the problem by an exact algorithm is time demanding and computationally intractable [11]. Therefore, this challenge has led to the need to develop algorithms that can produce near-optimal solutions in a reasonable amount of time [12]. During the years, different solution approaches have been developed including exact algorithms (such as branch-and-bound [13,14] and branch-and-cut [15]), heuristic algorithms (such as the Clarke-Wright savings algorithm [16]), and metaheuristic algorithms (such as simulated annealing [17,18], genetic algorithms [19], tabu search [20], and ant algorithms [21]. Earlier, conventional heuristic algorithms were designed as a response to limited computer processing power. However, in recent years, meta-heuristic algorithms, such as simulated annealing, genetic algorithms, and tabu search, have been developed to perform a more thorough search of solution space. These algorithms are more computationally expensive where the time needed for calculations and the degree of optimality depend on the specifications of the problem and the design of the meta-heuristics. Some coverage planning methods have been presented in the scientific literature [22][23][24]. However, the case of capacitated operations has not been covered extensively. Ali et al. (2009) [25] proposed a practical planning approach for harvesting a field with several capacitated combine harvesters and tractor-trailers.   [26] presented an algorithmic method for the optimization of capacitated field operations using the case of liquid fertilizing. The presented method in this study avoids turnings within the main cropping area of the field and limits all the maneuvers to the headland part. As a result, there is less soil compaction in the inner part of the field in comparison with the previously mentioned approaches.
The objective of this paper is to develop an approach for the optimization of capacitated field operations such as manure application in the field. The target users are farmers and agricultural advisers who are interested in reducing the operational cost and the environmental impacts whilst maximizing the field efficiency. The specific objectives included: • Decompose the fertilizing operation, by considering the operations elements (such as performing the main task) and unproductive elements (such as turnings in the headland part and idle transportation), to determine the operational performance of the machinery. Farmers' current practices are used for validating and benchmarking the proposed algorithm and model.

•
Develop an algorithm and tool to solve this specific field coverage problem with optimized application rates and minimized driving distances for all individual tracks in the field.

•
Develop an approach and tool to help farm managers to select the proper tank volume for the machinery system given specific constraints.

Characteristics of the Slurry Application
In capacitated operations such as slurry application, the depot will be visited several times for refilling to complete the task. The number of tracks covered in each tour (a path that starts and ends at the depot) depends on the capacity of the slurry tank and the application demand of those tracks. The demand of each work-area can be calculated by multiplying the targeted application rate (L/m 2 ) to the total area of that track (m 2 ). The target application rate is the dosage of applied manure which is predetermined to meet crop nutrient requirements. In the manure distribution process, farmers prefer to go to the field with a full tank and return to the depot with the empty tank. Therefore, the summation of applied material in the visited tracks for each tour should be equal to the volume of the slurry tank. Thus, by adjusting the demands of a group of tracks, it is possible to cover them by one load. In this study and based on the agronomical and environmental standards [27], a target application rate was considered for the field and an arbitrary tolerance of ± 30% was considered for adjusting the tracks' demand according to the information derived from a survey from farmers.

Mathematical Formulation
As mentioned before, the field route planning can be cast as a VRP problem which consists of determining a set of routes with minimum traveled distance for a vehicle, starting and ending at a single depot and satisfying the demand of the tracks with the constraints that each track is covered exactly once, and the total demands of the tracks in each route do not exceed the capacity of the vehicle. Bochtis D and Vougioukas S [7], formulated this problem as a weighted graph G = (V, E), where V = NU{0}, is the set of all vertices and N = {1, 2, 3 . . . n} is the set of nodes and E is the set of all edges in the graph. The depot is represented by vertex 0 and the nodes in set N are denoted as the customers. For each edge e ij ∈ E , i j, a non-negative cost C ij is considered as transit cost. Each node i ∈ N is associated with a non-negative demand d i .
The objectives or criteria considered in this study are: (I) each node is visited exactly once, (II) all routes start and end at the depot, (III) for each route, the total demand should be equal to the vehicle capacity VC in order to satisfy the criteria of returning to the depot with an empty tank, (IV) the application rate for each track should fit in the defined range (± 30% target application rate), (V) minimizing the non-working distance traveled by the vehicle.
All the applied variables and parameters in the following formulas were defined in Table 1. The objective function can be defined as follows: Min : The constraints are as follows: The Equations (3) and (4) state that only one of the sibling nodes appears in the solution. Equation (5) specifies that the total summation of demand related to each trip should be equal to the capacity of the vehicle in order to satisfy the criteria to return to the depot with an empty tank. Equation (6) represents that if the vehicle starts covering a track K, it will leave the track at the end of the process. Equation (7) ensures that the application rate for each edge will fall within the stipulated range. Lastly, Equation (8) states the type of decision variable X ij .

Solution Representation
The proposed method denotes a novel precision agriculture method, that combines organic and artificial fertilizers application on agricultural fields while at the same time optimize the in-field transport and subsequent soil compaction. According to this new method, the complete process of fertilization of a field is executed in two steps. Firstly, a slurry tanker is used to distribute the liquid organic manure on the field by adjusting the agronomic required target application rate for each track and minimizing the non-working distance as determined by the proposed method. As a result, an organic fertilizer's application map will be generated to be used in the next step. This next step comprises an artificial fertilizer spreader applying variable rate fertilizer based on the previously generated map to meet the full crop nutrient requirements and at the same time minimize the non-working distance. This novel method will minimize the overall combined in-field traffic and the non-working driving distances by performing all the optimal turnings and maneuvers in the headland area. Moreover, it will reduce the soil compaction in the main cropping area of the field by avoiding turnings in-field and empty transport in the main part of the field. This is achieved by adjusting the fertilizer's application rate to let the machine's operator cover the entire track to reach the headland part. According to the machine's variable rate application capability, the planned actual application rate for each track is going to be determined as a stipulated deviation from the target rate for a field and crop.

Field Representation
The first step in the field area coverage planning process is to generate a geometrical representation of the selected field. One of the common methods in this concept is using the spatial configuration planning approach, which applies the geometric primitives to represent the field [10]. This method consists of three tasks: First, representing the field with simple geometry shapes; second, determining the driving direction in the field and; third, generating the work areas (tracks) based on the driving direction. The output from this method is a set of line segments or polylines which depicts the fieldwork areas and headland passes that can be followed by the machine.

Headland Area Generation
The field headland area, as showed in Figure 1, is generated by inwardly offsetting the field's boundary equal to the multiplication of the operating width of machine (W) times to the number of headlands passes (h). The distance from the field boundaries to the first headland pass is equal to W/2, while the distance between subsequent passes of headland equals to the W. AgriEngineering 2020, 3 FOR PEER REVIEW 5 Figure 1. Headland area generation.

Track Generation and Edge Type
A set of straight tracks parallel to the defined driving direction is generated. Each track is represented by two ends (nodes) that are located on the inner field boundary. The distance between subsequent tracks is equal to (W). Figure 2 demonstrates the generated work rows with different types of edges in the field.
There are six types of edges between each pair of the vertices and each defined edge has a length: • Gate to Headland (G2H): two edges connect the gate to the headland path • Headland: the connection between two subsequent vertices in headland path • Track: the connection between two pairs of nodes (two ends of a fieldwork area). Once a vehicle selects one end as the track entry, it has to finish the operation in the current track and exits at the opposite end of the current track before moving to another track. • Track to Headland (T2H): the connection between track nodes and vertices in the headland path • Track to track (T2T): the connection between the end nodes of two adjacent tracks • Headland to headland (H2H): the connection between two headland paths For the first type, the edge cost is the traveled distance between the headlands' vertices and the field's gate. For the second type, the edge cost is the length of the headland edge when the vehicle is operating in unproductive mode. For the third type, the edge cost is the length of the edge when the vehicle is operating in an unproductive mode. For the last type, the edge cost is the distance corresponding to the headland turning from the exit point of the current track to vertices in the headland path.

Track Generation and Edge Type
A set of straight tracks parallel to the defined driving direction is generated. Each track is represented by two ends (nodes) that are located on the inner field boundary. The distance between subsequent tracks is equal to (W). Figure 2 demonstrates the generated work rows with different types of edges in the field.
There are six types of edges between each pair of the vertices and each defined edge has a length: • Gate to Headland (G2H): two edges connect the gate to the headland path • Headland: the connection between two subsequent vertices in headland path • Track: the connection between two pairs of nodes (two ends of a fieldwork area). Once a vehicle selects one end as the track entry, it has to finish the operation in the current track and exits at the opposite end of the current track before moving to another track. For the first type, the edge cost is the traveled distance between the headlands' vertices and the field's gate. For the second type, the edge cost is the length of the headland edge when the vehicle is operating in unproductive mode. For the third type, the edge cost is the length of the edge when the vehicle is operating in an unproductive mode. For the last type, the edge cost is the distance corresponding to the headland turning from the exit point of the current track to vertices in the headland path.

Track Generation and Edge Type
A set of straight tracks parallel to the defined driving direction is generated. Each track is represented by two ends (nodes) that are located on the inner field boundary. The distance between subsequent tracks is equal to (W). Figure 2 demonstrates the generated work rows with different types of edges in the field.
There are six types of edges between each pair of the vertices and each defined edge has a length: • Gate to Headland (G2H): two edges connect the gate to the headland path • Headland: the connection between two subsequent vertices in headland path • Track: the connection between two pairs of nodes (two ends of a fieldwork area). Once a vehicle selects one end as the track entry, it has to finish the operation in the current track and exits at the opposite end of the current track before moving to another track. • Track to Headland (T2H): the connection between track nodes and vertices in the headland path • Track to track (T2T): the connection between the end nodes of two adjacent tracks • Headland to headland (H2H): the connection between two headland paths For the first type, the edge cost is the traveled distance between the headlands' vertices and the field's gate. For the second type, the edge cost is the length of the headland edge when the vehicle is operating in unproductive mode. For the third type, the edge cost is the length of the edge when the vehicle is operating in an unproductive mode. For the last type, the edge cost is the distance corresponding to the headland turning from the exit point of the current track to vertices in the headland path.

Optimization Algorithm
The problem explored in this study is to find the optimal traversal sequence of fieldwork tracks with an adapted application rate for each track to minimize the non-working driving distance in the field. Since there is a large discrete search space in this problem, it can be classified as a well-known VRP NP-hard problem where tracks' nodes will be visited instead of customers and the vehicle will refill at the depot. The meta-heuristic algorithm simulated annealing (SA) is applied to find the global optimum solution for this problem. The optimization model was developed using the PYTHON programing language (version 3.7) [28] 2.4.1. SA Algorithm The simulated annealing as a stochastic algorithm is used in this study to investigate the near-optimal solution. This method implemented in several problems in various fields such as computer design, route planning, image processing, etc. The performance structure of this algorithm is based on the simulation of the cooling process, which called annealing where a solid is gradually cooled. The SA algorithm has some primary features such as the initial temperature, cooling rate, and termination policy, which should be defined carefully in advance since the performance of the algorithm depends on them. Table 2 represents the corresponding values for the primary features of the SA algorithm. For instance, at the beginning when we have a higher temperature, the probability of accepting newer worse solutions is high and by spending time, it will decrease. By this ability, the SA algorithm enabled to search of all solution space and escape from trapping in local minima [29]. The overview of the SA algorithm could be listed as the following Flowchart ( Figure 3): AgriEngineering 2020, 3 FOR PEER REVIEW 6

Optimization Algorithm
The problem explored in this study is to find the optimal traversal sequence of fieldwork tracks with an adapted application rate for each track to minimize the non-working driving distance in the field. Since there is a large discrete search space in this problem, it can be classified as a well-known VRP NP-hard problem where tracks' nodes will be visited instead of customers and the vehicle will refill at the depot. The meta-heuristic algorithm simulated annealing (SA) is applied to find the global optimum solution for this problem. The optimization model was developed using the PYTHON programing language (version 3.7) [28] 2.4.1. SA Algorithm The simulated annealing as a stochastic algorithm is used in this study to investigate the nearoptimal solution. This method implemented in several problems in various fields such as computer design, route planning, image processing, etc. The performance structure of this algorithm is based on the simulation of the cooling process, which called annealing where a solid is gradually cooled. The SA algorithm has some primary features such as the initial temperature, cooling rate, and termination policy, which should be defined carefully in advance since the performance of the algorithm depends on them. Table 2 represents the corresponding values for the primary features of the SA algorithm. For instance, at the beginning when we have a higher temperature, the probability of accepting newer worse solutions is high and by spending time, it will decrease. By this ability, the SA algorithm enabled to search of all solution space and escape from trapping in local minima [29]. The overview of the SA algorithm could be listed as the following Flowchart ( Figure 3):

Neighborhood Operators
A new solution in the SA algorithm can be generated by applying various neighborhood operators on an initial solution. Random swaps, random swaps of subsequences, random insertions, random insertions of subsequences, reversing a subsequence, and random swaps of reversed subsequences are among the neighborhood operators, which is applied in this study. All mentioned operators demonstrated with an example in the following table (Table 3).

Initial Solution
In the fertilization process to avoid crossing the already applied area, first the inner part of the field and later the headland parts are going to be covered. The following flowchart ( Figure 4) represents an overview of the process of generating an initial solution.

Application Rates
In order to show the process of calculating the application rate for the tracks, a sample route with two tracks ([0, A, B, 0]) were selected and demonstrated in Figure 5. The required data for calculating the application rate for the tracks in one load are the length of tracks, the demands, the capacity of the slurry tank, the working width of the machine, and the deviation tolerance of the stipulated target application rate. All the applied variables and parameters in the process of calculating the application rates are defined in Table 4.

Application Rates
In order to show the process of calculating the application rate for the tracks, a sample route with two tracks ([0, A, B, 0]) were selected and demonstrated in Figure 5. The required data for calculating the application rate for the tracks in one load are the length of tracks, the demands, the AgriEngineering 2020, 3 FOR PEER REVIEW 9 capacity of the slurry tank, the working width of the machine, and the deviation tolerance of the stipulated target application rate. All the applied variables and parameters in the process of calculating the application rates are defined in Table 4.

Cost Matrix Generation
The cost matrix represents the transition cost between every pair of vertices of the field graph. Let M be the n * n, n = |V| matrix with element C ij , where C ij is the transition cost from vertex i to j when i j and is equal to 0 when i = j. In order to calculate the C ij , the generated field representation in the previous section is transformed as an undirected graph G = (N, E), where N = {n 0 , n 1 , . . . , n n }, n ∈ Z is the set of nodes consisting of all ending points of headland passes, fieldwork tracks, and connection paths, where E = n i , n j , . . . , n i , n j ∈ N is the set of edges representing all headland passes, fieldwork tracks, and connection paths. Each edge n i , n j , i j is associated with a weight, which corresponds to the actual length from node n i to node n j . Based on the derived graph, the distance between every pair of vertices can be computed by using the shortest path search algorithm (Dijkstra algorithm). It needs to be clarified that the vertex set V is a subset of N, since V only consists of depots and ending points of all tracks, while N is the set of all points. In the cost matrix, only the distance between every pair of vertices in V is stored as every C ij element.

Simulation Algorithm
A simulation model for the slurry application system is developed to calculate the total estimated non-working traveled distance (includes turnings' length and reloading distance) and the total estimated non-working time (includes turnings' time and reloading time). The refill event happens at the depot when the slurry tank is emptied. Hence, the simulation model calculates the transport distance and time and reloading time. An overview of the model is presented in Figure 6. The simulation model was developed using the MATLAB ® technical programming language (version R2018b) [30]. In the simulation model, a tool is developed to geometrically represent the field. The boundary of the field and headland areas demonstrated by the polygons and the inner part of the field is partitioned as a series of parallel working-areas (tracks). The input data for the geometrical representation model are the field boundary, number of headlands passes, driving direction, working width, and turning radius of the applied machine. The optimal sequence of tracks and their proper application rates determined by the optimization algorithm. The input data for the optimization algorithm consist of coordinates of fieldwork tracks, coordinates of headland's passes, the capacity of the slurry machine, machine's traveling, working and turning speed, the target application rate for the field, and the acceptable tolerance from the target application rate. The output of the geometrical representation model and the optimization algorithm was used as an input for the simulation algorithm. The output of the simulation algorithm is the segmentation of the task time and traveled distance for each part of the operation (productive or non-productive).
to when ≠ and is equal to 0 when = . In order to calculate the , the generated field representation in the previous section is transformed as an undirected graph = ( , ), where = { , , … , }, ∈ is the set of nodes consisting of all ending points of headland passes, fieldwork tracks, and connection paths, where = { , , … }, , ∈ is the set of edges representing all headland passes, fieldwork tracks, and connection paths. Each edge , , ≠ is associated with a weight, which corresponds to the actual length from node to node . Based on the derived graph, the distance between every pair of vertices can be computed by using the shortest path search algorithm (Dijkstra algorithm). It needs to be clarified that the vertex set is a subset of , since only consists of depots and ending points of all tracks, while is the set of all points. In the cost matrix, only the distance between every pair of vertices in is stored as every element.

Simulation Algorithm
A simulation model for the slurry application system is developed to calculate the total estimated non-working traveled distance (includes turnings' length and reloading distance) and the total estimated non-working time (includes turnings' time and reloading time). The refill event happens at the depot when the slurry tank is emptied. Hence, the simulation model calculates the transport distance and time and reloading time. An overview of the model is presented in Figure 6. The simulation model was developed using the MATLAB ® technical programming language (version R2018b) [30]. In the simulation model, a tool is developed to geometrically represent the field. The boundary of the field and headland areas demonstrated by the polygons and the inner part of the field is partitioned as a series of parallel working-areas (tracks). The input data for the geometrical representation model are the field boundary, number of headlands passes, driving direction, working width, and turning radius of the applied machine. The optimal sequence of tracks and their proper application rates determined by the optimization algorithm. The input data for the optimization algorithm consist of coordinates of fieldwork tracks, coordinates of headland's passes, the capacity of the slurry machine, machine's traveling, working and turning speed, the target application rate for the field, and the acceptable tolerance from the target application rate. The output of the geometrical representation model and the optimization algorithm was used as an input for the simulation algorithm. The output of the simulation algorithm is the segmentation of the task time and traveled distance for each part of the operation (productive or non-productive).

Results
In order to validate the proposed algorithm and evaluate the enhancement in the operation efficiency caused by the optimized plan as compared to the operation executed based on conventional planning, a shallow injection fertilizing operation was performed and recorded. As a case study, a sample field was considered for comparison between the conventional and optimized plans. The conventional method indicates that farmers prefer to fill the slurry tank at the depot and after visiting the field, the distributor is returned to the depot with an empty tank. It means that all the capacity of the slurry tank should be applied in the field before returning to the depot. When adhering to this method, sometimes it is not possible to find a combination of tracks where the total required material for covering them is exactly equal to one load of the applicator. The ramifications of this are that farmers perform turns in the main part of the field or emptied the rest of the material into another field. Figure 7 demonstrates the target field and the location of slurry storage (Depot). The yellow points show the log file for the slurry application in that field. All the information regarding the field and characteristics of the applied machine presented in Table 5. Moreover, the information corresponding to the properties of the generated tracks showed in Table 6. These data are going to be used later in the optimization and simulation algorithm to calculate the optimal coverage plan.

Results
In order to validate the proposed algorithm and evaluate the enhancement in the operation efficiency caused by the optimized plan as compared to the operation executed based on conventional planning, a shallow injection fertilizing operation was performed and recorded. As a case study, a sample field was considered for comparison between the conventional and optimized plans. The conventional method indicates that farmers prefer to fill the slurry tank at the depot and after visiting the field, the distributor is returned to the depot with an empty tank. It means that all the capacity of the slurry tank should be applied in the field before returning to the depot. When adhering to this method, sometimes it is not possible to find a combination of tracks where the total required material for covering them is exactly equal to one load of the applicator. The ramifications of this are that farmers perform turns in the main part of the field or emptied the rest of the material into another field.  Figure 7 demonstrates the target field and the location of slurry storage (Depot). The yellow points show the log file for the slurry application in that field. All the information regarding the field and characteristics of the applied machine presented in Table 5. Moreover, the information corresponding to the properties of the generated tracks showed in Table 6. These data are going to be used later in the optimization and simulation algorithm to calculate the optimal coverage plan.     . In each path, the slurry tank starts from the depot (zero represents the depot) with full capacity (33,000 L) and it covers some tracks and when the tank is empty it returns to the depot for reloading. In the first path, tracks numbers 10, 6, 8, and part of track 9 were covered. In the second tour, the tracks' ID 7, 3, 5 and the rest of the track ID 9 were covered and the total required material for these tracks is 31,437 (liters) and there is 1563 (liters) of slurry remained in the tanker and the farmer decided to empty the tank into another field besides the target field (part a). The last route starts from the depot and it includes the tracks ID number 1, 4, 2, and remained material (8344 L) emptied in the adjacent field (parts b and c).
AgriEngineering 2020, 3 FOR PEER REVIEW 12 In the second tour, the tracks' ID 7, 3, 5 and the rest of the track ID 9 were covered and the total required material for these tracks is 31,437 (liters) and there is 1563 (liters) of slurry remained in the tanker and the farmer decided to empty the tank into another field besides the target field (part a). The last route starts from the depot and it includes the tracks ID number 1, 4, 2, and remained material (8344 L) emptied in the adjacent field (parts b and c). A log file was generated using recorded GPS data of the slurry distributor by considering a constant application rate for all the tracks. The information in the log file was used for calculation of non-working time and distance for a slurry tank in the field. The non-working time is the accumulation of GPS time when the distributer does not apply fertilizer in the field or when it is out of the field. The non-working traveled distance can be determined by calculation of the number of meters traveled by the machine in the non-working time. Table 7 demonstrates the results of the analysis for the conventional method. Based on the input data mentioned in Tables 5 and 6, the geometrical representation model divided the target field into 10 working areas with one headland part. Figure 9 demonstrates the target field with all the characteristics such as field boundaries, tracks, headland area, and the gate of the field. By applying the proposed method in this study, an optimal solution was generated for this field. Table 8 shows the results of the optimization and simulation model for the optimal solution. A log file was generated using recorded GPS data of the slurry distributor by considering a constant application rate for all the tracks. The information in the log file was used for calculation of non-working time and distance for a slurry tank in the field. The non-working time is the accumulation of GPS time when the distributer does not apply fertilizer in the field or when it is out of the field. The non-working traveled distance can be determined by calculation of the number of meters traveled by the machine in the non-working time. Table 7 demonstrates the results of the analysis for the conventional method. Based on the input data mentioned in Tables 5 and 6, the geometrical representation model divided the target field into 10 working areas with one headland part. Figure 9 demonstrates the target field with all the characteristics such as field boundaries, tracks, headland area, and the gate of the field. By applying the proposed method in this study, an optimal solution was generated for this field. Table 8 shows the results of the optimization and simulation model for the optimal solution. The solution is an optimal coverage plan which first covers the inner part of the field and then the headland area. The numbers without a letter (h) in the solution represent one of the track's sibling nodes. Moreover, headland passes divided into some edges which they are distinguished by letter (h). Zero represents the depot and in the solution, when there is a zero, it means a reloading event happened.
AgriEngineering 2020, 3 FOR PEER REVIEW 13 The solution is an optimal coverage plan which first covers the inner part of the field and then the headland area. The numbers without a letter (h) in the solution represent one of the track's sibling nodes. Moreover, headland passes divided into some edges which they are distinguished by letter (h). Zero represents the depot and in the solution, when there is a zero, it means a reloading event happened. To have a better comparison between the optimal solution and the conventional one, the time and distance from the field's gate to the depot was calculated and added to the optimal solution (traveled distance from the gate to the depot and back to the gate = 6861 m; time to go from gate to the depot and back to the gate = 1184 s). Since in the optimal solution there are three loads, therefore, the time and distance from the gate to the depot should be multiplied by 3 and added to the nonworking time and distance of the optimal solution.

Optimized Solution
<0, 12,13,16,17,0,20,9,8,5, 0, 4, 1, 3h, 4h, 5h, 6h, 7h, 8h, 9h, 10h, 11h, 12h, 13h, 14h, 15h, 16h, 17h, 18h, 19h, 20h, 21h, 22h, 23h, 24h, 25h, 26h, 27h, 28h, 29h, 30h, 31h, 32h, 33h, 34h, 35h, 36h, 37h, 38h, 39h, 40h, 41h, 42h, 43h, 44h, 45h, 46h, 47h, 48h, 49h, 50h, 51h, 52h, 53h, 54h, 55h, 56h, 57h, 58h, 59h, 60h, 61h, 62h, 63h, 64h, 65h, 66h, 67h, 68h, 69h, 70h, 71h, 72h, 73h, 74h, 75h, 76h, 77h, 78h, 79h, 80h, 81h, 82h, 83h, 84h, 0> The comparison showed that applying this new method can reduce the non-working driving distance by 18.6% and the non-working time by 28.1%. In the first tour of the conventional method, tracks' ID (10,6,8) are covered completely and due to the limited capacity, part of the track 9 is covered in that tour and machine-turned through the cropping part of the field to go back to the depot for refilling. These maneuvers inside the agricultural main field have negative effects on the soil and can reduce the crop's yield [7]. However, the proposed method in this study avoids turnings through the main field by adjusting the application rates to an optimal value to let the machine operate the entire track. Table A1 in the appendix part, demonstrates the values of application rates To have a better comparison between the optimal solution and the conventional one, the time and distance from the field's gate to the depot was calculated and added to the optimal solution (traveled distance from the gate to the depot and back to the gate = 6861 m; time to go from gate to the depot and back to the gate = 1184 s). Since in the optimal solution there are three loads, therefore, the time and distance from the gate to the depot should be multiplied by 3 and added to the non-working time and distance of the optimal solution. The comparison showed that applying this new method can reduce the non-working driving distance by 18.6% and the non-working time by 28.1%. In the first tour of the conventional method, tracks' ID (10,6,8) are covered completely and due to the limited capacity, part of the track 9 is covered in that tour and machine-turned through the cropping part of the field to go back to the depot for refilling. These maneuvers inside the agricultural main field have negative effects on the soil and can reduce the crop's yield [7]. However, the proposed method in this study avoids turnings through the main field by adjusting the application rates to an optimal value to let the machine operate the entire track. Table A1 in the Appendix A part, demonstrates the values of application rates for the generated optimized solution. In the second tour of the conventional method, to empty the slurry tanker, the farmer covered a small part of the adjacent field. The amount of applied slurry in the part (a) is equal to 1563 L and the machine traveled 329 m to reach the (a) part. Moreover, in the same situation, the vehicle traveled 342 m to arrive at the parts (b and c). The amount of slurry applied in these parts is equal to 8344 L. However, the proposed method in this study avoids these travels by adjusting the application rates to certain values that the machine's capacity fully used and at the end of each tour there is no remaining slurry in the tank.

The Strategic Decision for Machinery Selection
In order to make a strategic decision about the optimal capacity of the slurry distributor, the optimization and simulation model was applied to a range of capacities to find out in which capacity there is an optimal solution for the target field and the amount of non-working time/distance calculated for those solutions. Table 9 presents the average values of running the optimization and simulation model for five times for each possible capacity. The non-working time and distance were calculated based on the location of the field's gate. Table 9. The results of the optimal solution for possible capacities. The following diagram ( Figure 10) shows the values of non-productive time and distance related to each capacity that the optimization algorithm has found a solution for that. The comparison demonstrates that the capacity equal to 33 (m 3 ) has less non-working time and distance as compared to the other options.

Capacity of Slurry Tank (L) Weight (tons) Non-Productive Time (min)
Two more factors, the weight of the slurry tank and the required amount of tractor's horsepower for pulling the distributor were considered to make a better decision for choosing the proper capacity of the slurry tanker. Indeed, the weight of the distributer affects the soil compaction, and in turn, conceivably the crop yield [31][32][33]. To generalize this problem, the weighted sum model (WSM) [34] method was applied. According to this method, when there are m alternatives and n criteria, the best alternative is the one that satisfies (in the minimization case) the following expression: a ij w j , f or i = 1, 2, 3, . . . , m.
In this expression, the A * WSM−score is the WSM score of the best alternative, n is the number of decision criteria, a ij is the actual value of the i th alternative in terms of the j th criterion, and w j is the weight of importance of the j th criterion. AgriEngineering 2020, 3 FOR PEER REVIEW 15 Figure 10. Comparison of non-productive time and distances in different capacities.
Two more factors, the weight of the slurry tank and the required amount of tractor's horsepower for pulling the distributor were considered to make a better decision for choosing the proper capacity conceivably the crop yield [31][32][33]. To generalize this problem, the weighted sum model (WSM) [34] method was applied. According to this method, when there are alternatives and criteria, the In this expression, the − * is the WSM score of the best alternative, is the number of 354 decision criteria, is the actual value of the ℎ alternative in terms of the ℎ criterion, and is 355 the weight of importance of the ℎ criterion.

356
According to Table 9, there are 13 alternatives and 4 criteria and the only things that need to be 357 defined are the weight of each criterion. In this study, the significance of the criteria defined as the 358 following matrix (Table 10):  For simplicity, the values of Table 9 can be normalized before calculating the amount 361 of − * . Calculation shows that − * equal to 5.223 related to the capacity 33,000 L is the 362 minimum score among all the alternatives. Therefore, based on the defined weights in Table 10, the 363 proper capacity to cover the target field is 33,000 L.

364
Based on the defined significance of the criteria by the farmer, the selected capacity can change.

365
For instance, if the soil compaction get the highest weight and the weight matrix changes as follow:

Comparison of non-productive time and distance in different capacities
Non-productive time (minute) Non-productive distance (meter) Figure 10. Comparison of non-productive time and distances in different capacities.
According to Table 9, there are 13 alternatives and 4 criteria and the only things that need to be defined are the weight of each criterion. In this study, the significance of the criteria defined as the following matrix (Table 10): For simplicity, the values of Table 9 can be normalized before calculating the amount of A * WSM−score . Calculation shows that A * WSM−score equal to 5.223 related to the capacity 33,000 L is the minimum score among all the alternatives. Therefore, based on the defined weights in Table 10, the proper capacity to cover the target field is 33,000 L.
Based on the defined significance of the criteria by the farmer, the selected capacity can change. For instance, if the soil compaction get the highest weight and the weight matrix changes as follow: (W1 = 4, W2 = 1, W3 = 2, W4 = 3), as a result, the value of A * WSM−score would be equal to 6.104 which is the score of the capacity 15,000 L. Moreover, if we consider the same weight for all the criteria equal to one then, the proper capacity would be 31,000 L with the score equal to 2.603.

Discussion
Coverage planning for capacitated operation (such as fertilization) is presented in other studies [22,35]. In order to compare the presented method in this study with other techniques, another approach for coverage planning of capacitated field operations was considered to show the advantages or drawbacks of the proposed approach. According to the study by [26]  , they applied an algorithmic approach for the optimization of liquid fertilizing operation. Their method is based on the state-space search and the output solution is a sequence of pre-defined driving actions to change an initial state to a goal state by minimizing the non-working traveled distance in the target field. The soil injection fertilizing operation A was considered in this part for comparison. All the input data regarding the operation A were presented in the following table (Table 11).  Figure 11 demonstrates the optimized coverage plan for operation A. Based on the results, the total amount of non-productive distance can be determined by accumulating the non-working distance for each route, equal to 2730 m.

Discussion
Coverage planning for capacitated operation (such as fertilization) is presented in other studies [22,35]. In order to compare the presented method in this study with other techniques, another approach for coverage planning of capacitated field operations was considered to show the advantages or drawbacks of the proposed approach. According to the study by [26]  , they applied an algorithmic approach for the optimization of liquid fertilizing operation. Their method is based on the state-space search and the output solution is a sequence of pre-defined driving actions to change an initial state to a goal state by minimizing the non-working traveled distance in the target field. The soil injection fertilizing operation A was considered in this part for comparison. All the input data regarding the operation A were presented in the following table (Table 11).  Figure 11 demonstrates the optimized coverage plan for operation A. Based on the results, the total amount of non-productive distance can be determined by accumulating the non-working distance for each route, equal to 2730 m. The information presented in Table 11 was applied in the optimization and simulation model as input data and as a result, an optimized solution generated and provided in the following table (Table  12). Moreover, the list of optimized application rate for all the work rows generated and demonstrated in Table A2 in appendix part. The amount of traveled non-working distance based on the presented method in this study is slightly (58 m) more than the approach provided by [26]  . However, the presented solution in this study avoids all the turning through the main crop area of the field, and in this way, there is less soil compaction in comparison with the method provided by [26]  . Figure 12 shows the optimized coverage plan for operation A based on the proposed method in this study. The information presented in Table 11 was applied in the optimization and simulation model as input data and as a result, an optimized solution generated and provided in the following table (Table 12). Moreover, the list of optimized application rate for all the work rows generated and demonstrated in Table A2 in Appendix A part. The amount of traveled non-working distance based on the presented method in this study is slightly (58 m) more than the approach provided by [26]  . However, the presented solution in this study avoids all the turning through the main crop area of the field, and in this way, there is less soil compaction in comparison with the method provided by [26]  . Figure 12 shows the optimized coverage plan for operation A based on the proposed method in this study.  Figure 12. Field partition related to operation A.

Conclusions
A coverage-planning algorithm was developed for capacitated field operations such as slurry applications in agricultural fields. The meta-heuristic algorithm simulated annealing (SA) was applied to find the optimal traversal sequence of fieldwork tracks under the criterion of minimizing the non-working traveled distance during the headland turnings together with an adapted application rate for each track. Furthermore, a simulation model was developed to generate the segmentation of the task time and traveled distance for each element (productive/non-productive) of the operation.
The operations efficiency of the optimized plans generated by the proposed method in this paper were compared with conventional non-optimized methods used by farmers. The results showed that applying this new method could accomplish 18.6% and 28.1% reduction in the non-working traveled distance and the non-working time, respectively.
The proposed algorithm could be used as part of a tool to support strategic decisions about the capacity of the machine in capacitated operations such as slurry applications. Four factors (the weight of slurry tanker, non-productive time/distance, and required tractor's hp) were considered for comparison between different capacities to find the proper one.
Future research will focus on a method to find the optimal coverage plan for a target field by specific configurations of the machine (such as capacity). In this case, the algorithm will be able to find an optimal coverage plan for a sample field with any defined capacity of the slurry tank.

Conclusions
A coverage-planning algorithm was developed for capacitated field operations such as slurry applications in agricultural fields. The meta-heuristic algorithm simulated annealing (SA) was applied to find the optimal traversal sequence of fieldwork tracks under the criterion of minimizing the non-working traveled distance during the headland turnings together with an adapted application rate for each track. Furthermore, a simulation model was developed to generate the segmentation of the task time and traveled distance for each element (productive/non-productive) of the operation.
The operations efficiency of the optimized plans generated by the proposed method in this paper were compared with conventional non-optimized methods used by farmers. The results showed that applying this new method could accomplish 18.6% and 28.1% reduction in the non-working traveled distance and the non-working time, respectively.
The proposed algorithm could be used as part of a tool to support strategic decisions about the capacity of the machine in capacitated operations such as slurry applications. Four factors (the weight of slurry tanker, non-productive time/distance, and required tractor's hp) were considered for comparison between different capacities to find the proper one.
Future research will focus on a method to find the optimal coverage plan for a target field by specific configurations of the machine (such as capacity). In this case, the algorithm will be able to find an optimal coverage plan for a sample field with any defined capacity of the slurry tank.