An Improved MILP Formulation for Scheduling Multiproduct

Ind. Eng. Chem. Res. 2003, 42, 2365-2380

2365

An Improved MILP Formulation for Scheduling Multiproduct, Multistage Batch Plants S. Gupta and I. A. Karimi* Department of Chemical & Environmental Engineering, National University of Singapore, 4 Engineering Drive 4, Singapore 117576

Many multiproduct batch plants in the chemical industry employ multiple stages of nonidentical parallel units and operate on the basis of customer orders with different delivery dates. In this work, we present a new continuous-time mixed integer linear programming formulation without using time slots for the short-term scheduling of such plants. Our formulation allows both sequence-dependent and unit-dependent setup times and common operational considerations such as initial plant state and order/unit release times. We develop several novel constraints for the assignment of consecutive orders on a single unit and evaluate them thoroughly to identify the best constraints. Finally, we solve several examples to demonstrate the superiority of our proposed formulations. In comparison to existing works, our formulation requires roughly 30% fewer constraints, yields superior schedules, and reduces the computational times by 65%. We highlight the impact of M in big-M constraints on solution times and propose an industrially more realistic and computationally more efficient (99.8% reduction in solution time) scheduling objective of minimizing tardiness. Introduction Multiproduct batch plants with multiple stages of operation having parallel, nonidentical units in each stage are quite common in the batch chemical industry. Many of these plants operate on order-driven production with a unique due date for each order. In such plants, customer service (order delivery by due date) is the most important consideration, and scheduling of plant operations is a routine activity. However, because of the many alternate ways in which orders can be assigned to various units and produced in different sequences, the task of optimal scheduling is formidable. In fact, this study was motivated by exactly such a plant that produces electronic materials. General reviews on scheduling in batch chemical plants have been published by Reklaitis,1 Zentner and Pekny,2 and Shah.3 One key consideration in solving scheduling problems is the representation of the time domain. Two types of representation have been used: continuous and discrete. Kondili et al.4 developed a mixed integer linear programming (MILP) formulation using a discrete-time representation and introduced the concept of the statetask network (STN) to represent and schedule a general batch plant. Their formulation required significant computational times even for small problems. Shah et al.5 reformulated their model and reduced the computational time by a factor of 3. The main difficulty in the discrete-time formulation is the size of the resulting MILP problem. This limits their utility in solving largescale industrial problems. Recently, much work has been reported on continuous-time representations. Zhang and Sargent6 extended the concept of the resource-task network (RTN), proposed by Pantelides,7 and developed a continuous-time mixed integer nonlinear programming (MINLP) formu* Corresponding author. E-mail: [email protected].

lation for multipurpose batch/continuous plants. To solve this MINLP, they used linearization techniques to obtain a difficult-to-solve MILP. Schilling and Pantelides8 applied a special branch-and-bound algorithm that branches on both discrete and continuous variables to solve a continuous-time MILP formulation based on the RTN representation. Ierapetritou and Floudas9 presented a continuous-time MILP formulation for short-term scheduling in multipurpose batch plants represented by STNs. They reduced the number of binary variables as compared to the formulations of Zhang and Sargent6 and Schilling and Pantelides8 by decoupling the task events from the unit events. In addition to the above efforts for general batch plants, several works have addressed specific configurations to exploit problem-specific features for enhanced efficiency. For instance, Cerda et al.10 developed a continuous-time MILP formulation for the short-term scheduling of a single-stage batch plant with parallel units using a single set of variables for order-unit assignments and order sequencing. Later, for a structurally similar plant, Mendez et al.11 used separate variables for assignment and sequencing to achieve a reduction in the number of binary variables. Scheduling in a multistage plant is much more complex than that in a single-stage plant. Not much attention has been paid to scheduling in plants with multiple stages, despite its industrial significance. Pinto and Grossmann12 developed a continuous-time MILP model for short-term scheduling in multistage batch plants with parallel nonidentical units. They used parallel time axes for units and orders, and using different time slots, defined a set of four-index binary variables (stage-order-slot-unit) for unit-order assignments and order sequencing. They used minimizing earliness as the problem objective and did not consider sequence-dependent changeovers. Because their model could not solve industrial-size problems in reasonable

10.1021/ie020180g CCC: $25.00 © 2003 American Chemical Society Published on Web 05/02/2003

2366

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003

Figure 1. Schematic of a multistage batch plant with parallel units

amounts of time, they developed13 an alternate MILP model with a preordering heuristic to reduce the solution times. Following this work, Hui and Gupta14 presented a continuous-time formulation for multistage, multiproduct batch plants. They used three-index binary variables (order-order-stage) for order sequencing and did not use time slots. They addressed sequence-dependent changeovers and sought schedules with minimum earliness. Even they were not able to solve large problems and could solve only problems with five orders to optimality in reasonable times. In a later work, Hui et al.15 used the same model with slight modifications to address unit-dependent setup times, but still minimized earliness. All of these works on multistage batch plants did not consider the initial state of the plant and assumed that all units and orders were available at the start of the scheduling period. In this paper, we propose a new and much improved continuous-time MILP formulation for the short-term scheduling of a multistage batch plant with nonidentical parallel units. In contrast to the previous works12,14,15 on multistage plant scheduling, we propose minimizing tardiness as the scheduling objective because of its several advantages, as discussed later. We consider a general scenario in which the order changeovers can depend on both the sequence of orders and the processing unit. Our model incorporates the release times of units and orders and takes into account the initial state of the plant. The main contribution of this work is that it presents several models with new approaches for assigning consecutive orders to a single unit. The best of these approaches uses significantly fewer constraints and gives better objective values than the previous works12,14,15 with significantly less effort. We begin the paper with a detailed problem description and then develop several MILP formulations. Next, we test these alternate MILP formulations to identify the best among them by means of an exhaustive numerical evaluation. Finally, we compare our best model with the existing previous models12,14,15 to demonstrate the superiority of our model. Problem Description Figure 1 shows a schematic diagram of a multistage, multiproduct plant with S stages (s ) 1, ..., S). Let Js denote the set of parallel, nonidentical, batch or semicontinuous processing units in stage s and J denote the set of all units (j ) 1, ..., J) in the plant. The units within a stage are nonidentical in the sense that their suitabilities for processing products and their processing rates are different. The plant production is driven by product orders. In this paper, we address the problem of scheduling plant

operation to satisfy I distinct customer orders (i ) 1, ..., I) with known deterministic due dates di (i ) 1, ..., I). Let I denote the set of all orders, Ij denote the subset (Ij ⊂ I) of orders that a unit j (j ) 1, ..., J) can process, and Jis denote the subset (Jis ⊂ Js) of units that can process order i in stage s. Every changeover of orders on a unit requires some setup. This might involve cleaning and/or physical/ operational changes in the unit. In this paper, we address the most general scenario in which setup times are both sequence-dependent (sequence of orders) and unit-dependent. Because of product incompatibilities, some restrictions on the sequencing of orders on units exist. Let FIi denote the set of orders that cannot follow order i on any unit. Because the units in a stage are nonidentical, it is possible that two orders i and i′ might need completely different subsets of units in a stage s, i.e., Jis ∩ Ji′s ) L or the null set. Let NCis denote the subset of orders that cannot use any unit suitable for processing order i in stage s. In addition to the above, we make the following assumptions: (1) A unit cannot process more than one order at a time, and an order cannot be processed by more than one unit at a time. (2) Processing is nonpreemptive. (3) Unlimited intermediate storage (UIS) is available after each stage. (4) Batch transfer times are negligible and are included in the processing times. (5) Processing units do not fail, and processed batches are always satisfactory. (6) Due dates di (i ) 1, ..., I), processing times, and setup times are deterministic and known a priori. (7) Time zero denotes the start of the current scheduling period. (8) The states of all units (orders they are processing and the extents of their processing) in the plant at time zero are known. (9) A unit j can start processing an order i ∈ Ij only after time URTj, called the unit release time. Similarly, the processing of an order i ∈ I can begin only after time ORTi, called the order release time. As stated earlier, the aim is to schedule the production of I orders in the plant subject to the above features and restrictions. We assume that the scheduling goal is to maximize customer satisfaction in terms of order delivery by the due dates. In the next section, we develop an improved MILP formulation for this problem. The formulation uses a continuous-time representation, but no time slots. Mathematical Formulation This scheduling problem entails three principal decisions: (1) assigning each order to a specific unit in each stage of its processing, (2) sequencing orders on each unit, and (3) determining the start times for processing of orders in each stage of its processing. We now define variables and constraints that enable these three decisions. Throughout this section, unless a specific range is specified for an index in a constraint, it is assumed that the constraint is written for all possible values of that index. We begin with variables and constraints that assign orders to various processing units. Order Assignment. To assign orders to specific units for processing, we define a binary variable Zij as follows

Zij )

{

1 if order i is processed on unit j 0 otherwise

If an order i must be processed in stage s ∈ Si, then it must be assigned to exactly one unit j in stage s, i.e., j

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2367

NCis} as the set of feasible successors of i in stage s, we write

∈ Jis. Thus, we must have

∑ Zij ) 1

s ∈ Si

(1)

∑ Xii′s e 1

j∈Jis

i ∈ Is

(4)

i′∈FSis

Appendix A shows that Zij, although defined as a binary, can be treated as a continuous variable with bounds 0 e Zij e 1. Next, we sequence orders for processing on each unit. Order Sequencing. From all of the orders assigned to a unit j, one must be processed first. To select the first order, we define the binary variable ZFij as follows

ZFij )

{

1 if order i is processed first on unit j 0 otherwise

Because no more than one order can be first on unit j, we have

ZFij e 1 ∑ i∈I

(2)

j

The equality will hold true for units that process at least one order. If an order i′ is processed in stage s and it is not the first order on any unit in stage s, then it must follow another order i on some unit in stage s. To ensure this condition, we define another binary variable Xii′s as follows

Xii′s )

{

1 if order i′ is processed after order i on some unit in stage s 0 otherwise

2(Xii′s + Xi′is) +

∑ Xi'is + j∈J ∑ ZFij ) 1

i ∈ Is

is

for all pairs (i, i′) ∈ Is. When units cannot process all orders or when some orders are incompatible with each other, the above constraint need not be written for all pairs (i, i′) ∈ Is, as we explain next. First, an order i cannot follow order i′ in a stage s if no unit in stage s can process them both, i.e., Jis ∩ Ji′s ) L (null set). Earlier, we defined NCis as the set of all such infeasible orders for order i in stage s. Thus, we can say that i′ ∉ NCis in the last constraint. Furthermore, order i might also not follow order i′ in a stage as a result of product incompatibility. We defined FIi′- as the set of all such incompatible orders for order i′. Hence, we can also say that i ∉ FIi′- in the last constraint. This means that each order i has a set of feasible predecessors and feasible successors in each stage of its processing. Defining FPis ) {i′|i′ * i, i ∉ FIi′-, and i′ ∉ NCis} as the set of feasible predecessors of order i in stage s, we rewrite the last constraint as

∑ Xi′is + j∈J ∑ ZFij ) 1

i′∈FPis

i ∈ Is

∑

j∈(Jis-Ji's)

Now there are only two possible scenarios for an order i processed in stage s. Either it is processed first on some unit, or it follows another order i′ on that unit. In other words

i′∈Is

The above inequality allows an order i to process last on some unit in stage s, as last orders cannot have any successors. Although eqs 3 and 4 enable orders i and i′ to be consecutive in stage s, they do not assign them to a single unit. In the next subsection, we derive constraints that ensure that if i and i′ are consecutive orders in stage s, then they are processed on the same unit in stage s. Unit Assignment. To assign consecutive orders i and i′ to a single unit j in stage s, we first ensure that they cannot be assigned to units that cannot process them both. Such units are given by j ∈ (Jis - Ji′s) ∪ (Ji′s - Jis), and we call them uncommon units. This is same as ensuring that consecutive orders i and i′ can be assigned only to units that can process them both. Such units are given by j ∈ Jis ∩ Ji′s, and we call them common units. Even after ensuring that the consecutive orders are assigned to common units, it is possible that they might be processed on different common units. Therefore, we need additional constraints to force consecutive orders i and i′ to one common unit j. Note that two orders i and i′ are consecutive in stage s if Xii′s + Xi′is ) 1. One way to ensure that they cannot be assigned to uncommon units in stage s is to use

(3)

is

Now, an order i cannot have more than one successor i′ in a stage s. Also, i′ must be a feasible successor of i. Therefore, defining FSis ) {i′|i′ * i, i′ ∉ FIi , and i′ ∉

Zij +

∑

Zi′j e 2

j∈(Ji′s-Jis)

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is (5) Observe that Xii′s + Xi′is ) 1 w Zij ) Zi′j ) 0 for uncommon units from the above constraint. Because eq 5 is symmetric with respect to i and i′, we can write it for just combinations of (i, i′) and not permutations. By doing this, we reduce the number of constraints by onehalf. As mentioned before, we now need another constraint to ensure that consecutive orders i and i′ are processed on one common unit j, i.e., ZijZi′j ) 1 for some j ∈ Jis ∩ Ji′s. One way to achieve this is to use

Zij +

∑

Zij′ e Zi′j + 1 - Xii′s - Xi′is

j′∈(Jis-Ji′s)

j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is (6a) Observe that Xii′s + Xi′is ) 1 will make the second term on the LHS of eq 6a equal to 0 from eq 5, and then eq 6a will reduce to

Zij e Zi'j

j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is

(6b)

We now show that eq 6b forces ZijZi′j ) 1 for one common unit j, whenever one of the following two conditions is met:

1. Zij ) 1 for some unit j ∈ Jis ∩ Ji′s 2. Zi′j ) 1 for some unit j ∈ Jis ∩ Ji′s and Zij ) 0 ∀ j ∈ (Jis - Ji′s) When condition 1 is met, it is easy to see that Zij ) 1 w

2368


Table 1. Various Sets of Unit Assignment Constraints set 1

eq 5

unit assignment constraints

∑


Zij +

j∉(Jis-Ji′s)

6a

Zij +

∑

∑

index ranges Zi′j e 2

max no.a

i′ ∉ NCis, i′ > i, (i, i′) ∉ Is

IC

2

j ∉ Jis ∩ Ji′s, i′ > i, (i, i′) ∉ Is

IC

2

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is

IC

2

j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is

IC

2

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is

IC

2

j ∉ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is

IC

2

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is

IC

2

j ∈ Jis ∩ Ji′s, i′ ∉ FIi , i′ * i, (i, i′) ∈ Is

IP

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is

IC


IP

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is

IC

2

j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is

IC

2

i′ ∈ FSis, i ∈ Is

IP

2


IP

2

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is

IC

2

j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is

IC

2


IP

2


IP

2

j ∈ Jis, i′ ∉ NCis, i′ ∉ FIi , i′ * i, (i, i′) ∈ Is

IP

2

j∉(Ji′s-Jis)

Zij′ e Zi′j + 1 - Xii′s - Xi′is

j′∉(Jis-Ji′s)

2

5

∑


Zij +

j∉(Jis-Ji′s)

6c 3

5

4

5

∑


Zij +

Zi′j e 2

Zi′j e Zij + 1 - Xii′s - Xi′is

∑


Zij +

6e

Zij e Zi′j + 1 - Xii′s

5


∑

6f

Zi′j e Zij + 1 - Xii′s

7

Xii′s + Xi′is +

∑

∑

Zi′j e 2

j∉(Ji′s-Jis)

Zij +

j∉(Jis-Ji′s)

6

∑

j∉(Ji′s-Jis)

j∉(Jis-Ji′s)

5

Zi′j e 2

Zij e Zi′j + 1 - Xii′s - Xi′is j∉(Jis-Ji′s)

6d

∑

j∉(Ji′s-Jis)

∑

Zi′j e 2

2

2

j∉(Ji′s-Jis)

Zij e 1

2

j∉(Jis-Ji′s)

7

6c

Zij e Zi′j + 1 - Xii′s - Xi′is

8

Xii′s +

∑

Zij e 1

j∉(Jis-Ji′s)

8

6e


9

Xii′s + Xi′is +

∑

Zi′j e 1

j∉(Ji′s-Jis)

6d 9

10

Zi′j e Zij + 1 - Xii′s - Xi′is Xii′s +

∑

Zi′j e 1

j∉(Ji′s-Jis)

10

6f


22b

Xii′s + Zij +

∑Z

i′j′

e2

j′∉Ji′s j′*j

a Maximum numbers of constraints for each common unit. C refers to combinations and P to permutations of (i, i′). b Corrected form of eq 1 in Hui and Gupta,14 as discussed in Appendix C.

Zi′j ) 1 from eq 6b. Now, when condition 2 is met, Zi′j ) 1 w Zi′j′ ) 0, ∀ j′ ∈ Ji′s, j′ * j, from eq 1. This makes Zij′ ) 0, ∀ j′ ∈ Jis ∩ Ji′s, j′ * j, from eq 6b. Then, Zij ) 0, ∀ j ∈ (Jis - Ji′s), will make Zij ) 1 from eq 1. Thus, eq 6a ensures that ZijZi′j ) 1 for one common unit j, as long as one of the above two conditions is met. Because eq 5 meets both conditions, eq 5 and eq 6a form a valid constraint set to force consecutive orders i and i′ onto a single unit j. Note that eq 6a becomes relaxed for Xii′s + Xi′is ) 0 and it is enough to write eq 6a for combinations of i and i′ only and not permutations. The constraint set (eqs 5 and 6a) presented above is just one valid set of unit assignment constraints. It turns out that there are several, in fact nine, different ways of writing these constraints. For the sake of brevity, we derive the additional eight sets in Appendix B and list them in Table 1. Each set comprises two constraints. The first constraint forces Zij ) 0 or

Zi′j ) 0 or both for uncommon units, while the second forces Zij ) Zi′j ) 1 for one common unit j. The various sets differ in the numbers of actual constraints, the numbers of nonzeros in the constraints, and the way in which the Zij’s are set equal to 0 or 1. As all of these factors impact the computation times, the decision on which set is the best will be based on a detailed numerical evaluation in a later section. Having assigned and sequenced orders on specific units and having ensured that consecutive orders i and i′ are processed on one common unit, we proceed to set the start times for the processing of orders in each stage of their production. Order Timing. Let Tis denote the time at which order i starts processing in stage s and tij denote the duration of its processing on unit j. Furthermore, if i is being processed in stage s, then let nsis denote its next stage of processing. Because i must finish processing


in stage s before it can start processing in stage nsis, we have

Tis′ g Tis +

∑ Zijtij

s′ ) nsis, s ∈ Si

(11)

j∈Jis

M(1 - Xii's) + Ti′s g Tis +

∑ Zij(tij + cii′j) j∈J is

i′ ∈ FSis, i ∈ Is (12) where M is a large positive number that relaxes the constraint when Xii′s ) 0. Note that the above constraint is written for all permutations of i and i′. A unit j can start processing orders only after its release time URTj. Now, when a unit j is released, it might require a setup time to start processing an order i from a new production schedule. Let c0ij denote this time. Because the first order i on unit j (ZFij ) 1) can start processing only after the unit is released and set up to process i, we obtain

∑ ZFij(URTj + c0ij)

i ∈ Is

(13)

j∈Jis

Note that the above constraint gets relaxed for ZFij ) 0. Finally, the processing of an order i can begin only after its release time ORTi. If fsi denotes the first stage of processing of an order i, then we must have

Tis g ORTi

s ) fsi

(14)

Miscellaneous. It is clear that an order i can be processed first on a unit j only if it is assigned to that unit. In other words, ZFij ) 1 w Zij ) 1, or

Zij g ZFij

i ∈ Ij

(15)

This completes the development of constraints for our mathematical formulation. We now determine an objective function for the problem. Scheduling Objective. Production of an order i is complete when the order finishes being processed in its last stage. Let lsi denote the last stage of processing for order i. The completion time of order i in stage lsi is the sum of its start time in stage lsi and its processing time. To compute the delay in the completion of order i, we define Tdi by the constraint

Tdi g [Tis +

maximize

[ ∑ wis(Tis + ∑ Zijtij) - N × Tdi] ∑ i∈I s∈S j∈J i

As stated in the problem description, a unit j requires a setup to change over from order i to order i′. Let cii′j denote the setup time required for a changeover from order i to order i′ on unit j. Now, if an order i′ follows an order i on a unit j of stage s, i.e., Xii′s ) 1, then order i′ can start processing only after order i finishes processing and the unit is setup to process order i′. In other words

Tis g

completions. Indeed, this is the scheduling objective used in the previous works12,14,15 on this problem, i.e.

∑ Zij tij] - di j∈J

s ) lsi

(16)

(17)

is

where N is a large positive number and wis is defined by

0.2s max di i wis ) di

(

)

Equation 17 maximizes the completion times of orders in each stage of their processing and imposes a penalty on the tardy orders through the last term. This produces orders as close as possible to their due dates (just-intime scheduling), so it is equivalent to minimizing earliness. As we discuss later, there are several problems with this objective. First, in many cases, an order can be shipped as soon as it is finished, so it does not incur inventory cost. Second, the schedules obtained using eq 17 exhibit equipment idle times in the beginning of the scheduling period. Because most schedules invariably get revised after a few days because of changes in orders or other reasons, such idle times are not desirable because they reduce operating flexibility. Finally, as we see later, problems with this objective require more computation time than those with the objective of minimizing tardiness. Therefore, we propose the following objective function for our problem

minimize

pi × Tdi ∑ i∈I

(18)

where the parameter pi essentially reflects the “priority” or “importance” of order i. The above objective minimizes the total weighted tardiness of order completions. Note that the absolute values of the pi are not important at all. Only the relative values count, and in practice, the plant personnel can usually assign these in many different ways. In fact, the plant personnel can use a different set of pi from time to time. Although it is not advisable to suggest a strict guideline for selecting the values of the pi, one possible approach is as follows. Because each order has a dollar value (absolute value or marginal profit) and a customer, pi should reflect the relative value of the order, as well as the relative importance of the customer. Thus, one can assign a number (say, on a scale of 1-5) to each customer first. Let us call this customer priority, which indicates how valued the customer is to the company. Then, one can compute a value priority for each order by dividing the value of that order by the minimum value among all orders. Then, one can set

pi ) customer priority × value priority

is

In an industrial scheduling problem having varying due dates for customer orders, the goal is normally to finish processing of as many orders as possible before their due dates, i.e., minimize the tardiness of order completions. However, if orders are finished before their due dates, they might incur an inventory cost. To reduce this cost, one could also minimize the earliness of order

Equation 18 might produce a schedule with more earliness than eq 17, but it will eliminate costly equipment idle times in the beginning of the scheduling period. If orders complete processing much before their due dates, then either the customer lead time has been overestimated or the plant has excess production capacity. In any case, eq 18 will produce schedules with minimum tardiness, with reasonable earliness, and also

2370


Table 2. Processing Times (tij, h) for Orders on Various Units order i unit j

1

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25

18.1 18.1 18.1 18.1 18.1 18.1 14 5 5 12 12 12 12 12 12 12 12 -

23 23 23 23 23 23 14 5 5 12 12 12 12 12 12 12 12 12 12 100 48

18.1 18.1 18.1 18.1 18.1 18.1 14 5 5 24 24 24 24 24 24 24 24 9.3 9.3 24 24 -

20 20 20 20 20 20 11 5 5 12 12 12 12 12 12 12 12 7.9 7.9 24 24 -

17 17 17 17 17 17 14 5 5 12 12 12 12 12 12 12 12 12.5 12.5 24 24 -

15 14 13 12 18 15 15 7 7 13 12 15 17 17 18 19 16 13 13.5 14 14.5 12 12 23

31 30 34 32 31 31 31 31 14 15 16 41 15 81 16 16 21 12 13 11 11 11 11

12 12 14 15 16 16 15 12 13 14 14 14 14 14 14 14 10 10 9 8

13 7 8 9 16 15 15 11 13 7 8 9 16 15 15 11 12 15 17 17 22 21 -

12 4 23 12 14 13 13 13 12 4 23 12 14 13 13 13 6 12 23 12 22 12

8 9 13 -a 13 13 35 23 8 9 13 13 13 35 23 8 5 12 23 12 22 12

10 7 14 13 13 23 13 10 7 14 13 13 23 13 10 12 23 12 21 16

18.1 18.1 18.1 18.1 18.1 18.1 14 5 5 12 12 12 12 12 12 12 12 9.5 9.5 24 24 -

23 23 23 23 23 23 14 5 5 12 12 12 12 12 12 12 12 12 12 10 10

18.1 18.1 18.1 18.1 18.1 18.1 14 5 5 24 24 24 24 24 24 24 24 9.3 9.3 24 24 -

20 20 20 20 20 20 11 5 5 12 12 12 12 12 12 12 12 7.9 7.9 24 24 -

17 17 17 17 17 17 14 5 5 12 12 12 12 12 12 12 12 12.5 12.5 24 24 -

15 14 13 12 18 15 15 7 7 13 12 15 17 17 18 19 16 13 13.5 14 14.5 12 12 23

31 30 34 32 31 31 31 31 14 15 16 41 15 81 16 16 21 12 13 11 11 11 11

12 12 14 15 16 16 15 12 13 14 14 14 14 14 14 14 10 10 9 8 7 7

13 7 8 9 16 15 15 11 13 7 8 9 16 15 15 11 12 15 17 17 22 21 -

12 4 23 12 14 13 13 13 12 4 23 12 14 13 13 13 6 12 23 12 22 12

a

9.5 9.5 24 24 -

7 7

Order cannot be processed on this unit.

with minimum equipment idle times. Moreover, as shown later, it will also reduce the CPU time significantly. With this in mind, the MILP formulation for our scheduling problem can be stated as minimization of the weighted tardiness (eq 18) subject to eqs 1-4 and 1116 and one unit assignment constraint set from Table 1. Each set from Table 1 will result in a unique model. We call these models M1, M2, M3, ..., M9. For instance, M5 will use the set 5 constraints from Table 1 along with eqs 1-4 and 11-16. As stated earlier, it is not clear a priori which model is the best among M1-M9, as they use different numbers of constraints and nonzeros and require different CPU times to solve a given problem. Therefore, to identify the best model, we performed a thorough numerical evaluation of models M1-M9 using problems of varying sizes. Model Evaluation To evaluate the nine models, we used the example of Pinto and Grossmann,12 which Hui and Gupta14 also used. Whereas Pinto and Grossmann12 used unitdependent setup times, Hui and Gupta14 used sequencedependent setup times. The example involves a process with five stages (S ) 5, s ) 1-5) with parallel nonidentical units. The units are labeled 1-6 in stage 1, 7-9 in stage 2, 10-19 in stage 3, 20-22 in stage 4, and 2325 in stage 5. All orders pass through stages 1-2-34-5 in that sequence. Tables 2-5 present the processing times for orders, sequence-dependent setup times, unitdependent setup times and order due dates, respectively. All units and orders are available at the start of the new scheduling period, i.e., URTj ) ORTi ) 0 ∀ i, j. In the case of sequence-dependent setup times, a unit does not require any setup to process the first order, i.e., c0ij ) 0 ∀ i, j. To have more problems for model evaluation, we vary the numbers of orders and the scheduling objective in

the above example to form three sets of test problems. We use a notation including both the test problem set and the number of orders to identify these test problems. For instance, a problem in set 2 with 8 orders is called problem 2-8. Set 1. We use four test problems (labeled 1-5, 1-8, 1-10, and 1-12) with 5, 8, 10, and 12 orders, respectively, for this set. For these problems, we assume sequencedependent setups and minimize earliness of order completions (eq 17). Set 2. We again use four problems (labeled 2-10, 2-12, 2-17, and 2-22), but with more orders (10, 12, 17, and 22 respectively). We still use sequence-dependent setups, but we now minimize tardiness of order completions (eq 18). Because the previous work did not use 17 and 22 orders, we assumed additional data for problems 2-17 and 2-22. Moreover, we assumed all pi’s in eq 18 to be unity. Set 3. We use only two test problems (labeled 3-5 and 3-8) with 5 and 8 orders, respectively, for this set. In contrast to sets 1 and 2, we assume unit-dependent setups, but ew still minimize earliness of order completions (eq 17). To make a fair comparison with previous work,12,15 which solved problem 3-5, we also assume a due date of 500 h for all orders in problem 3-5. Because the above problems assume either sequencedependent or unit-dependent setups, which are special cases of the general models M1-M9, we first simplify M1-M9 slightly as follows. For sets 1 and 2, the setup times depend only on the sequence of orders, so cii′j becomes cii′, and eq 12 is replaced in all models by

M(1 - Xii′s) + Ti′s g Tis +

∑ Zijtij + cii′

j∈Jis

i′ ∈ FSis, i ∈ Is (19) Also, sets 1 and 2 have URTj ) ORTi ) c0ij ) 0 ∀ i, j, so eqs 13 and 14 disappear.

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2371 Table 3. Sequence-Dependent Setup Times (cii′, h) for Problems in Sets 1 and 2 order i′ order i 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 a

1 7 12 11 13 7 -a 8 4 5 5 4 5 7 12 11 13 7 8 4 5

2

3

4

5

6

7

8

9

10

11

12

13

14

15

16

17

18

19

20

21

22

7

8 6

9 5 7

10 8 6 2

5 4 4 8 2

7 5 5 6 1 9

8 6 7 6 6 7 8

9 5 6 7 5 8 4 8

7 4 5 7 8 7 4 9 7

8 7 6 6 7 7 5 8 7 8

9 8 6 7 8 8 6 8 9 8

5 4 4 8 2 7 8 4 5 5 4

7 5 5 6 1 2 3 7 5 6 6 5 7

8 6 7 6 6 5 4 8 6 7 6 6 8 6

9 5 6 7 5 6 5 9 7 8 7 7 9 5 6

7 4 5 7 8 7 6 4 8 9 6 8 7 4 5 7

5 4 4 8 2 7 6 5 9 4 7 7 5 4 4 8 2

7 5 5 6 1 9 4 5 6 5 8 5 7 5 5 6 1 9

8 6 7 6 6 7 8 5 8 6 8 6 7 6 6 7 8

9 5 6 7 5 8 4 8 6 8 9 7 9 5 6 7 5 8 4 8

7 4 5 7 8 7 4 9 7 5 8 7 7 4 5 7 8 7 4 9 7

8 9 4 2 3 7 5 6 6 5 7 5 8 9 4 2 3 7 5 6

5 4 5 4 8 6 7 6 6 8 6 1 5 4 5 4 8 6 7

6 6 5 9 7 8 7 7 9 5 5 5 6 6 5 9 7 8

7 6 4 8 9 6 8 10 8 6 2 3 7 6 4 8 9

6 5 9 4 7 7 5 4 4 8 2 5 6 5 9 4

5 6 5 8 5 7 5 5 6 1 9 5 5 6 5

8 6 8 6 7 6 6 7 8 5 -

8 9 7 9 5 6 7 5 8 19 8 5 8

8 7 7 4 5 7 8 7 4 9 7 5

8 8 7 6 6 7 7 4 8 7 8

9 8 6 7 8 8 5 6 8 9

4 4 8 2 7 8 4 5

5 6 1 2 3 7 5 6

6 6 5 4 8 6 7

5 6 5 9 7 8

7 6 4 8 9

6 5 9 4

5 6 5

-

8

Infeasible sequence of orders.

Table 4. Unit-Dependent Setup Times (cj, h) for Problems in Set 3 unit j

1-7

8-9

10-19

20-21

22

23-24

25

cj (h)

8

1

2.5

6

24

4

5

Table 5. Due Dates (di, h) of Orders order(s) i

di (h)

2, 13 1 3, 12, 15 4, 6, 16, 20 5, 7, 11 8, 9, 14, 19 21 17 10, 22 18

500 510 520 530 540 550 560 570 580 600

For set 3, the setup times depend only on the units, so cii′j becomes cj, and eq 12 is replaced in all models by

M(1 - Xii′s) + Ti′s g Tis +

∑ Zij(tij + cj)

j∈Jis

i′ ∈ FSis, i ∈ Is (20) Also, c0ij reduces to cj, URTj ) 0 ∀ j simplifies eq 13 as

Tis g

∑ ZFijcj

i ∈ Is

(21)

j∈Jis

and eq 14 disappears. In addition to our nine models (M1-M9), we also include the model of Hui and Gupta14 in our evaluation. We denote their model (as implemented by us) by M10 and use it in sets 1 and 3 only. However, note that we had to modify one constraint (eq 1 in Hui and Gupta14) in their formulation to match the number of reported constraints and more importantly to get a feasible solution (see Appendix C for details). The corrected equation is shown in Table 1. We implemented all 10 models (M1-M10) in GAMS 2.5 v19.6 (Brooke et al.16) and solved all problems using CPLEX 7.0 on a Sun

Enterprise 250 server running Sun OS 5.7 with a single Ultra SPARC II 400-MHz processor having 2 GB of RAM. Numerical Performance. In this subsection, we use the problems in sets 1-3 described earlier to evaluate models M1-M10. On the basis of these results, we select the best scheduling model for the present problem. Recall that M10 is our implementation of the corrected (Appendix C) model of Hui and Gupta.14 Table 6a,b lists the numbers of constraints, nonzeros, and variables (continuous and binary) for all problems. Because the number of variables depends only on the number of orders, all models have the same numbers of binary and continuous variables for a given problem. As the number of orders increases from 5 to 22, the number of binary variables increases 14-fold. This increase in problem size naturally affects the CPU time. For a given problem, the numbers of nonzeros and constraints differ among models. For instance, models M7 and M9 have the most constraints, whereas models M1, M2, M3, M6, and M8 have the fewest. In fact, the former have 40-60% more constraints than the latter. Similarly, models M6 and M8 have the fewest nonzeros, whereas model M10 has almost double that number. Thus, models M6 and M8 have the smallest numbers of both nonzeros and constraints. Also note that the pairs of models M2, M3 and M4, M5 have the same numbers of nonzeros and constraints. During our numerical work, we noticed surprisingly that solution times depended heavily on the value of M used in eq 12. It is worth noting that the solution times for the same problem can differ by as much as 10 h for two different values of M. Thus, even though M can be any large positive value, it has to be selected judiciously. The problem, however, is that the solution times do not follow any definite trend with the value of M, so that it is impossible to fix an optimal value of M a priori. We observed that models with very high and very low values of M required large solution times or failed to even locate an optimal solution. Moreover, this dramatic impact of M makes it imperative that a model be judged on the basis of its average performance over several values of M and not just one. Therefore, in our work,

2372


Table 6. Statistics for Various Models and Problems (a) Numbers of Nonzeros and Constraints for Various Models and Problems model problem

statistics

M1

M2

M3

M4

M5

M6

M7

M8

M9

M10

1-5

nonzeros constraints nonzeros constraints nonzeros constraints nonzeros constraints nonzeros constraints nonzeros constraints nonzeros constraints nonzeros constraints nonzeros constraints nonzeros constraints

2548 549 6549 1294 10 100 1948 14 546 2711 9832 1948 1 4227 2711 27 799 5210 46 558 8658 -

2500 549 6333 1294 9675 1948 13 701 2711 9407 1948 13 382 2711 26 343 5210 44 201 8658 -

2500 549 6333 1294 9675 1948 13 701 2711 9407 1948 13 382 2711 26 343 5210 44 201 8658 -

2888 743 7381 1818 11 299 2760 16 033 3877 11 031 2760 15 714 3877 31 285 7681 52 513 12 814 -

2888 743 7381 1818 11 299 2760 16 033 3877 11 031 2760 15 714 3877 31 285 7681 52 513 12 814 -

2498 549 6268 1294 9566 1948 13 534 2711 9298 1948 13 215 2711 26 069 5210 43 650 8658 2628 574 6482 1334

2888 785 7377 1945 11 291 2961 16 020 4178 11 023 2961 15 701 4178 31 268 8302 52 471 13 863 -

2494 549 6290 1294 9591 1948 13 532 2711 9323 1948 13 213 2711 26 029 5210 43 699 8658 2624 574 6504 1334

2888 785 7379 1945 11 293 2961 16 025 4178 11 025 2961 15 706 4178 31 275 8302 52 485 13 863 -

4440 709 12 057 1790 18 711 2729 26 567 3879 -a 4440 709 12 057 1790

1-8 1-10 1-12 2-10 2-12 2-17 2-22 3-5 3-8

(b) Numbers of Continuous and Binary Variables for Problems variables

a

number of orders

continuous

binary

5 8 10 12 17 22

136 223 279 332 467 610

189 433 635 881 1631 2645

Problems not solved in our study.

we used M ) 600, 700, 800, 900, and 1000 for problems in sets 1 and 3 and M ) 1000, 2000, 5000, 10 000, and 20 000 for set 2. For problem 1-10, we tried M ) 5000 and 10 000, but most models failed to give an optimal solution even after 6-7 h. For all problems in sets 1 and 3, we fixed the value of N in eq 17 equal to 10 000. Most published MILP models using big-M constraints have overlooked this important point. Two previous works14,15 on this problem also do not report M, whereas Pinto and Grossmann12 obtained it from the problem data by taking the value of M as the maximum of all of the due dates of orders. Thus, it becomes difficult to compare models if the value of M is not known. We begin our evaluation of models with problem set 1. Set 1. Recall that set 1 problems have minimizing earliness as their scheduling objective. Table 7 presents the results for set 1. Note that the solution time increases from a few seconds to a few hours as the number of orders increases from 5 to 10. In fact, all models fail to reach the optimal solution even in 5 h for problem 1-12. Thus, for problem 1-12, we terminate the solver after 5000 s and obtain a suboptimal solution for all models. To compare model performances for a given problem, we define average CPU time, relative CPU time, and average relative CPU time as follows

average CPU time ) average with respect to M of CPU times relative CPU time ) CPU time of a model lowest CPU time among all models

average relative CPU time ) average with respect to M of relative CPU times Furthermore, we define relative CPU time for a model as

model relative CPU time ) average with respect to all problems in a set of the average relative CPU times of a model In addition to these quantities, we define the model relative gap for a given problem as the average with respect to M of relative gaps. We use the last measure especially for problem 1-12, as no model can solve this problem to optimality. On the basis of the above criteria, we rank the models as shown in Table 7. Note that, for a small problem such as problem 1-5, it is difficult to say which model is the best, because the average CPU times are less than 1 s. For larger problems such as problems 1-8 and 1-10, we can make a clear distinction among models. For problem 1-8, models M6 and M8 perform the best. M6 takes about 33 s on an average, whereas M8 takes about 76 s. All other models need average CPU times of up to 15 times that of M6. The best (over M) CPU times for M6 and M8 are 23 and 38 s, respectively, which are reasonably close to each other. For problem 1-10, model M6 has an average CPU time of 3458 s, whereas model M8 has 3212 s. All other models need average CPU times of up to 7 times that of M8. In fact, some models (M2, M3, M5, M7, and M10) fail to locate an optimal solution even after 5-10 h (indicated in bold) for some values of M. Model M8 requires the lowest CPU time of 866 s for all 50 (10

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2373 Table 7. Computational Times (s) for Problems 1-5, 1-8, and 1-10 and Relative Gaps (%) for Problem 1-12 CPU Times (s) model M 600 700 800 900 1000 average avg relative

M1

M2

0.8 0.8 0.8 0.6 0.6 0.726 1.33

M3

0.7 0.6 0.7 0.7 0.6 0.638 1.17

M4

M5

M6

0.6 0.6 0.6 0.8 0.6 0.62 1.13

problem 1-5 1.0 0.7 1.0 0.7 1.0 0.8 0.9 1.0 0.8 0.9 0.936 0.82 1.71 1.50

0.8 0.7 0.6 0.6 0.5 0.63 1.16

M7 0.9 1.2 1.1 0.9 0.8 0.966 1.78

M8

M9

0.8 0.5 0.7 0.7 0.6 0.65 1.18

0.9 1.2 1.2 1.1 1.1 1.106 2.03

M10 1.5 1.7 2.1 2.1 1.4 1.744 3.18

600 700 800 900 1000 average avg relative

59.9 174.1 196.0 64.4 175.0 133.9 4.54

83.4 61.2 120.9 282.0 54.8 120.46 3.52

240.0 66.8 60.0 189.0 178.0 146.8 4.92

problem 1-8 102.6 242.0 44.6 110.0 1719.0 89.0 205.0 306.0 393.0 255.0 492.83 200.4 13.00 6.70

30.3 24.2 46.0 42.1 23.0 33.1 1.02

85.9 41.5 122.7 477.0 144.0 174.2 5.17

124.0 40.0 107.0 38.6 74.0 76.71 2.46

1298 79.4 133.0 66.8 198.0 355.03 11.87

59.8 95.5 68.9 199.0 292.0 143 5.05


4212 5739 10 497 23 034 5925 9881 7.86

8739 4731 3588 35 516b 6333 11 781 9.42

1937 6712 35 428 5543 1657 10 255 7.12

problem 1-10 1631 7445 17 463 10 342 1716 14 017 14 943 4613 997 14 722 7350 10 228 6.96 8.65

1910 3391 4682 4393 2916 3458 2.84

29 061 28 017 17 156 19 305 14 806 21 669 18.17

1558 866 4551 6382 2707 3213 2.47

6043 8424 2790 5199 1352 4761.6 4.13

3236 5233 21 695 1287 14 755 9241 7.31

600 700 800 900 1000 average (%) gap rank relative CPU timea time rank

4.99 4.48 8.28 5.42 5.26 5.68 4

5.62 4.72 5.25 4.11 3.58 4.65 1

Relative Gaps (%) for Problem 1-12 after 5000 s of CPU Time 5.50 5.92 4.91 4.92 6.71 6.05 5.97 5.46 2.95 6.56 4.85 5.09 7.13 7.72 6.91 6.88 5.49 4.54 6.65 12.21 7.10 4.81 6.72 6.84 14.09 6.08 5.46 5.75 5.82 9.30 7 3 5 6 110

6.46 3.15 4.67 4.13 5.35 4.75 2

Model Relative CPU Times and Corresponding Ranks Based on Problems 1-5, 1-8, and 1-10 458 470 439 7.23 5.61 1.67 8.37 2.03 4 5 3 9 7 1 10 2

7.20 4.86 11.99 8.88 6.41 7.87 8

6.53 8.12 10.02 5.74 12.59 8.60 9

6.01 8

15.18 6

a Average over problems 1-5, 1-8, and 1-10 of average relative CPU times. b CPU times in bold indicate that the model failed to give an optimal solution in the given time.

models each for five values of M) runs made for problem 1-10. Observe the unpredictable impact of M in that model M3 solves problem 1-10 in 1937 s for M ) 600 but fails to give an optimal solution even after 10 h for M ) 800. Comparing the overall performances (Table 7) of models M6 and M8 for problems 1-8 and 1-10, we see that modelM6 has the lowest model relative CPU time of 1.67 and model M8 has the next best model relative CPU time of 2.03. In comparison, the model relative CPU times of all other models exceed 4. Because we could not solve problem 1-12 to optimality, we compare the model relative gaps of the models after 5000 s. For problem 1-12, models M2 and M8 perform the best. In contrast to its performance for problems 1-8 and 1-10, model M6 ranks 6th out of 10 for problem 1-12. Note that, although model M2 performs the best here, it did badly for problems 1-8 and 1-10. Because model M8 ranks second consistently for all three problems (1-8, 1-10, and 1-12), we consider its overall performance as the best and conclude that M8 is the best model for problems in set 1, with the next best model being M6. Set 2. Recall that the problems in this set have minimizing tardiness as the objective. Table 8 shows the results for set 2. The average CPU times for problems 2-10 and 2-12 are around 5 and 20 s, respec-

tively. Comparing models with such small CPU times could be misleading. Even the model relative CPU times for average performance for problems 2-10 and 2-12 do not clearly point to the best model, so we use the larger problems 2-17 and 2-22. The model relative CPU times for average performance for problems 2-17 and 2-22 show that the four models M3, M6, M8, and M9 perform the best. They all rank first with almost the same model relative CPU time of 1.5. To identify the best among them, we compare the average CPU times for problems 2-17 and 2-22. On the basis of their average CPU times, we can clearly eliminate models M3 and M9, which never perform the best. Model M6 performs best for problem 2-17, and model M8 performs best for problem 2-22. Hence, if we consider the average and relative CPU times together, we can say that M8 and M6 are the best models for solving problems with tardiness as the objective. Set 3. As we mentioned earlier, set 3 problems have unit-dependent setup times and minimizing earliness as the objective. From our results for sets 1 and 2, models M8 and M6 are the best models so far. To distinguish between them, we solve problem 3-8 with models M6 and M8 using several values of M. Table 9 presents the results. As we can see, model M8 outperforms model M6. M8 requires an average CPU time of 1256 s, as compared to 2683 s for M6. Moreover,

2374


Table 8. Computational Times (s) for Problems 2-10, 2-12, 2-17, and 2-22 model M

M1

M2

M3

M5

M6

M7

M8

M9

problem 2-10 6.6 1.9 8.0 6.5 6.3 7.9 8.4 2.1 7.7 9.2 7.4 5.5 5.88 4.35

M4

11.3 1.4 2.0 1.3 10.7 5.4 4.21

12.2 1.1 11.6 8.2 11.1 8.8 6.91

8.9 10.9 12.7 11.1 1.5 9.0 7.31

1.2 1.6 9.8 10.2 10.4 6.6 5.02

11.9 12.4 9.4 8.7 9.4 10.4 8.33

1000 2000 5000 10 000 20 000 average avg relative

11.1 9.2 1.8 9.4 8.5 8.0 6.40

2.0 9.9 1.3 7.0 1.4 4.3 3.52


31.5 22.5 25.4 20.1 20.1 23.9 6.57

15.6 22.0 14.1 3.2 14.3 13.8 3.75

20.2 3.7 14.2 19.3 15.7 14.6 3.75

problem 2-12 31.7 30.3 28.3 17.2 25.6 26.6 7.23

4.3 3.3 3.6 21.9 21.4 10.9 2.37

16.7 16.9 2.5 13.9 14.6 12.9 3.18

21.0 21.0 18.4 20.4 21.8 20.5 5.45

22.4 2.9 25.5 16.6 15.0 16.5 4.54

21.7 20.6 20.2 20.0 17.0 19.9 5.51


87 94 105 124 87 99.4 1.97

66 81 90 107 66 82.0 1.63

97 55 65 114 56 77.4 1.51

problem 2-17 165 106 191 87 120 133.8 2.59

89 89 68 94 105 89.0 1.76

74 73 49 45 45 57.2 1.09

136 102 115 82 85 104.0 1.99

110 77 89 62 111 89.8 1.75

98 112 76 65 70 84.2 1.61


2797 6695 4418 33 902 3833 10 329 28.71

609 809 2479 468 1088 1091 2.71

514 391 793 633 647 596 1.49

problem 2-22 939 698 1857 885 1612 1198 2.95

969 852 713 648 925 821 2.03

1532 1002 402 342 923 840 1.99

4921 606 604 433 495 1412 3.04

719 424 452 546 433 515 1.27

578 351 400 953 618 580 1.47

model relative CPU times and corresponding ranks based on above problems 6.48 3.64 4.82 5.79 3.29 5.05 6.38 15.34 2.17 1.50 2.77 1.89 1.54 2.52 4 3 1 3 2 1 3

4.78 1.51 1

6.92 1.54 1

relative CPU timea relative CPU timeb time rankb

a Average over problems 2-10 and 2-12 of average relative CPU times. b Average over problems 2-17 and 2-22 of average relative CPU times.

Table 9. Computational Times (s) of Models M6, M8, and M10 for Problem 3-8 model M

M6

M8

M10


3351a 4687 2340 1959 1082 2684 2.96

2129 578 1117 973 1487 1257 1.07

3236 5347 3148 4527 5000 4252 4.57

a CPU times in bold indicate that the model failed to give an optimal solution in the given time.

model M6 fails to obtain optimal solutions for M ) 600 and 700 even after 1 h. Model M8 has an average relative CPU time of 1.07 as compared to 2.96 for model M6. This clearly shows that model M8 performs better than model M6 for the problems in set 3. On the basis of the above results and discussion of models M1-M10 for all three sets of test problems, we conclude that M8 is the best model overall. Remarks Our model evaluation exercise points to several useful observations. The first is the huge difference in the ease of solution between the problems in set 1 and those in set 2. Compare the average CPU times for problems

1-10 and 1-12 with those for problems 2-10 and 2-12. Model M8 requires an average (over M) CPU time of 3212 s for problem 1-10 as compared to 7 s for problem 2-10. None of the models M1-M10 could solve problem 1-12 optimally in even 5000 s, but they took only 17 s on an average (over M1-M9) to solve problem 2-12. Moreover, for set 2, we could solve even the larger problems 2-17 and 2-22 in reasonable CPU times. This clearly demonstrates that the problems with tardiness as the scheduling objective are much easier to solve than those with earliness as the objective. This result, combined with several practical considerations discussed later, makes tardiness the most preferred objective for scheduling in multistage batch plants. The second observation is the fact that even though the pairs of models M2, M3 and M4, M5 have the same numbers of nonzeros and constraints, they require very different CPU times for solving a given problem. This is because they differ in the sequencing of the orderd on the units. Finally, the value of M has a huge impact on the solution times of a problem. We can see this for problem 1-10 in Table 7 where the CPU times differ by up to a factor of 20 for two different values of M. Furthermore, there is no definite relation between CPU time and M for a given problem or model. In other words, an optimal a priori choice of M is almost impossible for models with big-M constraints. It is also important to note that M

Ind. Eng. Chem. Res., Vol. 42, No. 11, 2003 2375 Table 10. Comparison of Model M8 with the Model of Hui and Gupta14 for Problems 1-5, 1-8, 1-10, and 1-12

Table 11. Comparison of Model M8 with the Models of Pinto and Grossmann12 and Hui et al.15 for Problem 3-5

model performance measure

Hui and Guptaa

M10b

M8

number of constraints best objective number of nodes number of iterations average CPU time (s)c best CPU time (s)d

problem 1-5 709 7498 3357 8843 87.7

number of constraints best objective number of nodes number of iterations average CPU time (s) best CPU time (s)

problem 1-8 1790 12 329e 14 568 100 000 1268.2

1790 12 386.3 17 470 110 977 143.0 59.8

1294 12 386.3 4400 33 990 76.7 38.6

number of constraints best objective number of nodes number of iterations average CPU time (s) best CPU time (s)

problem 1-10 2729 16 120 8985 100 000 1489

2729 16 344.5 280 347 3 089 585 9241 1287

1948 16 344.5 68 550 680 495 3212 866

problem 1-12 number of constraints 3879 best objective 16 004 number of nodes 6561 number of iterations 100 000 CPU time (s) for best objective 1762 relative gap (%)g 19.84

3879 18 763f 85 100 1 371 595 1237 6.02

2711 19 185f 10 100 154 837 265 3.91

709 7498 204 820 1.7 1.4

549 7498 38 321 0.7 0.5

a As reported in ref 14. b Corrected form of Hui and Gupta’s14 model implemented and solved by us. c Average over five values of M. d Best among five values of M. e Objective values in bold indicate suboptimal values. f Best objective among five values of M. g Relative gap (%) with the relaxed objective of 19 966 for problem 1-12.

does not affect the RMIP objective and, thus, has no effect on the integrality gap. In this work, all models resulted in the same RMIP objectives, so we could not use the integrality gap to compare them. Comparison with Previous Models With model M8 as the best model, we proceed to compare its performance in detail with those of the models reported in the literature. We do this only for sets 1 and 3, as the previous works12,14,15 did not consider minimizing tardiness as a scheduling objective. Hui and Gupta14 solved the problems in set 1, whereas Pinto and Grossmann12 and Hui et al.15 solved the problems in set 3. Table 10 shows a comparison of model M8 with the model of Hui and Gupta14 for the problems in set 1. To be fair, we compare model M8 with model M10, which, you recall, is Hui and Gupta’s model implemented and solved by us on our processor. Although all models use the same numbers of variables (binary and continuous) for a given problem, our model M8 uses up to 30% fewer constraints. Model M8 solves problems 1-5, 1-8, and 1-10 in substantially less average (over M) CPU times as compared to model M10. Also, the best (with respect to M) CPU times for model M8 are far lower than the best for model M10. For problem 1-12, model M8 reaches an objective value of 19 185 with a relative gap of 3.9% in around 80% less CPU time (265 s) than model M10, which gives an objective of 18 763 with a relative gap

model

variables conbest binary continuous straints objective

M8 M10a Hui et al.15 Pinto and Grossmann12

189 189 189 161

problem 3-5 136 574 136 709 136 709 167 511

M8 M10

433 433

problem 3-8 223 1334 223 1790

6828.8 6720.7 6716.1 6151

CPU time (s) 0.63 5.06 24b 92.09b

10 986.4 973 10 804.0 4527

a Modified for unit-dependent setup times as done by Hui et al.15 b As reported in the respective works.

of 6% in 1237 s. In other words, model M8 achieved a solution better than that obtained by model M10 in much less time. Note here that the actual objective used in models M8 and M10 amounts to maximizing completion times for zero tardiness. This is the same as minimizing earliness, which is why a higher objective is a better objective. In contrast to models M8 and M10, Hui and Gupta14 could not solve problems 1-8 and 1-10 to optimality. Table 11 shows a comparison of model M8 with the models used in the previous works of Pinto and Grossmann12 and Hui et al.15 Our aim here is to show that our model M8 (modified for unit-dependent setup times as shown earlier) obtains a better optimal objective value than those of previous works.12,15 For problem 3-5, Hui et al.15 obtained an objective of 6716.1, which they showed was better than the objective of 6151 reported by Pinto and Grossmann.12 Using our model, we obtain an optimal objective of 6828.76, which is better than both of the preceding values. Both works12,15 considered the unit-dependent setup time as a part of the processing time, whereas we treat it separately. Our implementation of model M10 (modified for unit-dependent setup times as done by Hui et al.15) gives an objective of 6720.66 vs 6716.10 in their implementation. Interestingly, their heuristic objective is better than their15 reported exact optimum. In our comparison, we could not include problems larger than problem 3-5, because the previous works either did not solve such problems15 or did not report the data.12 However, as shown in Table 11, we do solve problem 3-8 (using the data of problem 3-5) with models M10 and M8 to reconfirm that model M8 gives a better optimal objective than model M10 (10 986.36 vs 10 803.96) in much less time (973 vs 4527 s). Note here that the actual objective used in models M8 and M10 is maximizing completion times. For zero tardiness, this is the same as minimizing earliness, which is why a higher objective is a better objective. In addition to demonstrating the superiority of model M8 over model M10 in terms of the optimal objective value, we now show that model M8 is far more efficient than model M10 for any value of M in problem 3-8. This is shown in Table 9, where the performances of models M8 and M10 are compared for five values of M. As we can see, model M8 has an average relative CPU time of 1.07 as compared to 4.57 for model M10. Moreover, model M10 fails to solve to optimality for three of five values of M (shown by CPU times in bold). Thus, we conclude that model M8 performs better than model M10 not only in terms of the objective value, but also in that it requires 70% less CPU time on an average.

2376


Table 12. Start Times (h) and Unit Assignments of Orders for Problems 1-10 and 2-10 via Model M8 Problem 1-10: Minimizing Earliness Using N ) 10 000 and M ) 700, CPU time ) 866 s stage 1

stage 2

order

unit

start time

1 2 3 4 5 6 7 8 9 10 average

6 1 5 3 2 4 1 1 6 1 -

419.4 312 439.6 433.1 457.5 461.5 442 496 480 527 446.8

unit

start time

8 8 8 9 9 8 7 8 9 8 -

437.5 335 457.7 453.1 474.5 473.5 473 508 496 539 464.7

stage 3

stage 4

unit

start time

14 12 10 12 13 11 10 19 11 11 -

442.5 340 462.7 458.1 479.5 480.5 504 524 507 552 475

stage 5

unit

start time

unit

start time

earliness (h)

21 22 20 20 21 22 22 21 20 21 -

454.5 352 486.7 470.1 491.5 492.5 518 534 514 556 486.9

24 25 24 23 23 25 23 25 24 25 -

464 452 496 478 504 507 529 543 529 568 507

22 0 0 28 12 0 0 0 0 0 6.2

Problem 2-10: Minimizing Tardiness Using M ) 1000, CPU time ) 1.18 s stage 1 order

unit

start time

1 2 3 4 5 6 7 8 9 10 average

2 5 5 1 5 1 5 6 6 2 -

0 0 29 21 53.1 0 71.1 0 24 25.1 22.3

stage 2 unit

start time

9 7 9 8 7 7 9 9 9 7 -

97.1 111.1 47.1 41 70.1 15 147.1 59.1 82.1 92.1 76.2

stage 3

stage 4

unit

start time

12 11 12 11 12 12 11 12 11 11 -

102.1 132.1 122.1 160.1 190.1 150.1 178.1 172.1 197.1 148.1 155.2

Finally, it is worth pointing out that model M10 performs better than the models of Hui and Gupta14 and Hui et al.15 in terms of both CPU times and objectives. The far lower CPU times for model M10 as compared to those of Hui and Gupta14 and Hui et al.15 speak volumes for the rapid advances in both computing power and optimization solvers. At this pace of technological advance, our ability to solve increasingly larger and practically more useful MILP problems should continue expanding tremendously.

stage 5

unit

start time

unit

start time

earliness (h)

21 22 21 20 20 21 21 20 20 21 -

231.4 144.1 210.1 226.1 236 263.9 193.1 186.1 204.1 247.9 214.3

23 25 23 23 23 25 23 24 23 23 -

245.1 244.1 277.1 322.1 348.1 296.1 306.1 196.1 219.1 380.1 283.4

240.9 207.9 218.9 183.9 167.9 210.9 222.9 346.9 308.9 187.9 229.7

Table 13. Start Times (h) and Unit Assignments of Orders for Problem 3-5 via Model M8 stage 1

stage 2

stage 3

stage 4

stage 5

start start start start start order unit time unit time unit time unit time unit time 1 2 3 4 5 avg

4 3 5 2 1 -

431.4 312 391.6 431.1 401.5 393.5

9 449.5 12 454.5 21 8 335 18 340 22 9 409.7 16 414.7 21 8 451.1 15 456.1 20 8 418.5 10 423.5 20 - 412.8 - 418 objective value ) 6828.7 h

466.5 352 438.7 468.1 435.5 432.2

24 25 24 23 23 -

476 452 448 476 448 460

Examples In this section, we first present solutions for problems 1-10 and 2-10 using model M8 and then for problem 3-5 using models M8 and M10. For problems 1-10 and 3-5, we use N ) 10 000 and M ) 700, and for problem 2-10, we use M ) 1000. Recall that problem 1-10 minimizes earliness, whereas problem 2-10 minimizes tardiness. Table 12 shows the unit assignments and start times of orders in all stages for problems 1-10 and 2-10. The schedule for problem 1-10 starts at 312 h and has an average start time of 447 h in stage 1. The plant remains idle before this time, which is not desirable. In contrast, the schedule for problem 2-10 starts at time zero and has an average start time of 22 h in stage 1. In fact, the schedule for problem 2-10 almost finishes processing all of the orders before the schedule for problem 1-10 even starts. As a result, there is no equipment idle time in problem 2-10, which is what one would desire in industrial practice. Moreover, problem 2-10 requires only 1.18 s to solve, whereas problem 1-10 requires 866 s. Also, observe that the average earliness for the schedule of problem 2-10 is 230 h as compared to 6 h for problem 1-10. This clearly indicates that the plant

is underutilized or the due dates are too conservative. It is easy to infer that, as the plant receives more orders to process or as the due dates or customer lead times are revised in view of the plant capacity, the average earliness for the schedule of problem 2-10 will decrease. Hence, the objective of minimizing tardiness is a realistic one for scheduling a fully utilized plant with reasonable due dates. We mentioned earlier that our approach to incorporate unit-dependent setup times gives a better objective value than the approaches used in the previous works of Pinto and Grossmann12 and Hui et al.15 We now solve problem 3-5 using models M8 and M10 to compare their schedules. Remember that problem 3-5 has unit-dependent setup times and an objective of maximizing completion times of orders in all stages. Model M8 is modified for unit-dependent setup times as mentioned earlier. Similarly, we modify model M10 for unit-dependent setup times as done by Hui et al.15 Tables 13 and 14 show the unit assignments and start times of orders in various stages for models M8 and M10, respectively. We see from the tables that, for stage 1, the average start time (393.5 h) for the schedule of model M8 is greater than that (368.2 h) for model M10.


Conclusions

Table 14. Start Times (h) and Unit Assignments of Orders for Problem 3-5 via Model M10a stage 1

stage 2

stage 3

stage 4

stage 5

start start start start start order unit time unit time unit time unit time unit time 1 2 3 4 5 avg

1 1 2 5 3 -

381.9 271.5 370.1 409.6 408 368.2

9 408 14 414 21 8 302.5 10 308.5 22 9 396.2 10 402.2 20 9 437.6 12 443.6 20 8 433 11 439 21 - 395.5 - 401 objective value ) 6720.66 h

428.5 323 428.7 458.1 453.5 418.4

24 25 23 24 23 -

444 447 444 472 472 455.8

a Model M10 is our corrected implementation (Appendix C) of Hui et al.’s15 model for unit-dependent setups.

Figure 2. Gantt chart for problem 3-5 from model M8 using N ) 10 000 and M ) 700. Each block represents processing on some unit, and the number below the block indicates the unit.

An improved continuous-time MILP formulation for the short-term scheduling of a multistage batch plant with nonidentical parallel units in each stage has been presented. In contrast to previous work12,14,15 on this problem, the proposed formulation incorporates realistic features such as release times for orders and units, initial state of the plant, and both sequence-dependent and unit-dependent setup times. It also presents new approaches for assigning consecutive orders to a single unit, from which the best was identified via a thorough numerical evaluation. The best formulation reduces constraints by roughly 30%, nonzeros by roughly 50%, and solution times by 65% as compared to the earlier works. This work also presents a new scheduling objective that eliminates equipment idle times and requires significantly less CPU time to be optimized as compared to the objective function used by previous work12,14,15 on this problem. Our formulation for unitdependent setup times yields schedules with better objective values than those of previous work.12,15 Finally, through this study, we highlight the tremendous impact of the value of M in big-M constraints on solution times and the need to evaluate models for several values of M before choosing the best. Acknowledgment The authors thank Professor David Hui and Mr. Avaneesh Gupta of Hong Kong University of Science and Technology for reconfirming the numbers of constraints in their paper (Hui et al.15) on our request. Notation Indices

Figure 3. Gantt chart for problem 3-5 from model M10 using N ) 10 000 and M ) 700. Each block represents processing on some unit, and the number below the block indicates the unit.

s, s′ ) stages j, j′ ) units i, i′ ) orders Sets

In fact, this is true for all stages. Thus, model M8 finds a schedule with greater completion times in all stages and hence gives a higher objective value (6828.76 h) as compared to model M10 (6720.66 h). A possible explanation for this is that the previous works12,15 considered the unit-dependent setup time as a part of the processing time, whereas we do not. Doing the former prevents an order i from starting to process in stage nsis until the unit in stage s (where order i was processed) is set up to process the next order in sequence in stage s. This is not a valid approach because an order i can start processing in stage nsis as soon as it finishes processing in stage s. Thus, the systems in the previous models require more time to finish processing orders. This is evident from Tables 13 and 14. The time span for the schedule of model M8 is 312-500 h as compared to 272500 h for model M10. Figures 2 and 3 show the Gantt charts for the schedules of problem 3-5 using models M8 and M10, respectively. Each order is represented by a bar with five sections, which indicate the five stages of processing of orders. The numbers below each section of the bar show the unit where the order is processed in that stage. Observe that the two schedules (for models M8 and M10) have different unit assignments for a given order in all stages.

S ) stages I ) orders J ) units Si ) stages in which order i must be processed Is ) orders that are to be processed in stage s Ij ) orders that can be processed on unit j Js ) units in stage s Jis ) units in stage s that can process order i NCis ) orders that cannot follow order i in stage s FIi ) orders that cannot follow order i FPis ) feasible predecessors of order i in stage s FSis ) feasible successors of order i in stage s Parameters fsi ) first processing stage of order i nsis ) next processing stage of order i being processed in stage s lsi ) last processing stage of order i di ) due date of order i tij ) processing time of order i on unit j cii′j ) time required for a changeover from order i to order i′ on unit j c0ij ) time required by unit j to start processing the first order i cii′ ) time required for a changeover from order i to order i′ on any unit

2378


cj ) time required for any changeover of orders on unit j wis ) weight to prioritize processing of order i in stage s pi ) weight to prioritize order i URTj ) release time of unit j ORTi ) release time of order i N, M ) large positive values Variables Zij ) binary variable for assigning order i to unit j ZFij ) binary variable for assigning order i to process first on unit j Xii′s ) binary variable for sequencing order i′ after order i in stage s Tis ) start time for processing order i in stage s Tdi ) tardiness in completion of order i

Appendix A In a noncyclic sequence of orders on unit j in stage s, all orders cannot follow each other. Exactly one order must be first in the sequence, i.e., ZFij ) 1 for some i. Now, ZFij ) 1 w Zij ) 1 from eq 15, and Zij ) 1 w Zij′ ) 0 ∀ j′ ∈ Jis, j′ * j, from eq 1. Thus, if unit j processes only one order, then Zij automatically becomes an integer. Now, if unit j processes multiple orders, then the first order i must have a successor, i.e., Xiks ) 1 for some order k. Then, any of the nine unit assignment constraints in Table 1 for assigning consecutive orders i and i′ to a single unit j will force Zkj ) 1 and Zkj ) 1 w Zkj′ ) 0 ∀ j′ ∈ Jks, j′ * j, from eq 1. In other words, Zkj ∀ j automatically becomes an integer. Continuing the argument further, it becomes clear that, as long as ZFij and Xii′s are taken as binary, all Zij’s will be binary and can be treated as continuous with 0 e Zij e 1. Appendix B The second term on the LHS of eq 6a serves no purpose other than to tighten constraint 6a in the continuous sense. Thus, it is possible to delete this term without affecting the feasibility. Removing this term reduces the number of nonzeros in the formulation, which might reduce the solution time. This idea provides an alternate constraint for eq 6a as

Zij e Zi′j + 1 - Xii′s - Xi′is j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is (6c) and eqs 5 and 6c form the second valid set to force consecutive orders i and i′ onto one common unit. Now, eq 6c is written for combinations only, so swapping i and i′ in eq 6c gives us the new constraint

Zi′j e Zij + 1 - Xii′s - Xi′is j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is (6d) When Xii′s + Xi′is ) 1, eq 6d reduces to Zi′j e Zij, j ∈ Jis ∩ Ji′s, i′ > i, (i, i′) ∈ Is. Similarly to eq 6b, this also requires Zi′j ) 1 for j ∈ Jis ∩ Ji′s (or Zij ) 1 for j ∈ Jis ∩ Ji′s and Zi′j ) 0 ∀ j ∈ Ji′s - Jis) to force ZijZi′j ) 1 for one common unit j. As stated earlier, eq 5 meets both of these conditions, so eqs 5 and 6d form the third valid constraint set for unit assignment. Now, if we use just Xii′s instead of Xii′s + Xi′is in eq 6c, then we must write eq 6c for all permutations of i and i′. This gives us the new constraint

Zij e Zi′j + 1 - Xii′s j ∈ Jis ∩ Ji′s, i′ ∉ FIi , i′ * i, (i, i′) ∈ Is (6e) Note that eq 6e gives the following two constraints for any two pair of orders i and i′


j ∈ Jis ∩ Ji′s

Zi′j e Zij + 1 - Xi′is

j ∈ Jis ∩ Ji′s

We can show that the above constratints force ZijZi′j ) 1 for one common unit j when one of the following two conditions is met:

1. (Zij ) 1 for one j ∈ Jis ∩ Ji′s when Xii′s ) 1) and (Zi′j ) 1 for one j ∈ Jis ∩ Ji′s when Xi′is ) 1) 2. (Zi′j ) 1 for one j ∈ Jis ∩ Ji′s and Zij ) 0 ∀ j ∈ Jis - Ji′s when Xii′s ) 1) and (Zij ) 1 for one j ∈ Jis ∩ Ji′s and Zi′j ) 0 ∀ j ∈ Ji′s - Jis when Xi′is ) 1) The above two conditions have separate scenarios for Xii′s ) 1 and Xi′is ) 1, as we do not know a priori which of Xii′s and Xi′is will be equal to 1. Again, eq 5 meets both conditions, so eqs 5 and 6e form the fourth valid constraint set to assign consecutive orders i and i′ to one common unit j. Similarly, we can rewrite eq 6d with only Xii′s and obtain the new alternative


j ∈ Jis ∩ Ji′s, i′ ∉ FIi , i′ * i, (i, i′) ∈ Is (6f)

As for eq 6e, eq 6f also requires one of the following two conditions to force ZijZi′j ) 1 for one common unit j

1. (Zi′j ) 1 for one j ∈ Jis ∩ Ji′s when Xii′s ) 1) and (Zij ) 1 for one j ∈ Jis ∩ Ji′s when Xi′is ) 1) 2. (Zij ) 1 for one j ∈ Jis ∩ Ji′s and Zi′j ) 0 ∀ j ∈ Ji′s - Jis when Xii′s ) 1) and (Zi′j ) 1 for one j ∈ Jis ∩ Ji′s and Zij ) 0 ∀ j ∈ Jis - Ji′s when Xi′is ) 1) Again, eq 5 meets both conditions, so eqs 5 and 6f form the fifth valid constraint set for assigning consecutive orders i and i′ to one common unit j. Examining the required conditions for eqs 6a and 6cf, we realize that it is not necessary to force both Zij and Zi′j to equal 0 for uncommon units. Even if we force one of the orders (say i) to process on a common unit, the other (i′) can be forced to process on the same unit using any of eqs 6a and 6c-f. Thus, we now present sets where we first force either Zij ) 0 or Zi′j ) 0 for uncommon units and then use any of eqs 6a and 6c-f to force ZijZi′j ) 1 for one common unit j. To this end, we first modify eq 5 to force only Zij ) 0 for uncommon units as follows

Xii′s + Xi′is +

∑

Zij e 1

j∈ (Jis-Ji′s)

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is (7) Note that Xii′s + Xi′is ) 1 w Zij ) 0 for uncommon units


j ∈ Jis - Ji′s. This forces Zij ) 1 for one j ∈ Jis ∩ Ji′s from eq 1. To get Zi′j ) 1 for that j, we use one of eqs 6a and 6c-f. Any of them will force Zi′j ) 1 if one of their respective conditions is met. We can see that eq 7 meets condition 1 of eqs 6a and 6c. As discussed earlier, however, eq 6a is equivalent to eq 6c, if we drop the second term on the LHS of eq 6a, so we will not consider eq 6a henceforth. Then, eqs 7 and 6c form the sixth valid constraint set for assigning consecutive orders i and i′ to one common unit j. Recall that eq 7, similarly to eq 5, is also written for combinations of i and i′ only and not permutations. Now, if we use only Xii′s in eq 7, then we must write it for all permutations of i and i′ to obtain the new constraint

∑

Xii′s +

Zij e 1


(8)

j∈(Jis-Ji′s)

This constraint, for any pair (i, i′), gives

Xii′s +

∑

Zij e 1

∑

Zi′j e 1

j∈(Jis-Ji′s)

Xi′is +

It can be shown that these inequalities satisfy condition 1 of eq 6e, so eqs 8 and 6e form the seventh valid constraint set for assigning consecutive orders i and i′ to one common unit j. Equation 7 forces Zij ) 0 for uncommon units j ∈ Jis - Ji′s. Alternatively, we can force Zi′j ) 0 for uncommon units j ∈ Ji′s - Jis by using

∑

Zi′j e 1

j∈(Ji′s-Jis)

i′ ∉ NCis, i′ > i, (i, i′) ∈ Is (9) Equation 9 meets condition 1 of eq 6d, so eqs 9 and 6d form the eighth valid constraint set. For eq 9 also, if we use only Xii′s, then the constraint must be written for all permutations i and i′ as follows

Xii′s +

∑

Zi′j e 1


(10)

j∈(Ji′s-Jis)

As argued for eq 8, eq 10 meets condition 1 of eq 6f, so eqs 10 and 6f form the last (ninth) valid constraint set. Appendix C In implementing the model of Hui and Gupta14 and Hui et al.,15 we could not reproduce their reported numbers of constraints. More importantly, we could not even obtain a feasible solution from their model. On closer scrutiny, we discovered that eq 1 as published in their papers14,15 had a probable misprint. Their eq 1 (as reported) is, in our notation

Xii′s + Zij +

Xii′s + Zij +

∑ Zi′j′ e 2

j′∈Ji′s j′*j

j ∈ Jis, i′ ∉ NCis, i′ ∉ FIi , i′ * i, (i, i′) ∈ Is (22) We now write this constraint for all units j where i can be processed, i.e., j ∈ Jis instead of j ∈ Jis ∩ Ji′s. This ensures that eq 22 will surely be forced for one j, no matter which Zij (common j or uncommon j) becomes equal to 1. Now, when Xii′s becomes equal to 1, Zi′j ) 1 will be forced when Zij ) 1. Thus, model M10 is Hui and Gupta’s14 model with their eq 1 replaced by our eq 22. Only when we used eq 22 in our implementation were we able to match our counts of constraints with those reported by Hui and Gupta.14 Literature Cited

j∈(Ji′s-Jis)

Xii′s + Xi′is +

only the common units between orders i and i′. Now, if Xii′s ) 1 and Zij ) 1 for some common unit j, then it will definitely force Zi′j ) 1. However, if Zij never becomes equal to 1 for some j ∈ Jis ∩ Ji′s, then this constraint is always relaxed. This would allow i and i′ to process on different units and still give Xii′s ) 1. To rectify this oversight, we modify the constraint as follows

∑ Zi′j′ e 2

j′∈Ji′s j′*j

j ∈ Jis ∩ Ji′s, i′ ∉ FIi , i′ * i, (i, i′) ∈ Is This constraint is supposed to assign two consecutive orders to the same unit. Observe that it is written for

(1) Reklaitis, G. V. Overview of Scheduling and Planning of Batch Process Operations. Presented at the NATO Advanced Study Institute on Batch Process Systems Engineering, Antalya, Turkey, 1992. (2) Zentner, M. G.; Pekny, J. F. Learning To Solve Process Scheduling Problems: The Role Of Rigorous Knowledge Acquisition Frameworks. In Foundations of Computer Aided Process Operations; Rippin, D. W. T., Hale J. C., Davis, J. F., Eds.; CACHE: Austin, TX, 1994. (3) Shah, N. Single and Multisite Planning and Scheduling: Current Status and Future Challenges. In Proceedings of the Third Conference on Foundations of Computer Aided Process Operations (FOCAPO′98); CACHE and CAST Division of AIChE: New York, 1998; pp 75-90. (4) Kondili, E.; Pantelides, C. C.; Sargent, R. W. H. A General Algorithm for Short-term Scheduling of Batch Operations. I. MILP Formulation. Comput. Chem. Eng. 1993, 17, 211-227. (5) Shah, N.; Pantelides, C. C.; Sargent, R. W. H. A General Algorithm for Short-term Scheduling of Batch Operations. II. Computational Issues. Comput. Chem. Eng. 1993, 17, 229-244. (6) Zhang, X.; Sargent, R. W. H. The Optimal Operation of Mixed Production FacilitiessA General Formulation and Some Approaches for the Solution. Comput. Chem. Eng. 1996, 20, 897904. (7) Pantelides, C. C. Unified Framework for Optimal Process Planning and Scheduling. In Proceedings of the Second Conference on Foundations of Computer Aided Process Operations (FOCAPO II); CACHE and CAST Division of AIChE: New York, 1994; pp 253-274. (8) Schilling, G.; Pantelides, C. C. A Simple Continuous-time Process Scheduling Formulation and a Novel Solution Algorithm. Comput. Chem. Eng. 1996, 20, S1221-S1226. (9) Ierapetritou, M. G.; Floudas, C. A. Effective ContinuousTime Formulation for Short-term Scheduling. 1. Multipurpose Batch Processes. Ind. Eng. Chem. Res. 1998, 37, 4341-4359. (10) Cerda, J.; Henning, G. P.; Grossmann, I. E. A MixedInteger Linear Programming Model for Short-Term Scheduling of Single-Stage Multiproduct Batch Plants with Parallel Lines. Ind. Eng. Chem. Res. 1997, 36, 1695-1707. (11) Mendez, C. A.; Henning, G. P.; Cerda, J. Optimal Scheduling of Batch Plants Satisfying Multiple Product Orders with Different Due Dates. Comput. Chem. Eng. 2000, 24, 22232245. (12) Pinto, J. M.; Grossmann, I. E. A Continuous Time Mixed Integer Linear Programming Model for Short-Term Scheduling of Multistage Batch Plants. Ind. Eng. Chem. Res. 1995, 34, 30373051.

2380


(13) Pinto, J. M.; Grossmann, I. E. An Alternate MILP Model for Short-Term Scheduling of Batch Plants with Preordering Constraints. Ind. Eng. Chem. Res. 1996, 35, 338-342. (14) Hui, C. W.; Gupta, A. A novel MILP formulation for shortterm scheduling of multistage multi-product batch plants. Comput. Chem. Eng. 2000, 24, 1611-1617. (15) Hui, C. W.; Gupta, A.; Meulen, H. A. J. A novel MILP formulation for short-term scheduling of multi-stage multi-product batch plants with sequence-dependent constraints. Comput. Chem. Eng. 2000, 24, 2705-2717.

(16) Brooke, A.; Kendrick, D.; Meeraus, A.; Raman, R. GAMS: A User’s Guide; GAMS Development Corporation: Washington, DC, 1998.

Received for review March 7, 2002 Revised manuscript received January 28, 2003 Accepted February 7, 2003 IE020180G

An Improved MILP Formulation for Scheduling Multiproduct

Recommend Documents