Research on Flexible Job Shop Scheduling Method for Agricultural Equipment Considering Multi-Resource Constraints

1. Introduction

As the main manifestation of agricultural productivity, agricultural mechanization equipment is important material foundations for ensuring food security in the context of reducing agricultural labor and achieving agricultural modernization in the context of complex crops. Agricultural mechanization equipment is an important material foundation for the development of modern agriculture. Promoting the development of agricultural mechanization equipment is an objective requirement for improving agricultural labor productivity, land output rate, and resource utilization rate. It is a practical need to support the development of agricultural mechanization, the transformation of agricultural development mode, and the quality and efficiency of agriculture. The production of agricultural mechanization equipment is greatly influenced by market demand, with fast updates and replacements of models and configurations. The demand is highly targeted, and different models need to be selected for different regions, environments, and operational needs[1]. At the same time, the agricultural operating environment varies greatly, and the busy farming hours also vary with climate change, resulting in rapid changes in market demand for various agricultural equipment and diverse demands for machine models. Due to the complexity of the market environment, compared with other mechanical equipment manufacturing industries, agricultural machinery and equipment enterprises have the characteristics of more model selection and unstable market demand cycles in production. Therefore, multi variety, small batch, and customized production methods have become the mainstream of agricultural machinery enterprises, and at present, in order to meet the diversified needs of the market, the production tasks have the characteristics of multi variety and small batch production[2].

In order to meet the diverse demands of the agricultural mechanization equipment production market and achieve flexible production goals of multiple varieties and small batches, manufacturing enterprises can adopt flexible production lines for production. Flexible Job Shop(FJS) is an advanced form of production organization that combines high flexibility and customized production capabilities to adapt to the changing market demands and customer orders in modern manufacturing. It can improve production efficiency, expand production capacity, and meet the assembly needs of multiple types of materials. It also has a certain degree of flexibility, and the products it produces can meet the different needs of customers, and can make corresponding adjustments according to the changes in different customer needs.

The flexible manufacturing shop of agricultural mechanization equipment belongs to the flexible job shop, and its corresponding workshop scheduling problem can be attributed to FJSP. FJSP was firstly introduced in reference[3], accompanied by a polynomial algorithm tailored for addressing small-scale problem instances. Since Brucker and Schlie[3] first delved into FJSP, it has garnered significant attention from numerous scholars. Consequently, many methods have been proposed to tackle this complex issue. The FJSP has been found widespread application in various industrial fields[4-8], including automotive assembly, agricultural machinery equipment production, semiconductor manufacturing and so on.

The extended FJSP which considers the worker-resources-constraint can be seen as a type of dual-resource constrained flexible job shop scheduling problem (DRCFJSP)[9]. Wu et al.[10] addressed a DRCFJSP that incorporated worker learning, with the goal of minimizing the makespan. They devised a hybrid genetic algorithm, integrating variable neighborhood search (VNS), to explore a wide range of potential solutions. Computational tests showed that their algorithm was effective in solving the problem. Kress, Müller, and Nossack[11] studied a DRCFJSP considering setup times and worker skills. They aimed to minimize makespan and tardiness using exact and heuristic methods by splitting the problem into routing and worker assignment subproblems. Gong et al.[12] studied a Double Flexible Job-Shop Scheduling Problem (DFJSP) with a goal to minimize makespan and maximal total worker cost. They proposed a new Hybrid Genetic Algorithm (NHGA) to solve this problem, which outperformed the Non-dominated Sorting Genetic Algorithm II (NSGA-II) in terms of accuracy and efficiency. Vital-Soto et al.[13] tackled the challenge of simultaneously determining machines, worker assignments, and sequencing flexibility in a flexible job-shop setting. Their objective was to minimize makespan, worker workload, and weighted tardiness. They developed a multi-objective optimization model and applied a unique variation of the elitist non-dominated sorting genetic algorithm (NSGA-II) to solve the complex problem they faced.

Collaborative scheduling in workshops involves simultaneously allocating processing machines and AGVs for tasks, aiming to optimize production and transportation resources in a specific workshop environment. Raman et al.[14] first introduced scheduling of AGVs in the flexible job scheduling problem (FJSPT), using a mixed integer programming model to minimize completion time. They treated transportation like operations, with AGVs as machines and travel times as processing times, applying traditional flexible scheduling methods to allocate both machines and AGVs. Yan et al.[15] studied the FJSP problem under limited transportation conditions in a digital twin workshop and designed an improved genetic algorithm to minimize makespan. Hu et al.[16] combined deep reinforcement learning with hybrid AGV scheduling rules for real-time scheduling in flexible workshops. Upon task completion or job arrival, a DQN-based scheduler selects a scheduling rule to allocate the most suitable AGV, aiming to minimize makespan and delay rates. Umar et al.[17] proposed a hybrid multi-objective genetic algorithm to obtain a Pareto solution to the problem of minimizing multiple conflicting objectives, including total makespan, AGV travel time, and delay cost. Similar studies include literature [18, 19].

In summary, many current researches on the dual resource constraints in FJSP mainly focuses on worker and machine constraints, (Automated Guided Vehicle)AGV and machine constraints. However, with the advent of intelligent manufacturing and Industry 4.0 era, the FJSP under multi-resource constraints i.e., Multi-Resource Constrained Flexible Job Shop Scheduling Problem(MRCFJSP) that considers workers, AGVs, and machines together is more in line with the reality of today's workshop production and has become increasingly important and complex. In this paper, a new mathematical model considering machine, worker and AVG in MRCFJSP is established. This model considers operation-sequence-constraint as well as resources constraint such as machine, worker and AVG. To effectively solve this problem, an effective MODGWO algorithm is proposed. The corresponding encoding/decoding method, population initialization strategy, crossover operators and local search technology are designed according to the features of this problem.

2. Materials and Methods

2.1. Problem description

The MRCFJSP considering multiple resource constraints of machines, workers, and AGVs can be described as follows:

There are n jobs J={J₁, J₂,……,J_n} to be processed, which need to be transported to m machines M={M₁, M₂,……,M_m} by v AGVs V={V₁, V₂,……,V_w} for operations by w workers W ={W₁, W₂,……,W_l}. The machining process of the job is determined, and these processes need to be processed one by one under specific process constraints. The machining process of different jobs is different, and the j-th operation(O_ij) of job J_i can be machined on one or more machines. The available set of processing machines for operation O_ij is M_ij⊆M. The processing time of the same operation on different machines is determined and different. Under the above conditions, with the goal of minimizing the maximum completion time(makespan) and minimizing the workers' continuous working time, allocate a machine M_k∈M_ij and a AGV V_v∈V to operation O_ij, allocate a worker W_w∈W to machine M_k, and determine the processing sequence of the same operation on the same machine and the transportation sequence of the same operation on the AGV, i.e., decide on the five sub-problems: machine assignment, worker assignment, and AGV assignment, operation sequence and transportation sequence. In addition, the following assumptions need to be satisfied:

(1) Each machine can process at most one operation at a time.

(2) Only one machine must be operated by a worker at a time.

(3) Only one machine and one AVG can be selected for each operation, and they cannot be interrupted during processing.

(4) All the jobs, machines and workers are available at zero time. All the jobs and AGVs are located in loading/unloading(LU) area at zero time.

(5) The processing operations for the same job shall meet the constraints of the operation sequence.

(6) Disregarding the transfer time of workers between machines, a worker can only operate one machine at a time.

(7) Machines and AVGs failures are not considered and the AVGs is sufficiently charged, and the loading and unloading time of the material on the machine is not considered.

(8) An AVG can only load one job at a time.

(9) After the current transportation task is completed, the AGV will not return to the loading and unloading area, but will go to the machine where the next transportation task is to be performed.

(10) Each machine has a buffer zone that can be used to park AGVs and store operations with sufficient buffer capacity.

Figure 1 illustrates a specific scenario. AGV, workers, and machines are regarded as three different resources in the production workshop. AGV is responsible for transporting raw materials from the loading area to various machines, thus establishing an important logistics chain. Starting from the loading area, the AGV picks up tasks one by one and transports them to the designated machine, as shown in loaded AGVs in Figure 1. Then, the machine is operated by a worker to start the processing task. Alternatively, the AGV can run empty to a certain machine and transport the semi-finished products to the machine to be processed in the next process, such as unloaded AVGs. Workers operate machines when they have processing tasks assigned, and they can take a break in the resting area when not assigned to any tasks. This continuous cycle process ensures efficient workflow and sustained productivity in the workshop.

2.2. Mathematical model

The notations required for the mathematical model are listed in Table 1.

Based on the above problem definition, constraints, assumptions and notations, the mathematical model of MRCFJSP is established as follows.

$$\min C_{\max }=\underset{1\le k\le m}{\max }(C_k)$$

(1)

(2)

This paper considers two optimization objectives. Eq. (1) represents minimizing the maximum completion time of last job(makespan). The second objective is to minimize the maximum continuous working hours of all workers, which is shown in Eq. (2).

Subject to:

$$CTij-STij={\displaystyle \sum _{\mathrm{k}=1}^{m}xijk}Pijk\mathrm{s},i=1,2,\mathrm{...}n;j=1,2,\mathrm{...}ri;s=1,2,\mathrm{...},l$$

(3)

(4)

$$CTij{\displaystyle \sum _{k=1}^{m}xijk}+L(1-{\displaystyle \sum _{k=1}^{m}ziji\text{'}j\text{'}k})\le STi\text{'}j\text{'}{\displaystyle \sum _{k=1}^{m}xi\text{'}j\text{'}k},i=1,2,\mathrm{...}n;j=1,2,\mathrm{...},ri;j\text{'}=1,2,\mathrm{...}ri\text{'}$$

(5)

(6)

$$CTijv{\displaystyle \sum _{v=1}^w\alpha ijv}+L(1-{\displaystyle \sum _{v=1}^w\beta iji\text{'}j\text{'}v})\le STi\text{'}j\text{'}v{\displaystyle \sum _{v=1}^w\beta iji\text{'}j\text{'}v},i=1,2,\mathrm{...}n;j=1,2,\mathrm{...},ri;j\text{'}=1,2,\mathrm{...},ri\text{'}$$

(7)

(8)

$$\begin{array}{l}[(CThg-CTij-SThg)xhgqsxijks\ge 0]\vee [(CTij-CThg-STij)xhgqsxijks\ge 0]\\ i,h=1,\mathrm{2...},n;j=1,2,\mathrm{...}ri;g=1,2,\mathrm{...}rh;k,q=1,2,\mathrm{...},m;s=1,2,\mathrm{...},l\end{array}$$

(9)

(10)

$${\displaystyle \sum _{s=1}^{l}xijks=1},i=1,2,\mathrm{...}n;j=1,2,\mathrm{...}ri;k=1,2,\mathrm{...}m$$

(11)

(12)

$$xijk\in \{0,1\},i=1,2,\mathrm{...}n;j=1,2,\mathrm{...}ri;k=1,2,\mathrm{...}m$$

(13)

(14)

$$yiji\text{'}j\text{'}s\in \{0,1\},i,i\text{'}=1,2,\mathrm{...}n;j,j\text{'}=1,2,\mathrm{...}ri,s=1,2,\mathrm{...}l$$

(15)

(16)

$$\alpha ijv\in \{0,1\},i=1,2,\mathrm{...}n;j=1,2,\mathrm{...}ri;v=1,2,\mathrm{...}w$$

(17)

(18)

where Eq. (3) indicates that the completion time of each job is after the start time. Eq. (4) represents the processing sequence constraint of a job. Eq. (5) to Eq. (7) represent the processing time constraint, which can guarantee that each operation cannot start until the assigned machine, worker and AVG are ready. Each machine can only process at most one operation and can only be operated by one worker at a time are shown in Eq. (8) and Eq. (9), respectively. Eq. (10) to Eq. (12) can ensure that each operation can be assigned to only one machine, one worker and one AVG respectively. Eq. (13) to Eq. (18) indicate that all decision variables are 0-1 variables.

2.3. Grey Wolf Optimization algorithm

The Grey Wolf Optimization(GWO)[20] algorithm was initially proposed by Australian scholars Seyedali Mirjalili et al. in 2014 for continuous real parameter optimization problems, and has been improved to handle some practical engineering problems[21-24]. As shown in Figure 2, the grey wolf pack has a very strict social hierarchy similar to the pyramid. The first three hierarchical categories are Alpha(α), Beta(β), and Delta(δ) as the best solutions to lead the rest of the wolves named Omega(ω) wolves toward promising areas in order to find the global solution.

A wolf(i) evaluates its proximity to the three best solutions based on Eq. (19)–(22) and subsequently updates its position using Eq. (23).

(19)

(20)

(21)

(22)

(23)

Where $\overrightarrow{X}$_i, $\overrightarrow{X}$_α, $\overrightarrow{X}$_β, and $\overrightarrow{X}$_δ, indicate the position vectors of wolf i, α, β and δ; $\overrightarrow{X}$_i(t+1) is the location of the wolf i in the next iteration; $\overrightarrow{r}$₁and $\overrightarrow{r}$₂are random vectors in [0,1]; $\overrightarrow{A}$ and $\overrightarrow{B}$ are two coefficient vectors where the vector $\overrightarrow{a}$ is linearly decrease from 2 to 0 over the course of iterations and can be calculated by Eq. (24).

(24)

Where where t is the present number of iteration and MaxIter represents the max number of iteration.

2.4. Multi-Object Discrete Grey Wolf Optimization algorithm for MRCFJSP

For continuous optimization problems, the grey wolf individual position can be updated by Eq.(19)-(23), however, MRCFJSP is a discrete problem and cannot be directly solved using the above formulas. This paper proposed a MODGWO algorithm for MRCFJSP based on the characteristics of the problem and combined with cross operation in genetic algorithm. First, an encoding/decoding method is introduced, followed by a method for initializing a mixed population. Then the Discrete Grey Wolf iteration operators are presented for global search. Finally, an effective VNS with four neighborhood structures is introduced to improve the local search ability.

To clarify the proposed algorithm, an 3×3×2×2 instance involving three jobs, three machines, two workers and two AVGs is provided in Table 2 and Table 3. Table 2 outlines the fundamental processing durations, while Table 3 displays the information of machines that workers can operate. The symbol ‘–’ signifies that a machine is incompatible with a particular operation or a worker is unable to operate a specific machine. In addition, by default, each AGV can transport all jobs.

2.4.1. Four-vector encoding

The efficiency of encoding and decoding techniques has a substantial influence on algorithm performance. Utilizing effective encoding and decoding methods can facilitate algorithm operations and enhance computational performance. The MRCFJSP needs to solve four sub-problems: operation scheduling, machine selection, worker selection, and AVG selection. Therefor in this subsection, we introduce an encoding scheme where each individual comprises four distinct vectors: the Operation Sequence(OS) vector, which denotes the order of operation execution; the Machine Selection(MS) vector, which specifies the assignment of machines; the Worker Selection(WS) vector, which outlines the allocation of workers; and the AVG Selection(AS) vector, which indicts the assignment of AVGs. Collectively, these four vectors constitute a viable solution for the MRCFJSP. In the four-layer coding method, we borrowed the idea of MOSO coding method proposed by Ho et al[25]. Figure 3 shows the four vector representation approach of the examples in Table 1 and Table2, where green indicates job 1, yellow indicates job 2 and blue indicates job 3.

(1) Operation sequence vector

An OS vector consists of an array of integers where the number of array elements equals the total number of operations in all jobs. Each operation is denoted by its corresponding job index. Specifically, the jth occurrence of index i signifies the jth operation within job i, and the count of index i repetitions is equivalent to the total number of operations in job i.

For example {3, 2, 1, 2, 3, 1, 1} is a legal vector or the aforementioned instance shown in Table 2 and Figure 3(a). The sequence of operational decoding for this OS vector is O₃₁→O₂₁→O₁₁→O₂₂→O₃₂→O₁₂→O₁₃.

(2) Machine selection vector

A MS vector consists of an array of integers, the length of which equals to the length of OS vector and each element indicates an available machine for the corresponding operation.

A MS vector {1, 3, 2, 1, 3, 1, 2} means that for the operation set {O₁₁, O₁₂, O₁₃, O₂₁, O₂₂, O₃₁, O₃₂}, the machine selection indicated by the MS vector is {M₁, M₃, M₂, M₁, M₃, M₁, M₂} correspondingly shown in Figure 3(b).

(3) Worker selection vector and AGV selection vector

The WS vector and AS vector are encoded in a similar way to MS vector. As shown in Figure 3(c), the first number 1 in the WS vector represents the O₁₁ is processed by worker W₁ on machine M₁. The WS vector {1, 2, 2, 1, 1, 1, 2} indicates that the operations {O₁₁, O₁₂, O₁₃, O₂₁, O₂₂, O₃₁, O₃₂} are processed by workers {W₁, W₂, W₂, W₁, W₁, W₁, W₂} on machines {M₁, M₃, M₂, M₁, M₃, M₁, M₂} respectively. Similarly, the AS vector {1, 1, 2, 1, 2, 2, 1} shown in Figure 3(d) represents that the operations{O₁₁, O₁₂, O₁₃, O₂₁, O₂₂, O₃₁, O₃₂} are transported by AGVs{V₁, V₁, V₂, V₁, V₂, V₂, V₁} respectively.

2.4.2. Decoding

There are mainly four approaches to decode a solution represented by the four-vector encoding scheme: non-delay, semi-active, active and hybrid schedules. In this paper, the active schedule is utilized to decode all the solutions. The decoding procedure is as follows: (1) Pick out the operations from the OS vector one by one from left to right, assign the corresponding machine, worker and AGV to the current job in sequence. (2) Determine the starting time of the current job from the earliest available start time of the corresponding machine, worker, and AVG, until the last job is processed. (3) Arrange the operations at their earliest feasible time. (4) Calculate the maximum completion time and maximum worker duration, which are the values of the two objective functions.

2.4.3. Population initialization

The standard of initial solutions significantly influences how well an algorithm performs in producing high-quality solutions. Recognizing the significant impact of population initialization on the algorithm's performance, we propose several adaptive strategies aimed at enhancing the quality of the initial population. To achieve a balance between the quality and diversity of initial solutions within the search space, we utilize a hybrid approach that combines random generation with strategic selection to form the initial population. Since the problem in this paper is divided into four sub-problems, operation sequencing, machine selection, worker selection and AVG selection, the population initialization is also carried out in four stages. Firstly, the initialization of OS vectors employs a hybrid approach[26] according to the following three rules:

(1) Most work remaining rule[27]: the job that has the maximum total processing time remaining is prioritized for execution first.

(2) Most number of operations remaining rule[26]: the job that has the most remaining operations is processed first.

(3) Random rule: operation sequencing is randomly generated.

The aforementioned three rules are utilized in the generation of OS vectors, and their respective proportions are 30%, 30%, and 40%. The initialization of MS vectors is initialized employing the following three rules:

(1) Global minimum processing time rule[28]: select the machine with the shortest processing time for the current operation.

(2) Workload considered rule[29]: Assign machines according to the machine workload.

(3) Random rule: Randomly assign machines.

The machine selection of the population is sequentially generated using the first and second rule in the proportion of 30% and the remaining 40% by random rule. The initialization process of the WS vectors and the AS vectors is similar to that of MS vectors, the only difference is that only workers and AVGs load, as well as random allocation, are used, with respective proportions of 50% and 50%.

2.4.4. Solution updating mechanism

Since the basic GWO cannot be directly applied to the discrete scheduling problem, this paper designs a discrete version of the solution updating mechanism for the machine selection, AGV selection and operation sequence. The grey wolf individual crosses with wolf α, wolf β, and wolfδaccording to a certain probability, and thus obtains a new individual, as shown in Eq. (25).

(25)

where X_k is the position vector of grey wolf t, X_α, X_β, and X_δ are the position vectors of grey wolf α, β, andδ, respectively, rand is a random number that follows a uniform distribution between 0 and 1. Cross is the crossover operator, and the specific crossover method for each vector is described in detail below.

In basic GWO, new individuals are generated based on the information of the three individuals α, β, andδ that are currently serving as the decision level individuals. To determine α, β, andδ individuals within a population, this paper adopts a method based on non dominated level and crowding distance to obtain decision level individuals, that is, sorting individuals according to their non dominated level and crowding distance in the population, and the top three individuals have the opportunity to become decision level individuals. For any two individuals, the individual with the lower level ranks is selected first, if two individuals have the same level, compare their crowding distance, and the individual with the larger crowding distance ranks first.

Meanwhile, we integrate the ideas of Jaya algorithm into the proposed MODGWO. The Jaya algorithm is a novel population-based metaheuristic recently proposed by Rao[30], which is based on the principle of continuous improvement, and improves the quality of the solution by moving the individuals closer to the best solution while moving away from the worst solution. Jaya algorithm iteratively evolves through Eq. (26). to obtain new solutions.

(26)

where X_k(t+1) is the updated solution corresponding to the initial solution X_k(t). X_best and X_worst represent the best and worst solutions of the current population, respectively. r₁ and r₂ are uniformly distributed. In the improved MODGWO, one of X_α, X_β, and X_δ is randomly chosen as X_bestusing Eq. (25), while the worst individual X_worst is selected from the last rank having the minimum crowding distance value. The detailed procedures for operation sequencing, machine selection, worker selection and AVG selection are summarized below.

(1) Operation sequencing

This paper draws on the method in literature [31] to delete sequences in the current solution that are identical to the worst solution and replace them with the best solution. The specific steps are as follows:

Step 1. Identify the worst individual (X_worst) and choose one form X_α, X_β, and X_δrandomly as X_best (assuming that X_β is selected) from the population.

Step2. Identify similar assignments of operations in the individual X_i and X_worst.

Step 3. Eliminate similar assignments from X_i and transfer the remaining assignments to a new solution X_inew.

Step 4. Beginning from the initial assignment of X_β, eliminate the assignments corresponding to the remaining assignments of individual X_inew.

Step 5. Beginning from the initial assignment of X_β, transfer the remaining assignments of X_β in the vacant elements of individual X_inew in the corresponding order.

Figure 4 illustrates the operation sequencing vector updating mechanism.

(2) Machine selection and AGV selection

The solution updating mechanism for the machine selection and AVG selection is the same and similar to that of operation sequencing.Only the detailed steps for updating machine selection are as follows:

Step 1. Identify the best individual X_best(Assuming it is X_α) and worst individual X_worst from the population.

Step 2. Identify similar assignments of machines assigned to the same machines in the individuals X_i and X_worst.

Step 3. Eliminate similar assignments from the individual X_iand transfer the remaining assignments to a new solution X_inew.

Step 4. Transfer the corresponding machine assignments of X_α to the vacant elements of individual X_inew.

Figure 5 illustrates the machine selection vector updating mechanism.

(3) Worker selection

The updating mechanism of worker selection is the same as the previous steps of machine selection and AVG selection, eliminating elements from the current individual that correspond to the position of the worst solution, and filling the eliminated elements with elements corresponding to the position of the best solution. Finally, it is necessary to conduct feasible checks on new individuals. Check for the feasibility of each assignment of X_inew based on the capable worker set of the corresponding machine in the updated machine selection vector. If infeasible randomly allocate an integer value within the capable worker set range.

2.4.5. Local search

To enhance the local search capability, variable neighborhood search operations are performed on the OS, MS, WS and AS of the optimal individuals after each iteration, respectively. Variable neighborhood search (VNS) can expand the search space of the algorithm by changing the domain structure based on problem features, avoiding the algorithm from falling into the local optimal solution. Firstly, the initial solution that requires neighborhood search is given. Secondly, the neighborhood structure is used for neighborhood search to obtain the neighborhood structure solution. If the target value of the neighborhood solution is better than the current solution, the current solution is updated and the weights are updated. The search returns to the first neighborhood structure and starts again until a better solution cannot be found. Then, the next neighborhood structure is reached for search until all neighborhood structures have been searched. Based on the characteristics of the MRCFJSP, the following 4 domain structures are constructed:

(1) VNS1(Insert, for OS vectors)Within the length range of OS vector, randomly generate two positions r1 and r2 (r1<r2), insert the element corresponding to r2 into position r1, and move the elements after position r1 backwards in sequence. For the OS vector in Figure 3, the two generated random numbers r1 and r2 are 2 and 5, respectively. The solution process of VNS1 is shown in Figure 6.

(2) VNS2(Choose the minimum processing time，for MS vectors)Randomly generate r ϵ [1, N] random numbers, where N is the total number of operations, and randomly select r positions from the MS vector. For each selected position corresponding to a operation, select the machine with the shortest processing time from its corresponding optional machine set for replacement. After replacement, check whether the workers corresponding to r positions in WS can operate the newly selected equipment. If it is not feasible, randomly select a worker from the operable worker set of the new equipment. For the solution generated in Table 2.1, MS is 1122312, and r=3 is randomly generated. Three positions are randomly selected, such as positions 2, 3, and 7. The solution process of VNS2 is shown in Figure 7.

(3) VNS3(Flip, for WS vectors)Within the length range of WS vector, randomly generate a position r1. In MS vector, position r1 should correspond to more than one available processing worker for the processing equipment, and then select a new worker from the set of available workers for equipment replacement.The solution process of VNS3 is shown in Figure 8.

(4) VNS4(swap, for AS vectors)Within the length of the AS vector, two randomly generated positions r1 and r2 (r1 ≠ r2) are exchanged for the elements at the r1 and r2 positions. The solution process of VNS4 is shown in Figure 9.

3. Results

The proposed approach is coded in Python 3.9 on an Intel i7 processor with a speed of 2.6 GHz on a Windows 11 operating system with 16 GB RAM.

3.1. Performance metrics

To assess the efficacy of the Pareto front derived from various perspectives, we employ two metrics: IGD[32], and the Hypervolume[33]. The aforementioned metrics encompass various facets concerning the quality of non-dominated solutions, particularly focusing on those derived from Pareto-based approaches. The detailed explanation of these two metrics is summarized below.

(1) Inverse Generational Distance (IGD): Measure the average distance from a predefined set of uniformly spaced points along the true or ideal Pareto Front (denoted as PF*) to the nearest solution in the obtained PF. Essentially, IGD assesses how well the obtained solutions approximate the true Pareto Front by calculating these distances and averaging them.

The formula for the IGD is typically expressed as follows:

(27)

where:

• |PF*| is the number of uniformly spaced points along the true Pareto Front.

• PF*_i​ represents the i-th point in PF*.

• PF is the set of obtained solutions (Pareto Front).

• PF_j represents the j-th solution in PF.

• d(PF*_i , PF_j)is the Euclidean distance between point PF*_i​and the closest solution PF_j​in PF.

A lower IGD value indicates that the obtained PF is closer to the true Pareto Front, suggesting better performance of the optimization algorithm in identifying high-quality solutions. Since the actual PF* of MRCFJSP in this paper is unknown, the non dominated solutions from multiple runs of all algorithms are first integrated into a hybrid solution set, and then the non dominated solutions in the hybrid solution set are used as PF*.

(2) Hypervolume(HV): HV is the volume of a hyper polyhedron enclosed by non dominated solutions and reference points in the search space, calculated as follows:

(28)

where P represents the Pareto Front approximate non dominated solution set obtained by the algorithm, F represents the reference points of the solution set, v is the volume of a hypercube enclosed by the p-th solution and the reference points. The target vectors corresponding to all Pareto solutions are normalized and the reference point for calculating HV is set to (1, 1).

3.2. Test instances generation

This paper extends the classic Brandimarte test case [27] in the flexible job shop scheduling problem by randomly adding worker and AVG information for solving the MRCFJSP. The number of workers in each workshop is calculated by multiplying the number of machine in the workshop by 0.6 and rounding up, and the probability that a worker can operate a certain machine is 0.5. For each test case, the ratio of the number of AGVs to the number of machines is considered to be 1:1, 1:2, and 1:3, respectively. Taking MK01 as an example, this example has 10 jobs and 6 machines. In the generated new example, the number of workers is 4, and the number of AVGs is 6, 3, and 2, respectively. The names of the extended examples are MK01_4_6, MK01_4_3, and MK01_4_2. In this way, 45 instances for MRCFJSP are used to test in this study, and Table 4 shows the relevant information, where n represents the number of jobs and m represents the number of machines..

3.3. Effectiveness verification for a hybrid initialization method

To verify the effectiveness of the hybrid initialization method and local search proposed procedure in this paper, it is compared with the random initialization method and method without local search, the comparison algorithms are as follows:

(1) MODGWO represents the method proposed in this paper.

(2) RI(Random Initialization) indicates that MODGWO uses a random initialization method to obtain the initial solution.

(3) WLS(Without Local Search) indicates that MODGWO does not use a local search procedure, but only a single tracking mode for searching.

Each method was independently run 20 times based on the expanded 45 test cases. The average values of IGD and HV were calculated based on the non dominated solution sets obtained by the three methods and calculated results are shown in Table 5. The results from Table 5 are shown as box plots in Figure 10.

The experimental results show that MODGWO can obtain better initial populations compared to RI which benefited from the effective diversity strategy introduced in their initialization. Compared with RI, MODGWO achieved a total of 40 minimum values on IGD and 38 maximum values on HV, indicating that the diversified initialization method proposed in this paper can obtain better quality and more Pareto solutions. And the 5 best IGD values obtained from RI method are all concentrated in smaller scale cases, indicating that for medium and large-scale cases, a diverse population initialization method is needed to obtain better solutions.

To assess the effectiveness of the proposed local search procedure, we compare WLS to MODGWO and RI for all problem instances. From Table 5 and Figure 10, it is seen that MODGWO and RI outperforms WLS for all data sets for IGD and HV metrics. The poor performance of WLS concerning the IGD and HV metrics can be attributed to the algorithm’s the absence of the local search procedure. On the whole, the introduced local search procedure plays a crucial role within the proposed methodology and has a substantial impact on enhancing the convergence performance of MODGWO. This indicates that the local search mechanism is integral to the overall effectiveness and efficiency of the MODGWO algorithm in finding optimal or near-optimal solutions to MRCFJSP.

3.4. Comparison with other algorithms

There have been many studies on flexible job shop scheduling under dual resource constraints for workers alone [34-36], as well as flexible job shop scheduling under AGV constraints alone [37-39]. However, there is a lack of research on flexible job shop scheduling under multi-resource constraints that comprehensively consider both workers and AGVs. Therefore, this paper chooses MOPSO proposed by Zhang et al. [40], MOJAYA proposed by Rylan H et al. [41], and NSGA II which is commonly used to solve multi-objective flexible job shop scheduling problems as comparative algorithms. These three algorithms have been often used to solve FJSP under resource constraints, which are most similar to the problem studied in this paper.

The experimental results for the two performance indicators have been organized in Table 6. Regarding the IGD metric detailed in Table 6, MODGWO demonstrates a notable superiority over its rival algorithms across all datasets. In 45 cases, MODGWO achieves 39 minimum values, while MOJAYA achieves 6 minimum values. It indicates that MODGWO exhibits superior exploration and exploitation capabilities. In comparison to MOPSO and NSGA-II, MOJAYA exhibits superior performance. Additionally, NSGA-II displays a competitive performance that is comparable to MOPSO. As for the HV metric, a similar trend is observed in Table 6, MODGWO demonstrates a significantly better performance as compared to other algorithms. In 45 cases, the MODGWO achieves 38 optimal values, and the MOJAYA achieves 7 optimal values. MOJAYA offers a better performance in terms of obtained metric values and the number of significantly better instances.

This indicates that both MODGWO and MOJAYA can obtain good solutions when solving MRCFJSP, with MODGWO yielding better results than MOJAYA. From Eq. (26), it can be seen that the Jaya algorithm drives the entire population to converge quickly to the global optimal solution through the optimal solution at each iteration, while from Eq. (25), the GWO algorithm drives the population to converge to the global optimal solution through the first three optimal solutions. This ensures the diversity of solutions and avoids the population from falling into local optimal solutions earlier. Therefore, MODGWO proposed in this paper combines the advantages of GWO algorithm and Jaya algorithm, which can effectively solve MRCFJSP.

To gain a deeper understanding of MODGWO's performance, especially regarding its convergence and diversity in tackling optimization challenges of varying scales, we compare its effectiveness based on small, medium, and large problem sizes. Benchmark instances are classified based on their jobs, machines, and operations. This detailed classification can more thoroughly evaluate the capabilities of the algorithms. The performance of the compared algorithms in terms of IGD, and HV metrics for selected instances are shown in Figure 11, Figure 12 and Figure 13.

Figure 11 depicts the comparative performance of the algorithms on smaller-scale instances. Specifically, as observed in Figure 11(a), MODGWO demonstrates superior performance in terms of the IGD metric for most instances, outperforming the other algorithms. Except for slightly lower performance than MOJAYA on MK_02_4_6 and MK05_3_2 instances, MODGWO achieves the best results on all other instances. Figure 11(b) shows the performance of the compared algorithms in terms of the HV metric. A similar trend is observed, as the IGD metric, with MODGWO performing poorly on MK02_4_6 and MK04_5_8 instances. MOJAYA surpasses NSGA-II and MOPSO in most instances while maintaining a competitive edge against MODGWO. Both NSGA-II and MOPSO exhibit comparable performance, suggesting a similar level of diversity and convergence in their respective Pareto solutions.

The comparison of medium scale instances is shown in Figure 12. As seen in Figure 12(a), in terms of IGD metrics, MODGWO performs well as compared to other algorithms for all instances except for MK06_6_5 and MK07_3_5 instances. MOJAYA is superior to NSGA-II and MOPSO across all instances. Regarding the HV metric depicted in Figure 12(b), it is evident that, aside from the instances of MK06_6_5, MK07_3_3, and MK10_9_5, MODGWO performs exceptionally well in all instances. NSGA-II and MOPSO also show good performance, and the difference between them was not very large, slightly worse than MOJAYA.

Figure 13 presents a comparison of algorithms on large-scale instances. Specifically, Figure 13(a) demonstrates that MODGWO surpasses its competitors in terms of the IGD metric across all instances except for MK11_3_3, and MK14_9_5, suggesting that the solutions in the Pareto Front are evenly distributed. In Figure 13(b) it is observed that the performance of MODGWO demonstrates a similar trend as compared to the IGD metric. MODGWO exhibits a superior performance compared to the other three algorithms, while the performance of MOJAYA, NSGA-II and MOPSO is not significantly different on the HV metric. This confirms that MODGWO efficiently balances maintaining diversity in the search process with improving the convergence of solutions along the Pareto Front.

The instance is solved using MODGWO proposed in this paper and the other three algorithms in Section 3.4, and the performance indicators are shown in Table 11. From Table 12, it can be seen that both indicator values of MODGWO are superior to the other three algorithms, indicating that MODGWO still has advantages compared to other comparative algorithms in solving enterprise instances.

In order to clearly display the results obtained by different algorithms, Table 13 lists the Pareto solution sets obtained by four algorithms. The first value in each Pareto solution represents the maximum completion time , i.e. makespan, and the second value represents the maximum continuous working hours of all workers(MCWH). It can be seen from the table that MODGWO obtained the highest quality solution.

The set of Pareto solutions in Table 13 is plotted as a scatter plot as shown in Figure 16, from which it can be clearly seen that the Pareto solutions obtained by MODGWO are more concentrated in the lower-left part of the coordinate system compared to those obtained by other comparative algorithms, which indicates that the quality of the results obtained by MODGWO is better. MODGWO does not obtain as many non-dominated solutions as other algorithms, but its solutions are of higher quality.

3.5. Scheduling Gantt Chart

To demonstrate the internal scheduling situation of the production workshop, select a scheduling result from MK05_3_4 and draw its scheduling Gantt chart. This case involves 15 jobs, 4 machines, 3 workers and 4 AGVs, where the set of machines that can be operated by the workers is shown in Table 7, and the transportation time of AGVs between each machine and the LU area is based on the Layout1 proposed by Bilge and Ulusoy[42]. The pareto optimal solution for MK05_3_4 obtained by MODGWO is [501, 32], [381, 54], [427, 37], [404, 48], [411, 44], [451, 36], [412, 40] and the Gantt chart obtained by selecting [501, 32] is shown in Figure 14. The horizontal axis in Figure 14 represents the processing time of the operations on the corresponding machine or the time of AVG transporting the jobs, in the vertical axis, M1-M4 represents the processing machine and corresponding processing operations, and AGV1-AGV4 represents the jobs transported by the corresponding AVG. If adjacent processes of the same job are processed on the same machine, AGV transportation is not required, so the transportation process is not shown in the figure. For example, operations O₄₃ and O₄₄ are both processed on machine M₄, AGV transportation of O₄₄ is not required. The numbers in the rectangular boxes corresponding to M1-M4 are in the form of 'x_y_z', where x represents the number of the job, y represents the operation number of the job, and z represents the number of the worker who processed the operation, e.g., ‘3_1_3’ in the first rectangular box of M₄ represents that the first operation of job 7 is processed by worker 3 on M₄. The numbers in the rectangular boxes corresponding to AGV1-AGV4 are in the form of 'x_y', which has the similar meaning as the numbers in the rectangular boxes corresponding to M₁-M₄, e.g., ‘4_1’ in the first rectangular box of AGV2 indicates that the first operation of job 4 is transported by AGV2. From Figure 14, it can be seen that scheduling needs to satisfy the multi-resource constraints of machines, operations, workers and AGVs simultaneously. For example, '5_3_1' on M₄ has to wait for AGV2 to transport '5_3' to M₄, while also satisfying the processing completion of operation O₅₂. At the same time, it also needs to wait for worker 1 to complete the previous operation ('10_2_1' on M₁).

At the same time, we conducted a scheduling scheme that only satisfies the single objective constraint of makespan, and compared it with the scheduling scheme that considers the continuous working hours of workers under the dual objective optimization, as shown in Figure 15. Figure 15 shows the scheduling Gantt chart from the perspective of workers and the text in the rectangular box is in the form of "x_y", which has the same meaning as the scheduling plan from the perspective of AGV in Figure 14. Figure 15(a) shows a single objective optimization scheduling Gantt chart considering only makespan, Figure 15(b) shows a dual objective optimization scheduling Gantt chart considering both makespan and worker continuous working hours, and Figure 15 (b) is also the scheduling Gantt chart of Figure 14 from the worker's perspective. By comparison, it can be found that although the makespan in Figure 15(a) is 359, each worker has a long duration of continuous work, which can lead to a decrease in work efficiency and an increase in production safety risks. Although the makespan in Figure 15(b) is 501, the longest continuous working time of workers is 32, which is from '3_7' to '5_6' of the worker W₂ in the figure. Moreover, the working time arrangement of each worker is relatively reasonable, which can ensure the completion time while allowing workers to have sufficient rest and adjustment when there are no additional work tasks, which is more in line with the modern people-oriented production concept.

3.6. Examples of Agricultural Machinery and Equipment Production

The products of agricultural machinery and equipment manufacturing enterprises involve in various aspects of agricultural production such as cultivation, planting, fertilization, irrigation, harvesting, transportation, and processing. The main processing techniques for its product components include "cutting", "turning", "boring", "milling", "drilling", "rolling", "inserting", and “heat treatment”. Due to its wide range of product models, many of which are personalized for customers, the product structure, standards, and process routes are complex and varied. As a result, it is very difficult to prepare the operation plan by manually operating the computer.

In this section, the actual production data of a large-scale seeding equipment manufacturing enterprise is used as an example, and the proposed MODGWO algorithm and other comparative algorithms are utilized for solving this problem, in order to verify and compare the effectiveness and performance of the proposed algorithm in solving practical problems. From the company's products, 5 jobs with 18 operations are selected to be processed on 11 machines, the number of AVGs is 5. The detailed information of the machines is shown in Table 8, and the detailed information of the operation processing is shown in Table 9. The transportation time of AGVs between each machine and the LU area is shown in Table 10, and the set of machines that can be operated by the workers is shown in Table 11.

4. Discussion

The research presented in this paper addresses the MRCFJSP by proposing a MODGWO algorithm. This problem is particularly complex due to the simultaneous consideration of machine, worker, and AGV constraints, aiming to optimize both the makespan and workers' continuous working time.

The MODGWO algorithm introduced in this study incorporates several innovations to effectively tackle the MRCFJSP. Firstly, an encoding/decoding method and a mixed population initialization strategy are designed to create a diverse set of initial solutions. This initialization process is crucial as it impacts the subsequent search performance of the algorithm. The results indicate that MODGWO, with its effective diversity strategy, outperforms a random initialization method, achieving better quality and more Pareto solutions.

Furthermore, the MODGWO algorithm integrates a solution updating mechanism tailored for discrete scheduling problems. This mechanism is inspired by the Grey Wolf Optimizer but adapted to include cross-operation principles from genetic algorithms. By incorporating the Jaya algorithm's principles of continuous improvement, MODGWO further enhances its solution quality. The integration of these algorithms allows MODGWO to effectively explore and exploit the search space, leading to superior performance compared to other metaheuristic algorithms such as MOJAYA, MOPSO, and NSGA-II.

The experimental results demonstrate MODGWO's superiority across various performance indicators, including the Inverted Generational Distance (IGD) and Hypervolume (HV). MODGWO achieves the best results in most instances, particularly in medium and large-scale problems, highlighting its robustness and scalability. In addition to quantitative results, this study also presents scheduling Gantt charts to visually demonstrate the internal scheduling situation of the production workshop. These charts illustrate the differences between single-objective and dual-objective optimization approaches, emphasizing the importance of considering workers' continuous working time alongside makespan. By doing so, MODGWO contributes to improving work efficiency and reducing production safety risks.

Despite the promising results, several areas for future research can be identified. First, exploring additional constraints and objectives within the MRCFJSP framework could further enrich the problem's complexity and practical relevance. For instance, incorporating energy consumption or machine maintenance schedules could provide a more holistic view of the scheduling problem. Second, hybridizing MODGWO with other advanced metaheuristic algorithms or machine learning techniques may yield even better performance. Lastly, implementing MODGWO in real-world industrial settings and evaluating its practical impact would be a valuable contribution to the field.

In conclusion, this paper proposes a novel MODGWO algorithm for solving the MRCFJSP, demonstrating its effectiveness and superiority over existing algorithms. By considering both makespan and workers' continuous working time, MODGWO contributes to more efficient and sustainable production scheduling. The research findings not only advance the theoretical understanding of complex scheduling problems but also offer practical insights for industrial applications.

5. Conclusions

Agricultural machinery and equipment production is greatly affected by market demand and seasons, and is characterized by fast production conversion and strong dynamics. The method of flexible job shop scheduling fully utilizes the flexible characteristics of various resources, which can greatly improve the production efficiency of enterprises when applied to agricultural machinery and equipment production enterprises. This paper proposes a MODGWO algorithm that integrates hybrid initialization strategy and multi neighborhood local search for flexible job shop scheduling problem under multi resource constraints of integrated workers and AGVs. On the basis of in-depth analysis of the actual operation of various resources, a flexible job shop scheduling model was constructed, which includes four production factors: jobs, machines, AGVs, and workers. The effectiveness of the hybrid initialization method and the multi-neighborhood local search is verified by comparison based on several instances. The results compared with various advanced algorithms show that the proposed algorithm can obtain better scheduling solutions when solving this type of scheduling problem.

The method proposed in this paper has been validated in the actual manufacturing of seeding locomotives. The scheduling scheme proposed in this paper not only improves production efficiency, but also balances the two goals of completion time and workers' continuous working hours. In the future, based on current research, we will further construct a multi-objective multi resource constrained scheduling model, study efficient multi-objective optimization methods, as well as utilize the current popular deep reinforcement learning method for solving.