Techno Blender
Digitally Yours.

Solving Two-Stage Stochastic Programs in Gurobi | by Nakul Upadhya | Oct, 2022

0 39


Formulating and solving a two-stage stochastic server farm problem

Stochastic programming (SP) is a framework for modeling optimization problems that involve uncertainty [1]. In many cases, SP models take the form of a two-stage problem. The first stage involves finding the optimal deterministic decisions. These decisions are based on information we know to be certain (AKA the here-and-now decisions). The second stage involves making decisions that rely on randomness (also called the recourse decision), given our first-stage decisions. The optimization problem aims to minimize the loss (or maximize the profit) resulting from the first stage decision plus the expected loss resulting from the second stage decision. Mathematically, this can be written in the following format:

In this article, I will demonstrate how to formulate and solve one of these problems using Gurobi.

Let’s put this framework into perspective with an example: a server farm problem. Let’s say we are trying to design a server farm for our company, and we need to install CPU, GPU, and TPU cores to help handle the computing resources of our company given the amount of space we have and our upfront budget for buying servers. We have a general estimation of the next-month demand for each resource and the distributions they follow. For the sake of simplicity, we will ignore future months after next month (though it is possible to do a multi-stage stochastic program to handle this).

Running the servers in the future will cost a small amount, but if we don’t meet the demand we will have to leverage cloud computing resources which will significantly increase our operational costs.
In this example, our first-stage decision (the here-and-now decision) is the number of each resource to buy. The second-stage (recourse) decision is the number of cloud-compute to buy. Our goal is to minimize the operational cost of the server farm. This optimization problem can be written as follows:

Where x represents how much of each core we will buy, y represents how much cloud compute we will leverage, i is our installation costs, B is our installation budget, s is the size of the cores, S is our available space, c is our cloud computing costs, and d(xi) is our random demand.

To solve the problem above, we will implement a technique called sample average approximation [2]. This method involves generating a large number (K) of potential realizations (scenarios) of the random variable through Monte Carlo sampling and treating each realization as an individual sub-problem [2]. This helps discretize the randomness and therefore makes it solvable without losing too much fidelity of information. This means that instead of just having 1 set of second-stage variables, we will have K variables and the expectation of our second-stage value is the average value of our realizations. With this, our formulation changes to:

While this formulation can then be implemented in a variety of programming languages and linear programming solvers (e.x. SCIP, CPLEX, CVXPY), we will focus on solving this problem using Gurobi — a commercial-grade, state-of-the-art mathematical programming solver. Gurobi has APIs for languages like Python, MATLAB, R, and Java, but most of its users prefer using the Python API [4], so we will use it here.

First, we will import our packages (numpy and gurobipy)as well as initialize (make up) the data for our problem.

In this problem i is the investment cost, sis the space each item takes up, Bis our budget, and Sis our maximum space. Another variable of note is the demand in each scenario d_xi. We can then use Gurobi to solve our problem. First, we will define our model with gp.model(), set the goal (minimize) as well as initialize our decision variables.

As you can see in the above examples, both our demand (d_xi) and our second-stage decision variable (y) are not single-dimensional vectors, but rather a matrix since both of these have a scenario index in our formulation. After this we can define our objective function:

One nice thing about Gurobi is that it supports numpy matrix multiplication operations when defining constraints and objectives, making it simple to code if you think in matrix multiplications (ex. i @ x). Gurobi also supports summation operations, making the barrier to entry lower as well. Here we calculate the first stage cost, then take the average second stage cost.

We can now move on to defining our constraints:

The first line of code ensures that the cost of installing in-house servers does not exceed our budget B. The second line ensures that we only buy as many servers as we have space for. The third line of code ensures the number of installed servers and amount of cloud compute resources in scenario k leveraged are enough to meet the demand in scenario k.

Now that the problem is fully defined, we can simply solve the problem and get our optimal decision by calling m.optimize():

The optimal variable values are then stored in each variable we defined and can be accessed via variable_name.x. While we could have gotten the optimal values for both the x and y variables, we only care about the decisions we can make now (x).

This approach gives us a high-quality solution, but there are drawbacks.

The first downside is that as the number of scenarios increases, the model will get exponentially harder to solve. For example, we only had 6 true decisions we needed to make, but using sample average approximation made our model make 153 decisions (3 first-stage, 3 second-stage for 50 scenarios). This means that decomposition algorithms (like Benders Decomposition for Linear SPs [3]) often need to be used to parallelize the computation.

Another consideration is that unfortunately, Gurobi is not free software. While Gurobi does offer an academic license that many students and researchers leverage, individuals not affiliated with a university may find it difficult to gain access to this software. One open-source alternative to Gurobi is SCIP [5]. Both share a similar modeling language so Gurobi users can easily adapt to SCIP and vice-versa. I chose to demonstrate this concept with Gurobi instead of SCIP due to personal preference along with Gurobi’s fast solver performance. That doesn’t mean that SCIP is inferior. In fact, many optimization researchers who focus on developing or improving optimization algorithms prefer SCIP since its open-source nature allows users to tweak every part of the solver allowing them to try new methodologies.

The full code for this problem can be found here.

Works Cited

[1] A. Philpott, What is Stochastic Programming. (2022), Stochastic Programming Society.

[2] A. J. Kleywegt, A. Shapiro, T, Homem-de-Mello. The Sample Average Approximation Method for Stochastic Discrete Optimization (2002), SIAM Journal on Optimization.

[3] S. S. Nielsen, S. A. Zenios. Scalable parallel Benders Decomposition for Stochastic linear programming (1997), Parallel Computing Volume 23, Issue 8.

[4] Gurobi. Starting with Gurobi (2022), Gurobi Website

[5] Bestuzheva et. al. The SCIP Optimization Suite 8.0 (2021), Optimization Online


Formulating and solving a two-stage stochastic server farm problem

Stochastic programming (SP) is a framework for modeling optimization problems that involve uncertainty [1]. In many cases, SP models take the form of a two-stage problem. The first stage involves finding the optimal deterministic decisions. These decisions are based on information we know to be certain (AKA the here-and-now decisions). The second stage involves making decisions that rely on randomness (also called the recourse decision), given our first-stage decisions. The optimization problem aims to minimize the loss (or maximize the profit) resulting from the first stage decision plus the expected loss resulting from the second stage decision. Mathematically, this can be written in the following format:

In this article, I will demonstrate how to formulate and solve one of these problems using Gurobi.

Let’s put this framework into perspective with an example: a server farm problem. Let’s say we are trying to design a server farm for our company, and we need to install CPU, GPU, and TPU cores to help handle the computing resources of our company given the amount of space we have and our upfront budget for buying servers. We have a general estimation of the next-month demand for each resource and the distributions they follow. For the sake of simplicity, we will ignore future months after next month (though it is possible to do a multi-stage stochastic program to handle this).

Running the servers in the future will cost a small amount, but if we don’t meet the demand we will have to leverage cloud computing resources which will significantly increase our operational costs.
In this example, our first-stage decision (the here-and-now decision) is the number of each resource to buy. The second-stage (recourse) decision is the number of cloud-compute to buy. Our goal is to minimize the operational cost of the server farm. This optimization problem can be written as follows:

Where x represents how much of each core we will buy, y represents how much cloud compute we will leverage, i is our installation costs, B is our installation budget, s is the size of the cores, S is our available space, c is our cloud computing costs, and d(xi) is our random demand.

To solve the problem above, we will implement a technique called sample average approximation [2]. This method involves generating a large number (K) of potential realizations (scenarios) of the random variable through Monte Carlo sampling and treating each realization as an individual sub-problem [2]. This helps discretize the randomness and therefore makes it solvable without losing too much fidelity of information. This means that instead of just having 1 set of second-stage variables, we will have K variables and the expectation of our second-stage value is the average value of our realizations. With this, our formulation changes to:

While this formulation can then be implemented in a variety of programming languages and linear programming solvers (e.x. SCIP, CPLEX, CVXPY), we will focus on solving this problem using Gurobi — a commercial-grade, state-of-the-art mathematical programming solver. Gurobi has APIs for languages like Python, MATLAB, R, and Java, but most of its users prefer using the Python API [4], so we will use it here.

First, we will import our packages (numpy and gurobipy)as well as initialize (make up) the data for our problem.

In this problem i is the investment cost, sis the space each item takes up, Bis our budget, and Sis our maximum space. Another variable of note is the demand in each scenario d_xi. We can then use Gurobi to solve our problem. First, we will define our model with gp.model(), set the goal (minimize) as well as initialize our decision variables.

As you can see in the above examples, both our demand (d_xi) and our second-stage decision variable (y) are not single-dimensional vectors, but rather a matrix since both of these have a scenario index in our formulation. After this we can define our objective function:

One nice thing about Gurobi is that it supports numpy matrix multiplication operations when defining constraints and objectives, making it simple to code if you think in matrix multiplications (ex. i @ x). Gurobi also supports summation operations, making the barrier to entry lower as well. Here we calculate the first stage cost, then take the average second stage cost.

We can now move on to defining our constraints:

The first line of code ensures that the cost of installing in-house servers does not exceed our budget B. The second line ensures that we only buy as many servers as we have space for. The third line of code ensures the number of installed servers and amount of cloud compute resources in scenario k leveraged are enough to meet the demand in scenario k.

Now that the problem is fully defined, we can simply solve the problem and get our optimal decision by calling m.optimize():

The optimal variable values are then stored in each variable we defined and can be accessed via variable_name.x. While we could have gotten the optimal values for both the x and y variables, we only care about the decisions we can make now (x).

This approach gives us a high-quality solution, but there are drawbacks.

The first downside is that as the number of scenarios increases, the model will get exponentially harder to solve. For example, we only had 6 true decisions we needed to make, but using sample average approximation made our model make 153 decisions (3 first-stage, 3 second-stage for 50 scenarios). This means that decomposition algorithms (like Benders Decomposition for Linear SPs [3]) often need to be used to parallelize the computation.

Another consideration is that unfortunately, Gurobi is not free software. While Gurobi does offer an academic license that many students and researchers leverage, individuals not affiliated with a university may find it difficult to gain access to this software. One open-source alternative to Gurobi is SCIP [5]. Both share a similar modeling language so Gurobi users can easily adapt to SCIP and vice-versa. I chose to demonstrate this concept with Gurobi instead of SCIP due to personal preference along with Gurobi’s fast solver performance. That doesn’t mean that SCIP is inferior. In fact, many optimization researchers who focus on developing or improving optimization algorithms prefer SCIP since its open-source nature allows users to tweak every part of the solver allowing them to try new methodologies.

The full code for this problem can be found here.

Works Cited

[1] A. Philpott, What is Stochastic Programming. (2022), Stochastic Programming Society.

[2] A. J. Kleywegt, A. Shapiro, T, Homem-de-Mello. The Sample Average Approximation Method for Stochastic Discrete Optimization (2002), SIAM Journal on Optimization.

[3] S. S. Nielsen, S. A. Zenios. Scalable parallel Benders Decomposition for Stochastic linear programming (1997), Parallel Computing Volume 23, Issue 8.

[4] Gurobi. Starting with Gurobi (2022), Gurobi Website

[5] Bestuzheva et. al. The SCIP Optimization Suite 8.0 (2021), Optimization Online

FOLLOW US ON GOOGLE NEWS

Read original article here

Denial of responsibility! Techno Blender is an automatic aggregator of the all world’s media. In each content, the hyperlink to the primary source is specified. All trademarks belong to their rightful owners, all materials to their authors. If you are the owner of the content and do not want us to publish your materials, please contact us by email – [email protected]. The content will be deleted within 24 hours.

Leave a comment