A global constraint for the capacitated single-item lot-sizing problem

TL;DR

This paper introduces a new global constraint for the capacitated single-item lot-sizing problem, enabling more efficient solution methods within a constraint programming framework by leveraging a novel lower bound and time decomposition.

Contribution

It formulates the lot-sizing problem as a global constraint with bound consistency, providing a new lower bound and filtering rules based on dynamic programming, improving solution efficiency.

Findings

Global constraint achieves bound consistency in pseudo-polynomial time.

Constraint programming outperforms decomposed models on difficult instances.

The proposed approach is competitive with mixed integer linear programming.

Abstract

The goal of this paper is to set a constraint programming framework to solve lot-sizing problems. More specifically, we consider a single-item lot-sizing problem with time-varying lower and upper bounds for production and inventory. The cost structure includes time-varying holding costs, unitary production costs and setup costs. We establish a new lower bound for this problem by using a subtle time decomposition. We formulate this NP-hard problem as a global constraint and show that bound consistency can be achieved in pseudo-polynomial time and when not including the costs, in polynomial time. We develop filtering rules based on existing dynamic programming algorithms, exploiting the above mentioned time decomposition for difficult instances. In a numerical study, we compare several formulations of the problem: mixed integer linear programming, constraint programming and dynamic programming. We show that our global constraint is able to find solutions, unlike the decomposed constraint programming model and that constraint programming can be competitive, in particular when adding combinatorial side constraints.

Tables

Table 1. Table 1: Indexing sub-problems

Table 2. Table 2: Three consistency levels

Consistency level	$D\left(x\right)$	$D\left(y\right)$
Initial domains	$\{1,2,4\}$	$\{4,5,6,7,8\}$
BC	$\{2,4\}$	$\{4,5,6,7,8\}$
RC	$\{2,4\}$	$\{4,6,8\}$
AC	$\{2,4\}$	$\{4,8\}$

Table 3. Table 3: Single-item lot-sizing - CP and DP

	CP_Basic				CP_LS				DP
Class	NODE	CPU	RNB	OPT	NODE	CPU	RNB	OPT	CPU	OPT
C1LS	8.3E+06	200	100%	0	1	1.1	0%	10	0.2	10
C2LS	1.7E+06	200	65%	0	1	1.0	0%	10	0.2	10
C3LS	6.5E+06	200	100%	0	28	2.1	0%	10	0.4	10
C4LS	3.4E+05	200	64%	0	35	2.2	0%	10	0.4	10
C5LS	8.5E+06	200	56%	0	1	1.1	0%	10	0.2	10

Table 4. Table 4: Single-item lot-sizing - MILP

	MILP_AGG					MILP_UFL
Class	NODE	CPU	RNB	LR	OPT	NODE	CPU	RNB	LR	OPT
C1LS	580	0.1	1%	10%	10	460	0.6	1%	3%	10
C2LS	1360	0.2	2%	10%	10	1643	1.4	2%	3%	10
C3LS	3213	1.7	2%	11%	10	13109	54.1	3%	3%	10
C4LS	2222	1.6	1%	10%	10	14366	61.6	2%	2%	10
C5LS	1691	0.3	2%	11%	10	9336	4.5	2%	3%	10

Table 5. Table 5: Scaling the global constraint

	CP_LS				MILP_AGG					DP
Class	NODE	CPU	RNB	OPT	NODE	CPU	RNB	LR	OPT	CPU	OPT
C1Peaks	125	1.1	2%	10	0	0.0	0%	13%	10	23.8	10
C2Peaks	468	5.7	2%	10	0	0.0	0%	6%	10	23.4	10
C3Peaks	2784	21.0	3%	10	1	0.0	0%	7%	10	22.6	10
C4Peaks	408	3.7	3%	10	0	0.0	0%	15%	10	23.4	10
C5Peaks	446	5.5	2%	10	0	0.0	0%	7%	10	24.0	10

Table 6. Table 6: Single-item lot-sizing with disjunctions

	CP_LS				MILP_AGG					DP
Class	NODE	CPU	RNB	OPT	NODE	CPU	RNB	LR	OPT	CPU	OPT
C1Disj	1	0.0	0%	10	6.4E+05	162.1	43%	52%	2	0.0	10
C2Disj	2	0.0	0%	10	1.4E+04	6.0	27%	38%	10	0.0	10
C3Disj	2	0.0	0%	10	1.6E+03	0.5	27%	38%	10	0.0	10
C4Disj	2	0.0	0%	10	2.2E+04	6.3	39%	48%	10	0.0	10
C5Disj	1	0.0	0%	10	9.6E+05	200.0	52%	61%	0	0.0	10
C6Disj	2	0.1	0%	10	3.2E+05	200.0	43%	52%	0	0.0	10
C7Disj	2	0.1	0%	10	5.0E+04	38.3	27%	39%	10	0.0	10
C8Disj	1	0.1	0%	10	5.7E+03	4.5	28%	39%	10	0.0	10
C9Disj	2	0.1	0%	10	3.6E+04	21.9	38%	48%	10	0.0	10
C10Disj	2	0.1	0%	10	3.3E+05	200.0	52%	60%	0	0.0	10

Table 7. Table 7: Single-item lot-sizing with Q/R

	CP_LS				MILP_AGG
Class	NODE	CPU	RNB	OPT	NODE	CPU	RNB	LR	OPT
C1QR	2	0.3	0%	10	17	0.0	1%	14%	10
C2QR	22	0.2	1%	10	56	0.1	1%	13%	10
C3QR	84	0.4	1%	10	27	0.0	1%	12%	10
C4QR	1	0.3	0%	10	20	0.0	1%	12%	10
C5QR	1	0.4	0%	10	13	0.0	1%	16%	10
C6QR	772	7.2	1%	10	943	0.5	2%	16%	10
C7QR	6488	42.8	0%	10	601	0.4	1%	14%	10
C8QR	26716	134.1	1%	5	392	0.3	1%	13%	10
C9QR	1	1.1	0%	10	510	0.3	1%	15%	10
C10QR	21	1.9	0%	10	1175	0.5	3%	18%	10

Table 8. Table 8: Single-item lot-sizing with disjunctives and Q/R

	CP_LS				MILP_AGG
Class	NODE	CPU	RNB	OPT	NODE	CPU	RNB	LR	OPT
C1DijsQR	1	0.3	0%	10	1808	0.8	7%	19%	10
C2DijsQR	6	0.2	0%	10	637	0.4	2%	14%	10
C3DijsQR	67	0.4	1%	10	710	0.4	2%	13%	10
C4DijsQR	3	0.4	0%	10	964	0.3	5%	16%	10
C5DijsQR	30	0.5	1%	10	209	0.2	12%	25%	10
C6DijsQR	669	4.7	0%	10	40150	53.4	5%	19%	10
C7DijsQR	3471	17.8	0%	10	4839	9.2	2%	15%	10
C8DijsQR	22386	94.0	1%	7	4618	7.4	3%	14%	10
C9DijsQR	53	1.5	0%	10	5066	8.8	4%	17%	10
C10DijsQR	7	2.6	1%	10	1663	1.7	11%	26%	10

Equations

$$\text{minimize}C$$

$$I_{t-1}+X_t$$

$$X_t$$

$$Cp$$

$$(MILP\_AGG)$$

$$Cs$$

$$X_t$$

$$I_t$$

$$Y_t$$

$$\text{eqref}objMILP1,\text{eqref}CpMILP1$$

$$X_t$$

$$(MILP\_UFL)$$

$$X_{tr}$$

$$\sum _{t=1}^r{X}_{tr}$$

$$X_{tr}$$

$$\underline{X_t^{\prime }}=0\text{ and }\overline{X_t^{\prime }}=\overline{X_t}-\underline{X_t}$$

$$\underline{I_t^{\prime }}=0\text{ and }\overline{I_t^{\prime }}=\overline{I_t}-\underline{I_t}$$

$$p_t^{\prime }=p_t$$

$$h_t^{\prime }=h_t$$

$$\displaystyle s^{\prime}_{t}=\left\{\begin{array}[]{rl}0&\text{if}\>\underline{X_{t}}>0\\ s_{t}&\text{otherwise}\end{array}\right.$$

$$d_t^{\prime }=d_t+\underline{I_t}-\underline{X_t}-\underline{I_{t-1}}$$

$${Cp}_{min}$$

$${Ch}_{min}$$

$${Cs}_{min}$$

$$C_{min}$$

$$X_t$$

$$I_t$$

$$Cp$$

$$Ch$$

$$Cs$$

$$C$$

$$\forall t\in [[1,T]]\text{ and }\forall I_t\in [[0,\overline{{\beta }_t}]]$$

$$\text{textsl}f(t,I_t)=\underset{I_{t-1}=a\dots b}{\min }\{\text{textsl}f(t-1,I_{t-1})+\text{mathds}{1}_{X_t>0}{s}_t+p_t{X}_t+h_t{I}_t\}$$

$$\forall t\in [[0,T-1]]\text{ and }\forall I_t\in [[0,\overline{{\beta }_t}]]$$

$$\text{textsl}{f}_r(t,I_t)=\underset{I_{t+1}=c\dots d}{\min }\{\text{textsl}{f}_r(t+1,I_{t+1})+\text{mathds}{1}_{X_{t+1}>0}{s}_{t+1}+p_{t+1}{X}_{t+1}+h_{t+1}{I}_{t+1}\}$$

$$d_t=0$$

$$s_t=0$$

$$[[u,v]]\cap [[u^{\prime },v^{\prime }]]=\varnothing$$

$$\sum _{i\in S}{w}_i\le C.$$

$$\text{textsl}wisp(0)=0$$

$$\text{textsl}wisp(i)=\max \{\text{textsl}wisp(i-1),\text{textsl}wisp(pre{c}_i)+w_i\},\forall i=1,\dots ,n$$

$$\text{textsl}{wisp}_r(n)=0$$

$$\text{textsl}{wisp}_r(i)=\max \{\text{textsl}{wisp}_r(i+1),w_i+\text{textsl}{wisp}_r(suc{c}_i)\},\forall i=1,\dots ,n$$

$$2x=y,x\in \{1,2,4\}\text{and}y\in \{4,5,6,7,8\}$$

$$Cp$$

$$Ch$$

Full text

1 Introduction

The field of production planning addresses numerous complex problems covered by operations research and combinatorial optimization. In particular, lot-sizing problems have been broadly studied. The core problem [30] and several variants have been solved by Dynamic Programming (DP) in polynomial time. Other variants (e.g. time varying production capacity and setup costs, multi-echelon) are NP-hard and are most of the time dealt with Mixed Integer Linear Programming (MILP) formulations (see e.g. [22, 5]).

State-of-the-art approaches for complex lot-sizing problems are currently based on polyhedral techniques such as cutting plane algorithms and can handle a large class of problems with side-constraints. Nonetheless theses techniques may eventually fail when facing combinatorial additional constraints. In this paper, we investigate alternative generic approaches based on combinatorial techniques and designed within the Constraint Programming (CP) framework. The rationale is that a lot of algorithmic results have been obtained on the fundamental problems in this field over the last sixty years. We propose to reuse them as filtering mechanisms and building blocks of a generic solver for lot-sizing. This paper is a first step in that direction: we introduce a new global constraint LotSizing embedding the single-item lot-sizing problem. LotSizing appears to be especially generic and suits well in the modeling of a great variety of lot-sizing problems. The problem being NP-hard, we prove several complexity results on achieving different consistency levels for the constraint. We use a time decomposition to propose a new lower bound for the single-item lot-sizing problem. This time decomposition combined with classical results, namely DP algorithms, enables us to derive interesting cost-based filtering algorithms for LotSizing.

The capacitated single-item lot-sizing problem

In this paper, we focus on the following single-item lot-sizing problem – denoted by $(L)$ – which is used as a building block to tackle more complex lot-sizing problems. The objective is to plan the production of a single product over a finite horizon of $T$ periods $\llbracket 1,T\rrbracket$ in order to satisfy a demand $d_{t}$ at each period $t$, and to minimize the total cost. The (per unit) production cost at $t$ is $p_{t}$ and a setup cost $s_{t}$ is paid if at least one unit is produced at $t$. A holding cost $h_{t}$ is paid for each unit stored at the end of period $t$. Furthermore the production (resp. the inventory) is bounded by minimal and maximal capacities $\underline{\alpha_{t}}$ and $\overline{\alpha_{t}}$ (resp. $\underline{\beta_{t}}$ and $\overline{\beta_{t}}$) at each period $t$.

Figure 1 shows the problem as a graph with the variables and parameters on each arc. For each period, the incoming arcs corresponds to the possible production (vertical arcs) and inventory from the previous period (horizontal arcs). The outgoing arcs correspond to the demand (vertical arcs) and inventory at the end of the period (horizontal arcs).

In the literature, one can find several models with upper bounds on either production or inventory. It is however unusual to include lower bounds. We make this assumption to be consistent with the CP framework that states domains for the variables.

Literature review

The CP literature is very limited in the field of lot-sizing problems. To the best of our knowledge, [17] is the only paper to study a lot-sizing-related global constraint. They consider a production planning problem in which a set of items has to be produced before their production deadline on a limited capacity machine, with the objective of minimizing stocking costs. This problem can be solved in polynomial time and is a special case of $(L)$ where production costs are set to zero ($p_{t}=s_{t}=0,h_{t}=1$), the production and inventory lower bounds are set to zero ($\underline{\alpha_{t}}=\underline{\beta_{t}}=0$), the production upper bound is constant and there is no inventory upper bound ($\overline{\alpha_{t}}=\alpha,\overline{\beta_{t}}=+\infty$). It can be seen as a scheduling problem with deadlines and the objective of minimizing the total earliness ($P|\tilde{d}_{j},p_{j}=1|\sum E_{j}$ with Graham notation). In their approach, a decision variable is associated to each item and specifies in which period the item has to be produced. This approach is suitable to deal with scheduling problems but seems less relevant to address lot-sizing problems for which large quantities of the same item can be produced in the same period. In [16], the authors extend their global constraint to varying production capacities and stocking costs. Note that CP solvers have been used in the past to solve lot-sizing problems (see e.g. [25] for a distribution multi-echelon system).

We now focus the literature review on some special cases of problem $(L)$. There is no paper, to our knowledge, that considers lower bounds on both production and inventory – see [21, 26] that consider inventory lower bounds only. [30] shows that the uncapacitated problem($\overline{\alpha_{t}}=\overline{\beta_{t}}=+\infty$) can be solved by DP in $O(T^{2})$. This complexity has later been improved to $O(T\log T)$ [9, 29, 1]. When adding a constant production capacity and a constant setup cost, ($s_{t}=s$, $\overline{\alpha_{t}}=\alpha$), the problem can be solved in $O(T^{4})$ with concave costs [11] and in $O(T^{3})$ with linear costs [27]. When the production capacity varies with time, the problem is NP-hard [7]. Note that when $p_{t}=h_{t}=0$, $(L)$ is equivalent to a knapsack problem. With time-varying inventory capacities, the problem can be solved in $O(T^{2})$ with production and inventory setup costs [21, 3]. Finally [4] is likely the only theoretical paper that studied the single-item problem with general capacities.

The rest of the paper is organized as follows. Section 2 presents algorithms from the literature that will be re-used later. Section 3 presents a new lower bound for this problem based on a time decomposition. Section 4 presents the LotSizing global constraint and states complexity results for achieving bound and range consistency. Section 5 presents cost-based filtering mechanisms for LotSizing. Section 6 compares numerically the performances of LotSizing with two MILP formulations, DP and a basic CP model. Section 7 considers two extensions with side constraints.

2 Preliminaries

This section presents classical MILP formulations and DP approaches, that will be used later in the paper. We also show that problem $(L)$ is equivalent to a problem without lower bounds on production and inventory.

2.1 MILP formulations

We list below a summary of the main notations.

Parameters

$T\in\mathbb{N}$: Number of periods.

• 2.

$p_{t}\in\mathbb{N}$: Unit production cost at $t$.

• 3.

$h_{t}\in\mathbb{N}$: Unit holding cost at $t$ (applied to the ending inventory).

• 4.

$s_{t}\in\mathbb{N}$: Setup cost at $t$ (paid if at least one item is produced at $t$).

• 5.

$d_{t}\in\mathbb{N}$: Demand at $t$.

• 6.

$\underline{\alpha_{t}},\overline{\alpha_{t}}\in\mathbb{N}$: Minimal and maximal production quantities at $t$.

• 7.

$\underline{\beta_{t}},\overline{\beta_{t}}\in\mathbb{N}$: Minimal and maximal inventory at the end of period $t$.

• 8.

$I_{0}\in\mathbb{N}$: Initial inventory.

Variables

$X_{t}\in\mathbb{N}$: Quantity produced at $t$.

• 2.

$Y_{t}\in\{0,1\}$: Setup variable that equals $1$ if at least one item is produced at $t$.

• 3.

$I_{t}\in\mathbb{N}$: Inventory at the end of period $t$.

• 4.

$C\in\mathbb{N}$: Total cost.

• 5.

$\mathit{Cp}\in\mathbb{N}$: Sum of production costs.

• 6.

$\mathit{Cs}\in\mathbb{N}$: Sum of setup costs.

• 7.

$\mathit{Ch}\in\mathbb{N}$: Sum of holding costs.

We assume that parameters and domains are integers. We denote by $X$, $I$ and $Y$ the vectors $\langle X_{1},\ldots,X_{T}\rangle$, $\langle I_{1},\ldots,I_{T}\rangle$ and $\langle Y_{1},\ldots,Y_{T}\rangle$. Without loss of generality we consider $I_{0}=0$. We also consider $I_{T}=0$. Indeed, we can compute the minimum mandatory quantity to store at the end of period $T$ from the production and inventory capacity constraints. If this quantity $q$ is strictly positive, we add a dummy period $T+1$ at the end of the time horizon with $p_{T+1}=h_{T+1}=0$, $\overline{\alpha_{T+1}}=\overline{\beta_{T+1}}=0$ and $d_{T+1}=q$.

Problem $(L)$ can be formulated as an aggregated MILP model (see e.g. [22]):

[TABLE]

where (2) are the flow balance constraints for each period and (3) are the setup constraints enforcing $Y_{t}$ to $1$ if a production is made at $t$. Finally, (4), (5) and (6) express the various costs. When $(L)$ is solved as a MILP, the variables $X$ and $I$ can be relaxed and considered real [22].

$(L)$ can also be modeled as a facility location problem [19]: the variables $I$ and $X$ are channeled to the variables $X_{tr}$, where $X_{tr}$ represents the proportion of demand $d_{r}$ produced in period $t$ and stored from $t$ to $r$. The model can be written as follows:

[TABLE]

Though the number of variables is increased, this model has the advantage to tighten the big M constraints (3) of the first formulation by stating constraints (12) and is known to provide an at least as good lower bound as the linear relaxation. Note that the shortest path reformulation (and $(l,S)$ inequalities) would provide an equivalent lower bound to UFL [22].

2.2 Linear relaxation

Solving the linear relaxation of MILP_AGG (i.e. $Y_{t}\in[0,1],\forall\>t\in\llbracket 1,T\rrbracket$) is equivalent to a minimum cost network flow problem [2]. The graph of this flow is presented in Figure 2. On each arc, $(u,c)$ represents the capacity ($u$) and unitary cost ($c$) of the arc. The units flow from source node $S$ to sink node $W$. On each production arc ($S$, $t$), the capacity is the production capacity and the cost is $p^{\prime}_{t}=p_{t}+\dfrac{s_{t}}{\overline{\alpha_{t}}}$. On each inventory arc ($t$, $t+1$), the capacity is the inventory capacity and the cost is $h_{t}$. Finally on each demand arc ($t$, $W$), there must be exactly $d_{t}$ units and the unitary cost is [math].

The flow problem can be solved in $O(T^{2})$ with the successive shortest path algorithm [2].

2.3 An equivalent problem without lower bounds

In this subsection, we show that problem $(L)$ is equivalent to a problem without lower bounds. It will allow us to re-use several classical lot-sizing algorithms and will also simplify the presentation of some results.

Solving a maximum flow problem on a network with lower bounds on flows is equivalent to solving a maximum flow problem on a transformed network without lower bounds as shown in [2]. We denote by $(L^{\prime})$ the resulting problem of that transformation applied to $(L)$ slightly adapted to take into account setup costs. The parameters of $(L^{\prime})$ are:

[TABLE]

$(X,I)$ is a solution of $(L)$ if and only if $(X^{\prime},I^{\prime})$ is a solution of $(L^{\prime})$. The intuition is that the production and inventory lower bounds are considered as mandatory quantities. As these quantities must be produced/stored at a precise period, no decisions have to be made about them and thus they can be removed from the problem. Note that the demands are also affected by the transformation.

The mandatory costs associated to the lower bounds are:

[TABLE]

and the variables of $(L)$ and $(L^{\prime})$ are linked as follows:

[TABLE]

Note that if the final demand $d^{\prime}_{t}$ is negative then $d_{t}+I_{t}<X_{t}+I_{t-1}$. This cannot be if constraints $\{(2),(3),(7),(8),(9)\}$ – which correspond to the feasibility of (L) – are bound consistent. We will show later in this paper (in section 5) that we can assume this property when removing the lower bounds.

2.4 Dynamic programming

$(L)$ can also be solved via DP [11]. We provide here the algorithm without lower bounds on production and inventory. The algorithm (called DPLS in the paper) iterates over the inventory levels. We denote $\textsl{f}\>(t,I_{t})$ as the minimum cost for producing the demands from $d_{1}$ to $d_{t}$ knowing that the stock level at $t$ is $I_{t}$:

[TABLE]

where $a=\max{\{0,d_{t}+I_{t}-\overline{\alpha_{t}}\}}$, $b=\min{\{\overline{\beta_{t-1}},d_{t}+I_{t}}\}$ and $X_{t}=I_{t}+d_{t}-I_{t-1}$. We define $I_{max}=\max\>\{\overline{\beta_{t}}\>|\>t\in\llbracket 1,T\rrbracket\}$. The initial states are $\textsl{f}\>(0,0)=0\text{ and }\forall\>I_{t}\in\llbracket 1,I_{max}\rrbracket\>,\>\textsl{f}\>(0,I_{t})=+\infty$. The value $\textsl{f}\>(T,0)$ gives the optimal cost of $(L)$. This dynamic programming algorithm runs in pseudo-polynomial time $O(TI_{max}^{2})$.

Note that DPLS consists in finding a shortest path in the graph for which there is a node for each inventory level at each period. The cost on an arc between two nodes $(t,I_{t})$ and $(t+1,I_{t+1})$ corresponds to the cost for satisfying demand $d_{t}$ and having an inventory level $I_{t+1}$ at the end of period $t+1$ knowing that there was an inventory level $I_{t}$ at the end of period $t$.

DPLS is a DP algorithm referred to as "forward" since it considers the periods in chronological order. We can also write the reverse (or "backward") DPLS. Let $\textsl{f}_{r}(t,I_{t})$ be the minimum cost for producing the demands from $d_{t+1}$ to $d_{T}$ knowing that the stock level at $t$ is $I_{t}$:

[TABLE]

where $c=\max{\{0,I_{t}-d_{t+1}\}}$, $d=\min{\{\overline{\beta_{t+1}},I_{t}-d_{t+1}+\overline{\alpha_{t+1}}\}}$ and $X_{t+1}=d_{t+1}+I_{t+1}-I_{t}$. The initial states are $\textsl{f}_{r}(T,0)=0$ and $\forall\>I_{t}\in\llbracket 1,I_{max}\rrbracket\>,\>\textsl{f}_{r}(T,I_{t})=+\infty$. As described above, $\textsl{f}_{r}(t,I_{t})$ can be seen as the shortest path from the node $(t,I_{t})$ to the node $(T,0)$. The value $\textsl{f}_{r}\>(0,0)$ gives the optimal cost of $(L)$.

3 A new lower bound for the single-item lot-sizing

In this section, we present a new lower bound for the total cost $C$ and how it can be adapted for the setup cost $C_{s}$. The general idea is to decompose $(L)$ into sub-problems, then to compute a lower bound on each of these sub-problems and finally combine them at best to find a global lower bound. We suppose here that $\underline{\alpha_{t}}=\underline{\beta_{t}}=0,\>\forall\>t\in\llbracket 1,T\rrbracket$. This assumption is not restrictive as production and inventory lower bounds can be easily removed in $(L)$ as shown in 2.3.

3.1 Lot-sizing sub-problem

A sub-problem $(L_{u,v})$, with $u<v$, is defined exactly as $(L)$ except that:

[TABLE]

The cost variables of $(L_{u,v})$ are denoted $C^{uv}$, $C_{p}^{uv}$, $C_{h}^{uv}$ and $C_{s}^{uv}$, corresponding to the total cost, sum of production costs, sum of holding costs and sum of setup costs of $(L_{u,v})$. Figure 3 illustrates the data used in sub-problem $(L_{u,v})$.

Sub-problem $(L_{1T})$ corresponds to the entire problem $(L)$. As there is no demand after period $v$ and no lower bounds of production, any solution of $(L_{u,v})$ is dominated by a solution with null inventory at the end of period $v$ ($I_{v}=0$). Note also that, in $(L_{u,v})$, some demands in $\{d_{u},d_{u+1},\ldots,d_{v}\}$ can be satisfied by a production made without setup cost before period $u$. Finally, an optimal solution of $(L_{u,v})$ provides a lower bound of the cost for satisfying the set of demands $\{d_{u},d_{u+1},\ldots,d_{v}\}$ in problem $(L)$. Indeed $(L_{u,v})$ is a relaxation of problem $(L)$.

There are $T(T-1)/2$ sub-problems and we order them by increasing end times first, then by increasing start times (see Table 1). Sub-problem $(L_{u_{i},v_{i}})$ will be referred to as sub-problem $i$.

Definition 1

Sub-problems $(L_{u,v})$ and $(L_{u^{\prime},v^{\prime}})$ are disjoint if

[TABLE]

3.2 Combining disjoint sub-problems provides a lower bound

For sub-problem $i$, we denote by $w_{i}$ a lower bound of its total cost $C^{u_{i}v_{i}}$. Disjoint sub-problems can be combined to obtain a lower bound for the total cost of $(L)$.

Theorem 1

For any set $S$ of disjoint sub-problems, we have

[TABLE]

Proof. Let $E^{*}$ be an optimal production plan for $(L)$ of cost $C^{*}$ and $S$ be a set of disjoint sub-problems. We will build from $E^{*}$ a feasible solution to each sub-problem $i$ of $S$ and prove that their costs add up to less than $C^{*}$.

Consider a sub-problem $i$ in $S$. For each demand $d_{t}$, $t\in\llbracket u_{i},v_{i}\rrbracket$ we produce $d_{t}$ at the same periods as it is produced in $E^{*}$ (with a First Come First Served policy). We obtain by this process a feasible solution to sub-problem $i$ and denote its cost by $K_{i}$.

As the sub-problems in $S$ are disjoint and we keep the same production orders, the sum of setup costs paid in all of these sub-problems is less than or equal to the sum of setup costs paid in $E^{*}$. The production and inventory costs are identical to the costs paid in $E^{*}$ for the demands included in $\cup_{i\in S}\llbracket u_{i},v_{i}\rrbracket$. It follows that $\sum_{i\in S}{K_{i}}\leq C^{*}$. Finally, as for each sub-problem $i$, $w_{i}$ is a lower bound of $C^{u_{i}v_{i}}$, we get: $\sum_{i\in S}{w_{i}}\leq\sum_{i\in S}{K_{i}}$. ◻

3.3 Combining lower bounds at best

Given a lower bound for each sub-problem, we wish to find the best lower bound of $C$, i.e. to determine the set $S$ of disjoint sub-problems that maximizes $\sum_{i\in S}{w_{i}}$.

This problem can be seen as a Weighted Interval Scheduling Problem (WISP) which can be solved in $O(n\log(n))$ where $n$ is the number of intervals [18]. The algorithm sorts the intervals in $O(n\log(n))$ and then applies a DP that runs in $O(n)$. In our case, there are $n=T(T-1)/2$ intervals (sub-problems) which are already sorted and the DP algorithm (called DPWisp in the paper) runs in $O(T^{2})$.

We use the indexing of intervals given in Table 1 and we denote by $\textsl{wisp}(i)$ the maximal weight that can be achieved using the $i$ first intervals. The forward DP writes as

[TABLE]

where $prec_{i}$ is the biggest integer, smaller than $i$ ($prec_{i}<i$), such that the intervals $prec_{i}$ and $i$ are disjoint. For each sub-problem $i$ such that $u=1$, we define $prec_{i}=0$. Hence $prec_{i}$ is the first interval before the $i^{th}$ one that is disjoint with it. For instance, $[3,4]$ is the $6^{th}$ interval and $prec_{6}=1$ since $[1,2]$ and $[3,4]$ are disjoint while $[2,3]$ and $[3,4]$ are not. The value $\textsl{wisp}(n)$ is then a lower bound of the global cost $C$.

We can also write the reverse version DPWisp (with the sub-problems considered backwards). The sub-problems are now sorted by decreasing start times first, then by decreasing end times.

[TABLE]

where $succ_{i}$ is the biggest integer, greater than $i$ ($succ_{i}>i$), such that the intervals $succ_{i}$ and $i$ are disjoint.

3.4 Computing lower bounds for sub-problems

In our algorithms, we will solve exactly the sub-problem by DP when the size is reasonable and solve the linear relaxation otherwise. The optimal cost of sub-problem $(L_{u,v})$ can be computed with DPLS (see Section 2.4) applied to the periods $u$ to $v$. We need to pre-compute $\textsl{f}\>(u-1,q_{u-1})$ for $q_{u-1}\in\llbracket 0,\min\>\{\sum_{t=u}^{v}{d_{t}},\>\overline{\beta_{u-1}}\}\>\rrbracket$. This can be done with a greedy algorithm which determines the cheapest periods in order to produce the requested quantity and to store it until $u$. The DPLS applied to a sub-problem $(L_{u,v})$ runs then in $O((v-u+1)(I_{max}^{uv})^{2})$ where $I_{max}^{uv}=\max\{\overline{\beta_{t}}\>|\>t\in\llbracket u,v\rrbracket\}$.

The linear relaxation can be seen as a minimum cost flow problem (see Section 2.2) and can be solved in $O(T^{2})$.

3.5 Adaptation to a lower bound on setup costs

We can re-use this previous approach to obtain a lower bound on the setup cost variable $\mathit{Cs}$. Sub-problems are defined similarly except that we remove unitary production costs and inventory costs ($p_{t}=h_{t}=0,\forall t$).

Note that we could use the same approach for $\mathit{Cp}$ and $\mathit{Ch}$. However, it is not necessary as $(L)$ is polynomial when removing setup costs.

4 The lot-sizing global constraint

In this section, we provide some CP background before presenting the LotSizing global constraint. We also study the complexity of achieving different consistency levels.

4.1 Constraint programming background

A Constraint Satisfaction Problem (CSP) [15] consists of a set of variables, with a finite domain of values for each variable, and a set of constraints on these variables. Upper cases are used for variables (e.g. $V_{i}$) and lower cases for values (e.g. $v_{i}$). We denote by $D(V_{i})$ the domain of variable $V_{i}$ and by $\overline{V_{i}}$ (resp. $\underline{V_{i}}$) the minimum (resp. maximum) value in $D(V_{i})$.

Let $c$ be a constraint on variables $\langle V_{1},\ldots,V_{n}\rangle$. A support for $c$ is a tuple $\langle v_{1},\ldots,v_{n}\rangle$ which satisfies $c$ and such that $v_{i}\in D(V_{i})$ for each variable $V_{i}$. A bound support is a tuple $\langle v_{1},\ldots,v_{n}\rangle$ which satisfies $c$ and such that $\underline{V_{i}}\leq v_{i}\leq\overline{V_{i}}$ for each $V_{i}$.

A variable $V_{i}$ is arc consistent (AC) for constraint $c$ if each value of $D(V_{i})$ belongs to a support for $c$. A variable $V_{i}$ is bound consistent (BC) for constraint $c$ if $\underline{V_{i}}$ and $\overline{V_{i}}$ belong to a bound support for $c$. A variable $V_{i}$ is range consistent (RC) for constraint $c$ if each value of $D(V_{i})$ belongs to a bound support for $c$. A constraint $c$ is AC (resp. BC, RC) if all its variables are AC (resp. BC, RC). A CSP problem is AC (resp. BC, RC) if each constraint is AC (resp. BC, RC).

The following example illustrates the three notions of AC, BC and RC. Consider the following linear constraint over two integer variables $x$ and $y$:

[TABLE]

The three levels of consistency are applied to the example and showed in table 2. The bound consistent domains are $\{2,4\}$ and $\{4,5,6,7,8\}$ since we just check if the bounds belong to a bound support: only $1$ can be removed from the domain of $x$ as $\langle 1,2\rangle$ is not a bound support. For instance, the bound support for $y=8$ is $\langle 4,8\rangle$. The range consistent domains are $\{2,4\}$ and $\{4,6,8\}$ since values $5$ and $7$ do not belong to a bound support. Indeed $\langle 2.5,5\rangle$ and $\langle 3.5,7\rangle$ are not bound supports as the values in a bound support must be in $\mathbb{Z}$. The value $6$ for $y$ is range consistent since $\langle 3,6\rangle$ is a bound support. Finally the arc consistent domains are $\{2,4\}$ and $\{4,8\}$ as we remove all inconsistent values.

4.2 Definition

We formally define here the global constraint LotSizing. This constraint is stated on the variable vectors $X=\langle X_{1},\ldots,X_{T}\rangle$, $I=\langle I_{1},\ldots,I_{T}\rangle$, $Y=\langle Y_{1},\ldots,Y_{T}\rangle$ and the four cost variables $\mathit{Cp}$, $\mathit{Ch}$, $\mathit{Cs}$, $C$ of $(L)$. The data of the problem is denoted by $data=\{(p_{t},h_{t},s_{t},d_{t},\underline{\alpha_{t}},\overline{\alpha_{t}},\underline{\beta_{t}},\overline{\beta_{t}})\>|\>t\in\llbracket 1,T\rrbracket\}$.

Definition 2

$\textsc{LotSizing}(X,I,Y,\mathit{Cp},\mathit{Ch},\mathit{Cs},C,data)$* has a solution if and only if there exists a production plan, solution of $(L)$ that satisfies:*

[TABLE]

The LotSizing global constraint has a solution if and only if the set of constraints $\{\eqref{flotMILP1},\ldots,\eqref{domYMILP1},\eqref{ubCp},\ldots,\eqref{ubC}\}$ has a solution.

4.3 Complexity

This subsection presents theorems on the complexity of achieving BC, RC or AC on LotSizing and one of its restrictions that does not take into account the costs and focuses on the flow equations. Let us first give a property on the complexity of $(L)$.

Property. *Problem $(L)$ with $h_{t}=p_{t}=\underline{\beta_{t}}=\underline{\alpha_{t}}=0$ and $\overline{\beta_{t}}=+\infty$ is NP-hard [12].

Theorem 2

Achieving BC for LotSizing can be done in pseudo-polynomial time.

Proof. BC for LotSizing can be achieved by solving a Shortest Path Problem with Resource Constraints (SPPRC) in the graph of DPLS (see Section 2.4) where the resources are the three intermediate costs of respective capacities $\overline{Cp},\overline{Ch},\overline{Cs}$ and the global cost $C$ is the objective. Finding the shortest path in this graph while respecting the three resources at the final node gives a bound support regarding the three intermediate costs. SPPRC is known to be weakly NP-hard [10]. ◻

We denote by $(L_{r})$ the feasibility problem associated to $(L)$.

[TABLE]

The constraints of $(L_{r})$ describes the dynamics of the problem, without considering costs.

Theorem 3

Achieving BC on $(L_{r})$ can be done in $O(T)$.

Proof. Figure 4 represents the constraint network $(L_{r})$ as well as the corresponding intersection graph. The rectangular-shaped constraints are the flow balance constraints and the dashed oval ones are the setup constraints. The intersection graph is built as follows: we set a vertex for each constraint and two vertices are linked if and only if the corresponding constraints have at least one variable in common. As the intersection graph is acyclic and each pair of constraints has at most one variable in common, the constraint network is Berge-acyclic. It is known that if we filter each constraint of a Berge-acyclic constraint network in an appropriate order then each constraint needs only to be woken twice in order to reach the fix-point [20]. Each constraint of the network can be filtered in $O(1)$, hence we can achieve BC on $(L_{r})$ in $O(T)$.

◻

In order to investigate RC for $(L_{r})$, we consider $(F)$ the more general problem of finding an integer flow in a directed graph $G=(V,E)$:

[TABLE]

where $x_{ij}$ is the flow going from node $i$ to node $j$, $\delta^{+}_{i}$ (resp. $\delta^{-}_{i}$) is the set of successor (resp. predecessor) nodes of node $i$ and $b_{i}\in\{-1,0,1\}$.

Theorem 4

BC and RC are equivalent for $(F)$.

Proof. By definition, RC implies BC for $(F)$. Conversely, assume BC for $(F)$.

In order to show the converse, we show that if $k\in\llbracket l_{ij},u_{ij}\rrbracket$, then $k$ belongs to a bound support for $(F)$.

Let $i_{0}$ be a direct predecessor of $j_{0}$ in the graph. BC implies that there exists $x^{\prime}$ (resp. $x^{\prime\prime}$) an integer solution of (F) such that $x^{\prime}_{i_{0}j_{0}}=l_{i_{0}j_{0}}$ (resp. $x^{\prime\prime}_{i_{0}j_{0}}=u_{i_{0}j_{0}}$). Let $\gamma\in[0,1]$ such as $k=\gamma\>l_{i_{0}j_{0}}+(1-\gamma)u_{i_{0}j_{0}}\in\mathbb{N}$. Let’s show that there exists an integer solution of $(F)$ where $x_{i_{0}j_{0}}=k$. We modify the bounds $\>\forall\>(i,j)\in E$:

[TABLE]

We can then show that $\forall\>i,j\>\>\widehat{l_{ij}}\geq l_{ij}$ and $\widehat{u_{ij}}\leq u_{ij}$. The constraint matrix of $(F)$ is totally unimodular since it is a flow problem, hence the application of the simplex algorithm to $(F)$ with the updated bounds gives an integer solution $x^{\prime\prime\prime}$ where $x_{i_{0}j_{0}}=k$. ◻

As constraints (2) are flow constraints and the setup variables $Y_{t}$ are binary, the following result follows.

Corollary 1

Achieving BC on $(L_{r})$ is equivalent to achieving RC on $(L_{r})$.

When there are no holes in the domains of $X$ and $I$, RC is equivalent to AC for $(L_{r})$. Note that there may exist holes in lot-sizing problems when considering batching constraints for instance.

5 Filtering the LotSizing constraint

This section describes the filtering of the LotSizing constraint. It also gives some implementation details to improve the incrementality of the global constraint.

Algorithm 1 gives an overview of the main filtering steps of LotSizing. Each step refers to the corresponding section for detailed explanations. When all setup variables are instantiated the problem amounts to a minimum cost flow problem (lines 1-2). If not, the general case is as follows. Firstly, the problem is transformed by removing all lower bounds (line 4) as LotSizing is defined with lower bounds and these can increase during the search. Secondly production and inventory costs lower bounds are computed (lines 5-6). Thirdly, when the overall problem is of reasonable size (lines 7-11) the remaining filtering is performed using dynamic programming. If not, the WISP relaxation is used and filtering is performed via the WISP support (lines 13-17).

5.1 Filtering when the setup variables are instantiated

When all the setup variables $Y$ are instantiated, problem $(L)$ becomes polynomial and amounts to a minimum cost flow problem. Solving a minimum cost flow problem on the flow graph presented in Figure 2 finds a solution to $(L)$ that minimizes $C$. This allows the user to branch only on the $Y$ variables since the solver can instantiate all the other variables in polynomial time when the setup variables are instantiated. Note that when using Lotsizing in a more complex model (with multiple Lotisizing or with additional constraints) the resulting problem, when all the $Y$ are instantiated, may not be polynomial. Therefore we let the user specify when stating the constraint if this property holds or not. Note also that a minimum cost flow dedicated filtering algorithm [24] can be used at this stage. Either of these two options allows the user to branch only on the $Y$ variables.

5.2 Filtering cost lower bounds

Lower bounds of the cost variables are computed as follows:

A lower bound on the production cost $\mathit{Cp}$ is computed by solving a minimum cost flow problem on the graph presented in §  2.2 considering only the variable production costs.

• 2.

A lower bound on the inventory cost $\mathit{Ch}$ is computed by solving a minimum cost flow problem on the graph presented in §  2.2 considering only the variable inventory costs.

• 3.

A lower bound on the global cost $\mathit{C}$ is computed using DPLS. The lower bound is given by the value $\textsl{f}\>(T,0)$.

• 4.

A lower bound on the setup cost $\mathit{Cs}$ is computed via dynamic programming as well. We consider here the problem $(L)$ without production or inventory costs. As mentioned in the literature review, this problem can be solved using the traditional knapsack dynamic programming algorithm slightly adapted to take into account the inventory upper bounds (we call this algorithm DPKnap).

5.3 Filtering X and I via dynamic programming

DPLS gives a lower bound of the global cost $C$. In order to filter the variables we use the tables created by DPLS and reverse DPLS. Remember that in the graph described in 2.4, $\textsl{f}\>(t,I_{t})$ can be seen as the shortest path from the node $(0,0)$ to $(t,I_{t})$ and $\textsl{f}_{r}(t,I_{t})$ is the shortest path from $(t,I_{t})$ to $(T,0)$. We filter each value $i_{t}$ in the domain of $I_{t}$ in $O(TI_{max})$:

[TABLE]

We filter each value in the domain of $X_{t}$ in $O(TI_{max}^{2})$:

[TABLE]

where $x_{t}=d_{t}+i_{t}-i_{t-1}$ and $cost(t,i_{t-1},i_{t})=\mathds{1}_{X_{t}>0}s_{t}+p_{t}x_{t}+h_{t}i_{t}$.

5.4 Scaling the filtering based on dynamic programming

In the case that DPLS has memory issues on the overall problem $(L)$, we solve a WISP (see Section 3.2) to find a lower bound on $C$. We can adapt the filtering rules (20) and (21) on the sub-problems of reasonable size. In order to compare the shortest paths to the global upper bound $\overline{C}$, we need to have a lower bound on the cost of the production outside the sub-problem. We use DPWisp and its reverse version to do so. We can then define:

$\textrm{lbBefore}(t)=\textsl{wisp}\>(\frac{t(t-1)}{2})$ which is the best bound we can get by combining the sub-problems ending by at most $t$. It is a lower bound on the satisfaction of the demands $d_{1}$ to $d_{t}$. We set $\textrm{lbBefore}(0)=0$.

• 2.

$\textrm{lbAfter}(t)=\textsl{wisp}_{r}(\frac{(T-t)(T-t+1)}{2})$ which is the best bound we can get by combining the sub-problems starting at $t$. It is a lower bound of the cost for satisfying the demands from $d_{t}$ to $d_{T}$. We set $\textrm{lbAfter}(T+1)=0$.

Figure 5 represents the different lower bounds computed while filtering the value $i_{t}$ for the variable $I_{t}$. $I_{t}$ belongs to the sub-problem $(L_{u,v})$, $lbBefore(u-1)$ and $lbAfter(v+1)$ are computed via DPWisp outside that sub-problem.

We filter the variables via the two following rules:

[TABLE]

For all the sub-problems of the support (or solution) of the WISP, we can compute DPLS and its reverse version. Note that here $\textsl{f}\>(t,i_{t})$ and $\textsl{f}_{r}(t,i_{t})$ come from DPLS applied to the sub-problems. Hence when considering sub-problem $(L_{u,v})$ and $u\leq t\leq v$, $\textsl{f}\>(t,i_{t})$ is a lower bound of the cost for satisfying demands $d_{u}$ to $d_{t}$ and having $I_{t}=i_{t}$. Similarly $\textsl{f}_{r}(t,i_{t})$ is a lower bound on the cost for satisfying the demands $d_{t+1}$ to $d_{v}$ and having $I_{t}=i_{t}$. Moreover, since the sub-problems consider only the demands $d_{u}$ to $d_{v}$, we cannot filter values greater than or equal to $\sum_{k=t+1}^{v}{d_{k}}$ for $I_{t}$ and $X_{t}$. Indeed, greater values might be used to satisfy demands outside the sub-problem and are not considered when computing DPLS.

Note that although the scaling may affect the quality of the filtering, it is a pragmatic rule applied to avoid wakening costly propagation when little filtering is expected.

5.5 Adaptation to take into account the setup cost

In order to adapt the filtering to take into account $\overline{\mathit{Cs}}$, we do not take into account the production and inventory costs. The resulting problem can then be tackled by dynamic programming with DPKnap. The variables can be filtered via the WISP support with the DPKnap computed on the sub-problems. The filtering rules are applied with the upper bound of $\mathit{Cs}$.

6 Numerical results on the single-item lot-sizing problem

This section validates our global constraint and the filtering mechanisms described above. We compare the performance of LotSizing to four other methods on the single-item lot-sizing problem.

Five methods

The five methods that solve the single-item lot-sizing are:

A basic CP model (CP_Basic), which is a decomposition of the single-item lot-sizing problem basically equivalent to the MILP model with the implication constraints $X_{t}>0\Rightarrow Y_{t}=1,\forall\>t\in\llbracket 1,T\rrbracket$ instead of the setup constraints (3)

• 2.

A CP model with our LotSizing global constraint (CP_LS)

• 3.

The dynamic programming algorithm presented in 2.4 (DP)

• 4.

The classical aggregated MILP model (MILP_AGG)

• 5.

The facility location MILP model (MILP_UFL)

The MILP models were implemented with CPLEX version 12.6 and the CP models in Choco 3.3 [23].

Branching heuristics for the CP models

A default branching heuristic is used to instantiate the variables to their lower bounds in a lexicographic order (chronological order here). The property described in 5.1 is valid for the single-item lot-sizing problem. For the sake of comparison, the same improvement is done for CP_Basic. The search space is thus restricted to the setup variables ($Y$) for both CP models.

Cost upper bound

Since we want to assess the quality of the filtering and LotSizing uses cost-based filtering, we choose to have the best possible upper bound on the global cost at the start of the resolution: the optimal cost. This means that our models still have to find the optimal solution. By setting the initial upper bound to the optimal value, we simply aim to avoid the issue of finding a good enough initial solution that would "activate the filtering" we want to assess. It gives a simple and identical set-up to all compared approaches and allows us to focus the analysis on the filtering we have been investigating. In a more realistic setting, the model CP_LS can be used to find upper bounds.

Instance parameters

The single-item instances are generated based on the parameters $d_{\text{avg}}$, $e$, $\delta$, $\theta$, $\lambda$, and $T$ as follows:

The inventory costs are constant and equal to $1$ (i.e. $h_{t}=h=1$).

• 2.

The setup and the unitary production costs are generated using two parameters: $e$ and $\theta$.

• (a)

$e$ represents the overall unitary production cost (i.e. the unitary production cost if the production capacity is saturated: $p_{t}\>\overline{\alpha_{t}}+s_{t}$ divided by $\overline{\alpha_{t}}$). We set $e=10$.

• (b)

$\theta\in[0,1]$ represents the portion of the setup cost to the unitary production cost. The overall production cost at $t$ (i.e. $e\>\overline{\alpha_{t}}$) will be imputable for $\theta$ to its setup cost and for $1-\theta$ to the unitary production cost at $t$. For each period, $\theta$ is uniformly randomly set in the interval $[0,1]$.

• 3.

The demand is uniformly randomly generated in the interval $[d_{\text{avg}}-\delta,d_{\text{avg}}+\delta]$.

• 4.

The production and inventory capacities are constant and equal to $\lambda\>d_{\text{avg}}$.

For each problem, we give the set of parameters that were used to generate the instances. Each class of instances contains 10 instances.

Experimentation setup

All the tests are run under Windows 8 on an Intel Core i5 @ 2.5 GHz with 12GB of RAM. We set a time limit of 200s and a memory limit of 4GB of RAM. The indicator NODE is the average number of nodes computed by each model on the class. CPU corresponds to the average CPU time used by the models. RNB is the average gap of the root node lower bound to the optimal value. LR is the average gap of the linear relaxation to optimal. Finally OPT is the number of solved instances in the class. The means are computed over all the instances of each class.

6.1 Single-item lot-sizing

The five instance classes are:

C1LS : $d_{\text{avg}}=1000$, $\delta=100$, $\theta\in[0.8,1]$, $\lambda=3$, $T=40$

• 2.

C2LS : $d_{\text{avg}}=1000$, $\delta=500$, $\theta\in[0.4,0.6]$, $\lambda=3$,$T=40$

• 3.

C3LS : $d_{\text{avg}}=1000$, $\delta=100$, $\theta\in[0.8,1]$, $\lambda=3$, $T=80$

• 4.

C4LS : $d_{\text{avg}}=1000$, $\delta=500$, $\theta\in[0.4,0.6]$, $\lambda=3$, $T=80$

• 5.

C5LS : $d_{\text{avg}}=1000$, $\delta=50$, $\theta=0.5$, $\lambda=3$, $T=40$

The results are presented in Tables 3 and 4.

These tables show that:

As expected, the basic CP model has a very large search space as it does not propagate any strong reasoning. We therefore did not use CP_Basic for the following results.

• 2.

As the upper bound provided is optimal and there is no upper bound on $\mathit{Cp}$,$\mathit{Ch}$ and $\mathit{Cs}$, CP_LS achieves AC at the root node and branches backtrack free towards an optimal solution.

• 3.

The linear relaxation of MILP_AGG is not as good as the MILP_UFL’s as it was expected due to the setup constraints (3). CPLEX however provides a better root node lower bound for MILP_AGG. MILP_UFL is not as competitive as MILP_AGG because of the number of variables and constraints. We therefore did not use MILP_UFL for the following results.

6.2 Scaling the global constraint

We then test the WISP support filtering described in 5.4 when the DP has memory issues. In order to generate memory issues for the DP, we add high consumption peaks in the instances. The peaks are added in periods 6 to 9, 12 to 15, 22 to 25 and 32 to 36 and correspond to demands of 50,000. When computing the global lower bound with the WISP, no sub-problem containing a demand peak is solved via dynamic programming since the peaks increase the complexity of the DP. The lower bound on these sub-problems is therefore their linear relaxation. The five instance classes have the following parameters:

C1Peaks : $d_{\text{avg}}=100$, $\delta=50$, $\theta\in[0.8,1]$, $\lambda=4$, $T=40$

• 2.

C2Peaks : $d_{\text{avg}}=100$, $\delta=50$, $\theta\in[0.4,0.6]$, $\lambda=4$,$T=40$

• 3.

C3Peaks : $d_{\text{avg}}=100$, $\delta=50$, $\theta=0.5$, $\lambda=4$, $T=40$

• 4.

C4Peaks : $d_{\text{avg}}=100$, $\delta=20$, $\theta\in[0.8,1]$, $\lambda=4$, $T=40$

• 5.

C5Peaks : $d_{\text{avg}}=100$, $\delta=20$, $\theta\in[0.4,0.6]$, $\lambda=4$, $T=40$

The branching heuristic is adapted to select first the setup variables of the high demand periods. Table 5 compares the three models CP_LS, MILP_AGG and DP on these big instances.

This table shows that:

The filtering is lighter, hence the root node lower bound gap increases as well as the number of nodes.

• 2.

The resolution is however faster than the DP.

• 3.

Although the linear relaxation degrades, CPLEX pre-processing behaves very well as shown by the root node lower bound.

7 Single-item lot-sizing with side constraints

We consider the single-item lot-sizing problem $(L)$ with three side constraints (domain disjunction, limited production rate and a combination of the two). The instances created here are generated the same way as before and we added leveled production and/or constrained production rate. For the following tests, we compared only CP_LS to MILP_AGG and to DP when it is relevant.

7.1 Disjunctive production constraints

We consider here that the production is leveled. The domains of each variable $X_{t}$ is defined by a disjunction of $n_{t}$ integer intervals $K_{k}=\llbracket\underline{K_{k}},\overline{K_{k}}\rrbracket,\>\forall\>k\in\llbracket 1,n_{t}\rrbracket$:

[TABLE]

DPLS can take into account the disjunctions without any loss of complexity. We add the following constraints to the MILP model MILP_AGG:

[TABLE]

The ten classes for this problem are:

C1Disj : $d_{\text{avg}}=100$, $\delta=50$, $\theta\in[0.8,1]$, $\lambda=5$, $T=40$

• 2.

C2Disj : $d_{\text{avg}}=100$, $\delta=60$, $\theta\in[0.4,0.6]$, $\lambda=5$,$T=40$

• 3.

C3Disj : $d_{\text{avg}}=100$, $\delta=70$, $\theta\in[0.3,0.8]$, $\lambda=5$, $T=40$

• 4.

C4Disj : $d_{\text{avg}}=100$, $\delta=30$, $\theta\in[0.6,1]$, $\lambda=5$, $T=40$

• 5.

C5Disj : $d_{\text{avg}}=100$, $\delta=50$, $\theta\in[0.9,1]$, $\lambda=5$, $T=40$

We generated C6Disj, C7Disj, C8Disj, C9Disj and C10Disj that have the same parameters than the five instances above, but with $T=80$. The disjunctions are added as follows: $D(X_{t})=\llbracket 0,30\rrbracket\cup\llbracket 100,150\rrbracket\cup\llbracket 200,240\rrbracket$. Table 6 gives the numerical results for the single-item lot-sizing with disjunctions.

Table 6 shows that the CP and DP models are very fast to solve these instances. The property described in 5.1 concerning the setup variables is not valid for this problem: indeed the flow with disjunctions is not polynomial. However as LotSizing’s filtering uses the DP, the global constraint can handle disjunctions on the domains of the production variables. Therefore CP_LS achieves AC at the root node and branches backtrack free towards an optimal solution. Unsurprisingly we note that the MILP model does not handle these disjunction constraints well.

7.2 Q/R constraints

Q/R constraints are interesting side constraints for single-item lot-sizing problems [13, 14]. They relate to the production rate and state that, given two integers $Q$ and $R$, there must be at least $Q$ and at most $R$ periods between two consecutive productions. Dynamic programming rapidly gets memory issues here, as the states should take into account what happened at least $R$ periods before. The Q/R constraints can be modeled by two Sequence constraints stated as follows:

[TABLE]

The Sequence constraint is defined as follows [6]: $\textsc{Sequence}(l,u,k,[Z_{1},\ldots,Z_{n}],v)$ holds if and only if:

[TABLE]

We add the following constraints to the model MILP_AGG:

[TABLE]

We add the $\#_{uv}$ variables that count the number of effective production periods between period $u$ and period $v$ included:

[TABLE]

These variables enable us to use the encoding of Sequence presented in [8] to propagate the Q/R constraints. We also add the useful following redundant constraints:

[TABLE]

The ten classes (C1QR, $\ldots$, C10QR) for this problem have the same parameters than C1Disj, $\ldots$, C10Disj to which we add (Q=2, R=6) for classes 1, 2, 3, 6, 7, 8 and (Q=3, R=7) for classes 4, 5, 9, 10. Table 7 compares CP_LS to MILP_AGG on the instances with Q/R.

The linear relaxation and root node lower bound of MILP_AGG has slightly worsened without degrading the performance of the model. CP_LS stays competitive on most of the instances.

7.3 Disjunctive with Q/R constraints

We add both Q/R and disjunctive production constraints. The problem cannot be tackled via DP due to the Q/R constraints, hence we compared CP_LS to MILP_AGG. The instances have the same parameters than before with both the disjunctions and the Q/R parameters presented for the latter problems. The results are shown in table 8.

On some classes, CP_LS does not solve all the instances yet is competitive compared to MILP_AGG and has a near optimal root node lower bound.

We summarize below the main conclusions of our numerical study that hold for the set of instances under consideration.

For the single-item problem:

CP_Basic fails to find the optimal solution in a reasonable time due to the size of the search space and the lack of pruning. As no information is given on the costs, its root node lower bound is very far from the optimal value.

• 2.

CP_LS achieves arc consistency at the root node since it is given the optimal upper bound and is a competitive approach to find the optimal solution. It still has to branch when there exist multiple optimal solutions (see for instance classes C3LS and C4LS). CP_LS can outscale DP on instances with demand peaks. The AC is not achieved at the root node since the scaling is based on a relaxation of the filtering mechanisms.

• 3.

DP and MILP_AGG are the most competitive approaches.

• 4.

MILP_UFL finds the best RNB of the two MILP models but is slower overall.

For the single-item problem with additional constraints:

Disjunctive constraints: CP_LS and DP outperform MILP which does not handle well combinatorial constraints.

• 2.

Q/R constraints: CP_LS stays competitive on most instances while DP is not a suitable approach.

• 3.

Q/R and Disjunctive constraints: CP_LS is the fastest on most instance classes. As for C7 and C8 classes, note that MILP_AGG was not very troubled by the disjunctive constraints (Table 5) nor by the Q/R constraints (Table 6) unlike CP_LS. Hence the difficulty of CP_LS to rapidly solve these instances (Table 7).

8 Conclusion

In this paper, we defined a global constraint LotSizing for a capacitated single-item lot-sizing problem. Firstly, we presented a new lower bound for this problem, based on a new decomposition of the problem into sub-problems. Secondly, we formally introduced our constraint and gave some complexity results. Thirdly, we developed filtering rules for the LotSizing global constraint based on dynamic programming. Finally, we presented a proof of concept for the filtering of the constraint via several numerical results. We can conclude that our approach based on constraint programming can yield interesting and competitive results for lot-sizing problems with side constraints. We however want to point out the limits of our numerical study that lie in the small size and variability of the set of instances that we chose for these first tests of the LotSizing global constraint.

The next step of this work will be to use the LotSizing global constraint as a building block to tackle multi-item and multi-echelon problems. In multi-item, each item can be modeled as a LotSizing constraint and each of them infers on the feasibility and costs of how to produce its item. They however share some variables if for instance we consider shared setup costs, the setup variables would be shared among the LotSizing constraints, communicating information. In multi-level problems, each retailer can be modeled as a LotSizing constraint. The global constraints are linked since the input of a retailer (production variables) is the output of the previous retailer (demands). The constraint programming framework built around the LotSizing global constraint might very well benefit from several work that has been done on the relaxation of these types of problems [28, 31].