Solver Max - Warehouse space for free: Linearized model

18 March 2024

In the first part of this article series, we formulated and implemented Model 1 for redesigning a pallet warehouse. Model 1, which in non-linear, works OK for small data sets, but it fails to find a good solution for the full data set.

In this article, we linearize our model using a standard technique. The idea is that, all going well, a linear model will be solvable at the scale we need.

Download the models

The models described in this series of articles are built in Python using the Pyomo library.

The files are available on GitHub.

Situation

The situation is described in the first part of this series: Warehouse space for free: Non-linear model.

Model 2 design

Model 1 includes a non-linear constraint, where we multiply the integer shelf height variables by the binary variables that allocate pallets to shelves, as shown in Figure 1. We want to linearize this constraint.

Figure 1. Equation (2) in Model 1

\begin{alignat}{1} &\text{Subject to} \\ &\quad \text{Pallets}_{p} \quad \le \quad \sum_{s=1}^{n} ShelfHeights_{s} \times Allocation_{p, s} \quad \forall \ p \in P \tag{2}\\ \end{alignat}

A method for linearizing this constraint is described in FICO Xpress Optimization Help, Product values. Specifically, after a bit of algebraic rearranging, the method is defined as shown in Figure 2.

Figure 2. Linearization of variable multiplication

\begin{alignat}{1} &\text{Given:}\\ &\quad \text{continuous variable } x, \text{where } x \text{ is between bounds }\text{L} \le x \le \text{U}\\ &\quad \text{binary variable } d\\ \\ &\text{The multiplication:}\\ &\quad y = x \times d\\ \\ &\text{Is equivalent to:}\\ \end{alignat} \begin{alignat}{1} &\text{Part 1 LB:} \quad &y \quad \ge \quad &\text{L} \times d \\ &\text{Part 2 UB:} \quad &y \quad \le \quad &\text{U} \times d \\ &\text{Part 3 UB:} \quad &y \quad \le \quad &x - \text{L} \times (1-d) \\ &\text{Part 4 LB:} \quad &y \quad \ge \quad &x - \text{U} \times (1-d) \\ \end{alignat}

In our situation, \(x\) is integer rather than continuous, but the method still works.

That is, we replace the non-linear constraint with an integer variable and four linear constraints.

Model 2 Formulation

We formulate Model 2 as an extension of Model 1, as shown in Figure 3.

Firstly, we remove the non-linear constraint from Model 1. Then we create new integer variables \(PalletShelf_{p, s}\), used for the pallet and shelf linearization, equivalent to the variable \(y\) in Figure 2 We also create a new constraint that performs the same purpose as Equation (2). Finally, we add four new constraints to ensure that our new constraint behaves as we want. That is:

Equation (2). We remove the existing non-linear constraint Equation (2).
Equation (24). Each pallet must be allocated to a shelf that is at least the height of the pallet. This constraint, in association with Equations (25) to (28), replaces Equation (2).
Equation (25). Part 1 LB, lower bound of allocation.
Equation (26). Part 1 UB, upper bound of allocation.
Equation (27). Part 2 UB, upper bound of not allocation.
Equation (28). Part 2 LB, lower bound of not allocation.

Figure 3. Model 2 formulation

\begin{alignat}{1} &\text{Subject to} \\ &\enclose{horizontalstrike}{ \quad \text{Pallets}_{p} \quad \le \quad \sum_{s=1}^{n} ShelfHeights_{s} \times Allocation_{p, s} \quad \forall \ p \in P \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \tag{2} }\\ \end{alignat} \begin{alignat}{1} &\quad \text{Pallets}_{p} \quad &\le &\quad \sum_{s=1}^{n} PalletShelf_{p, s} \quad &\forall \ p \in P \tag{24} \\ &\quad PalletShelf_{p, s} \quad &\ge &\quad 0 \times Allocation_{p, s} \quad &\forall \ p \in P, \ s \in S \tag{25} \\ &\quad PalletShelf_{p, s} \quad &\le &\quad \text{MaxShelfSize} \times Allocation_{p, s} \quad &\forall \ p \in P, \ s \in S \tag{26} \\ &\quad PalletShelf_{p, s} \quad &\le &\quad ShelfHeights_{s} - 0 \times (1 - Allocation_{p, s}) \quad &\forall \ p \in P, \ s \in S \tag{27} \\ &\quad PalletShelf_{p, s} \quad &\ge &\quad ShelfHeights_{s} - \text{MaxShelfSize} \times (1 - Allocation_{p, s}) \quad &\forall \ p \in P, \ s \in S \tag{28} \\ \\ &\text{Variables} \\ &\quad PalletShelf_{p, s} \quad &= &\quad \text{Pallet and shelf linearization, integer} \quad &\forall \ p \in P, \ s \in S \tag{29} \\ \end{alignat}

The bounds on the shelf heights are \(0\) (because some shelf positions may not be used) and \(\text{MaxShelfSize}\). For completeness, we include the multiplications by zero – on the assumption that the solver will remove those terms.

Note that we need to carefully translate the summation and indices from Equation (2) to the replacement constraints, to ensure that they are equivalent.

Model 2 implementation

Given the additional formulation for Model 2, the corresponding code is shown in Figure 4.

We define the new \(PalletShelf_{p, s}\) variable, then we use Pyomo's .deactivate() function to remove Equation (2), PalletFits, from the model. We also define the new PalletFitsLinear constraint, along with its four accompanying constraints that act as bounds on the new variable.

Figure 4. Model 2 code

Model.PalletShelf = pyo.Var(Model.P, Model.S, within = pyo.NonNegativeIntegers, bounds = (0, Model.MaxShelfSize), initialize = Model.MaxShelfSize)

Model.PalletFits.deactivate()   # This constraint is replaced by the following constraints

def rule_fitlinear(Model, p):   # Each pallet must be allocated to a shelf that is at least the height of the pallet
    return sum(Model.PalletShelf[p, s] for s in Model.S) >= Model.Pallets[p]
Model.PalletFitsLinear = pyo.Constraint(Model.P, rule = rule_fitlinear)

def rule_fitLB1(Model, p, s):   # Linearization of fit constraint, part 1
    return 0 * Model.Allocation[p, s] <= Model.PalletShelf[p, s]
Model.PalletFitsLB1 = pyo.Constraint(Model.P, Model.S, rule = rule_fitLB1)

def rule_fitUB1(Model, p, s):   # Linearization of fit constraint, part 2
    return Model.MaxShelfSize * Model.Allocation[p, s] >= Model.PalletShelf[p, s]
Model.PalletFitsUB1 = pyo.Constraint(Model.P, Model.S, rule = rule_fitUB1)

def rule_fitUB2(Model, p, s):   # Linearization of fit constraint, part 3
    return Model.ShelfHeights[s] - Model.MaxShelfSize * (1 - Model.Allocation[p, s]) <= Model.PalletShelf[p, s]
Model.PalletFitsUB2 = pyo.Constraint(Model.P, Model.S, rule = rule_fitUB2)
    
def rule_fitLB2(Model, p, s):   # Linearization of fit constraint, part 4
    return Model.ShelfHeights[s] - 0 * (1 - Model.Allocation[p, s]) >= Model.PalletShelf[p, s]
Model.PalletFitsLB2 = pyo.Constraint(Model.P, Model.S, rule = rule_fitLB2)

Model 2 solutions

Running Model 2 with the Gurobi solver with each of the data sets, we get the solutions shown in Figure 5.

Figure 5. Solutions for Model 2

That is:

The data with 200 pallets quickly solves to optimality, though the solve time is longer than Model 1 (which took 2 seconds).
The 2,000 pallet data performance is worse: Gurobi finds the optimal solution after 23 minutes, but takes almost 8 hours to prove optimality. Model 1 took about 5 minutes to prove optimality.
With the full data set of 20,000 pallets, after 8 hours Gurobi's solution is slightly worse than the current design, having 1,004 racks each with 5 shelves. Progress is very slow. Model 1 also struggled with this full data set.

Linearizing a model usually makes it easier and faster to solve. We have used this method successfully on many models, typically transforming a model from intractable to solvable. But not this model.

The issue is that, compared with Model 1, our linearized Model 2 has twice as many variables and four times as many constraints. Our linearization is not very "tight", meaning that many of our constraints are a long way from an optimal solution, so the solver has a lot of work to do. In contrast, Gurobi clearly has an efficient method for solving the quadratic Model 1, so it performs better than Model 2.

Since Model 2 is linear, we can also solve it using the free HiGHS solver, rather than using the commercial Gurobi solver. HiGHS takes several times longer to solve Model 2.

Next steps: Simplify our model

Model 2 works, but it is very slow. We failed to find a good solution for the full data set.

We need a new approach. In the first part of this series, we said that Model 1 has many moving parts. Linearizing Model 1, to form Model 2, added even more moving parts.

When we reconsider the situation, we note that there are relatively few potential values for some of the decisions being made by Models 1 and 2. Specifically, the number of shelves per rack is, realistically, an integer from 6 to 9, inclusive. We also suspect that the number of combinations of shelf sizes that sum to exactly 6.0m high is relatively small – hopefully small enough that we can enumerate them all. Therefore, we have an opportunity to simplify our model.

In the next article, we extract some of the decisions from the model and make them exogenous inputs. We enumerate all possible values for those decisions, then iterate over those values using a simplified model. Perhaps that will be sufficient to make the model solvable in a reasonable time?

Conclusion

In this article, we linearize our model with an expectation that it will be easier to solve than the non-linear Model 1.

But it isn't. Model 2's run time is substantially longer, and consequently we fail to find a good solution for the full data set.

Fortunately, we have some ideas for simplifying the model by making some decisions exogenous. So, in the next article, we build a simpler model. Hopefully, that model will perform better at a large scale, enabling us to find a solution to our warehousing problem.

If you would like to know more about this model, or you want help with your own models, then please contact us.

Warehouse space for free: Linearized model