Multi-Objective Parallel Machine Scheduling with Eligibility Constraints for the Kitting of Metal Structural Parts

: This paper studied a class of coupling problems of material assignment, part nesting, kit delivery and parallel machine scheduling. The aim of this paper was to solve the scheduling problem of metal structural part processing and welding assembly with eligibility constraints. A two-stage mixed-integer programming model was constructed. The eligibility constraints took into account the material type of parts and nesting. The objectives were to minimize the makespan, maximize material utilization and minimize the kit delivery metrics (kitting time and numbers of earliness and tardiness of kits). A hierarchical optimization approach was proposed. The scheduling model was solved by using the Gurobi solver in the ﬁrst stage, and the results were used to constrain the second stage. The second stage of the scheduling model was solved using an improved multi-objective genetic algorithm. Due to the strong coupling relationships among the sorting of parts, the sorting of each proﬁle and the sorting of each material, a hybrid encoding and decoding mode was designed for part sorting with eligibility constraints. Finally, the proposed scheduling approach was applied to actual production cases. The data showed that when the number of components exceeded 300 (the number of parts was about 1500), the material utilization reached 95%. Choosing a suitable number of machines, machine utilization reached 90%. The results demonstrated the effectiveness of the proposed scheduling model and algorithm.


Introduction
In the context of Industry 4.0, mass customization is developing rapidly, and customers can choose customized products with different configurations. Mass customization is a modern manufacturing mode that combines mass production with customization [1]. Mass customization can help balance the contradiction between order differentiation and manufacturing costs [2] and enables greater flexibility in manufacturing systems [3]. The key to achieving mass customization and rapid response is providing diverse product configurations to meet customer requirements and shortening the order delivery time [4]. It is an important way to attract customers and improve the competitiveness of manufacturing enterprises. However, it often adds complexity to the production scheduling process. Production scheduling is the core of modern manufacturing, and it aims at balancing customer requirements and manufacturing resources, as well as improving quality, increasing efficiency and reducing costs and consumption. Metal structural parts are widely used in aerospace products, ships, automobiles, transmission towers and other equipment. Generally, their production process includes part nesting, cutting processing and welding assembly. Due to the differentiation of orders, there are many types of metal structural parts and small batches, resulting in machine load imbalance, low material utilization, poor kitting and a backlog of work-in-progress in the production process.
The production process of metal structural parts could be regarded as a kind of parallel machine scheduling problem (PMSP). The classical PMSP assumes that there are n jobs, and each job must be processed by one of m parallel machines. The objectives of the PMSP are to minimize the makespan, the tardiness, etc. However, in the actual production process, part nesting makes the metal structural parts not independent. In addition, this is necessary for the close connection between machining and welding assembly. The absence of any one part leads to the retention of all the parts belonging to the same kit or component, and subsequent follow-ups cannot be carried out during assembly welding. Due to the variety of orders and the fluctuation in the price of metal materials, it is required that suitable scheduling must comprehensively consider the above factors to improve the machining and assembly efficiency of the metal structural parts and reduce the waste of materials. Therefore, eligibility constraints need to be considered, i.e., the material-type constraints of each structural part and the nesting constraints. Part nesting is necessary to optimize material utilization. When welding assembly is performed, it is necessary to ensure that all the parts of the same component or product are completed as close to the same time as possible, a process called kitting. Compared with the classical PMSP, this makes problem modeling and solution more complicated and puts forward higher requirements for the comprehensive performance of the scheduling approach.
Since McNaughton studied the PMSP for the first time [5], research on models and algorithms for the PMSP has achieved rich theoretical results [6,7]. Due to the wide range of practical applications, many researchers have focused on the PMSP in practical workshops [8][9][10]. Workshop scheduling needs to comprehensively consider the cooperation and combination of multiple production links, such as machining and assembly, making the scheduling problem more complex. Meng et al. analyzed a metal-cutting plant, considered processing times, setup times, pickup times, different delivery times and machine eligibility constraints, and addressed a bi-objective cutting parallel machine scheduling problem aiming to minimize the total makespan and the total tardiness [11]. Zhang et al. researched the integrated scheduling problem of production and delivery on parallel batch processing machines with non-identical capacities in different locations in cloud manufacturing and assumed the job release time, batch sizes, processing times and delivery times. A mixed-integer programming (MIP) model was presented to solve this problem, and the objective was to minimize the total service completion time [12]. Hatami et al. first described the distributed parallel machine and assembly scheduling problem (DPMASP) with eligibility constraints, aiming to optimize the scheduling of two manufacturing stages, namely, the production stage and the assembly stage. The production stage consisted of a set of identical factories with a set of unrelated parallel machines in one factory. In the assembly stage, the assembly machine was assigned to assemble all products. The proposed eligibility constraints indicated conditions whereby a job could only be assigned to an eligible factory. Minimizing the makespan of the products in the assembly stage was considered as the objective. The scheduling model was solved using CPLEX and Gurobi solver [13]. Amallynda et al. discussed an expansion of the DPMASP to a new parallel machine scenario and presented new meta-heuristics for scheduling problems. A set of products was manufactured through the assembly of several parts, each fabricated in the factories during the production stage. A product could be assembled in the assembly stage only after all the jobs in the production stage were completed. The objectives were to minimize the mean flow time and the number of tardy jobs [14]. The above studies addressed the PMSP of production and assembly or delivery with different optimization objectives. However, the scheduling objectives were mostly limited to part-level metrics; the focus was only on the point in time when the last job was finished, ignoring the actual duration of the jobs within the same component. The kit delivery efficiency was often overlooked, that is, the difference in the completion time of multiple parts of the same component was not considered. It has to be noted that the longer it takes to perform kitting, the more the work-in-process stock increases. For the case of kit delivery, Hsu et al. researched the customer order scheduling (COS) problem, where each order consisted of the set of jobs that had to be shipped as one batch at the same time. A new dispatching rule was proposed for COS in job shop scheduling, in order to reduce the total time it took to complete all jobs within the same order [15]. Chen et al. studied a synchronized scheduling problem of production and shipment with the objective to ensure that all products belonging to the same order were simultaneously completed (at least as close as possible). The proposed scheduling method was applied to a two-stage assembly flow shop scheduling [16]. Hu et al. developed a two-stage framework for the lot sizing and scheduling of the production activities at a kitting facility aiming to solve a kitting-specific production planning problem in a manufacturing setting [17]. According to the research background of this paper, it is also necessary to perform metal structural part nesting to improve the utilization rate of materials, which is a more realistic and complex PMSP. The part nesting problem [18,19], in general, has received increasing attention from researchers as an independent optimization problem [20,21]. However, the coupling problem of part nesting and parallel machine scheduling has rarely been considered. Some researchers studied the scheduling problem of cutting plates with parallel machines and solved it using branch and bound heuristic algorithms [22,23]. The above studies from the literature provide a basis for this paper, but there are few research studies on PMSP that simultaneously consider material assignment, material utilization and kit delivery.
In the production process of metal structural products, the complexity of the scheduling problem is determined by the differences in the ordered-product configuration and the resource constraints, making the scheduling problem a typical multi-objective optimization problem [24]. Various multi-objective optimization methods are reported in the existing literature, and we here summarize them [25,26]. (1) The weighted approach: The decision maker assigns different weights based on the importance of each objective. Multiple objectives are merged into a single-objective problem. (2) The hierarchical approach: The objectives are sorted based on their importance. Then, hierarchical solving is performed, and the optimal solution obtained in the previous level is used in the next level as a constraint. (3) The posterior approach: It aims to obtain a set of non-dominated solutions (also known as Pareto solutions). The decision maker selects an eligible optimal solution from this set of solutions. (4) The interactive approach: The decision maker expresses their preferred target value for each objective in the solution process. The overall computing mode extracts a feasible trade-off solution based on the objective equation and the decision maker's preference. Many scholars have combined the posterior approach with meta-heuristic algorithms to solve multi-objective optimization problems. Meta-heuristic algorithms characterized by simple structure, flexibility, gradient-free information, easy implementation and a strong capability of circumventing local constraints can achieve a global optimum [27]. Currently, multi-objective meta-heuristics include particle swarm optimization [28], the immune genetic algorithm [29], the bat algorithm [30], etc. Deb K. proposed a fast non-dominated sorting genetic algorithm with elite strategy (NSGA-II) [31], which is widely used to solve multi-objective workshop scheduling problems because of its less computing time and good solution quality. Chen et al. applied the non-dominated sorting genetic algorithm to the design of non-compact flow shop scheduling and successfully solved the multi-objective optimization problem considering process connections [32]. Wu et al. designed a flexible job shop scheduling problem (FJSP) model with energy considerations in line with the production field. They invented an adaptive population non-dominated sorting genetic algorithm to solve the scheduling problem. The results showed that the proposed algorithm was efficient and powerful when dealing with the multi-objective FJSP that included energy consumption [33]. Regarding multi-objective PMSP, Zheng et al. investigated the PMSP of matching orders and factories; the objectives included minimizing the total tardiness and the extra quality of service costs, using NSGA-II to efficiently solve large-scale instances [34]. Although the designed algorithm considers the eligibility constraint relationship between orders and machines of a factory, orders are independent and non-decomposable minimum job units, which cannot solve scenarios containing kitting requirements. Liu et al. studied a complex product assembly line scheduling problem considering multi-skilled worker assignment and parallel team scheduling and took the maximum completion time and the imbalance degree of team workload as the optimization objectives. The improved strategies were proposed on the basis of NSGA-II, and a hybrid coding method considered worker assignment and task order [35]. They focused the completion time due to the worker assignment scheme but did not consider the matching relationship between workers and jobs. Therefore, the proposed approaches were difficult to be applied to the optimization problem containing the coupled material assignment and the job scheduling studied in this paper. Bandyopadhyay et al. addressed a PMSP with three objectives by improved NSGA-II. Experimental results showed that the proposed algorithm also works effectively [36]. However, this algorithm was only applied to the classical PMSP scheduling problem and could not solve the scheduling problem constrained by the material assignment and the part nesting. Therefore, current research results on multi-objective optimization algorithms applied to the production scheduling did not yet have an effective algorithm to solve the multi-objective optimization problem of production scheduling in metal structural parts for kitting delivery considering both the material assignment and the part nesting constraints, as proposed in this paper.
Due to the particularity of the production process of metal structural parts, there is still a lack of effective scheduling models and algorithms to address such complex scheduling problems. Therefore, this paper studied a multi-objective scheduling problem of the kit delivery of metal structural parts with eligibility constraints. The scheduling process was divided into two stages to solve material assignment and part nesting, as well as kit delivery for welding assembly. Aiming at solving the two-stage scheduling problem, a mixed-integer programming scheduling model with eligibility constraints was established, and a hierarchical multi-objective optimization approach was designed. Briefly, the first stage of the scheduling model was solved using Gurobi, and the results of the first stage were used to constrain the second stage. The second stage of the scheduling model was solved using an improved multi-objective genetic algorithm. Furthermore, according to the scheduling characteristics of the metal structural parts, a hybrid encoding-decoding mode with eligibility constraints was designed. Finally, the validity of the model and algorithm was verified using the production scheduling data of the metal structural parts of bus frames.
The remainder of this article is arranged as follows: In Section 2, we describe the scheduling problem and propose a mixed-integer programming model. The design of an improved multi-objective genetic algorithm to solve the corresponding scheduling model is presented in Section 3. In Section 4, we conduct a numerical study. Finally, we summarize the contributions of this paper and suggest possible directions for future research in Section 5.

Problem Description
A typical product structure of metal structural parts is shown in Figure 1a. It can generally include multiple levels: parts, kits or components, and products. Welding is performed at each level to form sub-assemblies or assemblies and finally form products. Customized product configurations lead to the difference in the sub-assemblies at all the levels of the product and part processing. A machining workshop produces metal structural parts for multiple ordered products. Various types of materials are used in production. The type of a bundle of material corresponds to a cross-sectional shape and size. A bundle of material includes multiple identical profiles (the profile is a straight bar with a certain cross-sectional shape and size). All metal structural parts are cut from profiles one by one, and different parts could need to share the same material or even belong to the same profile due to part nesting. Therefore, the material type of the parts in the product and part nesting impose eligibility constraints on production scheduling. Before production, the material types of all parts are selected, and part nesting is performed to obtain the number of single profiles for each material, as shown in Figure 1b time of the part and the results of part nesting. The production process is divided into four steps, as shown in Figure 1c: (1) each material is sequentially allocated to buffers according to its total processing time and sequence. (2) Each profile in the buffer enters the machine in sequence according to profile sorting. (3) The single profile is cut into parts in sequence according to part sorting. (4) After the parts belonging to the same kit or product undergo kitting, they are delivered to the workshop.
Machines 2022, 10, x FOR PEER REVIEW 5 of 24 material types of all parts are selected, and part nesting is performed to obtain the number of single profiles for each material, as shown in Figure 1b. The processing time of a single profile and the total processing time of each material are obtained based on the known processing time of the part and the results of part nesting. The production process is divided into four steps, as shown in Figure 1c: (1) each material is sequentially allocated to buffers according to its total processing time and sequence. (2) Each profile in the buffer enters the machine in sequence according to profile sorting. (3) The single profile is cut into parts in sequence according to part sorting. (4) After the parts belonging to the same kit or product undergo kitting, they are delivered to the workshop.

Notation
According to the problem description, the mathematical symbols are as follows: Parameters: Processing time of the g-th profile of material h;

Notation
According to the problem description, the mathematical symbols are as follows: Parameters: i Index of machines: i ∈ M; r Index of orders: r ∈ O; k Index of kits (or components): k ∈ C r ; j Index of parts: j ∈ J; h Index of material type: h ∈ P; g Index of profiles that belong to material h : g ∈ P h ; M Set of machines: Set of kits (or components) that belong to order r : C r = {1, 2, · · · , N C r ; J Set of parts: J = {1, 2, · · · , N J ; P Set of material types: P = {1, 2, · · · , N P }; P h Set of profiles of material h : P h = {1, 2, · · · , N P h . z ihgg Equals 1 if the g-th profile, which belongs to material h, is adjacently processed before the g'-th one by machine i and 0 otherwise; w ihgg Equals 1 if part j, which belongs to the g-th profile of material h, is adjacently processed before j' by machine i and 0 otherwise.

Formulation
Based on the above analysis, the researched scheduling problem of manufacturing metal structural parts in this paper is detailed below.
Each material h is distributed to a machine i; the processing sequence of each profile g is determined, and each part j is assigned to a profile, so that all the parts that belong to one kit can undergo kitting and be delivered as closely to each other as possible. On the basis of satisfying the related constraints, multiple scheduling objectives are optimized. These objectives include the minimization of the makespan, the minimization of the material utilization ratio and objectives related to kitting.
The relative assumptions made are listed as follows: (1) Excluding the transport time of material and profiles; (2) Each material is only transported once, and its quantity can meet the cutting requirements of all parts using the material; (3) After each profile is processed, the handling time of off-cut is not counted; (4) The cutting edge width does not count; (5) The processing time of all parts is known and determined; (6) Each part can only be processed by one machine at most once; (7) The processing capacity of each machine is the same; (8) Each machine can only process at most one part at a time; (9) Preemption is not allowed; that is, any started operation cannot be stopped before completion; (10) All parts and machines are available at the beginning; machine failures and blocks are not considered.

Scheduling Model
In this section, an extended parallel machine scheduling model is established based on multi-objective optimization. The scheduling problem is separated into two stages. The first stage is to schedule all materials onto a set of parallel production machines M. Meanwhile, part nesting is carried out to determine the part set belonging to the same profile. In the second stage, all parts are scheduled and sorted according to the eligibility constraints to obtain the optimal processing sequence of materials, profiles and parts.
In actual production, the makespan and material utilization are given priority, as they directly affect the cost of the enterprise. Then, they need to meet the requirements of kitting to guarantee order delivery. Therefore, the optimization objectives corresponding to two-stage scheduling have priority. The optimization objectives in the first stage have higher priority than those in the second stage. We propose a hierarchical optimization approach. The results of the first stage are used to constrain the second stage. Hence, it is necessary to first achieve the objectives related to the makespan and part nesting and lastly optimize the objectives related to kitting.
Based on the notation defined, we propose the mathematical model in the next sections.

Material Assignment
The purpose of material assignment is to minimize the makespan. This is related to the number of machines and the processing time of each material. The mathematical model of material assignment is as shown below.
In the formulation, objective (1) is to minimize the maximum makespan of machines M in the first stage. Constraint (2) imposes that each material h is only processed by one machine and can only be processed once. Constraint (3) enforces that at least one material is processed on one machine. Equation (4) calculates the total processing time of material h. Constraint (5) specifies that decision variables x ih and ϕ hj are binary.

Part Nesting
The material utilization is expressed as a quantity. Here, only the influence of the part length on the nesting was considered, which was a one-dimensional nesting optimization problem. Hence, to optimize the quantity of each material, the mathematical model of material h is as shown below.
In the optimization model, objective (6) accounts for the minimization of the total quantity of material h required in the first stage. Constraint (7) imposes that each part j is only assigned to one profile g that belongs to material h and can only be assigned once. Constraint (8) imposes that the total length of parts distributed on material h cannot exceed the length of material h. Constraint (9) defines that decision variable ψ hgj is binary.
Formulas (1)-(9) above are the optimization objectives and constraints relative to the first stage. We carry out the distribution of each material and the part layout and obtain the optimal scheme. Then, we establish a mathematical model of optimization objectives related to kitting.

Optimization of Kitting
The second stage optimizes three evaluation indicators at the same time: minimizing the maximum kit completion time ( f 1 3 ), minimizing the number of earliness kits ( f 2 3 ) and minimizing the number of tardiness kits ( f 3 3 ). (1) Minimizing the maximum kit completion time: The completion time of kit k is the difference between the completion time of the last part of the kit and the starting time of the first part.
(2) Minimizing the number of earliness kits: Earliness kits are components whose last part's completion time is shorter than the starting time of the kit's delivery window. The number of earliness kits is the sum of the number of components completed ahead of schedule.
(3) Minimizing the number of tardiness kits: Tardiness kits are components whose last part's completion time is longer than the completion time of the kit's delivery window. The number of tardiness kits is the sum of the number of all delayed components.
Constraint (14) indicates that each part can only be processed by one machine and can only be processed once. Constraint (15) imposes that each machine needs to process parts. Constraints (16) and (17) indicate that material h is adjacently processed before h' by machine i. Constraints (18) and (19) denote that profile g belongs to material h and is adjacently processed before g' by machine i. Constraints (20) and (21) indicate that two parts on the g-th profile of material h are processed one after the other by a machine. Constraint (22) indicates that one machine can only process one part at a time. Constraint (23) is used to calculate the start of the processing time of each material. Constraint (24) is used to calculate the start of the processing time of each part. Constraint (25) indicates that every part is processed. Constraints (26) and (27) indicate that all processing is non-pre-emptive and that processing cannot be interrupted once it starts. Constraint (28) specifies that the decision variables, including φ ij , y ihh , z ihgg and w ihgjj , are binary.

Algorithm Design
The scheduling model consists of two stages: material assignment and part nesting constitute the first stage, and profile sorting and part processing sorting form the second stage. The two stages include three optimization models. From the perspective of the scale and computational complexity of the three optimization models, the first model is smaller than the second one, and the second one is smaller than the third one. Hence, we proposed a hierarchical optimization approach. Taking into account the accuracy and efficiency of the solution, model 1 and model 2 are solved with an accurate method that uses Gurobi as the solver in the first stage. Gurobi can solve many types of problems, including mixed-integer programming, and is the best and fastest linear programming/quadratic programming solver in the world today because of its wide range of problem types and high solution efficiency. At the same time, Gurobi can be easily integrated into third-party platforms such as MATLAB, Python or Net, and is more capable of handling a wide range in complex scientific research and applications, so Gurobi was chosen as the first-stage optimization solver. According to the solution results of the first stage, the number of profiles belonging to each material, the index of all profiles of each material and the part index of a profile can be obtained, and there is a fixed subordination relationship among them. The subordination relationships are eligibility constraints. The results of the first stage were used for eligibility constraints in the second stage. Model 3 belongs to the second stage, is a multi-objective optimization model, and the solution scale is large, so we proposed an improved multi-objective meta-heuristics algorithm based on the NSGA-II. Due to eligibility constraints, there are strong coupling relationships among the sorting of parts, the sorting of each profile and the sorting of each material. In order to express these sorting relationships, part-based coding is designed; in order to improve the problem of insufficient crossover due to too-long chromosomes, a segmented crossover method is designed. In order to improve the iterative quality of the population, the non-dominated sorting method is improved, and the elite retention strategy is used to improve the iterative quality of the population. The flow chart of the improved NSGA-II algorithm is shown in Figure 2.

Coding
This paper designs an integer coding method based on parts that can effectively express the sorting of materials and their assignment to each machine, the order of each material among the materials and the sorting of parts and their assignment to each profile. According to the optimization results of the model in the first stage, the subordination relationships among parts, profiles, materials and machines can be obtained. This kind of affiliation can be expressed with a tree diagram, as shown in Figure 3.
In order to improve the iterative quality of the population, the non-dominated sorting method is improved, and the elite retention strategy is used to improve the iterative quality of the population. The flow chart of the improved NSGA-II algorithm is shown in Figure

Coding
This paper designs an integer coding method based on parts that can effectively express the sorting of materials and their assignment to each machine, the order of each material among the materials and the sorting of parts and their assignment to each profile. According to the optimization results of the model in the first stage, the subordination relationships among parts, profiles, materials and machines can be obtained. This kind of affiliation can be expressed with a tree diagram, as shown in Figure 3.

Coding
This paper designs an integer coding method based on parts that can effectively express the sorting of materials and their assignment to each machine, the order of each material among the materials and the sorting of parts and their assignment to each profile. According to the optimization results of the model in the first stage, the subordination relationships among parts, profiles, materials and machines can be obtained. This kind of affiliation can be expressed with a tree diagram, as shown in Figure 3.  First, the initial reference code is designed; this is an integer sequence [1, 2, · · · , N J , where N J is the sum of the number of all parts, that is, the length of the chromosome. In the iterative process, the reference codes are reordered according to different chromosome sequences to form a feasible job code. A feasible encoding is [649251110121783]. Since all the permutations of the chromosome are merely the reordering of the original reference code, the number and type of certain elements in the reference code are not changed, nor are new elements formed. Therefore, this coding method does not generate illegal chromosomes. This coding method avoids reducing the algorithm search efficiency due to the identification of illegal chromosomes.

Decoding
The genetic-locus value in the designed chromosome represents the part number, and the parts are sorted according to the sequence of the genetic loci. The other parts on the profile are indexed according to the part number on the genetic locus to form the processing sequence of the parts on the profile. The other profiles of the material category are indexed according to the number of profiles to form the processing sequence of the profile. The other materials on the machine are indexed according to the number of materials to form the processing sequence of materials. The decoding process is shown in Figure 4. According to the above method, the sequencing of all the genetic loci is completed to form a sequencing scheme.
are new elements formed. Therefore, this coding method does not generate illegal chro-mosomes. This coding method avoids reducing the algorithm search efficiency due to the identification of illegal chromosomes.

Decoding
The genetic-locus value in the designed chromosome represents the part number, and the parts are sorted according to the sequence of the genetic loci. The other parts on the profile are indexed according to the part number on the genetic locus to form the processing sequence of the parts on the profile. The other profiles of the material category are indexed according to the number of profiles to form the processing sequence of the profile. The other materials on the machine are indexed according to the number of materials to form the processing sequence of materials. The decoding process is shown in Figure 4. According to the above method, the sequencing of all the genetic loci is completed to form a sequencing scheme.
[ 6 4 9 2 5 11 10 12 1 7 8   Combining the processing time of parts, the processing time of each profile and the processing time of each material, a feasible scheduling plan is finally formed. After decoding according to the chromosome generated above, Gantt charts are formed for material assignment scheduling and profile cutting scheduling. Moreover, the Gantt chart of part sorting is shown in Figure 5, assuming that the first kit contains parts 1, 2, 10 and 11, the second one contains parts 3, 5, 6 and 9 and the third one contains parts 4, 7, 8 and 12.
We can obtain the Gantt chart of kitting according to the Gantt chart of part sorting, as shown in Figure 6, and calculate the objective function related to the complete kit.  Combining the processing time of parts, the processing time of each profile and the processing time of each material, a feasible scheduling plan is finally formed. After decoding according to the chromosome generated above, Gantt charts are formed for material assignment scheduling and profile cutting scheduling. Moreover, the Gantt chart of part sorting is shown in Figure 5, assuming that the first kit contains parts 1, 2, 10 and 11, the second one contains parts 3, 5, 6 and 9 and the third one contains parts 4, 7, 8 and 12. We can obtain the Gantt chart of kitting according to the Gantt chart of part sorting, as shown in Figure 6, and calculate the objective function related to the complete kit.
mosomes. This coding method avoids reducing the algorithm search efficiency due to the identification of illegal chromosomes.

Decoding
The genetic-locus value in the designed chromosome represents the part number, and the parts are sorted according to the sequence of the genetic loci. The other parts on the profile are indexed according to the part number on the genetic locus to form the processing sequence of the parts on the profile. The other profiles of the material category are indexed according to the number of profiles to form the processing sequence of the profile. The other materials on the machine are indexed according to the number of materials to form the processing sequence of materials. The decoding process is shown in Figure 4. According to the above method, the sequencing of all the genetic loci is completed to form a sequencing scheme.  Combining the processing time of parts, the processing time of each profile and the processing time of each material, a feasible scheduling plan is finally formed. After decoding according to the chromosome generated above, Gantt charts are formed for material assignment scheduling and profile cutting scheduling. Moreover, the Gantt chart of part sorting is shown in Figure 5, assuming that the first kit contains parts 1, 2, 10 and 11, the second one contains parts 3, 5, 6 and 9 and the third one contains parts 4, 7, 8 and 12. We can obtain the Gantt chart of kitting according to the Gantt chart of part sorting, as shown in Figure 6, and calculate the objective function related to the complete kit.

Fast Non-Dominated Sorting
In a multi-objective optimization problem, in order to compare the pros and cons of two individuals, two measurement indicators-the non-dominated ranking level and the crowding distance-are commonly used. If none of the goals of individual p are worse than those of individual q and at least one goal of individual p is better than those of individual q, then individual p dominates individual q. First, all individuals are ranked according to non-dominated relationships. Obviously, the individuals at the first front end are dominant, and the individuals at the second front end are dominated by the first front-

Fast Non-Dominated Sorting
In a multi-objective optimization problem, in order to compare the pros and cons of two individuals, two measurement indicators-the non-dominated ranking level and the crowding distance-are commonly used. If none of the goals of individual p are worse than those of individual q and at least one goal of individual p is better than those of individual q, then individual p dominates individual q. First, all individuals are ranked according to non-dominated relationships. Obviously, the individuals at the first front end are dominant, and the individuals at the second front end are dominated by the first front-end individuals. Individuals are assigned to different front ends through sorting, and for individuals at the same front end, pros and cons can be determined by comparing the crowding distance until all individuals are sorted.
The evaluation of individual fitness and selection is performed using two methods: fast non-dominated sorting and crowding distance. Each individual in the population is assigned a non-dominated rank through rapid non-dominated sorting. For the population of pop individuals, it is assumed that n i individuals are dominated by the i-th individual and that the individuals dominating others form set s i . The process of rapid non-dominated sorting is described below.
Step 2: Determine the s k (k ∈ {1, 2, · · · , q}) of each individual on level 1; judge whether the n t−1 of each individual t in s k equals 0; set the individuals with n k−1 = 0 to level 2.
Step 3: Repeat Step 2 until all individuals are set to the corresponding rank. Then, the individuals on the same level are subjected to a fitness comparison based on the crowding distance.
The crowding distance of each individual is calculated as described below.
Step 1: Find the non-dominated individuals of the same front end and place them in a matrix W.
Step 2: Based on the c-th objective function, arrange matrix W in descending order, where c = 1.
Step 3: Let the crowding distance between the first and last individuals in matrix W be infinite.
Step 4: Calculate the crowding distance of all individuals in matrix W. The calculation formula is shown below.
where d(max) represents the maximum value of the c-th objective function; d(min) represents the minimum value of the c-th objective function; and d(f + 1) and d(f − 1) represent the c-th objective function value of the (f − 1)-th individual and the (f + 1)-th individual after sorting, respectively.
Step 5: c = c + 1, c ≤ T; obtain the crowding distance under each objective function.
Step 6: Add the crowding distance under each objective function to obtain the individual's final crowding distance. The greater the value of the crowding distance, the smaller the degree of individual crowding, and the better the diversity of the population.

Selection
In order to obtain better cross-offspring, a tournament selection strategy of two participating individuals is adopted. The size of the mating pool (pool = pop/2) is set; two different participating individuals, preferentially the participating individuals with the highest Pareto levels, are randomly selected to enter the parent population. When two individuals have the same Pareto level, the individual with a larger crowding distance is preferred to enter the parent population. When two contestants have the same Pareto level and the same crowding distance, the first contestant is selected. This process continues until the parent population of the pool size is selected.

Crossover/Mutation
Due to the large number of metal structural parts in the production workshop and the long chromosome length during encoding, the traditional cross-mutation method cannot obtain good results, so this paper proposes the segmented cross-mutation method.
The specific crossover process of segment crossover is as follows: First, parent chromosomes P1 and P2 are segmented, the length of the segment is set to h and the chromosome is divided into x segments. Then, two genetic loci, p1 and p2, in each chromosome of the parent are randomly selected. The chromosome between loci p1 and p2 of P1 is inherited by offspring O1. Finally, the genes in P2 that are different from O1 are inherited by O1, and offspring O2 is obtained in the same way. The schematic diagram of the cross-section algorithm is shown in Figure 7.
ticipating individuals is adopted. The size of the mating pool (pool = pop/2) is set; two dif-ferent participating individuals, preferentially the participating individuals with the highest Pareto levels, are randomly selected to enter the parent population. When two individuals have the same Pareto level, the individual with a larger crowding distance is preferred to enter the parent population. When two contestants have the same Pareto level and the same crowding distance, the first contestant is selected. This process continues until the parent population of the pool size is selected.

Crossover/Mutation
Due to the large number of metal structural parts in the production workshop and the long chromosome length during encoding, the traditional cross-mutation method cannot obtain good results, so this paper proposes the segmented cross-mutation method.
The specific crossover process of segment crossover is as follows: First, parent chromosomes P1 and P2 are segmented, the length of the segment is set to h and the chromosome is divided into x segments. Then, two genetic loci, p1 and p2, in each chromosome of the parent are randomly selected. The chromosome between loci p1 and p2 of P1 is inherited by offspring O1. Finally, the genes in P2 that are different from O1 are inherited by O1, and offspring O2 is obtained in the same way. The schematic diagram of the crosssection algorithm is shown in Figure 7. The mutation method and crossover method in this paper use the same principle. First, the chromosomes are segmented; then, two gene positions, p1 and p2, are randomly selected in each chromosome of the parent; finally, the values of the two gene positions are exchanged. At the beginning of evolution, a larger crossover probability should be selected. This search process is conducive to maintaining the diversity of the population. In the next stage, a careful search is required to prevent the optimal solution from being The mutation method and crossover method in this paper use the same principle. First, the chromosomes are segmented; then, two gene positions, p1 and p2, are randomly selected in each chromosome of the parent; finally, the values of the two gene positions are exchanged. At the beginning of evolution, a larger crossover probability should be selected. This search process is conducive to maintaining the diversity of the population. In the next stage, a careful search is required to prevent the optimal solution from being destroyed and accelerate the convergence rate, and a smaller crossover probability should be selected. The adaptive crossover probability used in this article is shown below.
where p c is the crossover probability; gen is the current number of iterations; maxgen is the maximum number of iterations; and p cmin is the minimum value of the crossover probability.

Elite Retention Strategy
The offspring population whose population size is pop after crossover and mutation and the parent population whose population size is pooled to form a population whose size is pop + pool, and the non-dominated sorting of the combined population is performed. The Pareto level and crowding distance of each individual in the merged population are calculated, and the first pop individuals are kept as the new population obtained after iteration.

Computational Results
In this section, we provide the results we obtained studying a simplified real-life instance of a bus manufacturing company in China to illustrate the proposed scheduling model and the proposed algorithm. Two sets of contrasting experiments with randomly generated instances were also conducted, aiming to further evaluate the performance of the proposed algorithm. The algorithm was coded in MATLAB R2016b, and the experiments were run on an Intel Core (TM) i5-3470 3.20 GHz processor with 6.00 GB of memory in the Windows 7 operating system.

Simplified Illustrative Real Case Study
The scheduling model and the proposed algorithm were applied to customized bus orders from a bus manufacturing company. Due to the large size of the custom order, we simplified the order, aiming to verify the feasibility of the model and algorithm. In this simplified real case, the company received orders for customized buses, as shown in Figure 1a; the simplified numbers of components and parts were lower than the actual number. There were three customized bus models, product 1 (H223), product 2 (S745) and product 3 (L645). Each bus had a customized frame, which was welded together starting from some kits. The first frame consisted of three kits, the second had two kits and the third consisted of four kits. Each component consisted of a variety of metal structural parts, which were machined by machines. We used three machine tools to complete those customized tasks. Figure 8 shows some of the metal structural parts. The complex structure, including the inclined planes, camber grooves and holes, had to be processed by a machine. The details of the corresponding input data are available in Table 1. As the part code and component code were long, they were reported as simplified codes. The data table contains 9 components, 36 parts and 4 materials, as well as the length and processing time of each part.
The offspring population whose population size is pop after crossover and mut and the parent population whose population size is pooled to form a population w size is pop + pool, and the non-dominated sorting of the combined population is perfor The Pareto level and crowding distance of each individual in the merged populatio calculated, and the first pop individuals are kept as the new population obtained iteration.

Computational Results
In this section, we provide the results we obtained studying a simplified real-li stance of a bus manufacturing company in China to illustrate the proposed schedu model and the proposed algorithm. Two sets of contrasting experiments with rand generated instances were also conducted, aiming to further evaluate the performan the proposed algorithm. The algorithm was coded in MATLAB R2016b, and the ex ments were run on an Intel Core (TM) i5-3470 3.20 GHz processor with 6.00 GB of mem in the Windows 7 operating system.

Simplified Illustrative Real Case Study
The scheduling model and the proposed algorithm were applied to customized orders from a bus manufacturing company. Due to the large size of the custom orde simplified the order, aiming to verify the feasibility of the model and algorithm. In simplified real case, the company received orders for customized buses, as shown in ure 1a; the simplified numbers of components and parts were lower than the actual n ber. There were three customized bus models, product 1 (H223), product 2 (S745) product 3 (L645). Each bus had a customized frame, which was welded together sta from some kits. The first frame consisted of three kits, the second had two kits and third consisted of four kits. Each component consisted of a variety of metal struc parts, which were machined by machines. We used three machine tools to complete t customized tasks. Figure 8 shows some of the metal structural parts. The complex structure, inclu the inclined planes, camber grooves and holes, had to be processed by a machine details of the corresponding input data are available in Table 1. As the part code and ponent code were long, they were reported as simplified codes. The data table conta components, 36 parts and 4 materials, as well as the length and processing time of part.  Due to the priority of the optimization objectives in the first stage being higher than that in the second stage, firstly, the scheduling scheme of the first stage was optimized. Through the analysis of the mathematical model, we found that the computational scale and complexity of the first stage were much smaller than those of the second stage. Therefore, in the first stage, we used the exact method to solve the problem, and the chosen solver was Gurobi.
First, according to the order information data, the parts belonging to the same material were clustered, and the total processing time of each kind of material was obtained. The first material (L1 = 2.7 m) included seven parts, that is, HJ1-1, HJ1-3, HJ2-2, HJ3-1, LJ1-6, LJ3-3 and LJ4-3. By adding up the processing time of these parts, the total processing time of this material was calculated to be 87.7 s. Similarly, the total processing time of the second material (L2 = 3.8 m) was 175.9 s, and this material consisted of 13 parts. Moreover, the total processing times of the third (L3 = 4.3 m) and fourth (L2 = 2.3 m) materials were 106.4 s and 95.9 s, respectively. Next, we distributed the four materials to three parallel machines. According to the principle of processing resource balance, all the kinds of materials were evenly allocated to the machines. Formulas (1)-(5) were solved with Gurobi, and the optimal result of material assignment was obtained. The scheduling Gantt chart of material distribution is shown in Figure 9. The optimal makespan was 183.6 s. The time reported in the scheduling scheme included only the part processing time, excluding other non-processing times, such as material handling, as stated in the model assumption.
Next, we distributed the four materials to three parallel machines. According to the principle of processing resource balance, all the kinds of materials were evenly allocated to the machines. Formulas (1)-(5) were solved with Gurobi, and the optimal result of material assignment was obtained. The scheduling Gantt chart of material distribution is shown in Figure 9. The optimal makespan was 183.6 s. The time reported in the scheduling scheme included only the part processing time, excluding other non-processing times, such as material handling, as stated in the model assumption. Finally, we completed the part nesting with the objective of maximizing the material utilization. Using Gurobi to solve Formulas (6)-(9), we generated the result of the optimal part cutting layout, and the corresponding Gantt chart is shown in Figure 10. Two-piece profiles were used for each material, and the parts to be machined were displayed on each profile. The optimization results tried to satisfy the maximum utilization ratio of each profile, unless the parts were not enough to fill the whole profile. We calculated the utilization rate of each profile according to the Gantt chart. The average utilization rate of materials was 80.2%.
We could infer that the optimization results of the first stage could be further improved with the increase in the number of parts, that is, the larger the scale of the scheduling is, the more obvious it is. On the basis of the scheduling scheme in the first stage, we optimized the scheduling in the second stage and used the improved multi-objective algorithm to solve Formulas (10)- (28). We set the appropriate algorithm parameters of improved NSGA-II. Through optimization, we obtained the iterative convergence curve, as shown in Figure 11. The iterative curve showed that the completion time of kits, the number of kits obtained in advance and the number of delayed kits were significantly Machines Makespan Finally, we completed the part nesting with the objective of maximizing the material utilization. Using Gurobi to solve Formulas (6)-(9), we generated the result of the optimal part cutting layout, and the corresponding Gantt chart is shown in Figure 10. Two-piece profiles were used for each material, and the parts to be machined were displayed on each profile. The optimization results tried to satisfy the maximum utilization ratio of each profile, unless the parts were not enough to fill the whole profile. We calculated the utilization rate of each profile according to the Gantt chart. The average utilization rate of materials was 80.2%. reduced. Finally, we converged to the optimal solution or near-optimal solution. The minimum completion time of kits was 131.2 s. None of the kits were completed in advance, and there were two kits that were delayed. We could infer that the optimization results of the first stage could be further improved with the increase in the number of parts, that is, the larger the scale of the scheduling is, the more obvious it is. On the basis of the scheduling scheme in the first stage, we optimized the scheduling in the second stage and used the improved multi-objective algorithm to solve Formulas (10)- (28). We set the appropriate algorithm parameters of improved NSGA-II. Through optimization, we obtained the iterative convergence curve, as shown in Figure 11.
The iterative curve showed that the completion time of kits, the number of kits obtained in advance and the number of delayed kits were significantly reduced. Finally, we converged to the optimal solution or near-optimal solution. The minimum completion time of kits was 131.2 s. None of the kits were completed in advance, and there were two kits that were delayed.  The part scheduling Gantt chart is shown in Figure 12, and the minimized makespan of the machine was 106.4, corresponding to machine M2. Taking machine M2 as an example, parts LJ3-2, LJ1-2, LJ1-4, SJ2-4, LJ2-1, HJ1-4, SJ1-2 and SJ2-2 were continuously processed. Combined with the results of material distribution and part layout, that is, Figures  9 and 10, we could infer that the first part (LJ3-2) belonged to the first profile (g1) of the third material and that the latter six parts belonged to the second profile (g2). The part scheduling Gantt chart is shown in Figure 12, and the minimized makespan of the machine was 106.4, corresponding to machine M2. Taking machine M2 as an example, parts LJ3-2, LJ1-2, LJ1-4, SJ2-4, LJ2-1, HJ1-4, SJ1-2 and SJ2-2 were continuously processed. Combined with the results of material distribution and part layout, that is, Figures 9 and 10, we could infer that the first part (LJ3-2) belonged to the first profile (g1) of the third material and that the latter six parts belonged to the second profile (g2). The kit scheduling Gantt chart is shown in Figure 13. It can be seen from the Gantt chart that the kit completion time distribution was between 115.3 and 183.6 s. The red boxes indicate the delivery time windows in this chart. It can be seen from the figure that no kits were completed ahead of schedule and that kits C1 and C3 were delayed. The rest The kit scheduling Gantt chart is shown in Figure 13. It can be seen from the Gantt chart that the kit completion time distribution was between 115.3 and 183.6 s. The red boxes indicate the delivery time windows in this chart. It can be seen from the figure that no kits were completed ahead of schedule and that kits C1 and C3 were delayed. The rest of the kits were completed within the due-date windows. The kit scheduling Gantt chart is shown in Figure 13. It can be seen from the Gantt chart that the kit completion time distribution was between 115.3 and 183.6 s. The red boxes indicate the delivery time windows in this chart. It can be seen from the figure that no kits were completed ahead of schedule and that kits C1 and C3 were delayed. The rest of the kits were completed within the due-date windows.

Randomly Generated Instance Studies
In order to further validate the effectiveness of the proposed approach for solving large-scale cases, the application validation of large-scale cases was required. Therefore, a large number of simulation cases in different scales were constructed and studied with reference to characteristics of the real case data. In these tests, we utilized the method of generating random instances designed by Olivier Holthaus [37]. There were two sets of randomly generated examples to evaluate the performance of the proposed algorithm. The input variables, including the number of parts, the processing time, the lengths of the parts and types of materials, were randomly generated with a uniform distribution. Each

Randomly Generated Instance Studies
In order to further validate the effectiveness of the proposed approach for solving large-scale cases, the application validation of large-scale cases was required. Therefore, a large number of simulation cases in different scales were constructed and studied with reference to characteristics of the real case data. In these tests, we utilized the method of generating random instances designed by Olivier Holthaus [37]. There were two sets of randomly generated examples to evaluate the performance of the proposed algorithm. The input variables, including the number of parts, the processing time, the lengths of the parts and types of materials, were randomly generated with a uniform distribution. Each instance was studied, and the detailed data are shown in Tables 2 and 3. For the first example, the set of small-size instances, the range of kits was (20,80). For the second set of large-size instances, the numbers of kits were set as 100, 200, 300, 400 and 500. The number of parts belonging to the same component was randomly generated from a uniform distribution of (2, 8) pieces. The processing time of each part was randomly generated in the range of (5, 20) s. The length of each part was randomly generated in the range of (0.2, 1) m. The number of types of materials was randomly generated in the ranges of (10, 15) (small sized) and (15,30) (large sized). The length of each type was randomly generated in the range of (2, 5) m. The material type of each part was randomly selected from the type set. The number of machines in each example is as shown in Tables 2 and 3.
For each group of kit and machine, 10 random instances were generated; then, the average values of 10 tests were calculated and recorded. When the numbers of kits were set as 300, 400 and 500, the algorithm had to run for a long time to find a feasible solution. Thus, considering the feasibility of the experiment, only five random instances were generated in such a group of tests. All the average values are reported in Tables 2 and 3.  Table 3. Computational results for large-size instances.

The Number of Components The Number of Machines
Objectives of Kitting effect of scale on production costs (that is, machine utilization and material utilization) in the first stage. The machine utilization and the material utilization presented in the experimental results were the optimization results obtained under the conditions of the used example data. The utilization ratio of materials gradually increased with the increase in the calculation scale. The utilization ratio of materials in all examples was more than 70%. When the number of kits was greater than 300, that is, the number of parts was greater than about 1500, the utilization ratio of materials was about 95%. This was due to the increase in the number of parts, making it easier for each profile to be filled and causing less off-cut. The material utilization was related to order size and part size. If the order size increases and the part size is appropriate, the proposed solution may result in the higher material utilization. From the point of view of machine task balance, the machine utilization was more than 60%. Meanwhile, we found that under the condition of the same number of components, for small-size instances (that is, the number of kits was smaller than 100 or the number of parts was smaller than 500), the machine utilization decreased with the increase in the number of machines. However, for large-size instances, the utilization increased at first and then decreased with the increase in the number of machines. Thus, the machine utilization was related to the order size and the number of machines. As shown in Figure 14, when there are 40 components and 5 machines, the machine utilization has reached 99.04%. However, the decision maker needs to consider the machine utilization and the material utilization together and choose a compromise set of optimized results. Through the analysis of the above trend, production managers can reasonably choose the number of machines according to the number of orders, so as to obtain appropriate utilization ratios of machines and materials.  Secondly, on the basis of the optimization results of each first example, we obtained the results of the second stage. The statistical results of small-size instances and large-size instances are listed in Tables 2 and 3, respectively. The data showed that the kit completion time and the numbers of earliness and tardiness of kits in each example were reduced. However, we also found that the optimization results were not obvious, especially the number of earliness and tardiness of kits. As the size in the example increased, they accounted for a maximum of nearly 50% of all processed kits. Because optimization was completed on the basis of the first-stage results, each material was scheduled in a single Secondly, on the basis of the optimization results of each first example, we obtained the results of the second stage. The statistical results of small-size instances and large-size instances are listed in Tables 2 and 3, respectively. The data showed that the kit completion time and the numbers of earliness and tardiness of kits in each example were reduced. However, we also found that the optimization results were not obvious, especially the number of earliness and tardiness of kits. As the size in the example increased, they accounted for a maximum of nearly 50% of all processed kits. Because optimization was completed on the basis of the first-stage results, each material was scheduled in a single batch, not in batches. As a result, the parts belonging to the same kits were distributed on different profiles or materials, resulting in not being processed closely together, resulting in a longer kit completion time and an unsatisfactory number of components.
Finally, we analyzed the solution time of each example in the two stages. As shown in Figure 15, it was found that with the increase in the size of the instance, the solution time of the first stage was correspondingly prolonged, but compared with that of the second stage, the solution time was significantly reduced. Therefore, from the point of view of solution efficiency, the feasibility of our scheduling model was verified, that is, the exact method was used in the first stage, and the approximate algorithm was used in the second stage. Although the accuracy of the solution was reduced in the second stage, the calculation time was improved, and satisfactory results were obtained at the same time.

Conclusions and Prospects
Aiming at the problem of the machining and delivery of metal structural parts with eligibility constraints, this paper proposes a complex multi-objective parallel machine scheduling approach that considers part nesting, material assignment and kit delivery. The proposed scheduling problem is divided into two stages: In the first stage, parts are laid out, and the materials are allocated to the machines. The optimization objectives are to maximize the utilization of materials and minimize the total makespan. The second stage sorts the profiles and parts and assigns them to each machine with the objectives of minimizing the kit delivery metrics (kitting time and the numbers of early and tardy kits). A mixed-integer programming model with eligibility constraints is established, and a hierarchical optimization method is proposed. The scheduling results of the first stage are used as the constraints of the second stage. Gurobi is used to obtain the first-stage optimization results; the multi-objective genetic algorithm is improved for the second-stage scheduling optimization. The algorithm includes encoding, decoding, fast non-dominated sorting and genetic operations, and a design considering eligibility constraints is employed. To verify the model and algorithm, we used the actual production data of a metal structure workshop of bus frames. In order to further verify the scalability of the proposed scheduling method, a larger-scale scheduling optimization was conducted through the derived simulation examples. The data showed that with the increase in the number of components, the utilization rate of materials gradually increased. When the number of components exceeded 300-that is, when the number of parts exceeded about 1500-the material utilization rate reached 95%. By selecting the appropriate number of machines,

Conclusions and Prospects
Aiming at the problem of the machining and delivery of metal structural parts with eligibility constraints, this paper proposes a complex multi-objective parallel machine scheduling approach that considers part nesting, material assignment and kit delivery. The proposed scheduling problem is divided into two stages: In the first stage, parts are laid out, and the materials are allocated to the machines. The optimization objectives are to maximize the utilization of materials and minimize the total makespan. The second stage sorts the profiles and parts and assigns them to each machine with the objectives of minimizing the kit delivery metrics (kitting time and the numbers of early and tardy kits). A mixed-integer programming model with eligibility constraints is established, and a hierarchical optimization method is proposed. The scheduling results of the first stage are used as the constraints of the second stage. Gurobi is used to obtain the first-stage optimization results; the multi-objective genetic algorithm is improved for the secondstage scheduling optimization. The algorithm includes encoding, decoding, fast nondominated sorting and genetic operations, and a design considering eligibility constraints is employed. To verify the model and algorithm, we used the actual production data of a metal structure workshop of bus frames. In order to further verify the scalability of the proposed scheduling method, a larger-scale scheduling optimization was conducted through the derived simulation examples. The data showed that with the increase in the number of components, the utilization rate of materials gradually increased. When the number of components exceeded 300-that is, when the number of parts exceeded about 1500-the material utilization rate reached 95%. By selecting the appropriate number of machines, the machine utilization rate reached 90%. The results verified the effectiveness of the scheduling model and algorithm proposed in this paper.
In the future, we aim to focus on extending the method proposed in this paper to more complex scheduling problems, continuing to add some practical constraints not considered in this paper, such as transit time, setup time, uncorrelated parallel machines, batch setup, integrating geometrical shape constrained for part nesting, etc., to formulate a new mathematical model for the scheduling problem of metal structural parts. Dynamic events of machine failures and order changes should also be considered. In addition, the multi-objective genetic algorithm should continue to be improved for scheduling problems, and the performance of the algorithm should be further improved by incorporating a new genetic operation mechanism.
Author Contributions: Conceptualization, X.Z. and J.G.; Data curation, X.Z.; Investigation, X.Z. and J.X.; Methodology X.Z. and J.G.; Project administration, J.X. and Y.W.; Resources and Software, J.X. and Z.X.; Funding acquisition, Z.X.; Writing-original draft, X.Z.; Writing-review and editing, X.Z. and J.G. All authors have read and agreed to the published version of the manuscript.