Stochastic Packing¶
Principles learned¶
Use list decision variables
Use a lambda expression to sum over a set
Use the sort operator
Problem¶
The stochastic packing problem is defined as a set of items that need to be grouped into bins. Each bin can contain any item, but the items have random weights. Here, randomness is represented by scenarios, where one scenario defines the weights of all items. Given a distribution of items in bins, each scenario yields a bin of maximum weight.
Now, which objective is most appropriate to minimize this stochastic maximum weight?
Minimizing the average can hide risky scenarios while minimizing the worst-case might be too pessimistic. A usual compromise to build a robust scenario is to optimize on a given percentile. It is what is done, minimizing the 90 th percentile of makespans. Thanks to the sort operator, such a nonlinear criterion is straightforward to model using Hexaly Optimizer.
Download the exampleData¶
For this example, we generate instances at runtime: first a uniform distribution is picked for each item, then for each scenario, the weight of each item is independently sampled from the corresponding uniform distribution.
Program¶
We use a set decision variable to represent the set of items contained in a bin. Those sets are constrained to form a partition as each item must be present in exactly one bin.
We then compute and store in the scenarioMaxWeight
array the maximum
weight corresponding to each scenario. To do so, we first need to compute the
total weight of each bin as a sum over a set and therefore need to
define a lambda expression.
We can then sort the array of maximum weights and access our objective function: its 9 th decile.
- Execution:
- hexaly stochastic_packing.hxm [hxTimeLimit] [solFileName=solution.txt]
use random;
/* Generate instance data */
function input() {
nbItems = 10;
nbBins = 2;
nbScenarios = 3;
rngSeed = 42;
// Pick random parameters for each item distribution
rng = random.create(rngSeed);
itemsMin[i in 0...nbItems] = rng.next(10, 101);
itemsMax[i in 0...nbItems] = itemsMin[i] + rng.next(0, 51);
// Sample the distributions to generate the scenarios
scenarioItemWeights[i in 0...nbScenarios][j in 0...nbItems] =
rng.next(itemsMin[j], itemsMax[j] + 1);
}
/* Declare the optimization model */
function model() {
// Set decisions: bins[k] represents the items in bin k
bins[k in 0...nbBins] <- set(nbItems);
// Each item must be in one bin and one bin only
constraint partition[k in 0...nbBins](bins[k]);
// Compute max weight for each scenario
scenarioMaxWeight[m in 0...nbScenarios] <- max[k in 0...nbBins](
sum(bins[k], i => scenarioItemWeights[m][i]));
// Compute the 9th decile of scenario max weights
stochasticMaxWeight <- sort(scenarioMaxWeight)[ceil(0.9 * (nbScenarios - 1))];
minimize stochasticMaxWeight;
}
// Parametrize the solver
function param() {
if (hxTimeLimit == nil) hxTimeLimit = 2;
}
/* Write the solution */
function output() {
println();
println("Scenario item weights:");
for [i in 0...nbScenarios] {
print(i + ": [");
for [j in 0...nbItems]
print(scenarioItemWeights[i][j] + (j == nbItems - 1 ? "" : ", "));
println("]");
}
println();
println("Bins:");
for [k in 0...nbBins]
println(k + ": " + bins[k].value);
}
- Execution (Windows)
- set PYTHONPATH=%HX_HOME%\bin\pythonpython stochastic_packing.py
- Execution (Linux)
- export PYTHONPATH=/opt/hexaly_13_0/bin/pythonpython stochastic_packing.py
from __future__ import print_function
import random
import math
import hexaly.optimizer
def generate_scenarios(nb_items, nb_scenarios, rng_seed):
random.seed(rng_seed)
# Pick random parameters for each item distribution
items_dist = []
for _ in range(nb_items):
item_min = random.randint(10, 100)
item_max = item_min + random.randint(0, 50)
items_dist.append((item_min, item_max))
# Sample the distributions to generate the scenarios
scenario_item_weights = [[random.randint(*dist) for dist in items_dist]
for _ in range(nb_scenarios)]
return scenario_item_weights
def main(nb_items, nb_bins, nb_scenarios, seed, time_limit):
# Generate instance data
scenario_item_weights_data = generate_scenarios(nb_items, nb_scenarios, seed)
with hexaly.optimizer.HexalyOptimizer() as optimizer:
#
# Declare the optimization model
#
model = optimizer.model
# Set decisions: bins[k] represents the items in bin k
bins = [model.set(nb_items) for _ in range(nb_bins)]
# Each item must be in one bin and one bin only
model.constraint(model.partition(bins))
scenarios_item_weights = model.array(scenario_item_weights_data)
# Compute max weight for each scenario
scenarios_max_weights = model.array(
model.max(
model.sum(bin,
model.lambda_function(
lambda i:
model.at(scenarios_item_weights, k, i)))
for bin in bins) for k in range(nb_scenarios))
# Compute the 9th decile of scenario max weights
stochastic_max_weight = \
model.sort(scenarios_max_weights)[int(math.ceil(0.9 * (nb_scenarios - 1)))]
model.minimize(stochastic_max_weight)
model.close()
# Parameterize the optimizer
optimizer.param.time_limit = time_limit
optimizer.solve()
#
# Write the solution
#
print()
print("Scenario item weights:")
for i, scenario in enumerate(scenario_item_weights_data):
print(i, ': ', scenario, sep='')
print()
print("Bins:")
for k, bin in enumerate(bins):
print(k, ': ', bin.value, sep='')
if __name__ == '__main__':
nb_items = 10
nb_bins = 2
nb_scenarios = 3
rng_seed = 42
time_limit = 2
main(
nb_items,
nb_bins,
nb_scenarios,
rng_seed,
time_limit
)
- Compilation / Execution (Windows)
- cl /EHsc stochastic_packing.cpp -I%HX_HOME%\include /link %HX_HOME%\bin\hexaly130.libstochastic_packing
- Compilation / Execution (Linux)
- g++ stochastic_packing.cpp -I/opt/hexaly_13_0/include -lhexaly130 -lpthread -o stochastic_packingstochastic_packing
#include "optimizer/hexalyoptimizer.h"
#include <cmath>
#include <iostream>
#include <random>
#include <vector>
using namespace hexaly;
class StochasticPacking {
private:
// Number of items
int nbItems;
// Number of bins
int nbBins;
// Number of scenarios
int nbScenarios;
// For each scenario, the weight of each item
std::vector<std::vector<int>> scenarioItemWeights;
// Hexaly Optimizer
HexalyOptimizer optimizer;
// Decision variable for the assignment of items
std::vector<HxExpression> bins;
// For each scenario, the corresponding maximum weight
std::vector<HxExpression> scenarioMaxWeight;
// Objective = minimize the 9th decile of all possible max weights
HxExpression stochasticMaxWeight;
void generateScenarios(unsigned int rngSeed) {
std::mt19937 rng(rngSeed);
std::uniform_int_distribution<int> distMin(10, 100);
std::uniform_int_distribution<int> distDelta(0, 50);
// Pick random parameters for each item distribution
std::vector<std::uniform_int_distribution<int>> itemsDists;
for (int i = 0; i < nbItems; ++i) {
int min = distMin(rng);
int max = min + distDelta(rng);
itemsDists.emplace_back(min, max);
}
// Sample the distributions to generate the scenarios
for (int i = 0; i < nbScenarios; ++i) {
for (int j = 0; j < nbItems; ++j) {
scenarioItemWeights[i][j] = itemsDists[j](rng);
}
}
}
public:
StochasticPacking(int nbItems, int nbBins, int nbScenarios, unsigned int seed)
: nbItems(nbItems), nbBins(nbBins), nbScenarios(nbScenarios),
scenarioItemWeights(nbScenarios, std::vector<int>(nbItems)), optimizer() {
generateScenarios(seed);
}
void solve(int timeLimit) {
// Declare the optimization model
HxModel model = optimizer.getModel();
bins.resize(nbBins);
scenarioMaxWeight.resize(nbScenarios);
// Set decisions: bins[k] represents the items in bin k
for (int k = 0; k < nbBins; ++k) {
bins[k] = model.setVar(nbItems);
}
// Each item must be in one bin and one bin only
model.constraint(model.partition(bins.begin(), bins.end()));
// Compute max weight for each scenario
for (int m = 0; m < nbScenarios; ++m) {
HxExpression scenario = model.array(scenarioItemWeights[m].begin(), scenarioItemWeights[m].end());
HxExpression weightLambda = model.createLambdaFunction([&](HxExpression i) { return scenario[i]; });
std::vector<HxExpression> binWeights(nbBins);
for (int k = 0; k < nbBins; ++k) {
binWeights[k] = model.sum(bins[k], weightLambda);
}
scenarioMaxWeight[m] = model.max(binWeights.begin(), binWeights.end());
}
// Compute the 9th decile of scenario max weights
HxExpression scenarioMaxWeightArray = model.array(scenarioMaxWeight.begin(), scenarioMaxWeight.end());
HxExpression sortedScenarioMaxWeight = model.sort(scenarioMaxWeightArray);
stochasticMaxWeight = sortedScenarioMaxWeight[(int)std::ceil(0.9 * (nbScenarios - 1))];
model.minimize(stochasticMaxWeight);
model.close();
// Parametrize the optimizer
optimizer.getParam().setTimeLimit(timeLimit);
optimizer.solve();
}
/* Write the solution */
void writeSolution(std::ostream& os) const {
os << "\nScenario item weights:\n";
for (int i = 0; i < nbScenarios; ++i) {
os << i << ": [";
for (int j = 0; j < scenarioItemWeights[i].size(); ++j) {
os << scenarioItemWeights[i][j] << (j == scenarioItemWeights[i].size() - 1 ? "" : ", ");
}
os << "]\n";
}
os << "\nBins:\n";
for (int m = 0; m < nbBins; ++m) {
os << m << ": { ";
HxCollection items = bins[m].getCollectionValue();
for (int i = 0; i < items.count(); ++i) {
os << items[i] << (i == items.count() - 1 ? " " : ", ");
}
os << "}\n";
}
}
};
int main(int argc, char** argv) {
int nbItems = 10;
int nbBins = 2;
int nbScenarios = 3;
int rngSeed = 42;
int timeLimit = 2;
try {
StochasticPacking model(nbItems, nbBins, nbScenarios, rngSeed);
model.solve(timeLimit);
model.writeSolution(std::cout);
return 0;
} catch (const std::exception& e) {
std::cerr << "An error occurred: " << e.what() << std::endl;
return 1;
}
}
- Compilation/Execution (Windows)
- copy %HX_HOME%\bin\Hexaly.NET.dll .csc StochasticPacking.cs /reference:Hexaly.NET.dllStochasticPacking
using System;
using Hexaly.Optimizer;
public class StochasticPacking : IDisposable
{
// Number of items
int nbItems;
// Number of bins
int nbBins;
// Number of scenarios
int nbScenarios;
// For each scenario, the weight of each item
int[][] scenarioItemWeights;
// Hexaly Optimizer
HexalyOptimizer optimizer;
// Decision variable for the assignment of items
HxExpression[] bins;
// For each scenario, the corresponding maximum weight
HxExpression[] scenarioMaxWeight;
// Objective = minimize the 9th decile of all possible max weights
HxExpression stochasticMaxWeight;
private void generateScenarios(int rngSeed)
{
Random rng = new Random(rngSeed);
// Pick random parameters for each item distribution
int[] itemsMin = new int[nbItems];
int[] itemsMax = new int[nbItems];
for (int i = 0; i < nbItems; ++i)
{
itemsMin[i] = rng.Next(10, 101);
itemsMax[i] = itemsMin[i] + rng.Next(51);
}
// Sample the distributions to generate the scenarios
scenarioItemWeights = new int[nbScenarios][];
for (int i = 0; i < nbScenarios; ++i)
{
scenarioItemWeights[i] = new int[nbItems];
for (int j = 0; j < nbItems; ++j)
scenarioItemWeights[i][j] = rng.Next(itemsMin[i], itemsMax[i] + 1);
}
}
public StochasticPacking(int nbItems, int nbBins, int nbScenarios, int rngSeed)
{
optimizer = new HexalyOptimizer();
this.nbItems = nbItems;
this.nbBins = nbBins;
this.nbScenarios = nbScenarios;
generateScenarios(rngSeed);
}
public void Dispose()
{
if (optimizer != null)
optimizer.Dispose();
}
void Solve(int limit)
{
// Declare the optimization model
HxModel model = optimizer.GetModel();
bins = new HxExpression[nbBins];
scenarioMaxWeight = new HxExpression[nbScenarios];
// Set decisions: bins[k] represents the items in bin k
for (int k = 0; k < nbBins; ++k)
bins[k] = model.Set(nbItems);
// Each item must be in one bin and one bin only
model.Constraint(model.Partition(bins));
// Compute max weight for each scenario
for (int m = 0; m < nbScenarios; ++m)
{
HxExpression scenario = model.Array(scenarioItemWeights[m]);
HxExpression weightLambda = model.LambdaFunction(i => scenario[i]);
HxExpression[] binWeights = new HxExpression[nbBins];
for (int k = 0; k < nbBins; ++k)
binWeights[k] = model.Sum(bins[k], weightLambda);
scenarioMaxWeight[m] = model.Max(binWeights);
}
// Compute the 9th decile of scenario max weights
HxExpression scenarioMaxWeightArray = model.Array(scenarioMaxWeight);
HxExpression sortedScenarioMaxWeight = model.Sort(scenarioMaxWeightArray);
stochasticMaxWeight = sortedScenarioMaxWeight[(int)Math.Ceiling(0.9 * (nbScenarios - 1))];
model.Minimize(stochasticMaxWeight);
model.Close();
// Parametrize the optimizer
optimizer.GetParam().SetTimeLimit(limit);
optimizer.Solve();
}
/* Write the solution */
private void WriteSolution()
{
Console.WriteLine();
Console.WriteLine("Scenario item weights:");
for (int i = 0; i < nbScenarios; ++i)
{
Console.Write(i + ": [");
for (int j = 0; j < nbItems; ++j)
Console.Write(scenarioItemWeights[i][j] + (j == nbItems - 1 ? "" : ", "));
Console.WriteLine("]");
}
Console.WriteLine();
Console.WriteLine("Bins:");
for (int m = 0; m < nbBins; ++m)
{
Console.Write(m + ": { ");
HxCollection items = bins[m].GetCollectionValue();
for (int i = 0; i < items.Count(); ++i)
Console.Write(items.Get(i) + (i == items.Count() - 1 ? " " : ", "));
Console.WriteLine("}");
}
}
public static void Main(string[] args)
{
int nbItems = 10;
int nbBins = 2;
int nbScenarios = 3;
int rngSeed = 43;
int timeLimit = 2;
using (
StochasticPacking model = new StochasticPacking(
nbItems,
nbBins,
nbScenarios,
rngSeed
)
)
{
model.Solve(timeLimit);
model.WriteSolution();
}
}
}
- Compilation / Execution (Windows)
- javac StochasticPacking.java -cp %HX_HOME%\bin\hexaly.jarjava -cp %HX_HOME%\bin\hexaly.jar;. StochasticPacking
- Compilation/Execution (Linux)
- javac StochasticPacking.java -cp /opt/hexaly_13_0/bin/hexaly.jarjava -cp /opt/hexaly_13_0/bin/hexaly.jar:. StochasticPacking
import java.util.Random;
import com.hexaly.optimizer.*;
public class StochasticPacking {
// Number of items
private int nbItems;
// Number of bins
private int nbBins;
// Number of scenarios
private int nbScenarios;
// For each scenario, the weight of each item
private int[][] scenarioItemWeights;
// Hexaly Optimizer
private final HexalyOptimizer optimizer;
// Decision variable for the assignment of items
private HxExpression[] bins;
// For each scenario, the corresponding max weight
private HxExpression[] scenarioMaxWeight;
// Objective = minimize the 9th decile of all possible max weights
private HxExpression stochasticMaxWeight;
private void generateScenarios(int rngSeed) {
Random rng = new Random(rngSeed);
// Pick random parameters for each item distribution
int[] itemsMin = new int[nbItems];
int[] itemsMax = new int[nbItems];
for (int i = 0; i < nbItems; ++i) {
itemsMin[i] = 10 + rng.nextInt(91);
itemsMax[i] = itemsMin[i] + rng.nextInt(51);
}
// Sample the distributions to generate the scenarios
scenarioItemWeights = new int[nbScenarios][nbItems];
for (int i = 0; i < nbScenarios; ++i) {
for (int j = 0; j < nbItems; ++j) {
scenarioItemWeights[i][j] = itemsMin[j] + rng.nextInt(itemsMax[i] - itemsMin[i] + 1);
}
}
}
private StochasticPacking(HexalyOptimizer optimizer, int nbItems, int nbBins, int nbScenarios, int rngSeed) {
this.optimizer = optimizer;
this.nbItems = nbItems;
this.nbBins = nbBins;
this.nbScenarios = nbScenarios;
generateScenarios(rngSeed);
}
private void solve(int limit) {
// Declare the optimization model
HxModel model = optimizer.getModel();
bins = new HxExpression[nbBins];
scenarioMaxWeight = new HxExpression[nbScenarios];
// Set decisions: bins[k] represents the items in bin k
for (int k = 0; k < nbBins; ++k) {
bins[k] = model.setVar(nbItems);
}
// Each item must be in one bin and one bin only
model.constraint(model.partition(bins));
// Compute max weight for each scenario
for (int m = 0; m < nbScenarios; ++m) {
HxExpression scenario = model.array(scenarioItemWeights[m]);
HxExpression weightLambda = model.lambdaFunction(i -> model.at(scenario, i));
HxExpression[] binWeights = new HxExpression[nbBins];
for (int k = 0; k < nbBins; ++k) {
binWeights[k] = model.sum(bins[k], weightLambda);
}
scenarioMaxWeight[m] = model.max(binWeights);
}
// Compute the 9th decile of scenario makespans
HxExpression scenarioMaxWeightArray = model.array(scenarioMaxWeight);
HxExpression sortedScenarioMaxWeight = model.sort(scenarioMaxWeightArray);
stochasticMaxWeight = model.at(sortedScenarioMaxWeight, (int) Math.ceil(0.9 * (nbScenarios - 1)));
model.minimize(stochasticMaxWeight);
model.close();
// Parametrize the optimizer
optimizer.getParam().setTimeLimit(limit);
optimizer.solve();
}
/* Write the solution */
private void writeSolution() {
System.out.println();
System.out.println("Scenario item weights:");
for (int i = 0; i < nbScenarios; ++i) {
System.out.print("" + i + ": [");
for (int j = 0; j < nbItems; ++j) {
System.out.print("" + scenarioItemWeights[i][j] + (j == nbItems - 1 ? "" : ", "));
}
System.out.println("]");
}
System.out.println();
System.out.println("Bins:");
for (int m = 0; m < nbBins; ++m) {
System.out.print("" + m + ": { ");
HxCollection items = bins[m].getCollectionValue();
for (int i = 0; i < items.count(); ++i) {
System.out.print("" + items.get(i) + (i == items.count() - 1 ? " " : ", "));
}
System.out.println("}");
}
}
public static void main(String[] args) {
try (HexalyOptimizer optimizer = new HexalyOptimizer()) {
int nbItems = 10;
int nbBins = 2;
int nbScenarios = 3;
int rngSeed = 42;
int timeLimit = 2;
StochasticPacking model = new StochasticPacking(optimizer, nbItems, nbBins, nbScenarios,
rngSeed);
model.solve(timeLimit);
model.writeSolution();
} catch (Exception ex) {
System.err.println(ex);
ex.printStackTrace();
System.exit(1);
}
}
};