Solving your first business problem in C++

The car sequencing problem is a well-known combinatorial optimization problem related to the organization of the production line of a car manufacturer. It involves scheduling a set of cars which are divided into classes by different sets of options. The assembly line has different stations where the various options (air-conditioning, sun-roof, etc.) are installed by teams. Each such station has a certain workload capacity given by the parameters P and Q: the station is able to handle P cars with the given option in each sequence of Q cars. The goal of the problem is to produce of a sequence of cars respecting the given demand for each class and such that overcapacities are minimized. In other words, for each option P/Q and for each window of length Q, a penalty max(0,N-P) is counted where N is the number of such options actually scheduled in this window.

In this section, we present a simple model for this problem. We illustrate how to read the input data from a file, we show how to model the Car Sequencing Problem by clearly defining the decisions, constraints, and objective of the problem, and we introduce some of the parameters that can be used to control the search. The full program for this problem is available in our example tour, together with its Hexaly Modeler, Python, Java, and C# versions, plus a certain number of instance files.

Reading input data

Data files for this problem are available in folder hexaly/examples/carsequencing and more files can be downloaded from http://www.csplib.org (problem number 001). The format of these files is:

  • 1st line: number of cars; number of options; number of classes.

  • 2nd line: for each option, the maximum number of cars with this option in a block (parameter P).

  • 3rd line: for each option, the block size to which the maximum number refers (parameter Q).

  • Then for each class: index of the class; number of cars in this class; for each option, whether or not this class requires it (0 or 1).

Below is the code to read the input data:

ifstream infile;
infile.exceptions(ifstream::failbit | ifstream::badbit);
infile.open(fileName.c_str());

infile >> nbPositions;
infile >> nbOptions;
infile >> nbClasses;

P.resize(nbOptions);
for (int o = 0; o < nbOptions; o++)
    infile >> P[o];

Q.resize(nbOptions);
for (int o = 0; o < nbOptions; o++)
    infile >> Q[o];

options.resize(nbClasses);
nbCars.resize(nbClasses);
for (int c = 0; c < nbClasses; c++) {
    int ignored;
    infile >> ignored;
    infile >> nbCars[c];
    options[c].resize(nbOptions);
    for (int o = 0; o < nbOptions; o++) {
        int v;
        infile >> v;
        options[c][o] = (v == 1);
    }
}

Modeling the problem

One of the most important steps when modeling an optimization problem with Hexaly is to define the decision variables. These decisions are such that instantiating them allows evaluating all other expressions of your model (in particular, the constraints and objectives). Hexaly Optimizer will try to find the best possible decisions with respect to the given objectives and constraints. Here, we want to decide the ordering of cars in the production line. In terms of 0-1 modeling, it amounts to defining nbClasses * nbPositions variables cp[c][p] such that cp[c][p] is equal to 1 when the car at position p belongs to class c and 0 otherwise:

for (int c = 0; c < nbClasses; c++) {
    for (int p = 0; p < nbPositions; p++) {
        cp[c][p] = model.boolVar();
    }
}

Defining the constraints of your model is another important modeling task. Generally speaking, constraints should be limited to strictly necessary requirements. In the car sequencing problem, it is physically impossible to have two cars at the same position. On the contrary, satisfying all P/Q ratios is not always possible, which is why we define it as a first-priority objective instead. Limiting the use of constraints is generally a good practice, as Hexaly Optimizer tends to benefit from denser search spaces. Therefore, symmetry-breaking constraints and other types of non-necessary constraints should also be avoided.

Here, we have two families of constraints expressing the assignment requirements: we must have exactly one car per position and the demand in cars for each class should be satisfied:

for (int c = 0; c < nbClasses; c++) {
    LSExpression nbCarsFromClass = model.sum(cp[c].begin(), cp[c].end());
    model.constraint(nbCarsFromClass == nbCars[c]);
}

for (int p = 0; p < nbPositions; p++) {
    LSExpression nbCarsOnPos = model.sum();
    for (int c = 0; c < nbClasses; c++) {
        nbCarsOnPos += cp[c][p];
    }
    model.constraint(nbCarsOnPos == 1);
}

Before writing the objective function, we introduce some intermediate expressions. First, we declare expressions op[o][p], which are equal to 1 when the car at position p has option o and 0 otherwise:

for (int o = 0; o < nbOptions; o++) {
    for (int p = 0; p < nbPositions; p++) {
        op[o][p] = model.or_();
        for (int c = 0; c < nbClasses; c++) {
            if (options[c][o]) op[o][p].addOperand(cp[c][p]);
        }
    }
}

Indeed, option o appears at position p if this position hosts one of the classes featuring this option. This ability to define intermediate expressions makes the model more readable and more efficient (if these intermediate expressions are used in several other expressions). Similarly, the number of times option o appears between position p and position p+Q[o]-1 is defined as:

for (int o = 0; o < nbOptions; o++) {
    for (int j = 0; j < nbPositions - Q[o] + 1; j++) {
        nbCarsInWindow[o][j] = model.sum();
        for (int k = 0; k < Q[o]; k++) {
            nbCarsInWindow[o][j] += op[o][j + k];
        }
    }
}

Finally, the penalty function on overcapacities can be written exactly as it was defined above in the specifications of the problem:

for (int o = 0; o < nbOptions; o++) {
    nbViolationsInWindow[o].resize(nbPositions - Q[o] + 1);
    for (int j = 0; j < nbPositions - Q[o] + 1; j++) {
        nbViolationsInWindow[o][j] = model.max(0, nbCarsInWindow[o][j] - P[o]);
    }
}

Note that we are using here a different way to express a sum. Indeed, you can define a sum with the sum operator, but you can also use arithmetic operators + and -. Ultimately, the last step consists in defining the objective function, which is the sum of all violations:

obj = model.sum();
for (int o = 0; o < nbOptions; o++) {
    obj.addOperands(nbViolationsInWindow[o].begin(), nbViolationsInWindow[o].end());
}
model.minimize(obj)

Here is the resulting code of the Hexaly model for the Car Sequencing Problem written using the C++ API:

LSModel model = localsolver.getModel();

// 0-1 decisions:
// cp[c][p] = 1 if class c is at position p, and 0 otherwise
cp.resize(nbClasses);
for (int c = 0; c < nbClasses; c++) {
    cp[c].resize(nbPositions);
    for (int p = 0; p < nbPositions; p++) {
        cp[c][p] = model.boolVar();
    }
}

// Constraints:
// For each class c, nbCars[c] are assigned to the different positions
for (int c = 0; c < nbClasses; c++) {
    LSExpression nbCarsFromClass = model.sum(classOnPos[c].begin(), classOnPos[c].end());
    model.constraint(nbCarsFromClass == nbCars[c]);
}

// Constraints: one car assigned to each position p
for (int p = 0; p < nbPositions; p++) {
    LSExpression nbCarsOnPos = model.sum();
    for (int c = 0; c < nbClasses; c++) {
        nbCarsOnPos += classOnPos[c][p];
    }
    model.constraint(nbCarsOnPos == 1);
}

// Intermediate expressions:
// op[o][p] = 1 if option o appears at position p, and 0 otherwise
vector<vector<LSExpression> > optionsOnPos;
optionsOnPos.resize(nbOptions);
for (int o = 0; o < nbOptions; o++) {
    optionsOnPos[o].resize(nbPositions);
    for (int p = 0; p < nbPositions; p++) {
        optionsOnPos[o][p] = model.or_();
        for (int c = 0; c < nbClasses; c++) {
            if (options[c][o]) optionsOnPos[o][p].addOperand(classOnPos[c][p]);
        }
    }
}

// Intermediate expressions: compute the number of cars in each window
vector<vector<LSExpression> > nbCarsInWindow;
nbCarsInWindow.resize(nbOptions);
for (int o = 0; o < nbOptions; o++) {
    nbCarsInWindow[o].resize(nbPositions - windowSize[o] + 1);
    for (int j = 0; j < nbPositions - windowSize[o] + 1; j++) {
        nbCarsInWindow[o][j] = model.sum();
        for (int k = 0; k < windowSize[o]; k++) {
            nbCarsInWindow[o][j] += optionsOnPos[o][j + k];
        }
    }
}

// Intermediate expressions: compute the number of violations in each window
vector<vector<LSExpression> > nbViolationsInWindow;
nbViolationsInWindow.resize(nbOptions);
for (int o = 0; o < nbOptions; o++) {
    nbViolationsInWindow[o].resize(nbPositions - windowSize[o] + 1);
    for (int j = 0; j < nbPositions - windowSize[o] + 1; j++) {
        nbViolationsInWindow[o][j] = model.max(0, nbCarsInWindow[o][j] - maxCarsPerWindow[o]);
    }
}

// Objective function: minimize the sum of violations for all options and all windows
obj = model.sum();
for (int o = 0; o < nbOptions; o++) {
    obj.addOperands(nbViolationsInWindow[o].begin(), nbViolationsInWindow[o].end());
}

model.minimize(obj);
model.close();

Parameterizing the optimizer

Some parameters can be defined after closing the model in order to control the search:

optimizer.getParam().setTimeLimit(30);
optimizer.getParam().setTimeBetweenDisplays(5);

All control parameters can be found in the HxParam class.

Having launched Hexaly Optimizer, you should observe the following output:

Hexaly 12.5.20240604-Win64. All rights reserved.
Load car_sequencing.hxm...
Run input...
Run model...
Preprocess model 100% ...
Run param...
Run solver...
Initialize solver 100% ...
Push initial solutions 100% ...

Model:   expressions = 27001,   decisions = 10000  constraints = 520, objectives = 1,
Param:   time limit = 30 sec, no iteration limit

[ optimization direction  ]          minimize

[  0 sec,       0 itr]: No feasible solution found (infeas = 1040)
[  5 sec,  217424 itr]: obj    =           33
[ 10 sec,  384536 itr]: obj    =           24
[ 15 sec,  584611 itr]: obj    =           19
[ 20 sec,  783454 itr]: obj    =           15
[ 25 sec, 1043836 itr]: obj    =           12
[ 30 sec, 1310253 itr]: obj    =           10
[ 30 sec, 1310253 itr]: obj    =           10

1310253 iterations performed in 31 seconds

Feasible solution:
  obj    =           10
  gap    =      100.00%
  bounds =            0
Run output...

Every five seconds (in accordance with the parameter in function setTimeBetweenDisplays), a brief report of the current state of the search is displayed: the number of elapsed seconds and the number of iterations performed, and the values of the objective functions given in lexicographic order (separated with |).

You can disable all displays by calling the setVerbosity method with parameter 0.

Writing solutions

Here, we define a function writing the solution of the problem in a file. It writes the solution in the following format:

  1. First line = the objective value

  2. Second line = the classes at positions 1 to nbPositions

You can retrieve the value of any expression (decisions, intermediate expressions) thanks to the getValue function:

ofstream outfile;
outfile.exceptions(ofstream::failbit | ofstream::badbit);
outfile.open(fileName.c_str());

outfile << obj.getValue() << endl;
for (int p = 0; p < nbPositions; p++) {
    for (int c = 0; c < nbClasses; c++) {
        if (cp[c][p].getValue() == 1) {
            outfile << c << " ";
            break;
        }
    }
}
outfile << endl;