The material allocation problem in the steel industry

Steelworks usually produce considerable surplus inventory, as when orders are canceled after materials have been produced or when some finished products are below the quality levels required for an order. To make use of surplus inventory, daily production planning is typically divided into two main steps. In the first step, the plan tries to satisfy orders as far as possible with existing materials (such as surplus inventory). In the second step, the plan designs production units for the manufacture of the remaining orders.

The material allocation problem (MAP) is a problem in the first step of the order-fulfillment process. In the MAP, the existing materials may include work-in-progress goods and surplus inventory. The work-in-progress goods are tentatively allocated to some orders, but they can be reallocated to new orders while solving the MAP.

The MAP is a kind of matching problem on a sparse bipartite graph, allocating a given set of materials to a given set of orders. Each order restricts its allocatable materials through the quality and other attributes of the order. In Figure 1(a), the dashed line from the material to an order indicates that the material can be allocated to that order. The material can be cut into any size and each order can be satisfied with multiple sources of material. For example, order A and order B respectively require five tons and four tons, and material batch V has a weight of ten tons, which can be allocated to both orders. When we allocate the material to these orders, the material will be split into three parts: five tons for order A, four tons for order B, and one ton remaining as surplus. [The trimming and yield losses shown in Figure 1(b) are discussed below in the section on matches]. In the following argument, we refer to surplus as the remainder of the material after allocations.

Figure 1

The MAP has several hard constraints that make it difficult to solve. One of the constraints is a unit weight constraint: When we allocate a batch of materials to an order, the materials must be cut into pieces whose sizes are within the range specified by the order. For example, suppose that an order requires 12 tons and that maximum and minimum unit weights specified in the order are five tons and four tons, respectively. (Customers place constraints on orders for reasons such as transport limitations and manufacturing conditions.) When we allocate a 12-ton batch of material to the order, the only solution that meets the constraint conditions is to cut it into three parts of four tons each. Other MAP constraints are described in the next section.

Each solution is evaluated from a broad perspective that includes the due dates of orders, preferences for matches, allocated weights, and the surplus weights of materials. It is evaluated primarily with respect to the total allocated weight, which also leads to reducing the total surplus. Because of the unit weight constraint, many small batches of surplus material remain after these allocations. Though a large surplus offers possibilities for allocation to orders in the near future, a small surplus (typically less than five tons) is hard to allocate because of the unit weight constraint. Thus, a small surplus is likely to be discarded, which is a pure loss for the steel company, so the minimization of the number of small surpluses is an important evaluation item.

It is easy to show that the MAP is NP-hard because the MAP includes a generalized assignment problem [1] as a subproblem. A typical generalized assignment problem is a problem assigning a set of jobs to a set of machines. The jobs and the machines correspond respectively to the orders and the materials in the MAP.

Because the MAP is an NP-hard problem and the sizes of the instances are very large, we use a local search to obtain a good solution. Our algorithm is based on a variable neighborhood search (VNS) [2] strategy. VNS combines local search with systematic changes of neighborhood and escape phases from local optima. The search first explores a small neighborhood around the current solution and allows an exchange with the current solution if a better one is found. This process continues until the search fails to find a better solution. When this happens, the search explores neighborhoods more distant from the current solution and continues the process. For the local search, we use nine types of moves, including ejection chain search [3]. In the local search algorithm, we adopt some acceleration techniques, such as the don't-look bit [4], and incorporate simple heuristics for surplus reduction into the local search algorithm.

We found that there are still a number of small surpluses produced by this basic algorithm, but that some of them can be removed by trying extensive searches. Therefore, we use a two-variable integer programming (TVIP) neighborhood search technique to reduce the number of small surpluses. TVIP defines a neighborhood for surplus reduction, and we can show that searching the neighborhoods is essentially equivalent to solving TVIP problems. In addition, we use an efficient algorithm to solve the TVIP problem that combines binary search and an algorithm for lattice-point counting in right triangles. For the lattice-point counting, an algorithm exists [5], but we have developed a simpler one.

The use of our algorithms by a Japanese steel company yielded considerable cost reductions.

The rest of this paper is organized as follows: We describe the details of the MAP, show our basic algorithm for the MAP, and propose a TVIP technique for surplus reduction. We then show our experimental results.

The MAP is similar to the generalized assignment problem [1] and the minimum cost flow problem, but these models cannot deal with the unit weight and other constraints.

Kalagnanam et al. [6, 7] modeled the MAP (also known as the surplus inventory matching problem) as a multiple knapsack problem with color constraints. Though the problems are similar, our models differ on many points: For example, the objective in the Kalagnanam model is to maximize the total allocated weights, but the objective in our model is to maximize the total profit. However, the crucial difference is that their model cannot deal with the minimization of the number of small surpluses. In [7], Kalagnanam et al. propose a heuristics based on the maximum flow algorithm for surplus reduction, but it is designed to reduce the total surplus, not to reduce the number of small surpluses.

In this section, we describe the formulation of the MAP. In the MAP, the set O of orders, the set S of (steel) materials, and the set M of the allocatable matches (pairs of orders and materials) are given. Formally, the MAP is the following optimization problem:

The objective is

subject to constraints (1) through (8):

(1)

(2)

nijPUimin ≤ awij ≤ nijPUimax     [∀(i, j) ∈ M],

(3)

where

where

and

(8)

where

awij ≥ 0, awij ∈ R [∀(i, j) ∈ M].

The packing constraint (6) is defined in the materials section below.

The variables are zo_i, zs_j, zm_ij, aw_ij, n_ij, and cw_j, and the other terms are labeled constants associated with orders, materials, and matches. The variables have the following meanings:

In the following subsections, we explain constraints (1)–(8).

Constraints (1)–(3) apply to orders. Each i ∈ O is associated with the following properties: maximum unit weight PU_i^max(>0), minimum unit weight PU_i^min(>0), target weight TW_i (>0), maximum total allocatable weight TW_i^max(>0), and profit per weight PF_i (>0), where PU_i^min ≤ PU_i^max and TW_i ≤ TW_i^max. The weight TW_i is the required weight for the order i, and TW_i^max is the maximum allocatable weight (so the total allocated weight for the order i cannot exceed TW_i^max).) This allowance with respect to the target weight is defined to avoid leaving many small surpluses after allocation. The profit per weight PF_i is a parameter for the objective function. An order i that is urgent is associated with a large PF_i.

Equation (1) defines the objective value zo_i for each order i. It is defined by the product of the profit per weight PF_i, the minimum value of the target weight TW_i, and the total allocated weight for the order i,

.

Since the allowance is only for satisfying the unit weight constraint, it is preferable to encourage allocating materials close to the target weight.

Constraint (2) ensures that the total allocated weight

does not exceed the maximum allocatable weight TW_i^max.

The unit weight constraint (3) restricts the allocated materials so that each batch of material that is allocated to the order i can be cut into pieces in the size range(PU_i^min, PU_i^max).

Constraints (4)–(6) apply to the materials. Each j ∈ S is associated with the following properties: weight SW_j (>0), profit per weight PF_j (>0), and cost per weight C_j (>0). The SW_j is the weight of the piece of material j ∈ S. The profit per weight PF_j and the cost per weight C_j are parameters for the objective value. A piece of material with a high profit per weight PF_j is favored for allocation. A piece of material with a high cost C_j should not be discarded.

Equation (4) defines the objective value zs_j. The cw_j is the consumed weight of the material j, which can be assumed to be an approximation of the sum of the allocated weights, so

[The precise definition of cw_j is given in Equation (8)]. Hence, the argument SW_j − cw_j for the function f represents the surplus weight of the material. The function f(x) is defined so that small surplus x (typically x < 5) leads to a large penalty. In our parameter settings, we set c₁ = 100, c₂ = 0.3, and c₃ = 0.05. Figure 2 shows the graph.

Figure 2

Constraint (5) ensures that consumed weight does not exceed the weight of the available material.

Constraint (6) is a packing constraint, which restricts the combinations of allocations of a material. Each piece of material is associated with a symmetric function g : M × M → {TRUE, FALSE}, as shown in Table 1. When we allocate a piece of material j to orders by using two matches m₁, m₂ ∈ M, m₁ and m₂ must satisfy g(m₁, m₂) = TRUE. When we allocate some material into three or more matches, each pair of matches must satisfy the function g. The packing constraint is due to the layout of the factory. Orders are associated with a unique route to fulfill the order according to its specifications. If orders with different routes are packed in the same batch of material, it may be impossible to cut the materials for the orders.

Equations (7) and (8) apply to matches. Each (i, j) ∈ M (⊆ O × S) is associated with the following properties: cost C_ij (≥0), trimming loss TL_ij (0 < TL_ij ≤ 1), and process yield PY_ij (0 < PY_ij ≤ 1). The cost C_ij is a parameter for the objective function. A match with a high C_ij is favored for allocation. TL_i,j and PY_i,j are for the yield loss calculation.

Equation (7) defines the objective value for the piece of material j. The objective value for each piece of material j is defined by the product of the cost C_ij and the allocated weight of the match aw_ij.

Equation (8) defines the consumed weight cw_j of a material j. When we allocate material to orders, two types of losses arise: Loss for process yield and trimming loss [Figure 1(b)]. The loss for process yield always arises when we allocate material to orders. (If a piece of material is not allocated to any orders, there is no loss for process yield.) The term SW_jPY_ij in Equation (8) represents the loss for the process yield. The process yield PY_ij is set so that PY_i1,j = PY_i2,j holds if g[(i₁, j), (i₂, j)] = TRUE for (i₁, j), (i₂, j) ∈ M. The trimming loss arises when the width of the material is larger than the required width of the order. The term aw_ij /TL_ij in Equation (8) represents the trimming loss for allocating match (i, j) with the weight of aw_ij.

In the MAP, the objective is to maximize the value of the objective function. Because both the profit PF_i of each order and the profit PF_j of each piece of material are relatively larger than the cost C_ij of each match, we obtain the highest score when the total allocated weight to an order i is equal to the target weight TW_i.

Because the formulation of the MAP contains nonlinear functions, we cannot use mixed-integer programming (MIP) solvers, such as ILOG CPLEX**. Therefore, we have designed a heuristic algorithm for this problem.

In this section we give an overview of our algorithm. In the first subsection, we describe our local search algorithm, in the second we show our acceleration techniques for the local search, and in the third we show our basic heuristics to reduce the surpluses.

As explained above in the Introduction, our algorithm is based on a VNS strategy. To construct an initial solution for the local search, we use a random fit strategy (a simple randomized variant of first fit). The random fit strategy permutes the matches in a random order, and our algorithm allocates the matches as far as possible one by one in that order.

For VNS, we constructed the following nine types of moves, which are divided into five groups. All of the moves are designed so that the objective score improves primarily by increasing the total weight of allocated materials.

Basic moves
In this group, we have two types of moves, a simple move and a flip move. In the simple move, shown in Figure 3(a), we try to allocate a match when there is an allocatable match that is not allocated in the current solution. In the flip move, we remove a match and replace it with another match. Figure 3(b) illustrates a flip: We remove the match (A, X) and replace it with the match (A, Y). The allocated weight of the match is set to 0 in the removing step, and the allocated weight in the replacing step is as close as possible to the target weight of the order. The changes of the allocated weights are done similarly in the following moves.

Figure 3

Cyclic moves
In this group, we have two types of moves, a two-cyclic move and a three-cyclic move. In both moves, we exchange matches cyclically—i.e., some allocated materials are unallocated and replaced with other allocations. For example, a three-cyclic move is shown in Figure 3(c): We remove the three matches (A, Z), (B, X), and (C, Y) and replace them with the three matches (A, X), (B, Y), and (C, Z). The two-cyclic move is defined similarly but can also be regarded (less formally) as a swap.

Shift moves
In this group, we have two types of moves, a one-shift move and a two-shift move. In both moves we remove a match and replace it with another match, and this move triggers additional allocations. For example, a two-shift move is shown in Figure 3(d): We remove the two matches (B, X) and (C, Y) and replace them with the three matches (A, X), (B, Y), and (C, Z). The one-shift move is defined in a similar way: We remove a match and replace it with two matches.

Combination moves
In this group, we have two types of moves, an order-optimal (order-opt) move and a material-optimal (material-opt) move. In the order-opt (material-opt) move, we remove all allocated matches with respect to a single order (material) and reallocate with the best allocations with respect to the order (material). Strictly speaking, we cannot find the best allocations among all of the candidates because there are too many candidates due to the complex objective function, even if the order or material has a few allocatable matches. Therefore, the reallocations are done using suboptimal allocations. For the order-opt moves, we simply assume that only one piece of material can have a surplus in the best allocation and find suboptimal allocations under that assumption. Suppose that the order has k allocatable matches. Then, by that assumption, k − 1 matches are assumed to be fully allocated or not allocated, and one batch of material can be partially allocated. Thus, there are at most k(2^k−1) candidates. Our algorithm finds the best allocation among them and replaces the old allocations with the best one. Similarly, for material-opt moves we assume that the orders, except for one, are allocated by their target weights, and the reallocations are done by using the suboptimal allocations.

Ejection chain search
The ejection chain search is the ninth move. This move is an extension of the shift moves; i.e., the ejection chain search performs a k-shift move for any integer k (see Figure 4). This type of move is used, for example, by Glover [3]. Because the ejection chain search is time-consuming, we first use one-shift and two-shift moves and then use the ejection chain search for k ≥ 3.

Figure 4

To find an ejection chain, we use the concept of an improvement graph (Figure 4). An improvement graph is constructed from a current solution. Each node in the improvement graph corresponds to an allocated match in the solution and is associated with a profit (the change in the objective value) when we remove the corresponding match. (Note that the profit is usually a negative value.) For example, node A shows that we reduce the objective value by 1 if we unallocate match (A, V). Each directed edge in the improvement graph corresponds to a match that is not allocated in the current solution. To illustrate, suppose that there are two allocated matches (A, V) and (B, W) and that a match (B, V) can be allocated if we remove the two matches (A, V) and (B, W). Then a directed edge corresponding to the match (B, V) is added to the improvement graph [the direction is from (A, V) to (B, W)]. The directed edge is associated with the profit when we allocate the match (B, V) in a situation in which the matches (A, V) and (B, W) are unallocated. For example, Figure 4 shows that we gain 2 in the objective score if we unallocate matches (A, V) and (B, W) and allocate match (B, V).

We search for a directed path or a cycle in the improvement graph that satisfies the following condition: Any two matches that correspond to two directed edges in the path or the cycle do not have a shared order or piece of material. If the path or the cycle satisfies this condition, the total profit associated with the path or the cycle is equal to the profit when we remove the matches that correspond to the nodes in the path or the cycle and replace the matches that correspond to the directed edges in the path or the cycle. For example, from the improvement graph in Figure 4, it can be seen that the allocations on the left are improved into the allocations on the right by using the path with solid arrows.

Therefore, all we have to do is to find a path or a cycle that satisfies the above condition in the improvement graph. In our implementation, this is done by a brute-force strategy. Though brute force is normally time-consuming, the algorithm finishes in a reasonable time when this search is applied to a current solution that all of the other moves cannot improve.

Variable neighborhood search
Overall, our algorithm applies the nine types of moves in the following order:

1. Simple move.
2. Flip move.
3. One-shift move.
4. Two-cyclic move.
5. Order-opt move.
6. Material-opt move.
7. Two-shift move.
8. Three-cyclic move.
9. Ejection chain search.

If a move successfully improves a current solution, our algorithm goes back to the first move (the simple move) and continues. If a move fails to improve a current solution, it tries the next move. Our algorithm finishes when all moves fail to improve the current solution. This order is approximately an increasing order of the neighborhood space of each move. Though this order may not be the best one, our experiments show that this order is the best among several candidates with respect to average execution time.

We now describe three techniques used to accelerate our local search algorithm.

Dealing with the unit weight constraints
Suppose there is an order whose minimum unit weight and maximum unit weight are a and b, respectively. When we check whether or not a given weight x satisfies the unit weight constraint of the order, the usual method is to check whether an ≤ x ≤ bn holds by repeatedly incrementing the integer n by 1. Since it is easy to prove that an ≤ x ≤ bn holds if and only if x ≤ ⌊x/b⌋ b holds,¹ we check the unit weight constraint with this inequality, which makes our algorithm efficient. If the weight x does not satisfy the unit weight constraint, it is also easy to show that ⌈x/a⌉ a is the minimum weight not less than x that satisfies the unit weight constraints, and that ⌊x/a⌋ b is the maximum weight not more than x that satisfies the unit weight constraint. These facts are also used in our algorithm.

Don't-look bit
We used the don't-look bit technique [4] to accelerate the local searches. This technique skips duplicates when trying to apply moves. In our implementation, we used this technique only to accelerate the local searches, though there are some implementations in which they accelerate local searches while sacrificing the quality of the solutions obtained in the local searches.

Grouping
Because of the packing constraint, some combinations of matches cannot be allocated at the same time. Therefore, to reduce the computation time for applying moves, we divide the matches associated with each batch of material into groups such that no pair of matches in different groups can be allocated at the same time. Since the function g for each batch of material can be regarded as an adjacency matrix for an undirected graph [i.e., if g(m₁, m₂) = TRUE, we can assume that the nodes m₁ and m₂ are connected by an edge], the grouping is done by identifying the connected components in the graph.

It is important to reduce the number of surpluses, especially the number of small surpluses. To avoid leaving a surplus, our basic algorithm uses a simple heuristic that tries to satisfy orders not only by their target weights but also by using the weights without surplus material. For example, suppose that a material with a weight of 5 tons is allocatable to an order that requires 4.5 tons. We normally try to allocate the material by the weight of 4.5 tons, but we also try to allocate the material with the weight of 5 tons if the order has enough allowance (i.e., the maximum allocatable weight for the order is at least 5 tons). Because a small surplus leads to a bad objective value in the objective function, the allocation of 5 tons is likely to be selected. This simple heuristic significantly reduces surplus materials.

The heuristic for surplus reduction in the previous section is not fully satisfactory. There still remain a certain number of surplus pieces of material that could be further reduced. Figure 5 illustrates an example in which the order A requires 14.3 tons and the order B requires 4.2 tons. Suppose that the minimum and maximum unit weights of the order A are 1.7 tons and 1.9 tons, and that the minimum and maximum unit weights of the order B are 2.0 tons and 2.3 tons. For this instance, the heuristic in the previous section fails to allocate the remaining surplus, while we can eliminate that surplus by allocating the material for order A with 5.6 tons and for order B with 4.4 tons. Therefore, we propose a move to deal with this problem which assuredly determines whether or not we can eliminate a surplus.

Figure 5

We designed a move called the TVIP neighborhood search. In the TVIP, we release all allocations with respect to a batch of material with a surplus and reallocate the material without the surplus. The reallocation is done by trying all of the combinations of allocatable matches, but it is not trivial, even when the number of candidate matches is small. This move is inserted after the three-exchange move in the VNS.

Usually only a few batches of materials are split into more than two orders, so we focus on allocating the batch of material to at most two orders in the TVIP. It is trivial to check whether the material can be allocated to a single order without any surplus. Therefore, we focus on determining whether or not a batch of material can be allocated to two orders.

Suppose that we are given two orders A and B as well as a piece of material that has a weight of w tons and is allocatable to both orders. Suppose that the maximum unit weights of the orders A and B are u_au and u_bu, the minimum unit weights of the orders A and B are u_al and u_bl, and the maximum allocatable weights of the orders A and B are a_u and b_u, respectively. Let a_l = 0 and b_l = 0. (To simplify the following arguments, we assume that the losses for process yield and trimming do not arise when we allocate matches. We can show that TVIP works, even if we remove this assumption.)

Determining whether the material can be allocated to the orders without a surplus is equivalent to determining whether there is a feasible solution subject to the following conditions:

where

x_a, x_b ∈ R⁺

and

n_a, n_b ∈ Z⁺.

We do not specify objective functions in the conditions (9) because we are concentrating only on reducing the number of surpluses in this section, and every feasible solution has no surplus.

In the conditions in (9), x_a, x_b, n_a, and n_b are variables that respectively represent a unit weight for order A, a unit weight for order B, the number of units for order A, and the number of units for order B, respectively. Hence, x_an_a and x_bn_b represent the weights of allocations for order A and order B.

Finding allocations without any surplus does not always lead to a better solution with respect to the objective function. However, because of the importance of reducing the number of small surpluses, finding allocations with no surplus usually leads to a better solution with respect to the objective function.

Though finding a feasible solution of the conditions in (9) is not trivial, there is an efficient algorithm for finding a feasible solution. Suppose that the conditions have a feasible solution x_a, x_b, n_a, n_b. Then there exists another solution x’_a, x’_b, n’_a, n’_b such that at least one of the eight equalities in (9a–d) holds (by perturbing x_a and x_b). When any of the four equalities in (9a) or (9b) holds, it is easy to find a feasible solution for all of the conditions in (9). For example, when we assume that there is a feasible solution such that x’_an’_a = a_l, we only have to check whether or not the order B can be allocated with the weight w − a_l, which is a trivial check. The other three cases can be tested in a similar way. When any of the four equalities in (9c) or (9d) holds, we can also find a feasible solution efficiently. To prove this, first we show that there is an equivalent TVIP formulation for each case. If x’_a = u_al, the conditions in (9) are equivalent to the following TVIP problem, subject to

(10)

For the other three cases, we find similar TVIP problems. Thus, all we have to do is solve a TVIP problem (10) efficiently.

The feasible region for the integer programming problem of (10) is depicted by the region ABCD in Figure 6, where

n_al = max {a_l /u_al, (w − b_u)/u_al}

and

n_au = min {a_u /u_al, (w − b_l)/u_al}.

Figure 6

Finding a feasible solution for the integer programming of (10) is equivalent to finding a lattice point (a point with integer coordinates) in that region.

It is easy to find a lattice point on line segments by using the following lemma (see Theorems 2–4 on p. 24 of [8]).

Lemma 1
The linear Diophantine equation ax + by = c has a solution if and only if c is divisible by d, where d = GCD(a, b) (the greatest common divisor of a and b). Furthermore, if (x₀, y₀) is a solution of this equation, the set of solutions of the equation consists of all integer pairs (x, y), where

x = x₀ + (b/d)k, y = y₀ − (a/d)k     (k = ···, −2, −1, 0, 1, 2, ···).

By using the extended Euclidean algorithm [9], GCD(a, b) as well as the integers x₀ and y₀ such that ax₀ + by₀ = GCD(a, b) can be computed. Thus, by using Lemma 1, it is easy to check whether or not there exists a lattice point on a given segment.

To find a lattice point inside the above region, we use a binary search (using vertical lines to divide the region). In order to do the binary search, we need an algorithm to determine whether or not a lattice point is in a given region. It is sufficient if there is an algorithm for counting the number of lattice points in a given region. The following is an efficient algorithm for counting lattice points in right triangles:

Algorithm calcN(a, b, c)
  a, b, and c are positive integers such that a ≥ b
begin
  m := ⌊c/a⌋;
  if a = b then
    return m(m − 1)/2;
  else
    k := ⌊(a − 1)/b⌋; h = ⌊(c − am)/b⌋;
    t := km(m − 1)/2 + mh
    return calcN[b, a − bk, c − b(km + h)] + t;
  endif
end

This algorithm counts the number of lattice points in triangle

T(a, b, c) = {(x, y) ∈ R²|ax + by ≤ c, x > 0, y > 0}.

Since the region can easily be divided into right triangles and rectangles, we can count the number of lattice points in the region efficiently. Thus, we can find a lattice point in the given region by using a binary search.

Lattice-point-counting algorithms have been studied intensively for high-dimensional cases. For example, Barvinok gives an efficient algorithm [10], but it is theoretical and too complicated. (The algorithm is not explicitly described. To the best of our knowledge, there exists only one implementation of the algorithm, named LattE [11] because of its complication). For a two-dimensional case, Beck and Robins give an algorithm [5] that is much simpler, but still more complex than ours (the time complexity of our algorithm is the same as that of the algorithm of Beck and Robins). The correctness of our algorithm is given in [12]. Note that our algorithm deals only with integers, which contributes to its simplicity. We also note that Eisenbrand and Rote give an algorithm for solving TVIP [13], but their algorithm is also theoretical and not practical. (Their algorithm is fast when integers require long bits to encode them; i.e., integers are encoded far longer than 32 bits or 64 bits.)

It might seem that lattice points could be found by using simpler heuristics. A heuristic may suffice if there is a lattice point in the region, but in most cases the region does not contain any lattice point. This is because the solutions we consider here were supplied by the simple surplus-reduction heuristics in the previous section. Our algorithm is especially efficient if the region contains no lattice points, because in those cases we avoid the binary search.

In this section, we show experimental results for our algorithm. All tests were executed on a PC workstation, an Intel Pentium** 4 running at 3.0 GHz, with 1 GB RAM using Microsoft Windows Server** 2003. The maximum time for execution was set to 3,600 seconds (one hour). This limitation is not too long for our client because the program runs at night.

To test our algorithm, we used six instances, sized as shown in Table 2. The differences among the six instances are the numbers of orders and materials. All of the instances were generated randomly in the following manner.

First, all matches were randomly generated. Every order was randomly divided into three groups (I, II, and III), and the packing constraints in the materials were determined so that orders in the same group could be allocated at the same time and orders in different groups could not be allocated at the same time.

The target weights of each order were randomly chosen from the range of 2.0–12.0 tons. The maximum total allocatable weights for each order were set to 1.2 times the target weight of the order. The minimum unit weight for each order was randomly chosen from 0.5 tons to 0.9 times the target weight of the order. The maximum unit weight for each order was randomly chosen from 0.5 tons to 1.5 times the minimum unit weight of the order. The profit per weight for each order was randomly chosen from the range [0–500].

The weight of each batch of material was randomly chosen from the range of 12.0–18.0 tons. The profit per weight for each material was randomly chosen from the range 0–500. The cost for each material was set to 1.

The cost of each match was set to 0. The trimming loss for each match was randomly chosen from the range [0.90–1.00]. The process yield of each match was set according to the group of the order of the match. The process yield was set respectively to 1.00, 0.98, and 0.96 for the group of orders I, II, and III.

We first tested the performance of our local search; Table 3 shows the results. As the initial solutions were made by the random fit algorithm explained above in the section on local search, we executed our algorithm five times for each instance. The results in Table 3 are averaged over the five executions.

The table shows that our algorithm runs within the specified times and improves the objective values considerably from the initial solutions. The gaps show that objective values differ within only 0.3% after local searches. These results show that our local search algorithm works well and produces near-optimal solutions. The results of Kalagnanam et al. are for reference only. Because their algorithm was not designed for our model, a direct comparison cannot be made with this data, and we cannot conclude that our algorithm is superior to theirs from these experiments alone.

We used a different experiment to show the performance of our TVIP technique for surplus reduction. We executed our algorithm with and without TVIP five times each. Table 4 shows the results. For each instance, it contains the average number of small surpluses without TVIP, the average number of small surpluses with TVIP, the average reduction ratio of surpluses, and the average ratio of increase of execution times. The numbers in parentheses are the average execution times.

Table 4 shows that TVIP succeeds in reducing the number of small surpluses without significantly increasing the computation time, which shows the effectiveness of applying TVIP in addition to the heuristics for surplus reduction in the basic algorithm. While the reduction is roughly 0.5–3.5 percent on random data in the experiments, TVIP reduces the number of small surpluses 5–10 percent on real data, which yields considerable profit to a steelworks.

We conducted another experiment to compare TVIP and the maximum-flow-based heuristics of Kalagnanam et al. [7], in which we replaced TVIP in our algorithm with their heuristics. However, the results showed that their heuristics did not improve either solution; that is, it produced the same solutions as our algorithm without TVIP. This shows that their heuristics may be effective for reducing the total number of surpluses but are not as effective for reducing the number of small surpluses.

We constructed an algorithm based on variable neighborhood search for the MAP, which has resulted in considerable cost savings in a real steelworks. In particular, our TVIP technique has contributed to reducing the number of small surpluses of material without significantly increasing the computation times.

TVIP can be applied to other problems in real applications that cannot be modeled as well-known simplified models. Such problems can essentially include TVIP problems as subproblems. Even if the integer programming is associated with an objective function, the binary search in TVIP can find the best solution.

Our future work is to generalize our TVIP technique in order to adapt it to three or more variable integer programming problems. This generalization will depend on a practical and efficient algorithm for solving fixed-variable integer programming problems or a practical and efficient algorithm for lattice-point counting in fixed-dimensional polytopes (higher-dimensional analogs of polygons and polyhedrons).

We thank the IBM workers who made a dedicated effort to deliver our algorithm to the client. We also thank our colleagues, especially Takayuki Osogami, for their useful comments on this paper.

**Trademark, service mark, or registered trademark of ILOG Inc., Intel Corporation, or Microsoft Corporation in the United States, other countries, or both.

¹In this paper, ⌊x⌋ denotes the greatest integer not more than x and ⌈x⌉ denotes the least integer not less than x.

Received September 19, 2006; accepted for publication November 28, 2006; Published online May 18, 2007.