# Benders decomposition for congested partial set covering location with
uncertain demand

Alice Calamita1 Ivana Ljubić2 Laura Palagi1

(1Department of Computer, Control and Management Engineering, Sapienza
University of Rome, Via Ariosto 25, 00185 Rome, Italy
{alice.calamita<EMAIL_ADDRESS>
2Department of Information Systems, Decision Sciences and Statistics, ESSEC
Business School of Paris, 3 Avenue Bernard Hirsch, 95000 Cergy, France
<EMAIL_ADDRESS>
###### Abstract

In this paper, we introduce a mixed integer quadratic formulation for the
congested variant of the partial set covering location problem, which involves
determining a subset of facility locations to open and efficiently allocating
customers to these facilities to minimize the combined costs of facility
opening and congestion while ensuring target coverage. To enhance the
resilience of the solution against demand fluctuations, we address the case
under uncertain customer demand using $\Gamma$-robustness. We formulate the
deterministic problem and its robust counterpart as mixed-integer quadratic
problems. We investigate the effect of the protection level in adapted
instances from the literature to provide critical insights into how sensitive
the planning is to the protection level. Moreover, since the size of the
robust counterpart grows with the number of customers, which could be
significant in real-world contexts, we propose the use of Benders
decomposition to effectively reduce the number of variables by projecting out
of the master problem all the variables dependent on the number of customers.
We illustrate how to incorporate our Benders approach within a mixed-integer
second-order cone programming (MISOCP) solver, addressing explicitly all the
ingredients that are instrumental for its success. We discuss single-tree and
multi-tree approaches and introduce a perturbation technique to deal with the
degeneracy of the Benders subproblem efficiently. Our tailored Benders
approaches outperform the perspective reformulation solved using the state-of-
the-art MISOCP solver Gurobi on adapted instances from the literature.

##  1 Introduction

The partial set covering location problem (PSCLP) belongs to the class of
facility location problems that incorporate a notion of coverage. Coverage is
defined by a proximity measure, often the distance or the travel time,
establishing whether a potential facility location can serve or cover a
certain demand point. The PSCLP aims to minimize the cost of opening the
facilities while ensuring the coverage of a specified (partial) amount of the
total demand. It derives from the renowned set covering location problem
(SCLP), which aims to locate a minimum-cost set of facilities such that all
demand points are covered at least once. The PSCLP emerged as a response to
the often expensive or impractical solutions produced by the SCLP, allowing
for partial coverage. Also PSCLP has some drawbacks: it does not take into
account the congestion at the facilities, which may arise from a sudden
increase in the customers’ demand. Studying congestion is crucial as it
directly impacts the performance (i.e., the quality of service) and efficiency
of this type of networks, and sometimes, depending on the application, it
could lead to additional cost due to diseconomies of scale (due to, e.g., the
employment of additional overtime workers, the use of more expensive
materials, or by neglecting or postponing equipment maintenance schedules
[32]).

We, therefore, introduce a novel problem that we denote as the congested
partial set covering location problem (CPSCLP), consisting of choosing where
to locate the facilities among a set of candidate sites, that can satisfy a
target (partial) demand in such a way as to minimize the cost of facility
opening and congestion. By seeking a more balanced solution, CPSCLP can
prevent facilities from being overloaded, thereby minimizing diseconomies of
scale and ensuring better resource allocation. Congestion in facility location
problems is typically modeled through a convex quadratic term in the objective
function (see, e.g., [16, 18]). This term acts similarly to a penalty function
associated with capacity constraints, as it penalizes every additional unit
being served by a given facility. Hence, taking into account the minimization
of congestion can be seen as imposing a limited capacity on the facilities, as
claimed in [16].

Another key aspect to consider is that in most real-world applications,
accurately estimating demand can be challenging due to its inherent
variability or lack of historical data. If the estimate is not correct and the
demand deviates from the expected value, it could lead to an optimal solution
that is actually impracticable or infeasible. To mitigate the sensitivity of
the solution with respect to changes in the problem parameters, we address the
CPSCLP under the assumption of uncertain customer demand. Specifically, we
propose to deal with data uncertainty using the approach known as
$\Gamma$-robustness, as introduced in [8]. This approach aims to create a
reliable and efficient network architecture that is robust against demand
changes. In particular, based on the assumption that nature is restricted in
its behaviour, we protect against the case in which deviations from expected
demand will occur for at most $\Gamma$ customers.

##### Literature overview

The problem we address emerges at the intersection of two optimization
problems: the PSCLP and the congested facility location problem (CFLP).
Despite their practical relevance, these problems have received little
attention in the scientific literature.

The minimum cost covering problem traces back to Hakimi’s work [31] in 1965,
where the problem of locating the minimum number of police officers so that
everyone is within a given distance from an officer is introduced; Hakimi
suggests a solution procedure based on enumeration. The first integer
programming formulation of the problem was proposed in [50] to solve the
problem of locating emergency service facilities in a discrete space. Then, in
[51], a heuristic is suggested to assign ladder trucks to fire stations and
the problem is formulated as a minimum cost covering problem.

PSCLP, specifically, was introduced by Daskin and Owen [15] in 1999. The
authors propose an approach based on a Lagrangian heuristic. In 2019, an exact
approach based on Benders decomposition was provided in [14], and several
large-scale instances were made available by the authors as a benchmark for
future works. More recently, in [11], five customized presolving methods have
been discussed to enhance the capability of employing mixed-integer
programming (MIP) solvers in solving PSCLPs.

Other studies have considered partial set covering problems in different
contexts than location, such as a mining application in [9], and a related
problem in which customers whose distance falls between a lower and an upper
bound from their nearest facility are only partially covered [7].

An alternative version to PSCLP is the maximal covering location problem
(MCLP), that aims for the subset of facilities maximizing the coverage while
respecting a budget constraint. The MCLP was formulated for the first time by
Church and ReVelle [13] in 1974 as a 0-1 linear programming and in [41] was
proved to be NP-hard. There are two versions of this problem in the
literature, one imposing an upper bound on the number of open facilities –
whose LP relaxation could be integral for relatively small size instances (see
[49])– and the other using a budget constraint – whose LP relaxation usually
leads to more fractional solutions. The MCLP received more attention than
PSCLP. However, despite the good quality of the LP relaxation, the exact
approaches proposed are few. Indeed, mostly heuristic and metaheuristic
algorithms have been proposed for this problem. Among them, we mention: greedy
procedures in [13, 44]; heuristics based on the Lagrangian relaxation of the
constraints in [17, 25, 26, 35]; a decomposition heuristic in [47];
metaheuristics in [45, 4, 40, 39, 38, 3]. As for the exact approaches, in
[17], a branch-and-bound framework using Lagrangian relaxations to dualize the
covering constraints is investigated. Finally, a Benders decomposition method
has also been provided for this problem in [14].

In this paper, we model congestion through convex quadratic functions, as it
was done for the congested facility location problem (CFLP). We therefore
review this class of problems. The CFLP was introduced in 1995 by Desrochers
et al. [16], inspired by [33], in which a similar formulation was proposed to
illustrate a brief example involving skiers waiting for chairlifts. The work
of [16] provides a column generation embedded in a branch-and-bound scheme,
and reports optimal solutions for very small instances of the problem. In
[32], different MIP formulations for the case of convex and piece-wise linear
production cost functions are compared. Other contributions include the two
master theses [36] and [46] investigating the problem, and the article [18],
where the authors propose a Benders decomposition method that is effective
even though the subproblem is not separable. A branch and bound algorithm
based on Lagrangian relaxation and subgradient optimization is suggested in
[12]. For a recent study on how off-the-shelf MIP solvers behave on CFLP
instances modeled throw mixed-integer second-order cone programs, see [5].

There are several other studies modeling facility congestion using queuing
theory (some of them also considering uncertainty in the data incorporated in
the models) that are not directly linked to our models. We therefore refer the
reader to the survey [10] and the literature overviews provided in [1, 52] for
further information.

For what concerns data uncertainty, there are several contributions in
stochastic or robust facility location (e.g., see the review in [48]).
However, none of them specifically addresses the PSCLP.

##### Main motivation and our contribution

A careful look into the literature suggests a lack of contributions in
problems of congested facility location and partial set covering location.
Recent works [14, 18, 19] showed promising results given by exact methods
based on Benders decomposition for the deterministic PSCLP and CFLP (both the
capacitated and uncapacitated cases). However, no consideration is given to
the inherent volatility of the parameters used to model these specific
problems. Motivated by the theoretical and practical relevance, we address the
robust and congested variant of PSCLP, which considers the minimization of the
congestion at the facilities and the protection against the changes in demand,
and investigate the use of Benders decomposition to solve large instances of
this problem.

The contributions of this paper can be summarized as follows:

  1. 1.

we state and formulate the congested partial set covering location problem
under the deterministic setting;

  2. 2.

we consider the uncertainty in the customer demand and formulate the robust
counterpart of the problem leveraging the $\Gamma$-robustness theory. The
demand changes affect both the congestion function and the covering
constraint. We show how the nonlinear protection functions can be linearized;

  3. 3.

we consider the case of quadratic convex congestion cost functions and apply
the perspective reformulation to the robust counterpart of the problem, aiming
for a tighter formulation. We then develop a tailored solution approach for
the resulting model based on Benders decomposition;

  4. 4.

we illustrate how to integrate our Benders approach within a MISOCP solver,
explicitly addressing all the ingredients that are instrumental for successful
implementation. Although Benders decomposition is a traditional technique, our
approach is innovative as it leverages the implementation of callback
functions in combination with quadratic constraints, a novel option provided
by modern MISOCP solvers. We assess both single-tree and multi-tree Benders
approaches through a comprehensive set of (adapted) instances from existing
literature, aiming to validate their efficiency and compare their performance
with a state-of-the-art MISOCP solver;

  5. 5.

we test the effect of the protection level in the coverage and load (together
and separately) to provide critical insights into how sensitive the planning
is to the protection level, offering a more comprehensive understanding of the
problem.

The paper is structured as follows. In Section 2, we introduce the
deterministic formulation of the problem and its robust counterpart,
explaining how to account for uncertainties in customer demand using the
$\Gamma$-robustness approach. Section 3 focuses on the Benders decomposition
method: in this section, we provide the master and subproblem formulations,
and derive the expression of the Benders cut. In Section 4 we show a cut-
strengthening technique to deal with the degeneracy of the Benders subproblem.
Section 5 provides insights into the integration of Benders approach in a
MISOCP solver. We describe in detail both the single-tree and multi-tree
approaches. In Section 6, we report the sensitivity analysis testing the
effect of $\Gamma$ on the optimal solutions and the computational experiments
comparing our tailored Benders to a state-of-the-art solver. Conclusions are
given in Section 7.

##  2 Problem Setting

We are given a set $I$ of potential facility locations with opening cost
$f_{i}\geq 0$ for $i\in I$, and a set $J$ of customer locations such that each
customer location $j\in J$ is associated with a demand $d_{j}\geq 0$. From now
on, we refer for short to facility and customer, omitting the word location.
For each customer $j$, we are also given a subset $I(j)\subseteq I$ of
facilities that can cover $j$, i.e., that can fully serve the demand $d_{j}$.
Similarly, for $i\in I$, let $J(i)$ be the subset of customers that can be
covered by $i$. More generally, for a subset of facilities $N\subseteq I$, let
$J(N)\subseteq J$ be the subset of customers that can be covered by an $i\in
N$.

Given a parameter $0<D\leq\sum_{j\in J}d_{j}$, we aim to identify a subset of
facilities to open in order to ensure that the total served customer demand is
at least $D$, while minimizing the overall costs, given by facility opening
expenses and congestion costs. To model congestion cost, we follow what is
done in the literature of CFLP [16, 18, 33], in which a nonlinear cost
function $F_{i}$ is used to penalize each additional unit of demand served by
a given facility $i$. For each facility $i\in I$, let
$F_{i}:\mathbb{R}\to\mathbb{R}_{+}$ be a convex and non-decreasing function
representing the congestion cost for each unit of served demand.

For each $i\in I$, the binary variable $y_{i}$ is set to one if facility $i$
is open and to zero otherwise. For each $i\in I$ and $j\in J$, the continuous
allocation variable $0\leq x_{ij}\leq 1$ denotes the fraction of the demand of
customer $j\in J$ served by facility $i\in I$.

###  2.1 The deterministic problem formulation

The most general deterministic CPSCLP can be formulated as a mixed-integer
nonlinear program (MINLP) as follows

$\displaystyle\displaystyle{\min}\quad$ $\displaystyle\sum_{i\in
I}f_{i}y_{i}+\sum_{i\in I}F_{{i}}\left(\sum_{j\in J}d_{j}x_{ij}\right)$ (1)
$\displaystyle{s.t.}$ $\displaystyle\quad\sum_{i\in I}\sum_{j\in
J}d_{j}x_{ij}\geq D$ (2) $\displaystyle\sum_{i\in I(j)}x_{ij}\leq 1$
$\displaystyle j\in J$ (3) $\displaystyle 0\leq x_{ij}\leq y_{i}$
$\displaystyle j\in J,i\in I(j)$ (4) $\displaystyle y_{i}\in\\{0,1\\}$
$\displaystyle i\in I.$

The objective (1) aims to minimize the sum of facility opening and congestion
costs. For the sake of simplicity, we assign equal weights to the two
conflicting objectives. However, different (non-negative) weights can be
easily incorporated without altering the structure of the problem. Constraint
(2) ensures the total customer demand served is at least $D$. Assignment
constraints (3) make sure that the fraction of covered demand of a customer
does not exceed the unit. Finally, constraints (4) ensure that allocation to a
facility is only possible if it is open.

This problem is NP-hard since it is a generalization of the traditional set
covering location problem, which is NP-hard [15].

###  2.2 The robust counterpart of the problem

In real-world scenarios, customer demand often varies or is difficult to
estimate. To capture this uncertainty, we assume that each entry $d_{j}$,
$j\in J$ of the vector of demand $d$ is modeled as an independent, symmetric
and bounded random variable (with unknown distribution) $\tilde{d}_{j}$, $j\in
J$ that takes values in $[d_{j}-\hat{d}_{j},d_{j}+\hat{d}_{j}]$. We allow the
possibility that the deviations from the nominal values could also be zeros,
i.e. that $\hat{d}_{j}=0$ for some $j\in J$. We adopt the notion of protection
introduced by Bertsimas and Sim in [8] known as $\Gamma$-robustness, which
assumes that only a subset of the components of $d$ will deviate from their
nominal values, adversely affecting the solution. Hence, we introduce an
integer number $\Gamma$, taking values in the interval $[0,|J|]$, that limits
the number of demand deviations. Parameter $\Gamma$ controls the level of
robustness against the solution: if $\Gamma=0$, we completely ignore the
uncertainty (deterministic setting), while if $\Gamma=|J|$, we are considering
all possible demand deviations (which is the most conservative strategy).

We note that the vector of demand $d$ is involved in the modeling of
congestion and coverage affecting both optimality and feasibility. By
protecting against the uncertainty, we mean that we are interested in finding
an optimal solution that:

  1. 1.

optimizes against all scenarios under which up to $\Gamma$ demand coefficients
can vary in such a way as to maximally influence the objective; the worst-case
scenario is given by the largest increase of demand;

  2. 2.

is protected against all cases in which up to $\Gamma$ demand coefficients
change affecting the feasibility; the worst-case scenario is given by the
largest decrease of demand.

We observe that the two worst-case realizations play against each other.

We now introduce a robust counterpart of the problem, that optimizes against
the worst-case realizations under demand uncertainty and reads

$\displaystyle{\min}$ $\displaystyle\sum_{i\in I}f_{i}y_{i}+\sum_{i\in
I}F_{i}\left(\sum_{j\in J}d_{j}x_{ij}+\max_{\\{S:\ S\subseteq
J,|S|\leq\Gamma\\}}\left\\{\sum_{j\in S}\hat{d}_{j}x_{ij}\right\\}\right)$ (5)
$\displaystyle{s.t.}$ $\displaystyle\sum_{i\in I}\sum_{j\in
J}d_{j}x_{ij}-\max_{\\{S:\ S\subseteq J,|S|\leq\Gamma\\}}\left\\{\sum_{i\in
I}\sum_{j\in S}\hat{d}_{j}x_{ij}\right\\}\geq D$ $\displaystyle\sum_{i\in
I(j)}x_{ij}\leq 1\hskip 199.73846ptj\in J$ $\displaystyle 0\leq x_{ij}\leq
y_{i}\hskip 203.43703ptj\in J,i\in I(j)$ $\displaystyle
y_{i}\in\\{0,1\\}\hskip 210.55022pti\in I$

where, for a given $\Gamma$ and allocation choice $x_{ij}$, the objective
function is taking into account the sum of the $\Gamma$ largest deviations in
case the demand is increasing from the nominal value, whereas in the covering
constraint, we are considering the sum of the $\Gamma$ largest deviations in
case the demand is decreasing from the nominal value.

We now introduce the auxiliary aggregated variables $v_{i}\geq 0$ for $i\in
I$, denoting the total demand served by facility $i$ (also known as facility
load). Then, Problem (5) becomes

$\displaystyle\displaystyle{\min_{x\geq 0,y\in\\{0,1\\}^{|I|}}}\quad$
$\displaystyle\displaystyle\sum_{i\in I}f_{i}y_{i}+\sum_{i\in I}F_{i}(v_{i})$
(6a) $\displaystyle{s.t.}\quad$ $\displaystyle\displaystyle
v_{i}\geq\sum_{j\in J}d_{j}x_{ij}+\max_{\\{S:\ S\subseteq
J,|S|\leq\Gamma\\}}\left\\{\sum_{j\in S}\hat{d}_{j}x_{ij}\right\\}$
$\displaystyle i\in I$ (6b) $\displaystyle\displaystyle\sum_{i\in I}\sum_{j\in
J}d_{j}x_{ij}-\max_{\\{S:\ S\subseteq J,|S|\leq\Gamma\\}}\left\\{\sum_{i\in
I}\sum_{j\in S}\hat{d}_{j}x_{ij}\right\\}\geq D$ (6c)
$\displaystyle\displaystyle\sum_{i\in I}x_{ij}\leq 1$ $\displaystyle j\in J$
(6d) $\displaystyle\displaystyle x_{ij}\leq y_{i}$ $\displaystyle j\in J,i\in
I(j)$ (6e)

Following [8], by applying strong duality to the inner maximization robust
terms, we can derive the following equivalent MINLP formulation of this
problem

$\displaystyle\displaystyle{\min_{(x,\tau,\rho,\pi,\sigma)\geq
0,\,y\in\\{0,1\\}^{|I|}}}\quad$ $\displaystyle\sum_{i\in
I}f_{i}y_{i}+\sum_{i\in I}F_{i}(v_{i})$ $\displaystyle{s.t.}\quad$
$\displaystyle v_{i}-\sum_{j\in J}d_{j}x_{ij}-\left(\Gamma\rho_{i}+\sum_{j\in
J}\sigma_{ij}\right)\geq 0$ $\displaystyle i\in I$ (7a)
$\displaystyle\sum_{i\in I}\sum_{j\in
J}d_{j}x_{ij}-\left(\Gamma\tau+\sum_{j\in J}\pi_{j}\right)\geq D$ (7b)
$\displaystyle\tau+\pi_{j}\geq\sum_{i\in I}\hat{d}_{j}x_{ij}$ $\displaystyle
j\in J$ (7c) $\displaystyle\rho_{i}+\sigma_{ij}\geq\hat{d}_{j}x_{ij}$
$\displaystyle j\in J,i\in I$ (7d) $\displaystyle\sum_{i\in I}x_{ij}\leq 1$
$\displaystyle j\in J$ (7e) $\displaystyle x_{ij}\leq y_{i}$ $\displaystyle
j\in J,i\in I(j)$ (7f)

We can state the following theorem, whose proof is in the Appendix A.

###### Theorem 2.1

Formulation (7) is equivalent to (6).

Problem (7) is a much more computationally tractable formulation. Indeed, the
problem is now formulated as a mixed-integer program with linear constraints
and a convex objective function.

####  2.2.1 The perspective reformulation

Following the work of [16, 18], we consider the specific case where the
congestion function is given by the convex function
$F_{i}(t)=a_{i}t^{2}+b_{i}t$ for $i\in I$, with $a_{i}$ and $b_{i}$ non-
negative input coefficients. Then, the robust counterpart of the problem,
formulated in (7), reads

$\displaystyle\displaystyle{\min_{(x,\tau,\rho,\pi,\sigma)\geq
0,\,y\in\\{0,1\\}^{|I|}}}\quad$ $\displaystyle\sum_{i\in
I}f_{i}y_{i}+\sum_{i\in I}b_{i}v_{i}+\sum_{i\in I}a_{i}v_{i}^{2}$
$\displaystyle{s.t.}\quad$ $\displaystyle\eqref{first}-\eqref{last}.$

which is a mixed-integer (convex) quadratic program (MIQP) with binary
variables. Therefore, we can derive the so-called perspective reformulation
[23, 28]

$\displaystyle\displaystyle{\min_{(x,u,v,\tau,\rho,\pi,\sigma)\geq
0,\,y\in\\{0,1\\}^{|I|}}}\quad$ $\displaystyle\sum_{i\in
I}f_{i}y_{i}+\sum_{i\in I}b_{i}v_{i}+\sum_{i\in I}a_{i}u_{i}$
$\displaystyle{s.t.}\quad$ $\displaystyle\eqref{first}-\eqref{last}$
$\displaystyle v_{i}^{2}\leq u_{i}y_{i}$ $\displaystyle i\in I$ (8a)

in which the quadratic term $v_{i}^{2}$ in the objective function is replaced
by a non-negative variable $u_{i}$. Constraints (8a) are rotated second-order
cone constraints (hence still convex) as the right-hand side is the product of
non-negative variables. They are obtained by the strengthening of the $u_{i}$
variables through constraints $v_{i}^{2}\leq u_{i}$, by replacing the convex
function $v_{i}^{2}$ with its perspective defined by $y_{i}(v_{i}/y_{i})^{2}$
if $y_{i}>0$, and zero if $y_{i}=0$; see [28] for details. It is well-known
(see e.g., [2, 21, 22, 24, 27, 28]) that the continuous relaxation of a
perspective reformulation produces stronger bounds than the bounds given by
the continuous relaxation of the original formulation. Consequently, we can
state that the original formulation (7), with $F_{i}(t)=a_{i}t^{2}+b_{i}t$,
lies in an extended space compared to the formulation (8), and for this
reason, we denote it as the extended formulation.

##  3 Benders Decomposition

Both the extended and the perspective formulations have a number of variables
and constraints that depend on the number of customers $|J|$ and facilities
$|I|$. In real world context, the number of customers can be quite large,
affecting the solution time of problems (7) or (8). In this regard, we propose
to use Benders decomposition, a classical method for mixed-integer linear
programming introduced by Benders [6] in 1962. Our idea is to reduce as much
as possible the number of variables in the problem, by projecting out of the
master problem (at least) all the variables that depend on the number of
customers. This will lead to a boosting in the solution time, as we will see
in Section 6. Specifically, we investigate the projecting out of the
continuous variables $(x,\tau,\rho,\pi,\sigma)$, leading to a convex master
problem containing the cones coming from the perspective reformulation and the
complete objective function, and an LP subproblem. The advantage of having the
cones in the master formulation is that the master problem has a very tight
relaxation, which is amenable when dealing with Benders decomposition as it
can substantially reduce solution times. We also observe that the subproblem –
which is solved several times to generate valid cuts for the master – is an LP
that can be therefore efficiently solved by any state-of-the-art solver.

###  3.1 The master problem

From the formulation (8), we project out the continuous variables
$(x,\rho,\tau,\pi,\sigma)$. The master problem is then defined on the
variables $(y,u,v)$ and is given by

$\begin{array}[]{rlrr}\displaystyle{\min_{(u,v)\geq
0,\,y\in\\{0,1\\}^{|I|}}}&\displaystyle\sum_{i\in I}f_{i}y_{i}+\sum_{i\in
I}b_{i}v_{i}+\sum_{i\in I}a_{i}u_{i}\\\ {s.t.}&\phi(y,v)\leq 0\\\
&v_{i}^{2}-u_{i}y_{i}\leq 0&&i\in I\end{array}$ (9)

where $\phi(y,v)$ is a convex function that evaluates the feasibility of the
vector $(y,v)$ for the constraints (7a)-(7f). It can be approximated by linear
cuts to be generated on the fly, known as Benders cuts, that are valid for any
given vector $(y,v)$. Usually, two types of Benders cuts are used: the
optimality cuts and the feasibility cuts. In our approach, there is no need
for Benders optimality cuts, as all the variables appearing in the objective
function $(y,v,u)$ belong to the master problem. Hence, we will only use
Benders feasibility cuts to discard infeasible points. The number of these
cuts is exponential, making enumerating them all at once impractical. Since
only some of them are necessary to find an optimal solution, they will be
dynamically separated by the decomposition approach. Consequently, the master
problem will contain only a subset of Benders cuts, belonging to a so-called
_relaxed master problem_.

The relaxed master problem is a mixed-integer second order cone program
(MISOCP), thus a convex MINLP. It can be solved as a mixed-integer linear
program by a branch-and-cut approach where the integrality requirement on $y$
is relaxed and linear outer-approximations of the conic constraints are
generated on the fly or a nonlinear programming solver is used for the
underlying NLP problem at every node [29].

###  3.2 The subproblem

Given a master solution $(\bar{y},\bar{v},\bar{u})$, the subproblem is then
given by

$\displaystyle\displaystyle\phi(\bar{y},\bar{v})={\min_{(x,\rho,\tau,\pi,\sigma)\geq
0}}\quad$ $\displaystyle 0$ $\displaystyle{s.t.}\quad$
$\displaystyle\sum_{j\in{J}}d_{j}x_{ij}+\Gamma\rho_{i}+\sum_{j\in{J}}\sigma_{ij}\leq\bar{v}_{i}$
$\displaystyle i\in I$ $\displaystyle\sum_{j\in
J}\sum_{i\in{I}}d_{j}x_{ij}-\Gamma\tau-\sum_{j\in J}\pi_{j}\geq D$
$\displaystyle\tau+\pi_{j}-\sum_{i\in{I}}\hat{d}_{j}x_{ij}\geq 0$
$\displaystyle j\in J$
$\displaystyle\sigma_{ij}+\rho_{i}-\hat{d}_{j}x_{ij}\geq 0$ $\displaystyle
j\in J,i\in I(j)$ $\displaystyle\sum_{i\in{I}}x_{ij}\leq 1$ $\displaystyle
j\in J$ $\displaystyle x_{ij}\leq\bar{y}_{i}$ $\displaystyle i\in I(j),j\in
J.$

which is a linear programming problem. The problem above can be feasible and
in this case there is no need to generate a Benders cut, or infeasible. In the
latter case a Benders feasibility cut is produced, i.e., a cutting plane
$\phi(\bar{y},\bar{v})\leq 0$ that discards the master solution
$(\bar{y},\bar{v},\bar{u})$ that leads to the infeasibilty of the subproblem.
Solving a subproblem that could potentially be infeasible may lead to
computational issues. An infeasible primal problem means that in the dual
space we are optimizing over an unbounded cone pointed at zero. Among the
successful strategies to overcome these difficulties are the normalization
techniques that consist of solving the dual LP over a bounded polyhedron.
There is abundant literature on different normalization techniques for Benders
feasibility cuts, see, e.g., [14, 20, 34, 37, 42]. We use a very natural
approach exploiting the fact that the solution $(\bar{y},\bar{v},\bar{u})$ of
the relaxed master problem is infeasible for the subproblem if and only if the
demand covered by $\bar{y}$ is strictly less than $D$. Indeed, in our case,
the covering constraint forms the irreducible infeasible subsystem, which is
the minimal subset of constraints whose removal makes the problem feasible.
Hence, instead of solving a feasibility LP given by the formulation above, we
search for the maximum coverage. The resulting subproblem is the following LP

$\displaystyle\displaystyle\phi^{\prime}(\bar{y},\bar{v})={\max_{(x,\rho,\tau,\pi,\sigma)\geq
0}}\quad$ $\displaystyle\sum_{j\in
J}\sum_{i\in{I}}d_{j}x_{ij}-\Gamma\tau-\sum_{j\in J}\pi_{j}$
$\displaystyle{s.t.}\quad$
$\displaystyle\sum_{j\in{J}}d_{j}x_{ij}+\Gamma\rho_{i}+\sum_{j\in{J}}\sigma_{ij}\leq\bar{v}_{i}$
$\displaystyle i\in I$
$\displaystyle\tau+\pi_{j}-\sum_{i\in{I}}\hat{d}_{j}x_{ij}\geq 0$
$\displaystyle j\in J$
$\displaystyle\sigma_{ij}+\rho_{i}-\hat{d}_{j}x_{ij}\geq 0$ $\displaystyle
j\in J,i\in I(j)$ $\displaystyle\sum_{i\in{I}}x_{ij}\leq 1$ $\displaystyle
j\in J$ $\displaystyle x_{ij}\leq\bar{y}_{i}$ $\displaystyle j\in J,i\in
I(j),$

in which we observe that for every given $(\bar{y},\bar{v})$, the vector
$(x,\rho,\tau,\pi,\sigma)=(0,0,0,0,0)$ is a feasible solution, i.e., this
problem is always feasible and bounded as $(x,\tau,\pi)$ are bounded. If the
optimal value of the subproblem $\phi^{\prime}(\bar{y},\bar{v})$ is greater or
equal to $D$, the master solution is feasible and we do not need to generate
any Benders cut. If, instead, $\phi^{\prime}(\bar{y},\bar{v})$ is strictly
less than $D$, the master solution is infeasible as it does not guarantee the
required coverage and we do generate a Benders feasibility cut. The Benders
cut is then given by $\phi^{\prime}(y,v)\geq D$ (and no more by $\phi(y,v)\leq
0$) for every solution $(y,v)$ of the relaxed master problem. We notice that
the extreme points of the polyhedron describing the problem
$\phi^{\prime}(\bar{y},\bar{v})$ correspond to extreme rays of the subproblem
associated to $\phi(\bar{y},\bar{v})$ in which the dual variable of the
covering constraint is fixed to one.

Because of concavity, $\phi^{\prime}(\cdot)$ can be overestimated by a
supporting hyperplane at $(\bar{y},\bar{v})$, so the following linear cut is
also valid

$\phi^{\prime}(\bar{y},\bar{v})+\xi_{y}(\bar{y},\bar{v})^{T}(y-\bar{y})+\xi_{v}(\bar{y},\bar{v})^{T}(v-\bar{v})\geq\phi^{\prime}(y,v)\geq
D$

where $\xi_{y}(\bar{y},\bar{v}),\,\xi_{v}(\bar{y},\bar{v})$ denote the
subgradients of $\phi^{\prime}$ with respect to $y,v$ in $(\bar{y},\bar{v})$.

Depending on the problem, the computation of the subgradients could be heavy.
Therefore, we introduce a simple reformulation of the Benders subproblem that
makes their calculation straightforward, following what was done, e.g., in
[18]. The reformulation of the subproblem reads

$\displaystyle\displaystyle\phi^{\prime}(\bar{y},\bar{v})=\max_{y,v,(x,\rho,\tau,\pi,\sigma)\geq
0}\quad$ $\displaystyle\sum_{j\in
J}\sum_{i\in{I}}d_{j}x_{ij}-\Gamma\tau-\sum_{j\in J}\pi_{j}$
$\displaystyle{s.t.}\quad$
$\displaystyle\sum_{j\in{J}}d_{j}x_{ij}+\Gamma\rho_{i}+\sum_{j\in{J}}\sigma_{ij}-v_{i}\leq
0$ $\displaystyle i\in I$
$\displaystyle\tau+\pi_{j}-\sum_{i\in{I}}\hat{d}_{j}x_{ij}\geq 0$
$\displaystyle j\in J$
$\displaystyle\sigma_{ij}+\rho_{i}-\hat{d}_{j}x_{ij}\geq 0$ $\displaystyle
j\in J,i\in I(j)$ $\displaystyle\sum_{i\in{I}}x_{ij}\leq 1$ $\displaystyle
j\in J$ $\displaystyle x_{ij}-y_{i}\leq 0$ $\displaystyle j\in J,i\in I(j)$
$\displaystyle y_{i}=\bar{y}_{i}$ $\displaystyle i\in I$ (10a) $\displaystyle
v_{i}=\bar{v}_{i}$ $\displaystyle i\in I$ (10b)

where we keep the master variables $y$ and $v$ as variables of the subproblem
as well, and we apply variable fixing through (10a)–(10b). By construction,
the subgradients are simply $\xi_{y}(\bar{y},\bar{v})=\bar{r}_{y}$ and
$\xi_{v}(\bar{y},\bar{v})=\bar{r}_{v}$, where $\bar{r}_{y}$ and $\bar{r}_{v}$
are the vectors of reduced costs associated with $y$ and $v$. The Benders cut
reads

$\phi^{\prime}(\bar{y},\bar{v})+\sum_{i\in
I}\bar{r}_{y_{i}}(y_{i}-\bar{y}_{i})+\sum_{i\in
I}\bar{r}_{v_{i}}(v_{i}-\bar{v}_{i})\geq D$ (11)

where each component of the reduced cost $\bar{r}_{y_{i}}$ (or
$\bar{r}_{v_{i}}$) defines an upper bound on the increase of the objective
function $\phi^{\prime}(\bar{y},\bar{v})$ when $\bar{y}_{i}$ (or
$\bar{v}_{i}$) increases.

##  4 Addressing degeneracy

When applying Benders decomposition, a highly degenerate subproblem admitting
several optimal solutions can result in the generation of shallow Benders
cuts. This slows down the convergence of the Benders decomposition, requiring
the addition of many cuts, which do not improve the bound that much. There are
many techniques proposed in the literature to address degeneracy (see, e.g.,
[43]). Inspired by the perturbation technique proposed in [19] to accelerate
the convergence of a cut loop (i.e., the cut separation at the root node of
the branching tree), we adapt this technique, that we will refer to as
$\epsilon$-technique, leading to stronger cuts. The result is an accelerated
convergence and a fewer number of generated cuts, as shown from the
computational experience provided later in Section 6.

We start by interpreting the Benders cut from the LP primal-dual perspective,
and then show how to apply this perturbation to derive valid Benders cuts for
the non-perturbed subproblem. Consider the non-perturbed subproblem formulated
as follows

$\displaystyle\phi^{\prime}(\bar{y},\bar{v})={\max_{y,v,(x,\rho,\tau,\pi,\sigma)\geq
0}}$ $\displaystyle\sum_{j\in
J}\sum_{i\in{I}}d_{j}x_{ij}-\Gamma\tau-\sum_{j\in J}\pi_{j}$ (12)
$\displaystyle{s.t.}$
$\displaystyle\sum_{j\in{J}}d_{j}x_{ij}+\Gamma\rho_{i}+\sum_{j\in{J}}\sigma_{ij}-v_{i}\leq
0$ $\displaystyle i\in I$
$\displaystyle\tau+\pi_{j}-\sum_{i\in{I}}\hat{d}_{j}x_{ij}\geq 0$
$\displaystyle j\in J$
$\displaystyle\sigma_{ij}+\rho_{i}-\hat{d}_{j}x_{ij}\geq 0$ $\displaystyle
j\in J,i\in I(j)$ $\displaystyle\sum_{i\in{I}}x_{ij}\leq 1$ $\displaystyle
j\in J$ $\displaystyle(\alpha_{j})$ $\displaystyle x_{ij}-y_{i}\leq 0$
$\displaystyle j\in J,i\in I(j)$ $\displaystyle y_{i}=\bar{y}_{i}$
$\displaystyle i\in I$ $\displaystyle(\beta_{i})$ $\displaystyle
v_{i}=\bar{v}_{i}$ $\displaystyle i\in I.$ $\displaystyle(\gamma_{i})$

where by $\alpha,\,\beta$ and $\gamma$ we denote the dual variables (or the
Lagrangian multipliers). Given the optimal dual solution
$(\bar{\alpha},\bar{\beta},\bar{\gamma})$, the optimal value of the dual of
(12) is

$\phi^{\prime}(\bar{y},\bar{v})=\sum_{j\in J}\bar{\alpha}_{j}+\sum_{i\in
I}\bar{y}_{i}\bar{\beta}_{i}+\sum_{i\in I}\bar{v}_{i}\bar{\gamma}_{i}.$

We can get the expression of the Benders cut from the LP dual by imposing that
the dual objective at the optimal solution is at least the target coverage,
namely

$\sum_{j\in J}\bar{\alpha}_{j}+\sum_{i\in I}\bar{\beta}_{i}y_{i}+\sum_{i\in
I}\bar{\gamma}_{i}v_{i}\geq D.$

The dual problem is highly degenerate as many objective coefficients, i.e.,
many components of $\bar{y}$ and $\bar{v}$ can be zero, resulting in several
free components of the dual variables $\beta$ and $\gamma$ and many equivalent
optimal solutions. To get the strongest Benders cut, we need the dual
multipliers appearing in the cut to take the smallest values possible. To get
this, we can apply the $\epsilon$-technique that consists of replacing the
zero objective coefficients of $\bar{y}$ and $\bar{v}$ with a sufficiently
small $\epsilon>0$ and solving the resulting problem. This induces the dual
model to minimize also the components of the dual variables $\beta$ and
$\gamma$ associated with the originally zero objective coefficients.

Unfortunately, solving the dual of the subproblem often leads to more
numerical issues than solving its primal, hence we explain how to apply this
technique to the primal subproblem (12). Specifically, given a sufficiently
small $\epsilon>0$, applying the $\epsilon$-technique to the primal problem
consists of replacing $\bar{y}$ and $\bar{v}$ with

$\bar{y}^{\epsilon}_{i}=\begin{cases}\bar{y}_{i}&\text{if }\bar{y}_{i}>0\\\
\epsilon&\text{if
}\bar{y}_{i}=0\end{cases}\qquad\text{and}\qquad\bar{v}^{\epsilon}_{i}=\begin{cases}\bar{v}_{i}&\text{if
}\bar{v}_{i}>0\\\ \epsilon&\text{if }\bar{v}_{i}=0\end{cases}$

for each $i\in I$. Then, the perturbed subproblem reads

$\displaystyle\phi^{\prime}(\bar{y}^{\epsilon},\bar{v}^{\epsilon})={\max_{y,v,(x,\rho,\tau,\pi,\sigma)\geq
0}}$ $\displaystyle\sum_{j\in
J}\sum_{i\in{I}}d_{j}x_{ij}-\Gamma\tau-\sum_{j\in J}\pi_{j}$ (13)
$\displaystyle{s.t.}$
$\displaystyle\sum_{j\in{J}}d_{j}x_{ij}+\Gamma\rho_{i}+\sum_{j\in{J}}\sigma_{ij}-v_{i}\leq
0$ $\displaystyle i\in I$
$\displaystyle\tau+\pi_{j}-\sum_{i\in{I}}\hat{d}_{j}x_{ij}\geq 0$
$\displaystyle j\in J$
$\displaystyle\sigma_{ij}+\rho_{i}-\hat{d}_{j}x_{ij}\geq 0$ $\displaystyle
j\in J,i\in I(j)$ $\displaystyle\sum_{i\in{I}}x_{ij}\leq 1$ $\displaystyle
j\in J$ $\displaystyle(\alpha_{j}^{\epsilon})$ $\displaystyle x_{ij}-y_{i}\leq
0$ $\displaystyle j\in J,i\in I(j)$ $\displaystyle
y_{i}=\bar{y}^{\epsilon}_{i}$ $\displaystyle i\in I$
$\displaystyle(\beta_{i}^{\epsilon})$ $\displaystyle
v_{i}=\bar{v}^{\epsilon}_{i}$ $\displaystyle i\in I.$
$\displaystyle(\gamma_{i}^{\epsilon})$

where by $\alpha^{\epsilon},\beta^{\epsilon}$ and $\gamma^{\epsilon}$ we
denote the dual variables. The Benders cut obtained from the Lagrangian dual
of (13) is

$\phi^{\prime}(\bar{y}^{\epsilon},\bar{v}^{\epsilon})+\sum_{i\in
I}\bar{r}_{y_{i}}^{\epsilon}(y_{i}-\bar{y}^{\epsilon}_{i})+\sum_{i\in
I}\bar{r}_{v_{i}}^{\epsilon}(v_{i}-\bar{v}^{\epsilon}_{i})\geq D$

where $\bar{r}_{y}^{\epsilon}$ and $\bar{r}_{v}^{\epsilon}$ are the vectors of
reduced costs associated with $y$ and $v$. We now state some properties of the
perturbed problem and clarify the relationship between the perturbed and the
non-perturbed problem.

###### Lemma 4.1

Let $(\bar{r}_{y},\bar{r}_{v})$ and
$(\bar{r}^{\epsilon}_{y},\bar{r}^{\epsilon}_{v})$ be the reduced costs
introduced for the Benders subproblem (12) and its $\epsilon$-perturbed
counterpart (13). The following equations are valid

$\bar{r}_{y}=\bar{\beta},\quad\bar{r}_{v}=\bar{\gamma},\quad\bar{r}^{\epsilon}_{y}=\bar{\beta}^{\epsilon},\quad\bar{r}_{v}^{\epsilon}=\bar{\gamma}^{\epsilon}.$

  * Proof

Consider the non-perturbed problem (12). Since (12) is an LP, from the
equivalence between LP duality and Lagrangian duality for LP problems we have

$\sum_{j\in J}\bar{\alpha}_{j}+\sum_{i\in I}\bar{\beta}_{i}y_{i}+\sum_{i\in
I}\bar{\gamma}_{i}v_{i}=\phi^{\prime}(\bar{y},\bar{v})+\sum_{i\in
I}\bar{r}_{y_{i}}(y_{i}-\bar{y}_{i})+\sum_{i\in
I}\bar{r}_{v_{i}}(v_{i}-\bar{v}_{i})$

which implies

$\sum_{i\in I}\bar{\beta}_{i}y_{i}=\sum_{i\in
I}\bar{r}_{y_{i}}y_{i},\quad\sum_{i\in I}\bar{\gamma}_{i}v_{i}=\sum_{i\in
I}\bar{r}_{v_{i}}v_{i},\quad\sum_{j\in
J}\bar{\alpha}_{j}=\phi^{\prime}(\bar{y},\bar{v})-\sum_{i\in
I}\bar{r}_{y_{i}}\bar{y}_{i}-\sum_{i\in I}\bar{r}_{v_{i}}\bar{v}_{i},$

meaning that $\bar{r}_{y_{i}}=\bar{\beta}_{i}$ and
$\bar{r}_{v_{i}}=\bar{\gamma}_{i}$ for each $i\in I$, namely the dual
variables associated with the fixing constraints are the reduced costs
associated with the fixed variables. By applying the same procedure to the
perturbed problem (13), we get
$\bar{r}^{\epsilon}_{y_{i}}=\bar{\beta}_{i}^{\epsilon}$ and
$\bar{r}_{v_{i}}^{\epsilon}=\bar{\gamma}_{i}^{\epsilon}$ for each $i\in I$.

$\square$

###### Proposition 4.2

If there exists a sufficiently small $\epsilon>0$ such that the dual solution
$(\alpha^{\epsilon},\beta^{\epsilon},\gamma^{\epsilon})$ of the perturbed
problem (13) is an optimal solution of the dual of the non-perturbed problem
(12), then

  * (1)

$\phi^{\prime}(\bar{y},\bar{v})=\displaystyle\phi^{\prime}(\bar{y}^{\epsilon},\bar{v}^{\epsilon})-\epsilon\sum_{i\in
I:\,\bar{y}_{i}=0}\bar{r}_{y_{i}}^{\epsilon}-\epsilon\sum_{i\in
I:\,\bar{v}_{i}=0}\bar{r}_{v_{i}}^{\epsilon}$, and

  * (2)

the following inequality

$\phi^{\prime}(\bar{y},\bar{v})+\sum_{i\in
I}\bar{r}_{y_{i}}^{\epsilon}(y_{i}-\bar{y}_{i})+\sum_{i\in
I}\bar{r}_{v_{i}}^{\epsilon}(v_{i}-\bar{v}_{i})\geq D$ (14)

is a valid Benders cut for the non-perturbed problem (12).

  * Proof

From the definition of $y^{\epsilon}$ and $v^{\epsilon}$ we have

$\displaystyle\sum_{i\in I}\bar{y}^{\epsilon}_{i}=\sum_{i\in
I}\bar{y}_{i}+\epsilon\sum_{i\in I:\,\bar{y}_{i}=0}1$ (15)
$\displaystyle\sum_{i\in I}\bar{v}^{\epsilon}_{i}=\sum_{i\in
I}\bar{v}_{i}+\epsilon\sum_{i\in I:\,\bar{v}_{i}=0}1.$ (16)

Consider the expression of the dual objective of the perturbed problem (13) at
optimal value

$\phi^{\prime}(\bar{y}^{\epsilon},\bar{v}^{\epsilon})=\sum_{j\in
J}\bar{\alpha}_{j}^{\epsilon}+\sum_{i\in
I}\bar{y}_{i}^{\epsilon}\bar{\beta}_{i}^{\epsilon}+\sum_{i\in
I}\bar{v}_{i}^{\epsilon}\bar{\gamma}_{i}^{\epsilon}.$

By using (15)-(16) we get

$\displaystyle\phi^{\prime}(\bar{y}^{\epsilon},\bar{v}^{\epsilon})=\underbrace{\sum_{j\in
J}\bar{\alpha}_{j}^{\epsilon}+\sum_{i\in
I}\bar{y}_{i}\bar{\beta}_{i}^{\epsilon}+\sum_{i\in
I}\bar{v}_{i}\bar{\gamma}_{i}^{\epsilon}}_{\phi^{\prime}(\bar{y},\bar{v})}+\epsilon\sum_{i\in
I:\,\bar{y}_{i}=0}\overbrace{\bar{\beta}_{i}^{\epsilon}}^{\bar{r}_{y_{i}}^{\epsilon}}+\epsilon\sum_{i\in
I:\,\bar{v}_{i}=0}\overbrace{\bar{\gamma}_{i}^{\epsilon}}^{\bar{r}_{v_{i}}^{\epsilon}}$

which proves (1) for the assumption we made on $\epsilon$ and using Lemma 4.1.

If $\epsilon$ is sufficiently small that
$(\alpha^{\epsilon},\beta^{\epsilon},\gamma^{\epsilon})$ is an optimal
solution of the dual of the non-perturbed problem (12), then for Lemma 4.1, we
can use $\bar{r}^{\epsilon}_{y_{i}}$ and $\bar{r}^{\epsilon}_{v_{i}}$ in (11)
and (2) is proved.

$\square$

##  5 Implementation Details

Having a correct and efficient generation of the Benders cuts is a crucial
point of Benders decomposition. In this section, we describe the steps for the
design of a Benders decomposition, revealing all the implementation
ingredients that play an important role in the design of an effective code. In
particular, we start by illustrating how to implement in practice Benders
decomposition in a MISOCP solver in Section 5.1. The discussion is based on
the solver we used, namely Gurobi [30]; however, our approach can be easily
extended to other solvers. Then, in Section 5.2 we report all the rounding
operations we made and the thresholds we used to mitigate numerical issues.

###  5.1 Implementing Benders decomposition in practice

There are two ways of implementing Benders decomposition: (1) using a cutting
plane procedure, also known as multi-tree approach, where each time a Benders
cut is generated, the cut is included in the master problem and the latter is
solved again as a MISOCP until optimality; or (2) using a branch-and-Benders-
cut approach, also known as single-tree approach, in which a single
enumeration tree for solving the MISOCP is created and Benders cuts are
separated at the nodes as in a classical branch-and-cut procedure. The
traditional approach is the multi-tree and it requires the solution of
possibly several master problems to optimality. The single-tree approach
instead requires the solution of only one master problem and has become very
popular in the recent literature (see, e.g., [18, 19, 34]). However, its
implementation was mainly applied in an MILP (rather than MISOCP) context. We
describe here the implementation of the two approaches.

Note that we do not use any specialized solution method to solve the master
problem or the subproblem as well. Indeed, using a general-purpose MISOCP
solver benefits many advantages, such as a better warm-start mechanism and
effective heuristics. This choice simplifies the implementation process
considerably.

####  5.1.1 Implementation of single-tree approach

For the single-tree approach, the common procedure when using a MIP solver is
to define a customized callback function to be called at every node of the
branch-and-cut applied to the master problem. The callback function, described
in Algorithm 1, defines how the separation is done. First, at each node of the
branching tree the continuous relaxation of the master problem is solved and
either a fractional or an integer solution $y$ is obtained. We restricted the
generation of a separating cut to integer solutions only for all nodes except
the root node. Specifically, at the root node we solve the subproblem to
possibly generate Benders cuts at any solution, integral and fractional, of
the master relaxation; at all the other nodes, we generate Benders cuts
exclusively when the solution of the master relaxation is integer. If at the
subproblem solution the target coverage is not reached, this solution must be
discarded. Therefore, we generate a normalized Benders cut, inject this cut
into the master problem and proceed through the branching tree until an
optimal solution is found or the maximum time is exceeded.

Algorithm 1 Single-tree Benders implementation

if (we are at the root node) OR (a MIP incumbent of the master problem is
found) then

extract the solution $(\bar{y},\bar{v})$ of the relaxed master problem

update the values of the master variables $(y,v)$ in the subproblem
constraints (10a)-(10b)

solve the subproblem

if (the target coverage is not reached) then

get the reduced costs associated with the master variables in the subproblem

add the normalized feasibility Benders cut (11) to the master problem

end if

end if

We observe that each variable-fixing equation (i.e., (10a)-(10b)) is meant to
be imposed by modifying the lower and upper bounds on the corresponding
variable, e.g., for the master variable $y$ we can impose $\bar{y}\leq
y\leq\bar{y}$, which can be handled very efficiently by the solver in a
preprocessing phase when the $\bar{y}$ is given.

We set the parameter timeLimit to stop the solution process when the maximum
time is reached. Moreover, in order to use a custom callback function, we set
the optimization parameters LazyConstraints to 1 and PreCrush to 1 for the
master problem. Note that not all commercial solvers allow you to apply a
callback function to a problem that contains quadratic constraints.

####  5.1.2 Implementation of multi-tree approach

The code of the multi-tree Benders decomposition method is reported in
Algorithm 2. It involves solving a master problem as an MISOCP iteratively,
updating the subproblem at every optimal solution of the master, and
generating Benders cuts until the master solution does not violate the
required coverage (we denote it as the termination condition in the algorithm)
or the process terminates due to exceeding of the time limit.

Algorithm 2 Multi-tree Benders implementation

while termination condition is not met or a maximum time is not exceeded do

solve the relaxed master problem as an MISOCP

extract the optimal solution $(\bar{y},\bar{v})$ of the master problem

update the values of the master variables $(y,v)$ in the subproblem
constraints (10a)-(10b)

solve the subproblem

if the target coverage is not reached then

get the reduced costs associated with the master variables in the subproblem

generate the normalized feasibility Benders cut

add it as a linear constraint to the master problem

else

termination condition is met

end if

end while

Also in this case, each variable-fixing equation (i.e., (10a)-(10b)) is
imposed by modifying the lower and upper bounds on the corresponding variable.

###  5.2 Thresholds and rounding operations

We now describe all the rounding operations we made and the thresholds we used
to deal with the numerical issues. When we get a MIP incumbent of the master
in the single-tree approach or an optimal solution of the master in the multi-
tree approach, we round to 0-1 the binary variables belonging to the master
solution. This is done even if the values should already be 0-1, as the solver
may lack precision. As for the fractional master solutions, we round to 0 all
the negative values of y and to 1 all the values of y exceeding 1 for the same
reason.

When applying the $\epsilon$-procedure, instead of replacing zero coefficients
with $\epsilon$, we replace the coefficients below a given small tolerance
$\text{tol}_{\beta}>0$ with $\epsilon$. This is because the solver is supposed
to treat values below some internal threshold as zeros.

Once the subproblem at a given master solution is solved, the condition we
should theoretically check to decide on the generation of the Benders cut is
if $\phi^{\prime}(\bar{y},\bar{v})<D$. To get a more precise check of this
condition from the numerical point of view, we implement this condition as
$\phi^{\prime}(\bar{y},\bar{v})/D<1-\text{tol}_{\alpha}$, where
$\text{tol}_{\alpha}>0$ is a given small tolerance taking two different
values, one for integer master solutions and another for fractional master
solutions. Indeed, for a correct implementation of the Benders decomposition
we only need to separate integral points. However, generating cuts also at
fractional solutions at the root node improves the quality of bounds and
reduces the size of the search space. Therefore, we decided to separate also
fractional solutions at the root node. Since this could sometimes slow down
the performance, to keep under control the number of Benders user cuts
generated, we set a higher value for the tolerance tolα, implying that their
generation is performed only if the coverage achieved at the fractional point
is far from the coverage required.

##  6 Computational Experiments

All the experiments run on a Ubuntu server with an Intel(R) Xeon(R) Gold 5218
CPU running at 2.30 GHz, with 96 GB of RAM and 8 cores. As for the optimizer,
we use Gurobi 10.0.3 with a time limit (TL) of 900 seconds. On Gurobi, to
guarantee better numerical precision, we set integrality tolerance IntFeasTol
to $10^{-9}$, feasibility tolerance FeasibilityTol to $10^{-9}$, the
optimality tolerances OptimalityTol to $10^{-9}$ and MIPGap to 0. As for the
setting of Benders decomposition methods in Gurobi, we enable the parameters
for lazy constraints (LazyConstraints=1) and user cuts (PreCrush=1) for the
single-tree approach to avoid any reductions and transformations that are
incompatible with a custom callback, and turn off presolve (Presolve=0) in all
the Benders approaches since during the presolve phase the model is not
complete as many constraints are missing and the optimizer could make
reductions without knowing the whole problem. We also set the following values
for the tolerances defined in the previous section:
$\text{tol}_{\alpha}=10^{-6}$ for integral points and
$\text{tol}_{\alpha}=0.5$ for fractional points, $\text{tol}_{\beta}=10^{-8}$,
$\epsilon=10^{-8}$. We also turn off presolve and use a dummy lazy callback on
the MISOCP model solved directly with Gurobi to have a fair comparison with
our Benders methods. To have a fairer comparison with our Benders, we turn off
presolve and use a dummy callback (i.e., an empty function declared as a
callback) to solve the MISOCP model via Gurobi blackbox. The use of a dummy
callback is typical in the MILP context to disable dynamic search and other
advanced features that solvers can apply when given a full model versus models
that require lazy cut separation.

The testbed is made of large-size instances derived from Cordeau et al. [14].
The original instances are tailored for realistic scenarios where the number
of customers is much larger than the number of potential facility locations
(i.e., $|J|>>|I|$). In particular, we selected instances with 1000 customers
and 100 potential facility locations. Customer demand in these instances was
generated by Cordeau et al. by uniformly sampling from the range [1,100] and
rounding to the nearest integer. The spatial coordinates for both customers
and facilities were randomly selected from the interval [0, 30].

We further adapted the testbed to suit our study on robust and congested
settings by injecting uncertainty in the customer demand and generating
congestion quadratic cost. Specifically, for each customer $j$, the deviation
$\hat{d}_{j}$ from the nominal demand $d_{j}$ was set as a random integer
drawn from the interval [$0,20\%\,d_{j}$], meaning that the maximum possible
deviation from the nominal value is 20%. Customers affected by uncertainty are
those for which $\hat{d}_{j}>0$. Protection from the uncertainty in customer
demand was introduced with varying levels of $\Gamma$. Congestion cost was
computed by setting $a_{i}=a$ and $b_{i}=b$ for each $i\in I$: the values of
$a$ and $b$ are chosen in such a way as to balance linear and quadratic costs,
and obtain optimal solutions requiring the opening of less than half of the
potential facilities.

###  6.1 Sensitivity analysis

In this paper we are studying a robust problem characterized by an unusual
feature: the worst realizations occur in the case of increasing demand for
what concerns the congestion, and in the opposite case of decreasing demand
for the coverage. Therefore, we find it interesting to investigate how the
optimal value, the coverage reached, the facility loads and the number of open
facilities vary based on the values of the protection level $\Gamma$ in each
of the following cases: the deterministic case, protection against load
uncertainty, protection against coverage uncertainty, and protection against
both load and coverage uncertainties.

The primary goal of this sensitivity analysis is to understand how changes in
the parameter $\Gamma$ influence the optimal solution of our problem.
Therefore, we solve eight different instances using Gurobi on the perspective
MISOCP reformulation of the problem, incrementally increasing the value of
$\Gamma$ in $\\{0,2.5\%,5\%,10\%,15\%,20\%,\\\ 30\%,40\%,$
$50\%,60\%,70\%,80\%,90\%,100\%\\}|J|$ and observing the resulting effects on
specific metrics. The instances are characterized by a different coverage
radius (R) and seed (s) used to generate the instance. In particular, for each
instance, starting from $\Gamma$ equal to zero, namely the deterministic
setting where we consider no protection towards the uncertainty, until
$\Gamma$ equal to the number of customers, namely the case with the maximum
protection towards the uncertainty, we report in Tables 1\- 2 the values of
the overall cost (ObjVal), the number of open facilities (#Fac), the facility
opening cost (OpenCost), the congestion cost (CongCost), the average load
(Load) and the coverage reached (Cov) at the optimal solution. Additionally,
Figures 1-3 graphically represent these metrics as functions of $\Gamma$ for a
representative instance.

From the results collected, we can claim that for increasing values of
$\Gamma$, we get increasing optimal values due to an increasing number of open
facilities (hence increasing opening cost) and an increasing average load at
the facilities (hence increasing congestion cost), until a certain value of
$\Gamma$ denoted in the figures as $\Gamma^{\star}$. For
$\Gamma\geq\Gamma^{\star}$, the values of these metrics remain stable, meaning
that over-conservative protection is not necessary as it does not affect the
optimal value. We can observe that $\Gamma^{\star}$ takes smaller values in
the case we only protect from the load uncertainty, higher values in the case
we only protect from the coverage uncertainty, intermediate values (around
$50\%|J|$) in the case we protect from both. Also note how the optimal
objective value in the double protection scenarios is closer to that of
coverage-only protection scenarios. As for the coverage reached, it remains
constant at the nominal coverage for the case of load uncertainty only, and it
is increasing with increasing values of $\Gamma$ for the other cases. We
remark that we are empirically overestimating $\Gamma$ (and $\Gamma^{\star}$
as well) as we do not consider all possible values of $\Gamma$.

Our findings indicate that the optimal solutions are highly sensitive to the
value of the protection level $\Gamma$, especially if we consider an uncertain
coverage. This understanding is crucial for making informed planning decisions
and assessing the reliability of a deterministic solution.

Table 1: Sensitivity analysis (part I). We denote by * the optimal value in
correspondence to $\Gamma^{\star}$.

|  | Uncertainty in Load & Coverage | Uncertainty in Coverage | Uncertainty in Load
---|---|---|---|---
(s-R) | $\bm{\Gamma}$ | ObjVal | #Fac | OpenCost | CongCost | Load | Cov | ObjVal | #Fac | OpenCost | CongCost | Load | Cov | ObjVal | #Fac | OpenCost | CongCost | Load | Cov
(1-5.5) | 0 | 26101.42 | 39 | 9880 | 16221.42 | 249.58 | 24958.00 | 26101.42 | 39 | 9880 | 16221.42 | 249.58 | 24958.00 | 26101.42 | 39 | 9880 | 16221.42 | 249.58 | 24958.00
25 | 27239.74 | 41 | 10740 | 16499.74 | 257.90 | 25097.36 | 26262.36 | 40 | 10280 | 15982.36 | 250.83 | 25082.95 | 27024.96 | 41 | 10740 | 16284.96 | 256.28 | 24958.00
50 | 27952.05 | 41 | 10820 | 17132.05 | 262.83 | 25176.00 | 26410.30 | 40 | 10310 | 16100.30 | 251.71 | 25170.58 | 27575.49 | 41 | 10740 | 16835.49 | 260.60 | 24958.00
100 | 28450.44 | 42 | 11190 | 17260.44 | 266.99 | 25367.02 | 26660.00 | 40 | 10370 | 16290.00 | 253.08 | 25308.27 | *27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
150 | 28680.39 | 42 | 11190 | 17490.39 | 268.76 | 25502.00 | 26849.04 | 40 | 10370 | 16479.04 | 254.51 | 25450.81 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
200 | 28898.38 | 42 | 11190 | 17708.38 | 270.39 | 25633.63 | 27021.09 | 40 | 10370 | 16651.09 | 255.81 | 25580.96 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
300 | 29228.79 | 42 | 11190 | 18038.79 | 272.83 | 25852.49 | 27297.09 | 41 | 10820 | 16477.09 | 257.56 | 25756.20 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
400 | 29477.27 | 43 | 11650 | 17827.27 | 274.43 | 26042.89 | 27540.93 | 41 | 10820 | 16720.93 | 259.57 | 25956.72 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
500 | 29634.08 | 43 | 11650 | 17984.08 | 275.64 | 26150.96 | 27746.09 | 41 | 10820 | 16926.09 | 261.25 | 26124.66 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
600 | *29675.66 | 43 | 11560 | 18115.66 | 276.59 | 26308.37 | 27889.45 | 41 | 10820 | 17069.45 | 262.37 | 26236.84 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
700 | 29675.66 | 43 | 11560 | 18115.66 | 276.59 | 26308.37 | *27954.82 | 42 | 11190 | 16764.82 | 263.15 | 26314.54 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
800 | 29675.66 | 43 | 11560 | 18115.66 | 276.59 | 26308.37 | 27954.82 | 42 | 11190 | 16764.82 | 263.15 | 26314.54 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
900 | 29675.66 | 43 | 11560 | 18115.66 | 276.59 | 26308.37 | 27954.82 | 42 | 11190 | 16764.82 | 263.15 | 26314.54 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
1000 | 29675.66 | 43 | 11560 | 18115.66 | 276.59 | 26308.37 | 27954.82 | 42 | 11190 | 16764.82 | 263.15 | 26314.54 | 27746.70 | 42 | 11110 | 16636.70 | 262.15 | 24958.00
(3-5.75) | 0 | 27046.64 | 34 | 9800 | 17246.64 | 240.46 | 24046.00 | 27046.64 | 34 | 9800 | 17246.64 | 240.46 | 24046.00 | 27046.64 | 34 | 9800 | 17246.64 | 240.46 | 24046.00
25 | 28018.51 | 35 | 10290 | 17728.51 | 247.32 | 24140.59 | 27142.45 | 35 | 10290 | 16852.45 | 241.12 | 24111.91 | 27860.37 | 35 | 10290 | 17570.37 | 246.22 | 24046.00
50 | 28683.83 | 36 | 10800 | 17883.83 | 251.90 | 24223.95 | 27234.65 | 35 | 10290 | 16944.65 | 241.77 | 24177.48 | 28355.61 | 36 | 10780 | 17575.61 | 249.73 | 24046.00
100 | 29264.81 | 36 | 10800 | 18464.81 | 255.98 | 24400.97 | 27413.54 | 35 | 10320 | 17093.54 | 242.81 | 24280.96 | *28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
150 | 29509.88 | 36 | 10800 | 18709.88 | 257.68 | 24552.56 | 27573.03 | 35 | 10320 | 17253.03 | 243.92 | 24391.61 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
200 | 29705.32 | 37 | 11290 | 18415.32 | 259.12 | 24661.42 | 27729.80 | 35 | 10320 | 17409.80 | 245.06 | 24505.88 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
300 | 30022.47 | 37 | 11290 | 18732.47 | 261.35 | 24866.54 | 28036.37 | 35 | 10320 | 17716.37 | 247.22 | 24722.26 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
400 | 30214.90 | 37 | 11290 | 18924.90 | 262.71 | 25007.61 | 28322.64 | 35 | 10320 | 18002.64 | 249.22 | 24922.34 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
500 | 30372.28 | 37 | 11290 | 19082.28 | 263.81 | 25136.47 | 28510.65 | 36 | 10800 | 17710.65 | 250.69 | 25068.79 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
600 | *30415.10 | 37 | 11290 | 19125.10 | 264.08 | 25226.96 | 28664.57 | 36 | 10800 | 17864.57 | 251.78 | 25178.33 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
700 | 30415.10 | 37 | 11290 | 19125.10 | 264.08 | 25226.96 | *28738.03 | 36 | 10800 | 17938.03 | 252.30 | 25230.16 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
800 | 30415.10 | 37 | 11290 | 19125.10 | 264.08 | 25226.96 | 28738.03 | 36 | 10800 | 17938.03 | 252.30 | 25230.16 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
900 | 30415.10 | 37 | 11290 | 19125.10 | 264.08 | 25226.96 | 28738.03 | 36 | 10800 | 17938.03 | 252.30 | 25230.16 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
1000 | 30415.10 | 37 | 11290 | 19125.10 | 264.08 | 25226.96 | 28738.03 | 36 | 10800 | 17938.03 | 252.30 | 25230.16 | 28554.81 | 36 | 10780 | 17774.81 | 251.14 | 24046.00
(4-6) | 0 | 27803.56 | 35 | 10350 | 17453.56 | 245.42 | 24541.50 | 27803.56 | 35 | 10350 | 17453.56 | 245.42 | 24541.50 | 27803.56 | 35 | 10350 | 17453.56 | 245.42 | 24541.50
25 | 28626.61 | 36 | 10840 | 17786.61 | 251.24 | 24628.51 | 27882.00 | 35 | 10350 | 17532.00 | 245.97 | 24596.95 | 28475.70 | 36 | 10840 | 17635.70 | 250.17 | 24541.50
50 | 29265.07 | 37 | 11330 | 17935.07 | 255.74 | 24697.81 | 27961.26 | 35 | 10350 | 17611.26 | 246.53 | 24652.80 | 28984.61 | 36 | 10840 | 18144.61 | 253.78 | 24541.50
100 | 29954.70 | 37 | 11330 | 18624.70 | 260.66 | 24873.50 | 28115.58 | 36 | 10840 | 17275.58 | 247.59 | 24758.61 | *29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
150 | 30204.67 | 37 | 11330 | 18874.67 | 262.42 | 25033.72 | 28270.32 | 36 | 10840 | 17430.32 | 248.70 | 24869.57 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
200 | 30403.79 | 38 | 11830 | 18573.79 | 263.78 | 25146.50 | 28426.84 | 36 | 10840 | 17586.84 | 249.81 | 24981.06 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
300 | 30723.60 | 38 | 11830 | 18893.60 | 266.05 | 25343.10 | 28726.17 | 36 | 10840 | 17886.17 | 251.93 | 25193.03 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
400 | 30923.64 | 38 | 11830 | 19093.64 | 267.47 | 25493.50 | 29015.64 | 36 | 10840 | 18175.64 | 253.98 | 25398.07 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
500 | 31076.66 | 38 | 11830 | 19246.66 | 268.55 | 25647.99 | 29231.95 | 36 | 10840 | 18391.95 | 255.52 | 25551.50 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
600 | *31088.18 | 38 | 11830 | 19258.18 | 268.63 | 25702.09 | 29383.90 | 37 | 11330 | 18053.90 | 256.61 | 25660.85 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
700 | 31088.18 | 38 | 11830 | 19258.18 | 268.63 | 25702.09 | *29442.15 | 37 | 11330 | 18112.15 | 257.03 | 25702.89 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
800 | 31088.18 | 38 | 11830 | 19258.18 | 268.63 | 25702.09 | 29442.15 | 37 | 11330 | 18112.15 | 257.03 | 25702.89 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
900 | 31088.18 | 38 | 11830 | 19258.18 | 268.63 | 25702.09 | 29442.15 | 37 | 11330 | 18112.15 | 257.03 | 25702.89 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
1000 | 31088.18 | 38 | 11830 | 19258.18 | 268.63 | 25702.09 | 29442.15 | 37 | 11330 | 18112.15 | 257.03 | 25702.89 | 29288.42 | 37 | 11330 | 17958.42 | 255.93 | 24541.50
(5-6.25) | 0 | 26409.11 | 39 | 9720 | 16689.11 | 253.18 | 25318.00 | 26409.11 | 39 | 9720 | 16689.11 | 253.18 | 25318.00 | 26409.11 | 39 | 9720 | 16689.11 | 253.18 | 25318.00
25 | 27346.55 | 40 | 10150 | 17196.55 | 260.24 | 25410.24 | 26492.39 | 39 | 9720 | 16772.39 | 253.82 | 25381.57 | 27203.98 | 39 | 9720 | 17483.98 | 259.16 | 25318.00
50 | 28034.17 | 40 | 10160 | 17874.17 | 265.35 | 25513.86 | 26576.66 | 39 | 9720 | 16856.66 | 254.46 | 25445.71 | 27703.10 | 40 | 10160 | 17543.10 | 262.88 | 25318.00
100 | 28716.61 | 41 | 10590 | 18126.61 | 270.55 | 25707.40 | 26748.20 | 39 | 9720 | 17028.20 | 255.76 | 25575.75 | 28008.73 | 40 | 10150 | 17858.73 | 265.27 | 25318.00
150 | 29000.45 | 41 | 10580 | 18420.45 | 272.73 | 25866.00 | 26923.89 | 39 | 9720 | 17203.89 | 257.08 | 25708.19 | *28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
200 | 29227.89 | 41 | 10580 | 18647.89 | 274.41 | 26005.37 | 27102.78 | 39 | 9720 | 17382.78 | 258.42 | 25842.25 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
300 | 29544.77 | 42 | 11130 | 18414.77 | 275.95 | 26176.28 | 27452.69 | 40 | 10150 | 17302.69 | 261.07 | 26107.28 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
400 | 29763.65 | 42 | 11140 | 18623.65 | 277.51 | 26352.13 | 27746.05 | 40 | 10150 | 17596.05 | 263.29 | 26328.77 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
500 | 29923.87 | 42 | 11140 | 18783.87 | 278.71 | 26467.47 | 27981.27 | 40 | 10150 | 17831.27 | 265.04 | 26503.68 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
600 | *30001.89 | 42 | 11140 | 18861.89 | 279.30 | 26623.75 | 28129.26 | 40 | 10150 | 17979.26 | 266.14 | 26614.09 | 28022.28 | 40 | 10150 | 17872.28 | 265.37 | 25318.00
700 | 30001.89 | 42 | 11140 | 18861.89 | 279.30 | 26623.75 | *28228.99 | 40 | 10150 | 18078.99 | 266.90 | 26690.21 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
800 | 30001.89 | 42 | 11140 | 18861.89 | 279.30 | 26623.75 | 28228.99 | 40 | 10150 | 18078.99 | 266.90 | 26690.21 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
900 | 30001.89 | 42 | 11140 | 18861.89 | 279.30 | 26623.75 | 28228.99 | 40 | 10150 | 18078.99 | 266.90 | 26690.21 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00
1000 | 30001.89 | 42 | 11140 | 18861.89 | 279.30 | 26623.75 | 28228.99 | 40 | 10150 | 18078.99 | 266.90 | 26690.21 | 28022.29 | 40 | 10150 | 17872.29 | 265.37 | 25318.00

Table 2: Sensitivity analysis (part II). We denote by * the optimal value in
correspondence to $\Gamma^{\star}$.

|  | Uncertainty in Load & Coverage | Uncertainty in Coverage | Uncertainty in Load
---|---|---|---|---
(s-R) | $\bm{\Gamma}$ | ObjVal | #Fac | OpenCost | CongCost | Load | Cov | ObjVal | #Fac | OpenCost | CongCost | Load | Cov | ObjVal | #Fac | OpenCost | CongCost | Load | Cov
(1-5.75) | 0 | 26101.42 | 39 | 9880 | 16221.42 | 249.58 | 24958.00 | 26101.42 | 39 | 9880 | 16221.42 | 249.58 | 24958.00 | 26101.42 | 39 | 9880 | 16221.42 | 249.58 | 24958.00
25 | 27162.26 | 41 | 10740 | 16422.26 | 257.30 | 25096.86 | 26254.34 | 40 | 10280 | 15974.34 | 250.77 | 25076.91 | 26953.31 | 40 | 10310 | 16643.31 | 255.95 | 24958.00
50 | 27856.05 | 41 | 10820 | 17036.05 | 262.10 | 25170.54 | 26396.61 | 40 | 10310 | 16086.61 | 251.60 | 25160.28 | 27490.00 | 41 | 10740 | 16750.00 | 259.94 | 24958.00
100 | 28406.33 | 42 | 11190 | 17216.33 | 266.64 | 25356.88 | 26642.23 | 40 | 10370 | 16272.23 | 252.97 | 25297.08 | *27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
150 | 28631.80 | 42 | 11190 | 17441.80 | 268.38 | 25498.00 | 26826.68 | 40 | 10370 | 16456.68 | 254.35 | 25435.08 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
200 | 28845.45 | 42 | 11190 | 17655.45 | 270.00 | 25621.57 | 26995.51 | 40 | 10370 | 16625.51 | 255.63 | 25563.49 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
300 | 29168.65 | 42 | 11190 | 17978.65 | 272.41 | 25841.94 | 27272.09 | 41 | 10820 | 16452.09 | 257.39 | 25738.81 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
400 | 29401.55 | 42 | 11190 | 18211.55 | 274.16 | 26036.76 | 27511.47 | 41 | 10820 | 16691.47 | 259.33 | 25932.75 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
500 | 29559.19 | 43 | 11560 | 17999.19 | 275.67 | 26177.64 | 27708.45 | 41 | 10820 | 16888.45 | 260.96 | 26096.40 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
600 | *29587.88 | 43 | 11560 | 18027.88 | 275.96 | 26277.24 | 27850.05 | 41 | 10820 | 17030.05 | 262.07 | 26207.05 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
700 | 29587.88 | 43 | 11560 | 18027.88 | 275.96 | 26277.25 | 27912.49 | 42 | 11190 | 16722.49 | 262.82 | 26282.42 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
800 | 29587.88 | 43 | 11560 | 18027.88 | 275.96 | 26277.24 | *27912.50 | 42 | 11190 | 16722.50 | 262.82 | 26282.42 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
900 | 29587.88 | 43 | 11560 | 18027.88 | 275.96 | 26277.25 | 27912.50 | 42 | 11190 | 16722.50 | 262.82 | 26282.42 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
1000 | 29587.88 | 43 | 11560 | 18027.88 | 275.96 | 26277.25 | 27912.50 | 42 | 11190 | 16722.50 | 262.82 | 26282.42 | 27710.68 | 42 | 11110 | 16600.68 | 261.88 | 24958.00
(3-5.5) | 0 | 27046.64 | 34 | 9800 | 17246.64 | 240.46 | 24046.00 | 27046.64 | 34 | 9800 | 17246.64 | 240.46 | 24046.00 | 27046.64 | 34 | 9800 | 17246.64 | 240.46 | 24046.00
25 | 28108.42 | 35 | 10290 | 17818.42 | 247.94 | 24148.01 | 27148.96 | 35 | 10290 | 16858.96 | 241.17 | 24116.58 | 27941.94 | 35 | 10290 | 17651.94 | 246.79 | 24046.00
50 | 28789.70 | 36 | 10800 | 17989.70 | 252.65 | 24238.32 | 27247.76 | 35 | 10290 | 16957.76 | 241.87 | 24186.80 | 28427.06 | 36 | 10780 | 17647.06 | 250.23 | 24046.00
100 | 29319.41 | 36 | 10810 | 18509.41 | 256.28 | 24413.69 | 27438.20 | 35 | 10320 | 17118.20 | 242.98 | 24298.11 | *28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
150 | 29563.69 | 37 | 11290 | 18273.69 | 258.11 | 24574.37 | 27609.19 | 35 | 10320 | 17289.19 | 244.14 | 24413.54 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
200 | 29766.99 | 37 | 11290 | 18476.99 | 259.54 | 24674.28 | 27769.77 | 35 | 10320 | 17449.77 | 245.32 | 24531.53 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
300 | 30089.87 | 37 | 11290 | 18799.87 | 261.80 | 24871.47 | 28084.11 | 35 | 10320 | 17764.11 | 247.55 | 24754.94 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
400 | 30285.64 | 37 | 11290 | 18995.64 | 263.18 | 25034.45 | 28374.23 | 35 | 10320 | 18054.23 | 249.57 | 24956.90 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
500 | 30451.85 | 37 | 11380 | 19071.85 | 263.75 | 25137.84 | 28560.08 | 36 | 10810 | 17750.08 | 250.95 | 25094.88 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
600 | *30486.15 | 37 | 11380 | 19106.15 | 263.97 | 25221.37 | 28713.89 | 36 | 10810 | 17903.89 | 252.05 | 25204.73 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
700 | 30486.15 | 37 | 11380 | 19106.15 | 263.97 | 25221.37 | *28782.36 | 36 | 10780 | 18002.36 | 252.74 | 25273.93 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
800 | 30486.15 | 37 | 11380 | 19106.15 | 263.97 | 25221.37 | 28782.36 | 36 | 10780 | 18002.36 | 252.74 | 25273.92 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
900 | 30486.15 | 37 | 11380 | 19106.15 | 263.97 | 25221.37 | 28782.36 | 36 | 10780 | 18002.36 | 252.74 | 25273.93 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
1000 | 30486.15 | 37 | 11380 | 19106.15 | 263.97 | 25221.37 | 28782.36 | 36 | 10780 | 18002.36 | 252.74 | 25273.92 | 28590.77 | 36 | 10780 | 17810.77 | 251.39 | 24046.00
(4-6.25) | 0 | 27803.56 | 35 | 10350 | 17453.56 | 245.42 | 24541.50 | 27803.56 | 35 | 10350 | 17453.56 | 245.42 | 24541.50 | 27803.56 | 35 | 10350 | 17453.56 | 245.42 | 24541.50
25 | 28566.13 | 36 | 10840 | 17726.13 | 250.81 | 24624.02 | 27880.53 | 35 | 10350 | 17530.53 | 245.96 | 24595.91 | 28423.85 | 36 | 10840 | 17583.85 | 249.80 | 24541.50
50 | 29186.63 | 37 | 11330 | 17856.63 | 255.18 | 24693.65 | 27958.29 | 35 | 10350 | 17608.29 | 246.51 | 24650.72 | 28921.03 | 36 | 10840 | 18081.03 | 253.33 | 24541.50
100 | 29926.45 | 37 | 11330 | 18596.45 | 260.46 | 24869.50 | 28109.22 | 36 | 10840 | 17269.22 | 247.54 | 24754.03 | 29279.70 | 37 | 11330 | 17949.70 | 255.87 | 24541.50
150 | 30190.79 | 37 | 11330 | 18860.79 | 262.32 | 25031.32 | 28260.64 | 36 | 10840 | 17420.64 | 248.63 | 24862.62 | *29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
200 | 30388.13 | 37 | 11330 | 19058.13 | 263.70 | 25145.50 | 28413.82 | 36 | 10840 | 17573.82 | 249.72 | 24972.25 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
300 | 30707.81 | 38 | 11830 | 18877.81 | 265.94 | 25339.55 | 28711.47 | 36 | 10840 | 17871.47 | 251.84 | 25184.28 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
400 | 30905.46 | 38 | 11830 | 19075.46 | 267.34 | 25489.77 | 29002.53 | 36 | 10840 | 18162.53 | 253.90 | 25389.67 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
500 | 31057.25 | 38 | 11830 | 19227.25 | 268.41 | 25644.22 | 29220.25 | 36 | 10840 | 18380.25 | 255.44 | 25543.67 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
600 | *31067.89 | 38 | 11830 | 19237.89 | 268.48 | 25694.96 | 29373.09 | 37 | 11330 | 18043.09 | 256.53 | 25653.11 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
700 | 31067.89 | 38 | 11830 | 19237.89 | 268.48 | 25694.96 | *29432.06 | 37 | 11330 | 18102.06 | 256.96 | 25695.68 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
800 | 31067.89 | 38 | 11830 | 19237.89 | 268.48 | 25694.96 | 29432.06 | 37 | 11330 | 18102.06 | 256.96 | 25695.68 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
900 | 31067.89 | 38 | 11830 | 19237.89 | 268.48 | 25694.96 | 29432.06 | 37 | 11330 | 18102.06 | 256.96 | 25695.68 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
1000 | 31067.89 | 38 | 11830 | 19237.89 | 268.48 | 25694.96 | 29432.06 | 37 | 11330 | 18102.06 | 256.96 | 25695.68 | 29279.72 | 37 | 11330 | 17949.72 | 255.87 | 24541.50
(5-6) | 0 | 26409.11 | 39 | 9720 | 16689.11 | 253.18 | 25318.00 | 26409.11 | 39 | 9720 | 16689.11 | 253.18 | 25318.00 | 26409.11 | 39 | 9720 | 16689.11 | 253.18 | 25318.00
25 | 27421.87 | 40 | 10150 | 17271.87 | 260.80 | 25414.51 | 26494.36 | 39 | 9720 | 16774.36 | 253.83 | 25383.07 | 27271.08 | 40 | 10150 | 17121.08 | 259.67 | 25318.00
50 | 28112.83 | 40 | 10160 | 17952.83 | 265.93 | 25523.00 | 26580.62 | 39 | 9720 | 16860.62 | 254.49 | 25448.73 | 27766.89 | 40 | 10160 | 17606.89 | 263.36 | 25318.00
100 | 28763.12 | 41 | 10590 | 18173.12 | 270.90 | 25726.17 | 26756.29 | 39 | 9720 | 17036.29 | 255.82 | 25581.87 | 28035.15 | 40 | 10150 | 17885.15 | 265.46 | 25318.00
150 | 29031.55 | 41 | 10690 | 18341.55 | 272.14 | 25858.38 | 26936.25 | 39 | 9720 | 17216.25 | 257.17 | 25717.48 | *28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
200 | 29253.31 | 41 | 10690 | 18563.31 | 273.79 | 25995.09 | 27119.67 | 39 | 9720 | 17399.67 | 258.55 | 25854.89 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
300 | 29577.26 | 41 | 10700 | 18877.26 | 276.09 | 26187.84 | 27478.95 | 40 | 10150 | 17328.95 | 261.27 | 26127.19 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
400 | 29797.02 | 42 | 11140 | 18657.02 | 277.75 | 26366.95 | 27775.56 | 40 | 10150 | 17625.56 | 263.51 | 26350.88 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
500 | 29957.95 | 42 | 11140 | 18817.95 | 278.95 | 26485.21 | 28012.94 | 40 | 10150 | 17862.94 | 265.27 | 26527.10 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
600 | *30036.75 | 42 | 11140 | 18896.75 | 279.49 | 26633.64 | 28161.54 | 40 | 10150 | 18011.54 | 266.38 | 26638.00 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
700 | 30036.75 | 42 | 11140 | 18896.75 | 279.49 | 26633.64 | *28250.86 | 40 | 10150 | 18100.86 | 267.06 | 26706.25 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
800 | 30036.75 | 42 | 11140 | 18896.75 | 279.49 | 26633.64 | 28250.86 | 40 | 10150 | 18100.86 | 267.06 | 26706.25 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
900 | 30036.75 | 42 | 11140 | 18896.75 | 279.49 | 26633.64 | 28250.86 | 40 | 10150 | 18100.86 | 267.06 | 26706.25 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00
1000 | 30036.75 | 42 | 11140 | 18896.75 | 279.49 | 26633.64 | 28250.86 | 40 | 10150 | 18100.86 | 267.06 | 26706.25 | 28040.16 | 40 | 10150 | 17890.16 | 265.50 | 25318.00

Figure 1: Plot of the variation of the overall objective value (on the left)
and the number of open facilities (on the right) based on $\Gamma$ values for
an instance having 1000 customers, seed 1 and coverage radius 5.5.

Figure 2: Plot of the variation of the facility opening cost (on the left) and
the congestion cost (on the right) based on $\Gamma$ values for an instance
having 1000 customers, seed 1 and coverage radius 5.5.

Figure 3: Plot of the variation of the average load (on the left) and the
coverage reached (on the right) based on $\Gamma$ values for an instance
having 1000 customers, seed 1 and coverage radius 5.5.

###  6.2 Performance analysis

We use a set of 128 instances of the CPSCLP to compare the performance of
different methods and formulations introduced in this paper. The
characteristics of each instance, including its optimal value and the number
of open facilities at the optimal solution, are reported in Table 3 of the
Appendix B. We use the following values for the protection level $\Gamma$:
$2.5\%|J|,5\%|J|,10\%|J|,15\%|J|,20\%|J|,30\%|J|,40\%|J|,50\%|J|$. We limit
$\Gamma$ at 50% of the number of customers to avoid overly conservative
models.

The purpose of this section is to evaluate the performance of our Benders
approaches as possible competitive exact algorithms to solve robust instances
of the CPSCLP. To the best of our knowledge, no previous computational study
appeared in the literature for such a problem, despite its potential
theoretical and practical relevance. For this reason, in this section we
compare the performance of our Benders algorithms against the direct use of
Gurobi applied as a black-box MISOCP solver to the perspective reformulation
of the robust counterpart of the problem. We also test Gurobi on the extended
MIQP formulation of the robust counterpart of the problem to show the
effectiveness of the perspective reformulation. We choose to compare our
Benders algorithms to Gurobi as it is considered one of the state-of-the-art
MIP solvers and is usually used as benchmark to compare the performance of
newly developed exact algorithms. Moreover, Gurobi allows the implementation
of callbacks on problems with quadratic constraints, a necessary feature for
the implementation of our single-tree approaches.

We discuss the results given by six exact algorithms: Gurobi on the extended
MIQP formulation (7) (denoted with MIQP from now on), Gurobi on the
perspective MISOCP reformulation (8) (denoted with MISOCP from now on), and
four versions of our Benders algorithms directly on the perspective
reformulation: the single-tree version (ST-BEN) and the single-tree version
using the $\epsilon$-technique (ST$\epsilon$-BEN), the multi-tree version (MT-
BEN) and the multi-tree version using the $\epsilon$-technique
(MT$\epsilon$-BEN). Numerical results are reported in Tables 4-6 of the
Appendix B. Using the values reported in these tables, we made a graphical
representation of the overall performance of the different exact algorithms in
Figure 4, illustrating the percentage of problems solved by each algorithm as
a function of the computational time. It is important to note that these
solution profiles do not represent the cumulative solution time but show the
percentage of problems the algorithm can solve within a certain amount of
time. The best performance are graphically represented by the curves in the
upper part of the plot. From this representation, we can also infer the
percentage of instances solved by each algorithm.

Figure 4: Solution profiles considering the whole benchmark set of 128
instances.

In particular, the profiles show that the MIQP model behaves poorly for this
problem as none of the instances can be solved by this method and the relative
gaps at the end of the optimization are pretty high (all over 50%). The second
and third worst performing methods for what concerns the number of solved
instances are the Benders approaches without the $\epsilon$-technique. Indeed,
ST-BEN solved 41 instances (around 32%) and MT-BEN solved 20 instances (around
16%). The Benders methods with $\epsilon$-technique, instead, solved most of
the instances: ST$\epsilon$-BEN solved 125 instances (around 98%) and
MT$\epsilon$-BEN solved 120 instances (around 94%). Finally, MISOCP solved all
128 instances.

As for resolution times, the algorithms running on the perspective
reformulation are much more faster than the one running on the extended
formulation: indeed, we cannot solve any instance within the time limit using
the MIQP model. The fastest methods are ST$\epsilon$-BEN and MT$\epsilon$-BEN
in a large percentage of the instances considered (around 97%), and MISOCP in
the remaining ones. With ST$\epsilon$-BEN, times of MISOCP are reduced on
average by almost 67%, with MT$\epsilon$-BEN, times of MISOCP are reduced on
average by almost 75%. We also remark that the two Benders methods using the
$\epsilon$-technique converge to the same optimal solution of MISOCP for the
chosen value of $\epsilon$.

For what concerns the direct comparison between ST$\epsilon$-BEN and
MT$\epsilon$-BEN, we can notice that both approaches are very effective, with
MT$\epsilon$-BEN being slightly superior in terms of solution times, and
ST$\epsilon$-BEN being slightly superior in terms of the number of solved
instances.

From Figure 5 we further compare the Benders approaches with respect to the
number of cuts generated (on integral and fractional points for the single-
tree approach, only at integral points for the multi-tree approach). From the
comparison of the box plots we can assess that the generation of cuts is
reduced in both single and multi-tree approaches when using the
$\epsilon$-technique, especially for the single-tree. The number of Benders
cuts generated on average by the multi-tree approach is reduced by the
$\epsilon$-technique by almost 75%, and in the single-tree approach by 83%.
Furthermore, using the $\epsilon$-technique the number of cuts generated
remains stable, while when it is not used the number of cuts reaches high
values (up to 522 cuts for ST-BEN) in several instances.

Figure 5: Boxplot of the number of Benders cuts generated by Benders
decomposition methods.

In Figures 6-7 we report the solution profiles for each fixed value of
$\Gamma$ to assess how $\Gamma$ affects the performance of each solution
algorithm.

Figure 6: Solution profiles for low values of $\Gamma$.

Figure 7: Solution profiles for medium to high values of $\Gamma$.

From these figures, ST$\epsilon$-BEN and MT$\epsilon$-BEN consistently emerge
as the best-performing methods across all tested values of $\Gamma$.
Specifically, MT$\epsilon$-BEN significantly outperforms ST$\epsilon$-BEN for
very low values of $\Gamma$ ($\Gamma=\\{2.5\%|J|,5\%|J|\\}$). As for Benders
methods that do not employ the $\epsilon$-technique, they exhibit good
performance at low $\Gamma$ values. Specifically, ST-BEN outperforms MISOCP in
more than 80% of instances when $\Gamma=2.5\%|J|$ and in over 40% of instances
when $\Gamma=5\%|J|$. Similarly, MT-BEN beats MISOCP in 70% of cases having
$\Gamma=2.5\%|J|$ and in nearly 30% when $\Gamma=5\%|J|$. However, their
performance drastically declines as $\Gamma$ increases, indicating a high
limitation in their scalability with respect to increasing levels of
protection against uncertainty.

##  7 Conclusions

We presented a novel investigation into solving a robust and convex quadratic
variant of the partial set covering location problem using modern Benders
decomposition. Our new model achieves a better resource allocation and more
balanced solutions to facility location problems affected by congestion and
demand uncertainty. We investigated the case with demand uncertainty affecting
both feasibility and optimality using $\Gamma$-robustness. We showed how a
simple implementation of Benders decomposition is successful on a perspective
reformulation of the robust counterpart of the problem, describing in detail
all implementation ingredients we used to design an effective code embedded in
a modern MISOCP solver. Even though Benders decomposition is a traditional
technique, the innovation of our proposal consists in being one of the first
to use callbacks in combination with quadratic constraints, a new option
provided by (few) state-of-the-art MISOCP solvers. Furthermore, we were able
to demonstrate the effectiveness of Benders decomposition, even though the
Benders subproblem is non-separable. We compared single and multi-tree
implementations of Benders decomposition with the default branch-and-cut of a
state-of-the-art MISOCP solver, also testing the use of a novel ad hoc cut-
strengthening technique for degenerate subproblems. Using this technique we
were able to reduce computational times of Gurobi on average by 75% with the
multi-tree and 67% with the single-tree procedure on 128 large-size instances
derived from the literature. We also conducted a sensitivity analysis to
understand how changes in the protection level $\Gamma$ influence the optimal
solutions of the instances. Our findings indicate that $\Gamma$ has a
significant impact on the overall cost, the number of facilities to open, the
coverage reached and the facility load, suggesting the unreliability of a
deterministic setting of our problem and the importance of planning decisions
that are robust towards the uncertainty.

## References

  * [1] R. Aboolian, O. Berman, and D. Krass.  Profit maximizing distributed service system design with congestion and elastic demand.  Transportation Science, 46(2):247–261, 2012.
  * [2] M. S. Aktürk, A. Atamtürk, and S. Gürel.  A strong conic quadratic reformulation for machine-job assignment with controllable processing times.  Operations Research Letters, 37(3):187–191, 2009.
  * [3] S. Atta.  An improved harmony search algorithm using opposition-based learning and local search for solving the maximal covering location problem.  Engineering Optimization, pages 1–20, 2023.
  * [4] S. Atta, P. R. Sinha Mahapatra, and A. Mukhopadhyay.  Solving maximal covering location problem using genetic algorithm with local refinement.  Soft Computing, 22:3891–3906, 2018.
  * [5] P. Avella, A. Calamita, and L. Palagi.  A computational study of off-the-shelf minlp solvers on a benchmark set of congested capacitated facility location problems.  arXiv preprint arXiv:2303.04216, 2023.
  * [6] J. F. Benders.  Partitioning procedures for solving mixed-variables programming problems.  Numerische Mathematik, 4(1):238–252, 1962.
  * [7] O. Berman, D. Krass, and Z. Drezner.  The gradual covering decay location problem on a network.  European Journal of Operational Research, 151(3):474–480, 2003.
  * [8] D. Bertsimas and M. Sim.  Robust discrete optimization and network flows.  Mathematical programming, 98(1):49–71, 2003.
  * [9] N. Bilal, P. Galinier, and F. Guibault.  An iterated-tabu-search heuristic for a variant of the partial set covering problem.  Journal of Heuristics, 20(2):143–164, 2014.
  * [10] B. Boffey, R. Galvao, and L. Espejo.  A review of congestion models in the location of facilities with immobile servers.  European Journal of Operational Research, 178(3):643–662, 2007.
  * [11] L. Chen, S.-J. Chen, W.-K. Chen, Y.-H. Dai, T. Quan, and J. Chen.  Efficient presolving methods for solving maximal covering and partial set covering location problems.  European Journal of Operational Research, 2023.
  * [12] T. R. L. Christensen and A. Klose.  A fast exact method for the capacitated facility location problem with differentiable convex production costs.  European Journal of Operational Research, 292(3):855–868, 2021.
  * [13] R. Church and C. ReVelle.  The maximal covering location problem.  In Papers of the regional science association, volume 32, pages 101–118. Springer-Verlag Berlin/Heidelberg, 1974.
  * [14] J.-F. Cordeau, F. Furini, and I. Ljubić.  Benders decomposition for very large scale partial set covering and maximal covering location problems.  European Journal of Operational Research, 275(3):882–896, 2019.
  * [15] M. S. Daskin and S. H. Owen.  Two new location covering problems: The partial p-center problem and the partial set covering problem.  Geographical Analysis, 31(3):217–235, 1999.
  * [16] M. Desrochers, P. Marcotte, and M. Stan.  The congested facility location problem.  Location Science, 3(1):9–23, 1995.
  * [17] B. T. Downs and J. D. Camm.  An exact algorithm for the maximal covering problem.  Naval Research Logistics (NRL), 43(3):435–461, 1996.
  * [18] M. Fischetti, I. Ljubić, and M. Sinnl.  Benders decomposition without separability: A computational study for capacitated facility location problems.  European Journal of Operational Research, 253(3):557–569, 2016.
  * [19] M. Fischetti, I. Ljubić, and M. Sinnl.  Redesigning Benders decomposition for large-scale facility location.  Management Science, 63(7):2146–2162, 2017.
  * [20] M. Fischetti, D. Salvagnin, and A. Zanette.  A note on the selection of Benders’ cuts.  Mathematical Programming, 124:175–182, 2010.
  * [21] A. Frangioni and C. Gentile.  Perspective cuts for a class of convex 0–1 mixed integer programs.  Mathematical Programming, 106:225–236, 2006.
  * [22] A. Frangioni and C. Gentile.  Sdp diagonalizations and perspective cuts for a class of nonseparable miqp.  Operations Research Letters, 35(2):181–185, 2007.
  * [23] A. Frangioni and C. Gentile.  A computational comparison of reformulations of the perspective relaxation: SOCP vs. cutting planes.  Operations Research Letters, 37(3):206–210, 2009.
  * [24] A. Frangioni, C. Gentile, and F. Lacalandra.  Tighter approximated milp formulations for unit commitment problems.  IEEE Transactions on Power Systems, 24(1):105–113, 2008.
  * [25] R. D. Galvão, L. G. A. Espejo, and B. Boffey.  A comparison of lagrangean and surrogate relaxations for the maximal covering location problem.  European Journal of Operational Research, 124(2):377–389, 2000.
  * [26] R. D. Galvão and C. ReVelle.  A lagrangean heuristic for the maximal covering location problem.  European Journal of Operational Research, 88(1):114–123, 1996.
  * [27] O. Günlük and J. Linderoth.  Perspective relaxation of mixed integer nonlinear programs with indicator variables.  In International Conference on Integer Programming and Combinatorial Optimization, pages 1–16. Springer, 2008.
  * [28] O. Günlük and J. Linderoth.  Perspective reformulation and applications.  In Mixed Integer Nonlinear Programming, pages 61–89. Springer, 2012.
  * [29] Gurobi Optimization.  MIQCPMethod, January 2023.  https://www.gurobi.com/documentation/10.0/refman/miqcpmethod.html.
  * [30] Gurobi Optimization, LLC.  Gurobi Optimizer Reference Manual, 2023.  https://www.gurobi.com.
  * [31] S. L. Hakimi.  Optimum distribution of switching centers in a communication network and some related graph theoretic problems.  Operations Research, 13(3):462–475, 1965.
  * [32] J. Harkness and C. ReVelle.  Facility location with increasing production costs.  European Journal of Operational Research, 145(1):1–13, 2003.
  * [33] J. Krarup, M. Labbé, and G. Rand.  Optimal location: minimum versus equilibrium allocation.  In Operational research’87, pages 718–729. North-Holland, 1987.
  * [34] I. Ljubić, P. Putz, and J.-J. Salazar-González.  Exact approaches to the single-source network loading problem.  Networks, 59(1):89–106, 2012.
  * [35] L. A. Lorena and M. A. Pereira.  A lagrangean/surrogate heuristic for the maximal covering location problem using hillman’s edition.  International Journal of Industrial Engineering, 9:57–67, 2002.
  * [36] D. Lu.  Facility location with economies of scale and congestion.  Master’s thesis, University of Waterloo, 2010.
  * [37] T. L. Magnanti and R. T. Wong.  Accelerating Benders decomposition: Algorithmic enhancement and model selection criteria.  Operations Research, 29(3):464–484, 1981.
  * [38] V. R. Máximo, J.-F. Cordeau, and M. C. Nascimento.  A hybrid adaptive iterated local search heuristic for the maximal covering location problem.  International Transactions in Operational Research, 2023.
  * [39] V. R. Máximo and M. C. Nascimento.  Intensification, learning and diversification in a hybrid metaheuristic: an efficient unification.  Journal of Heuristics, 25:539–564, 2019.
  * [40] V. R. Máximo, M. C. Nascimento, and A. C. Carvalho.  Intelligent-guided adaptive search for the maximum covering location problem.  Computers & Operations Research, 78:129–137, 2017.
  * [41] N. Megiddo, E. Zemel, and S. L. Hakimi.  The maximum coverage location problem.  SIAM Journal on Algebraic Discrete Methods, 4(2):253–261, 1983.
  * [42] N. Papadakos.  Practical enhancements to the Magnanti–Wong method.  Operations Research Letters, 36(4):444–449, 2008.
  * [43] R. Rahmaniani, T. G. Crainic, M. Gendreau, and W. Rei.  The benders decomposition algorithm: A literature review.  European Journal of Operational Research, 259(3):801–817, 2017.
  * [44] M. G. Resende.  Computing approximate solutions of the maximum covering problem with grasp.  Journal of Heuristics, 4(2):161–177, 1998.
  * [45] C. ReVelle, M. Scholssberg, and J. Williams.  Solving the maximal covering location problem with heuristic concentration.  Computers & Operations Research, 35(2):427–435, 2008.
  * [46] S. Şelfun.  Outer approximation algorithms for the congested p-median problem.  Master’s thesis, Bilkent Universitesi (Turkey), 2011.
  * [47] E. L. F. Senne, M. A. Pereira, L. A. N. Lorena, et al.  A decomposition heuristic for the maximal covering location problem.  Advances in Operations Research, 2010, 2010.
  * [48] L. V. Snyder.  Facility location under uncertainty: a review.  IIE transactions, 38(7):547–564, 2006.
  * [49] L. V. Snyder.  Covering problems.  Foundations of location analysis, pages 109–135, 2011.
  * [50] C. Toregas, R. Swain, C. ReVelle, and L. Bergman.  The location of emergency service facilities.  Operations Research, 19(6):1363–1373, 1971.
  * [51] W. Walker.  Using the set-covering problem to assign fire companies to fire houses.  Operations Research, 22(2):275–277, 1974.
  * [52] Y. Zhang, O. Berman, and V. Verter.  Incorporating congestion in preventive healthcare facility network design.  European Journal of Operational Research, 198(3):922–935, 2009.

## Acknowledgements

Alice Calamita acknowledges financial support from Progetto di Ricerca Medio
Sapienza Uniroma1 - n. RM12117A857F64F2. Laura Palagi acknowledges financial
support from Progetto di Ricerca Medio Sapienza Uniroma1 - n.
RM1221816BAE8A79.

##  Appendix A Appendix

We provide the proof to the Theorem 2.1, following what is done in Theorem 1
of [8].

###### Theorem A.1

The extended formulation (7) is equivalent to (6).

  * Proof

We show how to get formulation (7).

Given $i\in I$ and a vector $x_{i}^{\star}\in\mathbb{R}^{|J|}$, we define the
protection function $\alpha$ representing the sum of the $\Gamma$ largest
deviations for $i\in I$ as

$\alpha(x_{i}^{\star})=\max_{\\{S:\ S\subseteq
J,|S|\leq\Gamma\\}}\left\\{\sum_{j\in S}\hat{d}_{j}x^{\star}_{ij}\right\\}.$

By introducing the binary variables $w_{ij}$, the protection function $\alpha$
can be equivalently formulated as

$\begin{array}[]{rlr}\alpha(x_{i}^{\star})=\max&\displaystyle\sum_{j\in
J}\hat{d}_{j}x^{\star}_{ij}w_{ij}\\\ s.t.&\displaystyle\sum_{j\in
J}w_{ij}\leq\Gamma\\\ &w_{ij}\in\\{0,1\\}&\quad j\in J\end{array}$

which equals to

$\begin{array}[]{rlr}\alpha(x_{i}^{\star})=\max&\displaystyle\sum_{j\in
J}\hat{d}_{j}x^{\star}_{ij}w_{ij}\\\ s.t.&\displaystyle\sum_{j\in
J}w_{ij}\leq\Gamma\\\ &0\leq w_{ij}\leq 1&\quad j\in J\end{array}$ (17)

where we relaxed the integrality of variables $w_{ij}$ since the matrix of the
constraints is totally unimodular and $\Gamma$ is integral (linear programs of
this form have integral optima). Clearly the optimal solution of Problem (17)
consists of $\Gamma$ variables at 1. Associating multipliers $\rho_{i}$ to the
knapsack constraint and $\sigma_{ij}$ for $j\in J$ to the upper bound
constraints, we consider the dual of Problem (17)

$\begin{array}[]{rlr}\min&\Gamma\rho_{i}+\displaystyle\sum_{j\in
J}\sigma_{ij}\\\
s.t.&\rho_{i}+\sigma_{ij}\geq\displaystyle\hat{d}_{j}x_{ij}^{\star}&\quad j\in
J\\\ &\sigma_{ij}\geq 0&\quad j\in J\\\ &\rho_{i}\geq 0.\end{array}$ (18)

By strong duality, since Problem (17) is feasible and bounded for all integer
$\Gamma\in[0,|J|]$, then the dual problem (18) is also feasible and bounded
and their objective functions assume the same value in an optimal solution.
Hence, $\alpha(x_{i}^{\star})$ is equal to the optimal value of Problem (18)
for each $i\in I$. Thus constraints (6b) can be equivalently formulated as

$\begin{array}[]{cc}\displaystyle v_{i}-\sum_{j\in
J}d_{j}x_{ij}-\left(\Gamma\rho_{i}+\sum_{j\in J}\sigma_{ij}\right)\geq 0&i\in
I\\\ \rho_{i}+\sigma_{ij}\geq\displaystyle\hat{d}_{j}x_{ij}&i\in I,\ j\in J\\\
\sigma_{ij}\geq 0&i\in I,\ j\in J\\\ \rho_{i}\geq 0&i\in I.\end{array}$

Then, given vectors $x_{i}^{\star}$ for all $i\in I$, we define the total
protection function $\beta$ representing the sum of the $\Gamma$ largest
deviations as:

$\beta(x^{\star})=\max_{\\{S:\ S\subseteq
J,|S|\leq\Gamma\\}}\left\\{\sum_{i\in I}\sum_{j\in
S}\hat{d}_{j}x^{\star}_{ij}\right\\}.$

By introducing the binary variables $w_{ij}$, the protection function $\beta$
can be equivalently formulated as

$\begin{array}[]{rlr}\beta(x^{\star})=\max&\displaystyle\sum_{i\in
I}\sum_{j\in J}\hat{d}_{j}x^{\star}_{ij}w_{j}\\\ s.t.&\displaystyle\sum_{j\in
J}w_{j}\leq\Gamma\\\ &w_{j}\in\\{0,1\\}&\quad j\in J\end{array}$

which equals to

$\begin{array}[]{rlr}\beta(x^{\star})=\max&\displaystyle\sum_{i\in
I}\sum_{j\in J}\hat{d}_{j}x^{\star}_{ij}w_{j}\\\ s.t.&\displaystyle\sum_{j\in
J}w_{j}\leq\Gamma\\\ &0\leq w_{j}\leq 1&\quad j\in J\end{array}$ (19)

where we relaxed the integrality of variables $w_{ij}$ since the matrix of the
constraints is totally unimodular and $\Gamma$ is integral (linear programs of
this form have integral optima). Clearly the optimal solution of Problem (19)
consists of $\Gamma$ variables at 1. We next consider the dual of Problem
(19). Using the multipliers $\tau\in\mathbb{R}$ for the knapsack constraint
and $\pi\in\mathbb{R}^{|J|}$ for the upper bound constraints, we get

$\begin{array}[]{rlr}\min&\Gamma\tau+\displaystyle\sum_{j\in J}\pi_{j}\\\
s.t.&\tau+\pi_{j}\geq\displaystyle\sum_{i\in I}\hat{d}_{j}x_{ij}^{\star}&\quad
j\in J\\\ &\pi_{j}\geq 0&\quad j\in J\\\ &\tau\geq 0.\end{array}$ (20)

By strong duality, since Problem (19) is feasible and bounded for all
$\Gamma\in[0,|J|]$, then the dual problem (20) is also feasible and bounded
and their objective functions assume the same value in an optimal solution.
Hence, $\beta(x^{\star})$ is equal to the optimal value of Problem (20). Thus
constraint (6c) can be equivalently formulated as

$\begin{array}[]{lr}\displaystyle\sum_{i\in I}\sum_{j\in
J}d_{j}x_{ij}-\left(\Gamma\tau+\sum_{j\in J}\pi_{j}\right)\geq D&\\\
\tau+\pi_{j}\geq\displaystyle\sum_{i\in I}\hat{d}_{j}x_{ij}&\quad j\in J\\\
\pi_{j}\geq 0&\quad j\in J\\\ \tau\geq 0.\end{array}$

We finally get formulation (7).

$\square$

##  Appendix B Appendix

Tables 3 report the characteristics of the instances of congested partial set
covering location. Specifically, for each instance, the following data are
reported: a unique identifier ID, a number s linked to the seed of the
instance generator, the maximum number of customer deviations $\Gamma$ we
consider for the robust setting, the maximum distance of coverage R, the
optimal value for instance at hand ObjVal, and the number of open facilities
#Fac at the best feasible solution found on the instance.

Table 3: Characteristics of the testbed.

ID | s | $\bm{\Gamma}$ | R | ObjVal | #Fac |  | ID | s | $\bm{\Gamma}$ | R | ObjVal | #Fac |  | ID | s | $\bm{\Gamma}$ | R | ObjVal | #Fac
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
1 | 1 | 25 | 5.5 | 27239.74 | 41 |  | 43 | 3 | 100 | 6 | 29245.32 | 36 |  | 86 | 4 | 300 | 5.75 | 30739.45 | 38
2 | 1 | 25 | 5.75 | 27162.26 | 41 |  | 44 | 3 | 100 | 6.25 | 29197.75 | 36 |  | 87 | 4 | 300 | 6 | 30723.60 | 38
3 | 1 | 25 | 6 | 27080.35 | 41 |  | 45 | 3 | 150 | 5.5 | 29563.69 | 37 |  | 88 | 4 | 300 | 6.25 | 30707.81 | 38
4 | 1 | 25 | 6.25 | 27005.71 | 41 |  | 46 | 3 | 150 | 5.75 | 29509.88 | 36 |  | 89 | 4 | 400 | 5.5 | 30968.95 | 38
5 | 1 | 50 | 5.5 | 27952.05 | 41 |  | 47 | 3 | 150 | 6 | 29497.71 | 36 |  | 90 | 4 | 400 | 5.75 | 30941.75 | 38
6 | 1 | 50 | 5.75 | 27856.05 | 41 |  | 48 | 3 | 150 | 6.25 | 29466.02 | 37 |  | 91 | 4 | 400 | 6 | 30923.64 | 38
7 | 1 | 50 | 6 | 27769.04 | 41 |  | 49 | 3 | 200 | 5.5 | 29766.99 | 37 |  | 92 | 4 | 400 | 6.25 | 30905.46 | 38
8 | 1 | 50 | 6.25 | 27669.14 | 41 |  | 50 | 3 | 200 | 5.75 | 29705.32 | 37 |  | 93 | 4 | 500 | 5.5 | 31119.74 | 39
9 | 1 | 100 | 5.5 | 28450.44 | 42 |  | 51 | 3 | 200 | 6 | 29692.82 | 37 |  | 94 | 4 | 500 | 5.75 | 31094.42 | 38
10 | 1 | 100 | 5.75 | 28406.33 | 42 |  | 52 | 3 | 200 | 6.25 | 29660.42 | 37 |  | 95 | 4 | 500 | 6 | 31076.66 | 38
11 | 1 | 100 | 6 | 28368.57 | 42 |  | 53 | 3 | 300 | 5.5 | 30089.87 | 37 |  | 96 | 4 | 500 | 6.25 | 31057.25 | 38
12 | 1 | 100 | 6.25 | 28316.45 | 42 |  | 54 | 3 | 300 | 5.75 | 30022.47 | 37 |  | 97 | 5 | 25 | 5.5 | 27614.20 | 40
13 | 1 | 150 | 5.5 | 28680.39 | 42 |  | 55 | 3 | 300 | 6 | 30010.20 | 37 |  | 98 | 5 | 25 | 5.75 | 27515.77 | 40
14 | 1 | 150 | 5.75 | 28631.80 | 42 |  | 56 | 3 | 300 | 6.25 | 29979.06 | 37 |  | 99 | 5 | 25 | 6 | 27421.87 | 40
15 | 1 | 150 | 6 | 28597.00 | 42 |  | 57 | 3 | 400 | 5.5 | 30285.64 | 37 |  | 100 | 5 | 25 | 6.25 | 27346.55 | 40
16 | 1 | 150 | 6.25 | 28557.92 | 42 |  | 58 | 3 | 400 | 5.75 | 30214.90 | 37 |  | 101 | 5 | 50 | 5.5 | 28300.76 | 41
17 | 1 | 200 | 5.5 | 28898.38 | 42 |  | 59 | 3 | 400 | 6 | 30202.67 | 37 |  | 102 | 5 | 50 | 5.75 | 28207.67 | 40
18 | 1 | 200 | 5.75 | 28845.45 | 42 |  | 60 | 3 | 400 | 6.25 | 30169.95 | 37 |  | 103 | 5 | 50 | 6 | 28112.83 | 40
19 | 1 | 200 | 6 | 28807.26 | 42 |  | 61 | 3 | 500 | 5.5 | 30450.82 | 37 |  | 104 | 5 | 50 | 6.25 | 28034.17 | 40
20 | 1 | 200 | 6.25 | 28763.63 | 43 |  | 62 | 3 | 500 | 5.75 | 30372.28 | 37 |  | 105 | 5 | 100 | 5.5 | 28827.47 | 41
21 | 1 | 300 | 5.5 | 29228.79 | 42 |  | 63 | 3 | 500 | 6 | 30356.26 | 37 |  | 106 | 5 | 100 | 5.75 | 28801.14 | 41
22 | 1 | 300 | 5.75 | 29168.65 | 42 |  | 64 | 3 | 500 | 6.25 | 30321.74 | 37 |  | 107 | 5 | 100 | 6 | 28763.12 | 41
23 | 1 | 300 | 6 | 29123.02 | 42 |  | 65 | 4 | 25 | 5.5 | 28769.80 | 36 |  | 108 | 5 | 100 | 6.25 | 28716.61 | 41
24 | 1 | 300 | 6.25 | 29072.71 | 43 |  | 66 | 4 | 25 | 5.75 | 28688.22 | 36 |  | 109 | 5 | 150 | 5.5 | 29086.06 | 41
25 | 1 | 400 | 5.5 | 29477.27 | 43 |  | 67 | 4 | 25 | 6 | 28626.61 | 36 |  | 110 | 5 | 150 | 5.75 | 29060.88 | 41
26 | 1 | 400 | 5.75 | 29401.55 | 42 |  | 68 | 4 | 25 | 6.25 | 28566.13 | 36 |  | 111 | 5 | 150 | 6 | 29031.55 | 41
27 | 1 | 400 | 6 | 29350.02 | 43 |  | 69 | 4 | 50 | 5.5 | 29455.33 | 37 |  | 112 | 5 | 150 | 6.25 | 29000.45 | 41
28 | 1 | 400 | 6.25 | 29295.53 | 43 |  | 70 | 4 | 50 | 5.75 | 29349.76 | 37 |  | 113 | 5 | 200 | 5.5 | 29310.83 | 41
29 | 1 | 500 | 5.5 | 29634.08 | 43 |  | 71 | 4 | 50 | 6 | 29265.07 | 37 |  | 114 | 5 | 200 | 5.75 | 29284.08 | 41
30 | 1 | 500 | 5.75 | 29559.19 | 43 |  | 72 | 4 | 50 | 6.25 | 29186.63 | 37 |  | 115 | 5 | 200 | 6 | 29253.31 | 41
31 | 1 | 500 | 6 | 29504.86 | 43 |  | 73 | 4 | 100 | 5.5 | 30005.30 | 38 |  | 116 | 5 | 200 | 6.25 | 29227.89 | 41
32 | 1 | 500 | 6.25 | 29450.37 | 43 |  | 74 | 4 | 100 | 5.75 | 29978.57 | 37 |  | 117 | 5 | 300 | 5.5 | 29637.01 | 42
33 | 3 | 25 | 5.5 | 28108.43 | 35 |  | 75 | 4 | 100 | 6 | 29954.70 | 37 |  | 118 | 5 | 300 | 5.75 | 29608.42 | 41
34 | 3 | 25 | 5.75 | 28018.51 | 35 |  | 76 | 4 | 100 | 6.25 | 29926.45 | 37 |  | 119 | 5 | 300 | 6 | 29577.26 | 41
35 | 3 | 25 | 6 | 27951.82 | 35 |  | 77 | 4 | 150 | 5.5 | 30240.51 | 38 |  | 120 | 5 | 300 | 6.25 | 29544.77 | 42
36 | 3 | 25 | 6.25 | 27886.82 | 35 |  | 78 | 4 | 150 | 5.75 | 30219.40 | 37 |  | 121 | 5 | 400 | 5.5 | 29863.25 | 42
37 | 3 | 50 | 5.5 | 28789.70 | 36 |  | 79 | 4 | 150 | 6 | 30204.67 | 37 |  | 122 | 5 | 400 | 5.75 | 29832.83 | 42
38 | 3 | 50 | 5.75 | 28683.83 | 36 |  | 80 | 4 | 150 | 6.25 | 30190.79 | 37 |  | 123 | 5 | 400 | 6 | 29797.02 | 42
39 | 3 | 50 | 6 | 28608.38 | 36 |  | 81 | 4 | 200 | 5.5 | 30441.22 | 38 |  | 124 | 5 | 400 | 6.25 | 29763.65 | 42
40 | 3 | 50 | 6.25 | 28531.91 | 36 |  | 82 | 4 | 200 | 5.75 | 30418.71 | 38 |  | 125 | 5 | 500 | 5.5 | 30022.97 | 42
41 | 3 | 100 | 5.5 | 29319.41 | 36 |  | 83 | 4 | 200 | 6 | 30403.79 | 38 |  | 126 | 5 | 500 | 5.75 | 29992.16 | 42
42 | 3 | 100 | 5.75 | 29264.81 | 36 |  | 84 | 4 | 200 | 6.25 | 30388.13 | 37 |  | 127 | 5 | 500 | 6 | 29957.95 | 42
|  |  |  |  |  |  | 85 | 4 | 300 | 5.5 | 30763.43 | 38 |  | 128 | 5 | 500 | 6.25 | 29923.87 | 42

Tables 4-6 report the results of four exact algorithms: Gurobi on the extended
formulation (MIQP), Gurobi on the perspective reformulation (MISOCP), and two
versions of our Benders algorithms directly on the perspective reformulation:
a single-tree version using the $\epsilon$-technique (ST$\epsilon$-BEN), a
multi-tree version using the $\epsilon$-technique (MT$\epsilon$-BEN). The
metrics we use for the comparison are: i) the computational time (Time)
expressed in seconds, ii) the relative gap at the end of the optimization
(Gap), iii) the number (or the average number for the multi-tree approach) of
branching nodes explored (Nodes) and iv) the number of Benders cuts generated
on integral solutions (BCuts) and fractional solutions (FrBCuts). As for the
Gap, the internal relative gap of Gurobi was used where available. For the
single-tree method in the absence of feasible solutions and multi-tree method
running out of time, the gap was computed as
$100\frac{{\text{ObjVal}-\text{BestLB}}}{\text{ObjVal}}$, where ObjVal is the
optimal value of the instance (i.e., it is the value reported under column
ObjVal in Table 3), and BestLB is the best lower bound found by the method at
hand for the same instance.

Table 4: Performance of the exact procedures on the testbed (part I). Runs
reaching the time limit of 900 seconds are indicated by “TL”.

| MIQP | MISOCP | ST$\bm{\epsilon}$-BEN | MT$\bm{\epsilon}$-BEN
---|---|---|---|---
ID | Times[s] | Gap[%] | Nodes | Times[s] | Gap[%] | Nodes | Times[s] | Gap[%] | Nodes | BCuts | FrBCuts | Times[s] | Gap[%] | Nodes | BCuts
1 | TL | 59.20 | 1673 | 258.71 | 0 | 102 | 84.25 | 0 | 6769 | 23 | 3 | 22.16 | 0 | 3833.67 | 8
2 | TL | 60.06 | 1759 | 291.56 | 0 | 143 | 21.60 | 0 | 3811 | 10 | 3 | 20.79 | 0 | 3053.00 | 8
3 | TL | 60.48 | 1569 | 315.73 | 0 | 97 | 26.80 | 0 | 4260 | 11 | 3 | 15.00 | 0 | 2635.00 | 5
4 | TL | 61.09 | 1559 | 313.29 | 0 | 305 | 25.78 | 0 | 3718 | 14 | 3 | 18.74 | 0 | 1827.50 | 7
5 | TL | 58.08 | 2718 | 320.33 | 0 | 19 | 36.15 | 0 | 7225 | 25 | 3 | 42.61 | 0 | 5139.00 | 10
6 | TL | 58.79 | 2213 | 370.66 | 0 | 38 | 40.46 | 0 | 5630 | 23 | 3 | 33.36 | 0 | 4014.90 | 9
7 | TL | 59.94 | 1577 | 276.53 | 0 | 178 | 33.15 | 0 | 5989 | 14 | 3 | 44.57 | 0 | 4554.00 | 13
8 | TL | 60.68 | 1579 | 373.49 | 0 | 46 | 75.47 | 0 | 5124 | 20 | 3 | 31.46 | 0 | 5560.22 | 8
9 | TL | 55.28 | 8685 | 352.81 | 0 | 36 | 87.58 | 0 | 9931 | 28 | 3 | 41.13 | 0 | 4969.33 | 11
10 | TL | 58.32 | 2983 | 356.60 | 0 | 24 | 24.65 | 0 | 7948 | 13 | 3 | 29.45 | 0 | 5556.83 | 11
11 | TL | 59.16 | 2945 | 320.94 | 0 | 50 | 18.27 | 0 | 6430 | 13 | 3 | 23.58 | 0 | 9818.25 | 7
12 | TL | 60.05 | 1637 | 511.71 | 0 | 25 | 28.43 | 0 | 5388 | 17 | 3 | 37.83 | 0 | 5100.08 | 11
13 | TL | 56.32 | 5801 | 306.79 | 0 | 36 | 18.62 | 0 | 8331 | 21 | 3 | 47.97 | 0 | 5831.00 | 10
14 | TL | 58.32 | 4020 | 306.38 | 0 | 33 | 10.68 | 0 | 5142 | 11 | 3 | 50.39 | 0 | 7713.62 | 12
15 | TL | 58.55 | 2135 | 413.33 | 0 | 15 | 18.77 | 0 | 6575 | 19 | 3 | 31.77 | 0 | 5400.64 | 10
16 | TL | 60.24 | 1566 | 430.35 | 0 | 164 | 49.29 | 0 | 8958 | 31 | 3 | 26.90 | 0 | 6360.40 | 9
17 | TL | 55.80 | 7444 | 283.15 | 0 | 18 | 76.03 | 0 | 10677 | 23 | 3 | 48.53 | 0 | 5145.93 | 14
18 | TL | 57.75 | 2369 | 270.78 | 0 | 26 | 33.78 | 0 | 8335 | 25 | 3 | 42.39 | 0 | 4292.54 | 12
19 | TL | 58.46 | 1732 | 397.00 | 0 | 19 | 27.32 | 0 | 6780 | 14 | 3 | 37.01 | 0 | 5854.60 | 9
20 | TL | 59.43 | 1726 | 440.92 | 0 | 100 | 20.58 | 0 | 6107 | 14 | 3 | 33.88 | 0 | 3921.85 | 12
21 | TL | 55.84 | 5194 | 259.78 | 0 | 295 | 487.34 | 0 | 15980 | 99 | 3 | 52.23 | 0 | 5616.67 | 14
22 | TL | 56.94 | 2055 | 516.78 | 0 | 66 | 41.27 | 0 | 7614 | 22 | 3 | 104.70 | 0 | 6731.43 | 13
23 | TL | 58.02 | 1942 | 399.50 | 0 | 111 | 33.16 | 0 | 6466 | 27 | 3 | 68.73 | 0 | 5918.55 | 10
24 | TL | 59.50 | 1607 | 593.59 | 0 | 25 | 235.28 | 0 | 11487 | 17 | 3 | 26.46 | 0 | 3332.50 | 9
25 | TL | 56.43 | 4204 | 327.30 | 0 | 401 | 97.24 | 0 | 10956 | 31 | 3 | 53.12 | 0 | 4233.00 | 20
26 | TL | 57.19 | 2703 | 341.63 | 0 | 429 | 83.61 | 0 | 9859 | 26 | 3 | 46.07 | 0 | 4908.38 | 12
27 | TL | 57.45 | 1929 | 296.50 | 0 | 282 | 214.30 | 0 | 13550 | 25 | 3 | 36.24 | 0 | 5149.21 | 13
28 | TL | 58.58 | 1678 | 424.01 | 0 | 180 | 332.51 | 0 | 12828 | 34 | 3 | 45.39 | 0 | 6234.58 | 11
29 | TL | 56.02 | 4406 | 488.78 | 0 | 210 | 107.95 | 0 | 10360 | 26 | 3 | 58.53 | 0 | 4697.83 | 17
30 | TL | 56.12 | 2068 | 332.20 | 0 | 161 | 85.83 | 0 | 7897 | 24 | 3 | 49.61 | 0 | 3509.40 | 14
31 | TL | 58.11 | 1661 | 298.03 | 0 | 315 | 40.16 | 0 | 7321 | 25 | 3 | 34.36 | 0 | 4156.42 | 11
32 | TL | 59.28 | 1595 | 420.42 | 0 | 139 | 96.11 | 0 | 9181 | 20 | 3 | 40.30 | 0 | 4729.11 | 8
33 | TL | 59.30 | 6007 | 372.96 | 0 | 119 | 35.03 | 0 | 21349 | 14 | 3 | 27.82 | 0 | 5639.00 | 8
34 | TL | 60.61 | 2367 | 280.78 | 0 | 25 | 59.20 | 0 | 6700 | 30 | 3 | 28.43 | 0 | 6040.25 | 7
35 | TL | 61.23 | 4250 | 313.61 | 0 | 66 | 25.27 | 0 | 5505 | 15 | 3 | 29.35 | 0 | 4957.67 | 8
36 | TL | 63.20 | 1651 | 351.84 | 0 | 26 | 37.19 | 0 | 21822 | 11 | 3 | 24.58 | 0 | 10023.60 | 4
37 | TL | 58.59 | 3849 | 312.11 | 0 | 30 | 32.61 | 0 | 20932 | 11 | 3 | 46.55 | 0 | 15078.11 | 8
38 | TL | 60.26 | 1929 | 250.97 | 0 | 1 | 21.85 | 0 | 6045 | 14 | 3 | 30.44 | 0 | 11769.14 | 6
39 | TL | 60.93 | 2273 | 297.59 | 0 | 407 | 36.24 | 0 | 6368 | 24 | 3 | 38.55 | 0 | 10025.67 | 8
40 | TL | 61.84 | 2156 | 392.11 | 0 | 233 | 69.08 | 0 | 22147 | 12 | 3 | 35.39 | 0 | 9530.86 | 6
41 | TL | 55.89 | 10376 | 392.46 | 0 | 85 | 41.56 | 0 | 7610 | 28 | 3 | 84.04 | 0 | 15476.85 | 12
42 | TL | 58.26 | 6270 | 372.35 | 0 | 41 | 11.99 | 0 | 5234 | 11 | 3 | TL | 0.78 | 12335.80 | 4

Table 5: Performance of the exact procedures on the testbed (part II). Runs
reaching the time limit of 900 seconds are indicated by “TL”.

| MIQP | MISOCP | ST$\bm{\epsilon}$-BEN | MT$\bm{\epsilon}$-BEN
---|---|---|---|---
ID | Times[s] | Gap[%] | Nodes | Times[s] | Gap[%] | Nodes | Times[s] | Gap[%] | Nodes | BCuts | FrBCuts | Times[s] | Gap[%] | Nodes | BCuts
43 | TL | 60.52 | 1934 | 405.90 | 0 | 46 | 232.86 | 0 | 28824 | 41 | 3 | 35.81 | 0 | 9505.25 | 7
44 | TL | 61.43 | 2072 | 391.95 | 0 | 83 | 169.47 | 0 | 19606 | 23 | 3 | 46.80 | 0 | 15093.90 | 9
45 | TL | 56.00 | 9869 | 212.50 | 0 | 410 | 23.94 | 0 | 6312 | 19 | 3 | TL | 0.99 | 6104.20 | 4
46 | TL | 57.94 | 6300 | 411.84 | 0 | 165 | 16.75 | 0 | 5759 | 16 | 3 | TL | 0.12 | 12337.63 | 7
47 | TL | 59.74 | 4459 | 328.92 | 0 | 50 | 42.43 | 0 | 23302 | 20 | 3 | 43.61 | 0 | 14795.33 | 8
48 | TL | 61.32 | 3656 | 647.89 | 0 | 116 | 44.65 | 0 | 14049 | 20 | 3 | 110.68 | 0 | 9475.71 | 6
49 | TL | 56.93 | 6450 | 389.13 | 0 | 231 | 108.14 | 0 | 10941 | 36 | 3 | 89.88 | 0 | 11364.94 | 15
50 | TL | 59.05 | 3062 | 401.14 | 0 | 108 | 25.36 | 0 | 6453 | 20 | 3 | 143.52 | 0.07 | 9050.8 | 4
51 | TL | 59.54 | 2343 | 299.19 | 0 | 193 | 194.15 | 0 | 27922 | 75 | 3 | 64.25 | 0 | 22641.45 | 10
52 | TL | 61.48 | 1647 | 360.51 | 0 | 96 | 37.20 | 0 | 13995 | 22 | 3 | 63.37 | 0 | 13896.60 | 9
53 | TL | 56.71 | 6702 | 433.44 | 0 | 126 | 74.28 | 0 | 8888 | 31 | 3 | 92.55 | 0 | 9266.17 | 17
54 | TL | 57.90 | 4737 | 243.62 | 0 | 116 | 28.28 | 0 | 6666 | 21 | 3 | TL | 0.07 | 11656.89 | 8
55 | TL | 60.10 | 2657 | 433.56 | 0 | 15 | 70.24 | 0 | 24432 | 20 | 3 | 61.08 | 0 | 17704.44 | 8
56 | TL | 60.15 | 1972 | 394.15 | 0 | 1 | 25.40 | 0 | 13212 | 15 | 3 | TL | 1.87 | 12833.00 | 3
57 | TL | 57.12 | 4630 | 286.64 | 0 | 339 | 275.87 | 0 | 14930 | 71 | 3 | 821.38 | 0 | 14196.05 | 21
58 | TL | 58.41 | 2571 | 196.81 | 0 | 1 | 14.34 | 0 | 5737 | 13 | 3 | 47.41 | 0 | 10975.42 | 11
59 | TL | 58.86 | 1854 | 251.41 | 0 | 1 | 114.13 | 0 | 26320 | 27 | 3 | 49.38 | 0 | 16965.38 | 7
60 | TL | 60.45 | 1492 | 384.96 | 0 | 1 | 151.58 | 0 | 18180 | 23 | 3 | 392.73 | 0.85 | 10922.50 | 3
61 | TL | 56.05 | 7294 | 429.59 | 0 | 259 | 402.30 | 0 | 16159 | 95 | 3 | 172.82 | 0 | 11064.80 | 24
62 | TL | 57.57 | 4658 | 333.84 | 0 | 24 | 52.37 | 0 | 7056 | 39 | 3 | 49.12 | 0 | 12540.18 | 10
63 | TL | 58.72 | 2639 | 368.97 | 0 | 23 | 223.68 | 0 | 28817 | 36 | 3 | 68.52 | 0 | 13447.77 | 12
64 | TL | 59.76 | 1902 | 334.72 | 0 | 45 | 114.62 | 0 | 17733 | 27 | 3 | 59.32 | 0 | 13068.56 | 8
65 | TL | 62.93 | 532 | 419.54 | 0 | 43 | 41.47 | 0 | 7288 | 13 | 3 | 19.62 | 0 | 2953.29 | 6
66 | TL | - | 1 | 312.58 | 0 | 512 | 100.63 | 0 | 8131 | 36 | 3 | 23.13 | 0 | 4627.00 | 6
67 | TL | 62.46 | 1070 | 359.73 | 0 | 551 | 21.35 | 0 | 5798 | 11 | 3 | 27.68 | 0 | 4270.29 | 6
68 | TL | 63.26 | 1509 | 419.20 | 0 | 177 | 83.67 | 0 | 7753 | 23 | 3 | 30.00 | 0 | 3378.57 | 6
69 | TL | 59.76 | 1678 | 376.57 | 0 | 103 | 20.32 | 0 | 5519 | 14 | 3 | 32.42 | 0 | 8037.00 | 7
70 | TL | 60.88 | 2333 | 255.38 | 0 | 18 | 137.24 | 0 | 9423 | 26 | 3 | 41.67 | 0 | 9661.38 | 7
71 | TL | 61.73 | 1477 | 429.90 | 0 | 115 | 179.30 | 0 | 9549 | 36 | 3 | 38.96 | 0 | 7840.71 | 6
72 | TL | 62.61 | 1505 | 444.66 | 0 | 246 | 175.79 | 0 | 9660 | 43 | 3 | 40.49 | 0 | 5614.14 | 6
73 | TL | 57.82 | 6177 | 260.73 | 0 | 283 | 24.22 | 0 | 7934 | 27 | 3 | 42.26 | 0 | 7539.46 | 12
74 | TL | 60.24 | 3381 | 497.46 | 0 | 47 | 59.16 | 0 | 8485 | 23 | 3 | 40.09 | 0 | 6610.00 | 12
75 | TL | 61.31 | 1905 | 397.17 | 0 | 39 | 37.23 | 0 | 7684 | 24 | 3 | 34.85 | 0 | 8736.67 | 8
76 | TL | 62.33 | 1466 | 508.15 | 0 | 22 | 17.71 | 0 | 5004 | 17 | 3 | 43.92 | 0 | 13579.60 | 9
77 | TL | 58.44 | 6753 | 213.85 | 0 | 188 | 19.63 | 0 | 7078 | 17 | 3 | 57.25 | 0 | 7674.40 | 14
78 | TL | 59.87 | 2621 | 356.18 | 0 | 99 | 38.12 | 0 | 7721 | 19 | 3 | 46.82 | 0 | 8735.07 | 13
79 | TL | 60.90 | 1925 | 399.55 | 0 | 28 | 36.74 | 0 | 7446 | 18 | 3 | 40.36 | 0 | 13926.63 | 7
80 | TL | 61.92 | 2189 | 476.50 | 0 | 61 | 33.89 | 0 | 7132 | 19 | 3 | 52.57 | 0 | 10375.00 | 8
81 | TL | 58.75 | 3499 | 192.08 | 0 | 192 | 47.56 | 0 | 9630 | 30 | 3 | 59.45 | 0 | 8296.24 | 16
82 | TL | 60.21 | 1882 | 249.83 | 0 | 64 | 42.53 | 0 | 7947 | 16 | 3 | 39.39 | 0 | 10635.73 | 10
83 | TL | 61.30 | 1502 | 314.72 | 0 | 83 | 32.45 | 0 | 7126 | 25 | 3 | 41.31 | 0 | 7489.56 | 8
84 | TL | 62.03 | 1653 | 762.94 | 0 | 39 | 34.96 | 0 | 7870 | 17 | 3 | 43.79 | 0 | 11813.00 | 8
85 | TL | 57.88 | 3859 | 316.56 | 0 | 405 | 46.93 | 0 | 9319 | 33 | 3 | 51.13 | 0 | 6351.93 | 13

Table 6: Performance of the exact procedures on the testbed (part III). Runs
reaching the time limit of 900 seconds are indicated by “TL”.

| MIQP | MISOCP | ST$\bm{\epsilon}$-BEN | MT$\bm{\epsilon}$-BEN
---|---|---|---|---
ID | Times[s] | Gap[%] | Nodes | Times[s] | Gap[%] | Nodes | Times[s] | Gap[%] | Nodes | BCuts | FrBCuts | Times[s] | Gap[%] | Nodes | BCuts
86 | TL | 59.46 | 3217 | 345.93 | 0 | 127 | 46.74 | 0 | 8314 | 31 | 3 | 65.00 | 0 | 12021.62 | 12
87 | TL | 60.28 | 1855 | 358.78 | 0 | 204 | 20.30 | 0 | 5686 | 12 | 3 | 59.06 | 0 | 9716.25 | 11
88 | TL | 61.66 | 1534 | 647.00 | 0 | 62 | 28.32 | 0 | 6646 | 10 | 3 | 26.44 | 0 | 6924.60 | 4
89 | TL | 58.24 | 2277 | 321.66 | 0 | 58 | 25.74 | 0 | 7454 | 28 | 3 | 81.76 | 0 | 8646.84 | 18
90 | TL | 59.08 | 1951 | 339.09 | 0 | 183 | 20.19 | 0 | 6156 | 23 | 3 | 75.66 | 0 | 12627.75 | 15
91 | TL | 60.60 | 1453 | 387.58 | 0 | 108 | 51.42 | 0 | 7995 | 22 | 3 | 69.28 | 0 | 9662.94 | 15
92 | TL | 61.16 | 1497 | 400.76 | 0 | 146 | 34.59 | 0 | 7419 | 19 | 3 | 67.64 | 0 | 13801.82 | 10
93 | TL | 57.48 | 4943 | 310.48 | 0 | 174 | 28.87 | 0 | 6684 | 21 | 3 | 78.60 | 0 | 9130.50 | 15
94 | TL | 58.90 | 2538 | 309.35 | 0 | 73 | 22.32 | 0 | 5976 | 18 | 3 | 94.14 | 0 | 18118.19 | 15
95 | TL | 59.90 | 1671 | 359.47 | 0 | 54 | 20.39 | 0 | 5761 | 13 | 3 | 105.11 | 0 | 19957.29 | 13
96 | TL | 61.13 | 1512 | 349.15 | 0 | 75 | 27.34 | 0 | 6428 | 17 | 3 | 78.65 | 0 | 19785.73 | 10
97 | TL | 57.51 | 3068 | 283.46 | 0 | 168 | 32.86 | 0 | 2854 | 16 | 0 | 40.89 | 0 | 7548.92 | 11
98 | TL | 57.29 | 5126 | 266.62 | 0 | 155 | 68.27 | 0 | 2431 | 23 | 0 | 44.65 | 0 | 11081.40 | 9
99 | TL | 58.54 | 4201 | 288.46 | 0 | 57 | 59.54 | 0 | 2164 | 21 | 0 | 35.00 | 0 | 8192.78 | 8
100 | TL | 59.63 | 2818 | 414.16 | 0 | 164 | 29.32 | 0 | 1071 | 9 | 0 | 33.88 | 0 | 9048.38 | 7
101 | TL | 56.39 | 3915 | 338.86 | 0 | 133 | 114.89 | 0 | 6754 | 25 | 0 | 37.63 | 0 | 9021.27 | 10
102 | TL | 57.56 | 2040 | 408.32 | 0 | 167 | 843.48 | 0 | 7459 | 32 | 0 | 36.71 | 0 | 10459.30 | 9
103 | TL | 58.61 | 2707 | 291.86 | 0 | 188 | 23.38 | 0 | 1648 | 12 | 0 | 32.29 | 0 | 10258.88 | 7
104 | TL | 59.57 | 2012 | 380.54 | 0 | 105 | 42.42 | 0 | 2727 | 19 | 0 | 31.53 | 0 | 10005.25 | 7
105 | TL | 53.93 | 8614 | 279.69 | 0 | 191 | 332.45 | 0 | 3472 | 17 | 0 | 47.78 | 0 | 8094.64 | 13
106 | TL | 56.64 | 4641 | 496.08 | 0 | 77 | TL | 0.004 | 4847 | 33 | 0 | 40.98 | 0 | 10529.82 | 10
107 | TL | 58.16 | 1993 | 323.11 | 0 | 318 | 390.36 | 0 | 3045 | 21 | 0 | 42.82 | 0 | 10925.40 | 9
108 | TL | 58.83 | 1581 | 377.99 | 0 | 144 | 71.41 | 0 | 1825 | 13 | 0 | 34.13 | 0 | 9807.10 | 9
109 | TL | 53.71 | 8310 | 414.93 | 0 | 201 | 49.63 | 0 | 2336 | 32 | 0 | 41.39 | 0 | 10800.83 | 11
110 | TL | 55.34 | 6150 | 255.36 | 0 | 93 | 649.39 | 0 | 4035 | 21 | 0 | 37.78 | 0 | 10825.36 | 10
111 | TL | 58.01 | 3517 | 430.43 | 0 | 361 | 411.15 | 0 | 3351 | 21 | 0 | 48.69 | 0 | 9289.29 | 13
112 | TL | 58.65 | 2613 | 403.33 | 0 | 511 | 514.33 | 0 | 3654 | 16 | 0 | 39.25 | 0 | 13774.00 | 10
113 | TL | 54.84 | 5305 | 310.92 | 0 | 287 | TL | 0.10 | 3903 | 36 | 0 | 166.83 | 2.85 | 21814.67 | 2
114 | TL | 56.49 | 2951 | 326.95 | 0 | 386 | 54.76 | 0 | 2684 | 19 | 0 | 91.41 | 0 | 9101.88 | 16
115 | TL | 57.79 | 2306 | 479.59 | 0 | 421 | 151.00 | 0 | 3390 | 24 | 0 | 53.45 | 0 | 11614.50 | 13
116 | TL | 59.22 | 1670 | 358.71 | 0 | 430 | 779.97 | 0 | 4531 | 24 | 0 | 55.34 | 0 | 13458.54 | 12
117 | TL | 54.06 | 6569 | 329.25 | 0 | 557 | 107.68 | 0 | 2027 | 27 | 0 | 47.02 | 0 | 8207.93 | 14
118 | TL | 56.12 | 3915 | 348.10 | 0 | 724 | 239.35 | 0 | 2552 | 27 | 0 | 58.60 | 0 | 8993.50 | 15
119 | TL | 57.55 | 2369 | 540.82 | 0 | 230 | 65.39 | 0 | 3090 | 24 | 0 | 53.52 | 0 | 9119.93 | 14
120 | TL | 58.21 | 1793 | 455.66 | 0 | 551 | 30.74 | 0 | 2526 | 22 | 0 | 70.87 | 0 | 12007.00 | 15
121 | TL | 54.52 | 4205 | 324.42 | 0 | 412 | TL | 0.27 | 15924 | 27 | 0 | 74.74 | 0 | 10863.88 | 15
122 | TL | 56.40 | 2616 | 319.79 | 0 | 644 | 124.38 | 0 | 2072 | 28 | 0 | 65.27 | 0 | 6738.06 | 17
123 | TL | 57.20 | 1651 | 488.47 | 0 | 223 | 78.22 | 0 | 2464 | 18 | 0 | 58.73 | 0 | 9609.43 | 13
124 | TL | 58.68 | 1587 | 424.23 | 0 | 207 | 186.95 | 0 | 5358 | 23 | 0 | 58.83 | 0 | 11178.13 | 15
125 | TL | 54.07 | 6470 | 343.52 | 0 | 647 | 34.47 | 0 | 2536 | 37 | 0 | 70.27 | 0 | 8621.12 | 16
126 | TL | 55.21 | 4680 | 401.69 | 0 | 683 | 67.71 | 0 | 3553 | 35 | 0 | 68.55 | 0 | 8138.89 | 17
127 | TL | 56.30 | 2028 | 310.92 | 0 | 534 | 132.46 | 0 | 2597 | 21 | 0 | 57.51 | 0 | 10006.08 | 12
128 | TL | 58.04 | 1682 | 469.55 | 0 | 222 | 222.11 | 0 | 2520 | 21 | 0 | 61.26 | 0 | 11707.57 | 13