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

2_stoch_floc.ipynb Open In Colab Kaggle Gradient Open In SageMaker Studio Lab Hits

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.

previous

AMPL Development Tutorial 1/6 – Capacitated Facility Location Problem

next

AMPL Development Tutorial 3/6 – Benders Decomposition via AMPL scripting