# AMPL Development Tutorial 2/6 – Stochastic Capacitated Facility Location Problem¶

Description: This notebook continues our six-part series as the second installment.
In the first part, we tackled a basic facility location problem, showcasing the use of AMPL and the `amplpy`

module within a Jupyter Notebook to derive a solution.
Our analysis revealed that for low to medium customer demand scenarios, two facilities are adequate.
However, a spike in demand necessitates the activation of a third facility.
In this section, we delve into the stochastic programming formulation of our problem.
This approach aims to provide a ‘robust’ solution that accommodates any potential customer demand scenario in the future.
Moving on to the third notebook of the series, we will transition our focus from modeling to the development of algorithms within AMPL.
We plan to examine four unique methods for implementing the Benders decomposition, starting with AMPL scripting before moving the control flow over to the `amplpy`

module.
This shift will demonstrate the parallel solving of subproblems.
Finally, we will introduce the AMPLS library, designed to enhance the efficiency of Benders decomposition by allowing AMPL model instances to be exported as persistent solver representations.

Tags: amplpy, ampl, mip, stochastic, facility location

Notebook author: Gyorgy Matyasfalvi <gyorgy@ampl.com>

References:

AMPL a Modeling Language for Mathematical Programming – Robert Fourer et al.

SCIP Optimization Suite example projects – Stephen J. Maher (ZIB)

```
# Install dependencies
%pip install -q amplpy
```

```
# Google Colab & Kaggle integration
from amplpy import AMPL, ampl_notebook
ampl = ampl_notebook(
modules=["open"], # modules to install
license_uuid="default", # license to use
) # instantiate AMPL object and register magics
```

## Problem description¶

Facility location decisions are crucial and often involve significant investment for both public and private sector entities, bearing profound social, economic, and environmental implications. The strategic positioning of facilities, such as warehouses, factories, and service centers, can determine an organization’s operational efficiency, market reach, and overall sustainability.

Given the high stakes of these decisions, engineers and analysts have developed sophisticated models to aid organizations in identifying optimal locations. These models take into account a variety of factors, including but not limited to, transportation costs, proximity to customers and suppliers, labor availability, customer demand, and environmental regulations.

The challenge is compounded when considering the uncertainty inherent in future conditions. Factors such as fluctuating market demands, changes in infrastructure, and unpredictable socio-economic developments require a robust approach to facility location. Hence, engineers often employ stochastic models and robust optimization techniques that account for such uncertainties, ensuring that the chosen locations remain viable under a range of possible future scenarios.

## Mixed integer program¶

Below you can find the extensive form of the stochastic facility location problem as an explicit mixed integer program.

**Given:**

A set of facilities: $I$.

A set of customers: $J$.

Set of scenarios: $S$ (representing different customer demands).

**Task:**

Find the minimum cost facilities to open such that the customer demand can be satisfied in all scenarios.

### Variables¶

$x_i \in {0, 1} \quad \forall i \in I$

$x_i = 1$ if facility $i$ is opened.

$y_{ij}^s \geq 0 \quad \forall i \in I, \forall j \in J, \forall s \in S$

$y_{ij}^s$ is the level of demand for customer $j$ satisfied by facility $i$ in scenario $s$.

### Parameters:¶

$\alpha^s$: the probability of scenario $s$.

$f_i$: the fixed cost for opening facility $i$,

$q_{ij}$: the cost of servicing customer $j$ from facility $i$,

$\lambda_j^s$: the demand of customer $j$ in scenario $s$,

$k_i:$ the capacity of facility $i$.

### The extensive form¶

The extensive form of our stochastic program can be formulated as follows:

$ \begin{equation} \begin{array}{rll} \min \quad & \sum_{i \in I} f_i x_i + \sum_{s \in S} \sum_{i \in I} \sum_{j \in J} \alpha^s q_{ij} y_{ij}^s & \ & & \ \textrm{subject to} \quad & \sum_{i \in I} y_{ij}^s \geq \lambda_j^s & \forall j \in J, \forall s \in S \ & \sum_{j \in J} y_{ij}^s \leq k_i x_i & \forall i \in I, \forall s \in S \ & \sum_{i \in I} k_i x_i \geq \max_{s \in S} \sum_{j \in J} \lambda_j^s & \ & & \ & x_i \in {0, 1} & \forall i \in I \ & y_{ij}^s \geq 0 & \forall i \in I, \forall j \in J, \forall s \in S \end{array} \tag{1} \end{equation} $

## AMPL Implementation¶

Translating the mathematical formulation of our optimization problem into an AMPL model is a direct process. The AMPL code closely mirrors each inequality in the system (1), preserving the structure of the mathematical model.

AMPL’s expressive syntax allows for meaningful names for entities such as variables, parameters, and constraints, enhancing the model’s clarity.
For instance, we’ll represent the set $I$ as `FACILITIES`

, $J$ as `CUSTOMERS`

, and $S$ as `SCENARIOS`

.
Variables previously denoted as $x$ and $y$ will be named `facility_open`

and `production`

, respectively.

Similarly, we will rename our parameters for greater clarity: $f_i$ becomes `fixed_cost`

, $q_ij$ is now `variable_cost`

, $\lambda_j^s$ is referred to as `customer_demand`

, and $k_i$ is labeled as `facility_capacity`

.

Using descriptive names not only enhances the readability of the model but also its maintainability and the ease with which it can be shared and understood by others.

### Changes to the simple deterministic model¶

In this section, we will see that the AMPL model and data files bear a strong resemblance to those we examined in the previous section. The key differences include:

The introduction of a new indexing set,

`SCENARIOS`

, to accommodate multiple scenarios.The expansion of the

`production`

variable to include indexing over`SCENARIOS`

, allowing scenario-specific production levels.The enhancement of the

`customer_demand`

parameter to also be indexed over`SCENARIOS`

, reflecting varying demand across different scenarios.The addition of a

`prob`

parameter, which assigns specific probabilities to each scenario, enabling stochastic analysis.Modifications to the objective function and constraints to integrate these new elements and reflect the stochastic nature of the model.

These adjustments pave the way for a more robust model that can account for uncertainty in customer demand across various scenarios.
Notice that we did not have to make any changes to our `facility_open`

variable and related parameters.

```
%%writefile floc_ef.mod
# Sets
set FACILITIES; # set of facilities
set CUSTOMERS; # set of customers
set SCENARIOS; # set of scenarios
# Variables
var facility_open{FACILITIES} binary; # 1 if facility i is open, 0 otherwise
var production{FACILITIES, CUSTOMERS, SCENARIOS} >= 0; # production from facility i to satisfy customer demand j in scenario s
# Parameters
param fixed_cost{FACILITIES} >= 0; # fixed cost of opening facility_open i
param variable_cost{FACILITIES, CUSTOMERS} >= 0; # variable cost of satisfying customer_demand of customer j from facility_open i
param customer_demand{CUSTOMERS, SCENARIOS} >= 0; # customer_demand of customer j in scenario s
param facility_capacity{FACILITIES} >= 0; # facility_capacity of facility_open i
param prob{SCENARIOS} default 1/card(SCENARIOS); # probability of scenario s
# Objective
minimize TotalCost:
sum{i in FACILITIES} fixed_cost[i] * facility_open[i] + # Fixed cost of opening facility i
sum{s in SCENARIOS, i in FACILITIES, j in CUSTOMERS} prob[s] * variable_cost[i,j] * production[i,j,s]; # Variable cost of satisfying customer demand j from facility i in scenario s
# Constraints
s.t. satisfying_customer_demand{s in SCENARIOS, j in CUSTOMERS}:
sum{i in FACILITIES} production[i,j,s] >= customer_demand[j,s];
s.t. facility_capacity_limits{s in SCENARIOS, i in FACILITIES}:
sum{j in CUSTOMERS} production[i,j,s] <= facility_capacity[i] * facility_open[i];
s.t. sufficient_production_capacity:
sum{i in FACILITIES} facility_capacity[i]*facility_open[i] >= max{s in SCENARIOS} sum{j in CUSTOMERS} customer_demand[j,s];
```

```
Writing floc_ef.mod
```

```
%%writefile floc_ef.dat
set FACILITIES := Baytown_TX Beaumont_TX Baton_Rouge_LA;
set CUSTOMERS := San_Antonio_TX Dallas_TX Jackson_MS Birmingham_AL;
set SCENARIOS := Low Medium High;
param prob := Low 0.25 Medium 0.5 High 0.25;
param fixed_cost := Baytown_TX 400000 Beaumont_TX 200000 Baton_Rouge_LA 600000;
param facility_capacity := Baytown_TX 1550 Beaumont_TX 650 Baton_Rouge_LA 1750;
param variable_cost: San_Antonio_TX Dallas_TX Jackson_MS Birmingham_AL :=
Baytown_TX 5739.725 6539.725 8650.40 22372.1125
Beaumont_TX 6055.05 6739.055 8050.40 21014.225
Baton_Rouge_LA 8650.40 7539.055 4539.72 15024.325;
param customer_demand: Low Medium High :=
San_Antonio_TX 450 650 887
Dallas_TX 910 1134 1456
Jackson_MS 379 416 673
Birmingham_AL 91 113 207;
```

```
Overwriting floc_ef.dat
```

### Changes in complexity¶

Despite only minor alterations to the model, its complexity has likely increased significantly due to the addition of more variables and constraints. To assess the extent of these changes, we can utilize a simple AMPL script, similar to what we’ve used previously:

```
print "Number of variables: ", _nvars;
print "Number of constraints: ", _ncons;
option show_stats 1;
```

This script will help us understand the scale of the model’s evolution by providing the current counts of variables and constraints.

```
ampl = AMPL() # Instantiate an AMPL object
ampl.read("floc_ef.mod") # Read the model from file
ampl.read_data("floc_ef.dat") # Read the data from file
ampl.eval(
r"""
print "Number of variables: ", _nvars;
print "Number of constraints: ", _ncons;
option show_stats 1;
"""
) # Display total number of variables and constraints
ampl.option["solver"] = "highs" # Select the solver
ampl.option["display_precision"] = 0 # Turn off display rounding
ampl.solve() # Attempt to solve
print(
ampl.option["solve_result_table"]
) # Print the solve result table, this will inform us of the various solution codes.
# Check if the problem was solved if not print warning
srn = ampl.get_value("solve_result_num")
if srn != 0:
print(f"Warning:\tProblem not solved to optimality!\n\t\tsolve_result_num = {srn}")
else:
print("Success:\tProblem solved to optimality!")
# Print the solution
ampl.eval(
r"""
print;
display TotalCost;
display facility_open;
"""
)
```

```
Number of variables: 39
Number of constraints: 22
Presolve eliminates 1 constraint and 2 variables.
Adjusted problem:
37 variables:
1 binary variable
36 linear variables
21 constraints, all linear; 75 nonzeros
21 inequality constraints
1 linear objective; 37 nonzeros.
HiGHS 1.6.0:HiGHS 1.6.0: optimal solution; objective 16758018.6
18 simplex iterations
1 branching nodes
0 solved
100 solved?
200 infeasible
300 unbounded
400 limit
500 failure
Success: Problem solved to optimality!
TotalCost = 16758018.596250001
facility_open [*] :=
Baton_Rouge_LA 1
Baytown_TX 1
Beaumont_TX 1
;
```

## Conclusion¶

Our model, designed to meet demand across all scenarios, necessitates that all three facilities be operational in the stochastic programming solution.
Consequently, the `facility_open`

variable aligns with the high-demand solution of the simple deterministic model.
Moreover, the stochastic approach yields an estimate of the expected operational costs, offering a more accurate projection than the deterministic model’s estimate, which was **22,250,711.2**.

Transitioning from the simple deterministic model to the extensive form of a stochastic programming model proved to be relatively straightforward.

### Problem complexity¶

However, the transition has resulted in an increase in complexity: **the number of variables has expanded from 15 to 39, and the number of constraints has risen from 8 to 22**.
While this change is manageable for a small-scale problem, the complexity could grow exponentially with larger problems and additional scenarios.

### Dealing with complexity — decomposition¶

In the face of such complexity, what strategies are available? Fortunately, for two-stage stochastic programming problems like ours, a solution exists. The first-stage decisions involve determining which facilities to open, while the second-stage decisions focus on how to meet customer demand with the available facilities. This inherent structure allows for the problem to be broken down into a master problem and sub-problems, significantly reducing the complexity.

In the following section, we will explore the implementation of this decomposition strategy, known as Benders decomposition, within AMPL.