Antoine Perraudat
email: antoine.perraudat@emse.fr

Stéphane Dauzère-Pérès
email: dauzere-peres@emse.fr

Philippe Vialletelle
email: philippe.vialletelle@st.com

Robust tactical qualification decisions in flexible manufacturing systems

Keywords: Robust Optimization, Duality, Flexible Manufacturing Systems, Qualification Management

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Qualification management

In semiconductor manufacturing, Integrated Circuits (ICs) consist of transistors that are made in "front-end" factories. ICs are built on silicon wafers, and up to one thousand operations are required to complete the fabrication of one wafer. Each operation is performed in a work center grouping parallel machines. A work center ranges from a few machines to two hundred machines. As the number of operations is larger than the number of work centers, semiconductor factories are characterized by a high degree of re-entrant product flows in work centers. A wafer can visit more than forty times the same work center. To perform an operation, a "recipe" is run by a machine on a wafer. The recipe defines the pressure conditions, temperature conditions, chemicals and associated actions necessary to perform the operation. There might be more than one thousand different operations and recipes in a work center.

However, machines cannot simply run recipes once they are purchased and installed in the factory. They must be qualified to meet quality and yield requirements. In other words, when a machine is qualified, it is certified that the machine can run the recipe without deteriorating the product being manufactured on the wafer. If a machine is not qualified for the recipe, the machine cannot process the product at the operation. A qualification then corresponds to a couple (operation, machine). Note that not all recipes are "qualifiable" on a machine, i.e. only a subset of recipes can be qualified on a machine.

Satisfying the demand associated to each product is difficult in semiconductor manufacturing. Several hundred products compete for the same production machines in high mix manufacturing facilities. In addition, the demand by product is time-varying, often significantly from one month to another, and can be highly uncertain [START_REF] Ponsignon | Simulation-based performance assessment of master planning approaches in semiconductor manufacturing[END_REF]. There are also manufacturing risks (e.g. machine breakdowns, yield losses) that can prevent manufacturing facilities from satisfying the demand. When such conditions are met, the need for flexibility (the ability to respond effectively to changing circumstances, see [START_REF] Sethi | Flexibility in manufacturing: a survey[END_REF]) is imperative [START_REF] De Toni | Manufacturing flexibility: a literature review[END_REF]. Qualification management is closely related to the the notion of production flexibility, which is defined as all products a factory is able to produce without requiring additional major capital investment. Production flexibility is the result, among others, of process flexibility, which is defined as the ability of processing different products at the same time [START_REF] Boyer | Manufacturing flexibility at the plant level[END_REF][START_REF] Sethi | Flexibility in manufacturing: a survey[END_REF][START_REF] Jain | A review on manufacturing flexibility[END_REF]. Adding new qualifications improves the level of process flexibility of work centers and therefore improves the capacity of a factory to satisfy the demand.

In this paper, we are interested in the qualification optimization problem that typically arises at a tactical decision level where the planning horizon is between six and twelve months. The considered qualification optimization problem is a tactical capacity planning problem: The production capacity of a work center must be configured to satisfy the demand. There are existing machines in the work center, and new machines might be installed. Similarly, new products are being introduced in the factory, and new qualifications are necessary to increase the production capacity of new products and the production capacity of existing products with a ramp-up demand. This is because new qualifications enable operations associated to the product to be processed on more machines. More precisely, a set of new qualifications, i.e. new couples (operation, machine) to qualify, must be determined so that the demand for all products is satisfied while respecting production capacity constraints. The couple (operation, machine) must be either determined as to be qualified or not to be qualified.

Because qualifications can be expensive and time-consuming, between one week and several months mainly in the form of delay, the number of new qualifications to perform must be minimized and anticipated. Moreover, the demand by product, which is an external parameter to the company, is affected by uncertainty. In factories with a high product mix, i.e. many products, the uncertainty on the demand by product is particularly strong, as factories face frequent product mix changes with products that have short lifetimes. In other words, the set of qualifications determined to satisfy a nominal demand by product may be inappropriate if the realized demand by product is too different from the nominal demand by product. A significant change in the demand can significantly decrease the manufacturing performances. This is because the wafer of a product does not lead to the same workload of a wafer of another product due to different re-entrant flow factors and throughput rates [START_REF] Kotcher | Capacity planning in the face of product-mix uncertainty[END_REF]. Determining a "robust" set of new qualifications, which covers the uncertainty on the demand, is therefore also critical.

Related work 1.2.1. Process flexibility

Qualification management is closely related to the the notion of process flexibility. The scientific literature on process flexibility is mostly interested in measuring the performances of process flexibility designs (which could be called qualification configurations or designs in this paper) in terms of expected service levels using notably linear programming and max-flow models. The term "link" is preferred to the term qualification. In general, the literature on process flexibility deals with strategic problems at the supply chain level. Links are determined between products and factories. The quality of the links (the quality of the process flexibility design) between products and factories is evaluated. Link costs are constrained to a given budget. For instance, if n is the number of factories and products, then 2-chain designs considers at most 2n links. From a general point of view, the literature shows that a production system with limited process flexibility can achieve almost the same performances as a fully flexible system [START_REF] Wang | A review of flexible processes and operations[END_REF].

Under balanced (same number of factories and products) and symmetrical assumptions (each unit of product leads to the same amount of workload at any plant, and plants have the same production capacity), given a set of demand scenarios (demand is assumed to be independent and identically distributed), the seminal work of [START_REF] Jordan | Principles on the benefits of manufacturing process flexibility[END_REF] shows that effective sparse flexibility designs with at most 2n links can almost achieve the same benefits as full flexibility designs. In particular, they show that 2-chain designs (also referred as long chain designs in the literature) where each product is exactly linked to two factories (see Figure 1a) and where the design forms undirected cycle containing all machines and products, is almost as effective as the full flexibility designs (see Figure 1b) with much fewer links. They also show that there can exist multiple process flexibility designs with similar performances. Chain designs perform better than other sparse designs as they pool more products and factories, thus allowing to better face demand uncertainty. Based on this work, [START_REF] Boyer | Manufacturing flexibility at the plant level[END_REF], [START_REF] Graves | Process flexibility in supply chains[END_REF], [START_REF] Chou | Design for process flexibility: Efficiency of the long chain and sparse structure[END_REF], [START_REF] Simchi-Levi | Understanding the performance of the long chain and sparse designs in process flexibility[END_REF], [START_REF] Wang | Process flexibility: A distribution-free bound on the performance of k-chain[END_REF], [START_REF] Désir | Sparse process flexibility designs: is the long chain really optimal?[END_REF] and [START_REF] Bidkhori | Analyzing process flexibility: A distribution-free approach with partial expectations[END_REF] further study, validate and complement the benefits of sparse, chain and long chain flexibility designs. Interestingly, [START_REF] Chan | Process flexibility in baseball: The value of positional flexibility[END_REF] propose an analogy between flexibility in baseball, called positional flexibility, and process flexibility: "a baseball team can be viewed as a production network in which players (plants) produce innings-played to satisfy the demand for all positions (products) on a team." The authors show that positional flexibility, mostly through long (sub)chains, can contribute by itself to a few wins by season, notably by covering uncertainty sources such as injuries. Nevertheless, the main limits of chain flexibility designs for direct applications to qualification management in a manufacturing facility are:

• Most often, balanced systems (same number of factories and products) are studied, which is unrealistic in the work center of a semiconductor manufacturing facility.

• Most often, any factory can be linked to any product. This is impossible in qualification management in a manufacturing facility due to continuous investment. Machines belong to different generations, have different software and hardware restrictions and can be of different types.

They cannot be all qualified for the same fabrication operations. Consequently, chain designs are unlikely.

• Most often, only symmetrical systems. In such systems, products have the same demand distribution, each product unit leads to the same workload at any plant, all plants have the same production capacity, and the mean demand is equal to the total capacity. This is typically not true for instance in semiconductor manufacturing facilities, as two products can require different operations with different processing times. In addition, the mean demand in terms of workload is always smaller than the production capacity to control the fabrication time (see Hopp and Spearman 2011 and Section 3.1).

• Link delays are not considered. The qualification process of an operation on a machine may take several weeks to several months.

• Single period models are often considered. However, demands of products are highly dynamic (see e.g. Figure 2 in Section 2.2), which cannot be easily captured with single period models.

• Demand is often assumed to be independent and identically distributed. In high-mix (with a large portfolio of products) factories, demands are not independent and identically distributed. Typically, a few products are associated to most of the demand.

Other works contribute to the process flexibility literature by releasing some of these assumptions (but never all of them at the same time) to study more "general" manufacturing systems. [START_REF] Mak | Stochastic programming approach to process flexibility design[END_REF] propose a two-stage stochastic programming approach to determine process flexibility designs. The studied setting is a balanced system. Process flexibility costs are distinguished by factory and product and factories have different production capacities. They show that, when the demand by product is heterogeneous, the flexibility design determined with the stochastic programming approach generates a better profit than chain designs. For an unbalanced and unsymmetrical system, [START_REF] Chou | Design for process flexibility: Efficiency of the long chain and sparse structure[END_REF] identify underlying conditions such that sparse (not necessarily chained) flexibility designs achieve most of the benefits of the full flexibility design for an unbalanced and unsymmetrical system in a single period setting. They also show that adding a restricted number of links is often sufficient to significantly improve the ability of a production system to meet the demand. [START_REF] Deng | Process flexibility design in unbalanced networks[END_REF] formulate recommendations for process flexibility designs for unbalanced but symmetrical systems. [START_REF] Bidkhori | Analyzing process flexibility: A distribution-free approach with partial expectations[END_REF] derive a lower bound for chain designs when systems are unbalanced and factories have different production capacities. [START_REF] Chen | Optimal design of process flexibility for general production systems[END_REF] further study unbalanced and unsymmetrical systems by proposing a simple scheme to satisfy the expected demand with high probability in a single-period setting. [START_REF] Shi | Process flexibility for multiperiod production systems[END_REF] study flexibility designs in a multi-period setting for an unbalanced manufacturing system. The system is partially unsymmetrical as processing times are not differentiated between products. However, the demand is assumed to be identically distributed across time periods, which is not always realistic (see Section 2). [START_REF] Yan | On the design of sparse but efficient structures in operations[END_REF] propose a distributionally robust approach and a dual variable based heuristic to design effective sparse structures in a single period setting. From the full flexibility design, the dual variable based heuristic removes one link at a time while maximizing the total number of products. On numerical examples of balanced but unsymmetrical systems (processing times between products are not differentiated), they show that the heuristic proposes effective sparse structures, even long chain designs are proposed when enough links are removed. The heuristic performs also well for unbalanced systems. [START_REF] Wang | Robust optimization approach to process flexibility designs with contribution margin differentials[END_REF] propose a promising robust optimization approach to determine effective flexibility designs for unbalanced and unsymmetrical systems in a single period setting. Instead of modeling the demand with an uncertainty set and then determining the best flexibility design, the approach builds and ranks flexibility designs that are effective for a class of uncertainty sets. However, the studied class of uncertainty sets (partwise independently symmetric uncertainty sets) is not representative of demand uncertainty in our industrial context. [START_REF] Fiorotto | Process flexibility and the chaining principle in lot sizing problems[END_REF] present a deterministic lot-sizing problem motivated by the semiconductor industry. They propose two different lot-sizing optimization models to build the best long chain configuration or to find the best links (the total number of links is limited to a given number) while minimizing the setup, production, holding and backlogging costs. They analyze different flexibility designs and compare them to different long chain designs [START_REF] Jordan | Principles on the benefits of manufacturing process flexibility[END_REF]. They show that, when the capacity is tight or when inventory and backlogging costs are very different from one product to another, scenarios that are actually frequently encountered in high-mix factories, even the best long chain design is not satisfactory. Flexibility links can be misplaced because backlogging costs and setup times are not considered in the long chain principle. The authors show that the optimization obtains better cost effective designs with half the links used by the long chain design.

There are other features of high-mix semiconductor manufacturing that are not considered in the literature. For instance, the demand by product is also considered as uncertain, but another layer of complexity is added by the fact that qualifications are carried out for operations which may be common to some products. Probability distributions are not easy to obtain as new products are frequently introduced and are subject to product cannibalization (see Section 2). In addition, qualifying any operation on any machine, which is assumed in the literature, is not possible due to technological restrictions.

To improve the realism and for a relevant usability in semiconductor manufacturing, all of these assumptions should be ignored. As it is unlikely to determine analytic formulas under such conditions to help determine relevant process flexibility designs, and therefore relevant qualification configurations, solving complex combinatorial optimization problems is required as shown in [START_REF] Mak | Stochastic programming approach to process flexibility design[END_REF] and [START_REF] Fiorotto | Process flexibility and the chaining principle in lot sizing problems[END_REF].

Qualification management in semiconductor manufacturing

The literature is scarce on the design of qualification configurations in semiconductor manufacturing, in particular when the demand is uncertain. Stochastic programming has been a method of choice so far to deal with the uncertainty on the demand. [START_REF] Klemmt | A multistage mathematical programming based scheduling approach for the photolithography area in semiconductor manufacturing[END_REF] propose to design qualification configurations for a specific work center by covering a few scenarios on the demand by product, which is a common practice in the semiconductor industry. Nevertheless, the approach is not entirely detailed. [START_REF] Chang | Stochastic programming for qualification management of parallel machines in semiconductor manufacturing[END_REF] propose a stochastic programming optimization approach to maximize the weighted expected number of processed product quantities. The demand and the production capacities are subject to uncertainty. In addition, they consider that new qualifications lead to a stochastic capacity loss that can be described with a distribution probability. However, the approach proposed by [START_REF] Chang | Stochastic programming for qualification management of parallel machines in semiconductor manufacturing[END_REF] cannot be used at a tactical level. This is because their stochastic model does not ensure that the demand by operation has to be satisfied. Then, only a fraction of the operations corresponding to a product could be qualified, and the product could potentially never be delivered. [START_REF] Fu | Stochastic optimization of product-machine qualification in a semiconductor back-end facility[END_REF] also consider that the demand is uncertain in a qualification management optimization problem. Nevertheless, the problem is treated from an extended production planning standpoint and not from a capacity planning standpoint. Consequently, similarly to [START_REF] Chang | Stochastic programming for qualification management of parallel machines in semiconductor manufacturing[END_REF], the work of [START_REF] Fu | Stochastic optimization of product-machine qualification in a semiconductor back-end facility[END_REF] cannot be used at a tactical level. [START_REF] Liao | Optimizing site qualification across the supply network at western digital[END_REF] propose a two-stage stochastic programming optimization approach to maximize the total profit of a semiconductor company. The first stage problem consists in minimizing qualification costs while second stage problem consists in allocating product quantities to production sites to maximize revenue.

However, stochastic programming implies characterizing demand scenarios and associated probabilities. This is difficult as products tend to have dependent demands due to product cannibalization, which is not mentioned in the literature. Product cannibalization is particularly critical for manufacturers with a high product mix. Determining nominal demands and plausibility limits is a promising alternative: It is as natural as defining demand scenarios without requiring probabilities and can consider product cannibalization.

Contributions

Our contributions to the qualification management and robust optimization literature are as follows:

• We propose a new mixed integer linear programming mathematical model for the tactical qualification management problem when the demand is deterministic and the qualification lead times are considered. We show that the studied problem is NP-Hard.

• As the demand by product can be subject to uncertainty, we motivate the choice of robust optimization for the considered problem. We propose an uncertainty set based on the budget of uncertainty [START_REF] Bertsimas | The price of robustness[END_REF] to cover the demand uncertainty. A novelty of our approach is to take product cannibalization into account, rarely considered in the literature.

• We propose a new static robust reformulation of the deterministic model when the demand is considered as uncertain but can be described by D t .

• We propose a new decision-dependent uncertainty linear program to characterize the robustness of a set of qualifications. As the problem is NP-complete, a binary search solution approach is proposed when the uncertainty on the demand is symmetrical.

• In the computational study, we show on industrial data that the price of uncertainty is small, often less than a few qualifications, whereas the qualifications determined for the nominal demand often lead to capacity constraint violations.

The remainder of the paper is organized as follows. In Section 2, we describe and motivate the type of demand uncertainty faced in semiconductor manufacturing. We motive the use of robust optimization to cover demand uncertainty. In Section 3, the deterministic mathematical model is presented. Then, a mathematical robust optimization approach is proposed to cover demand uncertainty. In Section 4, we propose a mathematical model and discuss several approaches to determine the robustness of a given set of qualifications (e.g. the set of initial qualifications). In Section 5, a computational study on industrial data is conducted to evaluate the price of uncertainty [START_REF] Gorissen | A practical guide to robust optimization[END_REF], the practical tractability of the proposed optimization models, and possible the capacity constraint violations and consequences if the set of qualifications obtained by solving the deterministic optimization problem is used. In Section 6, we discuss how the proposed optimization models can be used for a practical use by capacity planners in a decision support system. Finally, in Section 7, we conclude and give some perspectives.

Uncertainty on the demand

Demand uncertainty and product cannibalization

Processing times, production capacities, qualification lead times and the demand by product can be subject to uncertainty. In this paper, only the demand uncertainty is considered, which is critical to a manufacturing company. The uncertainty on the demand is an external uncertainty, which is difficult, if not impossible, to control with discount prices and incentives even if the product is innovative. Considering the uncertainty on other parameters is left for future research.

Note that the uncertainty on the demand by operation is a consequence of the uncertainty on the demand by product. In the semiconductor industry, operations need to be run to process a product. However, all products do not share the same operations. Moreover, although two products share common operations, operations will not have the same processing times. This is due to differences in the re-entrant product flows. We are therefore interested in characterizing and modeling the demand uncertainty and linking it to the uncertainty on the demand by operation.

Although it is possible to accurately predict the total quantity of products that a manufacturing facility must complete in the future, it is often impossible to exactly know the quantity of each product. One important reason why is that high-tech companies such as semiconductor manufacturers with a large portfolio of products often face product cannibalization [START_REF] Moorthy | Market segmentation, cannibalization, and the timing of product introductions[END_REF][START_REF] Kim | Commonality in product design: Cost saving, valuation change and cannibalization[END_REF]. Product cannibalization occurs when a company manufactures different products that compete with each other on the same market. Consider the following example. A client that seeks to design an electronic system has the choice between several micro-controller among those that the company sells. A micro-controller is integrated circuit with essentially the same features as modern computers, i.e. computing unit, memory, input and output interfaces, but are dedicated to specific applications and require little energy. Several micro-controllers are suitable for a given application, and the final choice will be made based on cost, energy consumption and memory among other characteristics. The client will probably never buy all suitable micro-controllers. Therefore, selling one unit of a product may mean selling fewer units of other products. Nevertheless, a product cannot be replaced by any other product because all products are not used for the same application. Some products will be used in the automotive industry, whereas others will be used for industrial applications in factories, or telecommunication applications. Products are distinguished by their family. A product family is then a set of products that have similar characteristics, can be used for similar applications, and therefore compete on the same market segment.

Managing the demand uncertainty

To cover the demand uncertainty, two main methods exist: Stochastic optimization and robust optimization. Stochastic optimization assumes that the probability distribution of the demand uncertainty is known. Then, in general, the expected value of the objective function is optimized. In this paper, the objective would consist in minimizing the expected number of qualifications after generating, possible many, scenarios from the estimated probability distribution [START_REF] Birge | Introduction to stochastic programming[END_REF]. Robust optimization is different because the probability distribution of the uncertainty is not required. In robust optimization, the objective consists in minimizing the objective function while ensuring that the constraints are never violated [START_REF] Ben-Tal | Robust optimization-methodology and applications[END_REF][START_REF] Ben-Tal | Robust optimization[END_REF][START_REF] Bertsimas | Theory and applications of robust optimization[END_REF][START_REF] Gorissen | A practical guide to robust optimization[END_REF].

Robust optimization is more relevant when determining a set of qualifications at the tactical decision level. First, estimating the probability distribution of the demand of a product when it is correlated to the demands of other products is difficult. Furthermore, estimating the probability distribution of the demand for new products is difficult. This is because semiconductor manufacturers may not have enough data on the demands to derive relevant distribution probabilities as they experience frequent product mix changes (Bertsimas and Thiele, 2006a), i.e. the demand for a product strongly varies from one month to another. The demand for a product is therefore, in general, not identically distributed over time. Figure 2 provides an illustrative example using historical industrial data on the changes of the demand for one product over 12 months. For confidentially purposes, product names are not mentioned. In addition, the monthly demand is divided by the mean demand over the 12 months. The mean demand is several hundred units. The demand for this product is particularly interesting. There is a large, quick and intense ramp-up demand. Nevertheless, the demand quickly fades away. These demand fluctuations are critical, especially since their intensities are extremely difficult to predict in advance. This is when a robust optimization based approach is relevant to cover the demand uncertainty. Second, it is critical to anticipate relevant qualifications to cover the demand uncertainty. This is because, in general, it is important to perform the right qualifications and not all qualifications to respect capacity constraints and satisfy the demand [START_REF] Jordan | Principles on the benefits of manufacturing process flexibility[END_REF][START_REF] Benjaafar | The effect of routeing and machine flexibility on the performance of manufacturing systems[END_REF][START_REF] Graves | Process flexibility in supply chains[END_REF][START_REF] Chou | Design for process flexibility: Efficiency of the long chain and sparse structure[END_REF][START_REF] Johnzén | Flexibility measures for qualification management in wafer fabs[END_REF][START_REF] Fiorotto | Process flexibility and the chaining principle in lot sizing problems[END_REF][START_REF] Chen | Optimal design of process flexibility for general production systems[END_REF].

Furthermore, as qualification decisions are made at a tactical decision level, they have a major impact on all production planning and control management issues [START_REF] Hopp | Factory physics[END_REF]. For instance, if new qualifications are not properly determined, then effective robust production plans may not be found to satisfy the demand. Determining the right set of new qualifications is thus critical for manufacturing and financial performances.

Third, in practice, a way to deal with uncertainty is to frequently adjust the current set of qualifications by performing new qualifications when the nominal demand is updated. However, this is not always possible because the qualification process may sometimes take several weeks or months to validate the quality and the yield of the operation. Therefore, if the demand is updated late, it may be impossible to perform additional qualifications to satisfy the demand. Then, anticipating the right qualifications to cover the demand uncertainty is critical. Also, determining a set of robust qualifications could save critical time for capacity planners. This is because the set of qualifications would be determined in a less reactive manner but in a more proactive manner against demand changes. Capacity planners could therefore be assigned to other tasks. Note that the set of qualifications would still need to be adjusted when completely new products are introduced or old products are reintroduced because of unnoticed disqualifications.

Problem modeling

Problem description

Let us consider a work center of M unrelated parallel machines, both in terms of qualifications and throughput rates, which must process R different operations. Machines are unrelated because they are of different generations. A demand is associated to each operation on the considered horizon. The horizon consists of T periods. The work center is asymmetrical and unbalanced, i.e. the demand varies from one operation to another and the number of operations is much greater than the number of machines. A machine can only process qualified operations, and a "qualifiable" operation can only be processed on a machine if it is qualified. Qualifying an operation on a machine induces a qualification cost and is subject to a qualification lead time. The qualification matrix defines the initial set of active qualifications. A qualification is therefore a pair (operation, machine). The initial set of active qualifications is known and deterministic. Each machine has a finite production capacity that must be respected at each period on the considered horizon.

The objective is to minimize the total cost of the qualifications to perform, among the qualifiable pairs (operation, machine) not already qualified, while meeting demand and respecting capacity constraints.

This problem will be referred as the Minimum Cost Qualification Configuration Problem (MCQCP) in the remainder of the paper.

Deterministic modeling

Parameters: M : Number of machines, R: Number of operations, P : Number of products, T : Number of periods, q r,m : Is equal to 1 if machine m is initially qualified for operation r, to 2 if machine m is qualifiable for operation r, to 0 if machine m cannot be qualified for operation r, tp r,m : Throughput rate (per hour) of operation r on machine m, c t,m : Production availability (in hours) of machine m at period t, u max t,m : Maximum utilization rate allowed for machine m at period t, rf p,r : How many times (re-entrant flow factor) operation r needs to be run to produce one unit of product p, d t,p : Demand for product p at period t, l t,r,m : Lead time (in number of periods) when starting qualification procedure at period t of operation r on machine m, δ t : Discount factor at period t, cq r,m : Cost of qualifying operation r on machine m. Decision variables: OQ t,r,m ∈ {0, 1}: Is equal to 1 if there is qualification procedure to start for operation r at period t on machine m, and 0 otherwise, W IP t,r,m ∈ [0, 1]: Ratio of the demand for operation r processed by machine m at period t.

min t,r,m δ t cq r,m OQ t,r,m (1)
s. t. r ( p rf p,r d t,p )W IP t,r,m tp r,m ≤ c t,m u max t,m ∀t, ∀m (2)
m W IP t,r,m = 1 ∀t, ∀r | p rf p,r d t,p > 0 (3) W IP t,r,m ≤ q r,m ∀t, ∀r, ∀m | q r,m = 2 (4) W IP t,r,m ≤ t t =1|t-t ≥l t ,r,m OQ t ,r,m ∀t, ∀r, ∀m | q r,m = 2 (5) W IP t,r,m ≥ 0 ∀t, ∀r, ∀m (6) OQ t,r,m ∈ {0, 1} ∀t, ∀r, ∀m (7)
The objective function (1) minimizes the cost of performing qualifications on the planning horizon. The discount factor is used to decide if qualifications must be made as soon as possible or as late as possible. For instance, assume that qualification procedures must be started as late as possible. This is possible by ensuring that δ t ≥ δ t+1 ∀t ∈ {1, ..., T -1}. Constraints (2) ensure that the capacity constraint for each machine m and each period t is respected. Constraints (2) also limit the utilization rate of machine m at period t to a maximum of u max t,m . This controls the mean cycle time (fabrication time) in the work center as the mean cycle time increases exponentially with the utilization rate, and improves the responsiveness of the work center [START_REF] Hopp | Factory physics[END_REF]. Constraint (3) are the flow constraints. They ensure that the demand by operation must be satisfied. Constraints (3) are active only if there is demand for operation r at period t, ∀t, ∀r | p rf p,r d t,p > 0. For new operations or machines, if p rf p,r d t,p > 0 is not enforced for some periods, then new qualifications will be required because m W IP t,r,m = 1 has to be satisfied even if there is no demand. This condition is therefore used to avoid qualifying operations on machines in the early periods if the demand is only expected in the late periods. Constraints (4)-( 5) are the qualification constraints. They ensure that machine m is qualified for operation r at period t, if it has been newly qualified or was initially qualified while considering qualification lead times. Finally, Constraints ( 6) are the non-negativity constraints and Constraints (7) are the binary constraints.

Let us discuss below some important characteristics of our problem:

• The deterministic optimization model is relevant, although it does not consider demand uncertainty, because it considers essential features of qualifications which are qualification costs and delays, and models unbalanced and unsymmetrical systems. In the computational study on industrial data, we found that the deterministic model is easy to solve (see Section 5) for the considered work centers.

• MCQCP can also be solved factory-wide, i.e. by considering all work centers simultaneously. However, as two different work centers do not share operations, optimality is preserved when breaking down the problem by work center to reduce the size of the problem in terms of machines and operations.

• It is important to mention that MCQCP can be infeasible if the production capacities of machines are too small and if too few qualifiable pairs (operation, machine) exist to better balance the workload between the machines. Note that in the numerical experiments performed in Section 5, MCQCP is always feasible contrary to its robust counterpart.

• The deterministic model can still be used to determine a set of qualifications even if lead times are not modeled, i.e. if l t,r,m = 0 ∀t, ∀r, ∀m. In this case, decision variables OQ t,r,m should be interpreted as the period at which operation r must be qualified on machine m if OQ t,r,m = 1.

• Time-varying lead times are considered to better consider the fact that the periods can have different durations (as it is the case in the considered industrial context) and also different demand quantities, and therefore different loads, which influence both fabrication times [START_REF] Hopp | Factory physics[END_REF] and qualification times as test lots are used for qualifications. Time-varying lead times can also be used by capacity planners to simulate the acceleration or deceleration of qualification procedures during certain periods.

• In order to correctly consider new machines, it is sufficient to set c t,m to appropriate values until machine m is actually started-up in the factory. Start-up periods are notably used for qualification purposes. For instance, if the horizon is of 3 months with 3 periods of one month and the start-up period lasts one month, then c t,m must be equal to zero for the first two months.

• New qualifications can lead to capacity losses in the considered work center as it is required to run quality tasks on machines by using test products. [START_REF] Chang | Stochastic programming for qualification management of parallel machines in semiconductor manufacturing[END_REF] model this aspect by using a probability distribution. As it is complex to define relevant probability distributions, capacity losses due to new qualifications are modeled with available historical data as exogenous factors in the production capacity of each machine. Note that quality tasks are also frequently run even for existing qualifications, which is also considered in the production capacity of each machine.

• There are other sources of capacity loss, such as setups, which are considered in the value of the production capacity of each machine.

• In practice, operations are often subject to precedence constraints, i.e. operations must be performed in a precise order. In this paper, we are not interested in making detailed scheduling decisions but, at a tactical level, in verifying from a capacity planning standpoint that the demand can be met with the current and potentially new qualifications. Precedence constraints are therefore not explicitly considered as in most capacity planning models and for instance in [START_REF] Rowshannahad | Capacitated qualification management in semiconductor manufacturing[END_REF] for qualification management.

Computational complexity

Even for a single period, it is possible to show that MCQCP is a NP-Hard problem by reducing MCQCP to the Generalized Assignment Problem (GAP), known to be NP-Hard (see e.g. [START_REF] Nauss | Solving the generalized assignment problem: An optimizing and heuristic approach[END_REF].

Let us state GAP in terms of tasks and agents, where I is the number of tasks and J the number of agents. Let c i,j be the cost of assigning task i to agent j, r i,j the processing time required for task i by agent j, and b j the total capacity of agent j. Task i is assigned to agent j when X i,j = 1, and 0 otherwise. The question is "is there an assignment of tasks to agents such that the total assignment cost is equal to K, i.e.

i,j c i,j X i,j = K, and the constraints are satisfied, i.e. i r i,j X i,j ≤ b j ∀j, and j X i,j = 1 ∀i?".

Given a general instance of GAP, it is possible to build an instance of MCQCP as follows: Let T = 1, P = 1, the number of operations be equal to the number of tasks, i.e. R = I, and the number of agents be equal to the number of machines, i.e. M = J. Let us use the subscripts i and j in the remainder of this section. We can therefore set tp i,j = 1 r i,j ∀i, ∀j, c j = b j ∀j, d 1 = 1, q i,j = 2 ∀i, ∀j, l i,j = 0 ∀i, ∀j, and rf 1,i = 1 ∀i, u max j = 1 ∀j. Finally, let us set cq i,j = g + c i,j ∀i, ∀j, with g a large fixed cost such that g > i,j c i,j . The question is "is there a set of qualifications such that the total qualification cost is equal to K + Ig, and the capacity, flow and qualification constraints are respected?"

Assume GAP has a yes answer with a total assignment cost of K. The solution for GAP is also feasible for MCQCP because tasks and operations have the same throughput rates on agents and machines, and because the demand for each operation is equal to one unit. Therefore, MCQCP also has a yes answer with a total cost of K + Ig.

Assume that MCQCP has a yes answer. Because g is a large fixed cost such that g > i,j c i,j , each operation is qualified on one and only one machine, and therefore the total production flow of each operation is assigned to one and only one machine. In addition, because tasks and operations have the same throughput rates on agents and machines, GAP also has a yes answer with a total cost of K.

Robust modeling 3.4.1. Polyhedral uncertainty with budget of uncertainty

To consider demand uncertainty and product cannibalization, a polyhedral uncertainty set, based on budget uncertainty proposed by [START_REF] Bertsimas | The price of robustness[END_REF] are the plausibility limits for product p at period t. The uncertainty set D t that models the effect of product cannibalization by product family at period t is described below:

D t = {d t,p | d t,p ≥ d t,p -d t,p ∀p, d t,p ≤ d t,p + d t,p ∀p, p|α p,f =1 d t,p ≤ Γ t,f ∀f } (8)
In D t , the total demand by product family f at period t is limited to the budget of uncertainty Γ t,f , which is the maximum demand to cover for product family f at period t. Therefore, if the demand for a product in the product family increases above its nominal value, then the increase is made at the expense of another product in the product family, whose demand must decrease. In addition, if Γ t,f = p|α p,f =1 d t,p , then, for each product family f , the maximum overall quantity to produce is equal to the overall quantity in the nominal case, but the distribution of the demand between the products in the product family is unknown. Setting Γ t,f = p|α p,f =1 d t,p is a practical assumption. This ensures that qualifications are not determined to cover extreme cases where the quantity by product family would actually be much larger than the nominal quantity by product family, which is often unrealistic. Instead, qualifications are optimized to cover any demand realization given an overall quantity by product family. Note that, although the uncertainty set D t ensures that the total demand of all products in a family is not large, the total demand over all operations arriving in work centers can significantly increase as re-entrant flow factors significantly vary from one product to another.

Parameters d t,p and d t,p do not necessarily reflect the real uncertainty on the demand of product p at period t. They can be defined in a such way that they correspond to the uncertainty capacity planners want to manage if the real uncertainty is too expensive to cover [START_REF] Bertsimas | The price of robustness[END_REF].

Static reformulation

We investigate a static reformulation of the deterministic optimization problem. We follow Ben-Tal and Nemirovski (2002), [START_REF] Gorissen | A practical guide to robust optimization[END_REF] and [START_REF] Yanıkoglu | A survey of adjustable robust optimization[END_REF] to write the robust formu-lation of MCQCP. First, constraints with uncertain parameters, the demand, need to be identified, then the robust counterpart can be derived.

There are two constraints with uncertain parameters: The flow constraints (3), and the capacity constraints (2). Flow constraints (3). The demand is used to control when the flow constraint must be active. To make sure the flow constraints hold for any demand realization within D t , it is sufficient to replace the condition

m W IP t,r,m = 1 ∀t, ∀r | p rf p,r d t,p > 0 by m W IP t,r,m = 1 ∀t, ∀r | p rf p,r (d t,p + d t,p ) > 0.
Capacity constraints (2). If the demand uncertainty is row-wise and the uncertainty set is compact, then an optimal solution for the static reformulation problem is also an optimal solution for the adjustable robust reformulation problem [START_REF] Ben-Tal | Adjustable robust solutions of uncertain linear programs[END_REF][START_REF] Yanıkoglu | A survey of adjustable robust optimization[END_REF] Robust counterpart: The next step consists in determining the robust counterpart of the capacity constraints. The robust counterpart is independently determined from one capacity constraint to another. Consider one capacity constraint for a given machine m at period t:

Step 1 (worst-case reformulation):

max d∈Dt p d t,p ( r rf p,r W IP t,r,m tp r,m ) ≤ c t,m u max t,m
Intuitively, covering the worst-case realization in the uncertainty set D t will conduct to add qualifications to machines for operations common to many products, or for operations associated to products with large re-entrant flow factors, as they are the operations that will impact the most the utilization rate of machines.

Step 2 (duality):

The next step consists in taking the dual of the inner maximization problem. The inner maximization problem and its dual, which is a minimization problem, have the same objective value because the inner maximization problem is linear. For a given period t, the following optimization problem must be solved:

max p d t,p ( r rf p,r W IP t,
p|α p,f =1 d t,p ≤ Γ t,f ∀f (11)
The dual variables associated to Constraints ( 9)-( 11) are listed in Table 1. The dual of the inner maximization problem is a minimization problem. The minimization problem for a given capacity constraint for machine m at period t is modeled below: Step 3 (Robust Counterpart): The final step consists in omitting the minimization term to obtain the robust counterpart. Therefore, the robust counterpart of the capacity constraint for a given machine m and a given period t can be found below:

p (-(d t,p -d t,p )y min p ) + f (Γ f y gamma f ) + p ((d t,p + d t,p )y max p ) ≤ c t,m u max t,m -y min p + y max p + f |α p,f =1 y gamma f ≥ r rf p,r W IP t,
y gamma t,m,f ≥ 0 ∀t, ∀m, ∀f (16)
The objective function ( 12) minimizes the cost of performing qualifications, while Constraints ( 13)-( 14) are the "robustification" constraints. Constraints (15) ensure that the demand by operation must be satisfied if there is demand. Note that Constraints (15) are slightly different from Constraint (3) as it must be active when p rf p,r (d t,p + d t,p ) > 0 for operation r at period t instead of p rf p,r (d t,p ) > 0. Constraints ( 16)-( 17) correspond to the non-negativity constraints introduced by the "robustification" procedure.

Note that the robust optimization model ( 12)-( 15) can still be used when a product belongs to several product families. Only input parameters must be changed.

The robust optimization problem will referred as the Minimum Cost Robust Qualification Configuration Problem (MCRQCP) in the remainder of the paper.

Similarly to MCQCP, it is important to mention that MCRQCP can be infeasible if the production capacities of machines are too small and if too few qualifiable pairs (operation, machine) exist to better balance the workload between the machines. Note that in the numerical experiments performed in Section 5, MCRQCP is infeasible for some values of d t,p and d t,p .

Illustrative example on tractability

MCQCP and MCRQCP can be both modeled with mixed integer linear programs, and thus can be solved by standard solvers. Although no new binary variables are required in the robust reformulation of MCQCP, reformulating capacity constraints can modify the problem structure and introduce many more variables and constraints. The reformulation can also worsen the quality of linear relaxations, thus increasing the computational time required to reach an optimal solution in a branch and cut algorithm. It is then expected that MCRQCP requires more computational time to be solved than MCQCP. MCQCP has

T × M + T × R + 4 × T × R × M constraints, T × R × M continuous variables and T × R × M binary variables. MCRQCP has 2 × T × M × P + T × M × F more continuous variables and 3 × T × M × P + T × M × F more constraints than MCQCP.
Table 2 illustrates the additional computational effort required to solve MCRQCP by reformulating capacity constraints with respect to MCQCP in terms of number of decision variables and constraints.

The number of decision variables and constraints of MCQCP and MCRQCP are given for P = 238, R = 1208, F = 3, M = 20, T = 7. These values come from one work center (work center A) in the computational study. Assuming that the demand and worst-case demand are greater than 0 for all products and all periods, the increase of the number of continuous variables is equal to 16.6% and the increase of the number of constraints is equal to 12.8%. This makes the robust optimization problem more difficult to solve than the deterministic optimization problem as the robust optimization problem also tighten the capacity constraints. In practice, the robust optimization problem is much more difficult to solve than the deterministic optimization problem although most optimal solutions can be found in one hour (see Section 5.3.2). The first option consists in using values estimated by decision-makers in charge of defining and predicting future demands. However, determining relevant values can be difficult for some products, in particular for new products because data can be insufficient.

If it is too difficult to provide relevant d t,p for each product, another option is to propose initial values for d t,p . d t,p can first roughly initialized, e.g. initialized to d t,p , and then refined by characterizing the robustness of a set of qualifications (typically the set of initial qualifications) with respect to the demand uncertainty. More precisely, characterizing the robustness of a set of qualifications means determining to what extent a work center is able to correctly absorb the demand uncertainty. Characterizing a set of qualifications is similar to determining the largest d t,p for each product p at period t. Therefore, determining the robustness of a set of qualifications provides capacity planners with the tolerated changes on the demand by a work center. Then, if possible, demand changes should be made in the bounds defined by d t,p and d t,p to limit additional costs with outsourcing or new machines.

In the context of qualification management, [START_REF] Rossi | A robustness measure of the configuration of multi-purpose machines[END_REF] and [START_REF] Aubry | Maximizing the configuration robustness for parallel multi-purpose machines under setup cost constraints[END_REF] assume that satisfying the demand by product is a key issue to characterize the robustness of a set of qualifications. [START_REF] Rossi | A robustness measure of the configuration of multi-purpose machines[END_REF] seeks to characterize the robustness of a set of qualifications by determining the minimum additional quantity of products from the nominal demand that can be absorbed without the makespan exceeding a specified value. Robustness is defined as a distance in [START_REF] Rossi | A robustness measure of the configuration of multi-purpose machines[END_REF]. Similarly, [START_REF] Aubry | Maximizing the configuration robustness for parallel multi-purpose machines under setup cost constraints[END_REF] seeks to characterize the robustness of a set of qualifications by determining the largest perturbation from the nominal demand while ensuring that all machines have the same workload and that qualification costs do not exceed a predefined value. The L-1 norm is used. Similarly to [START_REF] Rossi | A robustness measure of the configuration of multi-purpose machines[END_REF] and [START_REF] Aubry | Maximizing the configuration robustness for parallel multi-purpose machines under setup cost constraints[END_REF], we assume that satisfying the demand by product is a key issue when characterizing the robustness of a set of qualifications. The major differences with [START_REF] Rossi | A robustness measure of the configuration of multi-purpose machines[END_REF] and [START_REF] Aubry | Maximizing the configuration robustness for parallel multi-purpose machines under setup cost constraints[END_REF] are that: (1) We do not assume that machines are uniform or related; (2) We consider large scale production systems with hundreds of products and thousands of operations; (3) Product cannibalization and correlated demands are considered. To characterize the robustness of a set of qualifications, we resort to robust optimization and the uncertainty set D t . More precisely, we seek to determine to what extent a set of qualifications is able to absorb the demand uncertainty when it is described by the uncertainty set D t . Assessing the robustness of a set of qualifications depends on the utility function used to evaluate it. First, we propose a generic mathematical model to model the robustness of a set of qualifications with respect to the demand uncertainty. Second, we propose a solution approach, based on a binary search approach, to determine the robustness of a set of qualifications.

Problem statement

The problem is mostly identical to the problem introduced in Section 3.1. The only difference is that the objective is to characterize the robustness of a set of given qualifications. This problem will be referred as the Maximum Robustness Budgeted Qualification Problem (MRBQP) in the remainder of the paper.

Problem modeling

Let us introduce a new decision variable θ t,p ≥ 0 ∀t, ∀p that is used to evaluate the robustness of a set of qualifications. Let us assume that d t,p is an uncertain parameter that depends on θ t,p : d t,p ∈ [d t,p -d t,p θ t,p , d t,p + d t,p θ t,p ] ∀t, ∀p. Let β t,p ≥ 0 ∀t, ∀p be a weight whose value corresponds the preferences of capacity planners when evaluating the robustness of a set of qualifications. The larger β t,p , the larger the emphasis on the robustness of product p at period t. Let f (θ) = t,p β t,p θ t,p be a utility function that evaluates the robustness of a set of qualifications, where θ = (θ 1,1 , ..., θ T,P ). Formally, the problem can be modeled as follows:

max f (θ) (18)
W IP t,r,m ≤ q r,m ∀t, ∀r, ∀m | q r,m = 2 (21) W IP t,r,m ≤ 0 ∀t, ∀r, ∀m | q r,m = 2 (22) θ t,p ≤ 1 ∀t, ∀p (23)
W IP t,r,m ≥ 0 ∀t, ∀r, ∀m (24)
θ t,p ≥ 0 ∀t, ∀p (25)
The objective function (18) maximizes the utility function f (θ). The capacity constraints (19) depend on θ, which is used to control the demand uncertainty. Constraints (20) model the flow constraints, while Constraints ( 21) and ( 24) model the qualification constraints. Constraints (23) ensure that the demand by product cannot be negative. Finally, Constraints (24) and ( 25) are the non-negativity constraints.

Solving MRBQP is equivalent to determining the robustness of the initial set of qualifications, or any set of qualifications as input parameter.

Robust counterpart

MRBQP is an optimization problem under decision-dependent uncertainty because the capacity constraints depend on the matrix θ used to control the demand uncertainty to cover. Optimization problems under decision-dependent uncertainty are known to be difficult to solve. When the uncertainty set is polyhedral, the objective function is linear and the constraints are linear, the optimization problem is NP-Complete [START_REF] Nohadani | Optimization under decision-dependent uncertainty[END_REF][START_REF] Lappas | Robust optimization for decision-making under endogenous uncertainty[END_REF].

Nevertheless, as in the classical robust optimization paradigm, it is possible to reformulate decisiondependent uncertainty constraints with duality. This is because θ is not a decision variable of the inner robust maximization problem. Let us consider the same uncertainty set as in Equations ( 8). The only difference stems from the fact that the plausibility limits of d t,p are now dependent on θ t,p . Similarly to Section 3.4.2, it is possible to "robustify" the capacity constraints (19).

We follow the same procedure as the one in Section 3.4.2, and the same notations for dual variables are used. Steps 1 and 2 are similar to Section 3.4.2. By deriving the robust counterpart of each capacity constraint, it is possible to write the robust reformulation of MRBQP below:

max f (θ) (26)
y gamma t,m,f ≥ 0 ∀t, ∀m, ∀f (29)
The objective function (26) maximizes the robustness of a set of qualifications. Constraints ( 27)-(30) correspond to the "robustification" constraints. They ensure that the capacity constraints must be respected for any realization in the uncertainty set D t .

Solving MRBQP leads to determining the largest θ t,p for product p at period t, and consequently to characterize the robustness of a set of qualifications. The main drawback of MRBQP is that it is computationally challenging to solve. This is because MRBQP contains products of variables, θ t,p and y max t,m,p , and θ t,p and y min t,m,p , which are introduced by the "robustification" procedure for capacity constraints. There are other possible MILP reformulations when one the variables is binary [START_REF] Nohadani | Optimization under decision-dependent uncertainty[END_REF][START_REF] Lappas | Robust optimization for decision-making under endogenous uncertainty[END_REF]. To determine an estimate of the robustness of a set of qualifications, a binary search solution approach is presented in Section 4.3.

Binary search approach

To characterize the robustness of a set of qualifications, it is possible, for each period, to maximize θ t,p assuming that θ t,p = θ t ∀p. For this objective, Algorithm 1, which a binary search like algorithm, can be used when θ 0 , an initial upper bound for θ, is provided. return θ min 18: end procedure The computational difficulty in Algorithm 1 comes from solving multiple large-scale linear programs. The computational burden can be lowered by warm-starts as only the coefficients y min t,m,p , y max t,p,m , and y gamma t,m,f variables in the "robustification" constraints must be changed. Note that if θ is assumed to be identical for all periods and products, Algorithm 1 returns the smallest θ min over all periods. From a practical standpoint, some products can be filtered out of Algorithm 1 if there is no uncertainty on the product, or if the uncertainty on the product does not need to be covered.

If Algorithm 1 is run when all new qualifications are started at t = 0, then an ideal value of θ is computed. This is an estimate of the largest value of θ for which the demand uncertainty can be covered in the work center. Reporting this value is interesting for capacity planners to assess the robustness of the work center against an ideal situation.

Computational study

The computational study is performed to answer the following questions: What is the price of uncertainty? Is it risky to use the set of qualifications determined by considering only the nominal demand? Is the robust optimization problem difficult to solve?

In Section 5.1, the instances used for the computational study and generated from industrial data are described. For confidentially purposes, raw values by product, by operation, by product family and by machine of parameters are not provided. Instead, means, minimums, maximums and standard deviations are presented. In Section 5.2, the design of experiments is presented, and the numerical results in Section 5.3. We show that the price of uncertainty, defined by comparing the number of qualifications determined for the robust optimization problem and for the deterministic optimization problem, when the demand is fully known (perfect handsight), is actually very small. Moreover, in a large number of experiments, the robustness of the set of qualifications determined by solving the deterministic optimization problem with the nominal demand is far from the robustness of the set of qualifications determined by solving the robust optimization problem whereas both qualification matrices have about the same number of qualifications. In addition, we show that only considering the nominal demand can lead to a large number of capacity constraint violations. The computational study highlights that selecting the right qualifications is more important for robustness than the number of qualifications.

Instance generation

In this section, the instances used in the computational study are described, and can be used to further generate instances from real industrial data. Note that the instances can be made available to the reader by contacting one of the authors. Work center: The computational study is performed by using industrial data from a semiconductor factory located at Crolles, France. Two critical work centers, work center A and work center B, of the factory are considered. Work center A has M = 20 machines. Work center B has M = 30 machines. Demand: A horizon of 7 periods, i.e. T = 7, is considered. Each period corresponds to one month. The nominal demand by product is given by internal forecasts for each period of the horizon. Exact demand values are not provided for confidentiality reasons. Instead, Table 3 illustrates the number of products with a non nul demand by period and the Coefficient of Variability (CV) of the demand by period. On the horizon, there are in total 238 products, i.e. P = 238. For work center A, these 238 products lead to 1,208 operations, i.e. R = 1, 208. For work center B, these 238 products lead to 401 operations. There is no uncertainty on the demand for the first month. Production capacities: For work center A, u max t,m = 0.95 ∀t, ∀m in the industrial data. Consider a given period t. For work center B, the mean of u max t,m is equal to 0.80, the minimum of u max t,m to 0.63, the maximum of u max t,m to 0.87, and the standard deviation of u max t,m to 0.079. Both work centers do not have the same values for u max t,m because machine types are completely different. Note that u max t,m is constant from one period to another. Similarly, the production capacity by machine c t,m is constant from one period to another, but is different from one machine to another. This is mainly because machines are non-identical and are of different ages and generations. Values for c t,m are given based on the length of period t. For work center A, the mean of c t,m is equal to 59% of the length of the period, the minimum of c t,m to 44%, the maximum to 66%, and the standard deviation to 6%. For work center B, the mean of c t,m is set to 75%, the minimum to 36%, the maximum to 85%, and the standard deviation to 9%. c t,m is not equal to 100% because machines have capacity losses, e.g. due to maintenance operations, engineering operations, setup times. Re-entrant flow factors: For work center A, the re-entrant flow factors vary between 14 and 72, with a mean of 41.2 and a standard deviation of 11.0. For work center B, the re-entrant flow factors vary between 1 and 28, with a mean of 16.0 and a standard deviation of 4.3. Product families: There are three product families, i.e. F = 3. Each product belongs to exactly one product family. The first product family contains 120 products. The second product family contains 64 products. The third product family contains 54 products. Qualification matrix: The initial set of qualifications is partially initialized, in particular because some machines are already qualified for existing operations. Consider work center A. The mean number of qualified machines by operation is equal to 4.2, and the standard deviation to 2.0. The minimum, respectively maximum, number of qualified machines for an operation is equal to 1, respectively 13. The mean number of qualified operations by machine is equal to 251.3, and the standard deviation to 188.0. The minimum, respectively maximum, number of qualified operations for a machine is equal to 25, respectively 645. Note that some operations cannot be qualified on some machines due to technological restrictions. In total, 2,843 new qualifications are possible in work center A. Consider work center B. The mean number of qualified machines by operation is equal to 3.5, and the standard deviation to 1.6. The minimum, respectively maximum, number of qualified machines for an operation is equal to 1, respectively 6. The mean number of qualified operations by machine is equal to 48.0, and the standard deviation to 43.8. The minimum, respectively maximum, number of qualified operations for a machine is equal to 0, respectively 130. Note that some operations cannot be qualified on some machines due to technological restrictions. In total, 1,266 new qualifications are possible in work center B. Some machines have no qualified operations because they are being started up. Qualification costs: We could not access to the qualification costs. Therefore, we assume that all qualification costs are identical and equal to one. This is a common assumption made by capacity planners in practice. Hence, in the computational study, the number of qualifications to perform must be minimized. Qualification lead times: Qualification lead times are rough estimates of the lead times to perform the qualification procedures. Qualification lead times vary between several days and two months. For qualification lead times that are smaller than 2 weeks, they are set to 0 because the considered period in the computational study is one month. Consider work center A. The minimum lead time for all operations and machines is equal to 0 period, the mean to 1.6, the standard deviation to 0.8, and the maximum to 2. Consider work center B. The minimum lead time for all operations and machines is equal to 0 period, the mean to 1.1, the standard deviation to 0.4, and the maximum to 2. Throughput rates: Throughput rates strongly vary from one machine to another and from one operation to another. Consider work center A. The minimum throughput rate for all operations and machines is equal to 11.4 wafers per hour, the mean to 221.6, the standard deviation to 126.7, and the maximum to 527.8. Consider work center B. The minimum throughput rate for all operations and machines is equal to 6.8 wafers per hour, the mean to 48.0, the standard deviation to 13.9, and the maximum to 83.3.

Design of experiments

Different values of θ are studied: 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7. Values of 0.2 and 0.3 for θ is not unusual even for early periods of the horizon for high mix factories. Larger values of θ are not considered because the robust optimization problem becomes infeasible from θ = 0.770 for work center A (see Section 5.3). In addition, the robust optimization problem becomes infeasible from θ = 0.294 for work center B. The budget of uncertainty Γ t,f is set to p|α p,f =1 d t,p , ∀t, ∀f . The discount factor δ t is set to 1 ∀t in numerical experiments. This means that there are no incentives on performing qualifications as soon as possible or as late as possible. In Algorithm 1, we consider that θ 0 t,p = 1 ∀t, ∀p.

In the experiments, MCQCP is solved once. The robustness of the optimized set of qualifications is evaluated with Algorithm 1. Then, for each possible value of θ, MCRQCP is solved. For each value of θ, 3,600 demand scenarios are generated to evaluate the capacity constraint violations if the nominal set of qualifications was considered, and the price of uncertainty. Because the true distribution of the demand is unknown and the demand between products is correlated, scenarios are randomly generated by using a linear program. The linear program, described in Appendix A. The linear program generates for a given θ a scenario on the demand by product and by period for a given demand level η t,f by product family and by period. In the experiments, it is assumed that η t,f is equal to the nominal demand by product family.

For the sake of presentation, in the remainder of the computational study, the set of qualifications determined by solving MCQCP for the nominal demand are called nominal qualifications, the set of qualifications determined by solving MCQCP for the perfect handsight demand scenario, perfect handsight qualifications, and the set of qualifications determined by solving MCRQCP, robust qualifications.

Note that the robust and nominal qualifications are not compared in a rolling horizon in the computational study, i.e. where qualifications could be updated at each period after demand realizations for the following reasons: (1) It is difficult to known the final practical decision when capacity constraint violations occur; (2) Qualification decisions must be anticipated due to long qualification processes; (3) Comparing the robust and nominal set of qualifications is possible and fair because both are computed from "static" optimization models.

Numerical results

Mathematical models and Algorithm 1 are implemented in Java 8 on a computer with an Intel Xeon CPU W3530 running at 2.80GHz with 8 threads and 12GB of RAM. Mathematical models are solved by using the solver IBM ILOG CPLEX 12.9 with default parameters. A computational time limit of one hour is given to the solver, is set to 0.0001 in Algorithm 1.

Section 5.3.1 answers the question "What is the price of uncertainty?", Section 5.3.2 the question "It the robust optimization problem difficult to solve?", and Section 5.3.3, the question "Is it risky to use the set of qualifications determined by considering only the nominal demand?"

What is the price of uncertainty?

The Price of Uncertainty (PoU) is computed by comparing the number of robust qualifications and the number of perfect handsight qualifications. [START_REF] Gorissen | A practical guide to robust optimization[END_REF] argue that a low mean PoU and standard deviation indicate a good robust solution. Table 4 shows the mean PoU, its standard deviation (std.) and its maximum value for each θ. Note that as for θ > 0.294, MCRQCP is infeasible for work center B, PoU is not presented.

Consider work center A. The mean PoU varies between 1.08 qualifications on average for θ = 0.1, with a standard deviation of 0.30, and 31.99 qualifications on average for θ = 0.7 with a standard deviation of 2.03. Note that the increase of PoU when θ increases is mainly due to the fact that the number of robust qualifications increases (see Figure 3). The standard deviation of PoU is small with respect to the mean PoU. To better put into perspective, the meaning of about 30 qualifications, consider θ = 0.7. In the worst case, PoU is equal to 35. Recall that the number of machines in work center A is equal to 20. In other words, to cover the demand uncertainty, it is required to add 35 20 = 1.75 qualifications on average to each machine, each having a few hundred qualifications on average, which seems acceptable in practice. Therefore, robust qualifications for work center A appear to be good solutions. In addition, a small number of additional qualifications, in the worst case 35, is required to cover the demand uncertainty. This is small compared to the 2,843 possible new qualifications. This suggests that it is possible to be robust by performing the right qualifications. Robust qualifications are also relevant because they can avoid capacity constraint violations contrary to nominal qualifications (see Section 5.3.3).

Similar observations can be observed for work center B (see Table 4 and Figure 4). The maximum PoU varies between 5 and 19 qualifications. Similarly to work center A, to better put into perspective the meaning of a PoU of 19 qualifications, recall that the number of machines in work center B is equal to 30. With respect to the perfect handsight qualifications, to cover the demand uncertainty, it is required to add 19 30 = 0.63 qualifications on average to each machine, each having a few tens of qualifications on average. This is also small compared to the 1,266 possible new qualifications. This again suggests that it is possible to be robust by performing the right qualifications.

Implementing perfect handsight qualifications is impossible because it is impossible to know in advance the demand realizations. A more practical price of uncertainty can be computed by comparing the number of robust qualifications and nominal qualifications. We found that the actual price of uncertainty is close to the PoU presented in Table 4. This is because, for both work centers, the number of nominal qualifications is equal to 4 and the mean number of perfect handsight scenarios is also close to 4 (see Figures 3 and4). Now consider the case where θ = θ max , where θ max is the largest possible value of θ for the considered work center. It can be computed by running Algorithm 1 when all new qualifications are started at t = 0. For work center A, this gives θ max = 0.77. When MCRQCP is solved for θ = θ max , 96 new qualifications are required (the set of new qualifications is optimal). 96 2,843 ×100 = 3.37% of all new possible qualifications are required to reach the same robustness than the one when all qualifications are performed. For work center B, θ max = 0.294. When MCRQCP is solved for θ = θ max , 135 new qualifications are required (optimality gap of 25.0% after 3,600 seconds). 135 1,266 ×100 = 10.6% of all new possible qualifications are required to reach the same robustness than the one when all qualifications are performed. In the best case, 135 -0.25 × 135 = 102 new qualifications are required, which corresponds to 8.05% of all possible new qualifications. This further suggests that it is possible to be robust by performing a limited number of qualifications. In other words, it can be ineffective to add many qualifications, if they are irrelevant. Similar observations can be found in other contributions on flexibility, e.g. on the long-chain and closed-chain principles [START_REF] Jordan | Principles on the benefits of manufacturing process flexibility[END_REF][START_REF] Chou | Design for process flexibility: Efficiency of the long chain and sparse structure[END_REF]. Thus, relevant qualifications must be carefully optimized and planned to immunize a work center against demand uncertainty. One of the reasons why PoU is small is that qualification costs are assumed identical in the computational study, which is in fact a common assumption in practice. PoU could potentially be larger if qualification cost profiles are different from one machine to another and from one operation to another. Nevertheless, PoU is not necessarily expected to be significantly larger since new qualifications must be paid in the nominal, perfect handsight and robust cases for the following reasons: (1) Qualifications are made for new operations or new machines, or existing operations that have never been qualified on existing machines and (2) A ramp-up demand for a product, even uncertain, implies adding new qualifications to machines to increase product capacity and balance the workload between the machines. If new qualifications are not performed, then it is impossible to satisfy the demand, and both MCRQCP and MCQCP are infeasible. If qualification cost profiles are very different, it may be possible to keep a small PoU by performing a lot of inexpensive qualifications and avoid performing expensive qualifications whenever possible. Finally, PoU is also small because of product cannibalization that limits the overall demand of products. Now assume that the manufacturer faces an extreme case where too many qualifications must be performed to cover the demand uncertainty with respect to the number of nominal qualifications. This information is still valuable for capacity planners because they will have to refine plausibility limits to limit additional outsourcing and machine purchasing costs. In this situation, MRBQP is relevant to help refining plausibility limits.

Finally, from a practical standpoint, as both work centers are located in the same factory, covering the demand uncertainty for θ larger than 0.3 in work center A is probably unnecessary as θ max is equal to 0.294 for work center B.

Is the robust optimization problem difficult to solve?

Consider work center A. For all values of θ, a set of optimal robust qualifications is determined. However, determining optimal robust qualifications is much more time consuming than determining optimal nominal qualifications (about 3 seconds). Similarly, determining optimal perfect handsight qualifications requires between 2 and 6 seconds in most cases, and never exceeds 13 seconds. Determining optimal robust qualifications requires between 46 seconds for θ = 0.1 and 1,551 seconds for θ = 0.7 (see Figure 5). For θ max , the optimal set of robust qualifications is determined in 656 seconds. It is also worth mentioning that all optimal nominal qualifications are determined at the root node by IBM ILOG CPLEX. Except for θ = 0.4, 0.5, 0.7 and θ = θ max , all robust qualifications are also determined at the root node by IBM ILOG CPLEX. This can be explained by the fact that modern solvers such as IBM ILOG CPLEX embed advanced preprocessing, probing, heuristic and cutting plane routines that are used to strengthen the linear relaxation of mixed integer linear problems (see e.g. [START_REF] Savelsbergh | Preprocessing and probing techniques for mixed integer programming problems[END_REF][START_REF] Atamtürk | Conflict graphs in solving integer programming problems[END_REF] and quickly to determine good solutions. It can also be observed that it is faster to get the optimal robust qualifications for θ = 0.6 than for θ = 0.5.

For work center B, determining nominal qualifications takes about 1 second, while, similarly to work center A, determining optimal robust qualifications is more difficult. For θ = 0.1, optimal robust qualifications are determined in 85 seconds (see Figure 5), and in 3,472 seconds for θ = 0.2. Branching in IBM ILOG CPLEX is required for both θ = 0.1 and θ = 0.2. 5.3.3. Is it risky to use the set of qualifications determined by only considering the nominal demand?

The numerical experiments show that it can be risky to implement nominal qualifications because it can lead to capacity constraint violations, which are computed with the following procedure:

1. A demand scenario is generated with the linear program in Appendix A. 2. Then, for the set of nominal qualifications and the generated demand, the Total Overtime (OT) is minimized with the linear program (B.1)-(B.6) in Appendix B. 3. If OT > 0, then there is at least one capacity constraint violation for the considered scenario.

In this case, to put into perspective what a positive overtime means, in particular in terms of machine utilization rates, we solve the nonlinear utilization balancing optimization problem proposed in [START_REF] Rowshannahad | Capacitated qualification management in semiconductor manufacturing[END_REF]. This avoids the problem where the total overtime for a period is set to a specific machine whereas, in practice, it would be balanbced with similarly qualified machines. The utilization balancing optimization problem is parameterized by an utilization balancing exponent γ, which is set to 20 in this paper.

For scenario i, the procedure enables us to determine what would be the utilization rate U i t,m of machine m at period t for a given demand by product and a given set of qualifications (here the nominal qualifications) if there is a capacity constraint violation. If U i t,m > u max t,m , then there is a capacity constraint violation for scenario i. Repeating this procedure for the 3,600 scenarios enables us to estimate the capacity constraint violations if only nominal qualifications were implemented.

Table 5 shows capacity constraint violations. Column "A" corresponds to the percentage of scenarios where there is at least one capacity constraint violation, Column "B" to the number of capacity constraint violations. Mathematically, the mean number of capacity constraint violations is computed as follows:

1 3,600×T ×M ×I 3,600 i=1 1(U i t,m -u max t,m )
, where 1(x) = 1 if x > 0, and 0 otherwise and I = 3, 600. The maximum (max.) number of capacity constraint violations is computed as follows: max i ( t,m 1(U i t,m -u max t,m )). Column "C" quantifies capacity constraint violations when there is at least one capacity constraint violation. Mathematically, the mean capacity constraint violation is computed as follows:

1 3,600×T ×M ×I 3,600 i=1 ( t,m max(0, U i t,m -u max t,m )), and the maximum capacity constraint violation is computed as follows max i,t,m max(0, U i t,m -u max t,m ). Columns "B" and "C" are computed only if there is at least one capacity constraint violation. 5), i.e. a relatively small number of scenarios. In addition, the number of capacity constraint violations is relatively small. In the worst case, 8 out of 140 (T = 7 and M = 20) capacity constraints are violated, capacity constraint violations are not very large, on average 0.004 and at most 0.010. This means that if the maximum utilization rate of a machine was set to 0.95, then on average, its real utilization rate would be equal to 0.9504, at most 0.96. Therefore, for θ = 0.1, using the nominal qualifications is probably acceptable. For θ = 0.2, 15.64% of the scenarios have a capacity constraint violation, which is significantly larger than for θ = 0.1. On average, the capacity constraint violation is equal to 0.010, and in the worst case to 0.029, which starts to be appreciable. For larger values of θ, using nominal qualifications is more risky. For instance, consider θ = 0.4 where 30.78% of the scenarios have at least one capacity constraint violation (see Table 5). In the worst case, 13 out of 140 capacity constraints are violated. In addition, the largest capacity constraint violation is equal to 0.066. This means, that if the maximum utilization rate of a machine was set to 0.95, then its real utilization rate would be equal to 1.003. The same observations can be made for larger values of θ. Utilization rates near 1.0 are not sustainable in terms of service levels. This is due to the fact that the cycle time increases almost exponentially with the utilization rate (queuing theory) and due to production variability [START_REF] Hopp | Factory physics[END_REF]. In other words, even small capacity constraint violations should be avoided.

Capacity constraint violations are more critical for work center B than for work center A. For θ = 0.1, 15.56% of the scenarios lead to at least one capacity constraint violation, and for θ = 0.2, 44.28% of the scenarios. In the worst case, there are 11 capacity constraint violations for θ = 0.1 and 22 capacity constraint violations for θ = 0.2. For θ = 0.1, the mean capacity constraint violation is equal to 0.010. u max t,m is set to low values (compared to work center A) because it is known that, in the industrial context, small increases of utilization rates can lead to much larger cycle times due to production variability.

Using nominal qualifications can lead to capacity constraint violations because nominal qualifications are not robust against demand uncertainty, and are in fact much less robust than robust qualifications. For work center A, Algorithm 1 for the nominal qualifications gives θ = 0.043, and θ = 0.024 for work center B. With a limited number of additional qualifications, robust qualifications lead to a much better robustness (see Section 5.3.1). Consider work center A and θ = 0.2, 7 robust qualifications are required instead of 4 nominal qualifications to avoid capacity constraint violations in 15.64% of the scenarios. For θ = 0.3, 9 robust qualifications are sufficient to avoid capacity constraint violations in 26.19% of the scenarios. Similar observations can be made for work center B. For instance, for θ = 0.1, 7 robust qualifications are required instead 4 of nominal qualifications to avoid capacity constraint violations in 15.56% of the scenarios. It is worth mentioning that robust qualifications are more robust against demand uncertainty because more qualifications are performed. Nevertheless, even by adding a large number of qualifications, the nominal qualifications are still outperformed by the robust qualifications in terms of demand uncertainty coverage. Let us consider the case where α-flexibility designs are enforced when nominal qualifications are determined by the optimization model ( 1)-( 7). An α-flexibility design enforces that at least α machines must be qualified by operation. For operations where it is not possible to have α qualified machines, the number of largest number of qualified machines is enforced. An α-flexibility design is enforced by adding the two following constraints: (1') t t =1|t +l t ,r,m ≤t OQ t ,r,m ≤ 1 ∀t, ∀r, ∀m, (2') m t t =1|t +l t ,r,m ≤t OQ t ,r,m ≥ min(α, α ) -α ∀t, ∀r. α is the the number of qualifiable and qualified couples (operation, machine) for a given period, and α is the number of qualified couples (operation, machine) for a given period. Constraints (1') are required, otherwise Constraints (2') could be satisfied by performing the new qualification at different periods. Table 6 shows that enforcing α-flexibility designs for the nominal qualifications does not lead to a better robustness against demand uncertainty than the robust qualifications even though many qualifications are performed. This is because there are many different ways to enforce an α-flexibility design. This reinforces the idea that if qualifications are not optimized, then even many qualifications may not be effective to tackle demand uncertainty. Practical consequences of capacity constraint violations are lower service levels, larger cycle time and larger inventory holding costs. Due to to capacity constraint violations, the number of products in the factory would have to be decreased so that the real utilization rates of machines violating their capacity constraint in the factory is at least lower than 1.0, and ideally lower than u max to control the cycle times. This can severely affect deliveries and the production objectives of the factory.

In practice, a method to deal with uncertainty is to continuously updating nominal qualifications each time the demand is updated. This should be avoided. This is because, as mentioned in Section 2.1, this does not guarantee to find feasible nominal qualifications because the qualification process may sometimes take several weeks or months to validate the quality and the yield of the operation. As the demand by product for the early months on the horizon is also subject to uncertainty, determining and planning robust qualifications is preferable for the whole horizon.

It is worth observing that, if θ is not adequately selected, there may also exist multiple sets of robust qualifications with the same number of qualifications. However, some sets of robust qualifications may actually be better to cover a larger demand uncertainty than other sets of robust qualifications, which is not captured by the robust optimization model because it only seeks to immunize the work center against the specified uncertainty. This is why Algorithm 1 is relevant to identify the most robust set of qualifications among all robust sets of qualifications. These observations are consistent with other observations in the literature: There may exist multiple robust solutions to an optimization problem. Although these robust solutions have the same worst-case objective value, they can have different performances for the nominal scenario [START_REF] Iancu | Pareto efficiency in robust optimization[END_REF][START_REF] Gorissen | A practical guide to robust optimization[END_REF][START_REF] De Ruiter | The impact of the existence of multiple adjustable robust solutions[END_REF][START_REF] De Ruiter | Robust optimization of uncertain multistage inventory systems with inexact data in decision rules[END_REF][START_REF] Yanıkoglu | A survey of adjustable robust optimization[END_REF].

Practical use of optimization models

Determining qualification decisions

A straightforward use of the robust optimization model ( 12)-( 14) is to determine new qualifications to perform to satisfy the demand while respecting capacity constraints and covering the demand uncertainty.

Further improving manufacturing performances

As illustrated on the industrial data in Section 5, a small number of qualifications among several hundreds of new qualifications is sufficient to cover the demand uncertainty. Consequently, it is likely that there are two different sets of robust qualifications that cover the demand uncertainty but lead to different performances, for instance in terms of utilization balance of the machines or production variability. It is necessary to distinguish them to further improve manufacturing performances. Differentiating identical sets of robust qualifications in terms of number of qualifications can be done by populating the solution pool after determining the minimum number of qualifications to perform. Modern solvers such as IBM ILOG CPLEX provide this functionality:

1. Two sets of robust qualifications may not be identical in terms of robustness. Algorithm 1 can be used to identify the most robust set of qualifications. 2. Two sets of robust qualifications may also be different in terms of real utilization rates although they all satisfy capacity constraints. [START_REF] Johnzén | Flexibility measures for qualification management in wafer fabs[END_REF] and [START_REF] Rowshannahad | Capacitated qualification management in semiconductor manufacturing[END_REF] propose a "time flexibility measure" to evaluate sets of new qualifications in terms of total utilization rate and utilization balance of the machines. This flexibility measure is interesting as maximizing the utilization balance contributes to further control and reduce cycle times. However, [START_REF] Johnzén | Flexibility measures for qualification management in wafer fabs[END_REF] and [START_REF] Rowshannahad | Capacitated qualification management in semiconductor manufacturing[END_REF] do not consider that demand uncertainty. Their model need to be robustifieds. 3. Robust qualifications can be differentiated in terms of production variability as a large production variability contributes to significantly increase cycle times [START_REF] Hopp | Factory physics[END_REF]. In semiconductor factories, partly due to re-entrant flow, it is unlikely that products arrive continuously in work centers. Work centers are often subject to large Work-In-Process (WIP) peaks leading to congestion. To better capture this phenomenon, [START_REF] Johnzén | Flexibility measures for qualification management in wafer fabs[END_REF] propose "a toolset" flexibility measure that captures the fact that operations with large demands must be more qualified than operations with low demands. [START_REF] Pianne | Ideal and potential flexibility measures for qualification management in semiconductor manufacturing[END_REF] argue that qualified process times should be balanced between machines in the work center. A machine should not be overqualified at the expense of other machines. This is because machines with few qualifications must process almost all their qualified products every qualified to meet the optimized utilization balance, which is difficult due to production variability. Associated flexibility measures are proposed in [START_REF] Pianne | Ideal and potential flexibility measures for qualification management in semiconductor manufacturing[END_REF]. They can be seen as ways to measure the quality of the balancing of the qualified process times, and not the quality of the utilization balance of the machines. 4. The principle of large closed chains or long chains can also be used to differentiate sets of qualifications. If one set of qualifications creates more closed chains or larger closed chains between machines and operations than other sets of qualifications, it is very likely that the former will deal better with WIP peaks than the latter [START_REF] Jordan | Principles on the benefits of manufacturing process flexibility[END_REF][START_REF] Graves | Process flexibility in supply chains[END_REF]. 5. Another straightforward way of differentiating sets of qualifications consists in enforcing αflexibility designs. However, note that enforcing α-flexibility designs without optimizing a criterion that helps to tackle WIP peaks, such as flexibility measures, may not necessarily lead to better performances (see Section 5.3.3).

Exploiting dual variables of robust reformulation

Bertsimas and Thiele (2006b) report that dual variables correspond to the sensitivity of the objective function to changes in parameters of the budget uncertainty set for an inventory management problem. Similarly, dual variables of the robust optimization model, namely y min t,m,p , y max t,m,p , y gamma t,m,f , can also be exploited: is the sensitivity of the number of qualifications to perform to changes in the parameter Γ t,f . In other words, if y gamma t,m,f decreases, Γ t,f indicates the potential reduction of the number of qualifications.

• y min t,
Exploiting the values of dual variables is particularly relevant from an industrial standpoint to identify if the demand uncertainty on some products or product families is very expensive in terms of number of qualifications. Reporting the values of dual variables can be used by capacity planners to refine the uncertainty set, i.e. by defining a smaller uncertainty set, and initiate a discussion with the departments in charge of defining future demands in the case where the number of qualifications to perform is overwhelming. Capacity planners can also initiate a discussion with the departments in charge of defining future demands that the demand uncertainty on some products or product families is not constraining for the production system. The departments can therefore consider new future potential product mixes, i.e. by defining a larger uncertainty set, that would have never been initially considered.

On infeasibilities

The optimization problems can be infeasible (see Section 5.3). For instance, this can be caused by large qualification lead times and too small production capacities to cover the demand uncertainty. Determining that optimization problems are infeasible is also valuable in practice.

If the nominal optimization problem is infeasible, it indicates to capacity planners that the demand must be changed. An option would then consist in adapting the product mix and production quantities, either by producing more products during some months or postponing production to make the best use of the installed process flexibility and still meet the demand [START_REF] Yayla-Küllü | Product line flexibility for agile and adaptable operations[END_REF]. However, it is difficult to anticipate how would be the new demand as it depends on different stakeholders (e.g. capacity planning, demand planning) within a company. For instance, if the nominal optimization problem is infeasible, the demand for products that generate a large workload at the work center can be decreased while the demand for products that generate a lesser workload can be increased. In this case, the total number of product units made may not decrease, backlogging costs may be acceptable, but lost sales may be incurred on critical products.

If both MCQCPLT or MCRQCPLT cannot be solved because capacity constraints cannot be respected, it is also possible to solve a utilization balancing problem where the demand is described by the uncertainty D t to highlight critical machines, i.e. machines for which U t,m > u max t,m . We refer the reader to [START_REF] Rowshannahad | Capacitated qualification management in semiconductor manufacturing[END_REF] and [START_REF] Christ | An iterated min-max procedure for practical workload balancing on non-identical parallel machines in manufacturing systems[END_REF] for existing utilization balancing approaches. These approaches also need to be robustified. In a decision support system, systematically solving a robust utilization balancing problem is relevant to either identify infeasibilities or most loaded and critical machines.

Conclusions and perspectives

In this paper, we first proposed a new mixed-integer linear programming mathematical model for a tactical qualification management problem, which is shown to be NP-Hard, when the demand is deterministic. We showed that the studied problem is NP-Hard. Second, we motivated the choice of robust optimization when the demand is uncertain, in particular for high mix factories. We proposed an uncertainty set based on the budget of uncertainty to describe product cannibalization and cover the demand uncertainty. Third, we proposed a new robust reformulation of the deterministic model when the demand is described by product cannibalization. Fourth, we proposed a linear program and a binary search approach to characterize the robustness of a set of qualifications when the demand is uncertain. Fifth, we performed a computational study by using industrial data from a high mix semiconductor manufacturer. In particular, we showed that, (1) The price of uncertainty is acceptable, often less than a few additional qualifications for each machine, (2) It is possible to achieve the same level of robustness as the case where all new qualifications are performed by only performing a restricted number of relevant qualifications, (3) Depending on the forecast uncertainty and the work center, the robust optimization problem can be difficult to solve, and (4) Using the nominal set of qualifications can lead to significant capacity constraint violations, although it can be used for some work centers when the forecast uncertainty is small. Finally, practical applications and implications of the developed models are discussed.

We believe the following perspectives are worth investigating in the future. First, other parameters can also be subject to uncertainty, e.g. production capacities, throughput rates of operations on machines, qualification costs and lead times. Studying the relevance and effect of uncertainty on these parameters can be valuable. Second, a large number of qualifications can be difficult to maintain at an operational level. Including disqualification decisions, e.g. constraining the number of qualifications by machines, or constraining the total number of qualifications in each period, could be relevant. Third, extending the static robust reformulation to adjustable robust reformulations may be valuable to further reduce qualification costs. Fourth, for work centers where the numbers of operations and machines are large, efficient solution approaches can be valuable. An option consists in using a cutting-plane solution approach with lazy constraints as proposed by [START_REF] Bertsimas | Reformulation versus cutting-planes for robust optimization[END_REF]. This might be a viable approach as the computational time required to solve MCQCP is small. Fifth, as there may exist several sets of robust qualifications in terms of number of qualifications given an immunization level, it would be interesting to use additional objective functions to select the most set of robust qualifications. This leads to considering a multi-objective optimization approach for the studied problem. Sixth, other solution approaches for MRBQP can be considered. Iterated max-min approaches are probably relevant not to restrict to the same value of θ for all products and periods. Seventh, studying the effect of different qualification cost profiles by machine or by machine and time dependent qualification decisions on the price of uncertainty can be interesting. Finally, as the ability of qualifications to cover the uncertainty on the demand strongly depends on the machines in the work center [START_REF] Hopp | Factory physics[END_REF], considering the investment decisions in terms of machines could also be investigated to cover the uncertainty on the demand. In the same vein than in [START_REF] Xu | Configuration of flexibility strategies under supply uncertainty[END_REF] with pricing flexibility for supply uncertainty or in [START_REF] Li | Supply chain coordination through capacity reservation contract and quantity flexibility contract[END_REF] with capacity reservation and quantity flexibility contracts, considering other types of flexibility should lead to interesting research avenues. Belgium, Spain, Greece, Sweden, Italy, Ireland, Poland, Hungary, Portugal, Denmark, Finland, Luxembourg, Norway, Turkey.

  (a) 2-chain design. (b) Full flexibility design.

Figure 1 :

 1 Figure 1: Visual comparison of different flexibility designs.

Figure 2 :

 2 Figure 2: Illustrative example on the demand profile for one product.

  (a) Number of robust qualifications by θ. (b) Mean number of qualifications by θ (PH).

Figure 3 :

 3 Figure 3: Work center A. Number of qualifications by θ.

  (a) Number of robust qualifications by θ. (b) Mean number of qualifications by θ (PH).

Figure 4 :

 4 Figure 4: Work center B. Number of qualifications by θ.

  (a) Work center A. (b) Work center B.

Figure 5 :

 5 Figure 5: Computational time (in seconds) required to determine the set of robust qualifications by θ.

  , is used. Let us introduce the new notations below: New parameters: F : Number of product families, d t,p : Nominal demand for product p at period t, d t,p ≤ d t,p : Maximum deviation from nominal demand for product p at period t, α p,f : Is equal to 1 if product p belongs to product family f , and 0 otherwise, Γ t,f : Budget of uncertainty for product family f at period t. The demand d t,p is assumed to be an uncertain parameter that takes values as follows: d t,p ∈ [d t,p -d t,p , d t,p + d t,p ] ∀t, ∀p. d t,p -d t,p and d t,p + d t,p

  . In this paper, the uncertainty set D t is compact: The uncertainty set D t is bounded, because 0 ≤ d t,p ≤ d t,p + d t,p ∀t, ∀p, and is closed because D t consists of a set of closed half spaces described by linear inequalities. However, the uncertainty is not row-wise because the uncertain parameter for period t, i.e. d t,p , is found in the capacity constraint of each machine. The uncertainty would be row-wise if the demand for a product also depended on the machine, which is impossible. Investigating adjustable robust reformulation is therefore interesting but left for future research.By considering the uncertainty set D t to model the demand uncertainty, capacity constraints become in a static reformulation: rf p,r d t,p )W IP r,m

	tp r,m	≤ c t,m u max t,m	∀t, ∀m, ∀d ∈ D t

r

( p

Table 1 :

 1 Dual variables associated to constraints in the uncertainty set Dt for a capacity constraint (2).

	Constraints Dual variables
	(9)	y min p
	(10) (11)	y max p y gamma f

  By deriving the robust counterpart for each capacity constraint and each time period and indexing the dual variables by period t and machine m, the overall robust optimization problem is:

	3.4.3. Robust optimization model
	min		δ t cq r,m OQ t,r,m	(12)
		t,r,m
	s. t.	(4) -(7)
			(-(d t,p -d t,p )y min t,m,p ) +	(Γ t,f y gamma t,m,f )
		p							f
		+		((d t,p + d t,p )y max t,m,p ) ≤ c t,m u max t,m	∀t, ∀m	(13)
			p
		-y min t,m,p + y max t,m,p
		+	f |α p,f =1	y gamma t,m,f	≥	r	rf p,r W IP t,r,m tp r,m	∀t, ∀m, ∀p	(14)
		y min t,m,p , y max t,m,p ≥ 0		∀t, ∀m, ∀p
									tp r,m	r,m	∀p
		y min p	, y max p	≥ 0		∀p
		y gamma f	≥ 0		∀f

m W IP t,r,m = 1 ∀t, ∀r | p rf p,r (d t,p + d t,p ) > 0 (15)

Table 2 :

 2 Comparison of the number of variables and constraints between MCQCP and MCRQCP. P = 238, R = 1208, F = 3, M = 20, T = 7.

		Optimization problem
	Number of	MCQCP MCRQCP Increase(%)
	Continuous variables 169,120	202,860	16.6
	Binary variables	169,120	169,120	0.0
	Constraints	685,076	785,456	12.8

4. Characterizing the robustness of a set of qualifications

4.1. Motivation

Determining intuitively relevant values for d t,p can be difficult.

Table 3 :

 3 Nominal demand by month.

				Month
	1	2	3	4	5	6	7
	CV 2.78 1.97 2.88 2.29 2.06 3.58 3.09

Table 4 :

 4 Price of Uncertainty (PoU).

		Work center A	Work center B
	θ	Mean Std. Max. Mean Std. Max.
	0.1 1.08 0.30	2	2.88 0.41	5
	0.2 3.10 0.63	4	16.66 0.94	19
	0.3 5.05 0.85	7	-	-	-
	0.4 7.96 1.11	10	-	-	-
	0.5 12.70 1.41	15	-	-	-
	0.6 18.44 1.68	21	-	-	-
	0.7 31.99 2.03	35	-	-	-

Table 5 :

 5 Capacity constraint violations.

		Work center A		Work center B
	A	B		C	A	B		C
	θ	Mean Max. Mean Max.		Mean Max. Mean Max.
	0.1 0.72% 0.058	8	0.004 0.010 15.56% 7.096	11	0.010 0.029
	0.2 15.64% 1.241	8	0.010 0.029 44.28% 8.947	22	0.014 0.048
	0.3 26.19% 2.080	12	0.015 0.046	-	-	-	-	-
	0.4 30.78% 2.670	13	0.020 0.066	-	-	-	-	-
	0.5 38.39% 3.693	24	0.024 0.086	-	-	-	-	-
	0.6 45.31% 5.058	36	0.027 0.106	-	-	-	-	-
	0.7 57.83% 7.704	51	0.030 0.126	-	-	-	-	-
	Consider work center A, θ = 0.1, 0.72% of the scenarios have a capacity constraint violation (see
	Table

Table 6 :

 6 Number of qualifications (NQ) and robustness (θ) of nominal qualifications when an α-flexibility design is enforced.

		Work center A Work center B
	α NQ	θ	NQ	θ
	1	4	0.043	4	0.024
	2	84	0.043	14	0.021
	3 251	0.061	77	0.012
	4 611	0.071	224	0.087
	5 1119	0.152	394	0.140

  m,p is the sensitivity of the number of qualifications to perform to changes in the parameter d t,p -d t,p . In other words, if d t,p -d t,p increases, y min t,m,p indicates the potential reduction of the number of qualifications. • y max t,m,p is the sensitivity of the number of qualifications to perform to changes in the parameter d t,p + d t,p . In other words, if d t,p + d t,p decreases, y max t,m,p indicates the potential reduction of the number of qualifications.

	• y gamma t,m,f

Acknowledgements

This work has been partly funded by the ANRT (Association Nationale de la Recherche et de la Technologie) through CIFRE funds (PhD number 2017/0750) and a cooperation contract between STMicroelectronics and Mines Saint-Etienne. This work has also received funding from the Electronic Component Systems for European Leadership Joint Undertaking under grant agreement No 737459 (project Productive4.0). This Joint Undertaking received support from the European Unions Horizon 2020 research and innovation program and Germany, Austria, France, Czech Republic, Netherlands,

Appendix A. Linear programming for scenario generation

The linear program (A.1)-(A.4) consists in simulating a (perfect handsight) scenario on the demand from a nominal demand and the uncertainty parameters defined in the uncertainty set D t . The w parameters are weights and can be randomly drawn to generate a scenario on the demand. Note that d t,p is a decision variable in the linear program (A.1)-(A.4) as a scenario on the demand must be generated. Here, q is the initial qualification matrix with new (nominal) qualifications.