Ujjwal Pratap

Leo Senique

Carlos Canudas-De-Wit
email: carlos.canudas-de-wit@gipsa-lab.fr

GTL-Healthmob: Simulation platform for urban mobility and epidemic control

Keywords: Simulation platform, Urban human mobility, Epidemic-mobility, Optimal lockdown, Grenoble metropolis

The GTL-Healthmob is an online simulation platform which aims at visualizing the population movement, simulating epidemic propagation and optimizing the mobility restrictions to limit the epidemic spread in the Grenoble metropolis. The platform takes user inputs for different epidemic and mobility control parameters and shows the results in terms of plots and heatmaps overlaid onto the map of the metropolis. The platform is being developed in the context of the ERC project Scale-FreeBack and the INRIA project Healthy-Mobility, to showcase the mobility and epidemic models, developed within the team by adapting and implementing them to the Grenoble metropolis, which can also enable policymakers to devise optimal mobility restriction policies.

Introduction

Although human mobility plays a vital role in a country's economy, it facilitates the spread of disease in the event of an epidemic by allowing contact between infected and susceptible populations. Therefore, controlling human mobility during an epidemic is a fundamental issue faced by policymakers. Such control can only be done optimally if human mobility is adequately modeled at the scale of a city or metropolis. Although, there is a large literature concerning epidemic spread and control of mobility, most of which focuses on control inter-city mobility or the epidemic spread in a single location. Although research exist for control of mobility for epidemic mitigation for urban cities, their implementation and extension on a large-scale urban mobility network is rarely available.

In this paper, we present and describe the GTL-Healthmob simulation platform, available through the link http://gtlville.inrialpes.fr/covid-19. It has been developed as an academic platform where users can a) simulate different scenarios of epidemic spread, b) simulate different scenarios of mobility restrictions and visualize its effect on epidemic spread c) simulate the optimal mobility restriction aimed at limiting epidemic spread and maximizing economy. This platform is part of the Grenoble Traffic Lab (GTL), which considers the Grenoble metropolis as a case study. Grenoble is metropolitan city in France located in the valley of Alpes mountain range. A metropolis is an administrative entity in France, in which several communes cooperate, and which has the right to levy local taxes. Grenoble metropolis is centered on the city of Grenoble. The considered area includes 55 communes including all the 49 communes in Grenoble metropolis and 6 communes in Gresivaudan, which are of economic importance to the Grenoble metropolis. The population of this area is about 500000 and has a surface area of about 600 km 2 . This paper is organized as follows: Section 2 presents an overview of the structure and main modules and submodules of the platform. Section 3 describes the methods and models used to calculate the evolution of epidemics and mobility, optimal policies based on the scientific work performed by the research team. Section 4 lists the datasets and their sources used to calibrate the models. Section 5 describes the user interface, features, and functionality available in the platform, and how they can be accessed by the general public. Finally, Section 6 concludes the paper. • The Web Interface module allows to visualise the mobility network, origin and destination nodes and also facilitates the users to choose input parameters and visualise the computed results.

Platform Overview

• The Computation module takes the parameters given as inputs by the users on the Web Interface module and static data from the Storage module to simulate the mobility and epidemic models and compute the evolution of mobility and epidemic spread. The results is then passed on to the Web interface for visualization. It is composed of three submodules:

-Mobility -Using different control and input parameters provided by the user, it computes the number of people present at each location at different times for each day during the simulation period.

-Epidemic -Using the epidemic parameters provided by the user and the number of people present at any location, given by mobility submodule, it computes the evolution of epidemic spread at each location.

-Optimization -Using the results from mobility and epidemic submodules, it computes the optimal mobility restriction that keeps maximum infection below the user defined threshold.

A more detailed description of the Computation submodules and the collected input data is presented in the following sections.

3 Methods and models

Mobility and epidemic spread

The mobility model which has been implemented on GTL-Healthmob to compute the mobility patterns and epidemic spread is presented in [START_REF] Pratap | Where, when and how people move in largescale urban networks: the Grenoble saga[END_REF] in discrete-time. Originally, this model was presented in [START_REF] Niazi | Optimal Control of Urban Human Mobility for Epidemic Mitigation[END_REF] in continuous-time. This model captures the mobility in an urban area between different locations and also integrates epidemic spread by coupling the mobility process with an SIR(Susceptible, Infected and Recovered) compartmental model at each location as depicted in Figure 2. The locations are of two types: origins and destinations. The origins correspond to the residential areas and the destinations correspond to the places of interest such as companies, research centers, schools, restaurants, cinemas, etc. The destinations have been grouped in different categories and subcategories. Every origin has a population and every destination has a capacity. The origins and destinations are connected by an origin-destination network. The urban human mobility is modeled on the network edges between different locations through flows and the epidemic spread is modeled locally at each location according to SIR model. Every day a certain number of people visit the destinations during specified time intervals and then return later on the same day and the epidemic spread at each location happens according to the SIR model [START_REF] Niazi | Optimal Control of Urban Human Mobility for Epidemic Mitigation[END_REF]. The daily movements of people are captured through time dependent functions called demand gating profiles and supply gating functions which represent the demand and supply of each location on a given day of the week. Moreover, the mobility between origins and destinations can be controlled by tuning either the capacity of destinations or by tuning their opening and closing times. In this work, these controls have been applied according to the destinations categories. This model can take as input the simulation duration, gating profiles, capacity control and schedule control(enforced closing time) for different destination categories, initial infected cases, infection rate, recovery rate.

Optimization

The mobility can be controlled by tuning the capacities of different destinations using capacity control inputs on different categories. The GTL-Healthmob shows the optimal control policies for epidemic mitigation in Grenoble metropolis. The optimization problem which has been implemented in the interface is presented in [START_REF] Niazi | Optimal Control of Urban Human Mobility for Epidemic Mitigation[END_REF] and an algorithm for solving it has been proposed in [START_REF] Pratap | Scalable socio-economic optimal mobility policies for epidemic mitigation under health constraints: application to Grenoble metropolis[END_REF]. This optimization problem aims at maximizing the economy of a city by allowing as much people to visit different destination as possible while keeping the total infection below a bound, which reflects the hospitals bed capacity.

For the implementation and adaption of the mobility model to the Grenoble metropolis, we consider to have 55 communes, and have grouped the destinations into 5 categories and 15 subcategories. For the platform, several types of data are needed which are as follows:

1. Map of the communes. (i) Time schedules (ii) Demand Gating Profiles for and Supply Gating Functions.

5. Average number of people visiting from an origin to a given destination.

The above-mentioned set of data is stored in Storage module which is imported by the Web Interface module for visualization of the nodes and the OD network and also for the simulation. A detailed description and methods to obtain these above mentioned datasets and construct the OD network using publicly available information have been presented in [START_REF] Pratap | Where, when and how people move in largescale urban networks: the Grenoble saga[END_REF]. Some of the most prominent sources are Openstreetmap, INSEE, Academy Grenoble, data.gouv.fr, EMD travel surveys, etc. These methods have also been presented in [START_REF] Pratap | From a toy example to the Large-Scale network of Grenoble," 1st Workshop "Epidemics: modeling, identification, control[END_REF].

User Interface

In this section we describe the interfaces and functionality of the web platform available to the general public. The welcome page of the GTL-Healthmob is shown in Figure 3. In the right panel, the map of Grenoble metropolis is shown with different buttons to visualise the network:

• Borders : It makes the border appear or disappear from the map.

• Nodes : A list of categories and subcategories appears. On selecting a category or a subcategory, the respective destination nodes are shown. When hovered over the nodes, their information pops up as in Figure 4. The 'aggregate' button aggregates the nodes in some subcategories which appear in italics. There are around 4000 destination nodes grouped in 5 categories and 15 subcategories but after aggregation, we have 374 destination nodes. For simulations, these aggregated nodes are used. • Population distribution : Shows the population distribution using a heatmap when origin is selected and distribution of capacities when destination subcategory is selected. See Figure 5a.

• O/D Graph : Shows the OD connections between different origins and the aggregated destination nodes. Figure 5b shows the connections for schools.

In the left panel of the welcome page, the user has different options for simulations. The options are : 1. Unrestricted mobility, 2. Manual restricted mobility 3. Optimal restricted mobility. On running the simulations, the user can visualise the results with the help of plots and heatmap animations which will be discussed in upcoming sections.

Unrestricted mobility

In this mode, user can visualise, how an epidemic, with different levels of infection rate and recovery rate, evolves in the Grenoble metropolis. The user has three epidemic scenarios to choose from. For all these scenarios, mobility is unrestricted but they differ in terms of infection level in various destination categories. The Parameters in the top left as can be seen in Figure 3 are:

• Simulation duration: The number of weeks the simulation runs. In this mode, it is fixed to be 15 weeks.

• Initial infected cases : The number of people initially infected.

• Select scenario : Choose the epidemic scenario. These scenarios mimic the behavior of different Covid-19 variants.

• Recovery period : This represents the number of days required for a person to recover if he is infected. It depends on the epidemic scenario.

• Infection level : This represents the rate at which an infected person might infect others. It ranges from low to high in different categories depending on the epidemic scenario.

Result visualization

When the run button is clicked, three plot visualization and animation of three heatmaps are shown.

The simulation uses the discrete-time SIR-mobility model presented in [START_REF] Niazi | Optimal Control of Urban Human Mobility for Epidemic Mitigation[END_REF]. In the left panel, user can see the following plots as can be seen in Figure 6.

• Epidemic evolution : Weekly average of proportion of infected people.

• mobility : Mobility patterns of people in each category for a week.

• ICU bed occupancy : Occupancy of the ICU beds in the hospitals of Grenoble metropolis.

Apart from the plots, user can choose to see the heatmap animations for mobility and infection evolution. For infections, heatmap is shown in absolute numbers and relative as well.

Manual restricted mobility:

In this mode, the user can choose to simulate mobility and epidemic scenario with parameters and controls of his choice. The parameters as can be seen in Figure 7 are:

• Simulation duration : The duration in weeks for which to run the simulation.

• Recovery period : Number of days, a person takes to recover after getting infected.

• Infection probability : The probability by which an infected person would infect a susceptible person.

• Initial infected cases : Initial number of infected persons.

• Allowed capacity : With its value in the interval [0, 1] for each category, user can impose restriction on mobility by restricting the capacity of destinations in different categories.

• Closing hours : With their value in the range [6, 24] for each destination category, user can simulate partial or full closure of the destinations by forcing them to close before the nominal closing hours.

• Contact rate : Average number of contacts people make in different categories.

Result visualization:

The results are visualised using the same plots and heatmaps described in Section 5.1.1. The simulation uses the discrete-time SIR-mobility model presented in [START_REF] Pratap | Where, when and how people move in largescale urban networks: the Grenoble saga[END_REF]. Plots and heatmap can be seen in Figure 7 for a simulation done for 15 weeks, when schools are running at half their capacities.

Optimal restricted mobility

This mode takes into account the implementation of the optimization problem presented in [START_REF] Pratap | Where, when and how people move in largescale urban networks: the Grenoble saga[END_REF] to the Grenoble metropolis. Here, the user can visualise optimal mobility restriction policies by restricting the capacities of destinations depending on epidemic status during the simulation period. The parameters are simulation duration, initial infected, recovery period, infection level and infection threshold. Infection threshold is a bound on maximum infection during the epidemic. This threshold reflects the hospital bed capacity. The rest of the parameters have been explained in previous sections. In this mode, the user has two options in terms of how the capacity control policies are applied :

• Category-free : Here, the same capacity control is applied to all the destinations irrespective of the category while maximising the economy and keeping the infection below a threshold. It is given by the algorithm developed in [START_REF] Pratap | Scalable socio-economic optimal mobility policies for epidemic mitigation under health constraints: application to Grenoble metropolis[END_REF].

• Category-dependent : Here, different capacity controls are applied to different the destinations in destination categories.

Result Visualization

The results are visualised through the following plots as shown in Figure 8:

• Capacity control : It shows the optimal capacity control inputs to be applied during different weeks for simulation to keep the infection in check and to protect the economy.

• Infection evolution : It shows the comparison between infections in uncontrolled case and the case when the optimal capacity control is applied.

• Mobility : It shows weekly mobility patterns. It is shown by an animation as the mobility pattern is subject to change with respect change in the capacity control. • Economic activity : It is depicted by a bar, which shows, in terms of percentage, how much economic activity under optimal mobility restriction was done with respect to that in case of unrestricted mobility.

• The animation for heatmap for mobility, absolute infection, and relative infection is shown as in the previous cases.

Conclusion and perspective

In this paper, we presented the main framework of the GTL-Healthmob simulation platform and functionality of its different components. The main objective of this platform is to facilitate other users to simulate different scenarios of epidemic spread or mobility restrictions to visualise the mobility patterns and its effect on the epidemic spread in Grenoble metropolis. Using it, an user can also get an overview of the current situation during an epidemic. Moreover, the platform uses different visualization techniques to make even a layman user to understand the trends of the mobility patterns and epidemic spread. It can also serve as a testbench for validating upcoming results.

Since, this platform is still under development, more features and improvements will appear in coming days. Finally, one of the perspective is to extend the mobility model to analyze the spatial distribution of the energy consumption of vehicles in a city.

Figure 1 :

 1 Figure 1: Modules and submodules of the GTL-Healthmob

Figure 2 :

 2 Figure 2: Depiction of mobility between different locations and epidemic spread inside each of them

2 .

 2 Origin and destination nodes.

  (i) Location (ii) Population of origins and capacities of destinations. 3. Origin-Destination network. 4. For destination in each subcategory for each day of the week.

Figure 3 :

 3 Figure 3: Welcome page of the GTL-Healthmob platform. In the left panel, there are different modes and parameters which an user can set to simulate different scenarios. In the right panel, users can visualise mobility-network related information and the heatmap animations for infection and mobility after the simulation.

Figure 4 :

 4 Figure 4: Node menu in GTL-Healthmob platform.

Figure 5 :

 5 Figure 5: (a) Population distribution. (b) Mobility network for schools. Origins are denoted in blue color and the destination in white color. The color of the connecting edge ranges from magenta to yellow depending on the amount of average daily flow along them.

Figure 6 :

 6 Figure 6: Unrestricted mobility: plots showing evolution of epidemic, mobility patterns and ICU occupancy and the right panel shows the heatmap for infection in absolute numbers at different locations on the day of infection peak during morning. Here, the selected scenario is similar to that of Covid19.

Figure 7 :

 7 Figure 7: Interface with manual restricted mobility mode. Here, we can see different parameters and control which an user can select for a simulation. It shows the results of a simulation when Schools are running at half their capacities.

Figure 8 :

 8 Figure8: Left panel of the Web Interface, when the mode "optimal restricted mobility" is selected. It shows the results in terms of plots and bar for category-free capacity control policy in case when the simulation duration is 15 weeks and infection threshold is set to 15% of the total population.

Acknowledgements

This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement 694209) and INRIA project Healthy-Mobility.

We also thank the help and participation of Federica Garin, M.U.B. Niazi, Alain Kibangou, P.A. Bliman, Vadim Bertrand in the theoretical development and practical implementation of this platform.


Kerem Akartunali
email: k.akartunali@ms.unimelb.edu.au

Andrew J Miller
email: andrew.miller@math.u-bordeaux1.fr

Kerem Akartunalı

A Computational Analysis of Lower Bounds for Big Bucket Production Planning Problems

Keywords: Production Planning, Integer Programming, Strong Formulations, Lagrangian Relaxation

In this paper, we analyze a variety of approaches to obtain lower bounds for multilevel production planning problems with big bucket capacities, i.e., problems in which multiple items compete for the same resources. We give an extensive survey of both known and new methods, and also establish relationships between some of these methods that, to our knowledge, have not been presented before. As will be highlighted, understanding the substructures of difficult problems provide valuable insights on why these problems are hard to solve. We conclude with computational results from widely used test sets and discussion of future research.

Introduction

Production planning problems have drawn considerable interest from both researchers and practitioners since the seminal paper of Wagner and Whitin [START_REF] Wagner | Dynamic version of the economic lot size model[END_REF]. These problems search for the production plan with the minimum total cost, which consists of fixed charges such as setup costs and linear charges such as inventory holding costs, that satisfies demand and follows restrictions of the production environment such as those imposed by capacities. The focus of this paper is on multi-level, multi-item production planning problems with "big bucket" capacities, i.e., each resource is shared by multiple items. These problems often have complicated BOM (Bills of Materials) structures, where the BOM details which items are required to produce another item. Due to these prerequisites, the BOM often has multiple levels, where the last level can be thought of as finished products, the next-to-last level can be thought of as items required to make finished products, and so forth.

The MRP (Materials Requirement Planning) approach and its successors MRP-II (Manufacturing Resource Planning) and ERP (Enterprise Resource Planning) have been widely used in the manufacturing industry to generate production plans. While they do provide accurate accounting for BOM structures, these approaches fail to account accurately for capacity restrictions and hence they do not consistently achieve feasible (let alone high quality) production plans. Realistic multi-level multi-item production planning problems are complicated and computationally challenging to solve, and therefore the development of computationally effective methods to tackle these problems is necessary.

On the theoretical side, the capacitated version of even the single-item lotsizing problem is N P-hard (see Florian et al. [START_REF] Florian | Deterministic production planning: Algorithms and complexity[END_REF] and Bitran and Yanasse [START_REF] Bitran | Computational complexity of the capacitated lot size problem[END_REF]). Because of problem complexity, dynamic programming algorithms have been proposed only for some special cases of the problem, see e.g. Zangwill [START_REF] Zangwill | A backlogging model and a multi-echelon model of a dynamic economic lot size production system-a network approach[END_REF], Florian and Klein [START_REF] Florian | Deterministic production planning with concave costs and capacity constraints[END_REF], Federgruen and Tzur [START_REF] Federgruen | A simple forward algorithm to solve general dynamic lot sizing models with n periods in O(nlogn) or O(n) time[END_REF], Aggarwal and Park [START_REF] Aggarwal | Improved algorithms for economic lot size problems[END_REF].

Heuristic algorithms have been employed for production planning problems by many researchers with the hope of obtaining good solutions in acceptable computational times. For a general review of earlier lot-sizing heuristics, refer to Maes and Van Wassenhove [START_REF] Maes | Multi item single level capacitated dynamic lotsizing heuristics: A computational comparison (part i: Static case; part ii: Rolling horizon)[END_REF]. Heuristic frameworks in general use some decomposition ideas, such as Lagrangian-based decomposition (e.g. Trigeiro et al. [START_REF] Trigeiro | Capacitated lot sizing with setup times[END_REF], Tempelmeier and Derstroff [START_REF] Tempelmeier | A lagrangean-based heuristic for dynamic multilevel multiitem constrained lotsizing with setup times[END_REF]), forward scheme and relax-and-fix (e.g. Afentakis and Gavish [START_REF] Afentakis | Optimal lot-sizing algorithms for complex product structures[END_REF], Belvaux and Wolsey [START_REF] Belvaux | bc-prod: A specialized branch-and-cut system for lot-sizing problems[END_REF], Stadtler [START_REF] Stadtler | Multilevel lot sizing with setup times and multiple constrained resources: Internally rolling schedules with lot-sizing windows[END_REF], Federgruen et al. [START_REF] Federgruen | Progressive interval heuristics for multi-item capacitated lot sizing problem[END_REF]) and coefficient modification (e.g. Katok et al. [START_REF] Katok | Lot sizing in general assembly systems with setup costs, setup times, and multiple constrained resources[END_REF], Van Vyve and Pochet [START_REF] Van Vyve | A general heuristic for production planning problems[END_REF]). The main disadvantage of the heuristic algorithms is the lack of guarantee of solution quality, and they also do not always provide useful insights about basic problem structures.

Mathematical programming results on production planning problems have usually focused on special cases such as single-item problems, and they have been limited for problems with big bucket capacities. We will briefly discuss these techniques in two subgroups: 1) Valid inequalities that are added into the original formulation using separation algorithms, and 2) Extended reformulations that solve the problem in a different variable space.

The first polyhedral study that defines problem-specific valid inequalities for production planning problems is the study of Barany et al. [START_REF] Barany | Strong formulations for multi-item capacitated lot-sizing[END_REF]. The authors propose the family of ( , S) inequalities for the single-item uncapacitated lotsizing problem, which describe the polytope of these problems. Some special cases of lot-sizing problems are investigated in Küçükyavuz and Pochet [START_REF] Küçükyavuz | Uncapacitated Lot-Sizing with Backlogging: The Convex Hull[END_REF] (uncapacitated problem with backlogging), Pochet and Wolsey [START_REF] Pochet | Polyhedra for lot-sizing with Wagner-Whitin costs[END_REF] (constant capacities), Loparic et al. [START_REF] Loparic | The uncapacitated lot-sizing problem with sales and safety stocks[END_REF] (uncapacitated problem with sales and safety stocks), and Constantino [START_REF] Constantino | A cutting plane approach to capacitated lot-sizing with startup costs[END_REF] (uncapacitated problem with start-up costs). Chan et al. [START_REF] Chan | Effective zero-inventory-ordering policies for the single-warehouse multi-retailer problem with piecewise linear cost structures[END_REF] study a warehouse problem that has a similar structure to a multi-item production planning problem having piecewise-linear costs associated with capacities. Atamtürk and Muñoz [START_REF] Atamtürk | A study of the lot-sizing polytope[END_REF] provide a recent polyhedral study that investigates the bottleneck cover structure in capacitated singleitem problems. Pochet and Wolsey [START_REF] Pochet | Solving multi-item lot-sizing problems using strong cutting planes[END_REF] study multi-item problems using valid inequalities, extending some single-item results to the multi-level case. On the other hand, Miller et al. [START_REF] Miller | Solving the multi-item capacitated lot-sizing problem with setup times by branch-and-cut[END_REF][START_REF] Miller | On the polyhedral structure of a multi-item production planning model with setup times[END_REF] provide rare results on multi-item problems with big-bucket capacities, where the authors study single-period relaxations and propose valid inequalities. In a recent study, Levi et al. [START_REF] Levi | Approximation algorithms for the capacitated multi-item lot-sizing problem via flow-cover inequalities[END_REF] study a version of the capacitated multi-item problem and they propose an approximation algorithm based on generating flow cover inequalities and randomized rounding.

Extended reformulations provide interesting results for production planning problems. A compact extended reformulation is the facility location reformulation of Krarup and Bilde [START_REF] Krarup | Plant location, set covering and economic lotsizes: an O(mn) algorithm for structured problems[END_REF], which defines the convex hull of the uncapacitated single-item problem when projected to original variable space. Eppen and Martin [START_REF] Eppen | Solving multi-item capacitated lot-sizing problems using variable redefinition[END_REF] study the shortest path reformulation, which is of smaller size compared to facility location reformulation. Rardin and Wolsey [START_REF] Rardin | Valid inequalities and projecting the multicommodity extended formulation for uncapacitated fixed charge network flow problems[END_REF] investigate the multi-commodity reformulation for fixed-charge network problems. Belvaux and Wolsey [START_REF] Belvaux | Modelling practical lot-sizing problems as mixedinteger programs[END_REF] and Wolsey [START_REF] Wolsey | Solving multi-item lot-sizing problems with an MIP solver using classification and reformulation[END_REF] are recent studies about reformulations and modeling issues. Anily et al. [START_REF] Anily | Multi-item lot-sizing with joint set-up costs[END_REF] provide tight reformulations for some special cases of the multi-item problem with joint setups.

In spite of this research, big bucket production planning problems remain hard to solve. Part of the reason for this is that most previous research focuses on developing and using results for single-item models, which are not sufficient to capture the fundamental sources of complexity of big bucket problems. The primary goals of this paper are to evaluate the strength of the relaxations defined by different mathematical programming techniques and to investigate why big bucket production planning problems are hard to solve in practice. More specifically, we are not primarily interested in extending single-item results to general production planning problems, but we want to discover relationships between different methods for generating lower bounds and the fundamental substructures that often make these methods insufficient to solve these problems well. We will consider all known methods for generating lower bounds of which we are aware, and we will investigate previously untried methods as well.

We can formulate the basic model that we study as follows:

min NT t=1 NI i=1 f i t y i t + NT t=1 NI i=1 h i t s i t (1) s.t. x i t + s i t-1 -s i t = d i t t ∈ [1, NT ], i ∈ endp (2)
x i t + s i t-1 -s i t = j∈δ(i) r ij x j t t ∈ [1, NT ], i ∈ [1, NI]\endp (3) NI i=1 (a i k x i t + ST i k y i t ) ≤ C k t t ∈ [1, NT ], k ∈ [1, NK] (4)
x i t ≤ M i t y i t t ∈ [1, NT ], i ∈ [1, NI] (5)
y ∈ {0, 1} N T xN I (6) x ≥ 0 (7) s ≥ 0 (8)
In this formulation, NT , NI and NK are the number of periods, items, and machines, respectively. The set endp includes all of the end-items, i.e. items with external demand; the other items are assumed to have only internal demand. (We lose no generality with this assumption, since any item that has both internal and external demand can be considered to be two distinct items, where the data not related to the demand and BOM is identical for these two items.) The variables x i t , y i t , and s i t represent production, setup, and inventory amounts for item i in period t, respectively. The setup and inventory cost coefficients are represented by f i t and h i t for each period t and item i. The parameter δ(i) represents the set of immediate successors of item i, and the parameter r ij represents the number of items required of i to produce one unit of j, not only for immediate but for all dependencies between i and j. The parameter d i t is the demand for end-product i in period t, and d i t,t is the total demand for i from period t to t , i.e., d i t,t = t t=t dt.

The parameter a i k represents the time necessary to produce one unit of i on machine k, and ST i k is the setup time for item i on machine k, which has a capacity of C k t in period t. Note that each item is processed by a preassigned machine, and we assume that each item is assigned only to one machine. (In many situations in both practice and the literature this assumption holds; when it does not, the formulation can be modified by including an additional index k on the x and y, updating the flow balance constraints, etc. In general, the results we discuss apply to these more general models as well-as has been previously observed by Miller [START_REF] Miller | Polyhedral Approaches to Capacitated Lot-Sizing Problems[END_REF], Stadtler [START_REF] Stadtler | Multilevel lot sizing with setup times and multiple constrained resources: Internally rolling schedules with lot-sizing windows[END_REF], and others).

The constraints (2) and (3) ensure production balance and demand satisfaction for end-items and intermediate items respectively, (4) are the big bucket capacity constraints, [START_REF] Florian | Deterministic production planning with concave costs and capacity constraints[END_REF] ensure that the setup variable is set to be 1 if there is positive production, and finally ( 6), [START_REF] Aggarwal | Improved algorithms for economic lot size problems[END_REF], and (8) provide the integrality and nonnegativity requirements. Note that we define M i t as follows:

M i t = min(d i t,N T , C k t -ST i k a i k ) i ∈ endp M i t = min( j∈endp r ij d j t,N T , C k t -ST i k a i k ) i ∈ [1, NI]\endp
We next define an echelon reformulation of the problem, see e.g. Pochet and Wolsey [START_REF] Pochet | Production Planning by Mixed Integer Programming[END_REF]. Our motivation for defining this reformulation is that it clearly exhibits the single-item structure that is present for each item, and it therefore enables us to apply results for single-item models to the multi-level model. We first define echelon demand parameters D i t and echelon stock variables E i t :

D i t = d i t + j∈δ(i) r ij D j t t ∈ [1, NT ], i ∈ [1, NI] (9)
E i t = s i t + j∈δ(i) r ij E j t t ∈ [1, NT ], i ∈ [1, NI] (10)
Substituting ( 10) into ( 2) and (3) for s i t , and using the definition (9), we obtain an equation that can replace (2) and (3) in the original formulation:

x i t + E i t-1 -E i t = D i t t ∈ [1, NT ], i ∈ [1, NI] (11)
To satisfy [START_REF] Maes | Multi item single level capacitated dynamic lotsizing heuristics: A computational comparison (part i: Static case; part ii: Rolling horizon)[END_REF], we add the following constraints:

E i t ≥ j∈δ(i) r ij E j t t ∈ [1, NT ], i ∈ [1, NI] (12)
E ≥ 0 (13)
Finally, to eliminate the original inventory variable s, we define echelon inventory holding costs

H j t = h j t -NI i=1 r ij h i t ,
H i t E i t ( 14
)
We can therefore define the feasible region of the production planning problem as X = {(x, y, E)|( 4)-( 7), (11)-( 13)}, which will be referred in the remainder of the paper as the "basic formulation". The production planning problem can be defined as min{( 14)|(x, y, E) ∈ X}. We could easily include overtime (i.e., extra capacity that can be bought with an additional cost) or backlogging (i.e., satisfying demand later than requested by the customer with a cost for customer dissatisfaction) variables to generalize this basic model, and some of the test problems we consider in Section 4 incorporate them.

For simplicity, we will sometimes use conv(a) to denote conv((x, y, E)|(a)), where (a) is a set of constraints. For example, {(x, y)|( 7) ∩conv(( 6))} represents {(x, y)|(7)} ∩ conv({(x, y)|(6)}) in our notation, or, equivalently {(x, y)| (7), 0 ≤ y ≤ 1}.

In Section 2, we provide a comprehensive survey of lower bounding methods presented in previous research, and we discuss previously untested methods as well. Section 3 is devoted to theoretical comparisons of different techniques, which can provide structural insight into multi-level big bucket problems. In Section 4, we present extensive computational comparisons obtained using widely known data sets. We conclude with future directions in Section 5.

Valid Inequalities, Reformulations, and Relaxations

In this section we discuss different approaches to obtain lower bounds. These methods vary from defining valid inequalities and reformulations to the use of Lagrangian relaxation.

The first technique we consider is the use of ( , S) inequalities of Barany et al. [START_REF] Barany | Strong formulations for multi-item capacitated lot-sizing[END_REF] defined for single-item problems, and generalized by Pochet and Wolsey [START_REF] Pochet | Solving multi-item lot-sizing problems using strong cutting planes[END_REF] to multi-level problems using the echelon reformulation. These can be defined as follows:

t∈S x i t ≤ t∈S D i t, y i t + E i ∈ [1, NT ], i ∈ [1, NI], S ⊆ [1, ] (15)
Since these inequalities are valid for the single-item submodels defined by each item, they are valid for the multi-item problem as well. Although there is an exponential number of these inequalities, a simple polynomial separation algorithm exists as shown in Barany et al. [START_REF] Barany | Uncapacitated lot sizing: The convex hull of solutions[END_REF], see Algorithm 1. As will be discussed later, there exist stronger formulations than that provided by using the ( , S) inequalities alone, but ( , S) inequalities have good practical use, especially when considering large problems.

The feasible region associated with this formulation can be defined as X LS = {(x, y, E)|( 4)-( 7), (11)-( 13), (15)}, and the problem can be defined as

Z LS = min{(14)|(x, y, E) ∈ X LS }.
The next technique is the facility location reformulation, originally defined by Krarup and Bilde [START_REF] Krarup | Plant location, set covering and economic lotsizes: an O(mn) algorithm for structured problems[END_REF] for the single-item problem. This reformulation divides

u i t,t ≤ D i t y i t t ∈ [1, NT ], t ∈ [t, NT ], i ∈ [1, NI] (16)
t t=1 u i t,t = D i t t ∈ [1, NT ], i ∈ [1, NI] (17)
x i t = NT t=t u i t ,t t ∈ [1, NT ], i ∈ [1, NI] (18)
u ≥ 0 [START_REF] Pochet | Polyhedra for lot-sizing with Wagner-Whitin costs[END_REF] This formulation adds O(NT 2 NI) variables and O(NT 2 NI) constraints to the problem.

One advantage of using the new variables u i t,t is that we can rewrite the capacity constraint (4) as follows:

NI i=1 (a i k ( NT t =t u i t,t ) + ST i k y i t ) ≤ C k t t ∈ [1, NT ], k ∈ [1, NK] (20)
This, along with constraints ( 16), can considerably help a state-of-the-art MIP solver generate knapsack cover cuts. Specifically, note that by adding NI i=1 a i k D i t,N T y i t on both sides and after rearranging the terms, (20) can be rewritten as

NI i=1 (a i k D i t,N T + ST i k )y i t ≤ C k t + NI i=1 NT t =t a i k (D i t y i t -u i t,t ) (21)
For each fixed pair of (t, k), and for any subsets I ⊆ {1, ..., NI} and T ⊆ {t, ..., NT }, we may generate cover cuts for each of the following continuous 0-1 knapsack constraints:

i∈I (a i k ( t ∈T D i t ) + ST i k )y i t ≤ C k t +   i∈I t ∈T a i k (D i t y i t -u i t,t )   (22)
Note that because of ( 16), the expression in the parenthesis on the righthand side of ( 21) or ( 22) can be considered as a single nonnegative continuous variable. Binary knapsack constraints with a single nonnegative continuous variable were studied by Marchand and Wolsey [START_REF] Marchand | The 0-1 knapsack problem with a single continuous variable[END_REF][START_REF] Marchand | Aggregation and mixed integer rounding to solve mips[END_REF] (see also Richard et al. [START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Basic theory and algorithms[END_REF][START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Superlinear lifting[END_REF]). Commercial solvers use the kinds of results they present to efficiently find subsets I and T and generate cover cuts that will approximate conv(X 6), ( 16), ( 19), ( 20)} is the feasible region of the intersection of these continuous 0-1 knapsack problems for a fixed (t, k) pair. Note that we can also define it as 6), ( 16), ( 19), ( 20), ( 18), ( 11)}, just for the convenience of having it in higher dimension. Related to X(t,k) KN , we will define

(t,k) KN ), where X (t,k) KN = {(y, u)|(
X (t,k) KN = proj y,u X(t,k) KN with X(t,k) KN = {(x, y, E, u)|(
X(t,k,{t(i)}) KN , for which we first choose a t(i) ∈ [t, NT ] for all i ∈ [1, NI],
for a given t. Then, we define

u i t,t 1 ≤ D i t 1 y i t t 1 ∈ [t, NT ], i ∈ [1, NI] (23)
u i t 1 ,t 2 ≤ D i t 2 y i t 1 t 1 ∈ [t + 1, t(i)], t 2 ∈ [t 1 , t(i)], ( 24
) i ∈ [1, NI] x i t = NT t 1 =t u i t,t 1 i ∈ [1, NI] (25)
E i t-1 = t-1 t 1 =1 NT t 2 =t u i t 1 ,t 2 i ∈ [1, NI] (26)
x i t + E i t-1 + t(i) t 1 =t+1 t(i) t 2 =t 1 u i t 1 ,t 2 ≥ D i t,t(i) i ∈ [1, NI] (27)
Then,

X(t,k,{t(i)}) KN = {(x, y, E, u)|(6), (19)
, ( 20), ( 23) -( 27)}. Note that we will use this explicit definition for the purposes of proving a key proposition in the next section.

On a separate note, basic continuous cover inequalities can also be generated as MIR inequalities, which are known to be effective for general mixed integer programs (see e.g. Günlük and Pochet [START_REF] Günlük | Mixing mixed integer inequalities[END_REF]). Of course, our approach will increase the problem size and it might easily become so large that it cannot be solved to optimality in an acceptable time. However, using this approach for the purpose of generating lower bounds can yield insights into the structure of our problems. This idea was initially suggested for single-level, single-machine problems by Van Vyve [42]. To the best of our knowledge, this approach has not been tested for multi-level problems before.

The feasible region associated with the facility location reformulation can be defined as X F L = {(x, y, E, u)|( 5) -( 7), ( 11) -( 13), ( 16) -(20)}, and the associated problem as Z F L = min{( 14)|(x, y, E, u) ∈ X F L }. On the other hand, generating all cover cuts approximates NT t=1 NK k=1 conv(X

(t,k) KN ), which is an ap- proximation for conv( NT t=1 NK k=1 X (t,k) KN )
. This leads us to define the polyhedron 5), ( 7), ( 11) -( 13), ( 17), (18

X KN F L = {(x, y, E, u)|(
)} ∩ conv( NT t=1 NK k=1 X (t,k)
KN ) and the associated problem

Z KN F L = min{(14)|(x, y, E, u) ∈ X KN F L }.
Next, we discuss the single-period relaxation of Miller et al. [START_REF] Miller | Solving the multi-item capacitated lot-sizing problem with setup times by branch-and-cut[END_REF][START_REF] Miller | On the polyhedral structure of a multi-item production planning model with setup times[END_REF], called as PI (Preceding Inventory). To describe the single-period formulation, for a given machine k ∈ [1, NK] and a given time period t ∈ [1, NT ], we choose a time period t(i) ≥ t for each i ∈ [1, NI]. Then we define

S i = E i t-1 + t(i) t=t+1 D i tt(i) y i t i ∈ [1, NI] D i = D i tt(i) i ∈ [1, NI]
Then, the single-period formulation can be written as follows:

x i t + S i ≥ D i i ∈ [1, NI] (28)
x i t ≤ M i t y i t i ∈ [1, NI] ( 29
) NI i=1 (a i k x i t + ST i k y i t ) ≤ C k t ( 30
)
x i t , S i ≥ 0 i ∈ [1, NI] (31)
y i t ∈ {0, 1} i ∈ [1, NI] (32)
We can define 28) -(32)} as the feasible region associated with a set of t(i) values, and

X (t,k,{t(i)}) P I = {(x, y, S)|(
X (t,k) P I = {t(i)} X (t,k,{t(i)}) P I
represents the feasible region for a given (t, k) pair. Note the similarity between this feasible region and X (t,k) KN we discussed earlier. Miller et al. [START_REF] Miller | Solving the multi-item capacitated lot-sizing problem with setup times by branch-and-cut[END_REF][START_REF] Miller | On the polyhedral structure of a multi-item production planning model with setup times[END_REF] define valid inequalities (namely cover and reverse cover inequalities) for PI, which are naturally valid for the original problem as well, and these inequalities can be seen as an approximation for conv(X (t,k)

P I ).
Next, we define the shortest path reformulation of Eppen and Martin [START_REF] Eppen | Solving multi-item capacitated lot-sizing problems using variable redefinition[END_REF]. In this formulation, which was originally defined for single-item uncapacitated models, we define new variables z i t,t , which are 1 if production of i in period t satisfies all the demand for i in periods t, ..., t , and 0 otherwise. Note the relationship between the new and original variables:

x i t = NT t =t D i t,t z i t,t t ∈ [1, NT ], i ∈ [1, NI] (33)
For the multi-level capacitated problem, we do not have the same optimality properties that we do for the single-item problem; we therefore let the z variables take fractional values. Also, using the echelon inventory holding costs

H i t , we define total inventory costs c i t,t = D i t,t NT j=t H i j .
Then the formulation is

min NT t=1 NI i=1 f i t y i t + NT t=1 NT t =t NI i=1 c i t,t z i t,t (34)
s.t. 1 = NT t=1 z i 1,t i ∈ [1, N I] (35) t -1 t=1 z i t,t -1 = NT t=t z i t ,t t ∈ [2, N T ], i ∈ [1, N I] ( 36
) NT t =t z i t,t ≤ y i t t ∈ [1, N T ], i ∈ [1, N I] (37)
NI i=1

(ST i k y i t + a i k NT t =t D i t,t z i t,t ) ≤ C k t t ∈ [1, N T ], k ∈ [1, N K] ( 38
) t t=1 NT t=t (D i t, tz i t, t - j∈δ(i) r ij D j t, tz j t, t) ≥ d i 1,t t ∈ [1, N T ], i ∈ [1, N I] (39)
z ≥ 0 (40)

y ∈ {0, 1} N T xN I (41)
The constraints [START_REF] Pochet | Production Planning by Mixed Integer Programming[END_REF] and ( 36) are the flow balance constraints, [START_REF] Marchand | The 0-1 knapsack problem with a single continuous variable[END_REF] provide the relationship between the linear and binary variables, [START_REF] Marchand | Aggregation and mixed integer rounding to solve mips[END_REF] is the capacity constraint, [START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Basic theory and algorithms[END_REF] ensures the relationship between different levels, and finally ( 40) and ( 41) provide the nonnegativity and integrality constraints. Note that for our multi-level problem, we derive the constraint (39) as follows: Using (11) and [START_REF] Belvaux | bc-prod: A specialized branch-and-cut system for lot-sizing problems[END_REF], and the assumption of zero initial inventory, we obtain

t t=1 (x i t -D i t ) ≥ t t=1 j∈δ(i) r ij (x j t -D j t ) (42)
Substituting ( 33) into (42) and rewriting results in [START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Basic theory and algorithms[END_REF]. Note that this formulation adds as many variables as the facility location reformulation, but number of constraints is only O(NT × NI). However, this formulation is not necessarily easier to solve, in part because the new constraints are comparatively dense and the coefficients on the new variables comparatively large.

The feasible region associated with this formulation can be defined as X SP = {(y, z)|( 35) -(41)}, and the problem can be defined as Z SP = min{(34)| (y, z) ∈ X SP }. Part of our motivation for completely substituting the x and E variables out of the formulation is that relaxing the constraints ( 35), [START_REF] Barany | Uncapacitated lot sizing: The convex hull of solutions[END_REF], and [START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Basic theory and algorithms[END_REF] decomposes the problem into NT distinct subproblems, one for each time period (an analogous observation was first made for single-level problems by Jans and Degraeve [START_REF] Jans | Improved lower bounds for the capacitated lot sizing problem with setup times[END_REF]). We will discuss this property in more detail later.

Next, we consider the multi-commodity reformulation proposed by Rardin and Wolsey [START_REF] Rardin | Valid inequalities and projecting the multicommodity extended formulation for uncapacitated fixed charge network flow problems[END_REF]. This approach is originally described for fixed-charge network flow problems. Like the facility location reformulation, it divides production using destination information, but since we have multiple levels, it also includes information about which end-item in the BOM it is produced for. Stock variables are also divided in a similar fashion. Thus, the new variables w i,j t,t indicate production of item i in period t to satisfy the demand of end-item j in period t , t ≥ t, and the new variables v i,j t,t indicate the inventory of item i held over at the end of period t to satisfy demand of end-item j in period t , t > t. The following constraints should be added to the basic formulation to finalize the reformulation:

x i t = NT t=t j∈endp w i,j t ,t t ∈ [1, NT ], i ∈ [1, NI] (43)
w i,j t,t ≤ r ij d j t y i t t ∈ [1, NT ], t ∈ [t, NT ], ( 44
) i ∈ [1, NI], j ∈ endp v i,i t-1,t + w i,i t,t = d i t t ∈ [1, NT ], i ∈ endp (45) v i,i t-1,t + w i,i t,t = v i,i t,t t ∈ [1, NT -1], t ∈ [t + 1, NT ], ( 46
) i ∈ endp v i,q t-1,t + w i,q t,t = j∈δ(i) r ij w j,q t,t t ∈ [1, NT ], i ∈ [1, NI]\endp, (47)
q ∈ endp v i,q t-1,t + w i,q t,t = v i,q t,t + j∈δ(i)

r ij w j,q t,t t ∈ [1, NT -1], t ∈ [t + 1, NT ], ( 48
) i ∈ [1, NI]\endp, q ∈ endp w, v ≥ 0 ( 49
)
This reformulation introduces O(NT 2 NI 2 ) additional variables and O(NT 2 NI 2 ) additional constraints. This is the main disadvantage of this reformula-tion, which can become computationally intractable as the problem size grows. However, it is the tightest compact, i.e., polynomial size, reformulation that we know for the problems in question.

The feasible region associated with this formulation can be defined as X MC = {(x, y, E, w, v)|( 4) -( 7), ( 11) -( 13), ( 43) -(49)}, and the problem can be defined as

Z MC = min{(14)| (x, y, E, w, v) ∈ X MC }.
Next, we discuss three approaches that employ Lagrangian relaxation to obtain structured subproblems and from those lower bounds for the original problem.

The first approach is to relax the capacity constraints (4), and obtain

LR 1 (λ) =min NT t=1 NI i=1 f i t y i t + NT t=1 NI i=1 H i E i t - NT t=1 NK k=1 λ k t C k t -( NI i=1 a i k x i t + ST i k y i t ) (50)
subject to (x, y, E) ∈ X LR1

where X LR1 = {(x, y, E)|( 5) -( 7), ( 11) -(13)}. Thus, the Lagrangian subproblem is a multi-item, multi-level uncapacitated lot-sizing problem. The Lagrangian dual problem is

LD 1 = max λ≥0 LR 1 (λ) (51)
The next Lagrangian relaxation approach relaxes the constraints linking separate levels, i.e. constraints [START_REF] Belvaux | bc-prod: A specialized branch-and-cut system for lot-sizing problems[END_REF], to obtain

LR 2 (µ) =min NT t=1 NI i=1 f i t y i t + NT t=1 NI i=1 H i E i t - NT t=1 NI i=1 µ i t   E i t - j∈δ(i) r ij E j t   (52) subject to (x, y, E) ∈ X LR2
where X LR2 = {(x, y, E)| (4) -( 7), ( 11), (13)}. The Lagrangian subproblem therefore decomposes into NK disjoint multi-item, big bucket single-machine problems, one for each machine. The Lagrangian dual problem becomes

LD 2 = max µ≥0 LR 2 (µ) (53)
Finally, the last Lagrangian approach extends the work of Jans and Degraeve [START_REF] Jans | Improved lower bounds for the capacitated lot sizing problem with setup times[END_REF] for single-level problems, which itself uses the shortest path reformulation of Eppen Martin [START_REF] Eppen | Solving multi-item capacitated lot-sizing problems using variable redefinition[END_REF]. Jans and Degraeve [START_REF] Jans | Improved lower bounds for the capacitated lot sizing problem with setup times[END_REF] simply relaxed the constraints linking time periods, yielding disjoint single-period subproblems. However, the problem in the multi-level case is that the constraints linking levels also involve multiple periods. Therefore, decomposing the problem into disjoint subproblems for each period is not possible, unless all constraints linking levels are also dualized. We dualize the constraints ( 35), ( 36) and [START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Basic theory and algorithms[END_REF] in the shortest path reformulation to obtain

LR 3 (β, γ) =min NT t=1 NI i=1 f i t y i t + NT t=1 NT t =t NI i=1 c i t,t z i t,t - NT i=1 β i 1 1 - NT t=1 z i 1,t - NI i=1 NT t =2 β i t   t -1 t=1 z i t,t -1 - NT t=t z i t ,t   (54)
- NI i=1 NT t =1 γ i t   t t=1 NT t=t (D i t, tz i t, t - j∈δ(i) r ij D j t, tz j t, t) -d i 1,t   subject to (y, z) ∈ X LR3
where X LR3 = {(y, z)| (37), ( 38), ( 40), ( 41)}. The Lagrangian subproblem decomposes into NKxNT disjoint capacitated multi-item, single-machine, singleperiod problems, and the Lagrangian dual problem is

LD 3 = max γ≥0,β LR 3 (β, γ) (55)
In the next section we provide theoretical comparisons for the various approaches we have described.

Exploring Relationships

Let the superscript LP indicate the LP relaxation of a problem, i.e., the binary variables y relaxed to be continuous with the bounds 0 ≤ y ≤ 1.

For example, Z LP LS is the problem Z LS with the integrality requirements for y variables relaxed. Similarly, X LP LS is the polyhedron of the LP relaxation of X LS .

Theorem 1 (Akartunalı and Miller [START_REF] Akartunalı | A heuristic approach for big bucket multilevel production planning problems[END_REF]) Z LP LS = Z LP F L = Z LP SP , i.e., the ( , S) inequalities, the facility location reformulation, and the shortest path reformulation all provide the same lower bound for the original problem.

For the proof of the theorem, please refer to Akartunalı [START_REF] Akartunalı | Computational methods for big bucket production planning problems: Feasible solutions and strong formulations[END_REF]. The proof uses Lagrangian duality and the fact that all these formulations provide equal lower bounds in the single-item case. See Krarup and Bilde [START_REF] Krarup | Plant location, set covering and economic lotsizes: an O(mn) algorithm for structured problems[END_REF], Eppen and Martin [START_REF] Eppen | Solving multi-item capacitated lot-sizing problems using variable redefinition[END_REF], and Barany et al. [START_REF] Barany | Uncapacitated lot sizing: The convex hull of solutions[END_REF] for the convex hull and integrality proofs in the single-item case.

Theorem 2 Z LP MC ≥ Z LP F L , i.e., the multi-commodity reformulation provides a lower bound that is at least as strong as that provided by the facility location reformulation. If the problem consists of a single level, then

Z LP MC = Z LP F L .
Although this result has been known by at least some researchers since the publication of Rardin and Wolsey [START_REF] Rardin | Valid inequalities and projecting the multicommodity extended formulation for uncapacitated fixed charge network flow problems[END_REF], it has never been formally stated and proven, to the best of our knowledge. We therefore provide a proof for the sake of completeness.

Proof. We will prove this by showing that proj x,y,E (X LP MC ) ⊆ proj x,y,E (X LP F L ) for the multi-level case. Let (v * , w * , x * , y * , E * ) ∈ X LP MC . First, observe that we can eliminate v * and rewrite ( 45)-( 48) in terms of w * , as follows:

t=t t=1 w * i,j t,t = r ij d j t t ∈ [1, NT ], i ∈ [1, NI], j ∈ endp (56)
Now, let

u * i tt = j∈endp w * ij tt ( 57
)
Obviously u * ≥ 0 since w * ≥ 0. Since w * satisfies (43), x * i t = NT t =t u * i tt . Similarly, summing (56) over j ∈ endp, we obtain t t=1 u * i tt = j∈endp r ij d j t = D i t , where the second equation follows from the definition of echelon demand [START_REF] Trigeiro | Capacitated lot sizing with setup times[END_REF]. Finally, using [START_REF] Akartunalı | A heuristic approach for big bucket multilevel production planning problems[END_REF] and (57), we obtain

u * i tt = j∈endp w * ij tt ≤ ( j∈endp r ij d j t )y * i t = D i t y * i t . This shows that (u * , x * , y * , E * ) ∈ X LP F L . Hence, proj x,y,E (X LP MC ) ⊆ proj x,y,E (X LP F L ). 2
The second part of the theorem can also be shown using the same technique as in the proof of first theorem, i.e., using Lagrangian duality and the fact that the multi-commodity reformulation and the facility location reformulation provide equivalent lower bounds in the single-item case (see Eppen and Martin [START_REF] Eppen | Solving multi-item capacitated lot-sizing problems using variable redefinition[END_REF] and Barany et al. [START_REF] Barany | Uncapacitated lot sizing: The convex hull of solutions[END_REF]).

This theorem shows us theoretically that the multi-commodity reformulation is stronger than the formulation defined by adding ( , S) inequalities, the facility location reformulation, and the shortest path reformulation. In the next section, we will computationally address the question of "how much stronger" for a variety of test problems.

So far we have made comparisons of different polyhedral approaches. Also interesting are the relationships between the Lagrangian approaches and these reformulations, as we investigate in the following results.

Theorem 3 Z LP MC ≤ LD 1 .
In words, the lower bound obtained by the Lagrangian that relaxes the capacity constraints is at least as strong as the lower bound obtained by multicommodity reformulation.

Proof. By the theorem related to the strength of the Lagrangian dual (see e.g. Theorem 10.3 of Wolsey [START_REF] Wolsey | Integer Programming[END_REF]), 5) -( 7), ( 11) -( 13))}

LD 1 = min{(14)|(x, y, E) ∈ (4) ∩ conv((
On the other hand, 5), ( 7), ( 11) -( 13), ( 43) -( 49)) ∩ conv(( 6))}}

Z LP MC = min{(14)|(x, y, E, w, v) ∈ (4) ∩ {((
Observe that {(x, y, E) ∈ conv(( 5) -( 7), ( 11) -( 13))} ⊆ proj x,y,E {(x, y, E, w, v) ∈ {(( 5), (7), ( 11) -( 13), ( 43) -(49)) ∩ conv(( 6))}}

This follows because conv(( 5) -( 7), ( 11) -( 13)) has integer extreme points because the polyhedron is the convex hull of an integer feasible region. On the other hand, {(( 5), (7), ( 11)-( 13), ( 43)-( 49))∩ conv((6))} does not necessarily have integral extreme points. Therefore,

Z LP MC ≤ LD 1 . 2 Theorem 4 Z LP F L ≤ Z KN F L ≤ LD 2 .
In words, the lower bound obtained by the Lagrangian that relaxes the level linking constraints is at least as strong as the lower bound obtained by the facility location reformulation strengthened to approximate the knapsack convex hulls.

Proof. The first relationship follows from the fact that Z KN F L is obtained by strengthening Z LP F L with additional constraints. For the second relationship, first observe that (using the same theorem as in the previous proof) 4) -( 7), ( 11), ( 13))} Observe also that conv(( 4) -( 7), ( 11), ( 13)) ⊆ proj x,y,E {(x, y, E, u)|( 5), (7), ( 11), ( 13), ( 17), [START_REF] Küçükyavuz | Uncapacitated Lot-Sizing with Backlogging: The Convex Hull[END_REF]

LD 2 = min{(14)|(x, y, E) ∈ (12) ∩ conv((
} ∩ conv( NT t=1 NK k=1 X (t,k) KN )
This concludes that Z KN F L is not as strong as LD 2 . 2

As mentioned before, generating cover cuts from [START_REF] Chan | Effective zero-inventory-ordering policies for the single-warehouse multi-retailer problem with piecewise linear cost structures[END_REF] only approximates the knapsack polyhedron and hence Z KN F L is the best possible bound that can be obtained by adding cover cuts to the LP relaxation of the facility location reformulation.

Theorem 5 Z KN F L = LD 3 .
We will use the following result for the proof of the theorem.

Lemma 6 (Pochet and Wolsey [START_REF] Pochet | Lot-size models with backlogging: Strong reformulations and cutting planes[END_REF]) All optimal solutions of the singleitem uncapacitated problem formulated using the facility location reformulation have the following property: 11) -( 13), (33)) ∩ conv(S 2 )}. Note that S 1 and S 2 are integer feasible regions whereas T 1 and T 2 are both polyhedra. Then, the proof of Theorem 5 follows.

u tt D t ≥ u tt +1 D t +1 t ∈ [1, NT ], t ≥ t
Proof. We will prove this by showing proj x,y,E (T 1 ) = proj x,y,E (T 2 ). First, let (x * , y * , E * , u * ) ∈ T 1 and hence (x * , y * , E * ) ∈ proj x,y,E (T 1 ). Therefore,

∃p j = (x j , y j , E j , u j ) ∈ S 1 , j ∈ [1, J], such that (x * , y * , E * , u * ) = J j=1 λ j p j for some λ ≥ 0, J j=1 λ j = 1. For all j ∈ [1, J], let {z i tN T } j = {u i tNT } j D i NT , where t ∈ [1, NT ] and i ∈ [1, NI].

Then, define recursively {z

i tt } j = {u i tt } j D i NT -NT t=t +1 {z i t t} j , for all t ∈ [1, NT ], t = NT -1, ..., t and i ∈ [1, NI]. Since NT t =t D i tt {z i tt } j = NT
t =t {u i tt } j and u j satisfies (20), z j satisfies [START_REF] Marchand | Aggregation and mixed integer rounding to solve mips[END_REF]. Next, note that

NT t =t {z i tt } j = {u i tt } j D i t ≤ {y i t } j
where the last inequality is essentially [START_REF] Van Vyve | A general heuristic for production planning problems[END_REF]. Finally, using Lemma 6, observe that

{z i tt } j = {u i tt } j D i t - {u * i tt +1 } j D i t +1 ≥ 0
Therefore, pj = (x j , y j , E j , z j ) ∈ S 2 , and using the same λ as before, (x * , y * , E * , z * ) = J j=1 λ j pj ∈ T 2 . Hence, (x * , y * , E * ) ∈ proj x,y,E (T 2 ). We conclude therefore that proj x,y,E (T 1 ) ⊆ proj x,y,E (T 2 ). Now, let (x * , y * , E * , z * ) ∈ T 2 and hence (x * , y * , E * ) ∈ proj x,y,E (T 2 ). Therefore,

∃q k = (x k , y k , E k , z k ) ∈ S 2 , k ∈ [1, K], such that (x * , y * , E * , z * ) = K k=1 µ k q k for some µ ≥ 0, K k=1 µ k = 1. For all k ∈ [1, K], let {u i tt } k = D i t NT t=t {z i t t} k ,
{u i tt } k = D i t NT t=t {z i t t} k ≤ D i t NT t=t {z i t t} k ≤ D i t {y i t } k
where the last inequality follows from [START_REF] Marchand | The 0-1 knapsack problem with a single continuous variable[END_REF].

Therefore, qk = (x k , y k , E k , u k ) ∈ S 1 , and using the same µ as before, (x * , y * , E * , u * ) = K k=1 µ k qk ∈ T 1 . Hence, (x * , y * , E * ) ∈ proj x,y,E (T 1 ). Therefore, proj x,y,E (T 2 ) ⊆ proj x,y,E (T 1 ). This concludes the proof. 2

Corollary 7 LD 3 ≤ LD 2 .
The proof for this corollary follows immediately from the Theorems 4 and 5. This result is our main motivation for skipping LD 3 in the computational tests discussed in the next section.

Proposition 8 For any given (t, k) pair and set of {t(i)} values,

proj x,y,E (conv(X (t,k,{t(i)}) P I )) = proj x,y,E (conv( X(t,k,{t(i)}) KN ))
This result, combined with Corollary 7, is our main motivation for omitting computationally testing the cover and reverse cover inequalities from Miller et al. [START_REF] Miller | Solving the multi-item capacitated lot-sizing problem with setup times by branch-and-cut[END_REF][START_REF] Miller | On the polyhedral structure of a multi-item production planning model with setup times[END_REF] in the next section.

Proof. We show first proj x,y,E (conv( X(t,k,{t(i)})

KN

)) ⊆ proj x,y,E (conv(X (t,k,{t(i)}) P I

)) for a given (t, k) pair and set of {t(i)} values. Let (x * , y * , E * , u * ) ∈ conv( X(t,k,{t(i)})

KN

). Then, we define

S * i = E * i t-1 + t(i) t=t+1 D i tt(i) y * i t.
It is easy to observe that (x * , y * , S * ) ∈ conv(X (t,k,{t(i)}) P I

).

Next, we prove that proj x,y,E (conv(X (t,k,{t(i)}) P I

)) ⊆ proj x,y,E (conv( X(t,k,{t(i)})

KN

)) for any given (t, k) pair and set of {t(i)} values. First, let (x * , y * , S * ) ∈ conv(X (t,k,{t(i)}) P I

). We define first

u * i t 1 ,t 2 = D i t 2 y * i t 1 for all t 1 ∈ [t + 1, t(i)] and t 2 ∈ [t 1 , t(i)]. Then, we define E * i t-1 = (S * i - t(i) t=t+1 D i tt(i) y * i t) + . Fi- nally, define u * i t,t = (min{D i t y * i t , x * i t -t -1 t=t u * i t, t})
+ for all t ∈ [t, t(i)], where they are calculated in the increasing order of t . Then, we can observe that (x * , y * , E * , u * ) ∈ conv( X(t,k,t(i))

KN

). 2 4 Computational Results

Overview

In order to provide diversified results, we used the following test instances for our computations:

• TDS instances: These test problems originate from Tempelmeier and Derstroff [START_REF] Tempelmeier | A lagrangean-based heuristic for dynamic multilevel multiitem constrained lotsizing with setup times[END_REF] and Stadtler [START_REF] Stadtler | Multilevel lot sizing with setup times and multiple constrained resources: Internally rolling schedules with lot-sizing windows[END_REF]. These include overtime variables in addition to the formulation in Section 2. Sets A+ and B+ involve problems with 10 items and 24 periods, and sets C and D involve problems with 40 items and 16 periods. Sets B+ and D include setup times. We chose the hardest instances from each data set for our computations, i.e., for each data set, we picked 10 assembly and 10 general instances with the highest duality gaps according to results of Stadtler [START_REF] Stadtler | Multilevel lot sizing with setup times and multiple constrained resources: Internally rolling schedules with lot-sizing windows[END_REF]. • LOTSIZELIB instances: These are the multi-level instances of LOT-SIZELIB [START_REF] Lotsizelib | Lot-sizing problems: A library of models and matrices[END_REF]. These include big bucket capacities, and backlogging is also allowed. The problems vary between 40 item, single end-item problems and 15 item, 3 end-item problems. All problems have 12 periods. • Multi-LSB instances: We have generated 4 sets of test problems based on the problem family described in Simpson and Erenguc [START_REF] Simpson | Modeling multiple stage manufacturing systems with generalized costs and capacity issues[END_REF], each set having 30 instances with low, medium and high variability of demand. From now on, we will call these sets SET1, SET2, SET3, and SET4. These instances are different from the previous sets in that they take component commonality into consideration and hence consider joint setup variables for each family, so setup times are defined for each family. While keeping the original BOM structures and holding costs, we removed the setup costs and added backlogging variables into the problem to obtain problems with a different nature from those of our other test instances. Except for the problems in SET2, which consider a horizon of 24 periods, all the instances have 16 periods. The main difference between SET1, SET2 and SET4 is that they have different resource utilization factors, which are all set over 100%, i.e., it is not possible to setup all families in a period and to produce that period's demand for all items. All problems have 78 items and an assembly BOM structure, and all instances allow backlogging to the last period. For more details about these instances, see Multi-LSB [START_REF]Multi-item lot-sizing problems with backlogging: A library of test instances[END_REF].

Note that average duality gaps after default times (see next section for more detail on "default times") for the test sets of TDS and Multi-LSB are provided in the Table 1 for an overview of problem complexity, where the basic formulation is strengthened with all violated ( , S) inequalities generated at the root node of the Branch&Bound tree using Algorithm 1. The main goal of this section is to computationally test the results we have theoretically proven and to observe how these strength relationships work in practice. This not only provides us with information about how strong the lower bounds actually are but also helps us to understand what prevents us from improving them. All the test instances are run on a PC with an Intel Pentium 4 2.53 GHz processor and 1 GB of RAM. All the formulations are implemented using Xpress Mosel (Xpress-MP 2004C, Mosel version 1.4.1).

In evaluating Lagrangians, we do not exactly solve any of the Lagrangian dual problems, which would require some method (such as a subgradient approach) to choose the optimal Lagranian multipliers. Instead, we first consider a strengthened LP formulation, i.e., the echelon formulation with all violated ( , S) inequalities generated at the root node, and then fix the the Lagrangian multipliers to the values of the optimal dual variables of the constraints to be relaxed in this formulation. We thus evaluate LR 1 (µ * ) and LR 2 (λ * ), respectively, for the optimal dual variables µ * of the capacity constraints and the optimal dual variables λ * of the level-linking constraints, respectively, in order to approximate LD 1 and LD 2 , respectively. These subproblems themselves are MIPs that, in general, are difficult to solve to optimality. Nevertheless, any lower bound on the optimal solution of the Lagrangian subproblem MIP is also a lower bound on the Lagrangian dual (and hence the original problem). Moreover, in every instance, for both LD 1 and LD 2 , the lower bound obtained computationally for the Lagrangian subproblem MIP is at least as strong as the lower bound provided by the original echelon formulation strengthened with ( , S) inequalities.

Similarly, as we discussed before, generating cover cuts on top of the facility location reformulation provides only an approximation of Z KN F L . Hence, the computational comparisons we provide for these relationships are all based on approximations. However, it seems that the approximations are often close. This gives us the chance to compare empirical results in addition to theoretically proven relationships.

Results

The detailed results for TDS instances can be found in the Appendix (as well as in Akartunalı [START_REF] Akartunalı | Computational methods for big bucket production planning problems: Feasible solutions and strong formulations[END_REF]). Note that we obtain the root node solution of the Branch&Bound tree for ( , S) inequalities, all generated through Algorithm 1, and for the multi-commodity reformulation (MC), without the effect of any solver cuts. For the facility location reformulation (FL), all the cover cuts generated by the solver are added at the root node and this strengthened formulation is used as FL lower bound. For comparison purposes, we also use the lower bound obtained by the heuristic in our companion paper (Akartunalı and Miller [44]), where the lower bound is based on the first iteration of a relaxand-fix framework, i.e., a partial LP relaxation of the original problem. For the Lagrangian relaxations that relax the capacity and level-linking constraints, we use the dual optimal values of the constraints from the strong LP relaxation as multipliers, and we set default times of 180 seconds for A+ and B+ instances, and 500 seconds for C and D instances. Note that if the Lagrangian relaxation subproblem is not solved to optimality in this preassigned time, the lower (LB) and upper (UB) bounds of this Lagrangian subproblem provide us the range where the actual lower bound (LB LD ) of the Lagrangian dual lies, i.e., LB ≤ LB LD ≤ UB. Hence, note that we use the lower and upper bounds of Lagrangian subproblem, i.e., LB and UB, in our discussions. Finally, note that due to Theorem 1 we omit the shortest path reformulation in our tests.

We review the results in pairwise comparisons, which are summarized in Table 2. One interesting computational comparison is the relationship we have proven in Theorem 2. As we can see from the detailed results, MC improves the ( , S) bound slightly, in general less than %1. The average improvements from the ( , S) inequalities bound to the MC bound, calculated as (MC bound -, S bound)/( , S bound) for each test instance, are provided in the column "MC vs. , S", and these values are around 0.20%. Considering the enormous size of the MC reformulation, these improvements are simply not worth the computational effort. The Lagrangian relaxation that relaxes the capacity constraints (1st LR) provides in general another slight improvement over the lower bounds of the MC reformulation, as can be seen in the second column of the same table (Column LB under "1st LR vs. MC"), which is calculated in a similar fashion, i.e., (1st LR bound -MC bound)/(MC bound). Note that we also provide averages calculated in the same way using the 1st LR's upper bounds instead of its lower bounds (Column UB under "1st LR vs. MC"). An interesting observation from the problems in set D, where 1st LR problems for all instances are solved to optimality, is that although in general 1st LR improves the MC bound, it is an approximation of LD 1 and it might result in a bound not as strong as the MC bound. However, as these results indicate, these two bounds are in general very close to each other. On the other hand, as the "FL vs. , S" column of Table 2 indicates, the facility location reformulation with cover cuts added (FL) improves in general the ( , S) bound more significantly compared to previous methods. These average percentages are calculated by (FL bound -, S bound)/( , S bound). Similar to our previous comparisons, we also provide the average improvements of the Lagrangian relaxation that relaxes level-linking constraints (2nd LR) over the FL bound in the column "2nd LR vs. FL", calculated by (2nd LR bound -FL bound)/(FL bound). Although one would expect the 2nd LR, the approximation of LD 2 , to improve the FL lower bounds, at first sight this does not seem to be the case for many problem instances, particularly due to negative averages in the LB column of Table 2. However, as can be seen from the UB column of the table, these problems are not close to optimality, particularly the bigger instances of test sets C and D, and the challenge here is that these problems need much more time than the assigned default times (or any reasonable amount of time) for optimality or even for an acceptable gap. For testing whether this is the case here, we experimented with a few of the small A+ and B+ instances that did not achieve the FL bounds earlier and ran them either until the lower bound was at least as strong as the FL bound or to optimality. However, this experiment failed due to memory problems for the instances from sets C and D.

Finally, the last two columns of Table 2 should also be addressed briefly. These columns indicate the duality gaps for the two Lagrangian problems, and as we mentioned before, the 1st LR problem is in general comparatively easier to solve than the 2nd LR problem. We had a total of 11 instances where the 1st LR could achieve the optimal solution in the assigned default times, compared to none for the 2nd LR.

Next, we present results for LOTSIZELIB instances in Table 3, where all values are shown explicitly, including the optimal solutions (OPT) in the last column. MC provides significant improvement over the ( , S) bound for some of these instances, whereas FL provides negligible improvement over MC. The 1st LR is comparatively more efficient on these instances than the 2nd LR. Note that 1st LR and 2nd LR do not necessarily improve MC and FL bounds respectively, similarly to the results for some TDS instances, since these are approximations for LD 1 and LD 2 . Also, note that all 2nd LR problems are at optimality or near, whereas 1st LR did not result in optimality in quite a few instances after the default time of 180 seconds. This indicates that these instances have the bottleneck not in capacity constraints but in the multilevel structure. This seems to be due in part to the fact that there is a single machine, and the capacity in these problems is comparatively loose. The detailed results on Multi-LSB instances can be seen in Akartunalı [START_REF] Akartunalı | Computational methods for big bucket production planning problems: Feasible solutions and strong formulations[END_REF], and the pairwise comparisons are summarized in Table 4, which is organized in the same fashion as Table 2. The default times for the first two sets are bottleneck lies in the capacities rather than the multi-level structure of these problems.

Summary

One of our main goals of this paper was to understand the structure of production planning problems and the underlying difficulties that make these problems very hard. In general, the Lagrangian relaxations we tested are helpful for this. First of all, recall that in general the Lagrangian relaxation that relaxes capacity constraints provides only slight improvement over the ( , S) bound. This bound is an approximation for the uncapacitated problem polyhedron, which indicates that removing capacities makes the problem much easier. This can also be observed by recalling that the final gaps after the default times were quite small for this Lagrangian relaxation in general.

On the other hand, the facility location reformulation with cover cuts and the Lagrangian that relaxes the level-linking constraints improve the lower bounds much more significantly. Recall that the cover cuts approximate the intersection of all knapsack sets included in the problem, and 2nd LR is an approximation for a single-level capacitated problem. Having higher duality gaps compared to the 1st LR, this Lagrangian relaxation problem is in general much harder to solve, indicating that the level-linking constraints are not the bottleneck of these problems. A similar comparison is achieved by Jans and Degraeve [START_REF] Jans | Improved lower bounds for the capacitated lot sizing problem with setup times[END_REF] for single-level problems, where their Lagrangian relaxation relaxing only period-linking constraints is a harder problem than the one that relaxes capacities. Recall that we did not report computational results on LD 3 , due to the result presented in Corollary 7.

Conclusion

In this paper, we have provided an extensive survey of different methodologies for obtaining lower bounds for big bucket production planning problems, and presented both theoretical and computational comparisons of them.

In summary, it seems that the multi-level structure by itself makes some of our problems challenging to solve. However, for most instances, and in particular for the most challenging, the single-level, capacitated substructures are clearly a much greater contributor to problem difficulty. It is this substructure for which the tools currently at our disposal are evidently not sufficient.

These observations indicate that the main bottleneck with these problems lies in the fact that there is no efficient polyhedral approximation of the multiitem, multi-period, single-level, single-machine capacitated problems. It seems that if we could solve these problems well or even adequately, our ability to solve multi-level bug bucket problems would increase dramatically. While initial efforts to find strong formulations for these problems have been made (e.g. see Miller et al. [START_REF] Miller | Solving the multi-item capacitated lot-sizing problem with setup times by branch-and-cut[END_REF]), this is a fundamental area in which it is crucial for the research community to improve the current state of the art. We will attempt to make contributions in this direction in future research.

Instance

Before starting the proof of Theorem 5 , let S 1

 51 KN = {(y, u)|(6),[START_REF] Van Vyve | A general heuristic for production planning problems[END_REF],[START_REF] Pochet | Polyhedra for lot-sizing with Wagner-Whitin costs[END_REF],[START_REF] Loparic | The uncapacitated lot-sizing problem with sales and safety stocks[END_REF]} and S 2 = {(y, z)|(37),[START_REF] Marchand | Aggregation and mixed integer rounding to solve mips[END_REF],[START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Superlinear lifting[END_REF], (41)}. Also let T 1 = {(x, y, E, u)|((11) -(13), (18)) ∩ conv(S 1 )} and T 2 = {(x, y, E, z)|((

  , which indicate the production of item i in period t to satisfy the demand of period t , where t ≥ t. The following constraints should be added into the basic formulation to finalize the reformulation:

	Algorithm 1. ( , S) separation
	Input: LP relaxation solution (x * , y * , E * )
	Output: Violated ( , S) inequalities
	for i=1 to NI
	for = 1 to NT
	Initialize S ← { }
	for t=1 to
	if x * i t > D i t, y * i t
	S ← S ∪ {t}
	if t∈S x * i t > t∈S D i t, y * i t + E * i
	Add the violated ( , S) inequality
	production according to which period it is intended for. This requires first
	defining new variables u i t,t

  NI]. Obviously, u k satisfies (19) since z k satisfies[START_REF] Richard | Lifted inequalities for 0-1 mixed integer programming: Superlinear lifting[END_REF]. Since

	NT t =t {u i tt } k = NT t =t D i tt {z i tt }

where t ∈ [1, NT ], t ∈ [t, NT ], and i ∈ [1, k and z k satisfies (38), u k satisfies (20). Finally, note that

Table 2

 2 Pairwise comparisons of lower bounds and LR gaps for TDS instances

	Test MC vs. 1st LR vs. MC	FL vs. 2nd LR vs. FL	Gaps
	Set	, S	LB	UB	, S	LB	UB	1st LR 2nd LR
	A+	0.29%	0.80% 2.99%	1.81% -0.05% 7.44%	2.09%	6.87%
	B+	0.28%	0.59% 3.06%	1.37% -0.35% 6.23%	2.38%	6.18%
	C	0.14%	0.20% 1.67%	0.86% -0.32% 6.25%	1.44%	6.14%
	D	0.21% -0.06% -0.06% 0.45% -0.43% 19.88%	0%	15.85%

Table 4

 4 Pairwise comparisons of lower bounds and LR gaps for Multi-LSB instances

	Test	MC vs. 1st LR vs. MC	FL vs. 2nd LR vs. FL	Gaps
	Set	, S	LB	UB	, S	LB	UB	1st LR 2nd LR
	SET1	0.02% -0.02% -0.02% 0.85% -0.29% -0.28%	0.00%	0.01%
	SET2	0.06% -0.06% -0.06% 0.28% -0.11% -0.05%	0.00%	0.06%
	SET3	6.28% -4.27% -4.27% 6.11% -5.14% 24.83% 0.00% 21.92%
	SET4	1.23% -1.14% -1.14% 3.40% -0.99% 4.34%	0.00%	4.76%

seconds, and for the last two sets 500 seconds. First of all, note that MC improves the ( , S) bound poorly in most of the instances. Also note that the 1st LR is solved to optimality for all these test problems, and as the table indicates, this approximation of LD 1 does not often provide an improvement over MC. This is due in part to poor multipliers generated from the ( , S) formulation.On the other hand, FL improves in general the ( , S) bound more significantly than MC, although the improvements are still minuscule. Note that 2nd LR does not solve to optimality for many test instances, particularly for the hard problems. Similar to the 1st LR, the 2nd LR does not provide necessarily an improvement over FL bound, due to poor multipliers. Compared to previous test problems, Multi-LSB instances are parallel to TDS problems, where the

Acknowledgement. The research carried out was supported in part by the National Science Foundation grant No. DMI 0323299. The authors are also thankful to two anonymous referees for their comments leading to an improvement of the presentation of the paper.