Car Sequencing (CARSEQ)#

`carseq.py`

```"""
## GAMSSOURCE: https://www.gams.com/latest/gamslib_ml/libhtml/gamslib_carseq.html
## MODELTYPE: MIP, MINLP
## KEYWORDS: mixed integer linear programming, mixed integer nonlinear programming, production planning, car manufacturing, line problem

Car Sequencing (CARSEQ)

A number of cars are to be produced; they are not identical, because
different options are available as variants on the basic model. The
assembly line has different stations which install the various options
(air-conditioning, sun-roof, etc.). These stations have been designed
to handle at most a certain percentage of the cars passing along the
assembly line. Furthermore, the cars requiring a certain option must not
be bunched together, otherwise the station will not be able to cope.
Consequently, the cars must be arranged in a sequence so that the capacity
of each station is never exceeded. For instance, if a particular station
can only cope with at most half of the cars passing along the line, the
sequence must be built so that at most 1 car in any 2 requires that option.
The problem has been shown to be NP-hard (Gent 1999).

This is the example given in Dincbas, et al. More instances can

Dincbas et al., Dincbas, M., Simonis, H., and Van Hentenryck, P.
Solving the car-sequencing problem in constraint logic programming.
In 8th European Conference on Artificial Intelligence (ECAI 88) ,
Y. Kodratoff, Ed. Pitmann Publishing, London, Munich, Germany, 290-295, 1988
"""

from __future__ import annotations

import os

import gamspy.math as gams_math
import numpy as np
import pandas as pd
from gamspy import (
Alias,
Card,
Container,
Equation,
Model,
Options,
Ord,
Parameter,
Sense,
Set,
Sum,
Variable,
)
from gamspy.math import ifthen

def main(mip=False):
classData_recs = np.array(
[
[1, 1, 0, 1, 1, 0],
[1, 0, 0, 0, 1, 0],
[2, 0, 1, 0, 0, 1],
[2, 0, 1, 0, 1, 0],
[2, 1, 0, 1, 0, 0],
[2, 1, 1, 0, 0, 0],
]
)
classData_recs = pd.DataFrame(
classData_recs,
columns=["numCars", "opt1", "opt2", "opt3", "opt4", "opt5"],
index=["class1", "class2", "class3", "class4", "class5", "class6"],
)
no_of_cars = classData_recs.numCars.sum()
classData_recs = classData_recs.reset_index().melt(
id_vars="index", var_name="Category", value_name="Value"
)

m = Container(
system_directory=os.getenv("SYSTEM_DIRECTORY", None),
)

# Sets
p = Set(
m,
name="p",
records=[f"pos{i}" for i in range(1, 11)],
description="position",
)
o = Set(
m,
name="o",
records=[f"opt{i}" for i in range(1, 6)],
description="options",
)
c = Set(
m,
name="c",
records=[f"class{i}" for i in range(1, 7)],
description="classes",
)

# Parameter
maxc = Parameter(
m,
name="maxc",
domain=o,
records=np.array([1, 2, 1, 2, 1]),
description="maximum number of cars with that option in a block",
)
bs = Parameter(
m,
name="bs",
domain=o,
records=np.array([2, 3, 3, 5, 5]),
description="block size to which the maximum number maxc refers",
)

classData = Parameter(
m,
name="classData",
domain=[c, "*"],
records=classData_recs,
description="class data",
)

if len(p) != no_of_cars:
raise Exception("inconsistent number of cars")

pp = Alias(m, name="pp", alias_with=p)

# Sets
blk = Set(
m,
name="blk",
domain=[o, p],
description="blocks of positions to monitor",
)
blkc = Set(
m,
name="blkc",
domain=[o, p, pp],
description="positions in the blocks",
)

blkc[o, p, pp].where[Ord(p) <= Card(p) - bs[o] + 1] = (
Ord(pp) >= Ord(p)
) & (Ord(pp) < Ord(p) + bs[o])
blk[o, p] = Sum(pp.where[blkc[o, p, pp]], 1)

# Variables
sumc = Variable(m, name="sumc", type="free", domain=[o, p])
cp = Variable(
m,
name="cp",
type="binary",
domain=[c, p],
description="class k is scheduled at position p",
)
op = Variable(
m,
name="op",
type="free",
domain=[o, p],
description="option o appears at position p",
)
v = Variable(
m,
name="v",
type="free",
domain=[o, p],
description="violations in a block",
)

if mip:
v.type = "positive"
op.type = "binary"

# Equations
defnumCars = Equation(
m,
name="defnumCars",
domain=c,
description="exactly numCars of class c assigned to positions",
)
defoneCar = Equation(
m,
name="defoneCar",
domain=p,
description="one car assigned to each position p",
)
defop = Equation(
m,
name="defop",
domain=[o, p],
description="option o appears at position p",
)
defopLS = Equation(
m,
name="defopLS",
domain=[o, p],
description="option o appears at position p",
)
defviol = Equation(
m, name="defviol", domain=[o, p], description="violations in a block"
)
defviolLS = Equation(
m, name="defviolLS", domain=[o, p], description="violations in a block"
)
defsumc = Equation(m, name="defsumc", domain=[o, p])

defnumCars[c] = Sum(p, cp[c, p]) == classData[c, "numCars"]
defoneCar[p] = Sum(c, cp[c, p]) == 1
defop[o, p] = Sum(c.where[classData[c, o]], cp[c, p]) <= op[o, p]
defsumc[o, p] = sumc[o, p] == Sum(c.where[classData[c, o]], cp[c, p])
defopLS[o, p] = op[o, p] == ifthen(sumc[o, p] >= 0.5, 1, 0)
defviol[blk[o, p]] = Sum(blkc[blk, pp], op[o, pp]) <= maxc[o] + v[o, p]
defviolLS[blk[o, p]] = v[o, p] == gams_math.Max(
Sum(blkc[blk, pp], op[o, pp]) - maxc[o], 0
)

obj = Sum(blk[o, p], v[o, p])

# Model
carseqMIP = Model(
m,
name="carseqMIP",
equations=[
defnumCars,
defoneCar,
defop,
defviol,
],
problem="mip",
sense=Sense.MIN,
objective=obj,
)
carseqLS = Model(
m,
name="carseqLS",
equations=[
defnumCars,
defoneCar,
defsumc,
defopLS,
defviolLS,
],
problem="minlp",
sense=Sense.MIN,
objective=obj,
)

if mip:
v.type = "positive"
op.type = "binary"
carseqMIP.solve(options=Options(relative_optimality_gap=0))
else:
carseqLS.solve(options=Options(relative_optimality_gap=0))

rep = Parameter(m, name="rep", domain=[p, c, o])
rep[p, c, o].where[(cp.l[c, p] > 0.5)] = classData[c, o]

print("Objective Function Value: ", carseqLS.objective_value)

if __name__ == "__main__":
main()
```