A Column Generation Algorithm for the Resource-Constrained Order Acceptance and Scheduling onUnrelated ParallelMachines

In this paper, we investigate the resource-constrained order acceptance and scheduling on unrelated parallel machines that arise in make-to-order systems. *e objective of this problem is to simultaneously select a subset of orders to be processed and schedule the accepted orders on unrelated machines in such a way that the resources are not overutilized at any time. We first propose two formulations for the problem: mixed integer linear programming formulation and set partitioning. In view of the complexity of the problem, we then develop a column generation approach based on the set partitioning formulation. In the proposed column generation approach, a differential evolution algorithm is designed to solve subproblems efficiently. Extensive numerical experiments on different-sized instances are conducted, and the results demonstrate that the proposed column generation algorithm reports optimal or near-optimal solutions that are evidently better than the solutions obtained by solving the mixed integer linear programming formulation.


Introduction
In many industries, the increasing demand for product variety and customization has intensified competitions among manufacturers to such an extent that attracting and retaining customers is more important than before. e fierce competition has stimulated manufacturers to switch from make-to-stock production to make-to-order production. In the make-to-order production system, manufacturing of the product begins after receiving a customer order. Undoubtedly, accepting all potential orders may cause overloads, order delays, and customer dissatisfaction, due to the limited production capacity and tight delivery requirements. erefore, the first and foremost decision for a make-to-order system is to decide simultaneously which orders to accept and how to schedule the accepted orders.
Due to its application in diverse industries, the order acceptance and scheduling problem has drawn considerate attention from academics as well as practitioners. It has been considered in single machine, parallel machine, and flow shop environments. Various formulations and solution methods for this problem have been proposed. Although there is abundant literature on the order acceptance and scheduling problem, most of these studies treat machine as the only required resource. However, in real-world make-toorder production systems, other resources such as machine operators, tools, fixtures, or industrial robots are also needed and restricted in quantity. An example can be found in the plastic injection molding industry. Plastic injection molding is a production method used to produce plastic products by injecting heated plastic into a preformed mold. e molten plastic is injected at a predetermined temperature and pressure and left to cool. When it is solidified and adapted the shape of the mold, it is ejected as a formed product and the molding cycle is repeated. Figure 1 illustrates the process of plastic injection molding. With no doubt, the injection machine can only start its processing when the mold suitable for the products to be produced is mounted on the injection machine. In general, the molds are quite expensive, and thus the number of molds available is limited. It is obvious that when a manufacturing system like this has limited additional resources besides machines, and thus is unable to meet all customers' demand, it has to decide which orders are to be accepted for processing as well as their associated schedule. erefore, it is of great practical significance to investigate the resource-constrained order acceptance and scheduling problem.
In this paper, we try to fill the gap and address the resource-constrained order acceptance and scheduling problem on unrelated parallel machines. e objective of this problem is to simultaneously select a subset of orders to be processed and schedule the accepted orders on unrelated machines so as to ensure that the resources are not overutilized at any time.
e main contributions of this research are threefold.
(1) We consider a new variant of the order acceptance and scheduling problem on unrelated machines, which takes the additional limited resources into consideration. We present an integer programming model and a set partitioning formulation for this problem. (2) We develop a column generation algorithm based on the set partitioning formulation. In order to solve the pricing subproblems efficiently, we devise a differential evolution algorithm to rapidly price out promising columns. (3) We conduct comprehensive numerical experiments on randomly generated instances to illustrate the efficacy of the proposed column generation algorithm.
e remainder of the paper is organized as follows. Section 2 reviews relevant literature. Problem definition and the integer linear programming model are provided in Section 3, followed by a description of the proposed column generation algorithm and its implementation in Section 4. In Section 5, we discuss our numerical results and validate the effectiveness and the efficacy of the proposed solution approach. Section 6 concludes the contribution and recommends directions for future research.

Literature Review
e scheduling problem addressed in this work is related to two fields of research: (i) order acceptance and scheduling and (ii) resource-constrained unrelated parallel machine scheduling. Existing studies on both topics is reviewed in this section.

2.1.
e Order Acceptance and Scheduling on Unrelated Parallel Machines. Due to its wide applications in the real world, the order acceptance and scheduling problem and its variants have been extensively studied. Interested readers can read Slotnick's study [1] for a comprehensive literature review of the order acceptance and scheduling problem. e following survey only includes a representative subset of papers addressing the order acceptance and scheduling problem on unrelated parallel machines.
Fanjul-Peyro and Ruiz [2] investigated two generalizations of the unrelated parallel machine scheduling problem, one of which is concerned with the integrated order selection and scheduling decision making problem. A mixed integer programming mathematical formulation is presented, along with a constructive method specifically tailored for such a generalization. Comprehensive numerical and statistical experiments show that the proposed method significantly improves the results yielded by the CPLEX solver. Hsu and Chang [3] dedicatedly study the unrelated parallel-machine scheduling problem Step 1: Install the mold suitable for the desired products and close it Step 2: Injection of the heated plastic Step 4: Ejection Step 3: Cooling

Mold
Injection machine Injection cycle with deteriorating jobs and rejection. e objective of this problem is to reject a subset of jobs and then schedule the accepted jobs in such a way to minimize the total cost. Computational experiments show that the problem can be solved in polynomial time when the quantity of machines is fixed. Jiang and Tan [4] study the unrelated parallel machine scheduling model with job acceptance option. Each job order is either accepted and processed without interruption, or is declined at a penalty cost. In an attempt to minimize the makespan of all accepted jobs, the penalty cost of rejected jobs, and the production cost of accepted jobs, the authors first develop a polynomial time algorithm to derive a lower bound and then propose a heuristic by introducing various relaxation techniques. e worst-case ratio bound of the proposed approach is proven to be 2.
Emami et al. [5] consider the simultaneous order acceptance and scheduling problem in an unrelated parallel machines environment with the assumption of uncertain revenue and job-processing times. e objective of this problem is to find a subset of orders and schedule them in such a way that the total profit is maximized. A Lagrangian relaxation algorithm enhanced by a cutting plane technique is implemented and capable of solving instances with up to 40 orders and 6 machines. e same problem is studied in Emami et al. [6] and solved by Benders decomposition approach. Wang and Wang [7] mainly focus on the biobjective optimization of order acceptance and scheduling on unrelated parallel machines with machine eligibility constraints, where a subset of jobs has to be selected to simultaneously maximize the profit level and to minimize the total completion time. e problem is formulated as a biobjective mixed integer linear programming. In view of its NP-hardness, a multiobjective parthenogenetic algorithm with parthenogenetic operators and Pareto-ranking and selection methods are implemented to obtain approximate Pareto front. Mor and Mosheiov [8] study the order acceptance and scheduling problem on unrelated parallel machines, which takes the position-dependent processing times into consideration. is problem aims at minimizing total absolute deviation of job completion times. It is shown that the problem is solved in polynomial time in the number of jobs, and a polynomial time solution approach based on solving a sequence of linear assignment problems is designed to solve this problem. Wang and Ye [9] investigate the order acceptance and scheduling problem on unrelated parallel machines. e purpose is to determine the orders to be accepted for processing and the processing sequence for the accepted orders to get the optimal profit. Two mixed integer programming formulations are presented, which are further enhanced by various valid inequalities. For the purpose of solving complicated instances, the authors propose a tailored branch-and-bound algorithm. e proposed branch-andbound algorithm which follows the paradigm of "divide and conquer" branches on the number of accepted orders and recursively solves the reduced unrelated parallel machine scheduling problem at each child node. Computational experiments confirm its effectiveness.

e Resource-Constrained Parallel Machine Scheduling.
Despite the plentiful achievements mentioned above, research on order acceptance and scheduling with explicit consideration of additional constrained resources is somewhat limited. Since the considered problem is closely related to the resource-constrained parallel machine scheduling, an insight into its recent advances has substantial implications in mathematical modeling and solution development for the studied problem. erefore, a brief literature review on it is presented as follows.
Ventura and Kim [10] investigate a resource-constrained parallel machine scheduling problem which arises in the burn-in process of semiconductor manufacturing. e problem is formulated as an integer programming model with the objective of minimizing the total absolute deviations of job completion times from due dates. e problem is proved to be solved in polynomial time when there is one single type of additional resource, and each job consumes one unit of additional resource at most. Ventura and Kim [11] extend the work of [10] by considering that the number of kinds of additional resources and the resource requirements per job are unrestricted. is more general problem is modelled as a 0-1 integer program, and the Lagrangian relaxation algorithm is applied to obtain high-quality solutions and tight lower bounds for large-scale problems. Chen [12] studies the resource-constrained unrelated parallel machine scheduling problem which originates from the job scheduling in the plastics-forming industry. Beside plastic moulding machines, scarce and costly dyes are also indispensable equipment and thus treated as secondary resource constraints for plastic forming operations. A heuristic is developed to minimize the makespan. e computational results show that it outperforms the simulation annealing algorithm in terms of efficiency and effectiveness. is work is further improved by Chen and Wu [13], who propose an effective tabu search algorithm. Ruiz-Torres et al. [14] investigate the problem of scheduling jobs on uniform parallel machines, where the speed of the machine is determined by the quantity of the assigned additional resource. e objective aim is to minimize the total number of delayed jobs. Two different versions of the problem are discussed. e first one makes the assumption that all jobs are preallocated to the machines, while the second one considers simultaneous optimization of job scheduling and resource allocation. An integer programming formulation is presented for the first case and five heuristics for the second. Edis et al. [15] comprehensively survey the related literature on resource-constrained parallel machine scheduling. ey classify the studies with respect to the characteristics of machine environments, objective functions, additional resource types, complexity of the problem and solution approaches, and other important issues. Strengths and limitations of the existing literature as well as open areas for future studies are also highlighted in this paper. Afzalirad and Rezaeian [16] deal with a resource-constrained unrelated parallel machine scheduling problem which considers various realistic features, such as release dates, sequence-Mathematical Problems in Engineering dependent setup times, machine eligibility, and precedence restrictions. An integer programming model is proposed to represent the complicated problem mathematically. Two different heuristic algorithms including genetic algorithm and artificial immune algorithm are designed to identify high-quality solutions. Comparative experiments are performed on a set of randomly generated test instances and demonstrate that both algorithms perform rather well for small-scale instances and the artificial immune algorithm is suggested to solve largescale instances. Li et al. [17] delicately address a jobscheduling problem arising in the manufacturing systems exposed to both machine and worker availability constraints. A branch population genetic algorithm is proposed to minimize the makespan and cost. e effectiveness of the proposed algorithm is validated by comparative experiments and an actual case study, with the conclusion that the proposed algorithm is very useful under the realistic conditions. Fanjul-Peyro et al. [18] investigate a resource-constrained parallel machine scheduling problem where the processing of jobs on the machines demands a number of units of a scarce resource. Two mixed integer linear programming (MILP) models are presented to solve small-sized instances. Besides, three meta-heuristic algorithms based on the MILP models are developed to solve large-sized instances. Vallada et al. [19] mainly focus on the solution development for the resource-constrained unrelated parallel machine scheduling problem. ey propose a basic scatter search algorithm, an enriched scatter search algorithm and an enriched iterated greedy algorithm, to solve this problem. In all the proposed methods, nonfeasible solutions with violation of resource limitation are allowed and further fixed to a feasible solution using a repairing mechanism. Different local search procedures based on insertion, swap, and restricted neighborhoods are implemented to enhance the performance of the proposed approaches. Nguyen et al. [20] deal with a resource-constrained job scheduling problem abstracted from mineral transportation. A parallel differential evolution algorithm is implemented in conjunction with iterated greedy search and column generation improvement strategies to solve this problem. ey conduct an exhaustive numerical evaluation of the proposed and existing methods using different-sized instances. e computational results confirm the efficiency and effectiveness of the proposed algorithm and advantage of the adopted improvement strategies. More recently, Yepes-Borrero et al. [21] address the unrelated parallel machine scheduling problem with setup times and additional resources in the setups. A MILP formulation and three meta-heuristics are proposed for this problem.
In summary, there are a number of studies on order acceptance and scheduling as well as resource-constrained parallel machine scheduling, but this field of research still requires further development due to its widespread applications. e key extension of our paper which differentiates ours from existing studies is the consideration of order acceptance and resource-constrained parallel machine scheduling in an integrated manner.

Problem Definition and Mathematical Model
In this section, we mainly elaborate on the problem statement and provide an integer programming formulation.

Problem Statement and Integer Programming Formulation.
is paper considers a setting where a manufacture receives orders from its customers but is confronted by additional limited resources. Specifically, let M � {1, ..., |M|} be the set of unrelated machines and N � {1, ..., |N|} be the set of orders. e maximum amount of resource k {k∈K} is constrained to be Q max k . Each machine can process only one order at any time and is always available within the planning horizon T � {0, 1, ..., |T|}. Each order can be processed on atmost one machine without interruption and requires a constant number of the additional resources during its entire process. us, each order can only be processed when one machine and sufficient additional resources are available during its processing time.
Associated with order i are as follows: (1) release date r i ; (2) processing time p im , corresponding to the time required to process order i on machine m; (3) due date d i ; (4) resource requirement of order i to resource k, q i k ; (5) revenue e i gained from its customer if it is accepted; and (6) penalty w i for each unit time of tardiness. Since tardiness penalties will incur a loss of revenue, the objective of this problem is to select a subset of orders and schedule them in such a way that all constraints are satisfied and the total net revenue is maximized. To formulate the problem, the following notations are defined.
Parameters: X i : a binary variable, equal to 1 if order i is accepted and 0 otherwise Y itm : a binary variable, equal to 1 if order i completes its processing at time t on machine m and 0 otherwise C i : completion time of order i T i : tardiness of order i e problem described above can be mathematically formulated as the following MILP model.
e objective function (1) seeks to maximize the total net revenue of accepted orders, which is computed by the total revenue minus the total tardiness penalty. Constraint (2) ensures that not more than one order can be assigned to any machine at any time period. Constraint (3) specifies that an order is either accepted or rejected. Constraint (4) restricts that the total amount of the additional resource k consumed at any time should not exceed the available amount Q max k . Constraint (5) restricts that an order cannot be processed before its release time. Constraint (6) computes the completion time for each order. Constraint (7) links variable Y itm and variable X i and is introduced to determine if an order is accepted or not. Constraint (8) calculates the tardiness for each order. Constraint (9) specifies the domains of decision variables.

An Illustrative Example.
To verify the proposed MILP formulation and intuitively understand this complicated scheduling problem, we consider a simple example involving 6 jobs, 2 machines, and one additional resource type. e maximum amount of the additional resource is 6 units. e parameters related to the example are provided in Table 1.
e mathematical formulation of this example is solved to optimality using CPLEX 12.10 software. e net revenue of the optimal solution is 48. In Figure 2, the optimal schedule is visualized in the Gantt chart. e horizontal axis of the chart shows the time when each order is under processing, while the vertical axis indicates the name of each machine.
Each block on the chart corresponds to an order, which is marked by its order name and the amount of resource required for processing.
It can be seen that all orders except order 3 are accepted. We can also observe that machine 2 is idle from time 4 to time 5 since order 5 and order 6 cannot be processed at the same time due to limited resources. From this example, we can know that the considered problem is difficult to solve because its solution requires not only selecting a subset of orders but also scheduling them on machines in such a way that the resources are not overused at any time.

Complexity of the Problem
e resource-constrained order acceptance and scheduling on unrelated parallel machines is strongly NPhard.
Proof. When we set the maximum amount of the additional resource to sufficiently large and introduce a dummy machine to which all rejected orders are assigned with zero processing time, zero revenue, and zero penalty, the studied problem can be reduced to the classical unrelated parallel machine scheduling problem which is known to be strongly NP-hard (Liaw et al. [22]). Because the classical unrelated parallel machine scheduling problem is a special case of the considered problem, it follows that the considered problem is also strongly NP-hard.

Column Generation Approach
e MILP formulation presented in Section 3 is capable of solving small instances to optimality in reasonable time, but it becomes computationally expensive with the increase of the problem size, as shown in Section 5. is fact necessitates the development of efficient solution methodologies. Column generation is an effective and efficient methodology which has been well studied and successfully applied to a wide range of problems. In this section, we provide an alternative set partitioning formulation and then develop a tailored column generation approach to solve the considered problem.

An Overview of the Proposed Column Generation
Algorithm. Column generation is an efficient approach for solving linear programming problems with a large number of variables (columns). Its main advantage is that not all  Machine 1  Machine 2   Order 1  1  2  5  10  2  3  2  Order 2  2  3  6  12  3  2  3  Order 3  3  4  6  10  4  3  3  Order 4  3  3  7  9  1  3  1  Order 5  4  4  8  11  2  1  2  Order 6  4  3  9  10  2  6  6 Mathematical Problems in Engineering 5 columns need to be enumerated. Instead, the problem is first formulated as a restricted linear master problem (RLMP). is RLMP considers only a subset of columns. And then, one or more subproblems with respect to current dual solution to the RLMP are solved to identify favorable new columns which can be incorporated into the RLMP to improve the objective value. e RLMP and subproblems are solved iteratively until no promising columns can be further added to the RLMP. is termination condition implies that the RLMP is optimal. e RLMP probably returns a fractional solution. In this paper, in order to obtain an integer solution, the RLMP with all columns which have ever been generated as decision variables is converted to an IP model and then solved by an IP solver.
In an attempt to enhance the performance of the column generation, we also propose a differential evolution (DE) algorithm to solve the subproblems efficiently. Besides, we develop a constructive heuristic to obtain the initial solution for the RLMP. e procedures of the proposed column generation approach are shown in Algorithm 1. Next, we present a detailed description of each step.

Master Problem: A Set Partitioning Formulation.
In this section, we reformulate the considered problem as a set partitioning formulation which serves as the master problem in the column generation framework. e set partitioning formulation which is a special class of linear integer programs has been utilized to formulate a wide range of practical problems in fields such as crew scheduling, vehicle routing, and timetable scheduling. Based on the special structure of the set partitioning formulation, efficient approaches techniques, for example, column generation, have been proposed to solve large-sized problems, making it possible to deal with many real-life applications. erefore, we attempt to propose an alternative mathematical formulation based on the set partitioning formulation.
In the set partitioning formulation, there is a set of columns (assignments) denoted by Ω that composes the set partitioning formulation matrix. Each column corresponds to a feasible assignment of a single order to a machine in accordance with the time-related constraints. is matrix contains three submatrices, denoted as A, B, and C, all containing |Ω| columns. Matrix A � (a iω ) contains a row for each order, and a iω � 1, if and only if column ω represents a feasible assignment of order i ∈ N. In matrix B � (b mtω ), m ∈ M; t ∈ T, a single row corresponds to a particular (machine, time) position. If machine m is in use at time t in assignment ω, then b mtω � 1 and 0 otherwise. In matrix C � (c tω k ), a single row corresponds to a time interval, each element of which represents the number of the additional resource assigned in time interval t to process the order which column ω represents. Each column ω is associated with a profit R ω which is defined as the net revenue of the respective assignment of a specific order to a machine. e profit R ω of column ω can be easily computed based on its corresponding order index and the (machine, time) position. e decision variable is denoted by x ω which equals 1 if column ω is used in the solution, and 0, otherwise. erefore, the set partitioning formulation is mathematically represented as follows: x e objective function (10) is to maximize the sum of revenue of selected assignments. Constraint (11) indicates that any order can be accepted atmost once. Constraint (12) limits that a given (machine, time) position can only process atmost one order. Constraint (13) ensures that the resources are not overutilized at any time. Constraint (14) defines the domains of variables.
In practice, the number of feasible assignments (columns) grows exponentially with the problem size. erefore, it is computationally intractable to completely enumerate all feasible columns, making the explicit resolution of the set partitioning formulation impractical for most real-world applications. As an alternative, we solve the set partitioning formulation problem in a column generation framework. e underlying idea of column generation is not to explicitly list all of the columns of the set partitioning formulation (10)-(14), but rather to dynamically generate promising columns. To be more specific, column generation starts with a "simple" formulation including a few columns ω ∈ Ω and relaxes the integrality restrictions on the variables, i.e., x ω ≤ 0, which is referred to as a restricted linear master problem (RLMP). Subproblems with respect to current dual variables of the RLMP are solved to identify favorable new columns which can be incorporated into the RLMP to improve the objective. e RLMP and subproblems are solved iteratively until no promising columns can be further identified.

4.
3. e Initial Solution to the RLMP. As mentioned before, column generation starts with only a small subset of columns and then generates necessary columns by an iterative procedure between the restricted master problem and subproblems. erefore, we should initialize the RLMP. In this section, we design a constructive heuristic to obtain the initial solution for the RLMP. e heuristic works as follows.
In the first step, the most urgent order is selected (i.e., the one with the earliest release time). en, in the second step, it is tentatively inserted into the end of each machine's current partial schedule. e order is rejected if its every possible insertion incurs a decrease in the net revenue. Otherwise, the machine which leads to the largest positive growth in the net revenue is selected to process the order, and the order is placed at the end of this machine's partial schedule. ese two steps iterate until all orders have been checked. e constructive heuristic can always yield a feasible scheduling plan. Afterwards, the scheduling plan is transformed into a series of assignments of orders to machines, which are mathematically represented by columns and construct our initial RLMP.

Mathematical Formulation of the Pricing Subproblem.
We solve the pricing subproblems to identify promising columns with positive reduced cost, whose inclusion in the restricted linear master problem will improve the objective value. Because machines are assumed to be unrelated, there are |M| pricing subproblems in total, one for each machine. at is, for each machine m, we try to find an assignment ω which has a positive reduced cost calculated by C mω � R ω − i∈N a iω β i − t∈T b mtω σ mt − t∈T c tω c t where β i , σ mt , and c kt denote the optimal dual solutions associated with constraints (11)-(13), respectively, in the current iteration of the restricted master problem. e net revenue R ω can also be represented in terms of the variables defined in the original integer programming model (1)-(9), i.e., where S itm is equal to 1 if order i starts its processing at time t on machine m and 0 otherwise. erefore, the pricing subproblem for machine m can be mathematically formulated as follows: Step 1: implement the proposed constructive heuristic to obtain a feasible solution Step 2: initialize the RLMP with the solution obtained by the constructive heuristic Step3: set the number of consecutive iterations of failing to find a column with positive reduced cost Iter non to be 0, i.e., Iter non � 0 Step 4: while Iter non is less than a specified number Iter nonCG , do Step 4.
In this model, the objective function (15) maximizes the reduced costs. Constraint (16) ensures that a column is associated with atmost one order. Constraint (17) ensures that an order can start its processing at time t on the given machine m only if the order is assigned to machine m. Constraint (18) guarantees that no order is processed before its release time. Constraint (19) computes the completion time of each order. Constraint (20) calculates the tardiness of each order. Constraints (21) and (22) define the regions of decision variables.
Ideally, the solution approach for solving the subproblems should be efficient because it has to be executed many times during the column generation procedure. Undoubtedly, the subproblems of large-sized instances will become computationally challenging for generic state-of-the-art solvers, such as CPLEX. Alternatively, we next introduce a DE algorithm for the pricing subproblems, which can solve the subproblems more efficiently.

DE Algorithm for Pricing Subproblems.
e DE algorithm is a simple yet effective metaheuristic algorithm which was first introduced by Storn and Price [23]. Similar to conventional evolutionary algorithms, DE maintains and evolves a population. Specifically, DE starts with a randomly initialized population and then iteratively proceeds from the current population to a new better-performing one via three evolutionary operators, including mutation, crossover, and selection operators. Mutation operator is responsible for generating new mutant individuals, while crossover operator is responsible for exchanging some subcomponents between parent and mutant individuals to produce new trial individuals. Based on a fitness metric, the selection operator compares the parent individuals with the trial ones and the better ones are admitted to the next generation. e three operators iterate until a termination condition is satisfied. e basic evolution process of the DE algorithm is shown in Algorithm 2.
Solution Representation. In order to apply the DE algorithm to solve the subproblems, we need to develop an appropriate encoding approach to represent a feasible solution of a given subproblem. As stated in Section 4.3, the purpose of solving the subproblems is to identify new columns with positive reduced cost, and each column corresponds to a feasible assignment of a single order to a machine in accordance with the time-related constraints. In view of this, we represent each individual of the population by a vector, whose length is 1 + |N|. e first element represents the candidate order, which takes value within [1, |N| ]. e remaining |N| elements corresponds to the starting times of |N| orders, each of which is sampled from its release time to the latest starting time. To better illustrate the solution representation scheme, a simple example for the case mentioned in Section 3.2 is presented in Figure 3.
In this example, the candidate order is order 3, and its starting time is 4. Suppose that the dual variable is given and we are solving the subproblem for machine 1, we can easily obtain a feasible assignment of order 3 to machine 1 and then compute its corresponding reduced cost (the objective function of this individual).
Initialization. e DE population maintains a population of NP vectors. e initial population is randomly generated in the feasible solution domain. For example, the i-th individual can be initialized as follows: where rand is a random number that is sampled from a uniform distribution within (0,1). LB and UB is the lower bound vector and upper bound vector of all elements, respectively.
Mutant Operator. At each generation, the mutation operator creates a new mutant individual V t i with respect to each individual (target individual X t r1 ) of the current population by the mutation operator.
ere are various mutation strategies. We adopt the DE/rand/1 scheme, which is expressed as follows: where X t r1 , X t r2 , and X t r3 are mutually different individuals randomly selected from the current population. e scaling factor F controls the mutation scale and is generally distributed within [0, 1].
If the components of the obtained mutant individual exceed the corresponding lower or upper bounds, they can be reinitialized randomly and uniformly within the given range. Since F is a fractional number, V t i may also contain fractional elements. erefore, if it happens, each fractional element in V t i is rounded to its nearest integer.
Crossover Operator. After mutation, in order to increase the diversity of the population, each pair of target individual X t i and its corresponding mutant individual V t i is crossed to produce a new trial individual U t i . In the DE algorithm, the most commonly used crossover strategy is the binomial crossover scheme, which is defined as follows: where rand is a function returning a randomly generated number within [0, 1]. CR is the crossover rate which is distributed within [0, 1] and mainly controls how many elements of the trial individual originates from the mutant individual. j rand is an integer number randomly selected from {1, 2, 3, ..., 1 + |N|} to enforce that the trial individual differs from its corresponding parent individual in at least one element.
Selection Operator. Finally, to keep a constant population size in the following generations, the selection operator is executed to decide whether the trial or the parent individual survives to the next generation. In this paper, we employ an one-to-one selection operator. To be specific, each trial individual U t i is evaluated, and then its objective function value (reduced cost) f(U t i ) is compared with the objective function value of the corresponding parent individual f(X t i ) in the current population. If the objective function value of the trial individual is greater than or equal to that of the corresponding parent individual, then the parent individual is replaced by the trial one. Otherwise, the parent individual is preserved for the next generation. e selection operator can be described as follows.
Stopping Criteria of the DE Algorithm. In the proposed DE algorithm, two stopping conditions are adopted. If any of these stopping conditions is met, the DE algorithm terminates. ese two stopping conditions are as follows: (1) Maximum number of iterations: Iter max iterations (2) Maximum number of nonimprovement iterations: Iter max−non iterations After termination of the DE algorithm, all columns with positive reduced cost that have been found in the final DE population are added to the RLMP.

Stopping Condition of the Column Generation Process.
e column generation process which is an iterative procedure between the restricted master problem and subproblems aims to identify columns with positive reduced cost. erefore, we terminate the column generation process when it fails to find columns with positive reduced cost for a specified number of consecutive iterations.

Finding an Integer Solution.
Column generation is used to solve the LP relaxation of the restricted master problem and thus may return a fractional solution. erefore, when column generation reports a fractional solution, we adopt a primal heuristic to obtain a feasible integer solution based on the columns already enumerated in the RLMP. e heuristic simply imposes integrality to the master variables, converts the RLMP model to an IP model, and then calls an IP solver to solve the IP model. Since we are solving the IP model on a subset of columns, this process does not necessarily obtain Step 1: initialize a population of NP individuals Step 2: evaluate the fitness values of the NP individuals Step 3: While the preset stopping condition is not satisfied DO Step 3.1: Mutation-create a mutant population using a mutation operator Step 3.2: Crossover-recombine the parent and the mutant to generate a trial population Step 3.3: Selection-compare the parent and trial individuals and select the better ones for the next generation End Step 4: output the best solution ALGORITHM 2: e procedure of the DE algorithm. Mathematical Problems in Engineering an optimal integer solution to the considered problem. However, it has been shown to be useful in practice, especially for the problem we address in this paper, as we will validate in Section 5.3. is heuristic has also been widely used in literature related to column generation-based heuristics (Hauge et al. [24]; Song et al. [25]; Liang et al. [26]).

Numerical Experiments
In this section, we conduct extensive computational experiments to verify the applicability and effectiveness of column generation based on randomly generated instances. Column generation is implemented in MATLAB 2016a and compared against the original MILP formulation solved by CPLEX 12.10. All runs are performed on a Windows 10 professional 64 bit operating system with an Intel Core i5-8250U CPU with 1.60 GHz and 8 GB RAM laptop. For each run, we enforce a time limit of 3600 s of CPU time.

Data Generation.
ere are no publicly available benchmark instances reflecting the setting of the studied problem.
erefore, we randomly generate instances in varying parameter and problem size via the generation schemes of Cesaret et al. [27] and Vallada et al. [19] with some adaptations to the characteristics of our problem. We Each |N| * |M| combination contains 9 subclasses, and each subclass comprises 5 randomly generated instances. Each subclass is characterized by parameters TF and R. More specifically, the parameter TF is the tardiness factor, which controls average due date values and parameter R determines the relative range of due dates. Parameters TF and R all take values from the set {0.2, 0.6, 1.0}, yielding the aforementioned 9 subclasses. Other parameters are generated as follows. Processing times are randomly generated from the interval [1,100], i.e., p im ∼ U [1,100], i ∈ N, m ∈ M. e due dates of jobs d i are generated from the formula where p � i∈N m∈M p im /|M| 2 , similar to the study by Cesaret et al. [27] and Wang and Ye [9]. e release time r i is generated from the interval [0, TF * p]. e number of resource types is set to be 2. For each type of resources, the resource consumption of each job is randomly generated in the interval [1,9]. And, the total number of each type of resources is set to be 6 * |M|. Penalties and revenues are randomly generated from a discrete uniform distribution in the interval [1,10].

Parameter Settings.
e parameter settings affect the quality of the solutions obtained by the proposed algorithm. In this paper, the Taguchi method is applied to tune the parameters.
e Taguchi method which is essentially a fractional factorial experiment approach uses orthogonal arrays to drastically reduce the number of experiments to be performed. e Taguchi method supposes two types of factors affect a process, i.e., controllable factors and uncontrollable noise factors. Since removal of the noise factors is impractical, the Taguchi method seeks to determine the optimal setting of controllable factors and to minimize the effect of noise, in order to achieve a robust performance. e responses at each setting of controllable parameters are considered as a measure, indicative of not only the mean of the performance but also its variance. e mean and the variance are transformed to a single performance measure known as the signal-to-noise (S/N) ratio. e term "signal" signifies the desirable value (mean response variable) and "noise" signifies the undesirable value (standard deviation), so the S/N ratio represents the amount of variation presents in the response variable. us, the goal is to maximize the S/ N ratio.
is paper uses the Taguchi method as a reliable method of setting the parameters for the proposed algorithm. e DE algorithm has five parameters: (1) population size NP; (2) scaling factor F; (3) crossover rate CR; (4) maximum number of iterations Iter max ; and (5) maximum number of nonimprovement iterations Iter max−non . And, the stopping condition of the column generation process involves one parameter Iter nonCG , i.e., the number of consecutive iterations of failing to find a column with positive reduced cost. Five levels are considered for each parameter, as shown in Table 2. e orthogonal array L 25 (6 5 ) shown in Table 3 is designed as the experimental layout, which displays 25 different combinations of factor levels for the experiments.
For each of three different-sized problems, five instances are randomly selected and then tested on the 25 parameter combinations.

Performance Metrics.
We use the gap between the revenue and the upper bound as the performance metric. e revenue is obtained from a given solution method, including either the proposed column generation approach or the integer programming formulation solved by CPLEX solver. e upper bound is obtained by solving the MILP formulation in Section 3 with a time limit of 3600 s. Furthermore, the number of optimal solutions out of 5 instances of each subclass and the average computation time are recorded as auxiliary metrics.
Specifically, columns "CH," "LP," "MILP," and "DE" report the time spent by the constructive heuristic, LP solver, MILP solver, and the proposed DE algorithm, respectively.
As indicated by the gap values in Tables 4 and 5, both the MILP formulation and column generation solve all instances to optimality. In comparison with the MILP formulation, column generation consumes significantly less running time to find the optimal solutions. at is, the MILP formulation needs computing time increasing from few seconds to hundreds of seconds, whereas column generation is able to find the optimal solution in seconds. For example, the average CPU time reduces from 119.14 s with the MILP formulation to 3.24 s with the proposed algorithm for instances with |N| � 20, |M| � 2, TF � 0.2, and R � 0.2.  To sum up, the numerical experiments on small-sized instances show that the proposed algorithm significantly outperforms the MILP formulation in terms of computational efficiency.

Results of Medium-Sized Instances.
We further examine the performances of the proposed method and the MILP formulation for solving medium-sized instances. Tables 6 and 7 present the experiment results.
As shown in Tables 6 and 7, we can observe that the proposed column generation algorithm is superior to the MILP formulation in terms of efficiency and effectiveness. For instances with |N| � 40, column generation finds 83 optimal solutions out of 90 instances with 0.13% deviation on average. In contrast, the MILP formulation identifies the  optimal solution in 80 out of 90 instances and reports 1.72% deviation on average. As expected, the performance of the MILP formulation degenerates since the solution space exponentially grows with the increase in problem size. When |N| equals 60, the MILP formulation only reports 36 instances out of 90 to optimality within 3600 s, whereas column generation finds 79 optimal solutions. e average CPU time of the MILP formulation increases to 2350.70 s, while column generation runs in 53.59 s on an average. Based on above comparison analysis, we can come to the conclusion that column generation still performs significantly better than the MILP formulation with regards to efficiency and effectiveness and is an effective approach for solving medium-sized instances.   Tables 8 and 9. ese instances become more challenging for the MILP formulation, which can be indicated by computation time. We can observe that column generation performs significantly better than the MILP formulation in terms of number of instances solved to optimality and computation times. Specifically, the average percentage deviation of the MILP formulation is consistently larger than that of column generation. We can also observe that the MILP formulation fails to find a feasible solution to instances with TF � 0.2, R � 1 and |N| � 100. In contrast, the proposed "-"indicates that the CPLEX fails to find feasible solutions for this subclass. " * " indicates that the CPLEX runs of memory when solving some instances of this subclass. column generation algorithm is able to prove the optimality of 177 instances and to provide very low gaps for the remaining, always below 1%. Column generation can generally solve each instance within 800 s, which further demonstrate its efficiency. e encouraging computational results of the proposed algorithm may be attributed to two aspects. Firstly, the column generation approach keeps most columns of the set partitioning formulation implicit and only generates necessary columns that can improve the objective value. erefore, a considerate number of times can be saved. Secondly, the DE algorithm is an efficient heuristic. It is used to solve the subproblems and helps to identify promising columns in a very short time.
In summary, the capability of the proposed algorithm to consistently find optimal or near-optimal solutions confirms its effectiveness and superiority against the MILP formulation. Our experimental results clearly show that column generation is a computationally viable solution approach for practical situations.

Effect of Parameters TF and R on the Performance of the
Methods. In this section, we further investigate the impact of the instance generation parameters on the performance of the methods. As stated in Section 5.1, parameter TF determines due date tightness while parameter R controls the due date range. Additionally, parameter TF also controls the dispersion degree of the release dates of the jobs, in a way that the higher the TF is, the more dispersed the release dates are.
As indicated by the gap and computation time values in Tables 4-9, we can notice that given that all other specifications are the same, instances with large values of TF are generally easier to solve than those with small values. A reasonable explanation for this behaviour could be as follows. A small TF value corresponds to a situation where the release dates are distributed in the very beginning of the time horizon and then are not very restrictive. And, a very large TF value leads to a situation where the release dates are very dispersed. Hence, as TF increases, we receive incoming orders released at dispersed time points with tight due dates. erefore, more orders are expected to be rejected, which in turn implies the reduction of the size of the solution space. Consequently, the instances with large TF values are less challenging. We can observe from the tables that the number of instances solved to optimality by the MILP formulation increases and the average CPU time spent by the MILP formulation tends to decrease as the TF gets large when it comes to medium-and large-sized instances. A similar observation can be found when instances are solved by the proposed column generation algorithm. at is, two subclasses are all exactly solved by the proposed algorithm, but the subclass with TF � 1 evidently consumes less CPU time than the subclass with TF � 0.2.
Given the TF value, the due date range increases with the increase of the parameter R. Results in Tables 4-9 show that given the TF value, some subclasses with large R value are more difficult to solve than those with small R value, but we can also find some exceptions. is coincides with observations reported by Wang and Ye [9]. erefore, there is no clear trend in correlation between the R value and the hardness of the problem.

Conclusion
In this paper, we study a resource-constrained order acceptance and scheduling on unrelated parallel machines. We present two different models for this problem, i.e., the MILP formulation and the set partitioning formulation. With the MILP formulation including all decision variables, we can directly solve problem instances via general-purpose commercial solvers without complicated programming. In order to solve large-scale instances, we develop a column generation algorithm based on the set partitioning formulation. In the proposed column generation algorithm, we propose a constructive heuristic to initialize the restricted linear master problem and design a differential evolution algorithm to solve the subproblems efficiently. Numerical experiments based on a set of randomly generated test instances are carried out to verify the performance of the proposed algorithm. e results demonstrate that the MILP formulation can identify optimal solutions for small-scale instances and report solutions with large optimality gaps at termination for medium-and large-sized instances. Column generation is advantageous over the MILP formulation in terms of effectiveness and efficiency and quite suitable for problems with different sizes. We also analyse the effect of parameters TF and R on the performance of the proposed algorithm and the MILP formulation.
Further research will concentrate on problem extensions and algorithmic enhancements. One possible extension of the current work is to consider stochasticities in processing times and develop a stochastic and robust counterpart for the column generation algorithm. Considering some practical assumptions including machine breakdown, parallel batch processing, and order interruption are complex but extremely valuable variations to the studied problem. Furthermore, we also attempt to improve the performance of the proposed column generation algorithm by solving the subproblems via a dynamic algorithm in conjunction with heuristics.

Data Availability
All data used to support the findings of the study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare no conflicts of interest.