C Staquet

P Bouruet-Aubertot

Mixing in weakly turbulent stably stratified flows

Keywords: Stably stratified flows, Turbulence, Mixing

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

Introduction

In a stably stratified inviscid fluid, the restoring action of the buoyancy force prevents vertical transport of mass from occuring: any fluid particle displaced vertically from its equilibrium position will oscillate around this position so that the net vertical displacement of the particle over a phase period will be zero. It follows that vertical transport of mass in a stably stratified fluid can only result from molecular diffusion, that is, from mixing. This does not mean that the vertical diffusivity of density is equal to the molecular diffusivity at all times: turbulent motions put into contact fluid particles with each other so that diffusive exchange of mass (which is roughly proportional to the density difference between the fluid particles) will be promoted. Fluid mixing is therefore enhanced, if the diffusive time scale of the motions which advect the fluid particles is not too large compared to the advective time scale.

The purpose of this note is to show that weakly turbulent stably stratified flows, in which the dissipation rate of kinetic energy is controlled by horizontal motions with a vertical shear, share a common feature when their mixing properties are compared. We focus upon three different flows, namely low amplitude breaking internal gravity waves, a three-dimensional homogeneous stably-stratified turbulence and a stably stratified shear layer, in two and three dimensions. The Boussinesq equations are solved for this purpose, using direct numerical simulations. Mixing is analysed using the method derived by [START_REF] Winters | Available potential energy and mixing in density-stratified fluids[END_REF] and [START_REF] Winters | Diascalar flux and the rate of fluid mixing[END_REF], which allows one to compute exactly the diffusive flux of density resulting in fluid mixing. This flux occurs across surfaces of constant density (rather than vertically) and is referred to as the diapycnal flux.

We shall show that in the three flows we consider, the diapycnal diffusivity displays the same dependency as a function of a dynamical parameter that characterizes the largest scales of the flow and we provide a simple expression of this law. The dynamical parameter is a Froude number when there is no background shear and a Richardson number when there is a mean shear.

The outline of the paper is the following. In Section 2, the method of mixing analysis is described, along with the numerical method we use to solve the equations of motions. Sections 3-5 successively report about the three flows we consider: in each section, we briefly specify the equations of motion that drive the flow dynamics, provide the striking features of those dynamics and focus upon the common mixing law shared by the three flows. Conclusions are drawn in Section 6.

A subset of the results presented in Sections 4 and 5 can be found in [START_REF] Godeferd | Statistical modelling and direct numerical simulations of decaying stably-stratified turbulence: Part 2. Large and small scale anisotropy[END_REF] and [START_REF] Staquet | Mixing in a stably stratified shear layer: two-and three-dimensional numerical experiments[END_REF], respectively, while Section 3 complements results obtained by [START_REF] Bouruet-Aubertot | Particle dispersion and mixing by breaking internal gravity waves: two-dimensional numerical experiments[END_REF]. The purpose of the present paper, in gathering these results and in further analysing them within a common frame, is to suggest that a characteristic law may exist between mixing and the large scale dynamics in weakly turbulent stably stratified flows.

Methodology

Analysis of mixing

The quantification of mixing in a stably stratified fluid relies upon the estimate of a mixing rate, which is the rate of increase per time unit of potential energy resulting from fluid mixing. This potential energy is denoted E b in the following. It is only a part of the total potential energy of the flow, E p , and was characterized by [START_REF] Lorenz | Available potential energy and the maintenance of the general circulation[END_REF] as being its instantaneous minimum value: E b is the potential energy of the fluid when all motions have ceased and the fluid has returned to rest. [START_REF] Lombard | Energetics of a stably stratified mixing layer[END_REF] showed that E b is the potential energy of a stable density profile introduced by [START_REF] Thorpe | Turbulence and mixing in a Scottish Loch[END_REF], obtained by sorting adiabatically the fluid particles so that the heaviest fluid particle has the lowest altitude. This sorted profile, which we shall refer to as ρ b (z), is stable by construction. The third step, realized by [START_REF] Lombard | Energetics of a stably stratified mixing layer[END_REF], was to compute the rate of change of E b (i.e. Ėb ) from the sorted profile. This problem was solved by [START_REF] Winters | Available potential energy and mixing in density-stratified fluids[END_REF] (1)

Among the three flows we consider, the surface fluxes are zero for the shear layer only.

A local diapycnal flux φ d (z, t) may also be defined, whose volume average is equal to Φ d /g, g being the modulus of the acceleration of gravity g. φ d (z, t) is expressed as [START_REF] Winters | Diascalar flux and the rate of fluid mixing[END_REF]

φ d (z, t) =-κ |∇ρ| 2 I dρ b /dz , (2)
where the average • I is performed on an isopycnal surface; z denotes the vertical coordinate and κ the molecular diffusivity of density.

When the fluid has returned to rest (if there is no shear) or to a hydrostatic balanced laminar state (if there is a horizontal mean shear), the diffusive flux of density results from purely molecular effects. We shall refer to this flux as Φ i . Its expression is

Φ i = κg ρ H , ( 3
)
where ρ is the density difference across the fluid domain and H is the domain height.

Φ d -Φ i is therefore the fluctuating diapycnal flux. Φ i is the minimum value of Φ d so that Φ d -Φ i is always positive.
Since the background potential energy E b results from the irreversible process of mixing, E p -E b defines the potential energy associated with reversible processes. The latter energy is usually referred to as the available potential energy and denoted E a . It is the potential energy of the density fluctuations about the sorted density profile. From the equation for the total potential energy

Ėp = Φ a + Φ i + Ėsurf adv p + Ėsurf dif p , (4)
and from Eq. ( 1), one infers immediately the equation that governs the rate of change of

E a Ėa = Φ a -(Φ d -Φ i ) + ( Ėsurf adv p -Ėsurf adv b ) + ( Ėsurf dif p -Ėsurf dif b ), (5)
where Φ a is the vertical advective flux of density; the last two terms in Eqs. ( 4) and ( 5) account for the rate of change of E p and E a , respectively, due to advective and surface mass fluxes across the bounding surface of the fluid domain. Eq. ( 5) shows that, in the fluid interior, a part of E a is reversibly converted into vertical kinetic energy, at a rate equal to Φ a while another part is irreversibly converted into fluctuating background potential energy, at a rate equal to Φ d -Φ i [START_REF] Winters | Available potential energy and mixing in density-stratified fluids[END_REF]. Several estimators of mixing have been proposed in the literature (see, for instance Fernando, 1991) and we shall focus in the present paper upon two of them. These are the diapycnal diffusivity and the mixing efficiency, which we now define.

As any diffusion coefficient, the diapycnal diffusivity K ρ is defined by the ratio of the flux of scalar along the background scalar gradient (φ d ) divided by this scalar gradient (-dρ b /dz). Hence [START_REF] Winters | Diascalar flux and the rate of fluid mixing[END_REF])

K ρ = φ d -dρ b /dz . ( 6
)
The efficiency of mixing (within the interior of the fluid domain) accounts for the amount of energy gained per time unit through the mixing process, Φ d , relative to the amount of energy brought to the fluid per time unit, Ṁ say. When the fluid system is not forced, which is the case in the three flows we consider, the dissipation rate of kinetic energy ǫ is used in place of Ṁ. Also, only mixing resulting from the flow dynamics, that is from density fluctuations about the sorted density profile, are of interest so that the fluctuating diapycnal flux Φ d -Φ i should be employed in place of Φ d . The definition of the mixing efficiency we shall use is thus

γ = Φ d -Φ i ǫ . (7)

Mathematical model and numerical method

The dynamics of the flows we consider are driven by the Navier-Stokes equations in the Boussinesq approximation: density fluctuations about a hydrostatic balanced state are supposed to be small enough (in magnitude and frequency) for the fluid motions to be incompressible. The effect of these density fluctuations upon the fluid motions comes into play through a buoyancy force (difference between the gravity and the Archimedean forces). These equations display a different form for each flow we consider so that their expression will be written in the next sections. The boundary conditions are either periodic or of the free slip type along the horizontal boundaries.

The Boussinesq equations are solved by a classical pseudo-spectral method (e.g. [START_REF] Canuto | Spectral Methods in Fluid Dynamics[END_REF]: the spatial derivatives are computed using Fourier transforms while the nonlinear terms are computed in physical space. The latter computation produces aliasing errors, which are eliminated in Fourier space by a standard truncation procedure. Time is advanced using a third order Adams-Bashforth scheme. When the boundary conditions are of the free slip type along a boundary, the numerical domain is doubled about that boundary assuming symmetry or antisymmetry properties of the fields about it. These fields become periodic upon the doubled domain and can be expanded in Fourier series. The specific properties of the associated Fourier coefficients are used to reduce the computation time, following a now classical method proposed by [START_REF] Orszag | Numerical simulation of incompressible flows within simple boundaries. I. Galerkin (spectral) representation[END_REF].

The analysis of mixing described in the previous section relies upon a sorted density profile, whose computation is performed as follows. The value of the density field at each grid point is identified with one fluid particle. These fluid particles are next sorted out with their density. The sorted density profile is constructed by filling out each horizontal plane with the heaviest fluid particles, starting with the bottom plane of the numerical domain. Because the density slightly varies from particle to particle, it is averaged over each plane to yield a one-dimensional profile. The word plane should be replaced by line when the flow evolves in two dimensions.

Breaking internal gravity waves of low amplitude

Overall flow dynamics

We consider a stably stratified fluid in which the Brunt-Väisälä frequency

N = - g ρ 0 dρ lin dz 1/2 (8)
has a constant value. ρ 0 + ρ lin is the density of a hydrostatic balanced state, ρ 0 being a constant reference density and ρ lin (z) a linear density profile. A monochromatic, high frequency (that is close to N) internal gravity wave of low amplitude, hereafter referred to as the primary wave, propagates in this fluid. Low amplitude means that the isopycnals are only slightly tilted with respect to the horizontal and far from overturning. We have studied the dynamics of this wave in the vertical plane formed by g and its wave vector k. The restriction to two-dimensional dynamics is motivated by a linear stability analysis by [START_REF] Klostermeyer | Two and three-dimensional parametric instabilities in finite amplitude internal gravity waves[END_REF], which predicts that such a low amplitude wave is most unstable to two-dimensional perturbations in the ( g, k) plane. Comparison of the two-dimensional wave dynamics with a three-dimensional computation [START_REF] Bouruet-Aubertot | Particle dispersion and mixing by breaking internal gravity waves: two-dimensional numerical experiments[END_REF] and with laboratory experiments [START_REF] Benielli | Excitation of internal waves and stratified turbulence by parametric instability[END_REF] confirms the validity of the two-dimensional approach when flow statistics are compared.

Let us consider a coordinate system (0,x,y,z), where (x, z) is the vertical plane of the primary wave with z pointing upwards. The Boussinesq equations in this plane are

∂∇ 2 ψ ∂t + J (∇ 2 ψ, ψ) = ∂ρ ′ ∂x + ν(∇ 2 ) 2 ψ (9) ∂ρ ′ ∂t + J (ρ ′ ,ψ) =-N 2 ∂ψ ∂x + ν Pr ∇ 2 ρ ′ . ( 10
)
As is customary, the use of a stream function ψ such that u =-∂ψ/∂z and w = ∂ψ/∂x implies that the incompressibility condition is satisfied implicitly. ρ ′ (x,z,t)is the density fluctuation about the hydrostatic density profile mentioned above, scaled as an acceleration:

ρ ′ (x,z,t) = (g/ρ 0 )(ρ(x, z, t) -(ρ 0 + ρ lin (z))
, where ρ(x,z,t) denotes the total density. J is the nonlinear Jacobian operator, ν the kinematic viscosity and

Pr = ν κ (11)
is the Prandtl number of the fluid (κ has been defined in the previous section).

The fluid domain is a square box and periodic boundary conditions are used along the two directions. Eqs. ( 9) and ( 10) are initialized with a propagating internal gravity wave (the primary wave) of wave vector k = (1, 0, 1), Froude number A(| k|/2π) 2 /N equal to 0.013 and Reynolds number Re = A/ν equal to 10700. A is the initial amplitude of the stream function of the wave and N is set to 1. Four computations will be presented in this section, with resolution 512 2 , which differ only in the value of the Prandtl number: Pr = 0.72, 1, 5 and 10.

Whether propagating [START_REF] Koudella | Two-dimensional instability mechanisms of a progressive internal gravity wave[END_REF] or standing [START_REF] Bouruet-Aubertot | Breaking of internal gravity waves trough two-dimensional instabilities[END_REF], the dynamics of the wave display the same behaviour: it becomes unstable through parametric instability among resonant triads, a result predicted theoretically by [START_REF] Drazin | On the instability of an internal gravity wave[END_REF] and [START_REF] Mied | The occurence of parametric instabilities in finite amplitude internal gravity waves[END_REF] for a propagating wave. The damping action of viscosity at small scales selects a range of unstable triads, among which the numerical grid operates a second selection. One particular triad eventually grows, which controls the whole instability. The amplification of the resonant waves within the triad, through extraction of energy from the primary wave, eventually makes the latter wave break. Breaking means here that the primary wave loses its coherence (but it does not grow in amplitude). A turbulent regime sets in, where the kinetic energy and potential energy spectra display a k -3 z law, k z being the vertical wavenumber. This regime lasts a few Brunt-Väisälä periods (2π/N) only because the primary wave is not forced (see [START_REF] Bouruet-Aubertot | Stratified turbulence produced by internal wave breaking: two-dimensional numerical experiments[END_REF] for standing waves; Koudella, 1999 for propagating waves). The breaking of the wave is illustrated in Fig. 1 for a propagating primary wave.

Analysis of mixing

Wavebreaking in an unforced uniformly stably stratified fluid appears to result in weak fluid mixing over the duration the (laboratory or numerical) experiments are conducted (e.g. [START_REF] Mcewan | Internal mixing in stratified fluids[END_REF][START_REF] Taylor | The energetics of breaking events in a resonantly forced internal wave field[END_REF][START_REF] Benielli | Excitation of internal waves and stratified turbulence by parametric instability[END_REF][START_REF] Bouruet-Aubertot | Particle dispersion and mixing by breaking internal gravity waves: two-dimensional numerical experiments[END_REF]. The reason is that the energy input to the fluid, a part of which only will be available for mixing, is much smaller than the potential energy of the background density profile. This point may be expressed alternately by noting that the vertical scale of the wave instability is much smaller than that of the background density profile. In laboratory experiments of an internal gravity wave propagating along an interface by contrast, the breaking of the wave may destroy the interface when the wave amplitude becomes comparable to the thickness of that interface (e.g. [START_REF] Hannoun | Turbulent mixing at a shear-free density interface[END_REF]. It follows that, when unforced wavebreaking occurs in a linearly stratified fluid, the sorted density profile only slightly departs from the linear profile: relatives changes are smaller than 1% and are associated with modulations of dρ b /dz over thin layers, whose thickness has been estimated in [START_REF] Bouruet-Aubertot | Particle dispersion and mixing by breaking internal gravity waves: two-dimensional numerical experiments[END_REF]. These changes, which would initiate the formation of density steps if the wave were forced, are ignored in the following. We shall thus consider that ρ b (z, t) ≃ ρ lin (z) (this approximation being further discussed at the end of this section). The diapycnal diffusivity (6) therefore takes the simpler expression

K ρ (z, t) ≃ g ρ 0 φ d (z, t) N 2 . ( 12
)
We shall focus in this section upon the relation between the diapycnal diffusivity (averaged over the height of the fluid domain) and a Froude number associated with the largest scales of the flow. The Froude number we consider is defined as follows:

Fr t = √ 2Z N = ǫ μN 2 1/2 , ( 13
)
where Z is half the domain averaged vorticity squared and ǫ = 2μZ denotes the dissipation rate of kinetic energy per unit area (μ = ρ 0 ν is the dynamical viscosity). In ordinary turbulence, such a Froude number is associated with the smallest scales of the flow (e.g. [START_REF] Imberger | Application of the Wigner-Ville distribution to temperature gradient microstructure: a new technique to study small-scale variations[END_REF]). In the present case where the kinetic energy spectrum evolves as k -3 z , however, the enstrophy spectrum evolves as k -1 z , meaning that the volume averaged enstrophy is controlled by the largest scales of the flow. Hence, in stably stratified turbulence resulting from internal gravity wave breaking, Fr t is a dynamical parameter associated with the largest scales.

This turbulent regime is statistically homogeneous along the vertical direction so that the diapycnal diffusivity K ρ may be averaged along the vertical direction. Normalizing this averaged value by the molecular diffusivity κ yields an estimate of the Cox number ( [START_REF] Osborn | Oceanic fine structure[END_REF])

Cox ≡ 1 H H 0 K ρ κ dz = Φ d Φ i , (14)
where Φ i has the simple expression κρ 0 N 2 in the present case of a uniformly stably stratified fluid and H is the domain height, as before.

The instantaneous values of the Cox number and of the turbulent Froude number have been computed at all times of the runs having Pr = 0.72, 5 and 10. Cox(t) is plotted as a function of Fr 2 t (t) in Fig. 2, once breaking has occurred. A linear law is displayed for the smallest values of Fr 2 t , which can be modelled by the following simple expression. From definition (7) for γ , one has: Φ d = ǫγ + Φ i . Expressing ǫ as a function of the turbulent Froude number (13), we find

Φ d = γμN 2 Fr 2 t + Φ i . Dividing by Φ i yields Cox = γ Pr Fr 2 t + 1, ( 15
)
where Pr is the Prandtl number. The coefficient of the linear law inferred from Fig. 2 is thus γ Pr. It yields a value for γ equal to ≃1.2 for Pr = 0.72, ≃0.5 for Pr = 5 and ≃0.35 for Pr = 10. The validity of expression (15) to account for the linear law displayed in Fig. 2 relies upon the mixing efficiency to have reached a constant value once breaking has occurred. γ is plotted as a function of time in Fig. 3: it seems indeed to relax toward an asymptotic value after breaking, with the approach faster as the Prandtl number increases. The estimated asymptotic values are lower than those inferred from Fig. 2 using the law (15) and an argument to account for this behavior is provided in the last paragraph of this section.

The existence of final steady regime also appears when the rate of production Φ a of the available potential energy is compared with the rate of destruction Φ d -Φ i (Fig. 4). Both We finally comment upon the initial behavior of γ displayed in Fig. 3: γ has a constant value, which is Prandtl number dependent. This is consistent with the fact that, for a low amplitude internal gravity wave propagating in a linearly inviscid stratified fluid in an infinite domain, γ is equal to 1/Pr; the sorted density profile and the linear density profile coincide at any time in this case because the flow is inviscid. Fig. 3 shows that the constant initial values of γ actually underestimate by 37% the theoretical prediction 1/Pr. The reason is that the sorted density profile does not coincide with the linear density profile near the horizontal boundaries of the fluid domain. A strong density gradient is observed instead at the horizontal boundaries, which modifies Φ d and, consequently, γ . The steepening of the density gradient near the boundaries is caused by the application of the sorting process in a bounded domain with periodic boundary conditions: at the lower boundary for instance, positive density fluctuations correspond to particles which would have been displaced upward in a infinite domain, thus stemming from a lower altitude, where the fluid is denser. The influence of this effect weakens when breaking occurs because strong gradients are created in the fluid interior as well.

Stably stratified homogeneous turbulence

Overall flow dynamics

Stably stratified homogeneous turbulence with a uniform density (or temperature) profile has been studied by many authors in the past, using laboratory and numerical experiments. In the latter case, a stably stratified turbulent flow is created by a method that has become standard [START_REF] Riley | Direct numerical simulations of homogeneous turbulence in density-stratified fluids[END_REF]: a velocity field is initialized in spectral (Fourier) space through a prescribed kinetic energy spectrum (usually peaked about a given wavenumber that sets the scale at which energy is injected) and evolves temporally up to the time t 0 the triple correlations have built. (In practice, this time is estimated by the skewness factor to have reached a value of order 0.45 or so, according to laboratory experiments, implying that a fully turbulent flow has developed.) The velocity field at t = t 0 is used as the initial condition of a stably stratified flow, created by setting the Brunt-Väisälä frequency to a non-zero value at t = t 0 . A density field may also be imposed at this time. In the simulations we present below, this density field is set so that it has the same spectral distribution as the velocity field at t = t 0 and the potential energy is half the kinetic energy of the flow at t = t 0 .

We solve the three-dimensional Boussinesq equations in a cubic numerical domain with periodic boundary conditions. These equations are

∂ u ∂t + ( u •∇) u =- 1 ρ 0 ∇P -ρ ′ i z + ν∇ 2 u, ( 16
) ∂ρ ′ ∂t + ( u •∇)ρ ′ = N 2 w + ν Pr ∇ 2 ρ ′ , (17)
∇• u = 0. ( 18
)
The notations are the same as in the previous section: u = (u, v, w) is the velocity field, P denotes the pressure field and ρ ′ (x,y,z,t)refers to the field of density fluctuations about a hydrostatic balanced state, scaled like an acceleration. ν and Pr, respectively, denote the kinematic viscosity and the Prandtl number of the fluid (defined by ( 11)). The numerical method has been described in Section 2.2.

A number of runs initialized by the method described above has been performed, at resolution 256 3 , in which N is equal to π or π/2 and 1/ν is varied from 150 to 400, the Prandtl number having the constant value of 1. Two runs will be presented in this section, corresponding to N = π/2 and N = π, the viscosity having the value 1/300 in both cases. The dynamics of these flows is best described in terms of three dimensional parameters, namely a Froude number Fr, a Reynolds number Re λ based upon the Taylor microscale λ and the Prandtl number Pr. In laboratory experiments (e.g. [START_REF] Lienhard | The decay of turbulence in thermally stratified flow[END_REF][START_REF] Yoon | The evolution of grid-generated turbulence under conditions of stable thermal stratification[END_REF], the Froude number is commonly defined by the ratio of two length scales, a buoyancy length scale L ρ = w ′ /N , which is the greatest distance a fluid particle can move against the density gradient and an integral length scale L = u ′3 /ǫ, which characterizes the scale of the largest energetic structures. u ′ and w ′ denote the rms value of the u and w velocity components, respectively. Thus

Fr = L ρ L = ǫ N w ′ u ′3 , Re λ = u ′ λ ν . ( 19
)
As is well-known, the buoyancy force stabilizes the flow, in converting vertical kinetic energy into potential energy, the latter energy being involved in irreversible processes such as nonlinear transfers toward small scales (at a faster rate than kinetic energy, e.g. [START_REF] Holloway | The buoyancy flux from internal gravity wave breaking[END_REF][START_REF] Lesieur | Turbulence in Fluids, 3rd Edition[END_REF][START_REF] Staquet | Statistical modelling and direct numerical simulations of decaying stably-stratified turbulence: I. Flow Energetics[END_REF] and molecular mixing. Consequently, as time evolves, the flow organizes into quasi-horizontal layers (e.g. Metais and Herring, 1989). This organization of the flow leads to a strong variability of the velocity components along the vertical direction, which is best put forward by considering the vertical derivative of a velocity component. An iso-surface of ∂v/∂z is thus displayed in Fig. 5 (see caption) at the end of the simulation with N = π. Since more energy is left in the horizontal plane, a horizontal velocity component is considered. This horizontal velocity component may be arbitrarily chosen because the flow is statistically isotropic about the vertical direction (so that u ≃ v).

The temporal evolution of Re λ is displayed in Fig. 6 for the two runs discussed in this section. Re λ starts from a value close to 45 and, from Nt/2π ≃ 3, relaxes toward a nearly constant value. This value depends upon N , which shows that the evolution of Re λ is not solely driven by buoyancy effects; the reason is that u ′ is mostly contributed by (large scale) horizontal motions. We found that other quantities relax toward a constant value from this time on, such as the ratio u ′ /w ′ [START_REF] Godeferd | Statistical modelling and direct numerical simulations of decaying stably-stratified turbulence: Part 2. Large and small scale anisotropy[END_REF]. This ratio remains close to the value of 1.4 from about Nt/2π ≃ 2, in good agreement with laboratory experiments [START_REF] Lienhard | The decay of turbulence in thermally stratified flow[END_REF][START_REF] Yoon | The evolution of grid-generated turbulence under conditions of stable thermal stratification[END_REF]. This suggests that a weakly turbulent, and even quasi-balanced, dynamical state may have been reached.

Analysis of mixing

As discussed in the previous section, the sorted density profile is close at all times to the linear density profile in such a homogeneous flow so that the available potential energy E a may be approximated by the variance of the density fluctuations about either profile: E a ≃ 0.5 ρ ′2 /N 2 , where • denotes a volume average (the validity of this approximation is discussed in [START_REF] Holliday | On potential energy density in an incompressible stratified fluid[END_REF]. It follows that a simple expression for the fluctuating diapycnal flux Φ d -Φ i is provided by the dissipation rate of the fluctuating density variance χ = κ |∇ρ ′ | 2 divided by N 2 : Φ d -Φ i ≃ χ/N 2 . This estimate is commonly used in the ocean (see, e.g. Toole, 1998 for a review). In this context, the mixing efficiency γ should be defined by γ h = χ/(N 2 ǫ), the volume averaged diapycnal Fig. 7a displays the rate of production Φ a and the rate of destruction Φ d -Φ i of the available potential energy. Both quantities nearly coincide during the final regime of the flow which attests that a quasi-steady regime has been reached.

The temporal evolution of γ h is plotted in Fig. 7b for the two runs we consider. As opposed to Re λ , the curves nearly coincide which shows that the evolution of γ h is driven by the accumulated effect of buoyancy. Apart from the very initial stage of the computation where a0.7 value is reached, γ h remains close to 0.45 at all times while slowly decaying.

In order to relate the mixing properties of the flow to the large scale dynamics, we express the Cox number Cox h as a function of the Froude number Fr defined by Eq. ( 19). When a dependency as a function of Fr 2 is sought, a simple expression is obtained

Cox h = Pr 20 u ′ w ′ 2 γ h (Fr Re λ ) 2 . ( 20
)
This expression involves Re λ and u ′ /w ′ , in addition to γ h (we recall that Pr = 1). The former three quantities relax toward a constant value as time elapses, implying that the Cox number should behave linearly upon Fr 2 when Fr is small enough. This is attested in Fig. 8, where Cox h is plotted as a function of (Fr Re λ ) 2 for the runs we consider. The slope of the linear dependency is about 0.045, in good agreement with the large time value of its theoretical expression (taking Pr = 1, γ h = 0.45 and u ′ /w ′ = 1.4 yields Pr/20(u ′ /w ′ ) 2 γ h = 0.044). One interest of this linear dependency is that, conversely, it provides an estimate of the mixing efficiency if the other parameters are known.

Stably stratified shear layer in two and three dimensions

Overall flow dynamics

The shear layer is the simplest case of an inhomogeneous flow. It is formed by two layers of fluids having different densities and horizontal velocities. The interface between the layers, which thickens by molecular diffusion, enlarges most efficiently if the interface is unstable. A necessary condition for instability is the Richardson number

J = min z g ρ 0 -dρ/dz dU/dz 2 (21)
to have a value smaller than 0.25 somewhere in the fluid (see, e.g. [START_REF] Drazin | Hydrodynamic Stability[END_REF]. The overbar in Eq. ( 21) refers to a horizontal average. If this condition is met, the Kelvin-Helmholtz instability may develop in the vertical plane of the unstable shear flow.

It is therefore a two-dimensional instability. As its amplitude grows, nonlinear effects gain importance and makes it saturate while the flow vorticity organises into quasi-horizontal billows. The billows are most unstable to a (two-dimensional) subharmonic instability which makes them pair (see, e.g. [START_REF] Ho | Perturbed free shear layers[END_REF] for a review). Secondary three-dimensional instabilities may also occur, at smaller scales than the two-dimensional instabilities. Further detail about the three-dimensional instabilities in a stably stratified shear layer may be found in, e.g. [START_REF] Schowalter | A study of streamwise vortex structure in a stratified shear layer[END_REF], [START_REF] Klaassen | The influence of stratification on secondary instability in free shear layers[END_REF] or [START_REF] Caulfield | Three-dimensionalization of the stratified mixing layer[END_REF].

To investigate the shear layer dynamics and resulting mixing properties, we solve the Boussinesq equations in two and three dimensions. Non-dimensional variables are used: we take half the thickness of the initial velocity profile as a length scale, which we denote as δ i , half the velocity difference across the flow as a velocity scale (U ) and half the density difference across the flow times √ Pr, as a density scale ( ρ √ Pr/2). In three dimensions, the resulting non dimensional equations are

∂ u ∂t + ( u •∇) u =-∇P + Jρ + 1 Re ∇ 2 u, (22)
∂ρ ∂t + ( u •∇)ρ = 1 Pr Re ∇ 2 ρ, (23)
∇ u = 0. ( 24
)
As in the previous section, u = (u, v, w) denotes the velocity field with u being the streamwise component and w the vertical one; P and ρ, respectively, are the dynamical pressure deviation and density deviation about a hydrostatic balanced state of reference density ρ 0 . J , Re and Pr are the three parameters, set by the scaling mentioned above, upon which the flow dynamics depend. J and Re are expressed as (Pr is defined by Eq. ( 11))

J = g ρ 0 ρ √ Pr 2 δ i U 2 , Re = Uδ i ν . ( 25
)
The initial mean (horizontally averaged) velocity and density fields are of the error function type, the thickness of the density profile being √ Pr larger than that of the velocity profile. The velocity field is destabilized by two-and three-dimensional perturbations. The boundary conditions are periodic along the streamwise (x) and spanwise (y) directions and of free slip type along the vertical. The same numerical method as in the previous sections is employed.

We have performed several runs, in two and three dimensions, in which J ranges between 0.04 and 0.16 and Pr has a value of 0.7 or 1.4. The resolution is (N x ,N z ) = (1024, 1025) in two dimensions implying that a high Reynolds number can be chosen (Re = 1000 or 2000); in three dimensions, the resolution is (N x ,N y ,N z ) = (128, 128, 257) and Re ≃ 300. (N i denotes the number of grid points along direction i).

An example of the flow development in two and three dimensions is displayed in Figs. 9 and 10, respectively.

Analysis of mixing

The mixing properties of a stably stratified shear layer have been studied by a few authors [START_REF] Koop | Instability and turbulence in a stratified fluid with shear[END_REF][START_REF] Scinocca | The mixing of mass and momentum by Kelvin-Helmholz billows[END_REF][START_REF] Cortesi | Numerical investigation of the entrainment and mixing processes in neutral and stably stratified shear layers[END_REF][START_REF] Staquet | Mixing in a stably-stratified shear layer[END_REF][START_REF] Staquet | Mixing in a stably stratified shear layer: two-and three-dimensional numerical experiments[END_REF][START_REF] Caulfield | The anatomy of the mixing transition in homogeneous and stratified free shear layers[END_REF]. Only in the three latter papers is the method derived by [START_REF] Winters | Available potential energy and mixing in density-stratified fluids[END_REF] used.

In such a vertically inhomogeneous flow, the sorted density profile strongly differs from the mean (horizontally averaged) density profile because the vertical scale of the instability is comparable to that of the vertical density gradient [START_REF] Staquet | Mixing in a stably stratified shear layer: two-and three-dimensional numerical experiments[END_REF]. Hence, as opposed to the previous sections, quantitative information into the mixing properties of the flow can only be gained by using the exact expression of the diapycnal flux, which relies upon the sorted density profile. The mixing efficiency γ defined by Eq. ( 7) is plotted in Fig. 11a for two-dimensional simulations and in Fig. 11b for three-dimensional simulations. The curves display similar features during the development of the Kelvin-Helmholtz instability (for t ≤ 45) and during the later stage of the flow (for t larger than ≃150). By contrast, during the intermediate turbulent regime, the two-dimensional simulations are unable to reproduce the mixing efficiency predicted the three-dimensional simulations because the conservation of enstrophy by nonlinear terms in two dimensions prevents kinetic energy from being efficiently transferred toward small scales. As a consequence, the dissipation rate of kinetic energy has a value too low relative to its three-dimensional counterpart, which overestimates γ .

As in the previous sections, we investigate the relation between the diapycnal diffusivity K ρ (normalized by the molecular diffusivity κ) and a large scale dynamical parameter. Note that, in the present inhomogeneous flow, the z-dependency of the parameters should be considered (for instance, a volume averaged Cox number would be meaningless). The presence of a mean shear implies that a Richardson number should be taken for the dynamical parameter. We make use of the monotonicity of the sorted density profile ρ b to introduce a Richardson number that remains one-signed at any time

Ri b (z) =-J dρ b /dz (dU/dz) 2 . ( 26
)
K ρ (z)/κ is plotted as a function of 1/Ri b (z) in Fig. 12 during the final regime of the two-dimensional runs displayed in Fig. 11a, when γ has reached a constant value. A linear dependency is found. It may be expressed by the same simple relation as in the previous sections, though a local formulation should be used to account for the vertical inhomogeneity. Since we deal with a sorted density profile, local means here that a given isopycnal is considered. We get where γ loc (z) = (φ d (z)φ i (z))/ǫ loc (z) is a local mixing efficiency (φ i (z) =-κdρ b /dz), which depends upon z a priori; we have also assumed that the local (in the above sense) dissipation rate of kinetic energy is mostly contributed by the vertical shear of the mean flow: ǫ loc (z) ≃ 2ν(dU/dz) 2 .

K ρ (z) κ = Pr γ loc (z) Ri b (z) + 1 (27)
The linear dependency we find in Fig. 12 implies that γ loc is actually a constant function of z. We found that this constant function has a value equal to 3γ . When the same analysis is conducted in a three-dimensional shear layer, this linear dependency is not recovered, even when Ri b (z) is defined using either a sorted or a horizontally averaged dissipation rate divided by 2ν (instead of (dU/dz) 2 ). Hence, in general, γ loc is a non trivial function of z.

The finding of a constant value in the two-dimensional shear layer but not in the threedimensional one suggests to analysing the main differences between the two-and threedimensional flow dynamics during the final regime. One important difference is that a quasi-steady state has been reached in two dimensions. For instance, the advective flux Φ a nearly matches the fluctuating diapycnal flux Φ d -Φ i during this later stage [START_REF] Staquet | Mixing in a stably stratified shear layer: two-and three-dimensional numerical experiments[END_REF]. The same result holds for the final regime of breaking gravity waves (Fig. 4) as well as of homogeneous stably stratified turbulence (Fig. 7a). The second difference is that the flow has organized into layers during this final regime, which implies that the dissipation rate of kinetic energy is mainly contributed by the vertical variability of horizontal motions. In homogeneous stably stratified turbulence for instance, 60% of ǫ is contributed by ν(∂u/∂z) 2 [START_REF] Godeferd | Statistical modelling and direct numerical simulations of decaying stably-stratified turbulence: Part 2. Large and small scale anisotropy[END_REF] and we found this ratio to be higher for low amplitude breaking gravity waves. This is also the case in the two-dimensional shear layer but does not hold at all during the final stage of the three-dimensional shear layer.

Conclusion

We have studied the mixing properties of three stably stratified flows: breaking gravity waves, a homogeneous stably stratified decaying turbulence and a stably stratified shear layer in two and three dimensions. We have found that, during the final regime of the homogeneous flows (stably stratified turbulence and breaking gravity waves) and of the two-dimensional stably stratified shear layer, the diapycnal diffusivity behaves linearly as a function of a Froude number squared when the flow is homogeneous and as a function of the inverse of a Richardson number in the two-dimensional shear layer case. This behaviour can be accounted for by a simple expression which shows that the coefficient of the linear law is the mixing efficiency of the flow, γ , times the Prandtl number of the fluid. Such a linear law thus amounts to the scaling often assumed by oceanographers (e.g. [START_REF] Munk | Abyseal recipes II: energetics of tidal mixing and wind mixing[END_REF]. The mixing efficiency is volume averaged in the homogeneous case and averaged over an isopycnal in the shear layer case. This linear law implies that γ has reached a constant value.

When analysing the final regime of the three flows sharing this common linear law, we find that a quasi-steady stable state has been reached, in which the production and dissipation rates of available potential energy balance; also, the dissipation rate of kinetic energy is mostly contributed by the vertical variability of horizontal motions. Thus, the flow has become weakly turbulent and organized into layers. Further work remains to be done to understand whether such a linear law is a generic mixing property of this type of flow.

  who provided analytical expressions for the fluxes that control Ėb . These fluxes are the diapycnal flux, denoted Φ d , as well as, if mass flux occurs across the bounding surface of the fluid domain, advective and diffusive mass fluxes across that surface. The latter fluxes are denoted Ėsurf adv b and Ėsurf dif b , respectively. The temporal evolution of E b thus obeys the equation Ėb = Φ d + Ėsurf adv b + Ėsurf dif b .

Fig. 1 .

 1 Fig.1. Contours of the density field (left column) and vorticity field (right column) in the (x, z) vertical plane, for an internal gravity wave of Froude number Fr t = 0.013 and Prandtl number Pr = 1 propagating toward the upper right corner of the domain. Upper row: views at Nt/2π = 32.6, shortly before breaking occurs: the parametric instability that grows through resonant interaction with the primary wave is of large enough amplitude to be visible. It is manifested as vorticity bands of alternate sign, organized into a wave packet (right frame). Lower row: the fields are displayed one Brunt-Väisälä period later, just after breaking has occurred. The vorticity layers have become unstable (the origin of this instability is discussed in[START_REF] Koudella | Two-dimensional instability mechanisms of a progressive internal gravity wave[END_REF]).

Fig. 2 .

 2 Fig. 2. Instantaneous Cox number Cox(t) defined by Eq. (14) as a function of the Froude number squared Fr 2 t (t) defined by Eq. (13) for the runs with Pr = 0.72 (-, bottom curve), Pr = 5(•••), Pr = 10 (-, top curve). All times after breaking are shown (Nt/2π ≥ 35).

Fig. 3 .

 3 Fig. 3. Mixing efficiency defined by Eq. (7) as a function of time (scaled by the Brunt-Väisälä period) for the same runs as in Fig. 2. The Pr = 1 run is also shown for comparison.

Fig. 4 .

 4 Fig. 4. Vertical advective flux of density Φ a and fluctuating diapycnal flux Φ d -Φ i vs. time for the Pr = 0.72 run. Time is scaled by the Brunt-Väisälä period.

Fig. 5 .

 5 Fig. 5. Run with N = π at 6 Brunt-Väisälä periods. Iso-surface of ∂v/∂z of value equal to 25% of the maximum value of ∂v/∂z in the fluid domain at that time (from Godeferd and Staquet, 2001).

Fig. 6 .

 6 Fig. 6. Temporal evolution of the Reynolds number based upon the Taylor microscale, defined by (19) for the two following runs, (-): ν = 1/300, N = π; (---): ν = 1/300, N = π/2. Time is scaled by the Brunt-Väisälä period 2π/N.

Fig. 7 .

 7 Fig. 7. (a) Comparison of the fluctuating diapycnal flux Φ d -Φ i (-) with the advective flux of density Φ a (---) for the N = π run. (b) Temporal evolution of the mixing efficiency γ h , defined in Section 4.2, (-): N = π ; (---): N = π/2. Time is scaled by the Brunt-Väisälä period 2π/N.

Fig. 8 .

 8 Fig. 8. Cox number Cox h (t) as a function of (Fr(t)Re λ (t)) 2 , where Fr denotes a Froude number and Re λ , a Reynolds number based upon the Taylor microscale, for all times of the two runs described in Fig. 6. (---): N = π/2; (-): N = π. All parameters are defined in Section 4.

Fig. 9 .

 9 Fig. 9. Vorticity contours of a shear layer evolving in a vertical plane (J = 0.16, Re = 2000, Pr = 0.7). (a) The Kelvin-Helmholtz instability has developed and saturated. (b) Turbulent regime, once two-dimensional secondary instabilities have developed. (c) Final time of the simulation: a steady regime (apart from diffusive effects) has been reached. The time unit is δ i /U (δ i and U are defined in Section 5.1).

Fig. 10 .

 10 Fig. 10. Temporal evolution of a three-dimensional shear layer simulation with J = 0.167, Re ≃ 300 and Pr = 1.4. A surface equal to 0.5 max D (|ω x |+|ω y |+|ω z |), where |ω i | is the i-component of the vorticity and D the fluid domain, is shown at successive times: (a) t = 34, once the Kelvin-Helmholtz instability has grown; the Kelvin-Helmholtz billow, already distorted by the growth of the three-dimensional perturbation, is visible in between vortex sheets (these are the usual braids, reinforced by the baroclinic torque); the streamwise axis is indicated by an X. (b) t = 108, when three-dimensional small scale secondary instabilities have developed. (c) t = 227, at the end of simulation; only 3/4 of the numerical domain along the vertical direction is shown for clarity. Time unit is defined in Fig. 9 (from Staquet, 2000).

Fig. 11 .

 11 Fig. 11. Mixing efficiency for a stably-stratified shear layer. (a) Two-dimensional shear layer: (-) Re = 2000, Pr = 0.5, (---) Re = 1000, Pr = 1.4, (•••) Re = 1000, Pr = 0.7( J = 0.167 for the three computations); the arrow indicates the time at which the Kelvin-Helmholtz instability saturates. (b) Three-dimensional shear layer; the mixing efficiency has been averaged over 7 computations with 0.04 ≤ J ≤ 0.167, Pr = 0.7 or 1.4 and Re ≃ 300. Time unit is defined in Fig. 9.

Fig. 12 .

 12 Fig. 12. Diapycnal diffusivity K ρ (z) defined by Eq. (6) as a function of the inverse of the Richardson number Ri b (z) defined by Eq. (26). The different curves correspond to successive times of the final regime of a two-dimensional simulation with J = 0.167, Re = 1000, Pr = 1.4.

Acknowledgements

Computations were performed at the Computer Center of CNRS (IDRIS).


Pierre-Alain Reynier
email: pierre-alain.reynier@lif.univ-mrs.fr

Arnaud Sangnier
email: sangnier@di.unito.it

Arnaud Sangnier Weak

Time Petri

Weak Time Petri Nets strike back!

come

Introduction

For verification purpose, e.g. in the development of embedded platforms, there is an obvious need for considering time features and the study of timed models has thus become increasingly important. For distributed systems, different timed extensions of Petri nets have been proposed which allow the combination of an unbounded discrete structure with dense-time variables.

There are several ways to express urgency in timed systems, as discussed in [START_REF] Sifakis | Compositional specification of timed systems[END_REF]. In timed extensions of Petri nets, two types of semantics are considered for time elapsing. In the weak semantics, all time delays are allowed whereas in the strong one, all transitions are urgent, i.e. time delays cannot disable transitions. While for models with finite discrete structure (such as timed extensions of bounded Petri nets or timed automata [START_REF] Alur | A theory of timed automata[END_REF]), standard verification problems are decidable for both semantics, for models with infinite discrete structure, the choice of the semantics has a deep influence on decidability issues. In this work, we consider the model of Time Petri Nets [START_REF] Merlin | Study of the Recoverability of Computing Systems[END_REF] (TPN) where clocks are associated with transitions, and which is commonly considered under a strong semantics. In this model, all the standard verification problems are known to be undecidable [START_REF] Jones | Complexity of some problems in Petri nets[END_REF]. On the other hand, in the model of timed-arc Petri nets [START_REF] Bolognesi | From timed Petri nets to timed LOTOS[END_REF], where clocks are associated with tokens and which is equipped with a weak semantics, many verification problems are decidable (coverability, boundedness...). Indeed, this semantics entails for this model monotonicity properties which allow the application of well-quasi-ordering techniques, see [START_REF] De Frutos-Escrig | Decidability of properties of timedarc Petri nets[END_REF][START_REF] Abdulla | Timed Petri nets and bqos[END_REF][START_REF] Abdulla | Dense-timed Petri nets: Checking zenoness, token liveness and boundedness[END_REF]. Note however that the reachability of a discrete marking is undecidable, as proven in [START_REF] Valero Ruiz | On non-decidability of reachability for timed-arc Petri nets[END_REF]. A natural question, which had surprisingly no answer until now, as mentioned in a recent survey on the topic [START_REF] Boyer | On the compared expressiveness of arc, place and transition time Petri nets[END_REF], is thus to study TPN under a weak semantics of time elapsing.

The time-elapsing policy states which delays are allowed in a configuration. The memory policy is concerned with the resets of clocks, and intuitively specifies, when firing a transition, which timing informations are preserved. The original model of Merlin [START_REF] Merlin | Study of the Recoverability of Computing Systems[END_REF] is equipped with an intermediate semantics which considers the intermediary marking bewteen consumption and production. Two others memory policies have been considered in [START_REF] Bérard | Comparison of different semantics for time Petri nets[END_REF] (the atomic and the persistent atomic) in which the firings of transitions are performed atomically. While these policies can be thought as cosmetic for the model of TPN, the results we obtain show this is not the case.

We are interested in the impact of the weak semantics on TPN, distinguishing between the different memory policies. We first study the decidability issues, and prove that for TPN with weak intermediate semantics, a discrete marking is reachable if and only if it is reachable in the underlying untimed Petri net. As a corollary, the problem of the marking reachability (and also coverability, boundedness) is decidable for this model. More surprisingly, we also prove that when changing the memory policy this result does not hold anymore and the verification problems become undecidable. In this work, we only consider untimed verification problems and we plan to study timed versions in future work. We then compare w.r.t. weak time bisimilarity (weak stands here for silent transitions) the expressive power of weak TPN looking at the different memory policies. We first prove that the persistent atomic semantics is strictly more expressive that the atomic semantics. Then, concerning atomic and intermediate memory policies, we provide a TPN which shows that the atomic semantics is not included in the intermediate one. Related works. As mentioned above, there are, up to our knowledge, only very few works considering TPN under a weak semantics. In [START_REF] Boyer | On the compared expressiveness of arc, place and transition time Petri nets[END_REF] the authors have proven that the weak intermediate semantics and the strong intermediate semantics are uncomparable. In another line of work, [START_REF] Haar | Equivalence of timed state machines and safe time Petri nets[END_REF] considers TPN under a semantics which is a kind of compromise between the standard strong and weak semantics. They provide translations between this model and timed state machines.

Definitions

Let Σ be a finite alphabet, Σ * is the set of finite words over Σ. We note Σ τ = Σ ∪ {τ } where τ / ∈ Σ represents internal actions. ǫ will represent the empty word. The sets N, Q, Q ≥0 , R and R ≥0 are respectively the sets of natural, rational, non-negative rational, real and non-negative real numbers. A valuation v over a finite set X is a mapping in R X ≥0 . For v ∈ R X ≥0 and d ∈ R ≥0 , v+d denotes the valuation defined by (v+d)(x) = v(x)+d. We note 0 the valuation which assigns to every x ∈ X the value 0.

As commonly in use for Time Petri Nets, we will associate rational intervals with transitions. Note that we could handle intervals with bounds given as real numbers if we abstract the problem of comparison of real numbers. We consider the set I(Q ≥0 ) of non-empty intervals (a, b) with non-negative rational bounds a, b ∈ Q ≥0 . We consider both open and closed bounds, and also allow a right open infinite bound as in [2, +∞[.

Petri Nets Definition 1 (Labeled Petri Net (PN)).

A Labeled Petri Net over the alphabet Σ τ is a tuple (P, T, Σ τ , • (.), (.)

• , M 0 , Λ) where:

-P is a finite set of places, -T is a finite set of transitions with P ∩ T = ∅, -• (.) ∈ (N P ) T is the backward incidence mapping, -(.)

• ∈ (N P ) T is the forward incidence mapping,

-M 0 ∈ N P is the initial marking, -Λ : T → Σ τ is the labeling function
As commonly in use in the literature, the vector

• (t) (resp. (t) • ) in N P is noted • t (resp. t • ).
The semantics of a PN N = (P, T, Σ τ , • (.), (.)

• , M 0 , Λ) is given by its associated labeled transition system N = (N P , M 0 , Σ τ , ⇒) where ⇒⊆ N P × Σ τ × N P is the transition relation defined by

M a =⇒ M ′ iff ∃t ∈ T s.t. Λ(t) = a ∧ M ≥ • t ∧ M ′ = M -• t + t • .
For convenience we will sometimes also write, for t ∈ T ,

M t =⇒ M ′ if M ≥ • t and M ′ = M -• t + t • . We also write M ⇒ M ′ if there exists a ∈ Σ τ such that M a =⇒ M ′ .
The relation ⇒ * represents the reflexive and transitive closure of ⇒. We denote by Reach(N ) the set of reachable markings defined by {M ∈

N P | M 0 ⇒ * M }.
It is well known that for PN the reachability problem which consists in determining whether a given marking M belongs to Reach(N ) is decidable; it has in fact been proved independently in [START_REF] Mayr | An algorithm for the general Petri net reachability problem[END_REF] and [START_REF] Kosaraju | Decidability of reachability in vector addition systems (preliminary version)[END_REF].

We introduce a last notation concerning Labeled Petri Nets. Given a PN N , a marking M of N and a multi-set ∆ = t 1 , . . . , t n of transitions of N , we write M ∆ |=⇒ M ′ if and only if the multi-set ∆ can be fired from M , meaning that there exists an ordering of transitions in ∆, represented as a permutation ϕ of {1, . . . , n}, such that the sequence of firings M

t ϕ(1) ==⇒ M 1 t ϕ(2) ==⇒ M 2 . . . t ϕ(n) ==⇒ M ′ exists in N .

Timed Transition Systems

Timed transition systems describe systems which combine discrete and continuous evolutions. They are used to define the behavior of timed systems such as Time Petri Nets [START_REF] Merlin | Study of the Recoverability of Computing Systems[END_REF] or Timed Automata [START_REF] Alur | A theory of timed automata[END_REF].

Definition 2 (Timed Transition System (TTS))

. A timed transition system over the alphabet Σ τ is a transition system S = (Q, q 0 , Σ τ , →), where the transition relation

→⊆ Q × (Σ τ ∪ R ≥0 ) × Q consists of discrete transitions q a -→ q ′ (with a ∈ Σ τ )
representing an instantaneous action, and continuous transitions q d -→ q ′ (with d ∈ R ≥0 ) representing the passage of d units of time.

Moreover, we require the following standard properties for TTS :

-TIME-DETERMINISM : if q d -→ q ′ and q d -→ q ′′ with d ∈ R ≥0 , then q ′ = q ′′ , -0-DELAY : q 0 -→ q, -ADDITIVITY : if q d -→ q ′ and q ′ d ′ -→ q ′′ with d, d ′ ∈ R ≥0 , then q d+d ′ ---→ q ′′ , -CONTINUITY : if q d -→ q ′ , then for every d ′ and d ′′ in R ≥0 such that d = d ′ + d ′′ ,
there exists q ′′ such that q

d ′ -→ q ′′ d ′′ -→ q ′ .
With these properties, a run of S can be defined as a finite sequence of moves

ρ = q 0 d0 -→ q ′ 0 a0 -→ q 1 d1 -→ q ′ 1 a1
-→ q 2 . . . an --→ q n+1 where discrete and continuous transitions alternate. To such a run corresponds the timed word w = (a i , η i ) 0≤i≤n over Σ τ where η i = i j=0 d j is the time at which a i happens. We then denote by Untimed(w) the projection of the word a 0 a 1 . . . a n over the alphabet Σ and by Duration(w) the duration η n . Note that in the word Untimed(w) the symbol τ does not appear. We will sometimes apply, without possible ambiguities, these notations to runs writing Untimed(ρ) and Duration(ρ). We might also describe the run writing directly q 0 w -→ q n+1 .

Time Petri Nets

Syntax. Introduced in [START_REF] Merlin | Study of the Recoverability of Computing Systems[END_REF], Time Petri Nets associate a time interval with each transition of a Petri net.

Definition 3 (Labeled Time Petri Net (TPN)).

A Labeled Time Petri Net over the alphabet Σ τ is a tuple (P, T, Σ τ , • (.), (.)

• , M 0 , Λ, I) where:

-(P, T, Σ τ , • (.), (.) • , M 0 , Λ) is a PN, -I : T → I(Q ≥0 )
associates with each transition a firing interval.

In the sequel, we associate with an interval its left bound and its right bound. More generally, given a transition t of a TPN, we will denote by eft(t) (resp. lft(t)) the left bound of I(t) (resp. the right bound of I(t)), standing for earliest firing time (resp. latest firing time). We have hence I(t) = (eft(t), lft(t)).

Semantics.

A configuration of a TPN is a pair (M, ν), where M is a marking over P , i.e. a mapping in N P , with M (p) the number of tokens in place p. A transition t is enabled in a marking M if M ≥ • t. We denote by En(M ) the set of enabled transitions in M . The second component of the pair (M, ν) is a valuation over En(M ), i.e. a mapping in R En(M ) ≥0

. Intuitively, for each enabled transition t in M , ν(t) represents the amount of time that has elapsed since t is enabled. An enabled transition t can be fired if ν(t) belongs to the interval I(t). The marking obtained after this firing is as usual the new marking M ′ = M -• t + t • . Moreover, some valuations are reset and we say that the corresponding transitions are newly enabled.

Different semantics can be chosen in order to realize these resets. This choice depends of what is called the memory policy. For M ∈ N P and t, t ′ ∈ T such that t ∈ En(M ) we define in different matters a predicate ↑ enabled s (t ′ , M, t) with s ∈ {I, A, P A} which is true if t ′ is newly enabled by the firing of transition t from marking M , and false otherwise. This predicate indicates whether we need to reset the clock of t ′ after firing the transition t at the marking M .

I:

The intermediate semantics considers that the firing of a transition is performed in two steps: consuming the tokens in • t, and then producing the tokens in t • . Intuitively, it resets the clocks of t and of the transitions that could not be fired in parallel with t from the marking M . Formally, the predicate ↑enabled I (t ′ , M, t) is defined by:

↑enabled I (t ′ , M, t) = t ′ ∈ En(M -• t + t • ) ∧ (t ′ / ∈ En(M -• t) ∨ t = t ′ ) A:
The atomic semantics considers that the firing of a transition is obtained by an atomic step. It resets the clocks of t and of the transitions t ′ which are not enabled at M . The corresponding predicate ↑enabled A (t ′ , M, t) is defined by:

↑enabled A (t ′ , M, t) = t ′ ∈ En(M -• t + t • ) ∧ (t ′ / ∈ En(M ) ∨ t = t ′ ) PA:
The persistent atomic semantics behaves as the atomic semantics except that it does not reset the clock of t.

↑enabled P A (t ′ , M, t) = t ′ ∈ En(M -• t + t • ) ∧ t ′ / ∈ En(M )
Finally, as recalled in the introduction, there are two ways of letting the time elapse in TPN. The first way, known as the strong semantics, is defined in such a matter that time elapsing cannot disable a transition. Hence, when the upper bound of a firing interval is reached then the transition must be fired or disabled. In contrast to that the weak semantics does not make any restriction on the elapsing of time. In this work, we focus on the weak semantics of TPN.

Definition 4 (Weak semantics of a TPN).

Let s ∈ {I, A, P A}. The weak s-semantics of a TPN N = (P, T, Σ τ , • (.), (.)

• , M 0 , Λ, I) is a timed transition system N s = (Q, q 0 , Σ τ , → s ) where Q = N P × R En(M ) ≥0
, q 0 = (M 0 , 0) and → s consists of discrete and continuous moves:

-the discrete transition relation is defined ∀a ∈ Σ τ by:

(M, ν) a -→ s (M ′ , ν ′ ) iff ∃t ∈ T s.t.                Λ(t) = a, and, t ∈ En(M ) ∧ M ′ = M -• t + t • , and, ν(t) ∈ I(t), and, ∀t ′ ∈ En(M ′ ), ν ′ (t ′ ) = 0 if ↑enabled s (t ′ , M, t) ν(t ′ ) otherwise -the continuous transition relation is defined ∀d ∈ R ≥0 by: (M, ν) d -→ s (M, ν ′ ) iff ν ′ = ν + d We also write a discrete transition (M, ν) t -→ s (M ′ , ν ′ ) to characterize the transi- tion t ∈ T which allows the firing (M, ν) Λ(t) ---→ s (M ′ , ν ′ ).
We extend this notation to words θ ∈ (T ∪ R ≥0 ) * , which correspond to sequences of transitions and delays and lead to a unique (if it exists) run ρ. We may write this run ρ : (M, ν) θ -→ s (M ′ , ν ′ ) and use Untimed(θ) (resp. Duration(θ)) to denote the word Untimed(ρ) (resp. to represent the delay Duration(ρ)). Finally, for s ∈ {I, A, P A}, we write

(M, ν) → s (M ′ , ν ′ ) if there exists a ∈ Σ τ ∪ R ≥0 such that (M, ν) a -→ s (M ′ , ν ′ ). The relation → *
s denotes the reflexive and transitive closure of → s . For a TPN N with an initial marking M 0 we define the following reachability sets according to the considered semantics:

Reach(N ) s = {(M, v) | (M 0 , 0) → * s (M, v)}.
Example 1. We illustrate the impact of the three memory policies in weak semantics.

Consider the net depicted on Figure 1, and the execution

t2, c, [0, 1] t1, a, [0, +∞[ • p Fig. 1. The TPN N 1 . (M, 0) 1 -→ s (M, 1) a -→ s (M, v) where M (p) = 1.
With the intermediate semantics, both clocks are reset as in the intermediate marking, the place p is empty. With the atomic semantics, the clock associated with t2 is not reset and the clock associated with t1 is reset because it corresponds to the fired transition. Finally, with the persistent atomic semantics no clock is reset.

Decidability

Considered problems and known results

Assume N = (P, T, Σ τ , • (.), (.)

• , M 0 , Λ, I) is a TPN. In this section, we will consider the following problems for s ∈ {I, A, P A}:

(1) The marking reachability problem : given M ∈ N P , does there exist ν ∈ R

En(M ) ≥0 such that (M, ν) ∈ Reach(N ) s ? (2)
The marking coverability problem : given M ∈ N P , does there exist M ′ ∈ N P and ν ∈ R It is well known that the "untimed" versions of these problems are decidable in the case of Petri nets. In fact, as mentioned before the marking reachability problem is decidable for Petri nets [START_REF] Kosaraju | Decidability of reachability in vector addition systems (preliminary version)[END_REF][START_REF] Mayr | An algorithm for the general Petri net reachability problem[END_REF] and the two other problems can be solved using the Karp and Miller tree whose construction is given in [11].

En(M ′ ) ≥0 such that M ′ ≥ M and (M ′ , ν) ∈ Reach(N ) s ? (
From [START_REF] Jones | Complexity of some problems in Petri nets[END_REF], we know that these problems are all undecidable when considering TPN with strong semantics no matter whether the semantics is intermediate, atomic or persistent atomic. In fact a TPN with strong semantics can simulate a Minsky machine. A Minsky machine manipulates two integer variables c 1 and c 2 and is composed of a finite number of instructions, each of these instructions being either an incrementation (q : c i := c i + 1) or a decrementation with a test to zero (q : if c i = 0 goto q ′ else c i := c i -1; goto q ′′ ), where i ∈ {1, 2} and q, q ′ , q ′′ are some labels preceding each instruction. There is also a special label q f from which the machine cannot do anything.

In [START_REF] Minsky | Computation: finite and infinite machines[END_REF], Minsky proved that the halting problem, which consists in determining if the instruction labeled with q f is reachable, is undecidable.

It is easy to encode an incrementation using a TPN (or even a PN), with a transition consuming a token in a place characterizing the current control state and producing a token in the next control state and in a place representing the incremented counter.

When encoding the decrementation with the test to zero, the strong semantics plays a crucial role. This encoding is represented on Figure 2. If there is a token in the place c i , there is no way for the TPN to produce a token in the place q ′ because time cannot elapse since the transition labeled with the interval [0, 0] is firable. The example of the Figure 2 shows that the strong time semantics allows to encode priorities (between transitions in conflict) and thus to encode inhibitor arcs. This construction obviously fails with the weak semantics.

• q q ′′ c i q ′ [0, 0] [1, 1]

The peculiar case of TPN with weak intermediate semantics

We prove here that the undecidability results we had before in the case of TPN with strong semantics do not hold anymore when considering the weak intermediate semantics. Before proving this we introduce some notations. For a TPN N = (P, T, Σ τ , • (.), (.)

• , M 0 , Λ, I), we denote by N U the untimed PN obtained by removing from N the component I. Furthermore given a set of configurations C ⊆ N P × R T ≥0 of N , we denote by Untime(C) the projection of C over the set N P . For s ∈ {I, A, P A}, we have by definition of the different semantics that Untime(Reach(N ) s ) ⊆ Reach(N U ) and as shown by the example given in Figure 2 this inclusion is strict in the case of the strong semantics. When considering the weak intermediate semantics, we prove that from any sequence of transitions ∆ firable in N U , we can effectively compute a reordering of ∆, and associated timestamps, leading to a correct run of N I . Theorem 5. For all TPN N , Untime(Reach(N ) I ) = Reach(N U ).

According to the previous remark, we only have to prove that Reach(N U ) ⊆ Untime(Reach(N ) I ). Therefore, we first state the following property expressing that if we reduce the intervals associated with transitions, this restricts the set of reachable configurations: Lemma 6. Let N , N ′ be two TPN identical except on their last component associating intervals to transitions, say respectively I and I ′ . If we have I ′ (t) ⊆ I(t) for any t ∈ T , then Reach(N ′ ) I ⊆ Reach(N ) I .

In the sequel, we will consider TPN in which intervals are reduced to singletons. That is we have I(t) = [eft(t), lft(t)] with eft(t) = lft(t) for all transitions t ∈ T . The proof of the result in this particular case thus entails the result in the general case. Before to proceed we introduce additional definitions for TPN.

Given a TPN N , a marking M of N and ∆ a multiset of transitions of N , we

define the set Candidate(M, ∆) = {t ∈ ∆ | M t ⇒ M ′ ∆\t |=⇒}.
We will then say that a configuration (M, ν) is compatible with a multiset ∆ iff:

M ∆ |=⇒ and ∀t ∈ Candidate(M, ∆), ν(t) ≤ lft(t).
We now prove the following proposition, which intuitively states how to turn a sequence of transitions in the untimed Petri net into a timed execution of the TPN. Proposition 7. Let N be a TPN with singleton intervals and (M, ν) be a configuration of N compatible with some multiset of transitions ∆. Then, for any transition t ∈ Candidate(M, ∆) such that δ(t) = lft(t)ν(t) is minimal (among the transitions of Candidate(M, ∆)), we have:

(i) (M, ν) δ(t) --→ I (M, ν + δ(t)) t -→ I (M ′ , ν ′ ), (ii) (M ′ , ν ′ ) is compatible with ∆ ′ = ∆ \ t, Proof. Let t ∈ Candidate(M, ∆) be such that for all t ′ ∈ Candidate(M, ∆), we have lft(t) -ν(t) = δ(t) ≤ δ(t ′ ) = lft(t ′ ) -ν(t ′ ). (i) First the time elpasing transition (M, ν) δ(t) --→ I (M, ν + δ(t)
) is possible as we consider the weak semantics. Second, the discrete transition (ii) To prove compatibility, first note that

(M, ν + δ(t)) t -→ I (M ′ , ν ′ ) is also possible as ν(t) + δ(t) = lft(t)
M ′ ∆ ′ |=⇒ because t ∈ Candidate(M, ∆). Second, let t ′ ∈ Candidate(M ′ , ∆ ′ ).
We distinguish two cases according to the value of the predicate ↑enabled I (t, M, t ′ ):

-If ↑enabled I (t, M, t ′ ) is true, then we have ν ′ (t ′ ) = 0 the result follows.

-Otherwise, the definition of ↑enabled

I (t, M, t ′ ) implies that M ≥ • t + • t ′ .
As a consequence, we have M

t ′ ⇒ t ⇒.
Then as t ′ ∈ Candidate(M ′ , ∆ \ t) we get that t ′ ∈ Candidate(M, ∆). Due to the minimality of δ(t) among the set Candidate(M, ∆), we obtain ν

′ (t ′ ) = ν(t ′ ) + δ(t) ≤ ν(t ′ ) + δ(t ′ ) = lft(t ′ ) as desired.
This concludes the proof.

⊓ ⊔

The inclusion Reach(N U ) ⊆ Untime(Reach(N ) I ) in the case of TPN with singleton intervals easily follows from this result. Indeed, consider some reachable marking M in Reach(N U ). There exists a sequence of transitions that leads to M from M 0 , we represent it through some multiset ∆. As initially all clock valuations are null in N I , the configuration (M 0 , 0) is thus compatible with ∆. An induction on the size of ∆, together with Proposition 7, thus gives the result. Note that Proposition 7 describes an effective procedure to compute a timed execution of N I : simply consider the transitions that are candidates, and choose one with the earliest deadline.

Using the decidability results in the case of PN, we obtain the following corollary:

Corollary 8. The marking reachability, marking coverability and boundedness problems are decidable in the case of TPN with weak intermediate semantics.

Undecidability for weak atomic and weak persistent atomic semantics

We consider now the case of the weak atomic and weak persistent atomic semantics. As for the strong semantics, but with a more involved construction, we will show that it is possible to encode the behavior of a Minsky machine into a TPN with weak (persistent) atomic semantics from which we will deduce the undecidability results. The TPN we build contains a place for each counter c i with i ∈ {1, 2} and a place for each label q of the considered Minsky machine. Furthermore, when executing the net, we will preserve the invariant there is a single place corresponding to a label q which is marked.

• q

[0, 0] t1 p1 [1, 1] t2 p2 ci [0, 0] t3 [0, 0] t5 2 p4 [0, 0] t4 q ′ p3 q ′′
Fig. 3. Encoding decrementation with weak atomic or persistent atomic semantics.

Encoding an incrementation can be done as in the strong semantics. Figure 3 shows how to encode the instruction (q : if c i = 0 goto q ′ else c i := c i -1; goto q ′′ ) using a TPN with weak atomic or persistent atomic semantics. We now explain the idea of this encoding. We consider the two following cases for the net shown in Figure 3: 1. Assume that the only place which contains a token is the place q, which means we are in the case where the value of c i is equal to 0 (no token in place c i ). The net then can only fire the sequence of transitions t1, t2, t3 and then t4 and finally it reaches a configuration where the only marked place is q ′ . 2. Assume now that there is a token in place q and that there is at least one token in place c i . We are in the case where the value of c i is different of 0. We have the following sequence of transitions:

-only the transition t1 is firable, so the net fires it; -afterwards the transition t2 and the transition t3 are firable. In fact, since we are considering weak semantics the deadline of t3 can be ignored thus making time passage in order to fire t2. Note that if the net chooses to fire t3, it will reach a deadlock state where no more transitions can be fired without having put a token in the place q ′ or q ′′ , therefore we assume that the transition t2 is first fired; -after having waiting one time unit and firing t2, the only transition which can be fired is t5. In fact since we are considering atomic (or persistent atomic) semantics, firing t2 does not make t3 newly enabled, whereas the weak intermediate semantics would have reset the clock associated to t3. So the net fires t5 consuming the token in p2, p3 and two tokens in c i (at least one was present from the initial configuration and the first firing of t2 added another one); -finally the net ends in a configuration with one token in q ′′ and the place c i contains one token less than in the initial configuration.

The above construction allows to reduce the halting problem for Minsky machine to the marking coverability problem for weak (persistent) atomic semantics. From this we can also deduce the undecidability for the marking reachability and boundedness problems. Hence: Theorem 9. The marking reachability, marking coverability and boundedness problems are undecidable for TPN with weak atomic or weak persistent atomic semantics.

In comparison with what occurs in the case of the strong semantics, this result is surprising, and it reveals the important role played by the memory policy when considering the weak semantics. Recall that as we have seen earlier, with the strong semantics, these problems are undecidable no matter which memory policy is chosen.

Finally, in the above construction, we can replace the edges between p2 and t2 by a read arc. Consequently, the considered problems are also undecidable for weak intermediate TPN with read arcs, unlike what happens for timed-arc Petri nets [START_REF] Bouyer | Timed Petri nets and timed automata: On the power of zeno sequences[END_REF].

Expressiveness

Preliminaries

Let S = (Q, q 0 , Σ τ , →) be a TTS. We define the relation

֒→⊆ Q × (Σ ∪ R ≥0 ) × Q by: -for d ∈ R ≥0 , q d ֒-→ q ′ iff
there is a run ρ in S such that ρ = q w -→ q ′ and U ntimed(w) = ε and Duration(w) = d, -for a ∈ Σ, q a ֒-→ q ′ iff there is a run ρ in S such that ρ = q w -→ q ′ and U ntimed(w) = a and Duration(w) = 0.

This allows us to define the following notion:

Definition 10 (Weak Timed Bisimulation). Let S 1 = (Q 1 , q 1 0 , Σ τ , → 1 ) and S 2 = (Q 2 , q 2 0 , Σ τ , → 2 ) be two TTS and ∼ be a binary relation over

Q 1 × Q 2 .

∼ is a weak timed bisimulation between S 1 and S 2 if and only if:

q 1 0 ∼ q 2 0 , and,

-for a ∈ Σ ∪ R ≥0 , if q 1 a ֒-→ 1 q ′ 1 and if q 1 ∼ q 2 then there exists q ′ 2 ∈ Q 2 such that q 2 a ֒-→ 2 q ′ 2 and q ′ 1 ∼ q ′ 2 ; conversely if q 2 a ֒-→ 2 q ′ 2 and if q 1 ∼ q 2 then there exists q ′ 1 ∈ Q 1 such that q 1 a ֒-→ 1 q ′ 1 and q ′ 1 ∼ q ′ 2 .
Two TTS S 1 and S 2 are weak timed bisimilar if there exists a weak timed bisimulation between S 1 and S 2 . We then write S 1 ≈ S 2 .

Definition 11 (Expressiveness w.r.t. Weak Timed Bisimilarity). The class C of TTS is less expressive than C ′ w.r.t. weak timed bisimilarity if for all TTS S ∈ C

there is a TTS S ′ ∈ C ′ such that S ≈ S ′ . We write C ⊑ C ′ . If moreover there is a S ′ ∈ C ′ such that there is no S ∈ C with S ≈ S ′ , then C is strictly less expressive than C ′ . We then write C ⊏ C ′ .

For s ∈ {I, A, P A}, we will denote by T PN s the class of TTS induced by TPN with s-semantics.

Atomic versus Persistent Atomic semantics

In [START_REF] Bérard | Comparison of different semantics for time Petri nets[END_REF], the authors prove that for TPN with strong semantics, the persistent atomic semantics is more expressive than the atomic semantics. We prove here that this result still holds in the case of the weak semantics. Intuitively, as it is shown on Figure 4, from a TPN with atomic semantics, we build another TPN in which we duplicate each transition. During an execution of this last TPN, at most one of the transitions obtained after duplication is enabled, and when it is fired it cannot be enabled again at the next step whereas the other one can. This trick allows us to simulate the atomic semantics with the persistent atomic one.

t Λ(t) I(t)

input places output places Proof. Let N = (P, T, Σ τ , • (.), (.)

t 1 Λ(t) I(t) t 2 Λ(t) I(t) • p 1 t p 2
• , M 0 , Λ, I) be a TPN over Σ τ . Figure 4 represents the construction of the TPN N ′ . Formally, its set of places P ′ is equal to P ∪ {p 1 t , p 2 t | t ∈ T } and its set of transitions T ′ contains two copies t 1 and t 2 of each transition t ∈ T . These copies are connected in the same way as the transition t is in N , plus additional edges to the places p 1 t and p 2 t , as depicted on Figure 4. Finally the initial marking of N ′ is M ′ 0 such that for all p ∈ P , M ′ 0 (p) = M 0 (p) and for all t ∈ T ,

M ′ 0 (p 1 t ) = 1 and M ′ 0 (p 2 t ) = 0. We now consider the relation ∼⊆ (N P × R T ≥0 ) × (N P ′ × R T ′ ≥0
) between the configurations of N A and the ones of

N ′ P A defined by (M, ν) ∼ (M ′ , ν ′ ) iff: -for all p ∈ P , M (p) = M ′ (p) and for all t ∈ T , M ′ (p 1 t ) + M ′ (p 2 t ) = 1, -for all t ∈ T , for all i ∈ {1, 2} if t ∈ En(M ) and t i ∈ En(M ′ ) then ν(t) = ν ′ (t i ).
It is then easy to verify that the relation ∼ is a weak timed bisimulation.

⊓ ⊔

We will now prove that the inclusion we obtain in the above proposition is strict. But before, we address a technical point which we will use to delay some sequences of transitions in weak TPN.

Lemma 13. Let s ∈ {I, A, P A} and consider a TPN N such that b is the smallest positive upper bound of the intervals of N . Let ρ be a run in N s of the form ρ :

(M, ν) δ>0 --→ s (M, ν + δ) t1 -→ s • • • tn -→ s ,
such that there exists a value η ≥ 0 verifying:

(i) ∀i ∈ {1, . . . , n}, t i ∈ En(M ) ⇒ ν(t i ) ≤ η, (ii) η + δ < b 2
Then the sequence ρ ′ : (M, ν)

δ+ b 2 ---→ s (M, ν ′ ) t1 -→ s • • • tn -→ s is firable in N s .
Proof. We introduce the following notations describing the run ρ:

(M, ν) δ>0 --→ s (M, ν + δ) = (M 1 , ν 1 ) t1 -→ s (M 2 , ν 2 ) t2 -→ s (M 3 , ν 3 ) • • • (M n , ν n ) tn -→ s (M n+1 , ν n+1 ).
We split the set of indices {1, . . . , n} into two disjoint subsets, I =0 = {i ∈ {1, . . . , n} | ν i (t i ) = 0}, and I =0 = {i ∈ {1, . . . , n} | ν i (t i ) = 0}. These two sets contain respectively the index of the transitions fired with a null (respectively non null) valuation.

We first show the following property:

∀i ∈ I =0 , 0 < ν i (t i ) < b 2 and [ν i (t i ), b[⊆ I(t i ) (1)
Consider the first part of this property. Note that these valuations are strictly positive by definition of I =0 , thus proving the left inequality. Consider now the right inequality, and let i ∈ I =0 . Since the sequence t 1 . . . t i is instantaneous and follows a non null delay step, t i is never newly enabled during t 1 . . . t i-1 (otherwise ν i (t i ) = 0). As a consequence, we have ν i (t i ) = ν(t i ) + δ, and by properties (i) and (ii) of ρ, we obtain ν i (t i ) ≤ η + δ < b 2 , as desired. By definition of b, and since t i is firable from (M i , ν i ), the inclusion [ν i (t i ), b[⊆ I(t i ) holds for any i ∈ I =0 .

We consider now the run ρ ′ : (M, ν)

δ+ b 2 ---→ s (M, ν + δ + b 2 ) = (M 1 , ν ′ 1 ) t1 -→ s (M 2 , ν ′ 2 ) t2 -→ s (M 3 , ν ′ 3 ) • • • (M n , ν ′ n ) tn -→ s (M n+1 , ν ′ n+1 ).
Note that the increasing of the delay is possible because of the weak semantics: in the strong one, the modification of the delay step could be impossible. To prove that this sequence is firable in N s , we proceed by contradiction. Assume there exists a position i ∈ {1, . . . , n} such that t i is not firable from (M i , ν ′ i ), and pick the smallest position verifying this property. We distinguish two cases:

1. If i ∈ I =0 . Then we have ν i (t i ) = 0. Since the instantaneous sequence t 1 . . . t i immediately follows in ρ the delay step δ > 0, t i is newly enabled by the firing of t 1 . . . t i-1 . Since this property only depends on discrete markings, which are preserved in ρ ′ , t i is also newly enabled by the firing of t 1 . . . t i-1 in ρ ′ . As a consequence, we have ν ′ i (t i ) = 0 = ν i (t i ), thus proving that t i is firable in ρ ′ , yielding a contradiction.

2. If i ∈ I =0 . Then we have ν i (t i ) = 0. As already mentioned above, we have in this case that t i is never newly enabled during t 1 . . . t i-1 in ρ. Since the discrete markings are preserved, this conclusion holds also in ρ ′ . As a consequence, we have

ν i (t i ) = ν(t i ) + δ and ν ′ i (t i ) = ν(t i ) + δ + b 2 , what yields ν ′ i (t i ) = ν(t i ) + b 2 .
Using Property (1), we obtain ν ′ i (t i ) ∈ I(t i ), and then t i is firable from (M i , ν ′ i ), contradicting the assumption.

Finally, we have obtained a contradiction in both cases, thus proving that ρ ′ is firable in N s .

⊓ ⊔

We now consider the TPN N 2 represented on Figure 5. Equipped a, [0, 1]

Fig. 5. The TPN N 2 .
with persistent atomic semantics, it accepts the set of timed words composed of letters a occurring before time 1. We will prove that this timed language cannot be accepted by any TPN equipped with the weak atomic semantics.

Proposition 14. There exists no TPN N (even unbounded) s.t. N A ≈ N 2 P A .

Proof. Assume there exists a TPN N such that N A ≈ N 2 P A . Denote by N the number of transitions of N , by b the smallest positive upper bound of the intervals of N , and consider a timed word w = (a, η 1 )(a, η 2 ) . . .

(a, η k ) such that ∀i, 1 -b 2 < η i < η i+1 < 1, and k ≥ N + 1.
This timed word w is recognized by N 1 P A and there exists thus a run of N A along w. We denote it by ρ and decompose it as follows :

ρ : θ0 -→ A d1 -→ A θ1 -→ A t 1 a -→ A θ ′ 1 -→ A • • • di -→ A θi -→ A t i a -→ A θ ′ i -→ A • • • d k -→ A θ k -→ A t k a -→ A
To obtain this decomposition we proceed as follows. We denote by t i a the i-th transition labelled by a. Then for each position i, we isolate the last delay step occuring before the transition t i a (it exists since η i > η i-1 ) and denote it by d i . Then we gather all the internal transitions occuring between this delay step and the transition t i a , and denote this sequence by θ i . The transitions between t i-1 a and the delay step constitute the sequence θ ′ i-1 . In particular, the following properties hold for any position i : andt i a occurs at time η i .

Λ(t i a ) = a, Untimed(θ i ) = Untimed(θ ′ i ) = ε, d i > 0, Duration(θ i ) = 0,
We claim there exists an index i ∈ {1, . . . , k} such that each transition t appearing in θ i t i a has already been fired since θ 0 , i.e. t also appears in

θ 1 t 1 a θ ′ 1 • • • θ i-1 t i-1 a θ ′ i-1 .
By contradiction, if it is not the case, then we can find, for each index i ∈ {1, . . . , k}, a transition, denoted t i , that never appears before. The choice of k verifying k ≥ N + 1 then implies that there exist two positions i = j such that t i = t j , thus yielding a contradiction. We can now fix an index i verifying the above described property.

We now show that Lemma 13 can be applied to the part of ρ associated with the sequence d i θ i t i a . More precisely, (M, ν) is the configuration reached after firing

θ 0 • • • t i-1 a θ ′ i-1 , the delay δ is equal to d i , the sequence t 1 • • • t n corresponds to θ i t i a , and η is defined as (η i -d i ) -(1 -b
2 ). In the atomic semantics, when a transition is fired, its clock is reset if it is still enabled. This property allows, together with timing constraints on the word w, to verify hypotheses (i) and (ii) the Lemma 13. Indeed, since each transition in θ i t i a has been reset along

θ 1 t 1 a θ ′ 1 • • • θ ′ i-1
, it has been reset since time η 1 . Since the global time associated with (M, ν) is equal to η id i , these valuations are bounded by above by the value

(η i -d i ) -τ 1 ≤ (η i -d i ) -(1 -b 2 ) = η. Second, we have η +δ = η i -(1-b 2 ) < b 2 , as desired (this follows from the inequalities 1 -b 2 < η i < 1).
Finally, Lemma 13 thus allows to delay of b 2 the firing of the sequence θ i t i a . In particular, this will produce a letter a at time η i + b 2 > 1. The TTS N A thus accepts a timed word not recognized by N 2 P A , providing a contradiction.

⊓ ⊔

Using the results of Propositions 12 and 14, we deduce that:

Theorem 15. T PN A ⊏ T PN P A .

About Atomic and Intermediate policies in weak and strong semantics

In this subsection, we discuss the comparison of the intermediate and atomic policies. As we will see, the situation is more complex than in the previous comparison.

On the inclusion of T PN I into T PN A . For the strong semantics, a construction has been proposed in [START_REF] Bérard | Comparison of different semantics for time Petri nets[END_REF] to transform any TPN with intermediate policy into an equivalent (w.r.t. weak timed bisimilarity) TPN with atomic semantics. A first attempt was thus to adapt this construction for the weak semantics. But studying this construction, we noticed that it is erroneous (even for the strong semantics). We present below an example exhibiting the error.

t ′ , c, [0, 1] t, a, [0, 1] •• p N3 t ′-, c, [0, 1] t ′+ , τ, [0, 0] t -, a, [0, 1] t + , τ, [0, 0] •• p pt p t ′ N ′
3 Fig. 6. A counter example to the construction of [START_REF] Bérard | Comparison of different semantics for time Petri nets[END_REF].

Example 2. Consider the net N 3 depicted on the left of Figure 6. The application of the construction proposed in [START_REF] Bérard | Comparison of different semantics for time Petri nets[END_REF] leads to the net N ′ 3 depicted on the right of Figure 6. According to [START_REF] Bérard | Comparison of different semantics for time Petri nets[END_REF], we should have, under the strong semantics, the relation N 3 I ≈ N ′ 3 A . However, it is easy to verify that in the TTS N ′ 3 A the letter c can be read after 2 times units (with the timed word (a, 1)(a, 1)(τ, 1)(τ, 1)(c, 2)) whereas it is not possible in N 3 I , thus proving that the construction proposed in [START_REF] Bérard | Comparison of different semantics for time Petri nets[END_REF] is erroneous.

This example leaves open the question of the inclusion of T PN I into T PN A for the strong semantics, and then for this semantics both inclusions are left open. For weak semantics, this inclusion is also open, but we show below that the converse inclusion is false.

Non inclusion of T PN

A into T PN I . We exhibit a TPN with atomic semantics which cannot be expressed in an equivalent way by any TPN with intermediate semantics (with weak elapsing of time). This is formally stated in the Proposition below. We consider the TPN N 1 represented on Figure 1. Interpreted in weak atomic semantics, the firing of the a-labelled transition does not newly enable transition labelled by c. This transition thus shares a token with transition a while preserving a time reference to the origin of global time, what is impossible in intermediate semantics.

Proposition 16. There exists no TPN N (even unbounded) such that N I ≈ N 1 A .

Proof. We proceed by contradiction and assume there exists such a TPN N , and denote by N its number of transitions, and b the smallest positive upper bound of its intervals. As in the proof of Proposition 14, we first exhibit a particular execution ρ of N I : Lemma 17. Let (η i ) 1≤i≤k be a set of timestamps such that for any

1 ≤ i ≤ k, 1 -b 2 < η i < η i+1 < 1 and k ≥ N + 1.
There exists a run ρ in N I of the following form:

ρ : 1-b 2 ---→ I θ1 -→ I d1 -→ I θ ′ 1 -→ I t 1 a -→ I θ ′′ 1 -→ I • • • θi -→ I di -→ I θ ′ i -→ I t i a -→ I θ ′′ i -→ I • • • θn -→ I dn -→ I θ ′ n -→ I t n a -→ I θ ′′ n -→ I
such that for any position i, Λ(t i a ) = a, the transition t i a occurs at time η i , Proof. We present here the details of the construction of the run ρ. This construction proceeds in three steps.

d i > 0, Untimed(θ i ) = Untimed(θ ′ i ) = Untimed(θ ′′ i ) = ε, Duration(θ ′ i ) = Duration(θ ′′ i ) = 0,
First step : construction of the structure of ρ. This step is the most involved one. We deeply use the bisimulation property between N 1 A and N I to build a first sequence. Figure 7 illustrates this construction. A dashed arrow between two configurations meens that these configurations are bisimilar. The direction of this arrow indicates which implication is used to obtain the bisimilarity (from N 1 A to N I , or conversely). In this figure, we omit the index I and A which should be associated to each of the step. We now detail step by step how this 2 ) time units 3 and then fire the a transition. By bisimulation, there exists a path in N I , written as u 1 t 1 a u ′ 1 , leading to a bisimilar configuration and such that u 1 only contains silent transitions and is of duration exactly η 1η 0 , t 1 a is labelled by a, and u ′ 1 is labelled with internal actions and has a null duration. The bisimilarity is indicated by a top-down dashed arrow since it is the existence of a path in N 1 A that implies the existence of a path in N I . Next, one can fire instantaneously c in N 1 A what implies the existence of a path σ 1 t 1 c σ ′ 1 , leading to a bisimilar configuration, and such that σ 1 and σ ′ 1 are labelled by internal actions and of null duration, and t i c is labelled by c. Then, we use the bisimulation property in the converse direction: the configuration reached after the prefix σ 1 , which is labelled by internal actions and of null duration, is still bisimilar with the configruation of N 1 A reached after the a transition, as indicated by the bottom-up dashed arrow (the only path in N 1 A of null duration and labelled by internal actions is the empty path). Then the same reasoning is applied from these two bisimilar configurations, and can be repeated arbitrarily many times.

N I : N1 A : (M0, 0) (q, 0) (M0, 1 -b 2 ) (q, 1 -b 2 ) 1 -b 2 1 -b 2 η1 -η0 a u1 t 1 a u ′ 1 c σ1 t 1 c σ ′ 1 η2 -η1 a u2 t 2 a u ′ 2 c σ2 t 2 c σ ′ 2
Fig. 7. Using bisimulation to build ρ.

Second step : definition of θ i , θ ′ i and d i . This step is much simpler. We simply split the paths u i in three (eventually empty) parts. The duration of u i is equal to η iη i-1 > 0, thus we can consider the last positive delay step appearing in u i , and denote it by d i . Then we simply write

u i = θ i • d i • θ ′ i .

Third step : definition of θ ′′

i . In this last step, we will define the path θ ′′ i as a prefix of the path u ′ i • σ i . This path possesses the following properties : its duration is null, it is labelled by internal actions, and there exists a c labelled transition t i c immediately firable after it. The only property missing is that its last transition newly enables transition t i c . To obtain this property, we will show that we can find a prefix of u ′ i • σ i having this additional property. Therefore we prove the following lemma, specific to the intermediate semantics:

Lemma 18. Let N be a TPN, and consider an execution of N I , denoted by ρ :

(M 1 , v 1 ) t1 -→ I (M 2 , v 2 ) t2 -→ I • • • (M n , v n ) tn -→ I (M n+1 , v n+1 ) t -→ I , such that for any i ∈ {1, . . . , n}, ↑enabled I (t, M i , t i ) = f alse.
Then, for any i ∈ {1, . . . , n}, the execution ρ

i : (M i , v i ) t -→ I ti -→ I exists in N I .
Proof (of Lemma 18). By definition of the predicate ↑enabled I , and since t is newly enabled by no t i and firable from (M n+1 , v n+1 ), we have that t is enabled in every M i and, because the sequence t 1 . . . t n is instantaneous (no delay step), v n+1 (t) = v i (t) for any i. Thus, t is firable from the configuration (M i , v i ). Moreover, we can notice, as t i is enabled in M i , that ↑enabled I (t, M i , t i ) = f alse implies that t i is still enabled after the firing of t from M i , and has not been newly enabled: ↑enabled I (t i , M i , t) = f alse. As a consequence, t i is firable from the configuration reached after firing t, thus proving the Lemma.

⊓ ⊔

Lemma 18 entails that there exists a transition in the sequence t i a u ′ i σ i which newly enables the transition t i c . Indeed, if it is not the case, the result of Lemma 18 implies that after the firing of u i , one can first fire the transition t i c , and immediately after fire the transition t i a . This leads to the production of a letter a after the letter c, which in impossible in N 1 A , thus leading to a contradiction. Finally, we define θ ′′ i as the (eventually empty if t i c is newly enabled by t i a ) prefix of u ′ i • σ i up to the first transition that newly enables t i c . This concludes the proof of the existence of the sequence ρ. ⊓ ⊔

We now return on the proof of Proposition 16. First, we claim that there exists an index i ∈ {1, . . . , k} such that each transition t appearing in θ ′ i t i a θ ′′ i has already been fired since time η 1 , i.e. t also appears in θ

′ 1 t 1 a θ ′′ 1 • • • θ i-1 θ ′ i-1 t i-1 a θ ′′ i-1 θ i .
The reasonning is similar to the one of the proof of Proposition 14: by contradiction, if there exists no such position, then we can find, for each index i ∈ {1, . . . , k}, a transition, denoted t i , that never appears before. The choice of k verifying k ≥ N + 1 then implies that there exist two positions i = j such that t i = t j , thus yielding a contradiction. We can now fix an index i verifying the above described property.

We now show that Lemma 13 can be applied to the portion of ρ associated with the sequence di -→

θ ′ i -→ t i a -→ θ ′′ i -→.
More precisely, let (M, ν) be the configuration reached after firing (1 -b

2 )θ 1 . . . t i-1 a θ ′′ i-1 θ i , the delay δ is equal to d i , the sequence t 1 • • • t n corresponds to θ ′ i t i a θ ′′ i , and η is defined as

(η i -d i ) -(1 -b 2
). In the intermediate semantics, when a transition is fired, its clock is reset if it is still enabled. This property allows, together with timing constraints on the run ρ, to verify hypotheses (i) and (ii) of the Lemma 13. Indeed, we first have that each transition in θ ′ i t i a θ ′′ i has been reset since time τ 1 . Since the global time associated with (M, ν) is equal to η id i , these valuations are bounded by above by the value

(η i -d i ) -η 1 ≤ (η i -d i ) -(1 -b 2 ) = η. Second, we have η +δ = η i -(1-b
2 ) < b 2 , as desired (this follows from the inequalities 1 -b 2 < τ i < 1). Lemma 13 thus allows to delay of b 2 time units the firing of thesequence θ ′ i t i a θ ′′ i . Moreover, as the transition t i c is newly enabled by the last transition of t i a θ ′′ i , and can be immediately firable after θ ′′ i in ρ, deduce that this immediate firing can also be performed in the delayed execution. We thus obtain an execution in N I with a c action following an a action after time 1, which is impossible in N 1 A , thus yielding a contradiction.

⊓ ⊔

Conclusion

We have studied in this paper the model of Time Petri Nets under a weak semantics of time elapsing, allowing any delay transition. We have first proven that for the intermediate memory policy, the set of reachable markings coincides with the reachability set of the underlying untimed Petri net. As a consequence, many verification problems are decidable for weak intermediate TPN. On the other hand, we have proven that the two other memory policies, namely atomic and persistent atomic, allow to simulate Minsky machines and thus are undecidable. Finally, we have studied expressiveness and have proven that (i) the atomic semantics is strictly less expressive than the persistent atomic one and (ii) the atomic semantics is not included in the intermediate one.

In further work, we plan to investigate properties concerning executions of weak intermediate TPN; such as time-optimal reachability, or LTL model checking. Indeed, while discrete markings are the same, the executions are different from those accepted by the underlying Petri net. Concerning expressiveness, we conjecture that intermediate and atomic semantics are uncomparable in general, and that bounded weak TPN are strictly less expressive than timed automata (without invariants).

Aknowledgments. We would like to thank Fabrice Chevalier for fruitful discussions and the anonymous reviewers for their insightful comments.

3 )

 3 The boundedness problem : does there exist b ∈ N such that for all (M, ν) ∈ Reach(N ) s and for all p ∈ P , M (p) ≤ b ?

Fig. 2 .

 2 Fig. 2. Encoding decrementation with strong semantics.

  by definition of δ(t), and since the intervals associated with transitions are all singletons.

Fig. 4 .Proposition 12 .

 412 Fig. 4. From atomic to persistent atomic semantics.

  Figure reads. Initially, due to weak semantics, the net N in intermediate semantics can choose to delay 1-b 2 time units. By bisimulation (and because there are no silent transitions in N 1 ), the same delay leads to a bisimilar configuration in N 1 A , thus indicated by a bottom-up dashed arrow. From this configuration one can in N 1 A delay η 1 -(1 -b

  and there exists a transition t i c , labelled by c, newly enabled by the last transition of t i a θ ′′ i and (immediately) firable from the configuration reached after θ ′′ i .

We denote by η0 the value 1 -b 2 .

⋆ Partly supported by the French project DOTS (ANR-06-SETI-003) and by the European project QUASIMODO (FP7-ICT-STREP-214755). ⋆⋆ Supported by a post-doctoral scolarship from DGA/ENS Cachan and by the MIUR-PRIN project PaCo -Performability-Aware Computing: Logics, Models and Languages.