How to Solve an Asset Storage Problem with Mathematical Programming

Moving out is hard and, to be honest, very annoying – it kind of sucks. It’s especially challenging if you don’t have a new place to go. In such cases, you’ll need to rent temporary storage place for your belongings while you stay with friends, which can be pricy since storage places charge by the area used (And before you tell me that these kind of situations never happen, allow me to point out that this exact situation happened to one of my friends when she moved from France to Poland).

If moving out is already difficult for individuals, imagine the complexity of moving out and temporarily storing an entire factory. You might think, "Wow, Luis, this is crazy," but such situations happen in real life. It happened to me, and unfortunately, I didn’t have the analytical tools to handle it at the time.

I used to work in the oil and gas service industry. During my last years there, business conditions in the country deteriorated rapidly and became volatile. For some services, the situation became untenable, leading management to cut losses and shut them down. One of the significant costs these services incurred was the rent for their operational facilities, so it was logical to stop using them to save on rent. However, shutting down a facility meant closing entire workshops and moving all the assets to a storage location. Remember, this is oil and gas – these assets are highly sophisticated tools worth millions of dollars. Just some of the spare parts of these assets were worth thousands and were easy targets for theft and black market sales. So, the challenge was not only to move the assets but also to store them temporarily and, most importantly, safeguard them.

My manager, a brilliant guy, approached me with this problem. He told me we were going to bring all the assets of an entire service division into the main facilities and store them in a secure, restricted area. Since these assets were so valuable, we wanted to avoid any possibility of theft or damage. For instance, some specialized trucks had parts that cost thousands of dollars, and given the country’s conditions, the temptation for theft was high. Thus, no one should enter the safe area unless authorized by management. To enforce this rule, we planned to enclose the area with an electric fence ⚡ and infrared sensors connected to an alarm 🚨 . The sensors would trigger the alarm if someone breached the perimeter (I am not kidding; this is a true story). The issue with this plan was that all that protection is expensive, and the company that rents it charges based on the total enclosed area. My manager wanted me to figure out the spatial arrangement of these assets to minimize the safe zone area, and thus the cost. The electric fence company required the enclosed area to be rectangular or square. Additionally, we couldn’t stack the assets to avoid damage (you can’t put a truck on top of another truck), so we were dealing with a two-dimensional storage problem, known formally in the literature as the assortment problem.

I wish I could tell you that I solved this problem using my mathematical programming superpowers, but that didn’t happen. This was long before I had my PhD, and at that time, I had no clue how to approach this problem with mathematical Optimization – in fact, I wouldn’t have even known how to start googling it (like how on earth I was supposed to known that the problem was called the assortment problem?). Sometimes in life, you just don’t know what you don’t know.

"There are known knowns. These are things we know that we know. There are known unknowns. That is to say, there are things that we know we don’t know. But there are also unknown unknowns. There are things we don’t know we don’t know." – Donald Rumsfeld

In fact learning how to solve some of these problems was my main motivation to continue my studies. So in this article I am going to teach you how to formulate this problem, how to solve it using Python and gurobipy and at the end I will briefly talk about some heuristic procedures that can be applied to tackle it.

If you have never used Gurobi before, especially from Google Colab, I invite you to read one of my previous articles where I explain how to set up your Gurobi web license in your Colab environment. You can find the link to the article below:

Optimizing Project Schedules: Using Gurobipy and Cheche_pm warm-starts for Efficient RCPSP…

If you also need a tutorial and more information about linear programming (LP) or integer programming (IP) techniques I strongly recommend you to check these excellent articles written by Bruno, I am leaving the links down below. Another great source is William’s HP book about Mathematical Programming [1].

Linear programming: Theory and applications

A Comprehensive Guide to Modeling Techniques in Mixed-Integer Linear Programming

Article index:

• The assortment problem
• Installing the libraries and setting up the Colab environment
• Input data & pre-processing
• Solving the problem with Gurobipy
• Constructive Heuristics solutions
• Solving the problem with a warm start from a heuristic
• Conclusions
• References

The assortment problem

A general version of the assortment problem involves selecting a subset of items from a larger set to maximize a certain objective, often revenue or profit, under constraints such as shelf space or budget. It’s a common problem in retail and operations management, involving optimization techniques and often consumer choice modeling. The specific assortment problem we are dealing with in this article is also known as the 2D rectangle packing problem, which frequently appears in logistics and production contexts [2–4].

In this problem, the goal is to determine the minimum size of the object such that the individual smaller units can be cut out from it. Formally, it can be defined as the problem of placing a given set of rectangles of the same or different sizes, without overlap, within a rectangle of minimum area [5]. This problem has multiple applications in facility layout situations, and there is an entire class of heuristic approaches for tackling it [6–8].

The problem requires the following variables:

• Variables x_i and y_i correspond to the (x,y) coordinates of each one of the assets of our problem that are represented by "Rectangles". We have a set I with all the assets or rectangles.
• Variables b_i_j_k are auxiliary binary variables that we are going to need for stablishing some "OR" condition on the model, their usage will be explained later with more detail.
• Variables X and Y correspond to the Total Width and Total Height, resulting from the arrangement of all the assets, with these two variables we will calculate the total area used to store all the assets.
• finally each rectangle i is going to have a variable r_i, that is going to be used to determine if the asset i is going to be rotated or not.

The problem can be modeled using the formulation below:

The formulation above might look daunting at first sight, but it is actually quite simple. The objective function corresponds to the total area, which is the product of X and Y, the total width and height, respectively.

Constraint (1) ensures that the total width value X is greater than the x-coordinate plus the width of each item. Similarly, constraint (2) ensures that the total height value Y is greater than the y-coordinate plus the height of each item.

Constraints (3–7) ensure there is no overlap between the positions of the different items. These constraints work as an OR condition, using the variables b_i_j_1, b_i_j_2, b_i_j_3, and b_i_j_4 for each item, along with the Big-M method. In summary, an item can be to the left, to the right, on top, or below another item, but never overlapping.

Constraints (8–11) ensure that the variables belong to the appropriate domains.

Installing the libraries and setting up the Colab environment

To use Gurobi in Google Colab, we first need to install it using the following code:

!pip install gurobipy

Once installed we can proceed to import the required libraries for this project.

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle
from itertools import chain, combinations, permutations, product
import gurobipy as gp
from gurobipy import GRB
from copy import deepcopy

We also need to initialize a gurobi session, so we need to create a params dictionary with all the relevant information of our gurobi license.

WLSACCESSID = '< Copy your WLSACCESSID here >'
WLSSECRET = '< Copy your WLSSECRET here >'
LICENSEID = '< Copy your LICENSEID here >'

params = {
"WLSACCESSID": WLSACCESSID,
"WLSSECRET": WLSSECRET,
"LICENSEID": LICENSEID,
}
env = gp.Env(params=params) #this line initializes the gurobi environment

Input data & pre-processing

For the sake of simplicity, we are going to work with only 34 items of 7 different dimensions, which are kept as small values for illustration purposes. Please note that in my previous experience, I had to deal with hundreds of assets of massive sizes, some of which were literally heavy load trucks. The dimensions of the example items are shown in Table 1.

To ensure that the data is available for use, we will store the different values in lists and then prepare an extended dataframe with each item as a separate row. We can achieve this with the following code:

widths = [4,1,2,3,2,6,10]
heights =[3,1,1,2,4,2,4]
quants = [4,10,8,5,2,3,2]

WIDTHS = []
HEIGHTS = []
for q,w,h in zip(quants,widths,heights):
  for x in range(q):
    HEIGHTS.append(w)
    WIDTHS.append(h)

data_df = pd.DataFrame()
data_df['HEIGHTS'] = HEIGHTS
data_df['WIDTHS'] = WIDTHS

N = len(data_df) # number of items
top = max(data_df['HEIGHTS'].sum(),data_df['WIDTHS'].sum()) # maximum value of X or Y
M = top # to be used as big M for the OR constraints

I = range(N) # for the indexes of each asset "i"
K = range(4) # for the index of the OR variables "b"

Solving the problem with Gurobipy

Now that we have all the required libraries and input data, we can proceed to create the mathematical model described in Figure 3. To implement this model using GurobiPy, we just need to follow the code snippet below:

model = gp.Model("Assortment",env=env)

# (x,y) Coordinate variables
x = model.addVars(I,lb = 0,ub = top,vtype=GRB.CONTINUOUS, name="x")
y = model.addVars(I,lb = 0,ub = top,vtype=GRB.CONTINUOUS, name="y")

# Rotation variables 
R = model.addVars(I,vtype=GRB.BINARY,name = 'R')

X = model.addVar(lb=0,ub = top,vtype = GRB.CONTINUOUS,name = "X")
Y = model.addVar(lb=0,ub = top,vtype = GRB.CONTINUOUS,name = "Y")

# b variables for OR condition
b_vars = [(i,j,k) for i in I for j in I if j!=i for k in K]
B = model.addVars(b_vars,vtype = GRB.BINARY,name = "B")

# Objective function
model.setObjective(X*Y,GRB.MINIMIZE);

# constraints (1) and (2)
for i in I:

  model.addConstr(X >= x[i] + WIDTHS[i]*R[i] + (1-R[i])*HEIGHTS[i])
  model.addConstr(Y >= y[i] + HEIGHTS[i]*R[i] + (1-R[i])*WIDTHS[i])

# Constraints (3-7)
for i in I:
  for j in I:
    if i == j:
      continue
    else:
      #constraint (3)
      model.addConstr(x[i] + WIDTHS[i]*R[i] + (1-R[i])*HEIGHTS[i] <= x[j] + M*(1-B[i,j,0]))
      #constraint (4)
      model.addConstr(x[j] + WIDTHS[j]*R[j] + (1-R[j])*HEIGHTS[j] <= x[i] + M*(1-B[i,j,1]))
      #constraint (5)
      model.addConstr(y[i] + HEIGHTS[i]*R[i] + (1-R[i])*WIDTHS[i] <= y[j] + M*(1-B[i,j,2]))
      #constraint (6)
      model.addConstr(y[j] + HEIGHTS[j]*R[j] + (1-R[j])*WIDTHS[j] <= y[i] + M*(1-B[i,j,3]))
      #constraint (7)
      model.addConstr(B[i,j,0] + B[i,j,1] + B[i,j,2] + B[i,j,3] >= 1)

Please note that constraints (8–11) are established during the variable creation lines, so we don’t need to explicitly add them via the model.addConstr method. Now that the model is complete, we can proceed to solve it. We will let this model run for 600 seconds (10 minutes) since this is a complicated problem to solve. We will also set up a MIP gap of 5%, meaning if the optimization process is within 5% of the optimum value, we can terminate the process and get the best solution obtained so far. To solve the problem, follow the code below:

tl = 600
mip_gap = 0.05

model.setParam('TimeLimit', tl)
model.setParam('MIPGap', mip_gap)
model.optimize()

After 600 seconds, we are still far from the optimal solution, but the used area has decreased substantially. It has gone from an initial value of 620 (please remember this value, as it will be relevant in the next section) to a final value of 238. That’s almost one-third of the initial solution, which is not bad.

With the problem solved, we just need to extract the solution and plot it. We can easily accomplish this using the Rectangle object from matplotlib by following the code below:

all_vars = model.getVars()
values = model.getAttr("X", all_vars)
names = model.getAttr("VarName", all_vars)

obj = round(model.getObjective().getValue(),0)

total_X = int(round((X.x),0))
total_Y = int(round((Y.x),0))

fig, ax = plt.subplots()

for item in I:

  coords = (x[item].x,y[item].x)

  if R[item].x <= 0.01:
    wid = HEIGHTS[item]
    hig = WIDTHS[item]
  else:
    wid = WIDTHS[item]
    hig = HEIGHTS[item]

  ax.add_patch(Rectangle(coords, wid, hig,
            edgecolor = 'black',
            facecolor = "Grey",
            fill = True,
            alpha = 0.5,
            lw=2))
ax. set_xlim(0, total_X )
ax. set_ylim(0, total_Y )

ax.set_xticks(range(total_X+1))
ax.set_yticks(range(total_Y+1))
ax.grid()
ax.set_title(f" Total area {total_X} x {total_Y} = {int(obj)}")

plt.show()

We can visualize the obtained solution of the assortment from Figure 5 below.

Finally, we can save the solution into a dataframe, which can later be downloaded as a CSV or Excel file if we want to share it with a colleague. We can easily do this by following the code snippet below:

output_list = []
for i in I:

  print(f"item {i} x:{x[i].x}, y:{y[i].x}, Rotated:{R[i].x <= 0.01}")
  row = {'item':i,'x':round(x[i].x,2),'y':round(y[i].x,2),'Rotated':R[i].x <= 0.01}
  output_list.append(row)

output_df = pd.DataFrame(output_list)
output_df.to_csv("output_solution.csv") # to download the solution as a .csv file

Constructive Heuristics solutions

As I mentioned in the introduction, although I would have loved to solve this problem optimally using the mathematical programming formulation from the previous section, I wasn’t aware of these techniques when I first faced this challenge. However, that doesn’t mean I didn’t provide a solution. In fact, I did arrange all the assets, and they were enclosed in a safe area. The solution wasn’t optimal, but it was a solution nonetheless. How did I do it? Most likely, I used a heuristic or a combination of various heuristics.

Formally, a heuristic is "any approach to problem solving that employs a pragmatic method that is not fully optimized, perfected, or rationalized, but is nevertheless good enough as an approximation or attribute substitution." We might not realize it, but we constantly apply heuristics to problem-solving, and this assortment problem is no exception. We can build various constructive heuristics. In the next section, we will discuss two of them.

FFDH: First Fit Decreasing Height Heuristic

With this approach, we simply sort the rectangles by height and then place each asset in the first available position that fits, moving from left to right. To implement this heuristic, we will create a new class in Python and then develop the heuristic as a function that works with instances of that class. We can easily create the new class and instantiate each asset using the code below:

# definition of the new class
class Rectangle_class:
    # constructor
    def __init__(self, width, height,index):
        self.width = width
        self.height = height
        self.x = 0
        self.y = 0
        self.index = index

# create a list of all the assets
rectangles = [] 

# initialize a rectangle instances for each asset
for i in range(len(data_df)): 

  h,w = data_df.iloc[i,0], data_df.iloc[i,1]
  REC = Rectangle_class(w,h,i)
  rectangles.append(REC)

The heuristic function can be implemented using the code below:

def ffdh(rectangles):
    # Sort rectangles by height in descending order
    rectangles.sort(key=lambda rect: rect.height, reverse=True)

    # Initialize variables to keep track of the positions
    current_y = 0
    current_x = 0
    row_height = 0
    total_width = 0

    for rect in rectangles:
        if current_x + rect.width > total_width:
            total_width = current_x + rect.width
        # If rectangle fits in the current row
        if current_x + rect.width <= total_width:
            rect.x = current_x
            rect.y = current_y
            current_x += rect.width
            row_height = max(row_height, rect.height)
        else:
            # Move to the next row
            current_y += row_height
            rect.x = 0
            rect.y = current_y
            current_x = rect.width
            row_height = rect.height

    total_height = current_y + row_height
    return rectangles, total_width, total_height

This function returns the modified list of the assets with their new x and y coordinates, but it also returns the total width and height used to store the assets. Let’s see what we get when we use this heuristic with our problem of 34 assets.

# call the function and store the outputs
packed_rectangles, total_width, total_height = ffdh(rectangles)

total_X = total_width
total_Y = total_height

# calculate the final area used
obj = total_X*total_Y

# plot the solution
fig, ax = plt.subplots(figsize=(16, 6))

for rect in packed_rectangles:

  coords = (rect.x,rect.y)
  wid,hig = rect.width,rect.height

  ax.add_patch(Rectangle(coords, wid, hig,
            edgecolor = 'black',
            facecolor = "Grey",
            fill = True,
            alpha = 0.5,
            lw=2))
ax. set_xlim(0, total_X )
ax. set_ylim(0, total_Y )

ax.set_xticks(range(total_X+1))
ax.set_yticks(range(total_Y+1))
ax.grid()
ax.set_title(f" Total area {total_X} x {total_Y} = {int(obj)}")
plt.xticks(rotation=30)
plt.show()

The obtained output solution is shown in Figure 6 below:

The solution provided by this heuristic has a total area of 620, which coincidentally was also the starting point for Gurobi. Before the solver starts the branch-and-bound algorithm, it applies various heuristics to refine the initial solution, and most likely this heuristic was one of them.

Shelf Heuristic

This heuristic arranges the assets on a plane with a maximum width constraint (that in this case is provided by us), sorting them by height in descending order. It places the assets on a current shelf until the width limit is reached, then moves to a new shelf, updating the total width and height accordingly. The heuristic is implemented using the function below:

def shelf_heuristic(rectangles, max_width):
    # Sort rectangles by height in descending order
    rectangles.sort(key=lambda rect: rect.height, reverse=True)

    current_x = 0
    current_y = 0
    shelf_height = 0

    for rect in rectangles:
        if current_x + rect.width > max_width:
            # Move to a new shelf
            current_x = 0
            current_y += shelf_height
            shelf_height = 0

        # Place the rectangle
        rect.x = current_x
        rect.y = current_y
        current_x += rect.width
        shelf_height = max(shelf_height, rect.height)

    # total_width = max_width
    # total_height = current_y + shelf_height

    total_width = max([rec.x + rec.width for rec in rectangles])
    total_height = max([rec.y + rec.height for rec in rectangles])

    return rectangles, total_width, total_height

Let’s see what we get when we use this heuristic with our problem of 34 assets, and using a shelf width of 20.

container_width = 20
packed_rectangles,total_width, total_height = shelf_heuristic(rectangles, container_width)

total_X = total_width
total_Y = total_height

obj = total_X*total_Y

# generating the plot of the solution
fig, ax = plt.subplots(figsize=(16, 6))

for rect in packed_rectangles:

  coords = (rect.x,rect.y)
  wid,hig = rect.width,rect.height

  ax.add_patch(Rectangle(coords, wid, hig,
            edgecolor = 'black',
            facecolor = "Grey",
            fill = True,
            alpha = 0.5,
            lw=2))
ax. set_xlim(0, total_X )
ax. set_ylim(0, total_Y )

ax.set_xticks(range(total_X+1))
ax.set_yticks(range(total_Y+1))
ax.grid()
ax.set_title(f" Total area {total_X} x {total_Y} = {int(obj)}")
plt.xticks(rotation=30)
plt.show()

The obtained output solution is shown in Figure 7 below:

The solution provided by this heuristic has a total area of 340, which is almost half of the area obtained from FFDH. This is impressive, especially considering that heuristics like this one run very quickly and that we are not allowing for rotation of the assets, which would likely improve the quality of the heuristic significantly. Please note that there are literally dozens of heuristics for this problem, we just covered two popular ones; for more information, please check references [6–8].

You might be thinking, wow Luis, it would be nice if the solver could start from this 340 area instead of the previous 620, and you would be totally right. That’s why we will cover this in the next section.

Solving the problem with a warm start from a heuristic

Modern solvers like Gurobi allow users to provide a "warm start" solution when setting up their model. This initial solution helps make the optimization process more efficient. From the solution obtained from our shelf heuristic, we can extract the values and transform the solution into the format expected by the model. Remember, this includes not only the coordinates but also the total width, height, rotation variables, and the overlapping variables b_i_j_k. We can easily transform a solution from a heuristic by using the code snippet below:

heuristic_dict = dict()
for rect in packed_rectangles: 
  index = rect.index
  heuristic_dict[index] = rect

b_values = dict()
for i in I:
  for j in I:
    if i == j:
      continue
    else:
      rect_i = heuristic_dict[i]
      rect_j = heuristic_dict[j]

      # initialize values of b at 0
      b_values[(i,j,0)] = 0
      b_values[(i,j,1)] = 0
      b_values[(i,j,2)] = 0
      b_values[(i,j,3)] = 0

      # check the different conditions and stablish the real values of b
      if rect_i.x + rect_i.width <= rect_j.x:
        b_values[(i,j,0)] = 1

      if rect_j.x + rect_j.width <= rect_i.x:
        b_values[(i,j,1)] = 1

      if rect_i.y + rect_i.height <= rect_j.y:

        b_values[(i,j,2)] = 1

      if rect_j.y + rect_j.height <= rect_i.y:

        b_values[(i,j,3)] = 1   

# store the solution values of x,y and r
x_dict = dict()
y_dict = dict()
r_dict = dict()

for i in I:

  REC = heuristic_dict[i]

  x_dict[i] = REC.x
  y_dict[i] = REC.y
  r_dict[i] = 1 # all assets are not rotated

Now that we have extracted the heuristic solution and transformed it into the proper form expected by our model, we can proceed to re-run the model, adding this solution as a starting point. The code below shows how to create the new model and add the warm start solution.

# Same model as before
model = gp.Model("Assortment_warm_start",env=env)

x = model.addVars(I,lb = 0,ub = top,vtype=GRB.CONTINUOUS, name="x")
y = model.addVars(I,lb = 0,ub = top,vtype=GRB.CONTINUOUS, name="y")

R = model.addVars(I,vtype=GRB.BINARY,name = 'R')

X = model.addVar(lb=0,ub = top,vtype = GRB.CONTINUOUS,name = "X")
Y = model.addVar(lb=0,ub = top,vtype = GRB.CONTINUOUS,name = "Y")

b_vars = [(i,j,k) for i in I for j in I if j!=i for k in K]
B = model.addVars(b_vars,vtype = GRB.BINARY,name = "B")

model.setObjective(X*Y,GRB.MINIMIZE);

for i in I:

  model.addConstr(X >= x[i] + WIDTHS[i]*R[i] + (1-R[i])*HEIGHTS[i])
  model.addConstr(Y >= y[i] + HEIGHTS[i]*R[i] + (1-R[i])*WIDTHS[i])

for i in I:
  for j in I:
    if i == j:
      continue
    else:
      model.addConstr(x[i] + WIDTHS[i]*R[i] + (1-R[i])*HEIGHTS[i] <= x[j] + M*(1-B[i,j,0]))
      model.addConstr(x[j] + WIDTHS[j]*R[j] + (1-R[j])*HEIGHTS[j] <= x[i] + M*(1-B[i,j,1]))

      model.addConstr(y[i] + HEIGHTS[i]*R[i] + (1-R[i])*WIDTHS[i] <= y[j] + M*(1-B[i,j,2]))
      model.addConstr(y[j] + HEIGHTS[j]*R[j] + (1-R[j])*WIDTHS[j] <= y[i] + M*(1-B[i,j,3]))

      model.addConstr(B[i,j,0] + B[i,j,1] + B[i,j,2] + B[i,j,3] >= 1)

# adding the warm start solution from the heuristic
for var in x:
  x[var].Start = x_dict[var]
  y[var].Start = y_dict[var]

  R[var].Start = r_dict[var]

X.Start = total_X
Y.start = total_Y

for var in B:
  B[var] = b_values[var]

Now that the model was initialized with the warm start, we can proceed to solve it.

tl = 600
mip_gap = 0.05

model.setParam('TimeLimit', tl)
model.setParam('MIPGap', mip_gap)
model.optimize()

If everything goes according to plan you should see a message from the solver, saying that it is going to start the optimization process from the solution provided by the user. Like in Figure 8 below.

Conclusions

In conclusion, this article has taken you through a detailed journey on how to solve a Two-Dimensional assortment problem using mathematical programming. We discussed a Mixed Integer Linear Programming formulation for this problem, as well as how to tackle it using some classical constructive heuristics. Moreover we covered how to include the solutions obtained from these heuristics as warm start for our solver. I hope this article serves as a helpful resource for anyone facing similar challenges in their current jobs. I would have loved 💙 to have had this knowledge back then, but in any case it’s never too late to learn 😄 . You can find the entire notebook with the code developed for this article in the link to my GitHub repository just below.

Youtube_optimization/2D_Assortment_Problem_Medium_Article_for_github.ipynb at main ·…

I sincerely hope you found this article useful and entertaining. If so, I’d love to hear your thoughts! Please feel free to leave a comment or show your appreciation with a clap 👏 . And if you’re interested in staying updated on my latest articles, consider following me on Medium. Your support and feedback are what drive me to keep exploring and sharing in this fascinating field. Thank you for taking the time to read, and stay tuned for more insights in my next article!

---

References

• [1] Williams HP. Model Building in Mathematical Programming. fourth ed. Chichester: Wiley; 1999
• [2] Huang, Y.H., Lu, H.C., Wang, Y.C., Chang, Y.F. and Gao, C.K., 2020. A Global Method for a Two-Dimensional Cutting Stock Problem in the Manufacturing Industry. Application of Decision Science in Business and Management, p.167.
• [3] Li HL, Chang CT. An approximately global optimization method for assortment problems. European Journal of Operational Research. 1998;105:604–612
• [4] Martello S, Vigo D. Exact solution of the two-dimensional finite bin packing problem. Management Science. 1998;44:388–399
• [5] Chen CS, Sarin S, Balasubramanian R. A mixed-integer programming model for a class of assortment problems. European Journal of Operational Research. 1993;65:362–367

• [6] Apple, J.M., and Deisenroth, M.P., "A computerized plant layout and evaluation technique – PLANET", Proceedings of AIIE Annual Conference, May 1973, 121.

• [7] Buffa, E.S., Armour, G.C., and Vollman, T.E., "Allocating facilities with CRAFT", Harvard Business Review, March-April, 1964, 130.

• [8] Seehof, J.M., and Evans, W.O., "Automated layout design program (ALDEP)", Journal of Industrial Engineering De￾cember (1967) 690.