Optimization Models for Scheduling of Jobs

This work is motivated by a particular scheduling problem that is faced by logistics centers that perform aircraft maintenance and modification. Here we concentrate on a single facility (hangar) which is equipped with several work stations (bays). Specifically, a number of jobs have already been scheduled for processing at the facility; the starting times, durations, and work station assignments for these jobs are assumed to be known. We are interested in how best to schedule a number of new jobs that the facility will be processing in the near future. We first develop a mixed integer quadratic programming model (MIQP) for this problem. Since the exact solution of this MIQP formulation is time consuming, we develop a heuristic procedure, based on existing bin packing techniques. This heuristic is further enhanced by application of certain local optimality conditions.


Introduction
The present study was motivated by a problem encountered by an aircraft logistics center that houses a number of facilities (hangars), each of which contains multiple work stations (bays). The center offers its customers a wide variety of aircraft maintenance services (jobs), including structural repairs, passenger/freight conversions, and engine changes. In this paper we concentrate on a single hangar which is equipped with several bays.
It is assumed that the company has contractual obligations for a number of existing jobs, which have to be completed in the hangar. These existing jobs are characterized by their specified starting times, bay assignments, and job durations (spans). Thus, a given bay will be idle during the time periods between existing jobs. We are interested in how best to schedule a number of new jobs into these vacant spaces. Specifically, the company wishes to construct a manufacturing schedule that can be flexible in accommodating future incoming jobs.
First we present a mixed integer quadratic programming (MIQP) approach to model this scheduling problem. Solution of this MIQP gives a schedule that optimizes the overall facility utilization, while providing increased flexibility to accommodate future jobs. Next we develop local optimality conditions for the common types of changes that can occur to a particular schedule (job switches). These conditions enable us to better understand optimal solutions and also to conduct sensitivity analyses.
The exact MIQP formulation of the scheduling problem has mn binary variables, where m denotes the number of vacant spaces and n denotes the number of new jobs. Since there are 2 mn possible choices for the binary variables, exact solution of our MIQP model is possible only for very small problems. Consequently we devel-This work is motivated by a particular scheduling problem that is faced by logistics centers that perform aircraft maintenance and modification. Here we concentrate on a single facility (hangar) which is equipped with several work stations (bays). Specifically, a number of jobs have already been scheduled for processing at the facility; the starting times, durations, and work station assignments for these jobs are assumed to be known. We are interested in how best to schedule a number of new jobs that the facility will be processing in the near future. We first develop a mixed integer quadratic programming model (MIQP) for this problem.
op a heuristic procedure by considering the vacant spaces as bins with varying capacities and applying some modified bin packing techniques. We can also incorporate the local optimality conditions to improve this basic heuristic approach.

Optimization Model
First we consider a very simplified scenario in order to motivate our optimization model. Suppose there is only one bay and, relative to the jobs already scheduled for that bay, there are three vacant spaces (bins) B 1 , B 2 , B 3 of lengths (capacities) 8 d, 12 d, 15 d. Also suppose that three new jobs J 1 , J 2 , J 3 with spans of 6 d, 7 d, 10 d are to be assigned to the bay. Several ways of assigning the three new jobs to the three bins are depicted in Fig.  1. Intuitively, it seems that assignment A 3 is the most desirable of the three, because it provides the largest contiguous space remaining in a bin (9 d) after assigning the three jobs. For example, by adopting assignment A 3 , the center is capable of handling a future job that requires a span of 9 d or handling two additional jobs with spans of 6 d and 3 d, whereas the other two assignments cannot. Here we observe that creating a large amount of space in a bin (or in a few bins) after assigning the jobs is generally preferred over small amounts of space that are scattered over many bins.
To be more precise, let us define the residual capacity of a bin as the capacity that remains after assigning a job (or jobs) to that particular bin. Since the sum of the residuals is fixed (in this example totaling 12) for all assignments, we consider the sum of squared residuals for these assignments. Assignment A 1 has residual capacities 8, 2, 2 with sum of squares 64 + 4 + 4 = 72; similarly, A 2 has residual capacities 2, 5, 5 with sum of squares 54, while A 3 has residual capacities 1, 2, 9 with sum of squares 86. We observe that assignment A 3 has the largest sum of squares, as a result of the contribu-tion from the highest residual capacity term. In general, we argue that an assignment with the largest sum of squares should be preferred because it leaves a higher residual capacity in a bin (or bins) instead of smaller residual capacities that are scattered over several bins. In turn this gives added flexibility to the company in accommodating future incoming jobs. We are thus motivated to use the maximization of the sum of squared residuals as our objective criterion. Since the sum of the residuals is fixed, we are equivalently maximizing the variance of the residuals.

Formulation
To formulate this problem, suppose that n new jobs are to be scheduled using m existing bins, with c k being the capacity of bin B k . Also suppose that job J i has span s i . Define the binary variable x ik = 1 if job J i is assigned to bin B k , and x ik = 0 otherwise. Also let z k denote the residual capacity of bin B k . Then our optimization model for the scheduling problem has the form: The above optimization model has m + n equality constraints involving mn binary variables and m continuous variables. The first m equality constraints define the residual capacity variables z k ; the nonnegativity of z k ensures that the set of jobs assigned to bin B k should not exceed that bin's capacity c k . The remaining n equality constraints require that each job should be  ∑ ∑ ∑ assigned to exactly one bin. Since we are optimizing a quadratic objective function over the (linear) constraint set, this model can be classified as a mixed integer quadratic programming (MIQP) problem. Since the objective function is a sum of squares, this problem is a type of convex maximization problem, a very difficult type of mixed integer optimization problem. Thus far we have not considered the starting time t k associated with each bin. We can easily incorporate such starting times into our MIQP model by simply modifying the objective function: where β is a nonnegative valued parameter. When β = 0, we recover our original model in which maximizing the sum of squared residuals ("good fit") is the sole criterion. However, by increasing the value of β, we emphasize the contribution of the starting time and thus encourage new jobs to be assigned earlier in time. The parameter β can be interpreted as the increased importance that the decision maker would like to give to earlier completion of jobs. It is not hard to verify that z(β) is a piecewise quadratic, increasing convex function of β.
The above model can be easily applied to situations in which there are multiple bays within a hangar. Indeed for the purposes of the optimization model MIQP, it is only necessary to maintain a list of bins with their capacities and starting times, without regard to the particular bay associated with each bin. In the multiple bay context, the parameter β can still be interpreted as the relative importance of a good fit versus an assignment of jobs to bins occurring earlier in time.

Examples
In order to get some insights into the proposed MIQP model we present two numerical examples. We have used the global optimization package LINGO [4] to solve these MIQP problems exactly. Example 1 involves m = 5 bins and two bays: the first three bins B 1 , B 2 , B 3 are associated with bay 1 and the other two bins B 4 , B 5 with bay 2. Capacities and starting times for these bins are given in Table 1. Four new jobs J 1 , J 2 , J 3 , J 4 are to be assigned; these jobs have the span times s i given in Table 2. Using LINGO we obtained the exact solution of the MIQP model for various values of β. It turned out that there were only three sets of optimal solutions as the parameter β was varied. For 0 ≤ β ≤ 0.11, the optimal assignment A 1 allocates jobs J 1 , J 2 , J 3 , J 4 to bins B 2 , B 4 , B 3 , B 1 respectively. A schematic illustration of assignment A 1 is given in Fig. 2. For the range 0.11 ≤ β ≤ 0.4, the optimal assignment A 2 allocates jobs J 1 , J 2 , J 3 , J 4 to bins B 4 , B 2 , B 3 , B 1 respectively. Thus the only modification that occurs to the optimal job schedule, in changing from assignment A 1 to assignment A 2 , is that jobs J 1 , J 2 interchange positions between bins B 2 , B 4 . Over the final range 0.4 ≤ β < ∞, the optimal assignment A 3 allocates jobs J 1 , J 2 , J 3 , J 4 to bins B 4 , B 2 , B 5 , B 1 respectively. The only difference between assignment A 2 and assignment A 3 is the movement of job J 3 from bin B 3 to bin B 5 . This behavior is consistent with our earlier observation that increasing β  Example 2 involves m = 9 bins and three bays, each containing three bins. We are interested in assigning five new job J 1 , ...., J 5 to these bins. The data for this problem are given in Tables 3 and 4. Here we also investigated the changes in optimal solutions as the parameter β is varied. For this problem, there were eight sets of optimal solutions as the parameter β was varied over the ranges:  [1,3], [3, ∞]. As one illustration, Fig. 3 presents the optimal solution A 4 over the first of these ranges. In assignment A 4 , jobs J 1 , J 4 , J 5 are assigned to bins B 3 , B 4 , B 5 whereas jobs J 2 , J 3 are assigned to bin B 9 .

Local Optimality Conditions
In the previous section, we observed several transformations that occurred in changing from one optimal schedule to another optimal schedule. Specifically, in changing from A 2 to A 3 , a job (J 3 ) moves from one bin to another bin, a transformation we call a move. In changing from A 1 to A 2 , two jobs (J 1 , J 2 ) interchange positions between two bins, termed a swap. In addition, another simple transformation involves transferring jobs from several bins to a common bin (grouping) or vice versa (ungrouping). Collectively we call these types of changes, and combinations thereof, job switches. In this section we will first study the individual effects of each type of job switch upon a given assignment. This enables us to develop local optimality conditions for these types of job switches, relative to a given assignment. Using these local optimality conditions, it is possible to identify job switches that can improve a given feasible solution. Moreover, these local optimality conditions can aid us in carrying out a sensitivity analysis with respect to the parameter β.

Move
Consider the two assignments A and B shown in Fig.  4. In assignment A, job J i is currently assigned to bin B u , whereas in assignment B, job J i has been moved to bin B v . Otherwise the two assignments are identical. The starting times of bins B u and B v are t u and t v respectively, and the span time of job J i is s i . Define a u (resp. a v ) as the capacity utilized by the other jobs that are already assigned to B u (resp. B v ), and r u (resp. r v ) as the remaining capacity of bin B u (resp. B v ) before the Let Z A , Z B denote the objective function values in Eq. (1) for assignments A, B. Since the only difference between the two assignments is the movement of job J i from bin B u to bin B v , we can ignore the other contributions in the objective function values Z A and Z B , writing If we define y u = r u + βt u , y v = r v + βt v , then the above expressions simplify to Consequently ∆Z = Z A -Z B = 2(y v -y u )s i is the change in objective function (1). Assignment A will be locally optimal with respect to the move in Fig. 4 if and only if ∆Z ≥ 0; since s i > 0 this is equivalent to y v -y u ≥ 0. However, if y v -y u < 0 then it is beneficial to move job J i from bin B u to bin B v , assuming this move is feasible.
It will be convenient to use an alternative representation to depict the move shown in Fig. 4. Namely, we denote the bins B u and B v by vertices u and v, and we use a directed edge i to denote the movement of job J i from bin B u to bin B v . This is depicted in Fig. 5.
To illustrate how this local optimality condition can be applied, consider Example 1 and assignment A 2 , which allocates jobs J 1 , J 2 , J 3 , J 4 to bins B 4 , B 2 , B 3 , B 1 . The only feasible moves involve transferring a job from one of the bins B 1 ,..., B 4 to bin B 5 . For example, assignment A 2 remains locally optimal for moving job J 3 from bin B 3 to B 5 if and only if 0 ≤ y 5 -y 3 = (16 + 28β) -(10 + 43β) = 6 -15β, so that β ≤ 0.4. The other conditions 0 ≤ y 5 -y 1 , 0 ≤ y 5 -y 2 , 0 ≤ y 5 -y 4 hold automatically since β ≥ 0. The result is that assignment A 2 is locally optimal under moves for all β ≤ 0.4, consistent with the results obtained earlier from LINGO.

Swap
We now consider the swapping of jobs, relative to two assignments A and B. In assignment A, jobs J i and J k are assigned to bins B u and B v respectively, whereas in assignment B, jobs J i and J k are assigned to bins B v and B u . Therefore the only difference between the two assignments is that jobs J i and J k are swapped between bins B u and B v . This is illustrated in Fig. 6. Again define a u (resp. a v ) as the capacity utilized by the other jobs that are assigned to B u (resp. B v ), and r u (resp. r v ) as the remaining capacity of bin B u (resp. B v ) before the assignment of jobs J i and J k . The residual capacities of bins B u and B v in assignment A are defined as z u and z v respectively, whereas the corresponding residual capacities of bins B u and B v in assignment B are defined as z′ u and z′ v . Thus   As before, let Z A , Z B be the objective function values for assignments A, B. Since the only difference between the two assignments is the swapping of jobs J i and J k between bins B u and B v , we can ignore the other contributions in the objective function values Z A and Z B , writing and assignment A is locally optimal with respect to the indicated swap if and only if ∆Z ≥ 0. Otherwise it is beneficial to swap jobs J i and J k , assuming that such a swap is feasible.

Grouping and Ungrouping
Another simple job switch involves the grouping of jobs, in which p jobs are moved from separate bins to a new bin. Specifically, suppose that in assignment A, jobs J 1 , ..., J p are assigned to bins B u 1 , ..., B u p , whereas in assignment B, they are all grouped together in the single bin B v . This transformation is illustrated in Fig. 7. Here we can define the utilized capacity of all bins after the assignment of jobs other than J 1 , ..., J p , as well as the residual capacities of bins B u 1 , ..., B u p , in a manner similar to that done previously.
The associated objective function values for these assignments are given by Simplification then produces An analogous development produces the following expression when p jobs J 1 , ..., J p are moved from a single bin B u to p separate bins B v 1 , ..., B v p : Notice the change in sign in the last summation from negative to positive for the ungrouping ∆Z, compared with the grouping ∆Z.

More Complex Job Switches
This section considers combining the previous transformations (moves, swaps, groupings, and ungroupings) into more complex job switches. To motivate the general case, we first consider an example involving eight jobs and five bins. Namely, the transformations that occur, in changing from the current assignment A to the new assignment B are as follows: jobs J 1 , J 5 swap between bins B a , B c ; jobs J 2 , J 3 ungroup from bin B a to bins B d , B e ; job J 4 moves from B b to B c while job J 6 moves from B d to B b ; and jobs J 7 , J 8 ungroup from B e to B b , B c . See Fig. 8 for a graphical representation of this more complex rearrangement. Following the previous development, we obtain Notice that for each directed edge (i, j) representing the movement of job J k in Fig. 8 there is a term 2(y jy i )s k in ∆Z. There are also terms in ∆Z to represent the ungrouping and grouping of jobs. Jobs J 1 , J 2 , J 3 ungroup at vertex a, giving rise to the product terms s 1 s 2 + s 2 s 3 + s 1 s 3 ; likewise the ungrouping of jobs J 7 , J 8 at vertex e gives the term s 7 s 8 . On the other hand, jobs J 1 , J 4 , J 8 group at vertex c, producing the (negative) term -(s 1 s 4 + s 4 s 8 + s 1 s 8 ), while jobs J 6 , J 7 group at vertex b, giving -s 6 s 7 .
In general suppose that the directed graph G = (V, E) represents the specified changes to a current assignment, where the set V of vertices represents the bins and the set E of directed edges represents the movement of jobs between bins. Denote the span time of the job associated with edge (i, j) by s ij . Also let Γ + (i) denote the set of edges leaving vertex i and let Γ -(i) denote the set of edges entering vertex i. Then the expression for ∆Z becomes

Heuristic Procedures
In Sec. 2.1, we modeled the scheduling problem as an MIQP with mn binary variables and m continuous variables, where m is the number of bins and n is the number of new jobs. Since the number of possible choices 2 mn for the binary variables rapidly becomes large, even for small m and n, the exact mathematical solution of the MIQP model is very time consuming. Therefore we consider heuristic solution approaches in this section, rather than exact procedures.
As noted earlier the vacant time interval between any two existing jobs can be considered a bin. Since the intervals between existing jobs are of different length, the corresponding bins have variable size. The underlying problem is then to assign the new jobs to these bins in such a way that the set of new jobs assigned to a bin fits within that bin's capacity; this is to be done in an "optimal" fashion. Thus our scheduling problem can be viewed as a variable-sized bin packing problem [3,5] with an unusual type of objective function. In this section we will first develop a simple heuristic procedure based on existing bin packing techniques. Next we will apply the local optimality conditions developed in Sec. 3 to improve this bin packing heuristic.

Bin Packing Heuristic
We begin by reviewing two bin packing algorithms that are well known in the literature [2]. They are the First Fit (FF) algorithm and the Best Fit (BF) algorithm. The objective of such standard bin packing algorithms is to minimize the number of bins that are need- ed to pack a given set of items. The FF algorithm assigns the next job into the lowest indexed bin into which it will fit. If the next job does not fit into any existing bin then we open a new bin and place the next job in the new bin. The BF algorithm is similar to the FF algorithm, except that it assigns the next job into that bin which will leave the smallest residual capacity after the assignment.
In our heuristic algorithm we first order the jobs in order of nonincreasing span times. There is an intuitive appeal to ordering the span times in this manner: we assign the jobs with higher span times first and hope that we can accommodate the jobs with smaller span times using the spaces that remain. A similar strategy is adopted in existing bin packing algorithms [1,6]. Before assigning the next job (in order of nonincreasing span) to a bin, we first identify the bins into which the job can fit. Among these candidate bins, we select a bin, say bin B u , having the minimum value of z u + βt u . As shown in the proposition below, this procedure enables us to (locally) improve the current solution as much as possible at the next step. Next we assign the new job to the selected bin and update the residual capacity of that bin. We follow this procedure until all jobs have been assigned. We now summarize more formally the steps of this heuristic algorithm.

bin_packing
1. Order the jobs by nonincreasing span time.
2. For the next job J k in order, with span s k , consider only those bins with residual capacity at least s k . Among such bins, select bin B u to minimize z u + βt u . 3. Assign job J k to bin B u and update z u ← z us k . 4. If there are additional jobs to be processed, go to Step 2.
Here β ≥ 0 is the same parameter introduced in Sec. 2.1. When β = 0 we assign the next job J k into a bin (into which it fits) having the minimum current residual capacity. Since this will also leave the smallest residual capacity after the assignment, our heuristic algorithm is similar to the BF algorithm for β = 0. When β is large, we assign jobs to bins in order of increasing starting time. In other words, jobs are assigned to bins that occur earliest in time. This is analogous to assigning a job to the lowest indexed bin into which it will fit, when bins are ordered by starting times. So when β is large, our heuristic algorithm behaves similar to the FF algorithm. In this way we have blended both BF and FF into our particular bin packing heuristic. We now demon-strate that our heuristic performs the locally "best" assignment for the current job at each step of the algorithm.
Proposition. The bin packing heuristic locally improves the objective function by as much as possible at the next step. Proof. Suppose that job J k with span s k is to be assigned and that bins B 1 , ..., B p are the bins into which J k can fit. Let z 1 , ..., z p be the current residual capacities of these bins. Select bin B u such that z u + βt u = min{z 1 + βt 1 , ..., z p + βt p }. ( Let Z j be the objective function value (1) obtained by assigning job J k to bin B j at the next step. We claim that Z u ≥ Z j for all 1 ≤ j ≤ p. Let α denote the contribution to the objective function from all bins other than B u and B j in the current assignment. Thus giving Z u -Z j = -2z u s k -2βt u s k + 2z j s k + 2βt j s k where the final inequality follows from Eq. (2).
We next select (in order) job J 2 , which can fit into bins B 2 , B 5 , with residual capacities 12, 16. Since min{12 + 0.3(22), 16 + 0.3(28)} = 18.6 is achieved for bin B 2 , we assign job J 2 to bin B 2 and update the residual capacity of bin B 2 to z 2 = 1. Continuing in this fashion, job J 4 is assigned to bin B 1 and job J 3 is assigned to bin B 3 . Thus the heuristic assigns jobs J 1 , J 2 , J 3 , J 4 to bins B 4 , B 2 , B 3 , B 1 respectively. This assignment is identical to the optimal assignment A 2 found for the range 0.11 ≤ β ≤ 0.4. In fact for Example 1 the heuris-Volume 111, Number 2, March-April 2006 Journal of Research of the National Institute of Standards and Technology tic produces optimal assignments over all three ranges for the parameter β.
For Example 2 the MIQP model gives eight optimal assignments corresponding to eight ranges of the parameter β. The assignments obtained by the heuristic procedure are identical to the optimal assignments obtained by the MIQP model in seven of these eight ranges. The only difference occurs over the range 0 ≤ β ≤ 0.004 for which the optimal assignment A 4 is shown in Fig. 3; it assigns jobs J 1 , ..., J 5 to bins B 3 , B 9 , B 9 , B 4 , B 5 . By contrast, the heuristic procedure obtains a different assignment A′ 4 , in which J 1 , ..., J 5 are assigned to B 3 , B 7 , B 1 , B 4 , B 5 . If however we perform in A′ 4 a grouping of jobs J 2 , J 3 from bins B 7 , B 1 to bin B 9 , then the change in objective function value is ∆Z = 2[(y 9 -y 7 )s 2 + (y 9 -y 1 )s 3 -s 2 s 3 ] = 2(-4 + 988β).
For 0 ≤ β < 0.004, the term in parentheses above is negative so that it is advantageous to perform this grouping. In other words, first applying the heuristic to obtain A′ 4 and then using an improving step (grouping) does indeed yield the optimal assignment A 4 over the range.
To summarize, for Example 1 the heuristic procedure obtained the optimal solution for all ranges. In Example 2, there was one instance in which the heuristic gave a suboptimal solution. However in this case, a single job switch (grouping) was sufficient to produce the optimal solution. This encouraging success suggests a hybrid heuristic procedure that first carries out the bin packing algorithm, followed by local improvements using selected job switches. In particular, it is straightforward to check for improving moves and swaps; it is a bit more tedious to evaluate all groupings and ungroupings relative to a given assignment.

Conclusions
In this paper we have considered the scheduling of different types of aircraft maintenance programs. In this initial study, we concentrated on the simplified case where there is only a single hangar. First we developed an MIQP to model this scheduling problem. The MIQP model incorporates a parameter β that reflects the relative importance of a good fit versus assigning of jobs to bins occurring earlier in time. For small test problems, it is possible to obtain an exact schedule by solving the MIQP model. Results from these test problems indicated that there were relatively few optimal schedules over the range of all possible values of the specified parameter β. We also developed local optimality conditions for certain types of job switches, relative to a given assignment. The local optimality conditions enable us to improve a given feasible job schedule. Because exact solution of the MIQP is limited to fairly small instances, we developed a simplified heuristic procedure based on existing bin packing techniques. An area for future research is to combine the bin packing heuristic with the intelligent application of the local optimality conditions.