1. Introduction
Intelligent manufacturing involves intelligent manufacturing technology and intelligent manufacturing systems[1–3]. Intelligent process planning (PP) is one of the most important intelligent manufacturing technologies and plays a crucial role in intelligent manufacturing systems[4–5]. PP can effectively shorten production cycle and improve product quality, and further, reduce the consumption for resources and energy [6]. Hence, PP has gained a considerable amount of attention in industrial applications. Nevertheless, due to the extensive internal flexibility within PP, it is very difficult to solve for an optimal process plan in practical production scenarios [7].
Solving a traditional PP problem usually involves three steps: process selection, resource allocation, and operation sequencing [5]. The first step is to fix the process selection, afterwards, the length of a process route can be determined because the number of operations varies from different process approaches. Next, the allocation of resources, including machines and tools, should be configured for all the operations. Sequencing the operations is the last and most critical step for obtaining a feasible process route that obeys the precedence constraints [8]. For example, milling and grinding must be arranged before fine milling and fine grinding operations, and tapping for a hole thread must be carried out after the operation of drilling a hole [9].
PP problems have been proved to be nondeterministic polynomial-time (NP)-hard[6,10] meaning that it is difficult to solve by only relying on the traditional gradient descent methods, graph theory methods, or simulation-based methods. Therefore, most researchers attempt to introduce meta-heuristics [11] to study PP problems. The main research approaches for PP problems include genetic algorithm (GA)[7,12], tabu search (TS) [13], particle swarm optimization (PSO) [8], ant colony optimization (ACO) [14], and honey-bees mating optimization (HBMO) [10], which aim to find high-quality solutions within less computational cost, and attach more importance to the efficiency, rather than the optimality.
In this background, optimal solutions for many open PP problems[8,15] have not been found yet. One important reason is that recent studies have paid more attention to improving the performance of intelligent algorithms than to the modification of current mathematical models, especially mixed-integer linear programming (MILP) mathematical models [14]. An effective MILP model can be solved by a mathematical programming solver, such as CPLEX, Gurobi, and so forth, and the optimal solution of a PP problem may be stably obtained under certain conditions. However, the solving effect of an MILP model is highly related to the problem scale. Once the scale becomes larger, the solving effect of the model deteriorates rapidly, and the model may even fail to find a feasible solution within a long computation time [16]. Therefore, current research work for PP does not focus on MILP models. In fact, the research status of MILP models is in a state of ‘‘compromise” because almost all reported PP models are feature-based[8,14,17]. These models convert an original network graph into a tabular form to simplify process representations. This kind of conversion, to some extent, changes the original process information of the network. (Discussions on this situation are presented in Section 3.3.) In addition, although there are other modeling methods that can avoid the abovementioned situation, they require complicated preprocessing, such as the generation of combinations of all possible processing operations, in order to eliminate and substitute OR nodes that represent a type of OR logic inside the network [18]. This preprocessing procedure is complicated, in particular, the presence of too many OR nodes or a complex topological structure will greatly increase the complexity and computational time of a PP problem.
In order to fill in the gap of recently reported research work, a novel MILP model is proposed based on the PP network topological structure. The main contributions of this paper are as follows:
(1) This paper proposes a new OR-node-based MILP model for PP problems with no need for any feature conversion or other preprocessing.
(2) The precedence relationships of the operations are discussed in detail, and three types of precedence relationship matrices are presented to illustrate the precedence relationship constraints.
The rest of this paper is organized as follows. Section 2 introduces the related work, while Section 3 proposes the model and discusses the related analysis mentioned above in detail. Section 4 presents several groups of comparative experiments to verify the advantage of the proposed model. Section 5 provides conclusions and outlines some future work.
2. Related work
The related work on PP can be divided into two categories: algorithms and mathematical models. Xu et al. [6] and Leo Kumar [4] have both provided good reviews of PP problems. Intelligent algorithms such as GA, simulated annealing (SA) algorithms, TS algorithms, and PSO algorithms have shown sufficient advantages and have been widely applied in PP problems [6]. To solve PP problems with complex prismatic parts, Li et al. [19] proposed a good hybrid algorithm comprising a GA and a SA algorithm. The local search capability of the hybrid algorithm was enhanced by a strategy of searching and selecting solutions based on the Hamming distance between alternative routes. Hua et al. [20] proposed a synthesis algorithm based on a GA to search for the global or near-global optimal solution in PP problems. Li et al. [15] proposed an effective genetic programming (GP) algorithm, which performed genetic operations on OR-node branches in the network graph. Taking into account manufacturing resources such as tool selection and feed direction, Shin et al. [21] proposed a symbiotic evolutionary algorithm to optimize three objectives, including balancing the machine workload, minimizing part movements, and minimizing tool changes. Wang et al. [22] put forward a PSO combining two local search strategies to solve the PP problem. Subsequently, Li et al. [8] presented a modified PSO to solve PP problems considering the transmission time between machines. Liu et al. [14] combined ACO with the constraint matrix and the state matrix of the problem, and applied it to solve the PP problem of two prismatic parts.
Although some achievements have been made in applying intelligent algorithms to solve PP problems, the solution quality can be further improved in most of the existing PP instances, since metaheuristic algorithms cannot guarantee the optimality of the solutions. Besides, as a method of describing problems, a mathematical model can help researchers to understand and comprehend problems more deeply and thoroughly [23]. Therefore, research on PP mathematical models is very meaningful. Floudas and Lin [24] analyzed several mixed-integer programming (MIP) models of PP problems on the methods of time representation, and then proposed several effective optimization approaches to improve the computational efficiency of the model. In view of the complementarity of PP and scheduling, Li et al. [25] established a mathematical model to integrate a PP problem and shop scheduling problem. Xia et al. [26] proposed a feature-based mathematical model for a reconfigurable PP problem. Similarly, based on features, Jin and Zhang [16] established MILP models for PP considering the transmission time between machines.
To the best of our knowledge, mathematical models for PP problems are always feature-based[8,14,17]. Although this kind of modeling method can describe most kinds of jobs, it is still inevitable that some direct precedence constraints must be added, which do not exist in the original network [16]. Therefore, the network graph method is more capable of describing different manufacturing flexibilities than a feature-based method. This paper proposes a novel MILP model directly based on the topology of the network graph; it can describe all types of manufacturing flexibilities of PP without adding or omitting any constraints. By solving the proposed model, new optimal solutions of some famous benchmarks in previous publications are successfully obtained.
3. The proposed MILP model for PP
3.1. Problem description
The types of flexibilities in PP include process flexibility, machine selection flexibility, and operation sequencing flexibility. There are several approaches to describe PP problems, such as Petri nets [27], feature tables [28], AND/OR graphs, and networks[29,30]. In Fig. 1, the manufacturing flexibilities are represented in the form of a network graph. This network graph is composed of five types of nodes: The starting node, which is virtual, represents the start of a part’s production; the ending node, which is also virtual, indicates the end of a part’s production; intermediate nodes represent operations; and the OR node, combined with the fifth type of node, the JOIN node, represents the process flexibilities [15]. An intermediate node contains three pieces of information: the operation number in the solid circle, the alternative machine number in { }, and the corresponding processing time in [ ]. For example, the intermediate node 6 indicates that operation 6 can be processed on any machine out of the three alternative machines 3, 7, and 13, and the required processing time is 44, 48, and 49, respectively. The time unit in the paper is omitted just like original data. Arrows connecting pairs of nodes in the network indicate precedence relationship constraints between operations [21]. For example, the arrow between operations 2 and 3 declares that operation 2 must be processed before operation 3. There are no fixed precedence constraints between operations not connected by an arrow. Only one link connected to the OR node will be selected. The operations contained in the link with other operations compose a feasible operation combination for the part [29]. Take Fig. 1 as an example: If operations 2 → 3 are selected between OR node 1 and JOIN node 1, and if operation 7 is selected between OR node 2 and JOIN node 2, then one of the complete feasible process routes is 1→2→3→5→6→7→9→10.
Fig. 1. A flexible process plan network. OR1, OR2, JOIN1, and JOIN2 represent OR node 1, OR node 2, JOIN node 1, and JOIN node 2, respectively.
3.2. Precedence relationship between operations
According to the definitions of binary variables, there are three modeling methods for PP problems [31]:
For the large amount of variables and constraints, the binary variable 
based on time periods is rarely mentioned in PP problem modeling. In the existing literature, the definition 
describes the operation sequence by the position of operation j[16,32,33], which is widely used in the current literature[24,34,35]. In this paper, the proposed MILP model for PP is established based on the third method, qjjʹ , for the first time.
based on time periods is rarely mentioned in PP problem modeling. In the existing literature, the definition describes the operation sequence by the position of operation j[16,32,33], which is widely used in the current literature[24,34,35]. In this paper, the proposed MILP model for PP is established based on the third method, qjjʹ , for the first time.As shown in Fig. 2, a feasible operation sequence of this part is 1→2→3→4→5→6→7→9. Operation 1 is in front of all the other eight operations; thus, according to the definition, 
= 1, where jʹ = 2, ..., 9. Operation 2 is processed before all the other seven operations, so 
= 1, where jʹ = 3, ..., 9. The number of ‘‘1”s in the matrix Q = [qjjʹ ] can be obtained as (n – 1)n/2, where n stands for the total number of operations, which is also the length of the operation sequence. Fig. 3 shows the precedence relationship transformation from a sequence to a Q matrix.
= 1, where jʹ = 2, ..., 9. Operation 2 is processed before all the other seven operations, so = 1, where jʹ = 3, ..., 9. The number of ‘‘1”s in the matrix Q = [qjjʹ ] can be obtained as (n – 1)n/2, where n stands for the total number of operations, which is also the length of the operation sequence. Fig. 3 shows the precedence relationship transformation from a sequence to a Q matrix.
Fig. 2. A process plan network.
Fig. 3. The transformation of the precedence relationship.
For a sequence that has been ordered, every two operations have a sequential relationship. Therefore, the priority relation can be represented by an n ×n matrix. The Q matrix contains all the precedence relationships of an operation sequence. The Q matrix can be regarded as a completely expressing matrix (CEM) for the precedence relationships of an operation sequence. Through the observation and analysis, the characteristic constraints of CEM Q can be concluded as follows:
(1) The diagonal elements of CEM Q are equal to 0:
(2) The sum of the two elements symmetrical about the diagonal is equal to 1:
(3) The sum of the elements in any two different columns is not equal:
Eq. (1) indicates that precedence relationships only exist between different operations. In Eq. (2), qjjʹ = 1 makes qjʹj = 0, whereas qjjʹ = 0 makes qjʹj = 1, because there is only one precedence relationship between the two operations. Eq. (3) is established because the sum of the column elements corresponds to the position that is unique in the operation sequence. The Q matrix contains all the precedence relationships between operations. According to the corresponding Q matrix, it can be quickly and easily determined whether a sequence satisfies the precedence constraints in a network.
The notation sjjʹ is defined to represent the precedence constraints in Fig. 2:
The corresponding constraint matrix S = [sjjʹ ] is shown in Fig. 4, where each ‘‘1” in the S matrix corresponds to an arrow in the network, representing a precedence constraint. If a sequence satisfies all the precedence constraints, its CEM Q is supposed to contain all the precedence values shown in matrix S, which can be formulated as follows:
Matrix S can be generated from the network graph, but the matrix Q is unknown because the operation sequence has not been determined yet. Therefore, Eq. (4) can be regarded as the constraint.
Fig. 4. Precedence constraints in a network and the corresponding matrix S.
The number of ‘‘1”s in CEM Q is (n – 1)n/2. However, the least number of ‘‘1”s required to determine an operation sequence in a matrix is (n – 1). For example, to determine the sequence 1→2→ 3→4→5→6→7→8→9, it only needs to set eight variables as ‘‘1,” which are q12, q23, q34, q45, q56, q67, q78, and q89. From this point, it is helpful to define a notation vjjʹ whose corresponding matrix V = [vjjʹ ] contains the least number of variables equal to 1. The matrix V can be regarded as an exactly expressing matrix (EEM) of the precedence relationships of the sequence, which contains the exact and least number of ‘‘1”s. The EEM V of the above operation sequence is shown in Fig. 5.
Fig. 5. The EEM V of the sequence.
The EEM matrix can also be directly called the precedence matrix, since it just contains the precedence relationship of the two directly adjacent operations. According to the direct precedence relationships, it is simple to obtain an operation sequence by sequentially identifying the elements of the EEM V. The EEM matrix has several characteristic constraints as follows:
(1) The number of variables that equal to 1 in the EEM V is (n – 1):
(2) Each row or column of the EEM V has at most one element equal to 1:
(3) The relationship between matrix Q and matrix V is expressed as follows:
Matrix V is a simplified representation of matrix Q, and both of them can determine a unique operation sequence. However, matrix S cannot determine a unique operation sequence due to its incomplete representation of the precedence relationships of a determined sequence. Therefore, matrix S is named the partly expressing matrix (PEM) in this paper. More than one sequence can satisfy the precedence constraints expressed by an S matrix. CEM Q, EEM V, and PEM S are shown in Fig. 6.
Fig. 6. Three types of precedence relationship matrices: (a) CEM Q, (b) EEM V, and (c) PEM S.
3.3. A discussion on feature-based and network-based process representations
In the current literature, all the mathematical models for PP problems are feature-based. For example, an instance adopted in the literature[8,16] is shown in feature tabular form in Fig. 7(a). In fact, this instance is derived from the 18th example in Ref. [29], as shown in Fig. 7(b). Therefore, Figs. 7(a) and (b) illustrate the same instance expressed in two different representation forms.
Fig. 7. Two representation forms of the same instance: (a) feature tabular form and (b) network graph.
In the feature table[8,16], features F2, F5, and F9 have alternative operations or operation sets, corresponding to the three OR nodes OR1, OR2, and OR3 in the network graph. However, in the tabular form, the operation sets under the same feature are constrained by direct precedence relationships that do not exist in the corresponding network graph. For example, if O4–O5 is chosen for feature F2, then O5 must be processed directly after O4. However, O5 does not have to be processed immediately after O4 according to the original network graph Fig. 7(b). As a result, the solution space might be changed, leading to the failure to obtain the optimal solution. Table 1 shows two optimal solutions respectively obtained by solving a feature-based model and a network-based model with the GAMS/CPLEX solver. The letter M with number subscript in the bracket means the allocated machine, and it is the same in the follow-up tables.
Observed from Table 1, the main distinction between these two sequences is that the sequence obtained by the network-based method does not have the direct precedence constraint of O4–O5. In addition, the production time of 356 obtained by the networkbased method is superior to the value of 357 obtained by the feature-based method. Therefore, the method of building a model directly on the basis of the network is superior to the featurebased method.
Table 1 Two optimal solutions obtained by using the feature-based model and the network-based model.
The process flexibility of the feature-based representation is implemented through alternative operations or operation sets for features. For the network-based representation, the process flexibility is expressed by selecting the links of OR nodes. The links refer to the arrow connected to the OR nodes. Each link of the OR nodes corresponds to a choice for the process flexibility. In Fig. 8, according to the order from the top to bottom and the left to right, the OR nodes and their links are numbered as shown. If link 1 of OR node 1 is selected, then operations 2, 3, and 4 will be selected. Additionally, whether link 1 or 2 of OR node 2 is selected, operations 6 and 7 will not be selected. The reason is that, in addition to OR node 2, operations 5, 6, 7, and 8 are controlled by link 2 of OR node 1. The selection performed by the OR node is valid only if the link at which this OR node lies is selected.
A binary parameter wjrl is introduced to describe the control function of OR nodes on the operations. The definition of the parameter wjrl is stated as follows:
For the example in Fig. 8, the corresponding values of wjrl are shown in Table 2.
Hence, the controlling function of the OR node can be concluded as: Operation j will be selected only under the condition that all the controlling links of operation j are selected. On the contrary, operation j will not be selected as long as one of its controlling links is not selected. A binary variable url is introduced to describe the choice of the links, and another binary variable xj is used to describe the operation choosing state. The definitions of url and xj are given as follows:
The model established in Section 3.4 is based on these two variables, url and xj, and on a parameter, wjrl.
Fig. 8. An example of a controlling discussion.
Table 2 The values of parameter wjrl.
3.4. Mathematical model for PP
The proposed MILP model is precedence-based and OR-nodebased. Most of the process optimization objective functions are time related [8] or cost related[14,18]. In this paper, with the objective of minimizing the production time, the transmission time between machines is taken into account in the model. The sets, subscripts, parameters, and variables of the model are introduced below (Table 3).
Table 3 Definitions of the sets, subscripts, parameters, and variables of the mathematical model for PP.
The total production time as the objective can be formulated as follows:
On the right side, the first part of Eq. (9) refers to the total transmission time, and the second is the total processing time. The constraints of the model are displayed as follows:
(1) OR-node controlling constraints:
The constraint in Eq. (10) indicates the ‘‘unselected” condition: Operation j will not be selected as long as one of all operation j’s controlling links is not selected. The constraint in Eq. (11) is the ‘‘selected” condition. The constraint in Eq. (12) means that only one link of an OR node can be chosen.
(2) Precedence constraints:
Eqs.(13)–(18) refer to the constraints in the precedence relationships. Unlike Eqs. (1)–(4), the operation selecting condition is added to this group of constraints. Eq. (13) corresponds to Eq. (1). Eqs. (14) and (15) correspond to Eq.(2), on the condition that operation j and jʹ are selected. Under the same condition, the constraint in Eq. (16) corresponds to Eq. (3). As for Eq. (17), the meaning is that the precedence relationships qjjʹ and qjʹj should be set as 0 if operation j is unselected. The constraint in Eq. (18) guarantees that the selected operation sequence obeys the precedence constraints.
(3) EEM V and CEM Q constraints:
Since a matrix V contains the precedence relationships between the two adjacent operations, the transmission time can easily be calculated. The constraints in Eqs. (19)–(21) describe the property of matrix V. Because matrix V is derived from matrix Q, Eq. (22) is presented.
(4) Machine selection constraint:
Eq. (23) means that there is only one machine that can be assigned for the selected operation, and there is no need to assign any machine for an unselected operation.
(5) Transmission constraints:
Eqs. (24) and (25) formulate the transmission time of operation j from the current processing machine to the next machine.
4. Experiments and discussions
To verify the proposed model, five groups of comparative experiments are carried out based on famous benchmarks. All the experiments are directly compared with the results of other reported methods. On a personal computer (PC) with 3.7 GHz and 16 GB random-access memory (RAM), the proposed model is coded in the GAMS, and the solver CPLEX is used to solve the PP problems. In this paper, the parameter Gap (%) is also introduced to evaluate the proposed model and the computation results. The Gap value represents the relative tolerance of the obtained solution; its definition is (BF – BP)/BP, where BF is the current best solution of the objective function, and BP is the lower bound. The smaller the Gap value is, the closer the current solution is to the optimal solution. The computation time of GAMS/CPLEX is set as 3600 s. If the optimal solution is not found within the time limit, the computation will be terminated and the best known solution will be output. The transmission time between machines adopted by the cases in experiments 1, 2, 4, and 5 is shown in Table 4[8,18].
4.1. Experiment 1
The three cases in experiment 1 are adopted from Jin and Zhang [16], where a dynamic programming (DP)-like heuristic algorithm is applied to solve the PP problem. The results obtained by the proposed MILP model and the DP-like heuristic are presented in Table 5. The production time of case 1 obtained by the MILP method is 357, which is better than the 360 provided by the DP-like heuristic. Furthermore, all three optimal solutions to the cases are found by the MILP model.
Table 5 Comparative results of experiment 1.
a Indicates that the optimal solution is found.
4.2. Experiment 2
The four cases in experiment 2 are from different publications: Case 1 comes from Zhang and Nee [36], and cases 2–4 come from Li and McMahon [37]. The detailed information of these cases can be found in the corresponding papers. Computational results from a modified PSO algorithm are given in the work of Li et al. [8]. The modified PSO algorithm is one of the state-of-the-art algorithms that are used to solve combinatorial optimization problems. The comparison between PSO and MILP is given in Table 6. It can be observed that the MILP model obtains the better solution in case 1. Furthermore, the optimal solutions for cases 2–4 are found by the MILP model method within short computational time (less than one second).
Table 6 Comparative results of experiment 2.
a Indicates that the optimal solution is found.
4.3. Experiment 3
The two cases in experiment 3 are adopted from Li et al. [15], and the machine transmission time matrix is shown in Table 7. These two cases employed the same part from Li et al. [15]. The only distinction between the two cases is that machine 2 is assumed to be broken down in case 2. The comparative results are listed in Table 8. Both the GP algorithm and the MILP model are able to find the optimal solutions of the two instances.
Table 7 Transmission time between the machines.
Table 8 Comparative results of experiment 3.
a Indicates that the optimal solution is found.
4.4. Experiment 4
The 17 cases in experiment 4 are adopted from the well-known Kim dataset [29], which consists of 18 parts, and the comparative calculation results are from Ref. [8]. Because there are some problems regarding the data of part 4 from Ref. [8], it is not selected here. The results obtained by the modified PSO algorithm, the simple GA, and the simple SA are presented in Table 9. Within a reasonable computational time, the proposed MILP model can find 13 optimal solutions out of 17 cases. For the cases that the optimal solutions are not found, the MILP model still obtains better solutions compared to the other algorithms, such as cases 3, 6, 12, and 15.
Table 9 Comparative results of experiment 4.
a Indicates that the optimal solution is found.
4.5. Experiment 5
The 11 cases in experiment 5 are adopted from another famous Shin benchmark [21], and the comparative calculation results are from Li et al. [8]. Because there are some problems regarding the data of some parts in Ref. [8], 11 parts (omitting parts 9, 10, 12, 13, 15, 17, and 18) are selected in this experiment group. The solutions obtained by the modified PSO algorithm, the simple GA, and the simple SA are presented in Table 10. Within a reasonable computational time, the proposed MILP model can find nine optimal solutions out of 11 instances. For the cases that the optimal solutions are not found, the MILP model still obtains better solutions, such as case 3.
Table 10 Comparative results of experiment 5.
a Indicates that the optimal solution is found.
4.6. Discussion
The proposed MILP model obtained 28 optimal solutions out of 37 instances within acceptable calculation time. The solutions found by the MILP model are better than those obtained by the high-performance heuristic [16] and meta-heuristic algorithms[8,15]. Experiments 4 and 5 were carried out on two widely-used benchmarks[21,29], and the better results suggest the superiority of the proposed model.
As shown in Table 11, the proposed model contains four types of subscripts, that is, operation, machine, OR node, and links, while the model reported in Ref. [16] contains six types of subscripts, that is, feature, operation set, operation, machine, position, and place. Fewer subscripts in the proposed model make the computation more effective than that of the model developed by Ref. [16]. Furthermore, the OR-node-based modeling method makes the proposed model more universal for solving different types of PP problems. This is why most optimal solutions can be obtained for both of the Kim [29] and Shin [21] benchmarks.
Table 11 The subscripts of the models.
5. Conclusions and future work
Considering the topology of network graph, this paper proposed a new MILP mathematical model based on OR nodes. Firstly, for precedence relationships between operations, three precedence matrices were introduced. Secondly, for better generality, the notations wjrl, url, and xj were introduced to describe the controlling function of the OR nodes. Finally, the proposed MILP model were coded in the mathematical programming solver CPLEX and tested on public benchmarks. The extensive comparative results verified the correctness and superiority of the proposed model.
In this work, an OR-node-based modeling method was proposed for the first time, demonstrating a new perspective for PP problems and their extension research. The analyses of the three precedence matrices in the paper also revealed the essence of the operation sequencing sub-problem, which is beneficial for further comprehension of the PP problem. However, there are still some limitations in this PP model research. Optimal solutions cannot be found for the minority of instances, and the computational efficiency is not always satisfactory, which implies that the proposed approaches can be further improved. Some simplification and speed-up strategies are urgently required for the further research work.
Acknowledgements
This work is supported in part by the National Natural Science Foundation of China (51825502 and 51775216), and in part by the Program for Huazhong University of Science and Technology (HUST) Academic Frontier Youth Team (2017QYTD04).
Compliance with ethics guidelines
Qihao Liu, Xinyu Li, and Liang Gao declare that they have no conflict of interest or financial conflicts to disclose.