Philippe Drobinski
email: philippe.drobinski@aero.jussieu.fr

Christian Haeberli

Evelyne Richard

Marie Lothon

Alain M Dabas

Pierre H Flamant

Markus Furger

Reinhold Steinacker

Scale interaction processes during the MAP IOP 12 south föhn event in the Rhine Valley

Keywords:

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

INTRODUCTION

Föhn is a generic term for strong downslope winds experiencing warming and heating in the lee of a mountain ridge. F öhn is linked with a decrease of cloudiness in the lee and is strong and gusty and often channelled in the numerous gaps and valleys of the ridge [START_REF] Brinkmann | 1974 Strong downslope winds at Boulder[END_REF][START_REF] Seibert | South foehn studies since the ALPEX experiment[END_REF]). F öhn may cause damage due to severe storms (Brinkmann 1974), snow melting (Hoinka 1985), or high pollution levels [START_REF] Nkemdirim | The effect of chinook weather on urban heat islands and air pollution[END_REF][START_REF] Hoinka | The surface layer on the leeside of the Alps during foehn[END_REF] and it is of course of great importance to better understand and predict the structure of the föhn ow.

The favourable synoptic-scale feature leading to south föhn north of the Alps is usually a south/south-westerly basic ow, in general preceding the passage of a cold front in the Alpine region [START_REF] Hoinka | 1985 Observations of the air ow over the Alps during a föhn event[END_REF]. A distinction is made between shallow föhn, which systematically starts typical föhn episodes [START_REF] Seibert | South foehn studies since the ALPEX experiment[END_REF], and deep föhn. In contrast to the deep föhn, the cross-Alpine southerly ow is essentially con ned to the valley transects, whereas the mid-tropospheric ow can have a westerly or even a north-westerly direction, which can be associated with an approaching cold front which ultimately terminates the föhn storm [START_REF] Egger | 1990a Thermally induced ow in valleys with tributaries. Part I: Response to heating[END_REF][START_REF] Heimann | Three-dimensional modeling of synthetic cold fronts interacting with northern Alpine foehn[END_REF].

The mesoscale processes that de ect and channel the westerly upper-level ow with the valleys during the shallow föhn phase have been investigated by [START_REF] Sprenger | Rotational aspects of strati ed gap ows and shallow föhn[END_REF] and Zängl (2002) and are similar to the counter-current theory (e.g. [START_REF] Wippermann | Air ow over and in broad valleys: Channelling and countercurrent[END_REF]). [START_REF] Sprenger | Rotational aspects of strati ed gap ows and shallow föhn[END_REF] also show that during shallow föhn, the complex three-dimensional Alpine ridge can give rise to ow separation and splitting in the vicinity of the valleys and gaps, and these may induce anomalies of potential vorticity and wake features to the lee. During deep föhn, the mechanisms causing the air to descend and displace the colder and denser air in the lee are closely connected to the existence of large-amplitude lee waves and/or hydraulic jump [START_REF] Brinkmann | 1974 Strong downslope winds at Boulder[END_REF]Hoinka 1985) which can eventually intensify the föhn ow along the gaps and valleys transecting the ridge (e.g. [START_REF] Zängl | Three-dimensional idealized simulations of the foehn in the region of Innsbruck[END_REF][START_REF] Flamant | Gap ow in an Alpine valley during a shallow south föhn event: Observations, numerical simulations and hydraulic analogue[END_REF]. The föhn penetrating down to the surface erodes the cold surface layer or cold pool, transported from polar regions and formed when an absence of clouds permits the ground to radiate freely towards space, up to 50 km north of the baseline [START_REF] Hoinka | The surface layer on the leeside of the Alps during foehn[END_REF].

In most of the aforementioned studies, the instrumental set-up did not allow an overall three-dimensional documentation of the föhn from the synoptic scale to the scale of the valleys, with suf cient spatial and time resolution to analyse the processes governing the evolution of föhn development and its spatial structure. During the Mesoscale Alpine Programme (MAP) eld experiment [START_REF] Bougeault | The MAP special observing period[END_REF], the target area for föhn studies was located in the Rhine Valley between the city of Chur and Lake Constance (Fig. 1).

The FORM ('FOehn in the Rhine Valley during MAP') programme was designed to improve the understanding and forecasting of the life cycle of föhn-related phenomena, including their three-dimensional structure and associated boundary-layer processes. In the framework of the FORM programme, Drobinski et al. (2001a) investigated theoretically and experimentally, the small-scale dynamical processes that explain how air blows in the Rhine Valley (hereafter, RV1 is used for the lower part of the Rhine Valley and RV2 for the upper part, see Fig. 1(b)) and the Seez Valley (hereafter called SV, see Fig. 1(b)), during föhn and in the absence of föhn. [START_REF] Drobinski | Statistical characterization of the ow structure in the Rhine valley[END_REF] used the rawinsonde observations acquired during the whole MAP Special Observing Period (SOP) to connect the small-scale circulation in the FORM target area with mechanical or thermal forcing at synoptic or local scales.

The present study investigates the spatial scale interactions during a south föhn event in the FORM target area between 29 October 1999 and 31 October 1999, i.e. during Intensive Observing Period (IOP) 12 of the MAP eld experiment. This south föhn event had the advantage of presenting clear separate development stages, each lasting several hours, with no intermittent breakdowns (unlike during IOPs 5, 8, 9 and 10) and could thus be reliably and accurately documented with time and space resolutions never achieved before, from its onset until its breakdown. This paper focuses on several points that in fact motivated the deployment of this unique and very dense instrumental set-up during MAP, to study the processes related to föhn in this region.

The key issues addressed in this study are:

²

The diurnal cycle of a south föhn event, i.e. the mechanism of temporal and spatial evolution and cessation of föhn in this complex valley system on a local scale.

²

The three-dimensional structure of the föhn ow from synoptic to local scale.

²

The interaction between low-level and mid-tropospheric föhn ow on the scale of the large Alpine valleys. This study makes use of the data collected in the FORM target area by the dense rawinsonde and surface station network, scintillometers and a ground-based Doppler lidar. At a larger scale, the rawinsonde and surface station networks provide the synoptic setting for föhn occurence, and constant-volume balloons provide information on the wave activity and the intensity of turbulence associated with the föhn event. In addition, mesoscale numerical simulations are used to investigate the role of the tributary valleys in the föhn ow structure in the FORM target area, and the coupling between the föhn ow and the wave activity.

Following the introduction in section 1, the instrumental set-up and the numerical model Meso-NH are presented in sections 2 and 3, respectively. Scale interactions are studied using a downscaling approach from the synoptic environment favourable to south föhn in the Rhine Valley (section 4) and the mesoscale processes forcing the ow at the Rhine Valley and tributary scale (section 5) to the development of the föhn ow at local scale in the FORM target area (section 6). This study is conducted for the different stages of development of the IOP 12 föhn episode.

INSTRUMENT SET-UP

(a) At the scale of the valley (i) Rawinsonde and surface station network. The valley ow structure in the FORM target area was investigated using rawinsonde stations in Buchs-Grabs, Diepoldsau, Heiligkreuz, Malans and Masein.

Surface measurements of pressure, temperature, humidity and sometimes wind speed and direction were recorded at the stations in Altenrhein, Chur, Quinten and Vaduz (Fig. 1) and were used to investigate the dynamical processes driving the ow structure in the FORM target area and to compute local pressure gradients.

(ii) Scintillometers. Two scintillometers with three sensors arranged in the shape of the letter L [START_REF] Poggio | Scintillometer wind measurements over complex terrain[END_REF] were located in RV2. The scintillometer transmitters were located in Triesenberg, the receivers in Flusa and Ergellen (see Fig. 1(b)) at about 500 m above ground level (AGL), i.e. 1000 m above sea level (ASL). Assuming a homogeneous wind eld in the whole area between the light beams, this V-shaped arrangement of the light beams allowed the simultaneous measurement of the vertical and horizontal components of the wind eld. This method was used to document the strong föhn wind storms of MAP IOP 8 to 10 (i.e. 20 to 24 October 1999) [START_REF] Furger | Comparison of horizontal and vertical scintillometer crosswinds during strong foehn with lidar and aircraft measurements[END_REF].

(iii) Transportable wind lidar. The transportable wind lidar (TWL) was located near Bad Ragaz in Vilters (see Fig. 1). The TWL provides wind radial-velocity measurements along the line of sight. During the MAP SOP, scans were performed by sweeping the azimuth angle from 262 ± (in the direction of RV1) to 142 ± (in the direction of SV), and the elevation from 4 ± to 28 ± with increments in elevation and azimuth of 2 ± and 5 ± , respectively (see Drobinski et al. 2001a for more details).

(b) At the synoptic scale (i) Rawinsonde and surface station network. The synoptic pressure gradients parallel and perpendicular to the Alpine ridge are computed from the surface measurements collected hourly in Lugano, Vaduz and Zurich. In addition, the evolution of the surface wind in the Reuss Valley at the Guetsch station (Fig. 1(a)) is analysed and the data collected by the rawinsondes launched from the operational meteorological station in Milano-Linate are used to document the upstream ow and stability on the southern part of the Alpine ridge.

(ii) Constant-volume balloons. During IOP 12, constant-volume balloons (CVBs) were launched from Ispra (Italy) located near the Lago Maggiore to document the dynamics of the isopycnic air ow (pressure, wind intensity and direction, temperature and humidity content) and particularly the characteristics of orographic waves and trapped lee waves (location, amplitude, period, intensity) (e.g. [START_REF] Kof | Trapped lee waves observed during PYREX by constant volume balloons: Comparison with Meso-NH simulations[END_REF].

(c) Surface-data analysis Pressure and temperature data from the surface stations (including those collocated with the upper-air stations) were interpolated at surface level to a regular grid with a spacing of 2 km using the Vienna Enhanced Resolution Analysis (VERA) algorithm. This algorithm is a further development of the method described in Steinacker et al. (2000) and has many similarities to the thin-plate spline method.

MODEL DESCRIPTION

The Meso-NH model solves the non-hydrostatic and anelastic equation system [START_REF] Lafore | The Meso-Nh atmospheric simulation system. Part I: Adiabatic formulation and control simulation[END_REF]. The turbulence scheme is uni-dimensional in the vertical and is based upon the physical mixing length of [START_REF] Bougeault | Parametrization of orography-induced turbulence in a meso-beta scale model[END_REF]. The convection and radiation schemes, as well as the parametrization of land surface processes are described in [START_REF] Bechtold | A mass ux convection scheme for regional and global models[END_REF]), Mocrette (1991) and [START_REF] Noilhan | A simple parameterization of land surface processes for meteorological models[END_REF], respectively.

A 36 h simulation initialized at 0600 UTC on 29 October 1999 was completed using stationary 2 km and 8 km resolution domains nested within a 32 km domain using twoway interfaces [START_REF] Stein | High-resolution non-hydrostatic simulations of ash-ood episodes with grid-nesting and ice phase parameterization[END_REF]. The initial and coupling elds were generated by rst interpolating the French analysis Action des Recherche Petite Échelle et Grande Échelle (ARPEGE) data available every 6 hours to the model grid. The vertical grid is made of 50 levels with a mesh stretched between 60 and 600 m. To ensure a good description of the atmospheric boundary layer, 12 levels are used below 1000 m. The size of the rst nine meshes is less than 100 m. Above 6400 m, the mesh size is constant and set to 500 m. An absorbing layer is set above 14 000 m. The top of the domain is located at 21 400 m.

In the following, Meso-NH simulations are validated and used to provide the threedimensional environment necessary for the interpretation of the measurements.

IOP 12 SYNOPTIC ENVIRONMENT

The föhn case of IOP 12 (30 October-31 October 1999) is also known as the 'Sahara dust föhn' of the MAP SOP. Ahead of a trough west of the Iberian peninsula on 29 October 1999, very warm air from north-western Africa was advected to southern France (Pau, foothill of the French Pyrenees, reported 29 ± C). A branch of this warm air also reached the regions north of the Alps. In Geneva (south-western Switzerland) large amounts of dust from the Sahara were registered on 29 October 1999 during the day. Also, at the MAP operations centre in Innsbruck a probe of Sahara dust was collected. Temperature maxima in the FORM target area reached a remarkable 20 ± C. The Iberian trough moved north-eastward and joined the strong westerly jet over the north-east Atlantic (Fig. 2).

The surface pressure eld was initially dominated by a large anticyclone with its centre over eastern Europe. This anticyclone moved further to the south-east. This movement, together with the trough moving towards the European continent, caused the large-scale mid-tropospheric ow to back from the west on 29 October in the afternoon to the south-west during the evening. In the same time interval, the cross-Alpine pressure difference increased to about 6 hPa within 6 hours (Fig. 3). This strong ampli cation of the pressure contrast is called the 'föhn knee' [START_REF] Brinkmann | 1974 Strong downslope winds at Boulder[END_REF]. The situation did not change much during the second half of the night. Until noon on 30 October 1999 the mid-tropospheric ow was reinforced and the cross-Alpine pressure difference increased to 8.5 hPa around 1200 UTC. At about the same time, föhn started. It ended when the cold front associated with the Iberian trough reached the target area shortly after midnight of 30/31 October 1999.

Precipitation during the whole episode was weak and spatially inhomogeneous. Very light precipitation in the northern parts of the target area started after the passage of the cold front after midnight of 30/31 October 1999. The rst rain observations are from Diepoldsau and Heiligkreuz in the FORM target area (Fig. 1

(b)) at 0100 UTC.

There is a considerable annual cycle in the characteristics of föhn ow. The mean cross-Alpine pressure difference (Vaduz-Lugano) during the föhn cases in September and early October 1999 varied between ¡3:1 hPa and ¡6:6 hPa, but increased to values between ¡4:9 hPa and ¡9:3 hPa during the cases between mid October and mid November. For a comprehensive comparison of all föhn cases observed during the MAP SOP, see [START_REF] Steinacker | Unstationary aspects of foehn in a large valley: Scienti c objectives and description of the foehn cases during the MAP subprogram FORM[END_REF]. To take into account the annual cycle, the comparison of the IOP 12 föhn case with other föhn events of the MAP SOP is restricted to those ve cases which occurred between mid October and mid November. With respect to the cross-Alpine pressure difference, the IOP 12 case was not as strong as the three preceding cases. But it was the case with the highest mean potential-temperature difference across the Alps (9.0 K). Although there are only a few simultaneous observations, there is a clear indication of a pronounced negative height difference at 700 hPa between Diepoldsau (situated to the north of the FORM target area) and Milano. Comparably strong differences were observed during only one other case (IOP 5) during the SOP.

During IOP 12, the ow direction at the Alpine crest was very similar to during the other IOPs, but the wind speed was up to 10 m s ¡1 lower than during the strongest cases.

DYNAMICAL PROCESSES OF THE F ÖHN LIFE CYCLE AT THE MESOSCALE

The föhn touches the ground in the Rhine Valley when the cross-Alpine pressure difference is larger than about 6 hPa (Fig. 3), i.e. on 30 October 1999 at about 0900 UTC, and lasts until 31 October 1999 at 0310 UTC. This pressure difference value is comparable to that observed during IOPs 2, 8, 9 and 10 and indicates that the blocking effect and the heating due the synoptic situation and föhn are strong [START_REF] Steinacker | Unstationary aspects of foehn in a large valley: Scienti c objectives and description of the foehn cases during the MAP subprogram FORM[END_REF]. In fact, substantial heating and drying (exceeding föhn criteria) is observed in Vaduz (Figs. 4(c) and(d)) in the FORM target area (Fig. 1(b)). This penetrating föhn period is associated with a south/south-westerly upper-level ow, channelled in the north-south oriented parts of the Rhine Valley (Fig. 4(b)). Preceding the penetrating föhn, a shallow föhn phase is observed in the night between 29 and 30 October 1999, during which a strong cross-Alpine pressure gradient allows the de ection of the westerly upper-level ow within the Rhine Valley (Fig. 4(b), [START_REF] Seibert | South foehn studies since the ALPEX experiment[END_REF]). However, katabatic drainage ow from the longitudinal (i.e. north-south oriented) and transverse (i.e. east-west oriented) tributaries of the Rhine Valley, prevents the föhn penetrating down to the surface of the Rhine Valley (Fig. 4, see section 6). The term 'sandwich föhn' might thus be appropriate to describe this ow (Vergeiner and Mayr 2000).

One purpose of this section is to study the upper-level ow de ection mechanism during the shallow f öhn. In fact, the previous studies on the channelling effect of the föhn ow by major Alpine valleys considered north-south oriented portions of Alpine valleys. In the FORM target area, the Rhine Valley has a very complex shape: below 46.8 ± N, the Rhine Valley is oriented west-east while above 46.8 ± N, it has a more southnorth orientation up to Lake Constance. Several major valleys merge with the Rhine Valley, these are the Val Lumnezia, the Domleschg, the Praettigau, the Seez Valley and the Walgau. Also, at the bifurcation between the Rhine and Seez Valleys, the Tamina gorge and the Weisstannen Valley play an important role. In other words, the Rhine Valley cannot be considered as a transect in the Alpine ridge and all these tributaries should play a role in the circulation of the f öhn ow in the FORM target area. As a consequence, a second objective in this section is to examine the role of the Alpine valley network in directing the föhn ow towards the FORM target area from its early stage of development (shallow föhn) until its breakdown.

(a) Sandwich shallow föhn On the windward side, the Po Valley and the valleys in the southern Alps are lled with a stably strati ed, relatively cool air mass. The stability is concentrated in an inversion layer between 2000 m and 3000 m. In the Po Valley a pronounced easterly ow is observed below the inversion providing evidence of blocking by the barrier (Fig. 5). In fact, due to the topography of the Alps and the Appennines, a convergence is to be expected in the region of Milano [START_REF] Seibert | South foehn studies since the ALPEX experiment[END_REF]). The air mass below the inversion remains in the south and is not part of the shallow föhn ow. The controlling parameters of upstream blocking are the Froude number, F r D U=N h, and the cross-mountain Rossby number, Ro D U=f L h [START_REF] Pierrehumbert | Upstream effects of mesoscale mountains[END_REF] (U is the mean upstream wind, N the Brunt-Väisälä frequency, h the ridge height, L h the half-width of the mountain, and f the Coriolis parameter). The parameters F r and Ro are found to be equal to 0.26 and 1.62 for h D 2500 m and L h D 80 km, respectively (Table 1), which predict strong blocking, the Coriolis force limiting the upstream extent of the decelerated layer. In the northern Alpine valleys, the air shows little motion during the shallow föhn (see Figs. 4 and6). Figure 6 shows a strong decoupling between lower-and upperlevel ows. This is clearly visible with a highly stable strati cation downstream, a layered structure of the humidity eld (not shown) and very different wind speeds and directions between the higher vertical levels of the atmosphere and the channelled ow. The observations show that on 29 October 1999 after 1700 UTC, the synoptic ow blows from the west/south-west while the ow within the Rhine Valley blows from the south in RV1 (Fig. 6).

The reason for the de ection of the westerly upper-level ow is a deceleration of the ow within the main longitudinal tributaries (i.e. the Reuss Valley, Val Lumnezia and the Domleschg, see Fig. 1(a)), due to the cross-section through which the ow has to stream. The retardation disturbs the equilibrium between the horizontal pressure gradient and the Coriolis force by diminishing the latter. Consequently, the ow is de ected in the The mean lapse rate dµ =dz, the Brunt-Väisälä frequency N , and the mean upstream wind U are obtained from the upstream sounding in Milano.

to about 1200 UTC that veers slowly to the south and stabilizes at about 1800 UTC. The wind speed is nearly constant and equal to 5 m s ¡1 and starts to intensify at 1800 UTC. The Meso-NH simulations catch with good accuracy the development of the southerly ow in the Reuss Valley but underestimate the wind intensity and generate large discrepancies in the wind direction before 1800 UTC on 29 October 1999. Horizontal cross-sections of Meso-NH horizontal wind elds show that near the surface, the ow at Chur is fed by the katabatic ow from south-western RV1 (Fig. 7(a)). At higher levels (Figs. 7(b) and (c)), the air from the Domleschg partly goes straight to the north over the Kunkelspass (which is only 700 m AGL, i.e. is about the same height as the crest of the gorge near Rothenbrunnen inside the Domleschg), and partly towards Chur in RV1. The reason why Val Lumnezia does not de ect signi cantly the westerly upper-level ow is the short length of the valley (Table 1) since the maximum meridional velocity behaves as L 1=2 where L is the length of the valley [START_REF] Sprenger | Rotational aspects of strati ed gap ows and shallow föhn[END_REF]. At 1500 m ASL, katabatic drainage ow blows from the transverse Praettigau. One can note a strong intensi cation of the ow in RV2 near Feldkirch where the ow originating from the splitting of the ow at the bifurcation of the Rhine and Seez Valleys merges with ow blowing from the east-west oriented Walgau. At the bifurcation between the Rhine and Seez Valleys, the ow comes from the Tamina gorge and RV1 (the ow rst splits between these two valleys near Kunkelspass) and splits between SV and RV2.

The detailed three-dimensional ow structure shown in Fig. 7 is in good agreement with the TWL observations (Fig. 10). The upper-level ow backs slowly from the west (Fig. 10(h)) to the south-west (Fig. 10(l)). At 1700 UTC, evidence of ow splitting is observed at 1000 m ASL with air blowing away from the TWL (hot colours) in SV and RV2. A strong transverse radial-velocity shear is visible in SV with a near-wall jet (about 8 m s ¡1 ) at the northern wall of SV, while the wind speed is very weak at the southern wall of SV. This could be due to wall separation and drainage out ow from the Weisstannen Valley (Fig. 7) and could explain why people living along the southern slope of SV report unexpectedly weak winds during föhn. At 1500 m ASL, the air is almost fully de ected in RV2, while a north-westerly ow blows in SV with a very thin layer of weak south-easterly ow at the northern sidewall of SV. At 2000 UTC, Fig. 10 (third row) shows that the layer where ow splitting occurs between the Rhine and Seez Valleys deepens as ow splitting is visible at 1000 m and 1500 m ASL. The valley ow and the ow above the crest are now clearly decoupled. [START_REF] Drobinski | Statistical characterization of the ow structure in the Rhine valley[END_REF] showed that ow splitting between the Rhine and Seez Valleys occurs during föhn events or at night. Figure 7 shows that westerly/north-westerly ow blowing from the Lake of Zurich and from the Klön Valley near Glarus (strong thermal gradient at this altitude in these valleys) prevents the air coming from RV1 to penetrate further in SV.

To provide some insight into the dynamics of the ow in the longitudinal tributaries of the Rhine Valley, the down-valley Rossby number R l D v=f l and the cross-valley Rossby number R L D l 2 =L 2 R l (where v is the along-valley wind speed, l and L the cross-and along-valley length-scales and f is the Coriolis parameter) are computed from Table 1. In the three valleys, R L ¿ 1 suggesting that the ow is in geostrophic balance in the valley, i.e. that the pressure gradient normal to the valley axis is balancing the Coriolis force associated with the along-valley component of the winds. As the largescale pressure gradient is oriented along the axis of the three valleys, then as R l > 1, it will be balanced by the ageostrophic acceleration of the ow down the axis of the valley and gap winds will result [START_REF] Overland | Scale analysis of marine winds in straits and along mountainous coasts[END_REF]. As the three valleys are narrow, i.e. l is lower than the Rossby radius of deformation l R D N 1H =f (1H is the depth of the valley and N the Brunt-Väisälä frequency), they favour rapid transitions and hydraulictype responses near the valley exits. The dimensionless valley depths g 1H D N 1H =U , which discriminate the nature (linear versus nonlinear) of the ow regime in the valley are close to the critical dimensionless mountain height ( g 1H » O.1/) and thus substantial nonlinear effects are expected (such as the hydraulic jump) [START_REF] Sprenger | Rotational aspects of strati ed gap ows and shallow föhn[END_REF].

(b) Penetrating föhn in the Rhine Valley

On 30 October 1999 at 1200 UTC, when föhn reaches the Rhine Valley, the upstream upper-level ow veers to the south/south-west and stability decreases as wider spacing between the isentropes is visible in Fig. 5. The vertical extension of the blocking region also narrows as the Froude number increases (F r D 0:44, see Table 1) [START_REF] Pierrehumbert | Upstream effects of mesoscale mountains[END_REF]. The Rossby number Ro also increases to 2.12, revealing persistent transience after the upstream ow initially impinges on the ridge. The transient upstream ow blocked by the ridge nally veers to the west with the föhn breakdown (Fig. 5). Above the inversion and the blocking area, Fig. 5 displays a stronger south-westerly ow between 2000 and 3500 m, veering westerly aloft with altitude. The cold-air stagnation blocked by the mountain range on the windward side is essential as it explains part of the strong cross-Alpine pressure gradient that gives birth to the föhn (Hoinka 1985).

Figure 11 shows a descending air mass in the lee: the air reaching the ground at the northern edge of the Alps originates from a level of 2000 m to 3500 m in the south (see isentropes in Fig. 11(a)). The air mass accelerates as it ows over the ridge and the Rhine Valley and experiences 15 km wavelength mountain waves (Fig. 12) penetrating down to the Rhine Valley (Fig. 13) and associated strong downward motion (Fig. 11). The differences between Meso-NH simulations and the observations along the trajectory of the CVB launched from Ispra, Italy on 30 October 1999 at 0811 UTC (Fig. 12) are less than 2 m s ¡1 for the horizontal wind speed, 10 ± for the wind direction, 1.2 K for the potential temperature and 0.25 g kg ¡1 for the mixing ratio. The CVB observations show an increase of the wind speed near the mountain ridge, a constant ow direction, no large variation of the temperature (305.2 K § 1.2 K) or humidity (3.0 g kg ¡1 § 0.3 g kg ¡1 ) (as observed during IOP 2 and 8, for instance) and vertical oscillations when the CVB reaches the mountain ridge. The isentropes are almost horizontal on the southern side of the Alps, with a drop of up to 3 km on the lee side. The induced downslope winds on the northern side (e.g. [START_REF] Smith | On severe downslope winds[END_REF] generate turbulence in the Rhine Valley as strong turbulent kinetic energy is simulated (not shown).

A rough estimate of the downslope surface wind can be computed by applying the steady-state equation of motion, i.e. when the pressure gradient term balances the advection term (e.g. Hoinka 1985;[START_REF] Hoinka | The surface layer on the leeside of the Alps during foehn[END_REF]. Using the observed pressure gradient of 8.5 hPa across the Alps (see section 4), we get a downslope wind velocity of 19 m s ¡1 which is in good agreement with the observations and simulations (e.g. Fig. 13). In the middle troposphere, wide spacing between the isentropes indicates destabilization and strong mixing with eventual wave overturning.

As the upstream upper-level ow is quasi-aligned with the longitudinal valleys, the channelling ef ciency of these valleys in directing the f öhn ow towards the FORM target area increases. Val Lumnezia plays a signi cantly more important role than during the shallow föhn, whereas the Praettigau does not channel the föhn ow (not shown). The intensi cation of the föhn ow is visible in Fig. 10(m) at 1000 m ASL where the ow splits at the bifurcation of the Rhine and Seez Valleys. At higher levels (Figs. 10(n) and (o)), the radial-velocity eld is very weak in SV, due to a strong downslope transverse wind from the Weisstannen Valley (note the blue area at the exit of the Weisstannen Valley meaning that the ow blows towards the TWL in Fig. 10(n)) and due to downward motion induced by the coupling between lower-and upper-level ows and mountain waves. The vertical velocity measured by the scintillometers provides evidence of a strong interaction between the low-level föhn jet and wave activity aloft in the free troposphere (Fig. 13). The magnitude of the downward motion can reach 3 m s ¡1 , so air masses advected from the south-west penetrate the Rhine and Seez Valleys and, thus, squeeze the ow layer that splits between the Rhine and Seez Valleys.

Table 1 shows that similarly to the previous day, the major geophysical characteristics of the three longitudinal tributaries on 30 October 1999, favour geostrophic adjustment and nonlinear processes, at least in the Reuss Valley and the Domleschg. The dimensionless depth g 1H D 0:82 of Val Lumnezia tends to a quasi-linear value, i.e. where linear theory provides meaningful predictions. Of course, a limitation of the scale analyses shown in Table 1 is the exclusion of: (1) strong local cross-valley ows on the kilometre scale due to valleys or gaps which drain the main valleys; (2) downslope acceleration due to thermal winds (e.g. breezes); or (3) slope-induced gravitational acceleration which can have a non-negligible impact.

DEVELOPMENT OF THE F ÖHN FLOW IN THE RHINE VALLEY AND ITS NEARBY

TRIBUTARIES

The IOP 12 föhn episode starts on 30 October 1999 at 0900 UTC and ends on 31 October 1999 at 0310 UTC. The detection of föhn and the determination of the start and end times of the episodes are given by the multi-parameter algorithm of [START_REF] Gutermann | Vergleichende untersuchungen zur föhnh äu gkeit im Rheintal zwischen Chur und Bodensee[END_REF] applied to the numerous surface measurements acquired in the Rhine Valley and its main tributaries. The algorithm is based on a thresholding of mean and gust valley winds, temperature increase >3 K and low humidities during the whole period. The föhn episodes were determined for each station individually. A föhn event is de ned by the earliest onset and by the latest end at any of the considered stations on the valley oor.

In this section, we analyse the different phases of the life cycle of this south f öhn event. We focus particularly on the three-dimensional structure of the föhn ow as a function of time and:

²

the time of occurence and cessation of ow splitting in the Rhine and Seez Valleys and its vertical extent, ² the spatial inhomogeneity of the föhn penetration within the Rhine Valley.

(a) Phase I: Pre-event situation The wind eld at surface level was very heterogeneous during the afternoon of 29 October 1999. The wind speed in the three valleys was at most 5 m s ¡1 . The weak heating by the sun in late fall and the cloud cover (all stations in the target area reported more than 4 octas of high and middle clouds) prevented the development of thermally driven wind systems.

Figure 6 displays ow reversal in the Heiligkreuz sounding at 1200 UTC, with higher wind speeds at 1500 UTC. This layer is about 1000 m deep. Above Buchs-Grabs, the layer with very weak northerly winds has about the same depth. Above 2500 m we observe a strong westerly component (see also Fig. 10(d)).

The shallow ow reversal at about 1000 m ASL can be observed with more detail in SV from the TWL data (Fig. 10(a)) and in RV2 from the scintillometer measurements of horizontal wind speed and direction (Fig. 13). Figure 13 shows a weak westerly/northwesterly valley wind between 1200 UTC and 1700 UTC. A very good agreement is found between the Meso-NH numerical simulations (dashed lines in Fig. 13) and the scintillometer measurements, during this phase as well as during the other phases of the föhn life cycle.

(b) Phase II: Sandwich shallow föhn Katabatic drainage ow in the Rhine Valley and its tributaries was established during the night from 29 October to 30 October 1999. The evolution of temperature, humidity and wind at surface stations in RV2 (see Fig. 4 for the Vaduz example) during the course of this night is very similar to that during one of the clear nights with very weak pressure gradients in late summer at the beginning of the SOP. In contrast to these clear nights, the contribution of radiative cooling to the katabatic ow on the night of 29/30 October 1999 is only small since large parts of the target area were covered by clouds. The only exception was the area between Malans and Heiligkreuz, where the cloud cover reduced to 4 octas during a few hours. The ow out of the Rhine Valley and its tributaries was supported by the isallobaric eld. Whereas pressure decreased north of the Alps due to the low-pressure system approaching from the west (Fig. 2), it remained more or less constant south of the Alps. The second effect of this differential pressure change was the establishement of a weak ow across the gaps in the Alpine ridge into the valleys north of the Alps. Wind speeds at the Alpine ridge varied between 5 and 7 m s ¡1 at Guetsch, and 1 and 5 m s ¡1 at Piz Martegnas, whereas in the Rhine Valley at 1000 m ASL the wind speed was between 9 and 10 m s ¡1 (Fig. 13). This ow supported the katabatic ow but it was too weak and too shallow to replace the air owing out of the Rhine and Seez Valleys and to become a real (shallow) f öhn ow. The air owing out of the Rhine and Seez Valleys was instead replaced by subsidence from above the large catchment of the Rhine Valley. The subsidence is con rmed by the temperature inversion which is higher in the southern part of the target area than in the northern parts (Fig. 6).

The katabatic ow reached Vaduz in RV2 and further than Quinten in SV. The wind speed at Vaduz was 3 to 4 m s ¡1 throughout this phase (Fig. 4). This is a typical value for thermally induced mountain wind (Egger 1990a, b). The humidity in the katabatic ow decreased from 70% on 29 October 1999 at 1900 UTC to 55% at the end of phase II at 0330 UTC on 30 October 1999.

North of Vaduz, the RV2 was covered by a shallow cold air pool. The thickness of this pool grew from less than 40 m at 1800 UTC on 29 October 1999 to 60 m at 2100 UTC and 0000 UTC, to 100 m at 0300 UTC, and 300 m at 0600 UTC on 30 October 1999. The temperature in this cold air was 5-6 K lower than in the katabatic ow (Figs. 14 and15, the Altenrhein station being in the cold air pool). Humidity in the cold air pool was near to saturation.

Due to the higher wind speed caused by channelling in the much narrower SV no distinct cold air pool developed in this valley (Figs. 14 and15) (the simulated potentialtemperature eld, not shown in the paper con rms this point). The channelling of the katabatic ow is clearly visible in the wind pro le (Fig. 6). The maximum wind speed gradually came down from 1500 m ASL to 800 m ASL which is a further indication of the subsidence. At the surface, wind speeds remained low throughout the night and the temperature at Heiligkreuz uctuated between 13 ± C and 15 ± C (Fig. 15).

Between the valley ow described above and the large-scale synoptic southwesterly ow a shallow layer with southerly ow established during this phase. This layer was more pronounced in the south and central part of RV2 than in any other part of the target area. Above Buchs-Grabs it stretched from 1000 to 1800 m ASL at 1800 UTC on 29 October 1999, was squeezed to 900-1500 m ASL with increasing wind speed (10 m s ¡1 at 1200 m) at 2100 UTC, 0000 UTC and 0300 UTC, and weakened again but stretched from 800 to 3100 m at 0600 UTC on 30 October 1999 (Fig. 6). This layer was continuously recorded by the scintillometers (Fig. 13). The changes in depth and speed coincide well with the higher wind speeds at the pass stations between 2100 UTC and 0300 UTC. Also, since the potential temperature in this ow is the same as at the passes (Guetsch at 2287 m ASL) we can rmly assume that this ow is that part of the ow which crosses the gaps in the Alpine divide. Because this shallow föhn ow does not reach the valley oor (as it should in order to ful ll the conventional de nitions) and it is caught between the katabatic ow below and the large-scale synoptic ow above, the term 'sandwich shallow föhn' is appropriate to describe the ow.

(c) Phase III: Cold air pool covers RV2 The regime described in the previous section lasted until 0400 UTC on 30 October 1999. During phase II, the mesoscale pressure minimum was situated in the northern part of RV2 or even north of it. Between noon on 29 October 1999 and midnight, pressure decreased by 1-2 hPa in the south of the target area, but by 4 hPa in its northern parts. Starting at 0300 UTC, the pressure eld established during phase II was disturbed. A small zone of higher pressure was established in the very northern parts of the Rhine Valley. Due to the steadily falling pressure, the pressure minimum in RV2 moved to the area between Vaduz and Sargans, in the southern parts of RV2. This caused a ow reversal in the upper parts of RV2, and allowed the cold air pool to pour into the region where we had katabatic ow during phase II (see Fig. 16 and Fig. 4 for the Vaduz example). The wind in this area turned to a northerly direction. The moisture content in the cold air pool was near to saturation (Fig. 15), but Diepoldsau was the only station which reported fog at 0700 UTC on 30 October 1999. Phase III lasted until 0700 UTC when the temperature started to increase rapidly. At most stations in the southern part of RV2 the same katabatic ow as in phase II established for a short time and later the temperature further increased. (d ) Phase IV: Penetrating föhn in RV1, RV2 and SV The main reason for the intensi cation of the föhn ow in RV and SV was the increasing cross-Alpine pressure difference which grew from 6 hPa at 0600 UTC to 8.5 hPa at 1200 UTC on 30 October 1999 (Fig. 3), and the rotation of the upperlevel ow from west to south-west (Fig. 10(p), fourth row). The change in the ow regime around 0600 UTC is clearly visible in Fig. 13 when the wind direction at about 1000 m ASL changed from 180 ± to 140 ± . This jump in wind direction, as well as the intensi cation of the ow, is visible in the Meso-NH simulation outputs. It becomes clear from the scintillometer wind speed, that the gusty and strong föhn ow started only at noon (Fig. 13). Thus, the main effect, until about 1200 UTC, of the föhn ow was the enhancement of the diurnal cycle, i.e. the temperature increase during the morning hours was steeper than the previous day. The rst temperature maximum was reached about 1 hour earlier and was 1-2 K higher than during the previous day (Fig. 15). The cold air pool had almost completely disappeared by 0900 UTC in RV1 and by 1200 UTC in SV and RV2 (Fig. 17(a)). Almost all stations between Masein in the southern tributary of RV1, and Ruggel in RV2, and Heiligkreuz in SV, show a clear positive jump in temperature and a negative jump in humidity between 1100 UTC and 1400 UTC on 30 October 1999 (Fig. 6). This is a clear indication that the föhn ow touched the surface. The distribution of potential temperature shows a well established föhn ow at 1500 UTC (Fig. 17(b)) which propagates at about 4.5 km h ¡1 . This led to a short-term small-scale pressure fall in the basin of Buchs-Grabs. The pressure minimum in RV2 was situated in this area for about 2 hours. As a consequence, the cold air pool situated in the northern parts of RV2 oscillated back to the south into the basin of Buchs-Grabs, and the temperature gradient between the föhn air and the cold air pool increased. At 1800 UTC, the border between the föhn ow and the cold air pool was again further to the north (Fig. 17(c)) and at 1900 UTC the föhn touched the surface at Diepoldsau (Fig. 17(d)). During the following hour, the föhn ow at the surface was con ned to the western valley wall due to the katabatic ow coming from the Walgau, east of Feldkirch. It was not until 2200 UTC that the föhn also touched the surface at Feldkirch. The cold air ow out of the Walgau was stronger and cut the föhn ow into a southern patch (from Balzers to Buchs-Grabs) and a northern part (Diepoldsau-Götzis) before the föhn was completely stopped by the approaching cold front at 0200 UTC.

CONCLUSION

This paper examines the three-dimensional structure and dynamics of the IOP 12 south föhn ow during its entire life cycle from 29 October 1999 until 31 October 1999 in the Rhine Valley, and speci cally in the FORM target area at the bifurcation between the Rhine and Seez Valleys, near Bad Ragaz, Switzerland. This south föhn case presents several levels of complexity:

²

In the time domain, this f öhn life cycle has strong similarities to [START_REF] Seibert | South foehn studies since the ALPEX experiment[END_REF] description, but also reveals the higher complexity of the föhn ow in the Rhine Valley compared with the comparatively simple structure in the Wipp Valley, on which Seibert's föhn life cycle is mainly based. As would be expected for a large valley with an extensive catchment, katabatic ow and cold air pools play an important role by preventing the föhn ow touching the ground. In addition, the topographical con guration in the Wipptal-Innsbruck area enable the quite simple distinction between katabatic ow and (shallow) föhn, because there is a change in the wind direction of roughly 90 degrees. In the Rhine Valley, katabatic ow and shallow föhn have very similar, at some places even identical, wind direction, which makes the distinction very dif cult [START_REF] Drobinski | Statistical characterization of the ow structure in the Rhine valley[END_REF].

²

In the space domain, strong interactions are shown between the synoptic scales that force cross-Alpine temperature, humidity and pressure gradients, the mesoscales that drive the föhn circulation in the complex valley network, and the local scales that cause highly inhomogeneous thermodynamical patterns.

In the FORM target area, the thermodynamical pattern is extremely inhomogeneous (both horizontally and vertically) and varies with time. This paper investigates in detail the multiscale processes driving interactively the structure of the cold air pool and the warm föhn air. These are:

²

Local forcing: radiative cooling, local pressure and temperature gradients (cold pool, thermal winds; : : : ).

²

Mesoscale forcing by the valleys: complex valley network tributaries directing the ow towards the FORM target area, ow splitting at bifurcations between valleys, föhn ow splitting near Feldkirch by cold air owing out the Walgau; : : : ² Synoptic-scale forcing, such as large-scale pressure gradients driving the direction and intensity of the katabatic and föhn ows in the valley, and coupling between the valley ow and tropospheric mountain waves.

Future work will focus on the processes of interaction between the tropospheric waves and the föhn ow, or more generally the planetary boundary layer. Also, in the framework of MAP, the Wipp Valley was highly instrumented, and thorough comparisons can be made between the föhn structure, intensity and evolution in the Rhine and Wipp Valleys.

Figure 1 .

 1 Figure 1. Topography of the FORM (see text) target area (a), zoomed on the region of interest on a local scale (between 9.3 and 9.6 ± E and 46.9 and 47.2 ± N) (b), which is also the most instrumented area (see section 2). The isocontour increment is 500 m from 500 m altitude to 3000 m altitude. Italic writing indicates names of valleys, while legends in the boxes are the name of cities. SV, RV1 and RV2 correspond to the Seez Valley, the lower Rhine Valley and the upper Rhine Valley, respectively. The dashed straight lines indicate the scintillometer light beams. The scintillometer transmitters are located in Triesenberg (on the eastern side of RV2), the receivers in Flusa and Ergellen (on the western side of RV2).

Figure 2 .Figure 3 .

 23 Figure 2. Synoptic situation for Intensive Observing Period 12 in 6-hourly intervals from ARPEGE analyses. Solid lines and shaded area: geopotential height at 500 hPa, with an increment of 50 m; dashed lines: mean-sealevel pressure with an increment of 5 hPa. The rectangle indicates the location of the FORM target area (see text).

Figure 4 .Figure 5 .

 45 Figure 4. Time series of (a) surface wind speed, (b) wind direction, (c) temperature (solid line) and temperature growth (dashed line), and (d) relative humidity measured at Vaduz. The shaded areas show the range for wind speed (¸3 m s ¡1) and direction (¸125 ± and •215 ± ), temperature growth (¸3 K h ¡1 for the onset of föhn and ¸¡6 K hr ¡1 for the end of föhn), and relative humidity (•45%) when föhn occurs at Vaduz. These criteria must be met simultaneously.

TimeFigure 6 .

 6 Figure 6. Same as Fig. 5 at (a) Heiligkreuz, (b) Buchs-Grabs, and (c) Malans.

Figure 7 .Figure 9 .

 79 Figure 7. 2 km resolution horizontal wind eld at (a) 1000 m above sea level (ASL), (b) 1500 m ASL, and (c) 2000 m ASL, as obtained from Meso-NH numerical simulations on 29 October 1999 at 1800 UTC.

Figure 10 .

 10 Figure10. Transportable Wind Lidar (TWL) radial-velocity eld at 1000 m above sea level (ASL, rst column), 1500 m ASL (second column), 2000 m ASL (third column) and 2500 m ASL (fourth column) retrieved from raster scans on 29 October 1999 at about 1200 UTC ( rst row), 1700 UTC (second row), 2000 UTC (third row) and on 30 October 1999 at about 1200 UTC (fourth row). Positive (negative) radial velocities correspond to air blowing away from (towards) the TWL. Absence of colour means no reliable data (low signal to noise ratio, presence of topography). The data quality-control procedure is explained in[START_REF] Dabas | Adaptive lters for frequency estimates of heterodyne Doppler lidar returns: Recursive implementation and quality control[END_REF]. The marker £ indicates the location of the TWL. SV, RV1 and RV2 denote the Seez Valley, the lower Rhine Valley towards Chur, and the upper Rhine Valley towards Lake Constance, respectively.

Figure 11 .Figure 12 .Figure 13

 111213 Figure 11. Vertical cross-section of Meso-NH (a) vertical velocity (grey scale, m s ¡1 ) and potential temperature (isentropes in solid line, K), and (b) south-north wind component (grey scale, m s ¡1 ) on 30 October 1999 at 1200 UTC between (9.5 ± E, 46.0 ± N) and (9.5 ± E, 48.0 ± N). This cross-section transects the centre of RV2 (see text).

Figure 14 .Figure 15 .

 1415 Figure 14. Surface potential-temperature analyses in the FORM (see text) target area during phase II using the VERA algorithm (see text). The analyses are shown (a) for 29 October 1999 at 1800 UTC and (b) for 30 October 1999 at 0300 UTC (the thermodynamical pattern does not evolve in between). The ordinates correspond to the longitude and the abscissae to the latitude.

Figure 16 .

 16 Figure 16. Same as Fig. 14 for phase III on 30 October 1999 at 0600 UTC.

Figure 17 .

 17 Figure 17. Same as Fig. 14 for phase IV on 30 October 1999 at (a) 1200 UTC, (b) 1500 UTC, (c) 1800 UTC and (d) 2100 UTC.

TABLE 1 .

 1 SCALE ANALYSIS OF THE FLOW IN THE REUSS VALLEY, VAL LUMNEZIA AND DOMLESCHG ¡4 s ¡1 , dµ=dz D 0:007 K m ¡1 , N D 1:54 £ 10 ¡2 s ¡1 and U D 17 m s ¡1 Along-valley wind speed v (m s ¡1 )

		Reuss Valley Val Lumnezia Domleschg
	Geometrical settings
	Height of the crest H (m)	2500	2100	2100
	Depth of the valley 1H (m)	1100	900	1200
	Length of the valley L (km)	45	16	32
	Width of the valley l (km)	5	7	6
	29 October 1999 at 1800 UTC
	f D 10 ¡4 s ¡1 , dµ=dz D 0:013 K m ¡1 ,
	N D 2:10 £ 10 ¡2 s ¡1 and U D 13 m s ¡1
	Along-valley wind speed v (m s ¡1 )	3	1	3
	Along-valley Rossby number R l	6.00	1.42	5.00
	Cross-valley Rossby number RL	0.07	0.21	0.17
	Dimensionless valley depth g 1H	1.78	1.45	1.94
	Rossby radius of deformation lR (km)	273	231	315
	30 October 1999 at 1200 UTC
	f D 10 12	4	8
	Along-valley Rossby number R l	24.00	5.71	13.33
	Cross-valley Rossby number RL	0.29	0.86	0.47
	Dimensionless valley depth g 1H	1.00	0.82	1.09
	Rossby radius of deformation lR (km)	200	169	231

ACKNOWLEDGEMENTS

The authors would like to thank the scientists of the MAP FORM working group, coordinated by H. Richner, and particularly G. Jaubert and G. Beffrey for fruitful discussions; B. Bénech and H. Berger for providing the constant-volume balloon data and help in analysing them; M. Aupierre, C. Boitel, A. Delaval, P. Delville, J. M. Donier, C. Loth and B. Romand for assistance in operating the TWL during the MAP eld experiment; M. C. Lanceau for help in collecting the referenced papers; E. Haeller for the MAP data centre management. The effort in collecting radiosonde measurements in a dense network and at high temporal resolution by the Swiss Army joint Meteorological Services is gratefully acknowledged. These stations were signi cantly supported with data loggers and computer hardware by Meteolabor AG and Polyinformatik GmbH. This research was funded by the Centre National de Recherche Scienti que (CNRS) through the Programme ATmosphère Océan à Méso-échelle (PATOM), the Institut National des Sciences de l'Univers (INSU) and by the Centre National d' Études Spatiales (CNES).


Afif Chandoul
email: afif.chandoul@g-scop.inpg.fr

Van-Dat Cung

Fabien Mangione
email: fabien.mangione@g-scop.inpg.fr

Chandoul Afif

Reusable containers within reverse logistic context

Keywords: Reusable transport packages, network flow, optimization, container management, supply chain

published or not. The documents may come

Introduction

Since the mid-1980s, reusable containers have been adopted by many companies. In 1987, for example, the Bergen Brunswig Drug Company in California purchased 120,000 returnable plastic containers to replace one-use corrugated cartons. The company ships from its 37 distribution centers to its 10,000 pharmacies, located in 40 states (D. Saphire, 1994). Many other manufacturers in electronic goods and the automobile industry have also switched to reusable containers. This adoption is essentially due to potential economic benefits. It enables companies to save money by decreasing packaging material requirements: generally the cost of reusable containers is amortized over their lifecycle and the more the container is used the more the cost per trip is decreased. It reduces product damage due to shipping and handling because reusable containers are generally sturdier than one-way container and they are designed to withstand multiple uses thus providing better product protection. Other savings are possible. Warehouse utilization can be improved by reducing storage space requirements since reusable containers can be stacked higher than one-way containers. Improvements in worker safety can reduce costs since their ergonomic design reduces injuries from box cutters, staples, debris and stray packaging(D. Saphire, 1994). Furthermore, reusable containers are better for the environment since they reduce packaging waste. This can be a tool to meet the waste reduction requirements of government regulations, which are especially strict in EU countries (M. Kärkkäinen et al., 2004). However, reusing containers requires a more complex supply chain. The purchasing costs of reusable containers are significantly higher than those of one-way containers. A reusable container, such as a plastic box, can cost 10 times more than a one-way container, such as a corrugated box (D. Saphire, 1994). Therefore, for reusable containers to be of benefit, efficient container-management is a top priority. Despite the importance of specific container management, there are few academic studies on the subject. Those that do exist can be classified under 4 headings:

-Organisational: these studies deal with the different issues related to the implementation of reusable container systems, like collaboration and cooperation between supply chain members. They, generally, stress the role of information systems [START_REF] Chan | A systematic approach to manufacturing packaging logistics[END_REF][START_REF] Chan | A pro-active and collaborative approach to reverse logistics—a case study[END_REF][START_REF] Eelco De | Making Waves : RFID Adoption in Returnable Packaging[END_REF][START_REF] Kärkkäinen | Managing the rotation of reusable transport packaging a multiple case study[END_REF]) . -Environmental: These studies present the environmental advantages of reusable containers compared to one-way containers [START_REF] Gassel | Returnable packaging for non-specific building materials[END_REF][START_REF] González-Torre | Environmental and reverse logistics policies in European bottling and packaging firms[END_REF]S. Paul Singh et al., 2006b). -Economic: in these studies, the economic benefits of reusable containers are shown compared to one-way containers. Often, these studies present the environmental advantages linked to economic benefits(Lerpong [START_REF] Jarupan | Evaluation of Trade-offs in Costs and Environmental Impacts for Returnable Packaging Implementation[END_REF][START_REF] Martin | Costing Structures of Reusable Packaging Systems[END_REF]Orbis, 2004;[START_REF] Stopwaste | A cost comparison model for Reusable Transport Packaging[END_REF]). -Operational: these studies propose mathematical tools that can help companies optimize the different aspects of their reusable container systems. We can find studies that propose models which optimize the total cost generated by the use of reusable containers (transport storage, maintenance…). Others deal with the problem of container depot localisation (ED Castillo and Cochran, 1996;Erera Alan L et al., 2005a;[START_REF] Erera Alan | Robust Optimization for Empty Repositioning Problems[END_REF][START_REF] Karimi | Scheduling tank container movements for chemical logistics[END_REF][START_REF] Kroon | Returnable containers: an example of reverse logistics[END_REF]. In this paper we will deal with the operational side of the problem. We will consider a system that contains three kinds of site: customer sites, supplier sites and depots: A customer site is a facility that receives the loaded containers from a supplier site. After the containers have been emptied the customer can return them to depots or supplier sites for storage and reuse. A supplier site is a facility that receives empty containers from depots or customer sites in order to load and send them to customers. A depot is a place in which we store empty containers after they have been used by customers. Depots can supply suppliers with empty containers. We suppose that storage is also possible at customer and supplier sites. Our purpose in this work is to propose a model that optimizes the total cost resulting from transportation between sites and storage and to propose an original resolution method for some particular cases. This method is based on simple recurring formulas that can be easily solved for the different decisions variables. This can help us to better understand the system behaviour and, in turn, can help us deal with other related problems, such as system dimensioning and purchasing policies for reusable containers. The remainder of this paper is organized as follows. In section 2, we will present a literature review of the container management problem from the operational side. In this review we will concentrate on the work presented by [START_REF] Kroon | Returnable containers: an example of reverse logistics[END_REF]. In this work the authors present a step-by-step analysis of the operating policies of container management systems. This work will be a reference for our paper. In section 3, a description of our problem will be presented in detail and an optimization model, in the form of a network flow problem, will be developed. This model is a basic formulation of the case with several customer sites, several supplier sites and several depots under some simplifying assumptions. In section 4, we will simplify the model and treat the case of one customer site, one supplier site and one depot. In this section an original resolution method will be developed. In section 5, a result analysis and interpretation will be presented. In the final section, we will outline our future intentions for extension and enhancement of the presented model.

Literature Review

Most of the existing works have focused on sea container management. This is due to their large cost. One of the problems which has been studied is container repositioning. [START_REF] Karimi | Scheduling tank container movements for chemical logistics[END_REF] presented a new linear programming methodology based on a continuous-time approach which they called the event-based "pull" approach. In this continuous-time representation, the principle is to fix all possible event times a priori, taking into account some simplifications and assumptions. Initially, a chronologically ordered superlist of possible instances is generated identifying which container movements (events) may occur and the types of movements that may occur at each such instance. This superlist of times and events is then used to develop a linear programming (LP) formulation whose solution will define the events that minimize the total logistic (flow) cost. [START_REF] Erera Alan | Global intermodal tank container management for the chemical industry[END_REF] proposed two different formulations, one deterministic (Erera Alan L et al., 2005a), one stochastic [START_REF] Erera Alan | Robust Optimization for Empty Repositioning Problems[END_REF], for the same problem. In the deterministic formulation, they integrated container booking and routing decisions with repositioning decisions in the same model. They believe that a container operator who uses this formulation via a global system may indeed be able to reduce costs and improve equipment utilization as compared to the standard approaches which ignore container routings. In the Stochastic formulation they presented a stochastic model in which uncertain parameters are assumed to fall within an interval around a nominal value and they established some conditions for robust flow under certain definitions. Similar work has been developed in the soft drinks field for the returnable bottle. Del Castillo and Cochran (1996) modelled the reusable bottle production and distribution activities of a Soft Drink Company. They optimized the problem in terms of increasing the number of empty bottles available in order to make the comany more competitive. This is an alternative to the total cost optimization when cost calculation is not possible. It is important to note that the majority of texts deal with the problem in specific cases (sea containers or returnable bottles). The only text that we found that considers the problem in the general case is the paper of [START_REF] Kroon | Returnable containers: an example of reverse logistics[END_REF]. In this paper, the authors classified returnable container systems into three groups depending on their operating policies: switch pool systems, systems with return logistics, and systems without return logistic. In a switch pool system, participants have their own allotments of containers. Thus the cleaning, maintenance and storage are the responsibility of each pool-participant. Poolparticipants may be the senders and recipients, or the senders, carriers, and recipients of the goods. If only the sender and the recipients are allotted containers, a transfer will take place when the goods are delivered to the recipient. In this case, the role of the carrier is to transport the loaded container to the recipients and to return the empty containers to the sender. Alternatively, the carrier is also allotted containers and replaces the loaded containers with empty containers when it picks up a load from a sender or recipient. In a system with return logistics, the containers are owned by a central agency. It is responsible for returning the containers after they are unloaded. There are two variants of this system.

• Transfer system: the core of this system is that the sender always uses the same containers. He is responsible for tracking, tracing, maintaining, cleaning and storing them.

• Depot system: Under the depot system, the empty containers that are not in use are stored at container depots, which can supply the senders with empty containers on demand.

In the last system, the system without return logistics, the containers are owned by a central agency, and the sender rents them for fixed periods. The containers return to the agency after use. The sender is responsible for return logistics, cleaning, maintenance, storage and control. In the same paper, [START_REF] Kroon | Returnable containers: an example of reverse logistics[END_REF] proposed a mathematical model for the design of a return logistics system for returnable containers with a depot system variant. The main purpose of the model is to determine the suitable number of containers, the appropriate number of container depots and their locations in order to minimize the total cost. The model considered in our paper assumes the operator is using a system with return logistics with the depot system variant. We suppose that the central agency can be both the customer (recipient), suppliers (sender) or third parties.

Problem description

We assume that the transportation network consists of depots, customer sites and supplier sites (see figure 1). For a period t, in order to satisfy customer demand, suppliers must have a sufficient number of empty containers at their disposal. These empty containers can be obtained from different customer sites, different depots or their own stocks. Once the containers are loaded they will be sent to their destinations (customers) during the same period. At customer sites the loaded containers will be emptied and can be returned to either depots or supplier sites as the best place for storage and reuse in subsequent periods. The different movements involved in this system generate large costs. Our purpose here is to optimize these costs and to propose the best plan for transportation and storage in each period. We will model the problem as a network flow (see figure 2) but first we will make some assumptions for simplification purposes and present the different decision variables.

Decisions variables

Let D be the number of depots (d=1, …, D), I be the number of customers sites (i=1,…, I), F be the number of Suppliers sites (f=1,…, F) and T be the number of periods in the planning horizon (t=1,…, T). We define the different decisions variables as follows:

•

Constraints

The constraints of the system are essentially the flow conservation constraints. They can be classified into four categories:

• the first represents the conservation flow before sending the loaded containers to different sites and depots

t i d t id f t if t i a u y a 1 + + = ∑ ∑ ] , 0 [ , ] , 0 [ I i T t ∈ ∀ ∈ ∀ ( 1 ) ∑ ∑ ∑ - - + = d t df i t if t f i t fi t f z y b x b 1 ] , 0 [ , ] , 0 [ F f T t ∈ ∀ ∈ ∀ ( 2 ) ∑ ∑ - + = i t id t d f t df t d u d z d 1 ] , 0 [ , ] , 0 [ D d T t ∈ ∀ ∈ ∀ ( 3 )
For example, constraint ( 1 ) represents the conservation flow at customer site i at period t; the principal inventory t i a is equal to the sum of departures from site i to depots ∑ • The second category of constraints represents the conservation flow after the emptying of loaded containers and their return to different sites.

∑ ∑

- - + = + d t id f t if t i t i t i u y a x a 1 1 1 1 ] , 0 [ , ] , 0 [ I i T t ∈ ∀ ∈ ∀ ( 4 ) t f d i t if t df t f b y z b 1 1 1 1 1 1 + + = ∑ ∑ - - + ] , 0 [ , ] , 0 [ F f T t ∈ ∀ ∈ ∀ ( 5 ) ∑ ∑ - + = + f t df t d i t id t d z d u d 1 1 1 1 ] , 0 [ , ] , 0 [ D d T t ∈ ∀ ∈ ∀ ( 6 )
• The third category represents the sum of containers arriving at a customer site i at period t.

∑ = i f t i t fi x x , ] , 0 [ , ] , 0 [ I i T t ∈ ∀ ∈ ∀ ( 7 )
• The last category of constraint is the lower band constraint, that is the sum of arriving containers must be greater than or equal to the demand at site i.

t i t i x ≤ δ ] , 0 [ , ] , 0 [ I i T t ∈ ∀ ∈ ∀ ( 8 )
We suppose that there is no capacity constraint on storage and transportation. The objective function can be expressed as: In the following sections, we will simplify the model; we will deal with the case with only one customer site, one supplier site and one depot. We will show that it is possible to determine all the decision variables with a simple recurring formula.

t f f t st bf t d d t st dd t i i t st ai d i t t id t id Tr d f t t df t df Tr t if t if i f t Tr t fi i f t Tr b C d C a C u u C z z C y y C x C C uid zdf yif xfi ∑ ∑ ∑ ∑ ∑ ∑ ∑ + + + + + + + + + = , , , , , , , , , , , ) 1 (

Studied case

In this case we add the following assumptions to the system (see figure 3):

• We have one customer site, one supplier site and one depot • We assume that the depot has the lowest cost, followed by the customer site and the supplier site being the most expensive. The different constraint can be simplified as follows

t t t t a u y a 1 + + = t t t t t z y b x b - - + = 1 t t t t u d z d - - = 1 t t t t t u y a x a 1 1 1 1 - - + = + t t t t b y z b 1 1 1 1 + + = + t t t t z d u d 1 1 1 1 - + = + t t x δ ≥
The objective function is given by

b C d C a C u u C z z C y y C x C C u z y x ∑ ∑ ∑ ∑ ∑ ∑ ∑ + + + + + + + + + = ) 1 ( ) 1 ( ) 1 (
Now we model the problem graphically as a network R=(X, U, a, b, c) for use in the different proofs. X: Set of nodes that represents the three sites at each period t U: Set of arcs joining nodes. These arcs can be storage arcs or transportations arcs. Each arc represent a decision variable and takes the same notation but in the capital form. For example, the arc which correspond to t y take the notation of t Y . a (u) : represent the cost of an arc U u ∈ . For example,

, 0 [ ) ( ) ( ) /( 0 ) ( T t X u u b X U u u b t t t ∈ ⎪ ⎩ ⎪ ⎨ ⎧ = ∀ = ∈ ∀ = δ c(u)
: represents the capacity of an arc. We suppose that it is infinite on all arcs

∞ = ∈ ∀ ) (u c U u
The different pairs of variables and arcs are shown with the corresponding arrow on figure 4 below.

In the remainder of this work, the majority of the properties will be demonstrated using the minimum cost network flow theorem. This theorem states that a necessary and sufficient condition for a feasible solution f to be optimal on a network R=(X, U, a, b, c) is that for each cycle

- + Γ ∪ Γ = Γ such that { } 0 )] ( ) ( [ ; )] ( ) ( [ > - - = - + Γ ∈ Γ ∈ u b u f Min u f u c Min Min u u γ it is necessary that 0 ) ( ) ( ≥ -∑ ∑ - + Γ ∈ Γ ∈ u u u a u a ) 1 , 1 ( 1 1 - - t t u U ) 1 , 1 ( 1 1 - - t t z Z ) 1 , 1 ( 1 1 - - t t y Y ) , ( t t u U ) , ( t t z Z ) , ( t t y Y ) , ( t t a A ) , ( t t b B ) , ( t t c C ) 1 , ( t t b B ) 1 , 1 ( 1 t c C ) , ( t t x X ) , ( 1 1 + + t t x X ) 1 , 1 ( t t a A ) 1 , ( t t u U ) 1 , ( t t z Z ) 1 , ( t t y Y ) , ( t t u U ) , ( t t z Z ) , ( t t y Y ) , ( 1 1 + + t t a A ) , ( 1 1 + + t t b B ) , 1 ( 1 1 + + t t c C ) 1 , 1 ( 1 1 + + t t b B ) 1 , 1 ( 1 1 + + t t c C ) 1 , 1 ( 1 1 + + t t a A ) 1 , 1 ( 1 1 + + t t u U ) 1 , 1 ( 1 1 + + t t z Z ) 1 1 ( 1 1 + + t t y Y

Figure 4 : Network flow illustration

Before resolving this problem we will define some new integer parameters, which will help us later when demonstrating properties. The first parameter is the principal threshold, a d N . This parameter will help us to decide whether or not to send surplus empty containers from the customer to the depot for storage and reuse in the following a d N periods. The second parameter is the secondary threshold, a d M . This parameter will show whether sending the container from the customer to the depot for storage without reuse is better then keeping it at the customer site in the following a d M periods.

Definition 1. in the case when we have

Tr Tr Tr z u y C C C + <
, we define the principal-storingthreshold of the customer to the depot, the integer a d N , such that:

⎪ ⎩ ⎪ ⎨ ⎧ + × + + > + × + + × + < + ×
C C N C C C N C C N C C C N ) 1 ( ) 1 ( Remark 1.
- + - + + = Γ a d a d a d N t N t t t N t Y A A A , in
other words this cost is the cost of storing one unit of container during a d N periods and transporting it from the customer to the supplier. The second term,

Tr st c a d Tr z u C C N C + × + , corresponds to the cost of path ) , ,..., 1 , , 1 ( 1 1 1 , - + - + - - = Γ a d a d a d N t N t t t t N t Z D D D U
which represents the total cost of sending one container from the customer to the depot, storing it during a d N in the depot and sending it on to the supplier. The intention here is to note that storing containers more than a d N at the customer site and transporting them subsequently to suppliers can be more expensive than transporting them to the depot, storing during more than a d N , and sending them on. It is easy to remark that

+ - Γ ∪ Γ a d a d N t N t , , constitutes a cycle 1 1 - t u t u t y t z t a t a1 t b1 t b t d t d1 1 1 - t y t u1 t y1 1 - a d N + Γ a d N t , - Γ a d N t , 1 - + a d N t y 1 - + a d N t z 1 - + a d N t a 1 - + a d N t d Figure 5 : Cycle + - Γ ∪ Γ a d a d N t N t , ,
Remark 2. It is easy to verify that for all integers i, such that

a d N i ≤ , Tr st d Tr Tr st a z u y C C i C C C i + × + < + ×
and that for all

a d N i > , Tr st d Tr Tr st a z u y C C i C C C i + × + > + × .
This can be proved as follows:

Let i be an integer such that

a d N i ≤ . i can be written in the form j N i a d - =
where j is an integer. Since

st d st a st d st a C j C j C j C j × - < × - ⇒ × > × summing this with Tr st d a d Tr Tr st a a d z u y C C N C C C N + × + < + × gives Tr st d a d Tr Tr st a a d z u y C C j N C C C j N + × - + < + × - ) ( ) ( which is exactly Tr st d Tr Tr st a z u y C C i C C C i + × + < + × . Now suppose that 1 + > a d N i
. i can be written in the form

j N i a d + + = 1 where 0 > j . Since st d st a C j C j × > ×
summing this with

Tr st c a d Tr Tr st a a d z u y C C N C C C N + × + + > + × + ) 1 ( ) 1 (
gives

Tr st d a d Tr Tr st a a d z u y C C j N C C C j N + × + + + > + × + + ) 1 ( ) 1 ( which is exactly Tr st d Tr Tr st a z u y C C i C C C i + × + > + × .
Remark 3. The first inequality (

Tr st c a d Tr Tr st a a d z u y C C N C C C N + × + < + × ) is valid only if Tr Tr Tr z u y C C C + <
. In fact, the first inequality can be written as: . In this case it is simple to see that if we have stock at the customer site at the end of period t it will be better to send all of it to the depot for storage and reuse.

Tr Tr Tr st d st a a d y z u C C C C C N - + < - × ) ( . Given that 0 ) ( > - × st d st a a d C C N
Definition 2. We define the secondary-storing-threshold of the customer to the depot, the integer a d M , such that:

⎪ ⎩ ⎪ ⎨ ⎧ × + + > × + × + < ×
] , 0 [ T t ∈ , 0 1 1 = = t t z y Proof. Suppose there exists a ] , 0 [ T t ∈ such that 0 1 ≠ t z
, and hence

0 1 ≠ + t b . Let - + Γ Γ = Γ t t t U
be a cycle composed of 4 arcs (see Figure 7) such that:

- Γ t =( ) , 1 1 + t t B Z + Γ t =( ) , 1 1 + + t t Z D t y1 1 + t a 1 + t b 1 + t d t z1 t u1 1 + t y 1 + t z 1 + t u Figure 7 : Cycle - + Γ Γ t t U We have )] ( ) ( [ u b u f Min t u - - Γ ∈ because 0 1 > t z and 0 1 > + t b . In addition, 0 )] ( ) ( [ > - + Γ ∈ u f u c Min t u
because the capacities of arcs are unlimited. Therefore,

{ } 0 )] ( ) ( [ ; )] ( ) ( [ > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min t t u u
γ and according to the minimum cost network flow theorem we must have that 0

) ( ) ( ≥ -∑ ∑ - + Γ ∈ Γ ∈ t t u u u a u a
. However, in this

case ∑ + Γ ∈ + = t u tr z st d C C u a ) ( < st b tr z u C C u a t + = ∑ - Γ ∈ ) ( which is absurd, therefore t z1 must equal zero, i.e. for all ] , 0 [ T t ∈ ⇒ 0 1 = t z
.

The same reasoning can be applied to t y1 .

Suppose there exists a ] , 0

[ T t ∈ such that 0 1 ≠ t y
, and hence

0 1 ≠ + t b . Let - + Γ Γ = Γ t t t U
be a cycle composed of 4 arcs (see Figure 8) such that:

- Γ t =( ) , 1 1 + t t B Y + Γ t =( ) , 1 1 + + t t Y A t y1 1 + t a 1 + t b 1 + t d 1 1 - t z t u1 1 + t y 1 + t z 1 + t u Figure 8 : Cycle - + Γ Γ t t U We have )] ( ) ( [ u b u f Min t u - - Γ ∈ because 0 1 > t y et 0 1 > + t b
. In addition,

0 )] ( ) ( [ > - + Γ ∈ u f u c Min t u
because the capacities of arcs are unlimited. Therefore,

{ } 0 )] ( ) ( [ ; )] ( ) ( [ > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min t t u u
γ and according to the minimum cost network flow theorem we must have that 0

) ( ) ( ≥ -∑ ∑ - + Γ ∈ Γ ∈ t t u u u a u a
. However, in this case

tr y st a u C C u a t + = ∑ + Γ ∈ ) ( < st b tr y u C C u a t + = ∑ - Γ ∈ ) ( which is absurd, therefore 0 1 = t y for all ] , 0 [ T t ∈
. Remark 5. Property 1 implies that the stock at the supplier site will remain equal to zero if we begin with a stock equal to zero at t=0 .

Property 2. For ] , 0 [ T t ∈ ,if 0 1 ≠ t u then ] , 1 [ a d N k ∈ ∀ , 0 = +k t z Proof. Suppose there exists a ] , 0 [ T t ∈ such that 0 1 ≠ t u . Let k t , Γ be a cycle with ] , 0 [ T t ∈ and ] , 1 [ a d N k ∈
(see figure 9). Suppose that

- + Γ ∪ Γ = Γ k t k t k t , , , , where )
, ,..., ( 1 , k t k t t k t Y A A + + + + = Γ
, is the path composed of the storage arcs between t+1 and t+k and the transportation arc between the customer and the supplier. ) , ,.., , 1 (

1 , k t k t t t k t Z D D U + + + -= Γ
is the path composed of the transportation arc from the customer to the depot at period t, the storage arcs at depot between t+1 and t+k and the transportation arc from the depot to the supplier at period t+k .

+ Γ k t , - Γ k t , t u1 1 + t u 1 + t y 1 + t z 1 + t a 1 1 + t a 1 1 + t b 1 + t b 1 + t d 1 1 + t d k t u + -1 1 k t u + k t y + k t z + k t a + k t a + 1 k t b + 1 k t b + k t d + k t d + 1 k t y + -1 1 k t z + -1 1 t y1 t z1 Figure 9 : Cycle - + Γ ∪ Γ k t k t , ,
Suppose there exists a 0 ≠ +k t z which is the first non-zero occurrence, that is if

1 > k then for all ] 1 , 1 [ - ∈ k i we have 0 = +i t z
. For, the case when k=1,

+ + + -= Γ
will have non-zero values, and

thus 0 )] ( ) ( [ , > - - Γ ∈ u b u f Min k t u and since 0 )] ( ) ( [ , > - + Γ ∈ u f u c Min k t u it follows that: { } 0 )] ( ) ( [ ; )] ( ) ( [ , , > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min k t k t u u γ . We Know that that ∑ + Γ ∈ = k t u u a , ) ( Tr st a y C C k + ×
and

Tr st d Tr z u k t u C C k C u a + × + = ∑ - Γ ∈ , ) ( this implies that ∑ ∑ + Γ ∈ - Γ ∈ < - k t u k t u u a u a , , 0 ) ( ) (
which is absurd because according to definition 1 for all

1 , 0 [ - ∈ t i 0 = i b . According to property 1, 0 1 1 1 1 = = - - t t z y
and therefore, 0

1 1 ≠ - t b . Given that 0 1 ≠ - t δ and 0 1 = - t b , it follows that ) 0 , 0 ( ) , ( 1 1 ≠ - - t t z y 1 st case : if 0 1 ≠ - t y Let - + Γ Γ = Γ t y t y t y , , ,

U

(see figure 10) be a cycle such that:

) , , 1 ( ) , 1 , ( 1 , 1 1 , t t t t y t t t t y Y A A B B Y - + - - - = Γ = Γ 1 1 - t u t y t z t a t b t d 1 1 - t y 1 - t u 1 - t y 1 - t z 1 - t a 1 1 - t a 1 1 - t b 1 - t b 1 - t d 1 1 - t d Figure 10 : Cycle - + Γ Γ t y t y , , U The path ) , 1 , ( 1 1 , t t t t y B B Y - - -= Γ
does not contain an arc of zero value, given that 0 1 ≠ t y and

0 1 1 ≠ - t b . Consequently, 0 )] ( ) ( [ , > - - Γ u b u f Min t y and given that 0 )] ( ) ( [ , > - + Γ u f u c Min t y it follows that { } 0 )] ( ) ( [ ; )] ( ) ( [ , , > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min t y t y u u γ . Since ∑ ∑ - + Γ ∈ Γ ∈ + = < + = t y t y u st b tr y u tr y st a C C u a C C u a , , ) ( ) ( (given that st a st b C C > ) ⇒ 0 ) ( ) ( , , < -∑ ∑ - + Γ ∈ Γ ∈ t y t y u u u a u a
this is absurd, therefore t b can never non-zero. Consequently,

given that 0 1 1 1 1 = = - - t t z y it also follows that 0 1 1 = - t b
.

2 nd case: if 0 1 ≠ - t z : Let - + Γ Γ = Γ t z t z t z , , ,

U

(see figure 11) such that:

) , , 1 ( ) , 1 , ( 1 , 1 1 , t t t t z t t t t z Z D D B B Z - + - - - = Γ = Γ Given that 0 1 ≠ - t z , 0 1 1 ≠ - t b and 0 ≠ t b , it follows that 0 )] ( ) ( [ , > - - Γ u b u f Min t z on ) , 1 , ( 1 1 , t t t t z B B Z - - -= Γ . In addition, since 0 )] ( ) ( [ , > - + Γ u f u c Min t z it follows that { } 0 )] ( ) ( [ ; )] ( ) ( [ , , > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min t z t z u u γ . In this case ∑ ∑ - + Γ ∈ Γ ∈ + = < + = t y t z u st b tr z u tr z st d C C u a C C u a , , ) ( ) ( (given that st d st b C C > ) ⇒ 0 ) ( ) ( , , < -∑ ∑ - + Γ ∈ Γ ∈ t z t z u u u a u a
which is also absurd, therefore t b can never be null. Therefore,

given that 0 1 1 1 1 = = - - t t z y it follows that 0 1 1 = = - t t b b ( 1 1 1 1 1 1 - - - + + = t t t t b y z b
).

Finally

U

be a cycle (see figure 12) such that:

) , 1 ( ) , ( 1 , , t t t u t t t u D U U A - + - = Γ = Γ t u1 1 + t y 1 + t z 1 + t a 1 + t b 1 + t d t y1
- + Γ Γ t u t u , , U Since 0 ≠ t u
, it also follows that 0 ≠ t a

, and as a consequence that

0 )] ( ) ( [ , > - - Γ u b u f Min t u . In addition, since 0 )] ( ) ( [ , > - + Γ u f u c Min t u it follows that { } 0 )] ( ) ( [ ; )] ( ) ( [ , , > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min t u t u u u γ . However, since ∑ ∑ - + Γ ∈ Γ ∈ + = < + = t u t u u tr a tr u u st d tr u C C u a C C u a , , ) ( ) ( ⇒ ∑ ∑ + - Γ ∈ Γ ∈ < - t u t u u u u a u a , , 0 ) ( ) (
this is in

U

(see figure 13) be a cycle such that:

t u1 1 + t y 1 + t z 1 + t a 1 + t b 1 + t d t
) ( ) ( , , ≥ -∑ ∑ - + Γ ∈ Γ ∈ t y t y u u u a u a but given that ∑ ∑ - + Γ ∈ Γ ∈ + = < = i y t y u tr x tr y u C C u a u a , , ) ( 0 ) (
(the cost of arc t A1 is null) there is a contradiction. Therefore, in this case

t t x δ = . 2 nd case: if 0 ≠ t z and 0 = t y (if 0 ≠ t y
we return to case 1)

Let i t y + , ] , 1 [ t T i - ∈
, or

j t u + 1 , ] , 0 [ t T j - ∈
, and be the first non-zero variable. For example, if i t y + is the first non-zero variable then all m t y + and

n t u + 1 (with ] 1 , 1 [ - ∈ i m and ] 1 , 0 [ - ∈ i n ) are zero. 1 st sub-case: if i t y + the first non-zero variable: Let - + Γ Γ = Γ t z t z t z , , ,

U

(see figure 14) such that:

) , ,..., 1 ( ) , ,..., , , ( , 1 , i t i t t t z i t i t t t t t z Z D D Y A A X Z + + + + + + - = Γ = Γ ) 1 ( ) , ( , , t t y t t t y A X Y = Γ = Γ + - t z t d1 i t y + i t z + t x i t d + i t a + Figure 14 : Cycle - + Γ Γ t z t z , , U
Given that i t y + is the first non-zero variable, it follows that between t+1 and t+i no containers depart from the customer site, i.e. all arcs of - Γ t z, are nonzero. In addition,

t t x δ > which simply means that 0 )] ( ) ( [ , > - - Γ u b u f Min t z . Therefore, { } 0 )] ( ) ( [ ; )] ( ) ( [ , , > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min t z t z u u γ and given that ∑ ∑ - + Γ ∈ Γ ∈ + × + + = < + × = i z t z u tr y st a tr x tr z u tr z st d C C i C C u a C C i u a , , ) ( ) ( (Because st a st d C C < )
this is in contradiction with minimum cost network flow theorem. Therefore

+ + + + + + - = Γ = Γ
In the case where j=0 the tow paths are reduced to

) 1 ( ) 1 , , ( , , t t z t t t t z D U X Z = Γ = Γ + - t z t d1 j t y + j t z + t x j t d + j t a + j t u + 1 j t d + 1 Figure 15 : Cycle - + Γ Γ t z t z , , U
Because it is supposed that t t x δ > then, given that no containers depart from the customer until t+j, there are no zero arcs on the path - Γ t z, and the t X arcs have a value greater than the lower band (

t t x δ > ). Therefore, 0 )] ( ) ( [ , > - - Γ u b u f Min t z and since 0 )] ( ) ( [ , > - + Γ u f u c Min t z it follows that { } 0 )] ( ) ( [ ; )] ( ) ( [ , , > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min t z t z u u γ Since, according to minimum cost network flow theorem, 0 ) ( ) ( , , ≥ -∑ ∑ - + Γ ∈ Γ ∈ t z t z u u u a u a
is a condition for optimality, and given that

∑ ∑ - + Γ ∈ Γ ∈ + × + + = < × = i z t z u u a tr x tr z u d C C j C C u a C j u a , , ) ( ) (
, this is absurd. Therefore, in this case t t x δ = .

3 rd sub-case: if all i t y + and

j t u + 1 , ] , 1 [ t T i - ∈ and ] , 0 [ t T j - ∈
, are equal to zero, then after t all transportation variables will be equal to zero except

k k x and z , ] , [ T t k ∈
. Now, suppose that the feasible flow f found in this case is optimal. Let f' be another feasible flow obtained by transforming f by subtracting

t t x δ λ - = from ) 1 , ,..., 1 , , , ( T T t t t t A A A A X Z = Γ and adding it to ) 1 , ,..., 1 , ( ' T T t t D D D D = Γ
. This gives a gain of ) ) ( ) ((

d z x a C t T C C C t T × - - + + × - × λ which implies that f is not optimal. Therefore, t t x δ > is impossible ⇒ t t x δ = . Property 6. If for ] , 0 [ T t ∈ , ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ ≥ + δ δ then ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + ≤ δ δ Proof. Assume there exists a ] , 0 [ T t ∈ such that ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + > δ δ and let ] , 1 [ a d N i ∈
be the index corresponding to the position of maximum demand, in other words

i t k t N k a d + + ∈ = δ δ ) ( max ] , 1 [ . So, according to property 2, if 0 1 ≠ t u then ] , 1 [ a d N k ∈ ∀ k t z + =0. In other words, in the interval ] , 1 [ a d N t t + +
the system uses only the stock available at the customer site. However, this stock is strictly lower than

i t + δ because ) ( max 1 1 ] , 1 [ 1 k t N k t t t t a d u a a + ∈ + < - + = δ δ
and consequently the demand at t+i cannot be satisfied, which is absurd. Therefore,

) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + ≤ δ δ . Property 7. Let ] , 0 [ T t ∈ and i∈ ] , 1 [ a b N represent the position index of maximal demand: this means that i t k t N k a d + + ∈ = δ δ ) ( max ] , 1 [ with i>1. If ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ > + δ δ
then we have:

0 1 = + j t u for all ] 1 , 1 [ - ∈ i j , 0 ≠ +k t a for all ] , 1 [ i k ∈ and 0 1 ≠ +l t a for all ] 1 , 1 [ - ∈ i l Remark 6. If i=1 (i is at the position of maximum demand), then 0 1 ≠ + t a Proof. Let ] , 0 [ T t ∈ and i∈ ] , 1 [ a b N such that i t k t N k a d + + ∈ = δ δ ) ( max ] , 1 [ and assume that ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ > + δ δ . ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + ≤ δ δ
which means that the customer stock at period

t t t t u a a 1 1 1 - + = + δ is greater than or equal to ) ( max ] , 1 [ k t N k a d + ∈ δ , i.e. ) ( max 1 1 ] , 1 [ 1 k t N k t t t t a d u a a + ∈ + ≥ - + = δ δ
. Assume that there is a non-zero departure

j t u + 1 for ] 1 , 1 [ - ∈ i j
which is the first non-zero occurrence, i.e. all departures

m t u + 1 for ] 1 , 1 [ - ∈ j m are null. Let + - Γ ∪ Γ = Γ
be a cycle (see figure 16) such that

) 1 ,..., , 1 ( ) 1 ,..., 1 , ( 1 1 1 j t t t j t t t D D U U A A + + + + + + - = Γ = Γ j t y + j t z + t x t u1 j t u + 1 1 + t a 1 + t d j t a + j t d + j t a + 1 j t d + 1 1 1 + t a 1 1 + t d Figure 16 : Cycle - + Γ Γ U
Since all m t u + 1 are null, there are no departures from the customer site, this means that stock

l t a + for ] , 1 [ j l ∈
can never decrease and will remain greater than

) ( max ] , 1 [ k t N k a d + ∈ δ . In addition, for all ] , 1 [ j k ∈ k t y + can never exceed ) ( max ] , 1 [ k t N k a d + ∈ δ ⇒ 0 1 > = - + + + k t k t k t a y a .
Therefore, the path ) 1 ,..., 1 , (

1 1 j t t t U A A + + + -= Γ cannot contain a null arc, i.e. 0 )] ( ) ( [ > - - Γ u b u f Min , and because 0 )] ( ) ( [ > - + Γ u f u c Min
, it follows that:

{ } 0 )] ( ) ( [ ; )] ( ) ( [ > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min u u γ Since, ∑ ∑ - + Γ ∈ Γ ∈ + × = < × + = u tr u st a u st d tr u C C j u a C j C u a ) ( ) ( (because st a st d C C <
) this is absurd and therefore there is no departure

j t u + 1 for ] 1 , 1 [ - ∈ i j
, and the stock

k t a + for ] , 1 [ i k ∈ can never decrease, i.e. ) ( max ] , 1 [ n t N n k t a d a + ∈ + ≥ δ . In addition, given that l t l t n t N n y a d + + + ∈ ≥ > δ δ ) ( max ] , 1 [ for ] 1 , 1 [ - ∈ i l , it follows that 0 1 > - = + + + l t l t l t y a a .
Consequently, for

j t u + 1 , ] 1 , 1 [ - ∈ i j , 0 ≠ +k t a for all ] , 1 [ i k ∈ and 0 1 ≠ +k t a for all ] 1 , 1 [ - ∈ i l . Property 8. For a )] , max( , 0 [ a d a d M N T t - ∈ , if ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ > + δ δ , then ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + = δ δ Proof. Let )] , max( , 0 [ a d a d M N T t - ∈ such that ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ > + δ δ and let i∈ ] , 1 [ a b
N such that the maximum demand position be defined as

i t k t N k a d + + ∈ = δ δ ) ( max ] , 1 [ . Assume that ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + < δ δ
. In other words, for the customer stock at period t+1,

) ( max ] , 1 [ 1 k t N k t a d a + ∈ + > δ
and because there is no departure (according to property 8,

0 1 = + j t u for all ] 1 , 1 [ - ∈ i j
), the stock will stay the same or increase at the customer site in ] , 1

[ i t t + +
. Therefore, we will find that at period Given that the stock at the customer site will stay the same until the period t+j (the first departure is

i t + , ) ( max ] , 1 [ k t N k i t a d a + ∈ + >
j t u + 1 , ] , [ a d N i j ∈
), it follows that all m t a + and

m t a + 1 for ] , [ j i m ∈ are nonzero. Thus, 0 )] ( ) ( [ > - - Γ u b u f Min and since 0 )] ( ) ( [ > - + Γ u f u c Min this implies that { } 0 )] ( ) ( [ ; )] ( ) ( [ > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min u u γ . Given that: ∑ ∑ - + Γ ∈ Γ ∈ + × = < × + = u tr u st a u st d tr u C C j u a C j C u a ) ( ) ( (because st a st d C C <
) this is absurd.

Therefore, in this case

) ( max 1 1 ] , 1 [ k t N k t t
( ) ( [ > - - Γ u b u f Min . Since, 0 )] ( ) ( [ > - + Γ u f u c Min , it follows that { } 0 )] ( ) ( [ ; )] ( ) ( [ > - - = + - Γ ∈ Γ ∈ u f u c Min u b u f Min Min u u γ Given that ∑ ∑ - + Γ ∈ Γ ∈ + × = < × + = u tr u st a u st d tr u C C k u a C k C u a ) ( ) (
this is absurd. Therefore, the inequality

) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + < δ δ is impossible ⇒ ) ( max 1 1 ] , 1 [ k t N k t t
+ + + + + + + - = Γ = Γ Given that 0 1 ≠ +n t a for ] , 1 [ a d N n ∈
and that there are no departures between periods t+ a d N and t+k, it follows that stocks

m t a + , for ] , 1 [ k N m a d + ∈
, and . Let f be the flow found in this case and suppose that it is optimal with the assumption that

n t a + 1 , for all ] 1 , 1 [ - + ∈ k N n
) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + < δ δ . Let 0 1 ) ( max 1 ] , 1 [ > - - + = + ∈ t k t N
) ) ( ) (( ) ( ) ( < - - - - × ⇒ × - + > × - st d st u st a st d Tr st a C t T C C t T C t T C C t T u α . Therefore, the cost of f' is less than f and consequently if ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + < δ δ
the flow cannot be optimal and, therefore, the inequality

) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + < δ δ is impossible ⇒ ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + ≥ δ δ . According to property 7 ) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + ≤ δ δ
, and, therefore,

) ( max 1 1 ] , 1 [ k t N k t t t a d a u + ∈ - + = δ δ . Property 9. If, ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ ≤ + δ δ for a )] , max( , 0 [ a d a d M N T t - ∈ , then 0 1 = t u Proof. Let )] , max( , 0 [ a d a d M N T t - ∈ such that ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ ≤ + δ δ and assume that 0 1 ≠ t u . According, to property 2 if 0 1 ≠ t u then for all ] , 1 [ a d N k ∈ , 0 = + k t z . Given that ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ ≤ + δ δ
, the customer stock cannot satisfy the maximum demand

) ( max ] , 1 [ k t N k a d + ∈ δ . Consequently, there will exist a ] , 1 [ a d N k ∈ such that 0 ≠ +k t z
, which is absurd according to property 2. Therefore, for all )] , max( , 0 [

a d a d M N T t - ∈ if ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ ≤ + δ δ then 0 1 = t u . Proprerty 10. For a )] , max( , 0 [ a d a d M N T t - ∈ , if ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ ≤ + δ δ then )) 1 ( , 0 max( 1 1 t t t t a z δ δ + - = + + ) , 1 min( 1 1 + + + = t t t t a y δ δ Proof. Let )] , max( , 0 [ a d a d M N T t - ∈ such that ) ( max 1 ] , 1 [ k t N k t t a d a + ∈ ≤ + δ δ
, which, according to property 10, implies that 0 and that one of the two variables are non-zero occurrences.

1 = t u . Suppose that { } ) 1 ( , 0 max 1 1 t t t t a z δ δ + - > + + . Since 1 1 1 + + + + = t t t z y δ , it follows that 1 1 1 1 ) 1 ( ) 1 ( 1 + + + + + - + = - + = t t t t t t t t z a y a a δ δ δ and, given that { } ) 1 ( , 0 max 1 1 t t t t a z δ δ + - > + + , it follows that 0 ) 1 ( ) 1 ( ) 1 ( 1 1 1 1 1 1 = + - + - + > + - + = + + + + + t t
1 st sub-case: if 0 1 ≠ +k t u and ] , 1 [ T t k t + ∈ +
is the first to be non-zero, which implies that there are no departures between period t+1 and t+k (unless

Γ ∪ Γ = Γ be a cycle such that: ) 1 ,..., , 1 ( ) 1 , 1 ,..., 1 , ( 1 1 1 k t t t k t k t t t D D U U A A A + + + + + + + - = Γ = Γ 1 1 t t t t a z δ δ + - > + + is also impossible in this case. 2 nd case: if all k t u + 1 and k t y + are zero for T k t ≤ +
. Let f be the optimal flow found in this case and let 0

1 1 1 1 1 > - - + = = + + t t
C t T C C t T u × - + > × - ) ( ) ( ⇒ 0 ) ) ( ) ( ( < × - - × - + × st a st d tr u C t T C t T C β this implies that f' is better than f. Therefore, taking )) 1 ( , 0 max( 1 1 t t

Interpretation and Discussion

In our case, the supplier site has the largest storage cost. We have seen that storage at the supplier site should be avoided. This is obvious; in fact there is no interest in sending empty containers for storage at the supplier site, neither from the customer site nor the depot site. Keeping empty containers and sending them as needed is better then sending them before they are needed and suffering the additional costs involved. Concerning movement between the customer and the depot; at the end of each period there are always deliveries of loaded containers from supplier site t δ . After they have been emptied the customer is able to return them to the best place for storage. Although the depot has the lowest storage cost, we have seen that the quantity sent from the customer, to the depot in each period is different to the total quantity in customer possession ( ). In fact, the decision of how many containers to send from the customer site to the depot must not only take into account the different storage costs but also the different transportation costs. The transportation cost from the customer to the supplier Tr ) are greater than maximum demand). This is because if we send a quantity to the depot such that the remaining stock becomes ) ( max In order to satisfy the demand for empty containers at the supplier site, the system always uses the customer inventory first and uses depot inventory to make up any shortfall. This is because if we use the depot inventory such that ); in other words, we will obtain a surplus at the customer site which will create additional costs.

Conclusion

Throughout this paper we have demonstrated that in this cost configuration ( ) a simple recurring formula solves the model. This result can be extended to other configurations costs such as:

• the supplier has the lowest cost followed by the depot ( ).

• the customer site has the lowest cost independently from whether the depot or the supplier has the second lowest cost. • the supplier site has the lowest cost, independently from whether the depot or the customer has the second lowest cost.

The same reasoning can be applied to each of these situations. This shows that the problem can be divided into four categories and easily understood. This result is important for two reasons. Firstly, it simplifies the resolution of the container management problem. Secondly, it permits greater understanding of system behaviour and facilitates the interpretation of all movements. This can help us understand other situations, such as when there are several customer sites, depots and supplier sites. In addition, this new resolution method opens new horizons for the treatment of other problems related to flow optimization, such as system container dimensioning and container purchasing policies. The last two problems have been examined for certain particular cases. The cases have only looked at one site systems and have neglected the different movements in the transportation network flow (D.J. [START_REF] Buchanan | Optimal policy for a periodic review returnable inventory system[END_REF]. It is our intention in future studies to investigate this problem from the global aspect and demonstrate possible simplifications.

  storage cost in each site and depot.

  ) : represent the lower band of arc U u ∈ . All arcs have a lower band equal to zero except t X arcs which have t δ as lower bands.

]

  Figure 5 represent an illustration of definition 1. The first term of the first inequality

  This definition (see figure6) shows that for a number of periods, i.e. greater than a d M , storing stock at the customer site will be more expensive than sending it to the depot for storage without reuse. The two parameters a d M and a d N will help us to establish the different properties.

  Figure 6 : + Γ a d M t , and -Γ a d M t ,

  Figure 12 : Cycle

γ

  Figure 13 : Cycle -+ Γ Γ t y t y , , U

  Figure 18 : Cycle -+ Γ Γ U

yC

  the transportation cost from the customer to the depot and the transportation cost from the depot to the supplier ) ( the customer site, the system holds a quantity exactly equal to the maximum demand on the horizon of the next a d N periods for direct use (in the case when the entries (

Figure 1 : Flow of container

	Suppliers	Loaded container Empty container
	Customers	Depots
	Assumptions
	-The system is deterministic. The different costs for storing and moving the empty and
	loaded containers between sites and depots are known a priori. Similarly, the demands
	for each customer site are also known a priori.
	-There is no loss in the system
	-Container transition time between different sites and depots are negligible. Container
	loading and unloading times are also neglected.
	-There is one type of container. All containers are homogeneous
	-Two movements of empty containers can be done between each pair of sites/depots or
	customers sites/ suppliers sites during the same period. The first takes an empty
	container to be loaded and sent to customers. The second takes place when the
	containers are emptied and it attempts to return empty containers to suitable locations
	for the next usage.
	-There are two types of inventory in sites and depots. The principal inventory, which is
	the inventory before the first movement of empty containers for loading, and the
	temporary inventory, which is the inventory after the movement of empty containers
	for loading. Only the cost of the principal inventory is considered in our model.
	-The storage capacities of sites and depots are unlimited
	-All demands are nonzero and are different

Figure 2 : Network flow

  : the principal inventory of containers in customer site i at period t. •

		t fi x : Number of loaded containers transported from suppliers f to a customer i at
		period t
	•	t i x : Number of containers needed in customer site i at period t
	•	t if y : Number of empty containers transported from customer i to supplier f at
		period t to be loaded and sent to customers.
	•	t if y 1 : Number of empty containers transported from customer i to supplier f at
		period t to be returned to the suitable location.
	•	z	t df	: Number of empty containers transported from depot d to suppliers f at
		period t to be loaded and sent to customers.
	•	z 1	t df	: Number of empty containers transported from depot d to suppliers f at
		period t to be returned to the suitable location.
	•	t id u : Number of empty containers transported from customer i to depot d at
		period t to be loaded and sent to customers.
	•	t id u 1 : Number of empty containers transported from customer i to depot d at
		period t for them to be returned to the suitable location.
	•		t i
		t i a 1 : the temporary inventory of containers in customer site i at period t.
	•	d	t d		: the principal inventory of containers in depot d at period t.
	•	t d d 1 : the temporary inventory of containers in depot d at period t.
	•		t f
		a	t i							t i a1		a	t i	+	1
												x	t i
					b	t f	y	t if	x	t fi	t f b1	t if y1	b	t f	+	1
					z	t df							t df z1
									t id u				t id u1
		d	t d								t d d1	d	1 + d t

a b

: the principal inventory of containers in supplier site f at period t. • t f b1 : the temporary inventory of containers in supplier site f at period t.

  The configuration of cost is as follow

	C	st d	<	C	st a	<	C	st b
	• the initial stock of empty containers at the supplier site is equal to zero
									x	t	Suppliers	Loaded container Empty container
									z	t	+	t 1 z
		a	t					y	t	+	t 1 y
				Customers	u	t	+	t 1 u	Depots

t d t b Figure 3 : Decision variables

  consequently

	C	Tr y	<	C	Tr u	+	C	Tr z	. Otherwise,	C	Tr y	>	C	Tr u	+	C	Tr z	would give	C	Tr y	+	C	st a	>	C	Tr u	+	C	st d	+	C	Tr z

m t u + 1 are null. Let + -

+ t y ), then the stock at customer site will be non-zero between these two periods. Let + -

Given that there are no departures between period t+1 and t+k (unless 1 + t y ), and that 0 1 ≠ +k t u

, it follows that all arcs of path ) 1 , ,..., 1 , (

, this is absurd. Therefore, in this case the inequality

is the first to be non-zero, then there are no departures from the customer site between periods t+2 and t+l. Let

, this is absurd and consequently )) 1 ( , 0 max(