J.-R Dreuzy

T 1 Le Borgne

Temporal evolution of age data under transient pumping conditions

Keywords: groundwater age, residence time distribution, environmental tracers, chlorofluorocarbons, heterogeneous aquifers, transient flow conditions

While most age data derived from tracers have been analyzed in steady-state flow conditions, we determine their temporal evolution when starting a pumping. Our study is based on a model made up of a shallowly dipping aquifer overlain by a less permeable aquitard characteristic of the crystalline aquifer of Ploemeur (Brittany, France). Under a pseudo transient flow assumption (instantaneous shift between two steady-state flow fields), we solve the transport equation with a backward particle-tracking method and determine the temporal 2 evolution of the concentrations at the pumping well of CFC-11, CFC-12, CFC-113 and SF 6 .

Apparent ages evolve because of the modifications of the flow pattern and because of the nonlinear evolution of the tracer atmospheric concentrations. To identify the respective role of these two causes, we propose two successive analyses. We first convolute residence time distributions initially arising at different times at the same sampling time. We secondly convolute one residence time distribution at various sampling times. We show that flow pattern modifications control the apparent ages evolution in the first pumping year when the residence time distribution is modified from a piston-like distribution to a much broader distribution. In the first pumping year, the apparent age evolution contains transient information that can be used to better constrain hydrogeological systems and slightly compensate for the small number of tracers. Later, the residence time distribution hardly evolves and apparent ages only evolve because of the tracer atmospheric concentrations. In this phase, apparent age time-series do not reflect any evolution in the flow pattern.

Introduction

Groundwater flow is by nature transient, because of the temporal variations of boundary conditions such as the variations of recharge over different time scales (seasons, decades, centuries or more) and because of anthropogenic forcings such as pumping or artificial recharge. Pumping has a significant impact on the flow pattern and on solute transport. They induce more convergent flow pattern and, even in some cases, some extension of recharge areas [START_REF] Bredehoeft | The Water Budget Myth Revisited: Why Hydrologists Model[END_REF][START_REF] Frind | Delineation of Three-Diemnsional Well Capture Zones for Complex Multi-Aquifer Systems[END_REF]. It will as well speed up flows and modify the relative role of structures, hydrodynamic properties and boundary conditionsincreasing for instance the effective recharge rate of unconfined aquifers in close connection to the surface [START_REF] Leray | Contribution of age data to the characterization of complex aquifers[END_REF][START_REF] Sophocleous | Groundwater recharge and sustainibility in the high Plains aquifer in Kansas[END_REF].

Environmental tracers have been widely used for water sources identification, estimation of residence time distribution and model calibration amongst others [START_REF] Castro | Noble gases as natural tracers of water circulation in the Paris Basin, 2. Calibration of a groundwater flow model using noble gas isotope data[END_REF][START_REF] Cook | Groundwater ages in fractured rock aquifers[END_REF][START_REF] Long | Age-distribution estimation for karst groundwater : Issues of parameterization and complexity in inverse modeling by convolution[END_REF][START_REF] Mcmahon | Use of geochemical, isotopic, and age tracer data to develop models of groundwater flow for the purpose of water management, northern High Plains aquifer, USA[END_REF][START_REF] Stichler | Use of environmental isotopes to define the capture zone of a drinking water supply situated near a dredge lake[END_REF]. Because they integrate velocities along flow paths, they reflect flow conditions over various time scales in the past. They are sensitive to transient phenomena affecting the flow field. More precisely, they are sensitive to transient phenomena occurring over time scales comparable to their characteristic time [START_REF] Zuber | On some methodological problems in the use of environmental tracers to estimate hydrogeologic parameters and to calibrate flow and transport models[END_REF]).

Yet, the influence of transient flow conditions on environmental tracer concentration has hardly been addressed. [START_REF] Sanford | Hydrochemical tracers in the middle Rio Grande Basin, USA: 2. Calibration of a groundwater-flow model[END_REF] have reconstructed transient recharge rates using 14 C data in the regional alluvial middle Rio Grande Basin. [START_REF] Schwartz | Ambiguous hydraulic heads and 14 C activities in transient regional flow[END_REF] have noticed that the interpretation of 14 C age in transient flow models can be ambiguous in terms of flow pattern as data distributed over the aquifer reflect different flow conditions. [START_REF] Long | Age-distribution estimation for karst groundwater : Issues of parameterization and complexity in inverse modeling by convolution[END_REF] have incorporated CFCs and 3 H data from a karstic system in binary mixing model with dilution allowing parameters to vary with time. Fewer studies have analyzed the role of transient flow conditions on residence time distribution. Using numerical simulations, [START_REF] Troldborg | Using environmental tracers in modeling flow in a complex shallow aquifer system[END_REF] have showed the effect of recharge seasonality on residence time distribution and have noticed distinct behaviors. In the shallowest part of the studied heterogeneous aquifer, residence times tend to be smaller than in steady-state flow conditions while they tend to be higher in the deepest part of the aquifer. The effect in a fully-penetrating well is however negligible. [START_REF] Zinn | Potential effects of regional pumpage on groundwater age distribution[END_REF] have analyzed the effect of the start of pumping on a synthetic configuration composed of an aquifer overlain by an aquitard. Their study have revealed important changes of the mean residence time at the pumping well and of the residence time distribution over long periods of time. Changes only come from the modification of the flow pattern as they solely focused on the mean residence time and not on the apparent ages obtained from tracers.

In this study, we determine the influence of the transient groundwater flow pattern induced by anthropogenic forcing on environmental tracers concentrations and SF 6 ) interpreted in apparent ages. We consider that the transient flow pattern is induced by the instantaneous start of a pumping well. Tracer concentrations are reported at the pumping well. Our study is based on the hydrogeological setting of Ploemeur, which is a well documented aquifer where water has been produced for the last two decades for the water supply of the nearby city (Le Borgne et al., 2004;[START_REF] Le Borgne | Assessment of preferential flow path connectivity, and hydraulic properties at single-borehole and cross-borehole scales in a fractured aquifer[END_REF][START_REF] Ruelleu | Impact of gently dipping discontinuities on basement aquifer recharge: An example from Ploemeur (Brittany, France)[END_REF]Touchard, 1999). Although based on a specific site, the results of this study can be generalized to shallowly dipping aquifers overlain by a leaking layer. Such condition has been proved to be of importance for groundwater resources in hard-rock aquifers [START_REF] Leray | Numerical modeling of the productivity of vertical to shallowly dipping fractured zones in crystalline rocks[END_REF]. While our objective is more methodological than targeted to a specific site, the Ploemeur aquifer still offers a complex and yet realistic hydrodynamic context. We use hydrogeological models previously calibrated in steady-state flow conditions under pumping [START_REF] Leray | Contribution of age data to the characterization of complex aquifers[END_REF] and determine the effect of transient flow conditions on apparent ages. We first aim at determining the causes of the temporal evolution of the apparent ages and specifically when they rather come from the transient modifications of the flow pattern and when they are more linked to the specificities of the tracers, especially their atmospheric concentrations. We second aim at assessing the interest of age data time series for models segregation. After recalling in section 2 the hydrogeological, flow and transport models, we present the results in section 3 and discuss them in section 4.

Hydrogeological, flow and transport models

We successively describe the hydrogeological models of the Ploemeur site that will be used as a basis of this study, the flow and the transport models as well as the numerical methods used.

We finally comment in details the derivation of the tracer concentrations and the corresponding apparent ages to highlight the possible causes of apparent age temporal variations.

Hydrogeological model

The study of the effects of transient flow conditions, induced by pumping, on age data is based on the Ploemeur aquifer, a highly heterogeneous hard-rock aquifer located on the south coast of Brittany near the city of Lorient (France). A previous study based on the inversion of gravimetric data has established a geological conceptual model [START_REF] Ruelleu | Impact of gently dipping discontinuities on basement aquifer recharge: An example from Ploemeur (Brittany, France)[END_REF]. This conceptual model is composed of two transmissive structures at large scale, the dipping contact zone and a North 20° normal fault, besides the Ploemeur and Guidel granites and overlying micaschists acting as a typical aquitard. Local heterogeneities are not represented in the model. The supplying area to the pumping well which amounts to a few square kilometers is limited in the North-South direction by these two granites. The pumping rate thus has a strong impact on flow pattern within this heterogeneous aquifer.

Because the shape and the dip of the contact zone are only partially known, the overall thickness of the aquitard-aquifer system remains relatively uncertain. To account for this uncertainty, a few structural models with distinct thickness have been built (Figure 1 and Table 1). The hydraulic properties of the different rocks have been set either to common values as for the granites which are found almost impervious (10 -11 m/s), or to measured values as for the micaschists (10 -7 m/s -5x10 -6 m/s) and the contact zone (1.9x10 -3 m 2 /s -3x10 -3 m 2 /s). In addition, the potential recharge rate R has been estimated at 200 mm per year [START_REF] Carn | Mise en valeur des ressources en eau souterraine du socle breton[END_REF][START_REF] Leray | Contribution of age data to the characterization of complex aquifers[END_REF]Touchard, 1999). Following these constraints, each model has been calibrated against the mean piezometric level measured at the pumping well (-5.5 masl) in steady-state flow under pumping conditions [START_REF] Leray | Contribution of age data to the characterization of complex aquifers[END_REF] by slightly adjusting the contact zone transmissivity previously estimated from long-term pumping tests (Le Borgne et al., 2006). Uniform porosity has also been calibrated against the CFC-12 age (30 years ± 1 year in 2009). Note that the overall volume of the system is about 1.5x10 9 m 3 and the mean residence time of the modeli.e. the first moment of the residence time distributionis around 13 years in ambient conditions and 50 years in pumping conditions. .

Our study has been carried out on about ten representative hydrogeological models differing by their structure, their micaschists permeability and their porosity. The interest of considering different models is to investigate the potential influence of the hydrogeological structure on the apparent ages and their evolution. We discuss further in section 4 how this sensitivity might be useful as an additional way to characterize the flow pattern. Among this set of models, only two are used here to illustrate the methodology as they all lead to the same conclusions. Table 1 synthesizes the parameters of the two chosen hydrogeological models.

Flow model

Transient flow conditions are induced by starting a pumping well. The transient pumping rate Q w (t) is a step function going from zero before the starting date t switch , to a positive value Q p :

.

In the particular case of the site of Ploemeur, Q p is set at 3.36x10 -2 m 3 /s. Pumping started in 1991 and the most part of the evolution of the piezometric levels occurred only in a few years after the start of pumping. t switch has thus been set at 1994. We solve the 3D diffusivity equation for the hydraulic head h(x, t) with free surface boundary conditions under a pseudo transient flow approximation:

( 2 )

( 3 ) ( 4 )

( 4, Figure 2) do not have any impact on the recharge areas captured by the pumping zone located at depth in the pumping well. Imposed heads, set at depth z 0 (5 meters below the ground surface level), are applied to the South and the North boundaries without any significant influence on the system because of the almost impervious granites (Equation 5, Figure 2).

Transport model

Transport is considered only advective as the macro-scale dispersion from local dispersion and diffusion is assumed to have a much smaller effect compared to the macro-scale dispersion induced by structural heterogeneity and sampling [START_REF] Labolle | Role of molecular diffusion in contaminant migration and recovery in an alluvial aquifer system[END_REF].

Transport is modeled by the advection equation [START_REF] Bear | Modelling and applications of transport phenomena in porous media. Theory and applications of transport in porous media[END_REF][START_REF] De Marsily | Quantitative hydrogeology: Groundwater Hydrology for Engineers[END_REF]): ( 10 )

Numerical methods

The unconfined flow equations are solved using a computationally effective finite-volume approach with a local adaptation scheme [START_REF] Bresciani | A finite volume approach with local adaptation scheme for the simulation of free surface flow in porous media[END_REF]. The advection equation is solved in backward-time that consists in reversing the flow field and adapting the boundary conditions [START_REF] Neupauer | Adjoint method for obtaining backward-in-time location and travel time probabilities of a conservative groundwater contaminant[END_REF][START_REF] Neupauer | Adjoint-derived location and travel time probabilities for a multidimensional groundwater system[END_REF][START_REF] Neupauer | Backward probabilistic model of groundwater contamination in non-uniform and transient flow[END_REF]. To solve the transport equation, we use a Lagrangian random walk method well-suited to purely advective transport [START_REF] De Dreuzy | Asymptotic dispersion in 2D heterogeneous porous media determined by parallel numerical simulations Water[END_REF][START_REF] Kinzelbach | The random-walk method in pollutant transport simulation[END_REF] and adapted to the unconfined conditions and to the backward-time resolution. We inject particles proportionally to flow [START_REF] Kreft | On the physical meaning of the dispersion equation and its solution for different initial and boundary conditions[END_REF] at the pumping well cell and tracked them to the aquifer free surface. The chosen number of injected particles (5x10 6 ) is determined through a convergence test of the mean and the variance of the residence time distribution. Flow and transport simulations are carried out using the H2OLAB platform [START_REF] Bresciani | A finite volume approach with local adaptation scheme for the simulation of free surface flow in porous media[END_REF][START_REF] Erhel | Multi-parametric intensive stochastic simulations for hydrogeology on a computational grid[END_REF][START_REF] Erhel | A parallel scientific software for heterogeneous hydrogeoloy[END_REF].

Computation and temporal evolution of apparent age

We derive the residence time distribution p(t) from the equally likely residence times given by the particles. The mean concentration C w of a tracer at the pumping well and at the sampling date t w is determined by the convolution of the residence time distribution p(t) with the atmospheric concentration C in (t w -t) of the tracer [START_REF] Kreft | On the physical meaning of the dispersion equation and its solution for different initial and boundary conditions[END_REF][START_REF] Maloszewski | Determining the turnover time of groundwater systems with the aid of environmental tracers, 1. Models and their applicability[END_REF]:

. ( 11
)
Note that, to highlight that the residence time distribution p(t) is highly dependent on a given sampling date t w , we later use it as p(t wt) where t wt represents the date at which water recharged, rather than as a function of the residence time "t" alone. Throughout the paper, p(t wt) is named recharge date distribution. The apparent age A [START_REF] Nir | On the interpretation of tritium 'Age' measurements of groundwater[END_REF] of a tracer at the pumping well x w and at the sampling date t w is then obtained from the difference between t w and the date at which the concentration C w is equal to the atmospheric concentration:

, ( 12
)
with the reciprocal function of the atmospheric concentration.

Equation 12 shows the two possible causes of the variation of age A with the sampling date t w .

The first cause is the evolution of the recharge date distribution p(t w -t) because of the modification of flow conditions that can be traced back to equation 1 through equations 2-11.

The shift from ambient to pumping conditions does not only change the flow magnitude but also the flow pattern, the recharge and discharge locations as well as the water mixing within the radius of action of the well [START_REF] Bredehoeft | The Water Budget Myth Revisited: Why Hydrologists Model[END_REF]. This first cause is solely linked to the flow conditions and does not depend on the tracer characteristics.

The second cause of the age variation directly comes from the tracer atmospheric concentration C in (t) through the convolution of equation 11. The non-linear temporal evolution of C in (t) modifies the sampling of the recharge date distribution p(t w -t). Monotonic and steep evolutions of C in (t) characteristic of the SF 6 and CFCs approximately between the 1970s and 1990s results in a large range of weights of the residence times in this period.

Flatter evolutions characteristic of the CFCs after the early 1990s give a higher but more equilibrated contribution of the shorter residence times. Depending on the sampling date t w , this equilibrated contribution is more or less significant. This effect occurs both in transient and steady-state flow conditions [START_REF] Troldborg | Using environmental tracers in modeling flow in a complex shallow aquifer system[END_REF][START_REF] Waugh | Relationships among tracer ages[END_REF][START_REF] Zhang | Numerical simulations of dating young groundwater with multiple atmospheric tracers : CFC-11, CFC-12, SF 6[END_REF].

Even under steady-state flow conditions, the non-linear evolution of the atmospheric concentration C in (t) lets the sampled concentration evolve. Under simpler terms, the evolution of the sampled concentration does not only come from the evolution of the system but also from the modification of the "observation device" (C in ).

It should be noted that this second effect is irrelevant when mixing is minimal like within the framework of the piston-flow model. In such cases, the residence time distribution resumes to a Dirac and the dependency to the tracer atmospheric concentration vanishes in equation 12 because of the direct transformation of C in by C in -1 . This effect is, however, very relevant close to the aquifer discharge zones (springs, wells) where the mixing of flow paths is maximal. The mixing of sampled flow lines enhances the importance of the non-linear evolution of C in (t). The case studied here pertains more to this second situation, in which, finally, both the temporal evolutions of the recharge date distributions p(t w -t) and the tracer atmospheric concentration C in (t) can modify the sampled concentration and the derived age

A(t w ) (equation 12).
To analyze the temporal evolution of the apparent ages, we have determined them using equation 12 by convoluting the recharge date distributions p(t w -t) obtained under the pseudo transient conditions with the CFC-11, CFC-12, CFC-113 and SF 6 atmospheric concentrations at the evolving sampling dates from 1994 to 2009 (Table 2). These dates correspond to sampling 0 to 15 years after the change of flow conditions from ambient to pumping conditions occurring at t switch =1994.

Results

This section first reports the temporal evolution of the apparent ages at the pumping well for the hydrogeological model numbered 1 of the site of Ploemeur described in section 1 and in Table 1. It then analyses the respective effects of the modification of the flow pattern and of the evolution of the tracer atmospheric concentration. The implication of the temporal evolution of the apparent ages for models segregation is further discussed in section 4.1 by comparing the results of the models numbered 1 and numbered 2.

Temporal evolution of apparent age

Figure 3 shows the evolution of the apparent ages from CFC-11, CFC-12, CFC-113 and SF 6 concentrations with the sampling date t w ranging from 1994 (ambient conditions) to 2009 (Table 2). At first sight, the apparent ages increase first sharply after the start of pumping and then more smoothly whatever the tracer. We indeed expect that the modifications of flow are maximal just after the start of the pumping and later decrease. Though similar, the apparent ages derived from the CFCs and SF 6 still exhibit some differences induced by the tracer atmospheric concentration C in (t). The temporal evolution is minimal for SF 6 because of the almost linear increase of its atmospheric concentration and maximal for the CFC-11 and CFC-12 because of the strong and non-monotonic variations of their atmospheric concentration (IAEA, 2006). CFC-11 and CFC-12 display similar variations because of the similar shapes of their atmospheric concentration chronicles.

The increase of the apparent ages in the models considered here is somehow counterintuitive.

Indeed, in the exponential model [START_REF] Haitjema | On the residence time distribution in idealized groundwatesheds[END_REF], the residence time distribution is independent of the pumping rate. Besides, in the piston-flow model, as said in section 2.5, apparent ages remain constant whatever the sampling dates due to minimal mixing. When starting a pumping in those conditions, one would expect an increase of velocity that would directly reduce the residence time and then the apparent ages. In the more complex system modeled here, pumping has a nontrivial effect on the apparent ages. Further insight is given by the analysis of the recharge date distribution p(t w -t) and its temporal evolution. Figure 4 shows the recharge date distribution sampled at three different dates : (a) at under ambient flow conditions, (b) at approximately one month after the pumping started and (c) at one year after the pumping started. At later dates (t w 1995), distributions stop to evolve. Therefore, only the distribution at is displayed.

For , the system is strongly transient (Figure 4a andb). The start of pumping induces a modification of the recharge date distribution from a narrow piston-like distribution at ambient conditions ( ) (Figure 4a) to a much broader, more exponential-like, later (Figure 4c). Intermediary distribution shapes occur in the first pumping year (Figure 4b).

The broadening of the distribution comes from the shift of status of the sampling zone. Under ambient flow conditions, it is a standard zone within the aquifer broken through by just a few flow lines. The different tracers then lead to very close ages (Figure 3). Under pumping conditions, it becomes the major discharge zone of the aquifer focusing a dense net of flow lines that initially discharged in much larger areas such as wetlands (Figure 2) and apparent ages from the different tracers spread out. For , the recharge date distribution broadens including both shorter and larger residence times than the ambient distribution. In the case of the Ploemeur site, shorter residence times eventually control the apparent ages.

This is also consistent with the increase of the mean effective recharge rate and the induced circulation speed up (from 160 mm/year in ambient conditions to 200mm/year in pumping conditions). As previously noted, the distributions sampled in 1995 and later are almost identical. Thus, for , sampled concentrations are likely to be more influenced by the temporal evolution of the tracer atmospheric concentrations C in (t) and less by the evolution of the recharge date distribution.

3.2.

Effect of the temporal evolution of the recharge date distribution

Methodology

To assess the role of the temporal evolution of the recharge date distribution (p(t w -t) in equation 12), we filter out the evolution of the tracer atmospheric concentration (C in (t) in equation 12). To this end, we translate the distributions (Table 2) to the same sampling date , the exponent "t" standing for translated. The recharge date distributions are shifted along the date axis by a translation of without any modification of their shape. The shifted distributions are noted :

( 13 ) This is illustrated on Figure 5 6). The apparent ages from the distribution initially sampled at (i.e. one month after the pumping started) constitute an intermediate case. Indeed, unlike the ages from the ambient distribution, they significantly diverge for the four tracers. However, their spreading is not as broad as for the distributions sampled at later dates ( ). We will discuss in section 4.1 the interest of these differences for model segregation. We have checked that the conclusions were the same for the other sampling dates of Table 2.

Results

3.3.

Effect of the temporal evolution of the atmospheric concentration

Methodology

To assess the influence of the tracer atmospheric concentration C in (t) independently of the flow pattern modifications, we compute the apparent ages at six translated sampling dates (Table 2) for a fixed recharge date distribution.

Results

Discussion

We first discuss the potential interest of the two successive evolutions of apparent ages, i.e. a first phase dominated by the transient flow pattern and a second one by the transient atmospheric concentration for flow pattern characterization and models segregation. We secondly compare these modeling results to the available data on the site of Ploemeur.

Information contained in the temporal evolution of apparent ages

The temporal evolution of apparent ages contains information on the hydrogeological system, in addition to information contained in the tracer concentrations analyzed at a single time.

Yet, information is largely different between the two evolutions identified in sections 3.2 and 3.3. The first evolution is the rapid shift when starting pumping from a restricted recharge date distribution to an extended one. The restricted recharge date distribution can be typically fitted by an inverse Gaussian model characterizing mainly the length of the flow paths from the recharge zone divided by the recharge rate [START_REF] Ginn | Notes on Groundwater Age in Forward and Inverse Modeling[END_REF][START_REF] Woolfenden | Modeled ground water age distributions[END_REF]. The broader recharge date distribution characterizes more globally the aquifer volume and the overall recharge rate [START_REF] Leray | Contribution of age data to the characterization of complex aquifers[END_REF], so does typically the exponential model [START_REF] Haitjema | On the residence time distribution in idealized groundwatesheds[END_REF]. Therefore, almost independent and useful pieces of information, such as a characteristic length vs. a volume, can be obtained using the same well before and after the start of pumping. Using equation 12 as an illustration, it comes down to use the temporal evolution of apparent ages to estimate the transient controlling parameters of the function p(t w -t).

Figure 9 compares the temporal evolution of the CFC-12 age for the models numbered 1 and numbered 2 (Table 1) at sampling dates ranging from 1994 to 2009 (Table 2). It shows that the two models give distinct apparent ages before 1995 i.e. when the recharge date distribution is transient. Later, once the recharge date distribution is broad and reaches a steady-state, and because the two models have been calibrated on CFC-12 age in 2009, the apparent ages are identical for both models. Differences only occur at short term after the pumping started and result from different transient evolutions of the residence time distribution. This underlines the interest of a continuous sampling of environmental tracer concentrations, particularly at the early times of pumping, to characterize the transient behavior of the system. In that case, apparent age time series can be used as an additional tool Still, CFC-113 exhibits an opposite age variation in comparison with CFC-12. Considering the age uncertainty due to analytical and sampling errors [START_REF] Leray | Contribution of age data to the characterization of complex aquifers[END_REF] Before pumping, the recharge area is a small zone upstream of the well (blue line). It quickly expands to a much larger area around the well only after one year of pumping (green line).

Actually, the pumping well collects water that naturally discharged in wetlands next to it and transforms a wide discharge zone to a point-like one.

The 2 nd phase corresponds to a slight increase of the chloride concentration (Figure 10). This increase does not come from the changes of the recharge area induced by the evolution of the flow pattern, in this case marginal (from orange to purple curves in Figure 11). It more likely comes from the modification of the deep water fraction coming from the main aquifer -i.e.

the contact zone -that dips quasi-vertically at about 1,500m from the pumping well (Figure 11, [START_REF] Ruelleu | Impact of gently dipping discontinuities on basement aquifer recharge: An example from Ploemeur (Brittany, France)[END_REF]). Its evolution only affects recharge dates earlier than 1940 and cannot consequently be detected by atmospheric gases. Finally, the stabilization of the chloride concentration (3 rd phase of Figure 10) agrees with the apparition of a new steadystate regime.

Conclusion

We have analyzed the temporal evolution of apparent ages in a complex aquifer under We identify two well-differentiated phases in the evolution of apparent ages. Apparent ages first evolve because of the pumping-induced modifications of flows. At this stage, the temporal evolution of apparent ages is distinct between two hydrogeological models and may be advantageously used to segregate them using for instance a classical least square comparison procedure of data with modeling results. The transient measure of the tracer concentrations contain some additional information and slightly compensate for the small number of available tracers. After one year of pumping, residence time distributions hardly evolve and apparent ages become solely modified by the transient evolution of the atmospheric concentrations. Apparent ages still significantly evolve but do not contain any additional information on the flow patterns beyond those contained in the steady-state data.

Tables Table 1: Parameters of the two hydrogeological models of the site of Ploemeur used in this study with their reference work. H TOT represents the mean thickness of the aquifer system composed of the micaschists and of the contact zone (Figure 1). The two models have been calibrated on the piezometric level at the pumping well and on CFC-12 apparent age at the pumping well. and SF 6 concentrations as functions of the sampling date t w for the model numbered 1 (Table 1). (Table 1). The translated distributions are the ambient distribution (initially sampled at , blue curve) and two pseudo transient distributions initially sampled at (wine curve) and at (green curve). Note that the ordinate axis has been broken between 0.05 and 0.19 to display the peak of the ambient distribution. of the model numbered 1 (Table 1). CFC-12 ages are computed at the pumping well and at sampling dates between 1994 and 2009 (Table 2). Labels next to the crosses stand for the corresponding sampling dates. We recall that the two models have both been calibrated on the CFC-12 age in 2009. Results for the tracers CFC-11, CFC-113 and SF 6 are similar and thus are not shown.

Highlights:

 Pumping start is modeled by an immediate shift between two steady-state flow fields  Ages deduced from CFCs and SF 6 concentrations evolve in two distinct phases  Transient flow patterns affect ages just after the pumping start but quickly vanish  Atmospheric concentrations then transiently weight the residence time distribution  Similar in terms of magnitude, these two regimes are helpful for models segregation

  x,t) the solute concentration at the position x and at the time t; θ the effective porosity; Q w (t) the rate of the pumping well; C w (t) the solute concentration at the pumping well; C 0 (x) the initial condition; C 1 (t) is the boundary condition on a first-type boundary (Γ 1 ) -such as the tracer atmospheric concentration at the free surface boundary -and q(x,t) the pseudo transient Darcy's flux derived from the head field h(x,t):

  for the three recharge date distributions of Figure 4 translated to . The apparent age determination is modified accordingly: ( 14 ) Equation 14 replaces equation 12 and filters out most of the effect of the tracer atmospheric concentration C in (t) to highlight that of the evolution of the flow pattern. The apparent ages are noted recalling that they are obtained after the translation of the recharge date distributions.

Figure 6

 6 Figure 6 compares the apparent ages obtained with the various recharge date distributions

  Figure 7 illustrates this analysis for the pseudo transient distribution initially sampled at and translated to 2004 and 2009. This transformation should not be confounded with the previous one in which we studied the influence of the recharge date distributions at fixed atmospheric concentration C in (t).

Figure

  Figure 8a, b and c display the results for the three distributions of Figure 4a, b and c. For the

  , the acquisition of additional concentration data over the next few years would help to confirm the age trends and eventually conclude about the flow pattern stabilization. Earlier flow conditions, and particularly initial ones i.e. just before the pumping started, are more difficult to constrain on the site of Ploemeur because of the lack of environmental tracers monitoring before 2006. Contrary to CFCs and SF 6 data, chloride and nitrate concentrations have been monitored on a longer date range starting at 1991. Because of reactivity processes and their evolution with time, the nitrate concentration chronicle can only give partial information on the flow pattern and on mixing processes. However, its non-negligible initial concentration (roughly 15 mg/l) before the start of pumping indicates short residence times consistent with the local recharge area and the short flow paths obtained in the models.Chloride has the advantage of being conservative and being characteristic of deeper waters and thus a good indicator of the flow pattern and its evolution. Figure10represents the temporal evolution of chloride concentration at the pumping well. It displays a sharp increase during the first two years of pumping (1 st phase of Figure10), a very slight increase for the next ten years (2 nd phase of Figure10) and finally stabilizes after 2005 (3 rd phase of Figure10). The first phase is consistent with the previous modeling result of a quick transition from the piston-like recharge date distribution to the more extended one. The short duration of the first phase is consistent with the fast temporal evolution of the recharge area (Figure11).

  pumping based on CFCs and SF 6 concentrations measured at the pumping well. The evolution of apparent ages at the pumping well with time can come from (1) the transient nature of flow conditions and (2) the transient evolution of the tracer atmospheric concentrations. To identify the respective role of these two sources, we proposed two successive analyses of the residence time distributions: (1) convolution at the same sampling date (i.e. with the same atmospheric concentration chronicle) of residence time distributions initially arising at different dates to assess the effect of the distributions, (2) convolution of one residence time distribution at various sampling dates (i.e. with various atmospheric concentration chronicles) to assess the effect of the atmospheric concentration.

FiguresFigure 1 :

 1 Figures

Figure 2 :

 2 Figure 2: Boundary conditions applied to the model i.e. no flow on the East and the West

Figure 4 :

 4 Figure 4: Recharge date distributions p(t w -t) of the model numbered 1 (solid curves, left

Figure 5 :

 5 Figure 5: CFC-12 atmospheric concentration C in (t) (grey dashes) superimposed on some

Figure 6 :

 6 Figure 6: Apparent ages derived from CFCs and SF 6 concentrations for the translated

Figure 7 :Figure 8 :

 78 Figure 7: Recharge date distribution p 3 at its initial sampling date (green solid

Figure 9 :

 9 Figure 9: CFC-12 ages of the model numbered 2 (blue crosses) as a function of CFC-12 ages

Figure 10 :

 10 Figure 10: Chloride concentration as a function of the sampling date t w measured at the

Figure 11 :

 11 Figure 11: Surface origin of the CFCs and SF6 (corresponding to recharge dates posterior to

is the hydraulic conductivity; n is the outgoing normal to the saturated domain; R is

  the potential recharge rate; z ground is the ground surface elevation; Γ s is the top of the saturated domain; Γ west , Γ east , Γ north and Γ south are respectively the West, East, North and South Equation 2 is solved in unconfined conditions since seepage conditions are a priori not known. Unconfined conditions are satisfied through both conditions at the free surface boundary of equation 3. When the free surface level is below the ground surface level, the effective recharge rate is equal to the potential recharge rate R; anywhere else, the free surface

	transport evolution. This is a reasonable approximation of transient flow conditions valid at
	low specific storage and with the advantage of being less costly numerically. Practically, it
	consists in ignoring the transition between the two steady-state velocity fields under ambient
	and pumping conditions.
	5 )
	( 6 )
	where K(x) sides of the domain; z 0 is a reference height; n w is the ingoing normal to the well screen; Γ w is
	the well screen surface and Q w (t) is the transient pumping rate defined in equation 1 and
	located at x w . The pseudo transient approximation intervenes in equation 6 and consists in
	assuming that steady-state flow conditions are quickly established compared to the solute

level is set at z ground and the effective recharge rate continuously evolves from negative values in the discharge zone to positive potential recharge rate R. The model has been built to include the nearest watersheds both in ambient and pumping conditions in order to minimize the potential effects of the boundary conditions. No-flow conditions applied on the West and East boundaries (Equation

Table 2 :

 2 Sampling dates where the apparent ages are determined.

Sampling dates t w Acronym Name and color for corresponding p (t w -t )

	1994	t w 1	p 1 = p (t w 1 -t ), blue
	1994.1	t w 2	P 2 = p (t w 2 -t ), wine
	1995	t w 3	P 3 = p (t w 3 -t ), green
	1999	t w 4	p 4 = p (t w 4 -t ), not displayed
	2004	t w 5	p 5 = p (t w 5 -t ), not displayed
	2009	t w 6	p 6 = p (t w 6 -t ), not displayed

Acknowledgments Funding was provided by the French National Research Agency ANR through the H2MNO4 project for the development of parallel simulation methods (ANR-11-MN). Research and monitoring on the site of Ploemeur is funded by the Environmental research Observatory H+ (Network of hydrogeological sites) and by the European Interreg IV project Climawat. J.-R. de Dreuzy acknowledges the European Union for its additional funding through the IEF Marie-Curie fellowship (PIEF-GA-2009-251710). The authors acknowledge Arash Massoudieh and an anonymous reviewer for their fruitful reviews.

to segregate hydrogeological models. Using equation 12 as an illustration, it comes down to use the temporal signal in apparent ages to segregate two functions p(t w -t) that exhibit distinct transient behaviors.

4.2.

Comparison to field observations on the site of Ploemeur 2012), degradation of CFC-11 under anaerobic conditions may explain these differences [START_REF] Cook | Transport of atmospheric trace gases to the water table -Implications for groundwater dating with chlorofluorocarbons and krypton-85[END_REF]. SF 6 apparent ages are less than 10 years. Such discrepancies cannot be explained by hydrodynamic mixing alone even at the pumping well [START_REF] Leray | Contribution of age data to the characterization of complex aquifers[END_REF] but are rather coming some geogenic production from the neighboring granites [START_REF] Busenberg | Dating young groundwater with sulfur hexafluoride: Natural and anthropogenic sources of sulfur hexafluoride[END_REF][START_REF] Koh | Evidence for terrigenic SF 6 in groundwater from basaltic aquifers, Jeju Island, Korea : Implications for groundwater dating[END_REF].

The results of the models in section 3 suggest that with negligible temporal change of the recharge date distribution, the apparent age should continuously increase whatever the tracer and within the date interval 2000-2010. The CFC-12 age increase at the Ploemeur site is similar in terms of magnitude to that of the models. CFC-12 data thus suggest that the tracer age evolution in that discharge area mainly comes from the temporal evolution of the atmospheric concentration weighting and that the flow pattern and the resulting recharge date distribution in the pumping zone are in turn not expected to vary significantly after 2006.


Romain Guillaume
email: romain.guillaume@irit.fr

Caroline Thierry
email: thierry@univ-tlse2.fr

Bernard Grabot
email: bernard@enit.fr

Modelling of ill-known requirements and integration in production planning

Keywords: requirement, demand, imprecision, uncertainty, theory of possibility

Making decisions on the base of uncertain forecasts is one of the key challenges for efficient Supply Chain Management. This article suggests the use of the theory of possibility for building a procurement plan on the base of ill-known requirements. These requirements, expressed in quantities by date, may come from various sources: forecasts or orders for instance. The possible types of imperfection pervading requirement are analysed and a unified representation model is suggested. A method is then described for calculating a plausible demand by period without loss of information; it is illustrated on an example in the last section.

INTRODUCTION

Today's industrial context is characterised by an increased uncertainty on the final customer's demand.

In many supply chains, and especially in the aeronautical industry, the supply times are much longer than the manufacturing times. In this context, an important challenge is to manage the procurement process in spite of the uncertainty on the demand.

Collaborative processes aim at reducing this supply risk. Vendor Management Inventory (VMI ® ) [START_REF] Marques | Vendor Managed Inventory, from concept to processes, for an unified view[END_REF] and Collaborative Planning, Forecasting and Replenishment (CPFR ® ) [START_REF] Ireland | Supply chain Collaboration: how to implement CPFR ® and other best collaborative practices, Integrated Business Management series and APICS[END_REF] are for instance widely used in an industrial context. These techniques aim at creating short and reactive decision loops between customers and suppliers in order to cope with the growing uncertainty on demand forecasting, due to the shortening of the product life cycle and to customers' versatility. Nevertheless, they can only deal with quite standard requirements, which is for instance the case in the automotive industry. Taking into account the uncertainty of the demand in the planning processes can be an alternative solution for more specific requirements, e.g. in the spatial and aeronautical industry for which the problem is still critical.

In this study is considered a dyadic elementary link of a supply chain, involving two actors (customer/supplier). Within these two partners, it will be shown that an uncertain demand can be considered in order to build a "fuzzy" procurement plan, firstly internally processed then transmitted to the upstream suppliers. Like a classical one, this procurement plan is expressed in quantities per periods, but does not only describe one (most plausible) situation, but a set of situations, with a quantification of their possibility of occurrence. It should then allow to make more informed decisions concerning the risk taken by ordering a given quantity.

The article is organized as follows: section 2 presents a state of the art on the application of possibility theory to production planning problems. The basis of possibility theory and the problem of representation of imperfect information are introduced in the following section, whereas a way to represent "imperfection" (uncertainty, imprecision…) on the requirements is proposed in section 3. The integration in the Material Requirement Planning of imperfect requirements, using the suggested method, is described in section 4. The approach is illustrated on an example in section 5.

LITERATURE REVIEW

In a manufacturing context, production planning aims at optimizing the use of the resources to fulfil the customers' demand. Usually, forecasts are estimated on the base of historical data, and allow to build a Sales and Operation Plan, implicitly considered as certain but periodically revised (SOP step of the MRP (Manufacturing Resource Planning) method). A Master Production Schedule then details the SOP, listing quantities of required products per periods of time. Using the bill of materials, these "independent" requirements are then used to generate requirements on components, also by periods. In an uncertain context, the problem is to build a feasible requirement plan on the base of an ill-known demand.

In that purpose, a large part of the literature focuses on rolling horizon planning methods: an initial planning is periodically updated in order to take into account the most recent knowledge on the demand. Strategies like adjusting the re-planning horizon size [START_REF] Galasso | Decision Support for Supply Chain Planning under uncertainty[END_REF], or freezing a part of the plan [START_REF] Zhao | Multilevel lot-sizing heuristics and freezing the master production schedule in material requirements planning systems[END_REF]) may be used in that purpose. On the other hand, production planning under uncertainty aims at modelling this uncertainty in order to provide the decision maker with a richer view on the various possible situations which may arise. Within this strategy, stochastic and possibilistic models have been mainly used to represent imperfect data. Stochastic models need a huge amount of historical data, and require that the context in which the data have been collected remains still valid in the future. Since this is more and more rare, we focus here on the case when it is difficult to use such methods. In that case, uncertainty may often be modelled on the base of the human expertise. Fuzzy logic and theory of possibility have shown their interest for representing subjective information and some works have already been performed on the representation of a subjectively assessed uncertainty on the demand.

Within these works, it is important to distinguish production planning under flexible constraints from production planning with ill-known parameters [START_REF] Dubois | Fuzzy scheduling: Modelling flexible constraints vs. coping with incomplete knowledge[END_REF]. In the literature on production planning under flexible constraints (Mula et al., 2007, Jamalnia and[START_REF] Jamalnia | A hybrid fuzzy goal programming approach with different goal priorities[END_REF], fuzziness models a level of satisfaction of constraints expressed by imprecise terms (often in natural language), or a degree of tolerance [START_REF] Bellman | Decision making in a fuzzy environment[END_REF]. Several surveys [START_REF] Kacprzyk | Fuzzy optimization and mathematical program-ming: a brief introduction and survey[END_REF][START_REF] Delgado | Fuzzy Optimization: "Recent Advances[END_REF][START_REF] Baykasoglu | A review and classification of fuzzy mathematical programs[END_REF] show how fuzziness can in that case be accommodated with linear programming.

In the literature on production planning under uncertainty, fuzziness is used for modelling the level of possibility that an event may occur. Both approaches refer to the possibility theory (Zadeh, 1978, Dubois and[START_REF] Dubois | Possibility Theory[END_REF] and to decision under fuzzy uncertainty [START_REF] Dubois | Decision-theoretic foundations of qualitative possibility theory[END_REF].

Production planning under uncertainty based on possibility theory may allow to compute a supply plan considering an uncertain demand, modelled by possibility distributions. Some authors suggest to make decisions on the base of a fuzzy plan, while others use defuzzification techniques for obtaining a crisp plan. In the first category can be found [START_REF] Fargier | The use of Possibilistic Decision Theory in Manufacturing, Planning and Control: Recent Results in Fuzzy Master Production Scheduling[END_REF], in which the MPS (Master Production Schedule) considers an uncertain demand. A method is proposed to support subcontracting decisions or negotiation with the customer. In the same way, in [START_REF] Grabot | Integration of uncertain and imprecise orders in the MRP method[END_REF], the imperfection on the customer's demand is taken into account in the Manufacturing Resource Planning through a method called F-MRP (Fuzzy-MRP), resulting in a "fuzzy plan".

Another solution is to model first the uncertain plan, then to defuzzify in order to obtain a crisp equivalent demand (demand after defuzzification) [START_REF] Niu | Application of Fuzzy MRP II in fast moving consumer goods manufacturing industry[END_REF][START_REF] Liang | Fuzzy multi-objective production/distribution planning decisions with multi-product and multi-time period in a supply chain[END_REF][START_REF] Peidro | Fuzzy optimization for supply chain planning under supply, demand and process uncertainties[END_REF]).

An alternative is to build a plan which takes into account a level of possibility, necessity or credibility in order to guarantee or optimise a "level of chance" that the plan is feasible. In [START_REF] Lan | Modeling fuzzy multi-period production planning and sourcing problem with credibility service levels[END_REF], the chance that the inventories will be superior to zero is guaranteed by the constraints while in [START_REF] Sun | Fuzzy Minimum risk material procurement planning problem[END_REF], the objective to minimize is the chance than the cost function is inferior to a given threshold. In [START_REF] References Aliev | Fuzzy-genetic approach to aggregate production-distribution planning in supply chain management[END_REF], the objective function to minimize is the chance that the demand will not be satisfied, combined with the profit. [START_REF] Liu | Theory and practice of uncertain programming[END_REF] presents a set of other models for optimisation under uncertainty with fuzzy parameters, which could be used in supply chain planning.

When dates related to precise orders and durations of tasks are considered, scheduling under uncertainty may result in a fuzzy schedule [START_REF] Chen | Order-fulfillment ability analysis in the supply-chain system with fuzzy operation times[END_REF]Huang, 2006, Balasubramanian and[START_REF] Balasubramanian | Scheduling optimization under uncertainty an alternative approach[END_REF]. These references are summarised in Table 1. Even if possibility theory has been intensively used for modelling uncertainty in planning processes, it is interesting to notice that none of these methods addresses the basic industrial problem in its more generic form: if the imperfection of the demand concerns both the quantities (i.e. the level of a requirement may change) and the dates (the requirement may be advanced or postponed), finding how many components may be required in a given period in order to fulfil this demand becomes a complex task. In order to address this problem, this article first suggests a method aiming at characterising and representing "uncertainty" (imperfection, uncertainty and imprecision) on the independent requirements, then shows how to turn these "imperfectly known" requirements into dependent requirements (a supply plan composed of required quantities of components by periods), in a way which stays as close as possible to industrial practices. This supply plan will be considered for deciding the quantities ordered to the upstream partners of the supply chain.

BASIC FRAMEWORK

The objective of this section is to position the various concepts covering the vague notion of "imperfect information", and to present some models allowing to represent them. Besides are introduced the main results of the theory of possibility which will be used afterwards.

Definition of imperfection

Even if the word "uncertainty" is commonly used, a distinction has to be made between several concepts: uncertainty, imprecision and incompleteness. The term "imperfection" will be used to gather these three concepts, as suggested in [START_REF] Bouchon-Meunier | La logique floue et ses applications[END_REF].

-Uncertainty occurs when it is not possible to state whether a particular statement is definitely true or false.

-Imprecision denotes the impossibility to enounce a piece of knowledge including crisp values in a sure manner.

-Incompleteness means in most cases that an information source does not contain full information about the values of the attributes of the real-world it describes (some properties may be missing from the description).

As stated above, two types of mathematic models can represent imperfect data: stochastic and possibilistic models [START_REF] Dubois | Représentations formelles de l'incertain et de l'imprécis[END_REF]. Because of the highly dynamic context of today's manufacturing, we consider in this article the case when no previous observation is available for building a stochastic model. In that case, possibility theory has shown its interest for modelling uncertainty on the base of expert knowledge.

Representation of imperfection

3.2.1Representation of uncertainty

For each event A from a subset S, the uncertainty may be evaluated through two levels: the possibility П(A), denoting to what extent the occurrence of event A if possible, and the certainty N(A), denoting to what extend the occurrence of event A is certain, with:

) ( ) ( , A A S A ∏ ≤ Ν ⊆ ∀ (1) and ) ( 1 ) ( , A A S A Π - = Ν ⊆ ∀ (2)
If the occurrence of A is certain П (A)=1 and N(A)=1; if the occurrence of A is impossible: П (A)=0 and N(A)=0. If there is no knowledge available: П (A)=1 and N(A)=0 (the event is fully possible but not necessary at all).

Representation of imprecision

An imprecise information may be defined as v∈A where A is a subset of S which contains more than one element. The imprecision may be expressed by a disjunction of values [START_REF] Dubois | Représentations formelles de l'incertain et de l'imprécis[END_REF] defined by a possibility distribution on S. v∈A means that all values from v outside A are supposed to be impossible.

A possibility distribution π v of v quantifies the plausibility of the information v. π v is a function of S in L such as 1 ) ( , , ) ( , = ∃ ∈ ∈ ∀ s s and L s S s v v π π
with v denoting an ill-known value in S, and L the scale of plausibility ([0,1] for the theory of possibility).

Representation of uncertainty and imprecision

Information can be pervaded with both uncertainty and imprecision. Uncertainty and imprecision is denoted by v∈AUB where A and B are two disjoints subset of S, one expressing the possibility that the event does not occur, the other defining the (possibly imprecise) content of the event if it occurs. The possibility distribution associated to the information v∈AUB is noted π(v) with:

( ) ( ) ) ) ( max( ) ) ( max( B v v B A v v A ∈ = Π ∈ = Π π π (3)
These two events being mutually exclusive, they are characterised by relation ( 4

): 1 )) ( 1 ), ( max( )) ( ), ( max( )) ( ), ( max( = - Π = Π Π = Π Π A N A A A B A (4)

Basic Operations

In order to describe the events by possibility distributions, trapezoidal distributions (cf. Figure 1), denoted (a, b, c, d, h), can be used without an important loss of generality. Indeed, these sets intend to model an expertise which suggests a global shape rather than a precise function.

For performing a Material Requirement Planning, we need first to define an operator of addition, which is made according to [START_REF] Zadeh | Fuzzy sets as a basis for a theory of possibility[END_REF] and [START_REF] Dubois | Possibility Theory[END_REF]. The sum of two trapezoidal distributions A i and A j defined by the quintuplets A i (a i , b i , c i , d i , h i ) and A j (a j , b j , c j , d j , h j ) is defined by [START_REF] Dubois | Possibility Theory[END_REF]:

A i ⊕ A j = (a, b, c, d, h) (5) with: ( ) d d d b b b c c c a a a h d h d h d h c h c h c h h h j i j i j i j i j j i i j j i i j i - + + + = + - - + =         + ⋅ =         + ⋅ = = , min If n A A A A ∪ ∪ = ... 2 1 and p B B B B ∪ ∪ = ... 2 1
with A i and B j defined by the quintuplets (a i , b i , c i , d i , h i ) and (a j , b j , c j , d j , h j )

p j n i j i B A B A .. 1 , .. 1 = = ⊕ = ⊕ ∪ (6)
Calculating gross requirements consists in allocating quantities of components to periods. If the quantity is imprecise or uncertain, we need to compute the possibility and the necessity that a set belongs to an interval.

The possibility of the event "x∈F ", denoted Π(x∈F), is the degree of intersection between A and F , defined by equation [START_REF] References Aliev | Fuzzy-genetic approach to aggregate production-distribution planning in supply chain management[END_REF]. The dual measure of necessity of "x∈F ", denoted N(x∈F), evaluates to what extent "x∈F " is certainly true. It is defined by equation (8).

)

)) ( ), ( (min( sup ) ( u u F x F A u µ µ = ∈ Π (7)
)

) ( 1 )) ( ), ( 1 (max( inf ) ( c F A u F x u u F x ∈ Π - = - = ∈ Ν µ µ (8)
where F c denotes the fuzzy complement of

F ) 1 ( C F F µ µ - = Quantity Possibility a h=1 b c d
If x is a real variable and A a fuzzy interval, tmin and tmax represent respectively the minimum and the maximum limit of the interval (see Figure 2). The possibility that x is between the limits of the interval is given by equation ( 9):

))

( ( sup ) ( max) min ( max min u t x t x t A t u t µ ≤ ≤ = ∈ Π = ≤ ≤ Π (9)
On the opposite, the possibility that x does not belong to the interval is given by equation ( 10

): ))) ( ( sup )), ( ( sup max( )) max, [ ( min]), , ( ( max( ) ( max min A u t x t x t x u t A t u µ µ ≤ ≤ = +∞ ∈ Π -∞ ∈ Π = ∉ Π (10)
The necessity that x belongs to the interval is given by (11):

) ( 1 ) ( max) min ( t x t x t x t ∉ Π - = ∈ Ν = ≤ ≤ Ν (11)
The results of equations ( 9) to ( 11) are illustrated in Figure 2.

CALCULATION OF THE GROSS REQUIREMENTS

Figure 3 illustrates the process allowing to compute the gross requirements. Ill-known requirements (coming from forecasts or orders) are firstly considered as input of the process ((1) in Figure 3). These requirements are modelled through possibility distributions (step (2)). The result of this step is that each requirement is now modelled by a fuzzy quantity on a fuzzy date (3). The possibility and necessity that each quantity belongs to a period of the plan is then computed (4). Since a quantity may have non null possibilities to belong to several periods, it is then possible to define various scenarios for each requirement (depending on the period on which its occurrence is considered) (5). The decision maker will then suggest criteria for extracting one or several scenarios (most possible, optimistic, pessimistic, etc.) (6). These criteria are used to filter the possible scenarios [START_REF] References Aliev | Fuzzy-genetic approach to aggregate production-distribution planning in supply chain management[END_REF]. The decision maker can then make a decision on the base of the provided information.

) ( t x ∈ Π min t min t min t max t max t max t A µ - 1 A µ ) ( t x ∈ Ν A µ Figure 3.

Flowchart of the computation of gross requirements

These steps will now be explained with more details.

Representation of ill-known requirements

A requirement o is characterised by a quantity Qo and a date (due date) τo.

4.1.1

Types of imperfections

The imperfection on the requirements may be on the requirement itself, i.e. the requirement can be uncertain. For example, if the requirement comes from forecasts, it can be confirmed or cancelled through time. The imperfection may also affect the attributes of the requirement: the quantity or the

General model of imprecision on the requirements

It is possible to define the imprecision of a quantity by two intervals: one defining the most possible values, and the second all the possible values [START_REF] Grabot | Integration of uncertain and imprecise orders in the MRP method[END_REF] or to infer a possibility distribution from very few measurements [START_REF] Mauris | Inferring a Possibility Disdribution from Very few Measurements[END_REF], then to compute the trapezoidal fuzzy number which approximates this possibility distribution [START_REF] Grzegorzewski | Trapezoidal Approximations of Fuzzy Numbers preserving the expected interval: Algorithms and properties[END_REF]. In order to be consistent with the expert knowledge available in the companies, we suggest to represent the imperfection (uncertainty and imprecision) on the requirements using two fuzzy sets:

-the first fuzzy set represents the imperfection (uncertainty and imprecision) on the required quantity, and the uncertainty of the occurrence of the requirement itself (possibility that the required quantity is 0), -the second fuzzy set represents the imperfection on the due date.

The considered fuzzy sets are so defined as unions of trapezoidal fuzzy sets, which preserves the generality of the representation. In the following figures are introduced the various types of elementary sets used to represent the requirement pervaded by different types of imperfection. sets the constraint that max(π(A))=1 (at least one event has to be completely possible), which is not the case for these sets if h and h' ≠1. In order to define an imprecise requirement which cancellation is more plausible than its confirmation, the union of the fuzzy sets of Figure 5 and 6 has to be defined, with (for instance) h=0.8 and h'=1.

i i i i i i ∪ = ∪ 1 ) ' , ( 1 max = = h h n j j
Table 2 summarises the way the different models can be combined in order to define the imperfections on the quantities, qa, qb... being the values allowing to define the fuzzy number which represents the quantities.

= τ τ τ 1 ) ( 1 max = = j h n j τ ∪ n i i i i i i h d c b a 1 ) ; ; ; ; ( = τ τ τ τ τ 1 ) ( 1 max = = j h n j τ

Material requirement planning

Models of requirements based on unions of two sets have been introduced. One of the sets represents the imperfection on the quantity and the uncertainty on the requirement, whereas the other represents the imperfection on the date. The problem is now to compute a gross requirement, which is a set of quantities expressed for each period of the horizon.

Different scenarios can be deduced from the elementary requirements expressed as in Tables 2 and3, since a requirement may belong to different periods (according to the imprecision on its date) and may be composed of different quantities (according to its uncertainty and imprecision).

To compute the gross requirements, let us define the following notations:

4.2.1

Design of requirement scenarios

For each requirement, the number of elementary scenarios is equal to the number of periods to which the requirement may belong. As a first step, the possibility for each requirement o to belong to each period t is calculated (equation 12). This means that the requirement exists in a given period t with a given level of possibility: [0;[START_REF] References Aliev | Fuzzy-genetic approach to aggregate production-distribution planning in supply chain management[END_REF]period 2= [7;14[ and period 3= [14;21[. Figure 7. Representation

of A τ ( ) 1 [ 7 ; 0 [ ) ( max ) 1 ( = ∈ = ∈ Π A A A τ τ π τ ( ) 5 . 0 [ 14 ; 7 [ ) ( max ) 2 ( = ∈ = ∈ Π A A A τ τ π τ ( ) 0 [ 21 ; 14 [ ) ( max ) 3 ( = ∈ = ∈ Π A A A τ τ π τ So { } 2 ; 1 = A T
The requirement belongs to period 1 with a degree of possibility 1, and to period 2 with a degree of possibility 0.5 (see Figure 7).

In the third step, the required quantity of the different scenarios is calculated. It can be noticed that the possibility degree that a requirement o belongs to period t (

) ( o o u ∈ Π τ
) expresses than the requirement can exist in this period, but also that it cannot exist in other periods of o T .

U { } ao o ao o u v o ao o ao U v U B U with U o ao - ∈ ∀ ∈ Π = ⇓ ∈ Π ⊆ T )) ( ; 0 ; 0 ; 0 ; 0 ( ) ( , , τ τ 1 = o j )) ( ; ; ; ; ( 1 , 1 , 1 , 1 , 1 , , 1 , ao o A A A A A U o u h d c b a B o o o o o ao ∈ Π × = τ …. o o J j = )) ( ; ; ; ; ( , , , , , , , ao o A A A A A U J o u h d c b a B o j o o j o o j o o j o o j o ao o ∈ Π × = τ … zo U .. ) ( , , = ⇓ ∈ Π ⊆ zo U v o zo o zo B U with U o τ 1 = o j )) ( ; ; ; ; ( 1 , 1 , 1 , 1 , 1 , , 1 , ao o A A A A A U o u h d c b a B o o o o o ao ∈ Π × = τ …. o o J j = )) ( ; ; ; ; ( , , , , , , , ao o A A A A A U J o u h d c b a B o j o o j o o j o o j o o j o ao o ∈ Π × = τ
Example: let us consider the same requirement A

Computation of the gross requirements

From the various scenarios, a set of gross requirements, which is a set of required quantities for each period on a planning horizon, can be computed.

From Table 4, we generate all the combinations of the requirement scenarios. We have

∏ = = O o ox S C 1 combinations. Let C be the set of these C combinations.
For a combination c of the requirement scenarios, the total required quantity D c,t for a period t is the sum of the required quantities

} )) ( ; ; ; ; ( , ] , 1 [ )) ( ; 0 ; 0 ; 0 ; 0 ( , , , , , , , , , , , , , , , , , , , , o o
A A A A A u j o o o o o o o u v o o u t o t u j o o u j o t u j o o o u h d c b a B J j u v u B with t u B E t u B E t w o j o o j o o j o o j o o j o o o o o o o o o o o ∈ Π × = ∈ ∀ - ∈ ∀ ∈ Π = ≠ ∀ = = ∀ = ∈ ∀ ∀ τ τ T T (14)
In order to manage the combinatory explosion of elementary scenarios, the decision maker may select a plausible set of combinations and thus may only consider the combinations for which the possibility is higher than a given value p.

The set of combination is then reduced to C p (see equation ( 15)). We compute the gross requirement of the combination c and the total required quantities for the period t of this combination:

{ } { } { } ] , 1
c t E D t u j o O A u j w o t c c o c o c c o c o c o ∀ ∀ = ⊕ ∈ = , ) ( , , , ) ; ( / , , , , , , (16)
))

( ( min _ ,t c t c D comb Pos Π = (17)
The combinatory is directly linked to the level of imperfection of the requirements, and to the threshold introduced by the decision maker.

Example: let us consider the same requirement A

The set of scenarios is { } { }

) 2 ; 2 ( ); 1 ; 2 ( ); 2 ; 1 ( ); 1 ; 1 ( = A w
We compute the quantity for period {1; 2}of the scenario quantity w A (in this case

E A,j,u,1 and E A,j,u,2 ) 2 , 2 , 2 , 2 , 2 , ; 2 , 1 , 2 , 2 , 1 , 1 , 2 , 1 , 1 , 2 , 1 , 1 , 1 , 1 , 1 , ; ; A A A A A A A A B E B E B E B E = = = = 2 , 1 , 1 , 2 , 2 , ; 2 , 1 , 1 , 2 , 1 , 1 , 2 , 2 , 1 , 2 , 1 , 2 , 2 , 1 , 1 , ; ; A A A A A A A A B E B E B E B E = = = =
Then, we compute the gross requirements:

) 1 ; 0 ; 0 ; 0 ; 0 ( ); 1 ; 2 ; 2 ; 12 ; 10 ( 2 , 1 , 1 , 2 , 1 1 , 1 , 1 , 1 , 1 = = = = A A E D E D ) 5
= = = = A A E D E D

ILLUSTRATIVE EXAMPLE

In this example, we consider four requirements concerning material C. Each requirement has its origin in four different products P1, P2, P3 and P4, C being one of their components. These requirements are imperfectly known. The time horizon is 28 days and starts at day 0. Each period is composed of 7 days.

The gross requirements of material C has to be computed.

The following qualitative information on each requirement is available (depending for instance on the schedule of assembly):

Requirement 1:

-Two quantities are possible 300 (very possible) and 600 (could be possible). These quantities are given at ± 20 units -The requirement date is 20. This date is given at ± 2 days Requirement 2:

-The most possible quantity is between 400 and 450, but the requirement is not totally sure. This quantity is given at ± 10 units -The requirement date is between 14 and 16 and this date is given at ±2 days Requirement 3:

-The most possible quantity required is between 350 and 360. This quantity is given at -The most possible requirement quantity is between 600 and 620. This quantity is given at ± 20 units -The most possible date is between 25 and 26; it is possible to need this quantity up to two days before, but not later than the most possible date.

These pieces of information are represented in Figures 8 and9.

Possibilistic model of the ill-known requirements

The first step of the method is to build the possibilistic models describing each requirement, on the base of their "expert" description. The result, obtained according to the principles described in section 4, are given for Requirements 1 and 2 in Table 6 and7 (300; 300; 20; 20; 1)U (600; 600; 20; 20; 0.5) (20; 20; 2; 2; 1) r 2 (400; 450; 10; 10; 1)U (0; 0; 0; 0; 0.5) (14; 16; 2; 2; 1) r 3 (350; 360; 10; 10;1) (5; 5; 1; 1; 1)U (15; 15; 3; 3; 0.5) r 4

(600; 620; 20; 20; 1) (25; 26; 2; 0; 1)

Material Requirement Planning

5.2.1

Design of requirement scenarios

The next step consists in building the scenarios according to the requirements belonging to the different time periods of the gross requirement horizon. o T being the set of the possible periods for requirement o, the level of possibility that each requirement belongs to the different periods has to be computed (Table 9). Now, let us consider that the four requirements have been considered as crisp for making possible a "classical" gross requirement calculation. In that case, the four requirements will be certainly considered as follows (see summary in Table 16):

Requirement 1: quantity 300, date 20; Requirement 2: quantity 425, date 15; Requirement 3: 355, date 5; Requirement 4: 610, date 25. Thus the crisp gross requirement is given in Table 16. It can be seen that the "crisp" scenario of Table 16, obtained by neglecting any knowledge on the imprecision of the input information, is close to combination 2. Therefore, a first result is that our method is consistent with a "good sense" defuzzification obtained by considering only the most plausible cases. What is more important is that another scenario, obtained with combination 1, has the same possibility (=1): as a matter of fact, there is a possibility equal to one that some requirements are also present in period 2, and as a consequence that the gross requirements are relatively balanced on the four considered periods. In addition, it can be seen that there are slight possibilities (0.33 and 0.25) that the gross requirements are completely unbalanced, with a quantity around 350 during period 2 and a huge quantity, around 1000, during period 4.

It is so clear that the method may allow to obtain a much more informative set of possible situations, together with an estimation of their possibility of occurrence. On that base, the decision maker is able to make a really "informed" decision, considering not only the most plausible case but also an estimation of the risk taken regarding to the possibility of occurrence of other situations.

CONCLUSION

Through a literature review, we have shown that the idea to model an imprecise demand on final products using possibility theory has sometimes been considered, but that the practical way for dealing with such a demand in the first steps of the MRP method (Material Requirement Planning) has not yet been fully explored. However, this step is the base for generating a supply plan for the upstream partners of the supply chain, which should take benefit from more precise information.

In this paper, uncertainty and imprecision on the requirements have been analysed and a model allowing to take them into account has been proposed. Based on this model, a method has been suggested for modelling uncertain gross requirements (describing plausible quantities by period). The result stresses that the resulting uncertain quantity which is required is not a simple fuzzy set, but describes a set of situations, which is much more informative than the result obtained by considering only the most plausible requirements. Indeed, this model does not lose any information on the relation between period and quantity.

As a following step, two different options are possible: the decision maker may choose a supply plan after assessing its feasibility, and may compute the risk taken (shortage or excess of inventories) if other situations occur. In order to do this, he can defuzzify the gross requirement before sending it to his supplier. The second solution is to transmit the set of the most plausible supply plans to the supplier as an element for negotiation. Indeed, the described situations are a good starting point for assessing how the risk taken by preparing a situation more than another can be shared between the two partners.

In this article, the incompleteness of the requirements has not been considered. Considering this other type of imperfection is an interesting perspective for the study, since the knowledge on the requirements may be incomplete for far periods. In that case, any evaluation of the required quantity linked to forecasts could be false. An incompleteness factor could be added to the requirements already known through long term contracts in order to address this problem.

Finally, fuzzy demands can be considered in other cases, such as economic production quantity models (Lee and [START_REF] Yao | Economic production quantity for fuzzy demand quantity and fuzzy production quantity[END_REF], inventory control [START_REF] Petrovic | Coordinated control of distribution supply chains in the present of fuzzy customer demand[END_REF] or optimal inventory (Chiang and Wenkai, 2000), also requiring a demand defined by quantities by period. As a perspective, the integration of our model in those methods could be another interesting issue.

Figure 1 .

 1 Figure 1. Trapezoidal distribution of possibility

Figure 2 .

 2 Figure 2. Computation of possibility and necessity levels that an element belongs to an interval It will be shown in next section how these basic operations allow to build a possible gross requirement plan.

Figure 4

 4 Figure4represents a precise value (a), of possibility h. Figure5shows an imprecise quantity with a maximum possibility equal to h. Figure6represents the uncertainty on the occurrence of a requirement; the value zero (which means that the requirement is cancelled) has a possibility level equal to h'.

Figure 4 .Figure 5 .

 45 Figure 4. Representation of (a; a; 0; 0; h)

  Index o: requirement with o∈[1,O] t: period with t ∈[1,T] j o : index of fuzzy set of requirement o with j o ∈[1,J o ] u o : period in which the requirement o could exist with u o ∈[ U ao ,U zo ] c: number of combination with c∈[1,C] Variables Q o : possible required quantity of requirement o, τ o : possible date of requirement o, period t of the fuzzy set j o of requirement o which belongs to period u which the fuzzy set indexed by j o of the requirement o belongs to period u o D c,t : gross requirement quantity for combination c for period t.

  step, for each requirement o, the set o T of the periods on which the requirement may occur is selected (equation 13u o ∈[ U ao ,U zo ] (13) Example: let us consider a requirement A on )

  the requirement, with j o the index of the fuzzy set (imprecise value); j o ∈[1,J o ]. The different scenarios o w of requirement o can be

  ) concerning all the scenarios of requirements w o of this combination belonging to this period t (equation 16).Thus, for each scenario w o and each period t in o requirement o belongs to the considered period t, the quantity is equal to indexed by j o of requirement o weighted by the possibility to exist in the period u o , and o belongs to another period than t (quantity=0, possibility denoting that the requirement belong to another period u o ).

{

  combination is characterised by a level of possibility (equation 17). The combinations computed are thus classified according to their degree of possibility c comb Pos _ .

o

  The most possible: around 5 (given at ±1 days) o The less possible: around 15 (given at ±3 days) Requirement 4:

Figure 8 .

 8 Figure 8. Information on quantities from the expert

Figure 12 .

 12 Figure 10. Representation of the quantities for requirements 1 and 2

Table 1

 1

			. Classification of previous works
	References		Use of fuzziness	type of model
	(Mula et al., 2007)		Flexible constraints	planning
	(Jamalnia and Soukhakian,
	2009)
	(Fargier and Thierry, 2000)	Integration of uncertainty	planning
	(Grabot et al., 2005)
	(Niu and Dartnall, 2008)		Defuzzification	planning
	(Liang, 2008)		then optimisation
	(Peidro et al., 2009)
	(Lan et al., 2009)		Optimisation	planning
	(Sun and Liu, 2008)		under uncertainty
	(Aliev et al., 2007)
	(Chen and Huang, 2006)		Optimisation	scheduling
	(Balasubramanian	and	under uncertainty
	Grossmann, 2003)

Table 2 .

 2 Representation a of quantity

		Certain quantity															Uncertain Quantity
		Precise quantity	Imprecise quantity	Precise quantity	Imprecise quantity
	Certain requirement	(qa; qa; 0; 0; 1)	(qa; qb; qc; qd; 1)	n ∪ i 1 =	(	i qa ( 1 max i qa ; = n h ; j j	; 0 1 = ; 0 )	i qh	)	n ∪ i 1 =	(	i qb n 1 max i qa ; = j	; (	i qc ) ; j h	i 1 qd =	;	i qh	)
	Uncertain requirement	(qa; qa; 0; 0; qh)U (0; 0; 0; 0; h') max(h',qh)=1	max(h',qh)=1 (0; 0; 0; 0; h') (qa; qb; qc; qd; qh) U	; 0 ; qa ; 0 ( ( i ∪ 1 i ∪ = n	; 0 ; qa ; 0 i	' ; h 0	) 0 ;	qh i	)	n	1	; 0 ; qa ; 0 ( (	) qc ' h ; ; 0 qb ; 0	;	qd	;	qh	)
						n
				1 max = j	(	h	j	,	h	'	)	=	1

Table 3

 3 summarises the combinations of models which may represent the various imperfections on the dates. The values allowing to represent the fuzzy numbers modelling a date are denoted τa, τb...

	Possibility
	h
	a	b	Quantity or Date
	c	d
	Possibility
	h'
	Quantity or Date

Table 3 .

 3 Representation of the requirement date

		Certain requirement date								Uncertain requirement date
		Precise date	Imprecise date			Precise date		Imprecise date
	Certain requirement	(τa; τa; 0; 0; 1)	(τa; τb; τc; τd; 1)	n ∪ i 1 =	i h ; τ a τ ( 1 ; max i a τ ( = n j	; 0 ) j	i h τ 1 ; 0 =	)	n ∪ i 1 =	i ) c τ j h ; τ i b τ ( 1 max a ; τ ( i = n j	d τ 1 = ;	i	;	i h τ	)
	Uncertain requirement	(τa; τa; 0; 0; τh) max(h';τh)=1	(τa; τb; τc; τd; τh) max(h'; τh)=1	n 1 ∪ i	(	a	i	;	a	i	;	; 0	; 0	i h	)

Table 4

 4 represents the values of the required quantities for each period, by taking into account that the requirement can exist in one and only one period of the horizon.

	U o being the number of periods of o T , we have	S	o	=	U	o	×	J	o	scenarios for each requirements o.
	Let us denote	S	o	,	j	o	,	u	o	=	∪	o E	,	j	0	,	u	o	,	t
											t	o T ∈

Table 4 .

 4 Representation of the possible quantities Quantity of requirement in other periods than u o

	Quantity of requirement in period u o
	ao

Table 5 .

 5 Results for the requirement A

	A r	⊆	1										B	A	1 , 1 ,	=	12 ; 10 (	;	2	;	2	) 1 ;
	Q	A	∈	A T	≠	2					B	A	,	2	1 ,	=	(	; 20	; 24	; 4	; 4	0	) 8 .
	B	A	,	2	1 ,	=	(	; 0	; 0	; 0	) 1 ; 0
	A r		⊆	2										B	A	, 1 ,	2	=	; 10 (	; 12	; 2	; 2	. 0	) 5
	Q	A	∈	A T		≠	1					B	A	,	2	,	2	=	(	; 20	; 24	; 4	; 4	) 4 . 0
	B	A	, 1 ,	2	=	(	; 0	; 0	; 0	; 0	) 5 . 0
	4.2.2

Table 7

 7

		. Model for Requirement 2
	Expression	Model		Representation
	Requirement can be cancelled.	Possibility to					Requirement 2
	The quantity belongs to an interval	be cancelled = 0.5		1,2
		Possibility to		1
		need 400=1; Possibility to need 450=1	possibility	0,4 0,6 0,8
				0,2
				0
				0	50	100	150	200	250	300	350	400	450	500
								quantity
	The quantities are given with	Trapezoidal					Requirement 2
	a tolerance	fuzzy number
		is used.		0,9 1
				0,8
			Possibility	0,3 0,4 0,5 0,6 0,7
				0,2
				0,1
				0
				0	50	100	150	200	250	300	350	400	450	500
									quantity
	The date is belong a interval	Trapezoidal					Requirement 2
	and is given with a tolerance	fuzzy number
		is used.		0,9 1
				0,8
			Possibility	0,3 0,4 0,5 0,6 0,7
				0,2
				0,1
				0
				0	2	4	6	8	10	12	14	16	18	20
									Dates

  Table8summarizes the model whereas Figures 10 and 11 represent graphically the models of the required quantities. Figure12summarizes the models of the requirements dates (black lines represent the periods).

	Table 8. Requirements
	Quantity	Requirement date
	r 1

Table 16

 16

		. "Crisp" gross requirement
	period	1	2	3	4
	quantity	355	0	725	610

Acknowledgements: This study has been performed with the support of Région Midi-Pyrénées and of the University of Toulouse.

We compute for each scenario w A the quantity

From Table 9, we determine the set of possible periods for each requirement

Then, we can build the different scenarios for each requirement (Tables 10 to 13):

Table 10. Results for the requirement 1

Table 12. Results for the requirement 3

5.2.2

Computation of the gross requirements

From Tables 10 to 13, the different gross requirements are obtained for some combinations of the requirement scenarios (Figure 13). There are 48 possible gross requirements, but only the gross requirements with level of possibility 1 are computed, so that the less possible gross requirement which could happen in the worst case: possibility = 0.25 (Table 14 andTable 15).