1 Introduction

The unsplittable capacitated network design problem plays an important role in many applications such as telecommunication network design, production distribution, and express package delivery; see [7, 14, 19] and the references therein. Given a network, a demand set with its origin-destination pairs of nodes for commodities, and a facility set with different types of facilities (with varying capacities and installation costs), the unsplittable capacitated network design problem is to install integer multiples of facilities on each arc of the network and route the flow of each commodity on a single path such that the total flow cannot exceed the total capacity on each arc, and the sum of facility installation costs and flow routing costs is minimized while meeting the demands of the commodities.

Let \(G = (V, E)\) be a directed graph with the node set V and the arc set E. Denote Q and T be the sets of commodities and facilities, respectively. The demand of commodity \(q \in Q\) from the source node \(\zeta ^q\) to the destination node \(\eta ^q\) is known as \(a_q\). Suppose that \(w^{ij}_q\) is the routing cost for commodity q on arc \((i,j)\in E\). If one module of facility \(t\in T\) is installed on arc (i, j) , let \(b_t\) and \(p^{ij}_t\) represent the additional capacity and the installation cost on arc (i, j) , respectively. The existing capacity on arc (i, j) is known as \(c^{ij}\). We introduce the binary variable \( x^{ij}_q \) to denote whether or not commodity q goes through arc (i, j) . The variable \( y^{ij}_t \) denotes the number of facility t installed on arc (i, j) . With these notations and variables, the mathematical formulation of the unsplittable capacitated network design problem is

$$\begin{aligned} \displaystyle \min _{x,y}\&\displaystyle \sum _{(i,j)\in E}\sum _{q\in Q} w_q^{ij}x_q^{ij}+ \sum _{(i,j)\in E}\sum _{t\in T} p_t^{ij}y_t^{ij}, \end{aligned}$$
(1)
$$\begin{aligned} \hbox {s.t.}&\displaystyle \sum _{(i,j)\in E}x^{ij}_q - \sum _{(j,i)\in E}x^{ji}_q = \left\{ \begin{array}{rl} 1, &{} \text {if}~ i = \zeta ^q; \\ -1, &{} \text {if}~ i = \eta ^q; \\ 0, &{} \text {otherwise}, \end{array} \right.&\quad \forall \,q\in Q,~i\in V, \end{aligned}$$
(2)
$$\begin{aligned}&\displaystyle \sum _{q\in Q}a_qx_q^{ij}\le \sum _{t\in T}b_ty^{ij}_t + c^{ij},&\quad \forall \, (i,j)\in E, \end{aligned}$$
(3)
$$\begin{aligned}&x^{ij}_q\in \{0,1\},~y_t^{ij} \in {\mathbb {Z}}_+,&\quad \forall \,q\in Q,\,t\in T,\,(i,j)\in E. \end{aligned}$$
(4)

In the above formulation, we minimize the sum of facility installation costs and flow routing costs in the objective function (1). Constraint (2) is the flow balance constraint. Constraint (3) is the capacity constraint that requires that the total flow cannot exceed the total capacity on each arc.

Problem (1)–(4) is \({\mathcal {N}}{\mathcal {P}}\)-hard even for \(|Q|=1\) and \(|T|=1\) [15]. Hence there is little hope to develop a theoretically efficient algorithm for solving it. Nevertheless, several polyhedral studies of some special cases of the problem have been done in the literature [2, 4, 7, 8, 13, 14, 19, 28, 30], which suggests us that it is possible to develop a computationally efficient algorithm if the polyhedral structure is well understood. Inspired by this, in this paper, we consider the convex hull of the set related to the capacity constraint on each arc, i.e., the so-called unsplittable flow arc-set polyhedron \( P = \mathrm{conv }\,(X) \) where

$$\begin{aligned} X = \left\{ (x,y)\in \{0,1\}^{|Q|}\times {\mathbb {Z}}_+^{|T|}: \sum _{q\in Q}a_qx_q \le \sum _{t\in T}b_ty_t+ c\right\} . \end{aligned}$$

Here the arc subscripts on variables \(x^{ij}_t\) and \( y^{ij}_t\), and parameter \(c^{ij}_t\) are dropped.

There exist several works studying the unsplittable flow arc-set polyhedron. In particular, Brockmüller et al. [13, 14] developed the c-strong inequality for the unsplittable flow arc-set polyhedron when there are only two facilities and the capacity of the second facility is an integer multiple of that of the first one. For problems with a single facility, i.e., \(|T|=1\), Atamtürk and Rajan [4] proposed the k-split c-strong inequality and the lifted cover inequality. Van Hoesel et al. [30] studied the lower convex envelope inequality. Their computational experiments on the c-strong inequality, the k-split c-strong inequality, the lifted cover inequality, and the lower convex envelope inequality demonstrate the effectiveness of integrating these inequalities in a branch-and-cut framework. Benhamiche et al. [8] generalized the c-strong inequality to solve a variant of the unsplittable capacitated network design problem.

Unfortunately, most of these studies are restricted to the unsplittable flow arc-set polyhedron with a single facility or two facilities with divisible capacities. The valid inequalities developed under these assumptions cannot be applied in the context of an arbitrary number of facilities and arbitrary capacities.

In this study, we do not make assumptions on either the number of facilities or the structure of the capacities. Instead, our approach is to develop an exact separation algorithm to solve the separation problem of the unsplittable flow arc-set polyhedron P with an arbitrary number of facilities and arbitrary capacities. More precisely, given a point \( ({\bar{x}},{\bar{y}}) \in {\mathbb {R}}^{|Q|}\times {\mathbb {R}}^{|T|} \), we want to generate a hyperplane to separate point \( ({\bar{x}},{\bar{y}}) \) from P or prove that point \( ({\bar{x}},{\bar{y}}) \in P\). To do this, we first analyze the coefficients in the nontrivial facet-defining inequality of polyhedron P, which is employed in formulating the separation problem as an optimization problem. We prove that the solution of the optimization problem corresponds to a facet-defining inequality of polyhedron P. To relieve the computational burden, we show that, in some special cases, a closed form of the optimization problem can be derived. For the general case, the exact separation algorithm, which includes the four steps: preprocessing, row generation, numerical errors, and sequential lifting, is employed in solving the optimization problem. Furthermore, a new technique is presented to accelerate the exact separation algorithm, which significantly decreases the number of iterations in the row generation subroutine. Finally, a comprehensive computational study is presented to test the effectiveness of the proposed algorithm.

It is worth noting that the considered exact separation for the unsplittable flow arc-set polyhedron can be seen as an extension of the exact separation for the 0–1 knapsack polytope; see [5, 6, 9,10,11,12, 20, 31] and the references therein. The difference is that the exact separation for the 0–1 knapsack polytope cannot handle non-binary integer variables, whereas the approach in this paper takes the non-binary integer variables into consideration such that it can be customized to solve the unsplittable capacitated network design problem.

The organization of this paper is as follows. In Sect. 2, we analyze the properties of the nontrivial facet-defining inequalities of polyhedron P and formulate the separation problem as an optimization problem. In Sect. 3, We consider some special cases for which a closed form of the optimization problem can be derived. In Sect. 4, we give a framework of the exact separation algorithm including preprocessing in Sect. 4.1, row generation in Sect. 4.2, numerical errors in Sect. 4.3, and sequential lifting in Sect. 4.4. In Sect. 5, we present the numerical results. Finally, in Sect. 6, we give some conclusions and future works.

Throughout this paper, let \(\varvec{e^i}\in {\mathbb {R}}^{|Q|}\) and \( \varvec{f^j}\in {\mathbb {R}}^{|T|}\) be the i-th |Q| -dimensional unit vector and j-th |T| -dimensional unit vector, respectively. Denote \({\varvec{e}}=(1,1,\ldots ,1)^{\top }\in {\mathbb {R}}^{|Q|}\) and \({\varvec{f}} = (1,1,\ldots ,1)^{\top }\in {\mathbb {R}}^{|T|}\). We use \( X_{\mathrm{{LP}}} \) to denote the linear relaxation of set X obtained by relaxing the integer variables to continuous variables. We assume that \( T = \{1, \ldots , |T|\} \ne \varnothing \), \( 0 < b_1 \le \cdots \le b_{|T|} \), and \( a_q > 0\) for all \( q\in Q \). Without loss of generality, we assume \(a^{\top } {\varvec{e}} - c >0\) since otherwise the capacity constraint in the unsplittable flow arc-set X is redundant.

2 Separation problem for the unsplittable flow arc-set polyhedron

In this section, we first study the polyhedral properties of the unsplittable flow arc-set polyhedron P. Then we formulate the separation problem over polyhedron P as an optimization problem and prove that there exists an optimal solution which corresponds to a facet-defining inequality of polyhedron P.

2.1 Characteristics of the unsplittable flow arc-set polyhedron

We first note that polyhedron P is full dimensional.

Proposition 1

The dimension of polyhedron P is \(|Q|+|T|\).

Next, the following characterizations of some vertices and extreme rays of polyhedron P are straightforward.

Proposition 2

The extreme rays of polyhedron P are \(({\varvec{0}}, \varvec{f^1}),\ldots ,({\varvec{0}},\varvec{f^{|T|}})\).

Proposition 3

The point \((x,\rho _t(x)\varvec{f^t})\) is a vertex of polyhedron P for each \(x\in \{0,1\}^{|Q|}\) and \( t \in T \), where

$$\begin{aligned} \rho _t(x)= \max \left\{ \left\lceil (a^{\top }x-c)/b_t\right\rceil ,0\right\} . \end{aligned}$$
(5)

The initial constraints \(x_q\ge 0\), \(x_q\le 1\), \(y_t\ge 0\), and \(a^{\top }x\le b^{\top }y +c\) are called trivial inequalities of polyhedron P. We now present a necessary condition to guarantee the nontrivial inequality

$$\begin{aligned} \alpha ^{\top }x \le \beta ^{\top }y + \gamma \end{aligned}$$
(6)

to be facet-defining for polyhedron P.

Proposition 4

Let (6) be a nontrivial facet-defining inequality of polyhedron P. Then

  1. (i)

    \(\alpha _q\ge 0\) and \(\alpha _q\le \lceil a_q/b_t\rceil \beta _t\) for each \( q\in Q\) and \( t\in T\);

  2. (ii)

    \(0< \beta _t\le \lceil b_t/b_k\rceil \beta _k\) for each \(t,k\in T\) with \( t\ne k \);

  3. (iii)

    \(\gamma \ge {\left\{ \begin{array}{ll} 0, &{} \text {if}~c\ge 0\mathrm{{;}}\\ -\lceil -c/b_t\rceil \beta _t,~\forall ~t\in T, &{} \text {otherwise}. \end{array}\right. }\)

Proof

For notation convenience, denote \(F = P \cap \{(x,y) : \alpha ^{\top }x =\beta ^{\top }y+ \gamma \}\).

  1. (i)

    For each \(q\in Q\), since F is a nontrivial facet of polyhedron P and inequality (6) differs from \(x_q\ge 0\), there exists a point \((x^{(1)}, y^{(1)})\in F\) with \(x^{(1)}_q=1\) and

    $$\begin{aligned} \alpha ^{\top }x^{(1)} = \beta ^{\top }y^{(1)} + \gamma . \end{aligned}$$
    (7)

    Since the coefficient \( a_q > 0 \), we have \((x^{(1)}-\varvec{e^{q}},y^{(1)})\in P\), and hence point \((x^{(1)}-\varvec{e^{q}},y^{(1)})\) satisfies the valid inequality (6), i.e.,

    $$\begin{aligned} \alpha ^{\top }x^{(1)} - \alpha _q \le \beta ^{\top }y^{(1)} + \gamma . \end{aligned}$$
    (8)

    Subtracting (7) from (8), we obtain \(\alpha _q\ge 0\). On the other hand, as inequality (6) differs from \(x_q\le 1\), there exists a point \((x^{(2)}, y^{(2)})\in F\) such that \(x^{(2)}_q=0\). This, combined with the fact that \((x^{(2)}+ \varvec{e^q},y^{(2)}+ \lceil {a_q}/{b_t}\rceil \varvec{f^t})\in P\) for each \( t \in T \), indicates that \(\alpha _q \le \displaystyle \lceil a_q/b_t\rceil \beta _{t}\).

  2. (ii)

    For each \(t\in T\), since inequality (6) differs from \(y_t\ge 0\), there exists a point \((x^{(3)},y^{(3)})\in F\) such that \(y^{(3)}_t\ge 1\). Then for all \(k\in T\backslash \{ t\}\), we have \((x^{(3)}, y^{(3)}-\varvec{f^t}+ \lceil {b_t}/{b_{k}}\rceil \varvec{f^{k}})\in P\), which further implies that \(\beta _t\le \lceil b_t/b_k\rceil \beta _k\). For each \( t \in T \), as \(({\varvec{0}},\varvec{f^t})\) is an extreme ray of polyhedron P (see Proposition 2), we have \(\beta _t\ge 0\). If \(\beta _\tau = 0\) for some \( \tau \in T \), then \(0\le \beta _t \le \lceil b_t/b_\tau \rceil \beta _\tau =0\) for all \(t\in T\backslash \{\tau \} \). Hence, inequality (6) reduces to \( \alpha ^{\top }x \le \gamma \). This, together with the fact that \( ({\varvec{e}}, \lceil (a^{\top }{\varvec{e}}-c)/b_t\rceil \varvec{f^t})\in P\), implies \(\alpha ^{\top }{\varvec{e}}\le \gamma \). However, this means that inequality (6) is dominated by the bound constraints, and thus it cannot define a facet of polyhedron P. Therefore, \( \beta _t > 0 \) for all \( t \in T \).

  3. (iii)

    If \(c\ge 0\), since \(({\varvec{0}},{\varvec{0}})\in P\), then \(\gamma \ge 0\); otherwise, as \(({\varvec{0}},\lceil -c/b_t\rceil \varvec{f^t} )\in P\), it follows that \(\gamma \ge -\lceil -c/b_t\rceil \beta _t\) for all \(t \in T\).\(\square \)

2.2 Separation problem

Given a point \(({\bar{x}},{\bar{y}})\in {\mathbb {R}}^{|Q|} \times {\mathbb {R}}^{|T|}\), the separation problem of polyhedron P is to construct a hyperplane (induced by inequality (6)) separating \(({\bar{x}},{\bar{y}})\) from P strictly, i.e.,

$$\begin{aligned} \alpha ^{\top }x\le \beta ^{\top }y+\gamma ,~\forall ~(x,y) \in P, \end{aligned}$$

and

$$\begin{aligned} \alpha ^{\top }{\bar{x}}> \beta ^{\top }{\bar{y}}+ \gamma , \end{aligned}$$

or prove that no such hyperplane exists, i.e., point \(({\bar{x}},{\bar{y}}) \in P \). The separation problem is trivial to be solved if point \(({\bar{x}},{\bar{y}})\notin X_{\mathrm{{LP}}}\) since one of the inequalities \( 0 \le x_q \le 1 \), \( q \in Q \), \( y_t \ge 0\), \( t \in T \), and \( \sum _{q\in Q}a_qx_q \le \sum _{t\in T}b_ty_t+ c \) must be violated by this point. Hence, we assume that point \(({\bar{x}},{\bar{y}}) \in X_{\mathrm{{LP}}} \) throughout this paper. Solving the separation problem is equivalent to solving

$$\begin{aligned} \begin{aligned} v = \displaystyle \max _{\alpha , \beta ,\gamma }~&{\bar{x}}^{\top }\alpha -{\bar{y}}^{\top }\beta -\gamma ,\\ \hbox {s.t.}~&x^{\top }\alpha - y^{\top }\beta -\gamma \le 0, ~\forall ~(x,y) \in P,\\ ~&(\alpha ,\beta ,\gamma )\in S, \end{aligned} \end{aligned}$$
(9)

where S is a closed convex set which guarantees that problem (9) is bounded. If \(v\le 0\), we prove \(({\bar{x}},{\bar{y}}) \in P\); otherwise, we find the hyperplane \( \alpha ^{\top }x = \beta ^{\top }y+\gamma \) separating \(({\bar{x}},{\bar{y}})\) from P strictly. Let \((x^1,y^1),\ldots ,(x^u, y^u)\), \( u \in {\mathbb {Z}}_+ \), be the vertices of polyhedron P. From the well-known Minkowski-Weyl theorem [25, 32] and the description of the extreme rays of polyhedron P in Proposition 2, problem (9) can be reduced to

$$\begin{aligned} \begin{aligned} v = \displaystyle \max _{\alpha , \beta ,\gamma }~&{\bar{x}}^{\top }\alpha -{\bar{y}}^{\top }\beta -\gamma ,\\ \hbox {s.t.}~&(x^k)^{\top }\alpha - (y^k)^{\top }\beta -\gamma \le 0,~k = 1,\ldots ,u,\\ ~&\beta _t\ge 0,~\forall ~t\in T,~ (\alpha ,\beta ,\gamma )\in S. \end{aligned} \end{aligned}$$
(10)

Notice that from Proposition 4, we have \( \beta _t > 0 \) for all \( t \in T \) in the nontrivial facet-defining inequality (6) of polyhedron P. Hence, here we consider a special choice of S:

$$\begin{aligned} S= \{(\alpha ,\beta ,\gamma ): \beta _1 = 1\}, \end{aligned}$$
(11)

where \( 1 \in T \). Based on this selection, problem (10) further reduces to

$$\begin{aligned} \begin{aligned} \displaystyle \max _{\alpha ,\beta ,\gamma }~&{\bar{x}}^{\top }\alpha -{\bar{y}}^{\top }\beta - \gamma ,\\ \hbox {s.t.}~&(x^k)^{\top }\alpha - (y^k)^{\top }\beta - \gamma \le 0,~k = 1,\ldots ,u,\\&\beta _1 = 1,~\beta _t\ge 0,~\forall ~ t\in T\backslash \{1\}. \end{aligned} \end{aligned}$$
(12)

Below we shall show that problem (12) is feasible and bounded.

Lemma 1

Given a point \( ({\bar{x}}, {\bar{y}}) \in X_{\mathrm{{LP}}}\), problem (12) is feasible and bounded.

Proof

We only need to show that problem (12) is bounded since the point \((\alpha ,\beta ,\gamma )=(\varvec{0}, \varvec{f^1}, 0)\) is a feasible solution. We shall use the contradiction argument to prove the boundedness of problem (12). Suppose that problem (12) is unbounded. Then there exists a vector \((\Delta \alpha , \Delta \beta , \Delta \gamma ) \in {\mathbb {R}}^{|Q|} \times {\mathbb {R}}^{|T|} \times {\mathbb {R}}\) such that \( \Delta \beta _1=0 \), \(\Delta \beta _t \ge 0\) for all \( t \in T\backslash \{1\} \),

$$\begin{aligned} (x^k)^{\top }\Delta \alpha - (y^k)^{\top }\Delta \beta - \Delta \gamma \le 0~\text {for all}~k = 1,\ldots ,u , \end{aligned}$$
(13)

and

$$\begin{aligned} {\bar{x}}^{\top }\Delta \alpha - {\bar{y}}^{\top }\Delta \beta - \Delta \gamma > 0. \end{aligned}$$
(14)

Combining (13) and (14), we have

$$\begin{aligned} ({\bar{x}}-x^k)^{\top }\Delta \alpha > ({\bar{y}}-y^k)^{\top }\Delta \beta ~\text {for all}~k=1,\ldots ,u. \end{aligned}$$
(15)

Define a point \( {\hat{x}} \in {\mathbb {R}}^{|Q|} \) such that

$$\begin{aligned} {\hat{x}}_q = {\left\{ \begin{array}{ll} 1,&{} \text {if}~\Delta \alpha _q \ge 0;\\ 0,&{} \text {if}~\Delta \alpha _q < 0.\\ \end{array}\right. } \end{aligned}$$

From Proposition 3, point \(({\hat{x}},\rho _1({\hat{x}})\varvec{f^1})\) is a vertex of polyhedron P where \(\rho _1({\hat{x}})\) is defined in (5). By substituting this point into (15) and using \(0 \le {\bar{x}}_q \le 1\) for all \( q \in Q\), \( {\bar{y}}_t \ge 0\) for all \( t \in T \) (as \( ({\bar{x}}, {\bar{y}}) \in X_{\mathrm{{LP}}} \)), \(\Delta \beta _1 = 0\), and \(\Delta \beta _t \ge 0\) for all \(t \in T \backslash \{1\}\), we have that

$$\begin{aligned} \begin{aligned} \displaystyle 0 \ge \sum _{q\in Q,\,\Delta \alpha _q\ge 0}({\bar{x}}_q-1)\Delta \alpha _q + \sum _{q\in Q,\,\Delta \alpha _q< 0}({\bar{x}}_q-0)\Delta \alpha _q >&\\ \sum _{t\in T\backslash \{1\}}({\bar{y}}_t-0)\Delta \beta _t + ({\bar{y}}_1 - \rho _1({\hat{x}})) \Delta \beta _1 \ge 0, \end{aligned} \end{aligned}$$

which leads to a contradiction. Thus problem (12) is bounded. \(\square \)

Lemma 1 guarantees that problem (12) contains at least one optimal solution. Next, we shall show that using the simplex method to solve problem (12), we will obtain a facet-defining inequality of polyhedron P.

Theorem 1

Given a point \( (\bar{x}, \bar{y}) \in X_{\mathrm{{LP}}}\), the basic optimal solution \( (\alpha , \beta , \gamma ) \) of the linear programming problem (12) defines a facet-defining inequality (6) of polyhedron P.

Proof

Let \( (\alpha , \beta , \gamma ) \) be a basic optimal solution of problem (12) and (6) be the corresponding inequality. For notation purpose, denote \(h = |Q|+ |T|\). From the linear programming theory, except the equality \( \beta _1 = 1 \), there exist another h active constraints at point \( (\alpha , \beta , \gamma ) \) such that the vectors of the coefficients of the constraints are linearly independent. Without loss of generality, we assume that these constraints are \( \beta _t = 0 \) for \(t = \ell _1, \ldots , \ell _\tau \) (\(\ell _i \in T\backslash \{1\}\), \(i=1, \ldots , \tau \)), and \( (x^k)^{\top }\alpha - (y^k)^{\top }\beta - \gamma = 0 \) for \( k = 1, \ldots , h - \tau \). Then

$$\begin{aligned} \begin{bmatrix} (x^{1})^{\top }&{} (y^{1})^{\top } &{} 1 \\ \vdots &{} \vdots &{}\vdots \\ (x^{h-\tau })^{\top } &{}(y^{h-\tau })^{\top } &{}1\\ {\varvec{0}}^\top &{} (\varvec{f^{\ell _1}})^{\top } &{} 0 \\ \vdots &{} \vdots &{}\vdots \\ {\varvec{0}}^\top &{} (\varvec{f^{\ell _\tau }})^{\top } &{} 0 \\ \end{bmatrix}\begin{bmatrix} \alpha \\ -\beta \\ -\gamma \end{bmatrix}= \begin{bmatrix} 0\\ \\ \vdots \\ \\ \\ 0\\ \end{bmatrix}, \end{aligned}$$

where the rank of the coefficient matrix is h. By adding the first row to the last \(\tau \) rows, we have

$$\begin{aligned} \begin{bmatrix} (x^{1})^{\top }&{} (y^{1})^{\top } &{} 1 \\ \vdots &{} \vdots &{}\vdots \\ (x^{h-\tau })^{\top } &{}(y^{h-\tau })^{\top } &{}1\\ (x^{1})^{\top }&{} (y^{1}+\varvec{f^{\ell _1}})^{\top } &{} 1 \\ \vdots &{} \vdots &{}\vdots \\ (x^{1})^{\top }&{} (y^{1}+\varvec{f^{\ell _\tau }})^{\top } &{} 1 \\ \end{bmatrix}\begin{bmatrix} \alpha \\ -\beta \\ -\gamma \end{bmatrix}= \begin{bmatrix} 0\\ \\ \vdots \\ \\ \\ 0\\ \end{bmatrix}. \end{aligned}$$

It follows immediately that the rank of the new coefficient matrix is also h. Hence the points \((x^{1}, y^{1}), (x^{2}, y^{2}), \ldots , (x^{h-\tau }, y^{h-\tau })\), \((x^{1}, y^{1}+\varvec{f^{\ell _1}}), \ldots , (x^{1}, y^{1}+\varvec{f^{\ell _\tau }})\) are affinely independent. Furthermore, these points satisfy (6) at equality. This, together with the fact that the dimension of polyhedron P is h in Proposition 1, implies that (6) is a facet-defining inequality of polyhedron P. \(\square \)

3 A closed form of the separation problem: sufficient conditions

Given a point \( ({\bar{x}}, {\bar{y}}) \in X_{\mathrm{{LP}}} \), in this section, we consider some special cases for which a closed form of problem (12) can be derived. The analysis result of these special cases will be used as a preprocessing technique to reduce the computational time of solving problem (12) for the general case. We first consider the case with a single facility, i.e., \(|T|=1\), and then generalize it to the multifacility case, i.e., \(|T|\ge 2\). All proofs of the propositions in this section are given in the “Appendix”.

To simplify the notation, in what follows, we refer problem (12) as the separation problem of polyhedron P. We say inequality (6) solves the separation problem (12) for polyhedron P if point \( (\alpha , \beta , \gamma ) \) is one of its optimal solutions. For notation purpose, denote

$$\begin{aligned} {r} =\max \left\{ \left\lceil - c/b_1\right\rceil , 0\right\} \end{aligned}$$
(16)

and

$$\begin{aligned} d \in \mathrm{argmax}\,_q\left\{ {\bar{x}}_q: q\in Q\right\} . \end{aligned}$$
(17)

3.1 Single facility

In this subsection, we consider the single facility case, i.e., \( T = \{1\} \). If \( Q = \varnothing \), problem (12) is a single variable problem which can be trivially solved, and hence we assume \( Q \ne \varnothing \). Throughout this subsection, we restrict to consider the case where (i) the capacity of one module of this facility is larger than or equal to the demand of each commodity; and (ii) if any commodity q goes through this arc, we need to install one more module of this facility on the arc. Mathematically, this can be written as the following two assumptions:

  1. (i)

    \( a_q \le b_1 \) for all \( q \in Q \);

  2. (ii)

    \( b_1r + c < a_q \) for all \( q \in Q \).

We next give a closed form of the separation problem (12) under some conditions in Propositions 5 and 6 .

Proposition 5

Let \(({\bar{x}},{\bar{y}}) \in X_{\mathrm{{LP}}}\). Suppose that \(T =\{1\}\), (i), (ii), and

$$\begin{aligned} \sum _{q \in Q}a_q\le b_1(r +1) + c \end{aligned}$$
(18)

hold. The inequality \(x_d\le y_1-r\) solves the separation problem (12). Moreover, it defines a facet of polyhedron P.

Example 1

Let \(X_1 = \left\{ (x, y)\in \{0,1\}^4\times {\mathbb {Z}}_+: 11x_1+15x_2+24x_3+50x_4\le 100y\right\} \) and \(({\bar{x}}, {\bar{y}}) = (0.3,0.5,0.9,0.1,0.38)\). Since \( r = 0\), \( c= 0 \le a_q \) for \( q =1,2,3,4 \), and \( \sum _{q=1}^{4} a_q = 11+15+24+50 = 100 \le 100= b_1\), by Proposition 5, inequality \(x_3\le y\) is a solution of the separation problem (12) which is violated by point \( ({\bar{x}}, {\bar{y}}) \).

The condition (18) in Proposition 5 requires a large b. To see this, suppose \( c = 0\). Then \( r = 0 \) and the condition reduces to \( \sum _{q \in Q} a_q \le b_1 \), which means one module of the facility is enough to cover all the demands through this arc. Due to this, we derive a closed form under a condition with a smaller b in the following.

Proposition 6

Let \(({\bar{x}},{\bar{y}}) \in X_{\mathrm{{LP}}}\). Suppose that \(T =\{1\}\), (i), (ii), and

$$\begin{aligned} \sum _{q \in Q}a_q - a_{{\bar{q}}} \le b_1(r+1)+ c< \sum _{q \in Q}a_q, \ \forall ~ {\bar{q}}\in Q \end{aligned}$$
(19)

hold. Then we have the followings.

  1. (a)

    If \( |Q| \le 2\), the inequality \(\sum _{q\in Q}x_q\le y_1-r\) solves the separation problem (12).

  2. (b)

    If \( |Q| \ge 3 \), one of the following three inequalities solves the separation problem (12), respectively:

    $$\begin{aligned} \left\{ \begin{aligned}&x_{d}\le y_1-r,&\text {if}\ \frac{\sum _{q\in Q}{\bar{x}}_q}{|Q|-1} \le {\bar{x}}_d;\\&\frac{1}{|Q|-1}\sum _{q \in Q} x_q\le y_1 -r,&\text {if}\ {\bar{x}}_d < \frac{\sum _{q\in Q}{\bar{x}}_q}{|Q|-1}\le 1;\\&\sum _{q\in Q}x_q\le y_1-r+|Q|-2,&\text {if}\ \frac{\sum _{q\in Q}{\bar{x}}_q}{|Q|-1}> 1. \end{aligned} \right. \end{aligned}$$

Moreover, in both two cases, the inequalities define facets of polyhedron P, respectively.

Example 2

Let \(X_2 = \left\{ (x, y)\in \{0,1\}^4\times {\mathbb {Z}}_+: 11x_1+15x_2+24x_3+50x_4\le 90y\right\} \). It is easy to verify the conditions of Proposition 6 are satisfied. Suppose \(({\bar{x}}, {\bar{y}}) = (0.4,0.5,0.4,0.4,0.47)\). By simple calculation, it follows that \( d=2 \), \( \bar{x}_d = 0.5\), and \(0.5< \sum _{i=1}^{4}{\bar{x}}_i/(4-1) = {\frac{17}{30}} < 1\). From Proposition 6, this gives us the inequality

$$\begin{aligned} \frac{1}{3}(x_1+x_2+x_3+x_4)\le y, \end{aligned}$$

which cuts off point \( ({\bar{x}}, {\bar{y}}) \).

3.2 Multifacility

In this subsection, we consider the multifacility case, i.e., \( |T| \ge 2 \). Throughout this subsection, we restrict to consider the case with (i), (ii), and

  1. (iii)

    \( \sum _{q\in Q} a_q \le b_t + c\) for all \( t\in T\backslash \{1\}\).

Assumption (iii) means that except the facility \(1 \in T\), the capacity of one module of other facilities is large enough to carry out all commodities.

Given a point \( ({\bar{x}}, {\bar{y}}) \in X_{\mathrm{{LP}}} \), under assumption (iii), we observe that if \(\sum _{t\in T\backslash \{1\}}{\bar{y}}_t\ge 1\), it follows that point \( ({\bar{x}},{\bar{y}}) \in P \). Indeed, since point \( (\bar{x}, \bar{y}) \in X_{\text {LP}} \), we have \( {\varvec{0}} \le \bar{x} \le {\varvec{1}} \). Then, using Proposition 4(i) and the fact that \( a_q > 0 \), in order to prove \( (\bar{x}, \bar{y}) \in P \), it is enough to show that point \( ({\varvec{e}}, \bar{y}) \in P\). Similarly, by Proposition 4(ii) and \( b_t > 0 \) for all \( t \in T\), it suffices to show that \( ({\varvec{e}}, \bar{y}) \in P \) with \( \sum _{t \in T \backslash \{1\}} \bar{y}_t = 1\) and \( \bar{y}_1 = 0 \). The latter is true since \( ({\varvec{e}}, \bar{y}) = \sum _{t \in T \backslash \{1\}} \bar{y}_t ({\varvec{e}}, \varvec{f^t}) \), and by assumption (iii), \( ({\varvec{e}}, \varvec{f^t}) \in P \) for each \( t \in T \backslash \{1\} \). Therefore, in the remaining of this subsection, we only consider the case \(\sum _{t\in T\backslash \{1\}}{\bar{y}}_t< 1\).

Let

$$\begin{aligned} \tilde{Q}:=\bigg \{q \in Q:{\bar{x}}_q > \sum _{t\in T\backslash \{1\}}{\bar{y}}_t\bigg \}\ . \end{aligned}$$
(20)

Based on Propositions 5 and 6 , we can derive similar results under the additional assumption (iii) in the case \( |T| \ge 2 \). This is summarized in Propositions 7 and 8 .

Proposition 7

Let \(({\bar{x}},{\bar{y}}) \in X_{\mathrm{{LP}}}\). Suppose that \(|T| \ge 2\), (i), (ii), (iii), (18) and \(\sum _{t\in T\backslash \{1\}}{\bar{y}}_t < 1\) hold. Then we have the followings.

  1. (a)

    If \( \tilde{Q} = \varnothing \), the inequality \(0 \le y_1 + r \sum _{t\in T\backslash \{1\}} y_t -r\) solves the separation problem (12).

  2. (b)

    If \( \tilde{Q} \ne \varnothing \), the inequality \(x_d\le y_1 + (r+1)\sum _{t\in T \backslash \{1\}}y_t-r\) solves the separation problem (12).

Moreover, in both two cases, the inequalities define facets of polyhedron P, respectively.

Proposition 8

Let \(({\bar{x}},{\bar{y}}) \in X_{\mathrm{{LP}}}\). Suppose that \(|T| \ge 2\), (i), (ii), (iii), (19), and \(\sum _{t\in T\backslash \{1\}}{\bar{y}}_t < 1\) hold. Then we have the followings.

  1. (a)

    If \( \tilde{Q} =\varnothing \), the inequality \(0 \le y_1 + r \sum _{t\in T\backslash \{1\}} y_t -r\) solves the separation problem (12).

  2. (b)

    If \( \tilde{Q}\ne \varnothing \) and \(\tilde{Q} \ne Q\), the inequality \(x_d\le y_1 + (r+1)\sum _{t\in T \backslash \{1\}}y_t-r\) solves the separation problem (12).

  3. (c)

    If \(\tilde{Q} = Q\) with \(|Q|\le 2\), the inequality \(\sum _{q\in Q}x_q\le y_1 + (r+|Q|)\sum _{t\in T \backslash \{1\}}y_t-r\) solves the separation problem (12).

  4. (d)

    If \(\tilde{Q} = Q\) with \(|Q|\ge 3\), one of the following three inequalities solve the separation problem (12), respectively:

    $$\begin{aligned} \left\{ \begin{aligned}&x_{d}\le y_1 + (r+1)\sum _{t\in T\backslash \{1\}} y_t-r,&\text {if}\ \frac{\sum _{q\in Q}{\bar{x}}_q-\sum _{t\in T\backslash \{1\}}{\bar{y}}_t}{|Q|-1} \le {\bar{x}}_d;\\&\frac{1}{|Q|-1}\sum _{q \in Q} x_q\le y_1+ \bigg (r+\frac{|Q|}{|Q|-1}\bigg )\sum _{t\in T\backslash \{1\}} y_t-r,&\text {if}\ {\bar{x}}_d < \frac{\sum _{q\in Q}{\bar{x}}_q-\sum _{t\in T\backslash \{1\}}{\bar{y}}_t}{|Q|-1}\le 1;\\&\sum _{q\in Q}x_q\le y_1+(r+2)\sum _{t\in T\backslash \{1\}} y_t-r+|Q|-2,&\text {if}\ \frac{\sum _{q\in Q}{\bar{x}}_q-\sum _{t\in T\backslash \{1\}}{\bar{y}}_t}{|Q|-1}> 1. \end{aligned} \right. \end{aligned}$$

Moreover, in all four cases, the inequalities define facets of polyhedron P, respectively.

Remark 1

Together with the trivial inequalities, the inequalities listed in Propositions 5–8 describe polyhedron P, respectively. Otherwise, suppose that there exists a facet-defining inequality (6) with \( \beta _1 = 1 \) differing from any of the inequalities in the list of Proposition 5 (Proposition 6, 7, or  8). Then there exists a point \( ({\bar{x}}, {\bar{y}}) \in P \) fulfilling inequality (6) at equality and the inequalities listed in the proposition are all inactive at this point. Considering the separation problem (12) with point \( ({\bar{x}}, {\bar{y}}) \), we know that the inequalities listed in the proposition cannot solve this problem, which leads to a contradiction.

4 An exact separation algorithm

Unlike Sect. 3, in this section, we focus on solving problem (12) without any assumptions on the data. To begin with, we note that since the number of constraints in problem (12) may be exponential, from a computational perspective, it is impractical to solve problem (12) when all the constraints are expressed explicitly. For this reason, we follow [31] to solve problem (12) by decomposing into the following four steps.

  1. 1.

    Preprocessing. In order to save computational time, we implement some preprocessing methods before solving problem (12). Firstly, instead of solving problem (12) directly, we solve a lower dimensional problem to find a violated inequality for the lower dimensional polyhedron. Besides, we avoid solving problem (12) if a most violated inequality is known or point \(({\bar{x}},{\bar{y}}) \in P\) (e.g., by using Propositions 5–8).

  2. 2.

    Row generation. We find a violated inequality or report that no violated one exists by solving problem (12) over the lower dimensional polyhedron with a row generation subroutine.

  3. 3.

    Numerical errors. To avoid numerical instabilities, the constructed inequality is scaled to obtain integral coefficients and the right hand side is recomputed to guarantee its validity.

  4. 4.

    Sequential lifting. The variables that are fixed in the preprocessing step are sequentially lifted to obtain a strong valid inequality for polyhedron P.

4.1 Preprocessing

Given a point \( ({\bar{x}}, {\bar{y}}) \in X_{\mathrm{{LP}}} \), let us fix the variables which take values on their bounds and consider the lower dimensional polyhedron \( P({\bar{x}}, {\bar{y}}) = \mathrm{conv }\,(X({\bar{x}},{\bar{y}}))\) where

$$\begin{aligned} X(\bar{x},\bar{y}) = \left\{ (x,y)\in \{0,1\}^{|{\bar{Q}}|}\times {\mathbb {Z}}_+^{|{\bar{T}}|}:\sum _{q\in {\bar{Q}}}a_qx_q\le \sum _{t\in {\bar{T}}}b_ty_t +\bar{c} \right\} , \end{aligned}$$

\({\bar{Q}}=\{q\in Q: 0<{\bar{x}}_q<1 \}\), \({\bar{T}}= \{t\in T: {\bar{y}}_t>0 \}\), and \(\bar{c} = c-\sum _{q\in Q,\,{\bar{x}}_q=1}a_q\). Similar to problem (12), we may solve the separation problem over polyhedron \( P({\bar{x}}, {\bar{y}}) \) by considering

$$\begin{aligned} \begin{aligned} \displaystyle \max _{\alpha , \beta ,\gamma }&\displaystyle \sum _{q\in {\bar{Q}}}{\bar{x}}_q\alpha _q - \sum _{t\in {\bar{T}}}{\bar{y}}_t\beta _t - \gamma ,\\ \hbox {s.t.}&\displaystyle \sum _{q\in {\bar{Q}}} x^k_q\alpha _q - \displaystyle \sum _{t\in {\bar{T}}}y^k_t\beta _t -\gamma \le 0,\ k =1 ,\ldots , \bar{u},\\&\beta _1 = 1,~\beta _t\ge 0,\ \forall \ t\in {\bar{T}}\backslash \{1\}, \end{aligned} \end{aligned}$$
(21)

where \((x^1,y^1),\ldots , (x^{\bar{u}},y^{\bar{u}})\), \( \bar{u} \in {\mathbb {Z}}_+ \), are the vertices of polyhedron \(P({\bar{x}}, {\bar{y}})\). Notice that here we, without loss of generality, assume that \( {\bar{y}}_1 >0 \). Both the numbers of variables and constraints in problem (21) are less than or equal to those in problem (12). Therefore, it can be expected that problem (21) is easier to be solved than problem (12), especially when the number of fixed variables is large.

We next present a preprocessing procedure based on the following observations.

  1. (i)

    If \( {\bar{Q}}= \varnothing \) and \( {\bar{T}}= \varnothing \), then point \( ({\bar{x}}, {\bar{y}}) \in P({\bar{x}},{\bar{y}}) \), and hence there does not exist a violated inequality.

  2. (ii)

    If \( {\bar{Q}}= \varnothing \) and \( {\bar{T}}= \{1\} \), there is only a single facet-defining inequality \( y_1 \ge \lceil -\bar{c}/b_1 \rceil \) of polyhedron \( P({\bar{x}},{\bar{y}}) \). If it is violated by point \( ({\bar{x}},{\bar{y}}) \), we directly move to the sequential lifting step in Sect. 4.4. Otherwise, we have \( ({\bar{x}}, {\bar{y}}) \in P ({\bar{x}}, {\bar{y}}) \).

  3. (iii)

    If

    $$\begin{aligned} \sum _{q\in {\bar{Q}}}a_q\lceil {\bar{x}}_q\rceil \le \sum _{t\in {\bar{T}}}b_t\lfloor {\bar{y}}_t \rfloor +\bar{c}, \end{aligned}$$
    (22)

    then we can conclude that \( ({\bar{x}}, {\bar{y}}) \in P ({\bar{x}}, {\bar{y}}) \). To see this, let

    $$\begin{aligned} K = \left\{ ({\bar{x}}, {\bar{y}}): x_q \in \{\lfloor {\bar{x}}_q\rfloor , \lceil {\bar{x}}_q\rceil \},\ \forall q\in {\bar{Q}},~y_t\in \{\lfloor {\bar{y}}_t\rfloor , \lceil {\bar{y}}_t\rceil \},\ \forall t\in {\bar{T}}\right\} . \end{aligned}$$

    It follows from (22) that \( K \subseteq X({\bar{x}}, {\bar{y}}) \). Therefore, we have \( ({\bar{x}}, {\bar{y}}) \in \mathrm{conv }\,(K) \subseteq P({\bar{x}}, {\bar{y}})\).

  4. (iv)

    For polyhedron \( P({\bar{x}}, {\bar{y}}) \), if one of the conditions in Propositions 5–8 is satisfied, we move to the sequential lifting step in Sect. 4.4 or conclude that point \( (\bar{x}, \bar{y}) \in P({\bar{x}}, {\bar{y}}) \) depending on whether or not one of the inequalities listed in the corresponding proposition is violated by point \( ({\bar{x}}, {\bar{y}}) \).

4.2 Row generation

We now describe the row generation subroutine for solving problem (21). Instead of solving the whole problem (21) with the potentially exponential many constraints, the row generation subroutine solves a problem with a subset of constraints, i.e.,

$$\begin{aligned} \begin{aligned} v(U) =&\displaystyle \max _{\alpha , \beta ,\gamma } \displaystyle \sum _{q\in {\bar{Q}}}{\bar{x}}_q\alpha _q - \sum _{t\in {\bar{T}}}{\bar{y}}_t\beta _t -\gamma ,\\&\hbox {s.t.}\displaystyle \sum _{q\in {\bar{Q}}}x_q\alpha _q -\sum _{t\in {\bar{T}}}y_t\beta _t - \gamma \le 0,\ \forall ~(x,y)\in U ,\\&\quad \beta _t\ge 0,\ \forall \ t\in {\bar{T}}\backslash \{1\}, \ \beta _1 = 1, \end{aligned} \end{aligned}$$
(23)

in each iteration, where \( U \subseteq P({\bar{x}}, {\bar{y}})\). We call (23) the partial separation problem. Apparently, the partial separation problem (23) is a relaxation of problem (21). Therefore, if \(v(U) \le 0\), we conclude that point \(({\bar{x}},{\bar{y}})\in P({\bar{x}}, {\bar{y}})\); otherwise, the solution \( ({\bar{\alpha }}, {\bar{\beta }}, {\bar{\gamma }}) \) of problem (23) corresponds to the inequality

$$\begin{aligned} \sum _{q\in {\bar{Q}}}{\bar{\alpha }}_qx_q \le \sum _{t\in {\bar{T}}}{\bar{\beta }}_ty_t+ {\bar{\gamma }}, \end{aligned}$$
(24)

which is violated by point \( ({\bar{x}}, {\bar{y}}) \). To further test whether or not inequality (24) is valid for \( X({\bar{x}}, {\bar{y}}) \), we solve the following unbounded integer knapsack problem:

$$\begin{aligned} z= \displaystyle \max _{x,y}\left\{ \displaystyle \sum _{q\in {\bar{Q}}}{\bar{\alpha }}_qx_q -\sum _{t\in {\bar{T}}}{\bar{\beta }}_ty_t- {\bar{\gamma }}: (x,y) \in X({\bar{x}}, {\bar{y}})\right\} . \end{aligned}$$
(25)

If \(z\le 0\), inequality (24) is valid for \( X({\bar{x}}, {\bar{y}}) \); otherwise the optimal solution of (25) violates inequality (24), and hence, we add this solution into set U and the procedure continues.

In contrast to the separation problem (12), for which it is shown to be bounded in Lemma 1, the partial separation problem (23) can be unbounded.

Example 3

Let \(X^0= \left\{ ( x, y)\in \{0,1\}\times {\mathbb {Z}}_+: 3x\le 5y\right\} \) and \(( {\bar{x}}, {\bar{y}}) = (0.5,0.3)\). Initializing \( U = \varnothing \), the partial separation problem (23) reduces to \(v(\varnothing ) = \max _{\alpha ,\gamma } \left\{ 0.5\alpha -0.3-\gamma \right\} \), which is unbounded.

In order to avoid the partial separation problem (23) to be unbounded, we shall provide some bounds on the variables. As it has been shown in Theorem 1, the basic optimal solution of problem (21) corresponds to a facet-defining inequality of polyhedron \( P({\bar{x}}, {\bar{y}}) \). Combining it with Proposition 4 and the fact that \( \beta _1=1 \), we can add the bound constraints

$$\begin{aligned} \begin{aligned}&0 \le \alpha _q\le \lceil a_q/ b_1\rceil ,&\forall ~q\in {\bar{Q}},\\&1 \le \beta _t \le \displaystyle \lceil b_t/b_{1}\rceil ,&\forall ~t\in {\bar{T}},\\&\min \left\{ -\lceil -\bar{c}/b_1\rceil , 0\right\} \le \gamma ,&\end{aligned} \end{aligned}$$
(26)

to problem (21) and then derive a stronger partial separation problem:

$$\begin{aligned} \begin{aligned} v(U) =&\displaystyle \max _{\alpha , \beta ,\gamma } \displaystyle \sum _{q\in {\bar{Q}}}{\bar{x}}_q\alpha _q - \sum _{t\in {\bar{T}}}{\bar{y}}_t\beta _t -\gamma ,\\&\hbox {s.t.}\displaystyle \sum _{q\in {\bar{Q}}}x_q\alpha _q -\sum _{t\in {\bar{T}}}y_t\beta _t - \gamma \le 0, \forall \ (x,y)\in U ,\\&\quad 0\le \alpha _q \le \lceil {a_q}/{b_1}\rceil ,\forall \ q\in {\bar{Q}},\\&\quad 1 \le \beta _t\le \lceil {b_t}/{b_1}\rceil ,\forall \ t\in {\bar{T}},\\&\quad \gamma \ge \min \{ -\lceil -\bar{c}/b_1\rceil , 0\}. \end{aligned} \end{aligned}$$
(27)

Clearly, problem (27) is bounded.

Another issue needed to be addressed is that the optimal solution of the unbounded knapsack problem (25) may not be unique. The inequality derived by some optimal solution may be stronger than those derived by other optimal solutions.

Example 4

Let \(X_3 = \left\{ (x, y)\in \{0,1\}^4\times {\mathbb {Z}}_+: 11x_1+15x_2+24x_3+50x_4\le 60y\right\} \). Consider the point \(({\bar{x}},{\bar{y}}) = (0.9, 0.5, 0.7, 0.1, 0.7)\). After solving the partial separation problem (27) with \(U = \{(1,1,1,1,2)\}\), we obtain the solution \( (\alpha , \beta , \gamma ) = (1,0,1,0,1, 0) \) corresponding to the inequality \(x_1+x_3\le y\). Considering the associated unbounded integer knapsack problem (25), we have two optimal solutions (1, 0, 1, 0, 1) and (1, 1, 1, 0, 1) which correspond to the inequalities

$$\begin{aligned} \alpha _1+\alpha _3 -1 -\gamma \le 0 \end{aligned}$$
(28)

and

$$\begin{aligned} \alpha _1+ \alpha _2+\alpha _3 -1 -\gamma \le 0, \end{aligned}$$
(29)

respectively. Obviously, inequality (29) is stronger than inequality (28). If we add inequality (29) to problem (27), in the next iteration, we will obtain the optimal solution \( (\alpha , \beta , \gamma ) = (1,0,0,1,1,0) \) corresponding to the inequality \(x_1+x_4\le y\), which is violated by point \( ({\bar{x}}, {\bar{y}}) \). Moreover, by solving unbounded integer knapsack problem (25), we know that \(x_1+x_4\le y\) is valid for \( P({\bar{x}},{\bar{y}}) \). However, if we add inequality (28) to problem (27), it can be checked that it needs more iterations to solve the separation problem (21).

We see from Example 4 that, to speed up the row generation procedure, among (possible) multiple optimal solutions of the unbounded integer knapsack problem, it is crucial to select one which corresponds to a stronger inequality for problem (27). Due to this, we next describe an iterative approach to modify an existing optimal solution for problem (25) such that the new optimal solution corresponds to a (possibly) stronger constraint in problem (27).

Let \((\tilde{x},\tilde{y})\) be an optimal solution of the unbounded integer knapsack problem (25). The corresponding inequality in problem (27) is

$$\begin{aligned} \sum _{q\in {\bar{Q}}}\tilde{x}_q\alpha _q -\sum _{t\in {\bar{T}}}\tilde{y}_t\beta _t - \gamma \le 0. \end{aligned}$$
(30)

Now suppose that \(\tilde{x}_{q'} = 0\) for some \( q' \in {\bar{Q}}\) and point \((\tilde{x}+\varvec{e^{q'}},\tilde{y})\in X({\bar{x}},{\bar{y}})\). Notice that as point \(({\bar{\alpha }}, {\bar{\beta }}, {\bar{\gamma }})\) is feasible solution of problem (27), we have \( {\bar{\alpha }}_{q'} \ge 0 \). Hence \((\tilde{x}+\varvec{e^{q'}},\tilde{y})\) is also an optimal solution of problem (25) corresponding to the inequality

$$\begin{aligned} \sum _{q\in {\bar{Q}}}\tilde{x}_q\alpha _q + \alpha _{q'}-\sum _{t\in {\bar{T}}}\tilde{y}_t\beta _t - \gamma \le 0, \end{aligned}$$
(31)

which is obviously stronger than inequality (30). Furthermore, we can recursively use this argument to strengthen an inequality based on the new optimal solution. We describe this iterative approach in Algorithm 1.

figure a

In Algorithm 1, we first sort the variables such that \(a_1\ge a_2 \ge \cdots \ge a_{|{\bar{Q}}|} \) in Step 1. We then recursively modify a point to a new point which corresponds to a (possibly) stronger constraint for problem (27) in Steps 3–5.

To conclude this subsection, we present our row generation procedure in Algorithm 2.

figure b

In Algorithm 2, we first initialize the set \(U = \varnothing \) in Step 1. We then solve problem (27) in Step 2 using the dual simplex method with the warm start information in the last step, see for example [21]. In Step 3, the optimal value of problem (27) is nonpositive, and hence point \((\bar{x},\bar{y}) \in P({\bar{x}}, {\bar{y}})\). In Step 4, we modify the code BOUKNAP, written by Pisinger [27], to solve the unbounded integer knapsack problem (25) with fractional value costs. In Step 5, we know the current inequality is invalid for \(P({\bar{x}},{\bar{y}})\). Therefore, we first modify the solution (obtained in Step 4) using Algorithm 1 such that it corresponds to a (possibly) stronger inequality. Then we add the corresponding inequality to problem (27) and go to Step 2. Finally, in Step 6, if \(\tilde{z} \le 0\), we obtain the valid inequality (24), which is violated by point \(({\bar{x}},{\bar{y}})\).

4.3 Numerical errors

Due to the numerical errors incurred in solving problem (27), we may get an invalid inequality in Algorithm 2. To avoid this, we scale the inequality to obtain integral coefficients by solving the integer programming problem

$$\begin{aligned} \begin{aligned} \displaystyle \min _{\alpha ,\beta ,\theta } ~&\theta ,\\ \hbox {s.t.}~&{\bar{\beta }} \theta = \beta ,\\ ~&{\bar{\alpha }} \theta = \alpha ,\\ ~&\beta _t\in {\mathbb {Z}},\, ~\forall ~t\in {\bar{T}},\\ ~&\alpha _q\in {\mathbb {Z}}, ~\forall ~q\in {\bar{Q}},\\ ~&\theta \ge 1. \end{aligned} \end{aligned}$$

Vasilyev et al. [31] suggested to use the enumeration of the multiplier \(\theta \) from 2 to \(10^4\) and checked the integrality within the tolerance \(10^{-5}\). To further avoid too much computational efforts, here we use the approach described in [1]. More formally, let \(\frac{\lambda _q}{\mu _q}\) be the rational representation of coefficient \({\bar{\alpha }}_q\). Note that \(\frac{\lambda _q}{\mu _q}\) can be obtained by using the Euclidean algorithm within a small tolerance. To avoid too large numbers, we give the following requirements:

$$\begin{aligned} \bigg |\frac{\lambda _q}{\mu _q} - {\bar{\alpha }}_q \bigg |\le \epsilon := 10^{-9},~|\lambda _q|\le \lambda _{\max } := 10^6,~\text {and}~\mu _q\le \mu _{\max } := 10^3. \end{aligned}$$

The representation of \({\bar{\beta }}_t\) for each \( t \in T \) is also required to be satisfied with the same restrictions. Let \( \mu \) be the least common multiple of all denominators. We also require \(\mu \le \mu _{\max }\). We scale \( {\bar{\alpha }} \) and \( {\bar{\beta }} \) by setting \(\alpha _q = \mu {\bar{\alpha }}_q\) for all \(q\in {\bar{Q}}\) and \(\beta _t = \mu {\bar{\beta }}_t\) for all \( t\in {\bar{T}}\) with the requirement that \(|\mu \alpha _{q}|\le \mu _{\max }\) and \(|\mu \beta _{t}|\le \mu _{\max }\). If all the requirements are satisfied, we accept this inequality; otherwise, we drop it.

After scaling \( \alpha _q \) for all \( q \in Q \) and \( \beta _t \) for all \( t \in T \), we recompute the right hand side \(\gamma \) by solving the unbounded integer knapsack problem

$$\begin{aligned} \gamma = \max _{x,y}\left\{ \sum _{q\in {\bar{Q}}}\alpha _qx_q - \sum _{t\in {\bar{T}}}\beta _ty_t : \displaystyle (x,y) \in X({\bar{x}}, {\bar{y}}) \right\} . \end{aligned}$$
(32)

This leads to the inequality

$$\begin{aligned} \sum _{q \in {\bar{Q}}} \alpha _q x_q \le \sum _{t \in {\bar{T}}} \beta _t y_t + {\gamma }, \end{aligned}$$
(33)

which is valid for \(P({\bar{x}}, {\bar{y}})\).

4.4 Sequential lifting

The inequality (33) is valid for \( P({\bar{x}}, {\bar{y}}) \). However, in general, it may be invalid for P. To resolve this problem, the variables, which are fixed in the preprocessing step, are sequentially lifted according to a given lifting order \( \Pi \), i.e., a permutation of \( (Q \backslash {\bar{Q}}) \cup (T \backslash {\bar{T}})\). We now illustrate the procedure to lift the first variable with index k based on inequality (33). Denote

$$\begin{aligned} W(C) = \max _{(x,y) \in \{0,1\}^{|{\bar{Q}}|} \times {\mathbb {Z}}_+^{|{\bar{T}}|} } \left\{ \sum _{q\in {\bar{Q}}}\alpha _qx_q - \sum _{t\in {\bar{T}}}\beta _ty_t: \displaystyle \sum _{q \in {\bar{Q}}}a_q x_q \le \sum _{t \in {\bar{T}}} b_t y_t + C\right\} . \end{aligned}$$

We have the following three cases.

  1. 1.

    If \( k \in Q \backslash {\bar{Q}}\) and \({\bar{x}}_k = 0\), then the lifted inequality is

    $$\begin{aligned} \sum _{q\in {\bar{Q}}}\alpha _qx_q + \alpha _kx_k \le \sum _{t\in {\bar{T}}}\beta _ty_t+ {\gamma }, \end{aligned}$$
    (34)

    where \(\alpha _k ={\gamma } - W(\bar{c}-a_k).\)

  2. 2.

    If \( k \in Q \backslash {\bar{Q}}\) and \({\bar{x}}_k = 1\), then the lifted inequality is

    $$\begin{aligned} \sum _{q\in {\bar{Q}}}\alpha _qx_q + \alpha _k(x_k-1) \le \sum _{t\in {\bar{T}}}\beta _ty_t +{\gamma }, \end{aligned}$$
    (35)

    where \(\alpha _k = W(\bar{c}+a_k)-{\gamma }.\)

  3. 3.

    If \( k \in T \backslash {\bar{T}}\) with \({\bar{y}}_k = 0\), then the lifted inequality is

    $$\begin{aligned} \sum _{q\in {\bar{Q}}}\alpha _qx_q \le \sum _{t\in {\bar{T}}}\beta _ty_t+\beta _ky_k+{\gamma }, \end{aligned}$$
    (36)

    where

    $$\begin{aligned} \beta _k=\max \left\{ (W(\bar{c}+\ell b_k)- {\gamma })/\ell :\ell = 1,\ldots ,{\bar{\ell }}\right\} . \end{aligned}$$

    Here we set \({\bar{\ell }}=\lceil (\sum _{q\in {\bar{Q}}}a_q-\bar{c})/b_k\rceil \) since \(W(\bar{c}+\ell b_k) = \sum _{q\in {\bar{Q}}}\alpha _q\) for all \(\ell \ge {\bar{\ell }}\).

Similarly, we may continue to lift the other variables, which are fixed in the preprocessing step in Sect. 4.1, to obtain a valid inequality of polyhedron P. During the whole lifting process, we need to solve several integer knapsack problems. This can be done via the dynamic programming algorithm. For more details, we refer to [31].

5 Numerical results

In order to test the effectiveness of the exact separation algorithm for solving the unsplittable capacitated network design problem, we implement it in C++ linked with IBM ILOG CPLEX optimizer 12.7.1 [16] library. Following [31], to avoid changing the problem structure, the presolving features are turned off in our experiments. Moreover, the dual simplex method is used to reoptimize the linear programming problem (27) after adding cutting planes. To eliminate the effect of multithreads, the computations are implemented in a single thread. The time limit is set to 7200 seconds. Except where explicitly stated, the other parameters in CPLEX are set to the default ones. The exact separation procedure stops if the optimal value of the LP relaxation problem of the unsplittable capacitated network design problem improves by less than \(0.01\%\) between two adjacent calls.

5.1 Testsets

We conduct our computational study on three testsets of the unsplittable capacitated network design problem (1)–(4). The first testset NDP1, studied in Atamtürk et al. [4], contains 20 instances with a single facility. We use this testset to compare the performance effect of the cuts constructed in our exact separation procedure with that of existing cuts studied in the literature. The second testset NDP2 includes 26 realistic network instances generated by the Survivable Network Design Library (SNDlib 1.0) [26]. 9 of them are instances with a single facility and 17 of them are instances with multifacility. To possibly reduce the unstable behavior of integer programming solvers (see, e.g., [18, 22]), we solve each instance using 10 different random seeds in 7200 seconds. We treat every pair of instance and seed as an individual model, which results in 200 models for testset NDP1 and 260 models for testset NDP2 .

The third testset NDP3 is randomly generated based on [23]. We use this testset to evaluate the performance effects of different capacity module sizes and module costs on the unsplittable capacitated network design problem using standard integer programming solver or our exact separation procedure. Table 1 lists different capacity module sizes and module costs studied in [23]. For the same facility, we assume that its module costs on all arcs are the same. In total, we study problems with 27 different capacities and costs structures. We generate the underlying graphs with 50 vertices using the procedure described in [23]; see also [24, 29]. We generate 10 graphs with 20 commodities with random source nodes and destination nodes. The demand of each commodity is chosen uniformly in \( \{10,11,\ldots , 190 \} \). For each graph, 27 problems are generated based on each item of capacity module sizes and costs in Table 1. Thus, in total, we have 270 models for testset NDP3.

Table 1 Capacity module sizes and module costs of testset NDP3
Full size table

5.2 Different lifting orders

In the sequential lifting step, different lifting orders lead to different inequalities [31]. Hence the first experiment is conducted to test performance effect of different lifting orders. We consider the following four lifting orders.

  • LIFT1: Variables, fixed to one, are lifted first in decreasing order of their coefficients. Then variables, fixed to zero, are lifted also in decreasing order of their coefficients.

  • LIFT2: Variables are lifted in decreasing order of their coefficients.

  • LIFT3: Variables are lifted in decreasing order of their reduced costs.

  • LIFT4: Variables are lifted in increasing order of their reduced costs.

Table 2 Performance comparison of different lifting orders
Full size table

Table 2 presents the arithmetic means of the gap closeds [33] of all models in the corresponding testsets. The gap closed is defined as

$$\begin{aligned} 100 \cdot \displaystyle \frac{z_{\mathrm{{root}}}-z_{\mathrm{{LP}}}}{z_{\mathrm{{ub}}} - z_{\mathrm{{LP}}}}, \end{aligned}$$

where \(z_{\mathrm{{root}}}\) is the objective value of the LP relaxation at root node after adding cuts, \(z_{\mathrm{{LP}}}\) is the value of the LP relaxation before adding cuts, and \(z_{\mathrm{{ub}}}\) is the value of the best known value of the model. Table 2 shows that using the lifting order LIFT4, the gap closed is slightly better than those of other lifting orders in all the three testsets. Therefore, in further computational studies, the lifting order LIFT4 is used in our exact separation algorithm.

Table 3 Performance comparison of the exact separation procedure with the default setting for testset NDP1
Full size table
Table 4 Performance comparison of the exact separation procedure with the default setting for testset NDP2
Full size table
Table 5 Performance comparison of the exact separation procedure with the default setting for testset NDP3
Full size table

5.3 Performance effect of the exact separation procedure

In this subsection, we evaluate the performance effect of adding cuts generated by the exact separation algorithm into the solver. Tables 3, 4, and 5 compare the computational results obtained by using CPLEX (CPX) and our exact separation procedure (EXACT). We report the number of solved models, the average running time, and the average number of explored nodesFootnote 1. Besides, columns “Faster” and “Slower” report the number of models that get at least 10% faster and slower, respectively. We group the three testsets into several brackets. The bracket “All” contains the models which can be solved by at least one of the settings. The bracket [n, 7200] contains the models which can be solved by the slower setting in at least n seconds. The larger n is, the harder of the model is. For each bracket, we report the number of considered models in column “total”.

As it can be observed in Tables 3, 4, and 5, our exact separation algorithm has a positive effect on all these three testsets, especially on the hard models. In particular, using the exact separation algorithm, we can solve 6, 8, and 71 more models than those using the default setting on testsets NDP1, NDP2, and NDP3, respectively. This clearly shows that our exact separation algorithm can improve the performance of the solver on solving the unsplittable capacitated network design problem.

Specifically, for testset NDP1, using the exact separation algorithm, we have 79 models for which the running times are slower than that of the default setting while only 34 models are solved faster. This is due to the fact that these models are easier than models in testsets NDP2 and NDP3, and, as a result, the benefit of the exact separation algorithm cannot compensate for its additional overhead on these easy models.

For testset NDP2, we notice that these models are harder than those in testset NDP1 or NDP3. Among these 260 models, CPX and EXACT only solve 55 and 63 of them, respectively. Nevertheless, for the solved models, using EXACT, the average running time decreases from 291s to 269s.

Table 6 Results of different capacity module sizes and module costs on testset NDP3
Full size table

For testset NDP3, Table 5 shows a significant improvement of using our exact separation algorithm. The average running time decreases from 1180s to 145s. To further see where the improvement comes from, we report the computational results of different capacity module sizes and costs of facilities independently in Table 6. To be more specific, we list the results in 9 different groups corresponding to Table 1. In each group, the number of facilities is different in each item (i.e. \( i\_j\_k \)). There are 10 different models in each item with different network structures and commodities. We report the number of average nodes, the average running time, and the number of solved models in each item. In the same group, we can observe that with our exact separation algorithm, the models with more facilities are easier to be solved than those with fewer facilities. However, the same behavior cannot be observed in the computational results of CPX. This shows that compared with problems with a single facility, our exact separation algorithm works better in problems with more facilities.

Fig. 1
figure 1

Performance profile for testsets NDP1, NDP2, and NDP3

Full size image

To end of this subsection, we plot the performance profiles of CPX and EXACT in Fig. 1 to further compare their performance. For each model in the testsets NDP1, NDP2, and NDP3, we compute a factor \(\tau \) as the ratio of the running time to solve to optimality of the considered setting to the minimum running time of two settings CPX and EXACT. Each point \((\tau , \sigma )\) of the curves in Fig. 1 represents that in \(\sigma \) percentage of the models, this particular setting is at most \(\tau \) times slower than the faster setting; for more details, see [17]. Consequently, the higher the curve is, the better the setting performs. Fig. 1 clearly shows that EXACT performs much better than CPX. In particular, EXACT is able to solve \(58.1\%\) of the models to optimality while CPX is only able to solve \(46.3\%\) of the models to optimality.

5.4 Comparison with c-strong inequality [13, 14]

According to [3], the c-strong inequality [13, 14] is quite effective in solving the unsplittable capacitated network design problem. Other proposed inequalities [4, 30] can provide additional improvement but the marginal effect on top of the c-strong inequality is limited. Therefore, we only compare the performance effects of the inequality generated by the exact separation procedure (EXACT) with the c-strong inequality (CSTRONG). Since the c-strong inequality can only be applied in the case of one facility or two facilities with divisible capacities, we only conduct experiments on testset NDP1. The results are reported in Table 7.

Table 7 Performance comparison with the c-strong inequalities on testset NDP1
Full size table

Table 7 shows that the performance of EXACT is better than that of CSTRONG, especially on the hard models. In total, EXACT solves 4 more models than CSTRONG. This shows that even compared with the existing polyhedral studies on the unsplittable flow arc-set polyhedron, our exact separation algorithm is more effective in solving the unsplittable capacitated network design problem.

5.5 Performance effect of using Propositions 5–8

We now report the performance effect of Propositions 5–8 in our procedure in Table 8. The average results over all the models in each testset are presented. Compared with EXACT, NOPROS refers to the setting of calling the row generation subroutine to solve the separation problem (21) even if one of the conditions in Propositions 5–8 is satisfied. For each setting, we report the running time spent in the row generation subroutine in column “Rgtime”. For EXACT, we additionally report the numbers of times that fulfilling the conditions in Propositions 5–8 in “P5”–“P8”, respectively. In column “Ncalls”, we report the number of calling the exact separation algorithm. In column “Rate”, we list the successful rate computed by the number of times fulfilling at least one of the conditions in Propositions 5–8 over the number of calling the exact separation algorithm.

Table 8 Performance effect of using Propositions 5–8
Full size table

As it can be seen in Table 8, using Propositions 5–8, the running time of the separation algorithm decreases considerably. This is due to the fact that, among the total number of the calls of the exact separation algorithm, 66.42%, 31.72%, and 46.17% of them can be computed directly using Propositions 5–8 for testsets NDP1, NDP2, and NDP3 without calling the time-consuming row generation subroutine. Since the models in testset NDP1 contain only a single facility, the conditions in Propositions 5 and 6 occur frequently while those in Propositions 7 and 8 never occur. Compared with those in Propositions 5 and 6, the conditions in Propositions 7 and 8 occur less frequently. This is not surprising since the conditions in Propositions 7 and 8 are much stricter than those in Propositions 5 and 6. Finally, we observe that the reduction on the running time is not consistent with that of the successful rate. This is because, using the row generation subroutine, the separation problems that fulfill one of the conditions in Propositions 5–8 are easier to be solved than those that do not fulfill any of them.

Table 9 Performance effect of using Algorithm 1
Full size table

5.6 Performance effect of using Algorithm 1

Finally, we report the performance effect of employing Algorithm 1 in the exact separation procedure. Compared with EXACT, NOALG2 refers to the setting of the exact separation procedure without using Algorithm 1. Notice that this may lead to a weaker constraint (30) in the row generation subroutine. For each setting, we report the running time of the row generation subroutine and the average iteration of a row generation call in columns “Rgtime” and “Iter”, respectively. We present the average results among all the models. As it can be seen in Table 9, with Algorithm 1, the iteration of row generation subroutine per call decreases significantly, which in turn, reduces 43.66%, 36.77%, and 34.44% of the row generation time in testsets NDP1, NDP2, and NDP3, respectively. This confirms that our proposed algorithm indeed works well in practice. Notice that it is reasonable to observe reductions on the average iteration exceeds reduction on the running time since Algorithm 1 also leads to a denser constraint for the linear programming problem (27).

6 Conclusion and future work

In this paper, we have considered the separation problem of the flow arc-set polyhedron in the unsplittable capacitated network design problem. By solving the separation problem, we generated the facet-defining inequality for the considered polyhedron. We showed that in some special cases, a closed form of the separation problem can be derived. For the general case, we used the exact separation algorithm to solve the separation problem. Moreover, a new technique was proposed to reduce the computational time in the row generation subroutine of the exact separation algorithm. The numerical experiments showed the effectiveness of the exact separation algorithm in solving the unsplittable capacitated network design problem and the advantage of the proposed technique in reducing the exact separation time.

There still exist some ideas to be explored in this study. For the flow arc-set polyhedron, we have proposed a new technique to speed up the row generation subroutine; see Algorithm 1. It deserves to test whether or not the same improvement can be observed for the knapsack polyhedron. In this study, we only implemented the exact separation algorithm on the flow arc-set polyhedron, but it can be extended on the cut-set polyhedron, see for example [2]. We are currently conducting this topic to see the performance effect of generalizing this to the cut-set polyhedron.