BIM production with perturbed data#

# install dependencies and select solver
%pip install -q amplpy

SOLVER = "highs"

from amplpy import AMPL, ampl_notebook

ampl = ampl_notebook(
    modules=["highs"],  # modules to install
    license_uuid="default",  # license to use
)  # instantiate AMPL object and register magics

Using default Community Edition License for Colab. Get yours at: https://ampl.com/ce
Licensed to AMPL Community Edition License for the AMPL Model Colaboratory (https://colab.ampl.com).

Problem description#

The company BIM realizes that a \(1\%\) fraction of the copper always gets wasted while producing both types of microchips, more specifically \(1\%\) of the required amount. This means that it actually takes \(4.04\) gr of copper to produce a logic chip and \(2.02\) gr of copper to produce a memory chip. If we rewrite the linear problem of the basic BIM problem and modify accordingly the coefficients in the corresponding constraints, we obtain the following problem

\[\begin{split} \begin{align*} \max \quad & 12x_1+9x_2 \\ \text{s.t.} \quad & x_1 \leq 1000 & \text{(silicon)}\\ & x_2 \leq 1500 & \text{(germanium)}\\ & x_1+x_2 \leq 1750 & \text{(plastic)}\\ & 4.04 x_1+2.02 x_2 \leq 4800 & \text{(copper with waste)}\\ &x_1, x_2 \geq 0. \end{align*} \end{split}\]

If we solve it, we obtain a different optimal solution than the original one, namely \((x_1,x_2) \approx (626.238,1123.762)\) and an optimal value of roughly \(17628.713\). Note, in particular, that this new optimal solution is not integer, but on the other hand in the linear optimization problem above there is no constraint requiring \(x_1\) and \(x_2\) to be such.

%%writefile BIM_perturbed_LO.mod

var x1 >= 0;
var x2 >= 0;

maximize profit: 12 * x1 + 9 * x2;

s.t. silicon: x1 <= 1000;
s.t. germanium: x2 <= 1500;
s.t. plastic: x1 + x2 <= 1750;
s.t. copper: 4.04 * x1 + 2.02 * x2 <= 4800;

Overwriting BIM_perturbed_LO.mod

m = AMPL()
m.read("BIM_perturbed_LO.mod")

m.option["solver"] = SOLVER
m.solve()

print(
    "x = ({:.3f}, {:.3f}), optimal value = {:.3f}".format(
        m.var["x1"].value(), m.var["x2"].value(), m.obj["profit"].value()
    )
)

HiGHS 1.5.3: HiGHS 1.5.3: optimal solution; objective 17628.71287
2 simplex iterations
0 barrier iterations
 
x = (626.238, 1123.762), optimal value = 17628.713

In terms of production, of course, we would simply produce entire chips but it is not clear how to implement the fractional solution \((x_1,x_2) \approx (626.238,1123.762)\). Rounding down to \((x_1,x_2) = (626,1123)\) will intuitively yield a feasible solution, but we might be giving away some profit and/or not using efficiently the available material. Rounding up to \((x_1,x_2) = (627,1124)\) could possibly lead to an unfeasible solution for which the available material is not enough. We can of course manually inspect by hand all these candidate integer solutions, but if the problem involved many more decision variables or had a more complex structure, this would become much harder and possibly not lead to the true optimal solution.

A much safer approach is to explicitly require the two decision variables to be non-negative integers, thus transforming the original into the following mixed-integer linear optimization (MILO) problem:

\[\begin{split} \begin{align*} \max \quad & 12x_1+9x_2 \\ \text{s.t.} \quad & x_1 \leq 1000 & \text{(silicon)}\\ & x_2 \leq 1500 & \text{(germanium)}\\ & x_1+x_2 \leq 1750 & \text{(plastic)}\\ & 4.04 x_1+2.02 x_2 \leq 4800 & \text{(copper with waste)}\\ &x_1, x_2 \in \mathbb{N}. \end{align*} \end{split}\]

The optimal solution is \((x_1,x_2) = (626,1124)\) with a profit of \(17628\). Note that for this specific problem both the naive rounding strategies outlined above would have not yielded the true optimal solution. The Python code for obtaining the optimal solution using MILO solvers is given below.

%%writefile BIM_perturbed_MILO.mod

var x1 integer >= 0;
var x2 integer >= 0;

maximize profit: 12 * x1 + 9 * x2;

s.t. silicon: x1 <= 1000;
s.t. germanium: x2 <= 1500;
s.t. plastic: x1 + x2 <= 1750;
s.t. copper: 4.04 * x1 + 2.02 * x2 <= 4800;

Overwriting BIM_perturbed_MILO.mod

m = AMPL()
m.read("BIM_perturbed_MILO.mod")

m.option["solver"] = SOLVER
m.solve()

print(
    "x = ({:.3f}, {:.3f}), optimal value = {:.3f}".format(
        m.var["x1"].value(), m.var["x2"].value(), m.obj["profit"].value()
    )
)

HiGHS 1.5.3: HiGHS 1.5.3: optimal solution; objective 17628
2 simplex iterations
1 branching nodes
 
x = (626.000, 1124.000), optimal value = 17628.000