Marco Baldovin

Grégoire Vallet

Gaetan Hagel

Emmanuel Trizac

Caroline Champenois

Gaetan Hagel

Gaëtan Hagel

Observation of self-diffusion in a strongly coupled non-neutral plasma

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.

I. INTRODUCTION

Laser-cooled clouds of atomic ions stored in a radiofrequency trap are practical realisations of a finite-size One Component Plasma (OCP) in the strongly coupled regime. The OCP is a reference model in the study of strongly coupled Coulomb systems [1]. By tuning the density n i and the temperature T of the sample, different regimes can be explored from gas to liquid and crystals. Standard kinetic theories [2] fail to describe transport plasma properties under conditions of strong Coulomb coupling because they neglect effects of spatial and temporal correlations induced by nonbinary collisions [3]. This fundamental problem needs to be solved to accurately model the transport properties, and equations of state of dense laboratory and astrophysical plasmas.

Even if measurements are important to benchmark potential models and test plasma theories out of the conventional plasma regimes [4][5][6], few experiments can give access to relevant diffusivity parameters like the ion selfdiffusion constant. The main contributions so far have come from experiments based on ultra-cold neutral plasmas created by photo-ionisation of an ultra-cold atomic cloud [3]. Despite the short lifetime of these neutral plasmas, the experimental data allow for studying the effects of strong coupling on collisional processes, which is * marco.baldovin@cnr.it of interest for dense laboratory and astrophysical plasmas [3,7,8]. Here, we propose to use another model system for a strongly correlated plasma, which benefits from an infinitely long trapping lifetime. The presence of the external confining potential induces a finite size for the system, while Coulomb repulsion forces tend to maximize the mutual distance between ions. The interplay between these two effects typically results in a shell structure [9,10], which can be fairly modeled by a regular lattice geometry. This will be the starting point for the building of our model.

The next section introduces the experimental set-up and the characteristics of the system used. In section III, we propose a model for the diffusion process, to identify the parameters that are accessible to experiments. Section IV compares the predictions of the model for the stationary regime to the results of the experiments. Section V proposes experimental strategies to measure the diffusivity of the ions, while conclusions are drawn in section VI.

II. EXPERIMENTAL SETUP

In the experiment considered as a support for studying diffusion properties, clouds of few hundreds to thousands Ca + ions are stored in a linear radio-frequency (rf) quadrupole trap where the role of the neutralising particles is played by the confining potential [11]. The technical details concerning the set-up can be found in [START_REF] Champenois | AIP Conference Proceedings[END_REF][START_REF] Kamsap | [END_REF] and we recall here the useful facts. In the pseudo-potential approximation [14,[START_REF] Gerlich | State-selected and state-to-state ionmolecule reaction dynamics, Part I[END_REF] which is relevant in the context of these experiments, the effective trapping potential can be described by

V trap (x, y, z) = 1 2 mω 2 r (x 2 + y 2 ) + 1 2 mω 2 z z 2 (1)
with m the mass of a single ion, and (x, y, z) its Cartesian coordinates. The potential depth is of several eV and its cylindrical symmetry is related to the quadrupole geometry which is built on four electrodes along the direction Oz. By means of Doppler laser-cooling, temperatures T of the order of 1 to 100 mK can be reached [16,17]. The density n i of the cloud is controlled by the strength of the rf trapping field, and scales with ω 2 r . Through the tuning of the density and of the temperature, the plasma parameter

G p = q 2 4π 0 ak B T (2)
can be tuned over several orders of magnitude, where a is the Wigner-Seitz cell radius defined as (3/4πn i ) (1/3) , when n i has reached the cold limit for the density [17], q the charge of the Ca + ion and k B the Boltzmann constant. It is possible to assign a gas, liquid and crystal state to such a sample, using the two body correlation function [18,19]. With the control on Doppler lasercooling and on the steepness of the trapping potential, the plasma parameter of a trapped-ion based finite OCP can span from gas (G p lower than 0.1) to liquid (G p of the order of 1 to 100) and crystal phases (G p larger than 200). For temperatures lower than 1 K, the thermal kinetic energy is small compared to the trapping and Coulomb repulsion potential energies and ions arrange in a stationary structure that minimise the total potential energy, to form what is called a Coulomb crystal, of an ellipsoidal shape, characterised by a radius R and length L, with an aspect ratio R/L controlled by the trapping potential aspect ratio ω 2 r /ω 2 z [17]. One can show [11,17,20,21] that starting from a gaseous ("high temperature") region and cooling to the liquid phase, the ion density is uniform in the cloud, except for an outside layer, of a width of the order of few µm, small compared to the sizes of the cloud which are of the order of several hundreds of µm (see Fig. 1). When cooled further to the crystal phase, the outer shape of the cloud does not change and the mean density remains equal as in the liquid phase. It is of the order of 10 8 cm -3 for the experiments mentioned in the following. An example of this structure formed by 1240 ±50 ions is shown on Fig. 1a, obtained by the image of their laser induced fluorescence on an intensified CCD camera. The fluorescence is driven by the laser excitation used for Doppler cooling. In the case of such 3D-trapped sample, a single laser beam propagating along the Oz direction is sufficient for Doppler laser cooling. Here it propagates toward the positive z. The pixel signal is proportional to the number of emitted photons integrated along the line of sight, which is one of the direction perpendicular to the trap axis Oz. Such systems can be considered as a finite size realisation of a one-component plasma and we can get an estimate of the self-diffusion coefficient D from the work of Daligault, based on molecular dynamics simulations [22] or computation by practical model [23]. For the typical range of parameters covered by the experiments detailed in the following, these results lead to a value for D of the order of 1 to 10 × 10 -6 m 2 /s.

Introducing extra lasers, the internal structure of Ca + allows to shelve the electronic state in a metastable state [24,25], or to trap it in a dark state involving a coherent superposition of two [26,27] or three [28,29] stable and metastable states. In these three situations, the ions cannot be excited by the cooling laser and thus do not scatter photons. These dark states have a lifetime in the millisecond to second range, which is far longer than the 6.9 ns lifetime of the excited state involved in the Doppler laser cooling. Because of several experimental imperfection like Doppler effect or collisions, the shelving or coherent trapping process does not involve all the ions at the same time and some ions are still enrolled in the laser cooling process, inducing a cycling of absorption and spontaneous emission. The net recoil induced by each cycle is responsible for a radiation pressure [START_REF] Metcalf | Laser cooling and trapping[END_REF] that is applied only on these bright ions. This state selective force is responsible for the spatial segregation between the bright and dark ions that is visible on figure 1b. In spite of this, thermalisation between the dark ions and the laser cooled ones allows to keep the cloud in a steady state, where segregation by the state-selective radiation pressure offers a unique tool to measure diffusion properties within a strongly correlated non-neutral plasma.

III. ANALYTICAL MODEL

In this section we propose a simplified model for the considered system, assuming that the dynamics can be mapped, at a microscopical level, to an exclusion process on a lattice. Such a regular geometry is inspired by the ordered structure at short scale of OCP [9]. The dynamics does not allow the simultaneous presence of more than one ion in a "site" of the lattice, hence the restriction to an exclusion dynamics, in the spirit of the celebrated Asymmetric Exclusion Process [START_REF] Derrida | [END_REF]. At first we will neglect the internal state long lived dynamics, i.e. we will assume that the ions are not allowed to pass from the bright to the dark state and vice-versa. We will say that the states of the ions are "frozen". This simplifying hypothesis accounts for the fact that the typical life-times of the shelved states are much longer than the characteristic times for the spatial displacement of the ions. While the main features of the model will be already present in the one-dimensional setting, an accurate functional form to be compared with experimental density profiles will be provided by considering a more realistic ellipsoidal geometry. In subsection III B the consequences of a nonfrozen regime will be discussed, as to provide strategies for future experiments in this regime also.

A. Limit of infinite lifetime of the dark state

Derivation of the model

We consider a one-dimensional lattice, extending along the z axis, composed by N cells with size a. Each site of the lattice is occupied by an ion, and each ion can be found either in a bright or in a dark state. For the moment, let us assume that transitions to and from the metastable state do not occur on the time-scale of the observed dynamics (frozen states). The only allowed evolution is the swapping of neighbour particles along the z direction. We will denote by γ u the rate at which a bright particle placed on site n exchanges its position with a dark particle at site n + 1, by γ d the exchange rate toward site n-1. The two rates are determined by thermal fluctuations and by the effect of radiation pressure, as it will be discussed in the following. In particular, we expect γ u and γ d to be identical in the absence of external forces; the cooling laser along the z axis, acting only on the bright ions, leads instead to unbalance.

Let us denote by p n (t) the probability that, at time t, the site n is occupied by a bright particle. We can write down an evolution equation for the {p n }, 1 < n < N by recalling the swapping rules introduced before:

dp n dt = -γ u p n (1 -p n+1 ) -γ d p n (1 -p n-1 ) + γ u p n-1 (1 -p n ) + γ d p n+1 (1 -p n ) . (3)
The first two terms on the right hand side are loss terms, as they account for the cases in which a bright particle initially present in the nth site leaves it and goes to a neighbour site, occupied by a dark ion; the two gain terms stand for the opposite transitions. The evolution equations are completed by the boundary conditions

dp 1 dt = -γ u p 1 (1 -p 2 ) + γ d p 2 (1 -p 1 ) dp N dt = -γ d p N (1 -p N -1 ) + γ u p N -1 (1 -p N ) . (4)
This scenario relies on the approximation that p n and p n+1 are independent probabilities. A more accurate description of the system would involve the conditional probabilities of finding bright particles in the neighbor sites n -1 and n + 1, given the occupation in the site n. The equations ruling the time evolution of such two-site probabilities would requires three-site terms, leading to a hierarchy of coupled equations hardly addressable. The factorization hypothesis, which is reminiscent of Boltzmann's Stosszahlansatz ("molecular chaos" hypothesis), allows to close the equations at the first level of the hierarchy (see, e.g., Chapter 3.3 of Ref. [START_REF] Huang | Statistical Mechanics[END_REF]). The quality of the approximation is checked a posteriori, by comparing the predictions of the model with the experimental results.

Physical interpretation

It is useful to switch to a continuous description, as it is typically done when considering the hydrodynamic behaviour of a gas of particles. Passing to such a coarsegrained model yields the twofold benefit of (i) getting insight about the physical interpretation of the terms ruling the evolution and (ii) finding explicitly -at least in some cases -the stationary state. The price to pay is a lower accuracy on the short length-scales. To this end, we define a density of bright particles as

z 0 dz ρ(z ) = z/a n=1 p n ,
where x denotes the largest integer smaller than x. We are interested in the limit N 1, with the total length of the lattice, L = N a, finite. We can then derive an evolution equation for the density from Eq. ( 3), by substituting

p n → aρ(z) p n±1 → aρ(z ± a) aρ(z) ± a 2 ∂ z ρ(z) + a 3 2 ∂ 2 z ρ(z) + ...

One gets

∂ t ρ = -V ∂ z ρ + 2aV ρ∂ z ρ + D∂ 2 z ρ , (5)
where

V = a(γ u -γ d ) D = a 2 (γ u + γ d ) 2 . ( 6
)
The evolution equation ( 5) is a Burgers' equation with viscosity [START_REF] Ablowitz | Nonlinear dispersive waves: Asymptotic analysis and solitons[END_REF], often encountered when taking the hydrodynamic limit of asymmetric exclusion processes [START_REF] Mallick | [END_REF]35]. The first term on the r.h.s. is a systematic drift due to the presence of the external forcing by the radiation pressure; the second one, nonlinear in the density, accounts for the exclusion processes that favor the occupation of a site by a single particle; the last term, proportional to the second derivative of the density, accounts for diffusion. Equation ( 5) can be written in the form of a conservation law

∂ t ρ = -∂ z J(z) , (7)
where

J(z) = V ρ(z) -aV ρ 2 (z) -D∂ z ρ(z) (8)
plays the role of a density current. A physical interpretation of the above scenario can be obtained by considering an effective potential U (z) felt by a bright particle moving in the lattice, as sketched in Fig. 2. This potential accounts for the repulsion between adjacent particles and for the effect of the external force F pointing along z > 0 which tilts the potential profile.

Let us call ∆U the potential barrier that a bright particle placed in z needs to overcome in order to swap with a neighbour dark particle placed in z + a. Because of the tilt, the potential barrier associated with a swapping with a dark particle in z -a is thus ∆U + F a. An Arrhenius law can be used to relate the transition rates with these energies like

γ u ∝ exp [-β∆U ] γ d ∝ exp [-β(∆U + aF )] , (9)
where β = (k B T ) -1 is the inverse temperature. From Eq. ( 6) follows

V = a(γ u + γ d ) γ u -γ d γ u + γ d = a(γ u + γ d ) tanh(βaF/2) a 2 2 β(γ u + γ d )F , (10)
the last approximation holding if aF is small with respect to the typical thermal fluctuations. In this case, we have an explicit expression for the mobility

µ = V F a 2 2 β(γ u + γ d ) . (11)
Comparing this result to Eq. ( 6) one gets

µ = βD , (12)
i.e. Einstein's relation, which is indeed expected to hold at equilibrium. The characterization of the coefficients µ and D is the goal of the measurement strategies proposed in the following sections.

Stationary profile

The stationary solutions to Eq. ( 5) are of the form

ρ(z) = 1 2a + D σV a tanh z -z 0 σ , (13)
where z 0 and σ are free parameters that are fixed by the boundary conditions of the problem: z 0 represents the center of the transition front between the bright and dark region of the ion cloud and depends on the ratio of bright and dark ions; σ is a typical width of this transition front. By substituting Eq. ( 13) in the density current given by Eq. 8 one finds that in the stationary state, the density current is

J(z) = V 4a - D 2 aV σ 2 .
Because the stationary state is an equilibrium state where no density currents are present, we impose the additional condition J = 0, leading to

σ = 2D V (14)
hence

ρ(z) = 1 2a 1 + tanh V 2D (z -z 0 ) . ( 15
)
The constant z 0 is fixed by the normalization on the ion density. Denoting by N b the number of bright ions, the relation

L 0 dzρ(z) = N b (16)
leads to

L 2a + D aV ln cosh LV 2D -tanh V z 0 2D sinh LV 2D = N b ,
hence, z 0 depends on σ = 2D/V as

z 0 = σ tanh -1 cosh (L/σ) -e (2aN b -L)/σ sinh (L/σ) . ( 17
)
If the front width |σ| is much smaller than the total length L, we can approximate the argument of the tanh -1 as

cosh (L/σ) -e (2aN b -L)/σ sinh (L/σ) 1 -2e 2(aN b -L)/σ if σ > 0 2e 2aN b /σ -1 if σ < 0 .
By recalling the small-x expansion

tanh -1 (1 -x) - 1 2 log x 2 - x 4 + O(x 2 ) ,
z 0 takes the simple form:

z 0 L 1 -N b N + o (σ) if V > 0 N b L N + o (σ) if V < 0 . ( 18
)

Ellipsoidal geometry

The experimental setup described in Section II confines the ion cloud in an ellipsoidal shape. It is experimentally verified that the action of the forcing laser does not alter this geometry. It is thus natural to model the dynamics as taking place in a 3d lattice with reflecting boundaries on an ellipsoidal domain. Along the z axis, the dynamics is the one described in the previous paragraphs. Along the xy plane, since no external forces are exerted, we expect self-diffusion.

The evolution equation for the density of ions is now given by

∂ t ρ = -V ∂ z ρ + 2a 3 V ρ∂ z ρ + D(∂ 2 x + ∂ 2 y + ∂ 2 z )ρ , ( 19
)
where ρ is a volume density, justifying the different dimensional factor in front of the nonlinear term, with respect to Eq. ( 5). Since, for every fixed value of z, the dynamics on the accessible domain of the xy plane is purely diffusive, the stationary solution is a density profile ρ(x, y, z), whose dependence on the x and y variables is only due to the constraint of the confining potential:

ρ(x, y, z) = Θ 4R 2 z L - z 2 L 2 -x 2 -y 2 ρ(z) a 2 , ( 20
)
where Θ(•) is the Heavyside step-function, and R the transversal semi-axis of the ellipsoid. The linear density ρ(z) is defined by Eq. ( 13), where the parameter z 0 needs to be fixed by taking into account the new geometry of the system. By recalling that the area of the section perpendicular to the z axis measures

S(z) = 4πR 2 z L - z 2 L 2 , (21)
the density of bright ions projected to the z axis is ρmeas (z) = dx dy ρ(x, y, z)

= S(z)ρ(z) = 2πR 2 a 3 z L - z 2 L 2 1 + tanh V 2D (z -z 0 ) . (22)
This is the quantity that is actually measured in experiments. The value of z 0 can be fixed again by imposing the normalization condition (16). One gets the relation

3 L 0 dz L z L - z 2 L 2 1 + tanh V 2D (z -z 0 ) = N b N , (23)
which can be inverted numerically to find the value of z 0 , given V /D and the density of bright ions (see Fig. 5 in the next Section).

B. Dynamics with a finite lifetime of the dark state

So far we have neglected the shelving dynamics and have assumed that the bright or dark state of each ion is fixed during the observational time. In this section we relax this hypothesis, assuming that the ions can switch from bright to dark, and vice-versa, during the dynamics. When variations of this state occur much more frequently than the typical displacements, a uniform distribution of bright ions along the z axis is expected to be observed. If, instead, the characteristic dark state lifetimes are of the same order of the ones responsible for the swapping of neighbour ions, nontrivial competing effects are expected to arise.

Let us consider again the simple one-dimensional model of Eq. ( 3) where we now include the possibility that a bright ion becomes dark with rate Γ d , and a dark ion becomes bright with rate Γ b . The shelving process is induced by laser excitation and thanks to a coherent three photon process, these two rates can be tuned independently [28,29]. The discrete lattice model becomes then

dp n dt = Γ b (1 -p n ) -Γ d p n -γ u p n (1 -p n+1 ) -γ d p n (1 -p n-1 ) + γ u p n-1 (1 -p n ) + γ d p n+1 (1 -p n ) , (24)
leading to the continuous-space evolution equation

∂ t ρ = Γ b a -(Γ b + Γ d )ρ -V ∂ z ρ + 2aV ρ∂ z ρ + D∂ 2 z ρ . ( 25
)
The stationary state for the density cannot be expressed in closed form. However, some quantitative predictions can still be made concerning this state. First of all, the average number of bright ions is controlled by the relation

N b = Γ b Γ d + Γ b N (26)
which results from a balance between bright and dark ions in the stationary state. For a force exerted by the laser on the bright ions oriented toward the positive direction of the z axis and for a transition region size which is smaller than the whole cloud, the density is expected to be close to a -1 (only bright ions) when z approaches L and to 0 (only dark ions) when z 0. In the former regime, the stationary state coming from Eq. ( 25) can be approximated as

∂ z ρ(L) = Γ d aV ; (27)
and in the latter case by

∂ z ρ(0) = Γ b aV . (28)
We thus get a normalization condition for the density profile and an approximation for its behavior close to the cloud boundaries. This information is useful when devising diffusivity measurement strategies in the nonnegligible shelving regime. Other useful information are also obtained by numerical simulations of the system as ruled by Eq. ( 25) as reported in Appendix B. A rich phenomenology of bright ion density profiles can be observed when varying two scaling parameters: the ratio θ = (γ u -γ d )/(γ u + γ d ) which scales the unbalance between the swapping probabilities and ε = (Γ b + Γ d )/(γ u + γ d ) which scales the probability of state exchange relative to the mean swapping probability.

IV. EXPERIMENTAL VALIDATION

In this section, we present an experimental validation of the model in the frozen limit (ε 1) described in section. III A. In practice, to go from the situation of Fig. 1a to Fig. 1b, we leave the Doppler cooling lasers on and we add an extra laser beam at 729 nm that drives part of the ions in a metastable state (shelving process) [24], like sketched in Fig. 3. The lifetime of the metastable state D 5/2 is 1.2 s [36], which is far longer than the swapping times 1/γ u and 1/γ d . The experimental data are extracted from the pictures of the laser induced fluorescence emitted by a cloud of ions, with and without shelving like explained in section II. The data processing allows to extract the characteristics of the ellipsoid formed by the picture of the cloud and to deduce the cloud sizes L and R based on the pixel dimension, which is 13 µm and the optical magnification measured to 12.6±0.1. By defining the boundary of cloud picture, it is possible to integrate the signal in the other direction transverse to the force direction and compute an integrated signal I(z). By taking care of setting the experimental conditions to have a linear response of the photon counting pixels and a uniform laser beam intensity over the whole cloud, the integrated signal I(z) can be used as a proxy for the density of bright ions ρ(z), with a scaling factor η, taking into account the probability for each ion to scatter photons and the detection efficiency and detector gain. The signal integrated over the whole ellipsoid is proportional to the number of bright ions in the cloud. The ratio of the two signals collected with and without shelving of part of the ions in a dark state gives access to

N b /(N b + N d ).
In the case of Fig. 1, this ratio is 0.58 ± 0.02, resulting in a ratio Γ b /Γ d = 1.38 ± 0.04.

In Fig. 4a, the integrated signal I(z) is plotted as a function of z/L, in the absence or presence of shelving of part of the ions in a dark state, for the same cloud. In both cases, I(z) is expected to obey Eq. 22 scaled with η, with V = 0 in the absence of shelving. In this case, a 1-parameter fit of the parabolic profile allows to defines the common prefactor η4πR 2 /(2a 3 ). Therefore, in the case of shelving, the fit of the functional form Eq (22) only includes two free parameters, namely z 0 and D/V . In Fig. 4b the inverse of the latter is plotted for different values of the radiation pressure force exerted by the cooling laser. To keep the temperature of the sample constant along the experiment, the laser cooling exci- The black line is the best fit with a parabolic profile, which is expected because the projected density of ions is just proportional to the transversal section of the cloud, Eq. ( 21). Red curve : number of collected photons I(z), as a function of the relative coordinate z/L, in the case where part of the ions are shelved in a metastable state, so that they do not scatter photons. The black curve is a fit which obeys Eq. ( 22), taking into account the scaling factor given by the previous fit. Panel (b) shows the value of V /D, obtained from this fit, for several signal profiles I(z) obtained with different values of the effective pressure force F .

tation process is set in the linear regime and the cooling laser beam is split in two counter-propagating beams and the balance between the two beam intensity is tuned to reach a variable effective pressure force with a conserved number of scattered photons. The method used to estimate the value of the effective force is explained in Appendix A. The linear behaviour observed in Fig 4b

is in agreement with the theoretical prediction, assuming a constant temperature. Indeed, from Eqs. ( 6) and the linear approximation of Eq. ( 10), one easily obtains 18) for the linear-lattice geometry. The green solid line is computed taking into account the ellispoidal geometry, see Eq. ( 23).

V D = F k B T . (29)
From a linear fit of the plot in Fig. 4b it is then possible to evaluate the temperature of the system, which in this case is

T 22 mK .
Given the Wigner-Seitz radius estimated from the density to 13.8 µm, this value gives a plasma parameter G p = 55, which defines a liquid phase. This temperature is consistent with the linear approximation of Eq. ( 10) as βaF/2 0.05 for a force of 2 × 10 -21 N. Furthermore, it is compatible with the typical temperatures reported for ion clouds where only a fraction of the ions are laser cooled [37,38]. The same method applied to other sets of measures (not shown here), with different values of the density, shows the same linear dependence and leads to values of the temperature of 19 mK and 16 mK. For a precise and accurate estimation of this temperature, all the cause of uncertainty concerning the detection efficiency and the force estimation must be evaluated and reduced.

As an additional check of the validity of our model, we plot in Fig. 5 the values of z 0 obtained from the fits of the density profiles against V /D, and we compare them with our theoretical prediction (23). The agreement between the two is a further confirmation of the relevance of our description, even for situations where a front is not visible on the picture, like for the smallest values of |V /D|. This comparison shows also that the simplified relation (18), which does not take into account the ellipsoidal geometry of the system, provides a fair approximation for large enough values of the external force.

The proposed model, based on a lattice-like description of the system, turns out to provide meaningful predictions in the considered cases, even if the ion cloud is in a liquid regime where a crystalline structure is not present. This is not completely surprising. The main ingredient for the exclusion-process model is the presence of a strong local order, and this is guaranteed as soon as G p exceeds a few tens [39]: in this case the pair correlation function g(r) is known to exhibit a well-defined first peak, from which a first layer of neighbors can be identified [40]. When deriving the model no assumption was made, instead, on the long-range order of the system, which is of course missing in the liquid phase.

V. STRATEGIES TO THE DIRECT MEASUREMENT OF ION DIFFUSIVITY

The analysis proposed in the previous section shows that it is not possible to extract definitive information about the diffusivity of the particles by only looking at static measurements in the regime where the probability for the particle to change state is negligible. While the ratio k B T between D and the mobility µ can be inferred and compared with the expected value of the local temperature, no conclusion can be drawn about D alone.

In order to circumvent this issue, hereafter we propose two possible strategies that may be pursued in future experiments. The first one still considers the frozen regime, but it focuses on dynamical measures of the density profile, acquired with fast frequency. The second one does not involve dynamical measurements, but it requires that the shelving rates are of the same order as the displacement rates.

Before proceeding with the discussion of methods to measure the diffusivity, we compare here the proposed model to the already available experimental results.

A. Dynamical measurements at high frequency

Let us assume that the dark state has an infinite lifetime (frozen regime) and that the system is prepared in the stationary state described by Eq. ( 15), with a radiation pressure force pointing toward positive z. The density profile of the bright ions is characterised by z, the mean value of their position defined as:

z = 1 N b L 0 dz ρ(z)z . (30)
At time t = 0 the two counter-propagating cooling lasers are tuned in such a way that γ u = γ d , i.e. V = 0, so that the following evolution is purely diffusive. The dynamics of z is obtained from Eq. ( 5) and reads

∂ t z = D L 0 dz ∂ 2 z ρ(z)z = DL N b ∂ z ρ(L) - D N b [ρ(L) -ρ(0)] . (31)
If the initial condition is characterized by a small front width compared to the total length L, the bright ion density profile is flat at the boundaries, and verifies

ρ(0) 0 ρ(L) 1 a ,
This means that at the beginning of the relaxation evolution, until the shape of the profile changes significantly, one has

∂ t z - D N b a . (32)
The diffusivity can be thus measured from the variation of z by repeating the same protocol with an increasing waiting time.

B. Static profile for finite lifetime of the dark state

As discussed in section III B, provided that the system can be brought to a regime where the rates of bright to dark are comparable with those of the ion displacements, the slope of the density profile close to the boundaries (Eqs. ( 27) and ( 28)) give access to V , provided Γ d and Γ b have been measured previously. In the corresponding case where Γ d and Γ b can be neglected (frozen regime), the fit of the density profile gives access to the front length σ = 2D/V from which one can deduce D. The control on Γ d and Γ b can be reached by tuning the atomlaser interaction parameters in the three-photon process detailed in [28].

It would be tempting to estimate also the ratio D/V from the fit in the regime where the lifetime of the dark state is finite, skipping the preliminary measurement of σ in the frozen one. However, one must take into account that the width of the transition region is given, in this case, by the non-trivial combination of two effects: (i) the self-diffusion of the ions, which we aim at measuring; (ii) the dynamical effect due to the fact that one ion can be suddenly shelved (or unshelved), and start travelling from one side of the system to the other one. When the shelving rates are high enough, the latter effect results in a broadening of the transition region between the two coexisting phases (dark and bright) of the cloud, which does not coincide with 2D/V anymore.

At a practical level, the measurement would require two steps. First, one should perform the experiment in a regime where the internal state of the ions is fixed. By repeating the analysis discussed in Sec. III A 3, it is possible to measure the value of σ, which relates the diffusivity D to the drift V . Then the shelving laser parameters should be switched to reach a regime where Γ b and Γ d are not negligible, to measure them by fitting the stationary profile of the bright ions distribution. It can be useful to fit it with the phenomenological law:

ρ(z) (a 1 + b 1 z) + (a 2 + b 2 z) tanh z -z 0 c , (33)
where it is understood that the profile has been already normalized by the local volume of the ellipsoid, in order to make it comparable with the 1-dimensional case. The stationary profile is reported for a 1dimensional lattice model composed of L = 128 sites. The blue solid curve is obtained as a fit on the phenomenological law [START_REF] Ablowitz | Nonlinear dispersive waves: Asymptotic analysis and solitons[END_REF]. Parameters: Γ d = 2 • 10 -3 s -1 , Γ b = 10 -3 s -1 , γu = 1.25s -1 , γ d = 0.75s -1 . Average over 5 • 10 5 measures of the system's histogram, at 1s time-intervals. Time-step for the evolution: ∆t = 10 -2 s.

Although Eq. ( 33) is not an exact solution of the stationary state for the considered model, it is expected to catch the essential qualitative features of the profile when the shelving effect is small. In particular it reproduces the expected linear behaviour far from the transition region, i.e. when the value of the hyperbolic tangent is almost constant (and equal to -1 or 1 depending on which boundary is considered). The frozen limit is recovered when b 1 and b 2 are equal to zero. An example from numerical simulations is shown in Fig. (6).

The values of b 1 and b 2 obtained by the fit are nothing but the slopes appearing in the l.h.s. of Eqs. ( 27) and (28). Since Γ d and Γ b can be measured independently by a spectroscopic method, the two equations provide an estimate for V . The ratio σ having been previously measured, the diffusivity follows as D = σV /2 .

VI. CONCLUSION

In this paper, we have shown how different aspects of the atom-laser interaction can be used to measure the self-diffusion coefficient of ions within a trapped lasercooled ion cloud. The measurements rely on a model developed to describe the external dynamics of the ions, when a spatial segregation is induced by the radiation pressure that is encountered only by ions that are not trapped in a dark state. The validation of this model allows to measure the temperature of the sample and to propose several protocols that should give access to other relevant parameters, like the self-diffusivity coefficient (not directly measured in the present work). These protocols rely on the control of the internal state dynamics that is permitted by a coherent multi-photon process. Together with the control gained on the sample temperature by Doppler laser cooling, these atom-laser interaction processes should allow to measure the self-diffusivity coefficient for strongly coupled non neutral-plasma, for a broad range of plasma parameters.

The probability P e can be estimated from the number N e of photons emitted in a given time interval τ meas by the whole cloud, and the total number of ions N :

P e = τ e N e τ meas N . ( A3
)
This quantity is basically the ratio between the time spent by a single ion in an excited state and the measurement time, and it represents therefore the probability searched for. Estimating N e requires to take into account the efficiency of the photon detector and this is certainly the largest source of uncertainty in the described protocol. In the considered experimental setup, based on an intensified CCD camera, an efficiency of 1 detected photons out of 10 has be evaluated. We obtain an estimate of about 2 • 10 9 emitted photons per second (slightly fluctuating) in the measurements without shelving, with a cloud of 1240 ions, resulting in P e 0.011.

FIG. 1 .

 1 FIG. 1. Panel (a): picture of the laser induced fluorescence of a trapped ion cloud made of 1240 ±50 Ca + ions. Given an optical magnification of 12.6 ± 0.1 and a pixel size of 13 µm, the scale of the picture in pxl is very close to the scale of the cloud in µm. The cloud dimensions are : R = 107(±2)µm for the largest radius and L = 576(±14)µm for the total length. Panel (b) : same cloud with an extra laser admitted, propagating toward z > 0 in real space and on the picture, this laser is shelving part of the ions in a metastable dark state.

FIG. 2 .

 2 FIG.2. Sketch of the potential felt by a bright particle moving along the lattice. The tilt of the periodic potential (green curve) is induced by a force F pointing toward z > 0 (red line).

2 FIG. 3 .

 23 FIG. 3. Scheme of the electronic levels of Ca + involved in the two regimes considered in the experiment. Panel (a) represents the cycle of excitation (black straight arrow) and spontaneous emission (black wavy arrow) in the laser-cooling process. Panel (b) refers to the regime in which shelving is present also (excitation to the shelved state is represented by the red arrow).

FIG. 4 .

 4 Photon counts (× 10 5 )

FIG. 5 .

 5 FIG.5. Transition front coordinate z0/L as a function of the drift. Red squares represent the values measured in experiments (i.e., from the fits of density profiles, as shown in Fig.4a. The blue solid line shows the approximate estimate provided by Eq. (18) for the linear-lattice geometry. The green solid line is computed taking into account the ellispoidal geometry, see Eq. (23).

FIG. 6 .

 6 FIG.6. Numerical simulation of a case with non-negligible shelving.The stationary profile is reported for a 1dimensional lattice model composed of L = 128 sites. The blue solid curve is obtained as a fit on the phenomenological law[START_REF] Ablowitz | Nonlinear dispersive waves: Asymptotic analysis and solitons[END_REF]. Parameters:Γ d = 2 • 10 -3 s -1 , Γ b = 10 -3 s -1 , γu = 1.25s -1 , γ d = 0.75s -1. Average over 5 • 10 5 measures of the system's histogram, at 1s time-intervals. Time-step for the evolution: ∆t = 10 -2 s.

ACKNOWLEDGMENTS

The authors acknowledge the contribution of Marie Houssin in the experimental preparation of the ion cloud and its interaction with the involved lasers. This work was supported by the LABEX Cluster of Excellence FIRST-TF (ANR-10-LABX-4801), within the Program "Investissements dAvenir" operated by the French National Research Agency (ANR), and by the ANR through grant ANR-18-CE30-0013. MB was supported by ERC Advanced Grant RG.BIO (Contract No. 785932).

Appendix A: Estimation of the laser's pressure force

The radiation pressure force encountered by the ions is due to the recoil associated with the absorption of photons from the laser beam k L , where k L is the wave vector of the laser (k L = 2π/λ L , λ L = 397 nm), pointing to the direction of the beam. In the limit of a low saturation of the excited transition, the stimulated emission following this absorption can be neglected and the only emission process is spontaneous. The averaged recoil over thousands of emission is then null and we only take into account the absorption induced recoil. The mean force depends on the number of absorption/emission cycles per unit time Γ e P e , with P e the probability for an ion to be in the excited state, and Γ e the probability for an ion in the excited state to decay to the ground state. This spontaneous decay rate is the inverse of the excited state lifetime, which is τ e = 6.9 • 10 -9 s for Ca + ions [41].

For a single laser beam, the mean pressure force exerted by the laser is thus given by [START_REF] Metcalf | Laser cooling and trapping[END_REF]:

(A1)

For the measurement of V /D like shown on Fig. 4, the effective pressure force needs to be tuned. Tuning P e would induce a modification of the laser cooling efficiency and thus of the temperature. To keep the temperature constant over the experimental run, the laser beam is split in two counter-propagating beams +k L and -k L with a shared intensity x + I L and x -I L with I L the total laser intensity, and x ± the tuning parameters. Again in the limit of low transition saturation, we can assume a linear response of the ions to the laser excitation and consider that P e = P e+ + P e-with P e± scaling with each laser beam intensity x ± . The average force on the z axis is thus equal to

Appendix B: Effect of the shelving on the stationary profile

The stationary density profile of a system described by Eq. ( 25) shows a quite rich phenomenology depending on the values of the dynamical parameters. Such a stationary solution does not admit a closed analytical form, but it can be explored by mean of numerical simulations. In Fig. 7 some examples are shown for varying values of the parameters of the experiment. If the shelving rate is small enough compared to the displacement rates, the functional form provided by Eq. ( 33) approximates pretty well the measured curve. We recall that the functional form [START_REF] Ablowitz | Nonlinear dispersive waves: Asymptotic analysis and solitons[END_REF] is purely phenomenological and it is not expected to fit the profile for any choice of the parameters. Fig. 7 shows indeed that the agreement gets worse when the ratio between the shelving and the displacement rates increases.


Pascal Cuoq

Florent Kirchner

Nikolaï Kosmatov

Virgile Prevosto

Julien Signoles

Boris Yakobowski

Patrick Baudin

Richard Bonichon

Bernard Botella

Loïc Correnson

Zaynah Dargaye

Philippe Herrmann

Benjamin Monate

Yannick Moy

Anne Pacalet

Armand Puccetti

Muriel Roger

Nicky Williams

Frama-C A Software Analysis Perspective *

Frama-C is a source code analysis platform that aims at conducting verification of industrial-size C programs. It provides its users with a collection of plug-ins that perform static analysis, deductive verification, and testing, for safety-and security-critical software. Collaborative verification across cooperating plug-ins is enabled by their integration on top of a shared kernel and datastructures, and their compliance to a common specification language. This foundational article presents a synthetic view of the platform, its main and composite analyses, and some of its industrial achievements.

Introduction

The past forty years have seen much of the groundwork of formal software analysis being laid. Several angles and theoretical avenues have been explored, from deductive reasoning to abstract interpretation to program transformation to concolic testing. While much remains to be done from an academic standpoint, some of the major advances in these fields are already being successfully implemented [START_REF] Cousot | The ASTRE É analyzer[END_REF][START_REF]Polyspace[END_REF][START_REF] Dahlweid | VCC: Contractbased modular verification of concurrent C[END_REF][START_REF] Randimbivololona | Applying formal proof techniques to avionics software: a pragmatic approach[END_REF][START_REF] Williams | PathCrawler: automatic generation of path tests by combining static and dynamic analysis[END_REF]] -and met with growing industrial interest. The ensuing push for mainstream diffusion of software analysis techniques has raised several challenges. Chief among them are: a. their scalability, b. their interoperability, and c. the soundness of their results.

Point a is predictably important from the point of vue of adoptability. Scaling to large problems is a prerequisite for the industrial diffusion of software analysis and verification techniques. It also represents a mean to better understand how language idioms (e.g. pointers, unions, or dynamic memory allocation) influence the underlying architecture of large software developments. Overall, achieving scalability in the design of software analyzers for a wide range of software patterns remains a difficult question.

Point b -interoperability -enables the design of elaborate program analyses. Consider indeed the interplay between program analyses and transformations [START_REF] Demange | A provably correct stackless intermediate representation for java bytecode[END_REF], the complementarity of forward and backward analyses [START_REF] Bardin | Refinement-based CFG reconstruction from unstructured programs[END_REF], or the precision gain afforded when combining static and dynamic approaches [START_REF] Beckman | Proofs from tests[END_REF]. Yet running multiple source code analyses and synthesizing their results in a coherent fashion requires carefully thought-out mechanisms.

Point c -soundness -is a strong differentiator for formal approaches. By using tools that over-approximate all program behaviors, industrial users are assured that none of the errors they are looking for remain undetected. This guarantee stands in stark contrast with the bug-finding capabilities of heuristic analyzers, and is paramount in the evaluation of critical software. But the design and implementation costs of such high-integrity solutions are hard to expend.

The Frama-C software analysis platform provides a collection of scalable, interoperable, and sound software analyses for the industrial analysis of ISO C99 source code. The platform is based on a common kernel, which hosts analyzers as collaborating plug-ins and uses the ACSL formal specification language as a lingua franca. Frama-C includes three fundamental plug-ins based on abstract interpretation, deductive verification, and concolic testing; and a series of derived plug-ins which build elaborate analyses upon the former. The extensibility of the overall platform, and its open-source licensing, have fostered the development of an ecosystem of independent third-party plug-ins [START_REF] Ceara | Taint dependency sequences: A characterization of insecure execution paths based on input-sensitive cause sequences[END_REF][START_REF] Cuoq | Fan-C, a Frama-C plug-in for data flow verification[END_REF][START_REF] Delmas | Taster, a Frama-C plug-in to encode Coding Standards[END_REF][START_REF] Demay | Sidan: A tool dedicated to software instrumentation for detecting attacks on non-control-data[END_REF][START_REF] Moy | Automatic Modular Static Safety Checking for C Programs[END_REF][START_REF] Pariente | Formal Verification of Industrial C Code using Frama-C: a Case Study[END_REF]. This article is intended as a foundational reference to the platform, its three main analyses, and its most salient derived plug-ins.

The Platform Kernel

Architecture

Figure 1 shows a functional view of the Frama-C architecture. Frama-C is based on CIL [START_REF] Necula | CIL: Intermediate Language and Tools for Analysis and Transformation of C Programs[END_REF], a front-end for C that parses ISO C99 programs into their normalized representation: loop constructs are given a single form, expressions have no sideeffects, etc. Frama-C extends CIL to support dedicated source code annotations expressed in ACSL (see § 2.2). This modified CIL front-end produces the C + ACSL AST, an abstract view of the program shared among all analyzers.

The Frama-C kernel provides several services, helping plug-in development [START_REF] Signoles | Frama-C Plug-in Development Guide[END_REF] and providing convenient features to the end-user.

-Messages, source code and annotations are uniformly displayed; parameters and command line options are homogeneously handled. -A journal of user actions can be synthesized, and be replayed afterwards, a feature of interest in debugging and qualification contexts. -A project system, presented in § 2.3, isolates unrelated program representations, and guarantees the integrity of their analyses. -Consistency mechanisms control the collaboration between analyzers ( § 2.4).

Analyzers are developed as separate plug-ins on top of the kernel. Plug-ins are dynamically linked against the kernel to offer new analyses, or to modify existing ones. Any plug-in can register new services in a plug-ins database stored in the kernel, making these services available to all plug-ins.

ACSL

Functional properties of C programs can be expressed within Frama-C as ACSL annotations [START_REF] Baudin | ACSL: ANSI/ISO C Specification Language Preliminary design, version 1.5[END_REF]. ACSL is a formal specification language inspired by Java's JML [START_REF] Burdy | An overview of JML tools and applications[END_REF], both being based on the notion of function contract introduced by Eiffel [START_REF] Meyer | Object-oriented Software Construction[END_REF]. In effect, the specification of a function states the pre-conditions it requires from its caller and the post-conditions it ensures when returning. Among these postconditions, one kind of clause plays a particular role by saying which memory locations the function assigns, i.e. which locations might have a different value between the pre-and the post-state. For instance, Fig. 2 provides a specification for a swap function. The first precondition states that the two arguments must be valid (int) pointers, i.e. that dereferencing a or b will not produce a run-time error. In addition, the second pre-condition asks that the two locations do not overlap. \valid and \separated are two built-in predicates: ACSL features various functions and predicates to describe memory states. However, it does not introduce any notion beyond the C standard, leaving each plug-in free to perform its own abstractions over the concrete memory layout. The assigns clause states that only the locations pointed to by a and b might be modified by a call to swap; any other memory location is untouched. Finally, the post-condition says that at the end of the function, * a contains the value that was in * b in the pre-state, and vice versa.

In addition to function specifications, ACSL offers the possibility of writing annotations in the code, in the form of assertions (properties that must be true at a given point) or loop invariants (properties that must be preserved across any number of loop steps). Annotations are written in first-order logic, and it is possible to define custom functions and predicates for use in annotations together with ACSL built-ins. Plug-ins can provide a validity status to any ACSL property and generate ACSL annotations. This allows ACSL annotations to play an important role in the communication between plug-ins, as explicited in § 2.4.

Projects

Frama-C allows a user to work on several programs in parallel thanks to the notion of project. A project consistently stores a program with all its required information, including results computed by analyzers and their parameters. Several projects may coexist in memory at the same time. A non-interference theorem guarantees project partitioning [START_REF] Signoles | Foncteurs impératifs et composés: la notion de projet dans Frama-C[END_REF]: any modification on a value of a project P does not impact a value of another project P ′ . Such a feature is of particular interest when a program transformer like Slicing ( § 6.1) or Aoraï ( § 6.2) is used. The result of the transformation is a fresh AST that coexists with the original, making backtracking and comparisons easy. Another use of projects is to process the same program in different ways -for instance with different analysis parameters.

Analyzers Collaboration

In Frama-C, analyzers can collaborate in two different ways: either sequentially, by chaining analysis results to perform complex operations; or in parallel, by combining partial analysis results into a full program verification.

The former consists in using the results of an analyzer as input to another one thanks to the plug-ins database stored by the Frama-C kernel. Refer to § 6.1 for a comprehensive illustration of a sequential analysis.

The parallel collaboration of analyzers consists in verifying a program by heterogeneous means. ACSL is used to this end as a collaborative language: plug-ins generate program annotations, which are then validated by other plug-ins. Partial results coming from various plug-ins are integrated by the kernel to provide a consolidated status of the validity of all ACSL properties. For instance, when the Value plug-in ( § 3) is unable to ensure the validity of a pointer p, it emits an unproved ACSL annotation assert \valid(p). In accordance with the underlying blocking semantics, it assumes that p is valid from this program point onwards. The WP plug-in ( § 4) may later be used to lift this hypothesis. The kernel automatically computes the validity status of each program property from the information provided by all analyzers and ensures the consistency of the entire verification process [START_REF] Correnson | Combining Analyses for C Program Verification[END_REF]: "if the consolidated status of a property is computed as valid [resp. invalid] by the kernel, then the property is valid [resp. invalid] with respect to ACSL's semantics".

Fundamental Analysis: Abstract Interpretation

The Value plugin (short for Value Analysis) is a forward dataflow analysis based on the principles of abstract interpretation [START_REF] Cousot | Abstract interpretation: a unified lattice model for static analysis of programs by construction or approximation of fixpoints[END_REF]. Abstract interpretation links a concrete semantics, typically the set of all possible executions of a program, to a more coarse-grained, abstract one. Any transformation in the concrete semantics must have an abstract counterpart that captures all possible outcomes of the concrete operation. This ensures that the abstract semantics is a sound approximation of the runtime behavior of the program.

Value, and abstract interpreters in general, proceed by symbolic execution of the program, translating all operations into the abstract semantics. Termination of looping constructs is ensured by widening operations. For function calls, Value proceeds essentially by recursive inlining of the function (recursive functions are currently not handled). This ensures that the analysis is fully context-sensitive. If needed, the user can abstract overly complex functions by an ACSL contract, verified by hand or discharged with another analysis.

Abstract domains

The domains currently used by Value to represent the abstract semantics are described below.

Integer computations. Small sets of integers are represented as sets, whereas large sets are represented as intervals with congruence information [START_REF] Granger | Static analysis of linear congruence equalities among variables of a program[END_REF]. For instance, x ∈ [3..255], 3%4 means that x is such that 3 ≤ x ≤ 255 and x ≡ 3 mod 4.

Floating-point computations. The results of floating-point computations are represented as IEEE 754 [START_REF]IEEE Std[END_REF] double-precision finite intervals. Operations on singleprecision floats are stored as doubles, but are rounded as necessary. Obtaining infinities or NaN is treated as undesirable errors.

Pointers and memory. To verify that invalid (e.g. out-of-bounds) array/pointer accesses cannot occur in the target program, Value assumes that the program does not purposely use buffer overflows to access neighboring variables [35, §6.5.6:8]. Abstract representation of memory states in a C program reflects this assumption: addresses are seen as offsets with respect to symbolic base addresses, and have no relation with actual locations in virtual memory space during execution.

Memory representation is untyped. It is thus straightforward to handle unions and heterogeneous pointer conversions during abstract interpretation. The abstract memory state maps each base address to a representation of a chunk of linear memory. Each such object itself maps ranges of bits to values. Given an array of 32-bit integers t, and reading from * (((char * )t) + 5), the analyzer determines that the relevant abstract value is to be found between bits [40..47] of the value bound to &t. Using bit as unit, instead of byte, allows to handle bit-fields [35, §6.2.6.1:4] by fixing a layout strategy (the C standard itself does not specify bit-fields layout, but then again, no more than any kind of data).

Finally, the content of some memory locations is deemed indeterminate by the C standard. Examples include uninitialized local variables, struct padding, and dereferencing pointers to variables outside their scope [35, §6.2.4:2]. Having indeterminate contents in memory is not an error, but accessing an indeterminate memory location is. To detect those, the values used to represent the contents of memory locations are taken, not directly from the abstract domain used for the values of an expression, but from the lattice product of this domain with two two-valued domains, one for initializedness and the other for danglingness.

Propagation of unjoined states Value's domains are non-relational. Instead, the datastructures representing the abstract semantics have been heavily optimized for speed and reduced memory footprint, to allow the independent propagation of k multiple states. Typically, the user can set a high enough k, so that finite loops are entirely unrolled. Successive conditionals are also handled more precisely. This parameter can be adjusted on a per-function basis, alleviating for a large part the need for relational domains, by implicitly encoding relations in the disjunction of abstract states.

Alarms Each time a statement is analyzed, any operation that can lead to an an undefined behavior (e.g. division by zero, out-of-bounds access, etc.) is checked, typically by verifying the range of the involved expression -the denominator of the division, the index of the array access, etc.

If the abstract semantics guarantees that no undesirable value can occur, one obtains a static guarantee that the operation always executes safely. Otherwise, Value reports the possible error by an alarm, expressed as an ACSL assertion. This alarm may signal a real error if the operation fails at runtime on at least one execution, or a false alarm, caused by the difference in precision between the concrete and abstract semantics. More precise state propagation typically results in fewer false alarms, but lenghten analysis time. Upon emitting an alarm, the analyzer reduces the propagated state accordingly, and proceeds onwards.

Fundamental Analysis: Deductive Verification

The WP plug-in is named after the Weakest Precondition calculus, a technique used to prove program properties initiated by Hoare [START_REF] Hoare | An axiomatic basis for computer programming[END_REF], Floyd [START_REF] Floyd | Assigning meanings to programs[END_REF] and Dijkstra [START_REF] Dijkstra | A constructive approach to program correctness[END_REF]. Recent tools implement this technique efficiently, for instance Boogie [START_REF] Leino | This is Boogie 2[END_REF] and Why [START_REF] Filliâtre | Why: a multi-language multi-prover verification tool[END_REF]. Jessie [START_REF] Marché | The Jessie plugin for Deduction Verification in Frama-C, version 2.30[END_REF], a Frama-C plug-in developed at INRIA, also implements this technique for C by compiling programs into the Why language. Frama-C's WP plug-in is a novel implementation of a Weakest Precondition calculus for generating verification conditions (VC) for C programs with ACSL annotations. It differs from other implementations in two respects. First, WP focuses on parametrization with respect to the memory model. Second, it aims at being fully integrated into the Frama-C platform, and to ease collaboration with the other verification plug-ins (especially Value) as outlined in § 2.4.

The choice of a memory model is a key ingredient of Hoare logic-based VC generators that target C programs (or more generally a language with memory references). A weakest precondition calculus operates over a language that only manipulates plain variables. In order to account for pointers, memory accesses (both for reading and writing) must be represented in the underlying logic. The simplest representation uses a single functional array for the whole memory. However, this has a drawback: any update to the array (the representation of an assignment * p=v) has a potential impact on the whole memory -any variable might have been modified. In practice, proof obligations quickly become intractable. Thus, various refinements have been proposed, in particular by Bornat [START_REF] Bornat | Proving Pointer Programs in Hoare Logic[END_REF], building upon earlier work by Burstall [START_REF] Burstall | Some techniques for proving correctness of programs which alter data structures[END_REF]. The idea of such memory models is to use distinct arrays to represent parts of the memory known to be separated, e.g. distinct fields of the same structure in the "component-as-array trick" of Burstall and Bornat. In this setting, an update to one of the arrays will not affect the properties of the others, leading to more manageable VC.

However, abstract memory models sometimes restrict the functions that can be analyzed. Indeed, a given model can only be used to verify code that does not create aliases between pointers that are considered a priori separated by the model. In particular, Burstall-Bornat models that rely on static type information to partition the memory are not able to cope with programs that use pointer casts or some form of union types.

In order to generate simpler VC when possible while still being able to verify low-level programs, WP provides different memory models that the user can choose for each ACSL property. The current version offers three main models:

-The most abstract model, hoare, roughly corresponds to Caveat's model [START_REF] Randimbivololona | Applying formal proof techniques to avionics software: a pragmatic approach[END_REF].

It can only be used over functions that do not explicitly assign pointers or take the address of a variable, but provides compact VC. -The default model is store. It is a classical Burstall-Bornat model, that supports pointer aliasing, but neither cast nor union types. -The runtime model is designed for code that perform low-level memory operations. In this model, the memory is seen as a single array of char, so that most C operations can be taken into account, at the expense of the complexity of the generated VC.

As a refinement, store and runtime can avoid converting assignments into array updates when the code falls in the subset supported by hoare. In particular, variables whose addresses are not taken and plain references -pointers that are neither assigned nor used in a pointer arithmetic operation -are translated as standard Hoare logic variables. This way, the overhead of other models with respect to hoare is kept to the places where it is really needed. Once a VC has been generated, it must be discharged. WP natively supports two theorem provers: the automated SMT solver Alt-ergo [START_REF] Bobot | The Alt-Ergo Automated Theorem Prover[END_REF], and the Coq proof assistant [START_REF]The Coq Proof Assistant Reference Manual[END_REF]. Other automated provers can also be used through the multiprover backend of Why. Advantages of using a dedicated back-end rather than relying completely on Why are twofold. First, it removes a dependency over an external tool, meaning that for verification of critical software, there is one component less that needs to be assessed. Second, WP can take advantage of specific features of Alt-ergo, most notably native support for arrays and records (that occur quite often in typical VC), that are not supported by Why yet.

In contrast to Jessie, that relies on an external tool for VC generation, WP operates entirely within Frama-C. In particular, WP fills the property status table described in § 2.4 for each annotation on which it is run. The dependencies of such a status are the annotations taken as hypothesis during the weakest-precondition calculus, the memory model that has been used, and the theorem prover that ultimately discharged the VC. The memory model has a direct impact on the validity of the result: an annotation can very well be valid under model store but not under runtime, as the former entails implicit separation hypotheses that are not present in the latter. In theory, the choice of a theorem prover is not relevant for the correctness of the status, but this information is important to fully determine a trusted toolchain.

Having WP properly embedded into Frama-C also allows for a fine-grained control over the annotations one wants to verify with the plug-in. WP provides the necessary interface at all levels (command-line option, programmatic API, and GUI) to verify targetted annotations (e.g. those yet unverified by other means in Frama-C, cf § 2.4) as well as to generate all the VC related to a C function.

Fundamental Analysis: Concolic Testing

Given a C program p under test and a precondition restricting its inputs, the PathCrawler plug-in generates test cases respecting various test coverage criteria. The all-path criterion requires covering all feasible program paths of p. Since the exhaustive exploration of all paths is usually impossible for real-life programs, the k-path criterion restricts exploration to paths with at most k consecutive iterations of each loop. The PathCrawler [START_REF] Williams | PathCrawler: automatic generation of path tests by combining static and dynamic analysis[END_REF][START_REF] Botella | Automating structural testing of C programs: Experience with PathCrawler[END_REF] method for test generation is similar to the so-called concolic (concrete+symbolic) approach and to Dynamic Symbolic Execution (DSE), implemented by other tools (e.g. DART, CUTE, PEX, SAGE, KLEE).

PathCrawler starts by: a. constructing an instrumented version of p that will trace the program path exercised by the execution of a test case, and b. generating the constraints which represent the semantics of each instruction in p. The next step, illustrated by Fig. 3, is the generation and resolution of constraint systems to produce the test cases for a set of paths Π that satisfy the coverage criterion. This is done in the ECLiPSe Prolog environment [START_REF] Schimpf | ECLiPSe -from LP to CLP[END_REF] and uses Constraint Logic Programming. Given a path prefix π, i.e. a partial program path in p, the main idea [START_REF] Kosmatov | Artificial Intelligence Applications for Improved Software Engineering Development[END_REF] is to solve the constraints corresponding to the symbolic execution of p along π. A constraint store is maintained during resolution, and aggregates the various constraints encountered during the symbolic execution of π. The test generation method follows the following steps: ). (G 5 ) Compute the next partial path, π, to cover. π is constructed by "taking another branch" in one of the complete paths already covered by a previous test case. This ensures that all feasible paths are covered (as long as the constraint solver can find a solution in a reasonable time) and that only the shortest infeasible prefix of each infeasible path is explored.

(G 1 )
PathCrawler uses Colibri, a specialized constraint solving library developed at CEA LIST and shared with other testing tools GATeL [START_REF] Marre | Test sequences generation from Lustre descriptions: GATeL[END_REF] and OSMOSE [START_REF] Bardin | OSMOSE: automatic structural testing of executables[END_REF]. Colibri provides a variety of types and constraints (including non-linear constraints), primitives for labelling procedures, support for floating point numbers and efficient constraint resolution. PathCrawler is a proprietary plug-in, also available in the form of a freely accessible test-case generation web service [START_REF] Kosmatov | Online version of PathCrawler[END_REF].

6 Derived Analyses

Distilling values

The outputs of the Value plugin are twofold. In addition to emitting alarms for statements it cannot guarantee are safe ( §3), Value automatically computes a per-statement over-approximation of the possible values for all memory locations. The derived analyses below reuse those synthetic results. In each case, the analysis is sound. Value's results are used to evaluate array indexes or resolve pointers, ensuring that e.g. pointer aliasing are always detected.

Outputs: over-approximates the locations a function may write to. Operational inputs: over-approximates the locations whose initial values are used by the function. Functional dependencies: computes a relation between outputs and inputs of a function; x FROM y, t [START_REF] Bardin | OSMOSE: automatic structural testing of executables[END_REF] (and SELF) means that output x is either unchanged (SELF), or that its new value can be computed exclusively from inputs y and t [START_REF] Bardin | OSMOSE: automatic structural testing of executables[END_REF]. Program Dependency Graph (PDG): produces an intra-procedural graph that expresses the data and control dependencies between the instructions of a function, used as a stepping stone for various analyzes [START_REF] Ferrante | The program dependence graph and its use in optimization[END_REF]. Defs: over-approximates which statements define a given memory location. Impact: computes the values and statements impacted (directly or transitively) by the side effects of a chosen statement. Slicing: returns a reduced program (a slice), equivalent to the original program according to a given slicing criterion [START_REF] Horwitz | Interprocedural slicing using dependence graphs[END_REF]. Possible criteria include preserving a given statement, all calls to a function, a given alarm, etc.

Analyses such as Defs or Impact make compelling code understanding tools, as they express in a very concise way the relationships between various parts of a program. Slicing goes one step further: while it is essentially dual to the impact analysis, it also builds reduced, self-contained programs, that can be re-analyzed independently. Those three analyses are fully inter-procedural. The analyses above illustrate sequential collaboration ( § 2.4). PDG makes heavy use of Functional dependencies, while Defs, Impact and Slicing leverage the information given by PDG. Some of those analyses can optionally compute callwise versions of their results, yielding one result per syntactic call, instead of one result per function. The improved precision automatically benefits the derived analyses. All the results are stored by the Frama-C kernel, and can be reused without being recomputed.

Annotation Generator

The Aoraï plug-in [START_REF] Stouls | Aoraï plug-in tutorial, version Nitrogen-20111001[END_REF][START_REF] Groslambert | Vérification de propriétés LTL sur des programmes C par génération d'annotations[END_REF] plays a particular role among the core Frama-C plug-ins. Indeed, it is one of the few whose primary aim is to generate ACSL annotations rather than attempting to verify them. Aoraï provides a way to specify that all possible executions of a program respect a given sequence of events, namely the call and return of functions, possibly with constraints on the program's state at each event. The specification itself can be given either as a Linear Temporal Logic (LTL, [START_REF] Pnueli | The temporal logic of programs[END_REF]) formula or in the form of an automaton. In the former case, Aoraï uses ltl2ba [START_REF] Gastin | Fast LTL to Büchi automata translation[END_REF] to obtain an equivalent Büchi automaton.

Given such an automaton, Aoraï provides ACSL specifications for each function transitions for the corresponding atomic event (call or return from f), and the specification of f itself. As the automaton is not deterministic in general, Aoraï uses a set of boolean variables to represent possible active states (aorai_state_S_ * ). Functions advancing the states of the automaton provide, for each of those variables, a complete set of behaviors indicating when they are set to 1 or 0.

The specification of f comprises various items. First, at least one state among a given set must be active before the call. This set is determined by a rough static analysis made by Aoraï beforehand. In addition, for each active state, at least one transition must be activated by the call event. The main post-condition is that when the function returns, at least one state is active among those deemed possible by Aoraï's static analysis. It is refined by additional clauses relating active initial states with active final states, and the state of the program itself with the active final states. Finally the main function has an additional postcondition stating that at least one of the acceptance states must be active at the end of the function (Aoraï does not consider infinite programs at the moment).

Annotation generation is geared towards the use of deductive verification plugins such as WP and Jessie for the verification of the specification. In particular, the refined post-conditions are mainly useful for propagating information to the callers of f in an Hoare-logic based setting. Likewise Aoraï also generates loop invariants for the same purpose. However it does not preclude the use of the Value plug-in to validate its specification, and therefore attempts to generate annotations that fit in the subset of ACSL that is understood by Value.

SANTE

The Sante Frama-C plug-in (Static ANalysis and TEsting) [START_REF] Chebaro | Program slicing enhances a verification technique combining static and dynamic analysis[END_REF] enhances static analysis results by testing. Given a C program p with a precondition, it detects possible runtime errors (currently divisions by zero and out-of-bound array accesses) in p and aims to classify them as real bugs or false alarms.

The Sante method contains three main steps illustrated in Fig. The number of slices generated, hence the number of test generation sessions, is influenced by various Sante options. The all option generates a unique slice p 1 including all alarms of p, while the converse option each generates a slice for each alarm. Options min and smart take advantage of alarm dependencies (as computed by the dependency analysis) to slice related alarms together. The smart option improves min by iteratively refining the slices as long as one can hope to classify more alarms running PathCrawler on a smaller slice.

Adoption

Adoption in the academic world has stemmed from a variety of partnerships. Foremost is the Jessie plug-in [START_REF] Marché | The Jessie plugin for Deduction Verification in Frama-C, version 2.30[END_REF] developed at Inria, which relies on a separation memory model but whose internal representation precludes its combination with other plug-ins. Verimag researchers have implemented a taint analysis [START_REF] Ceara | Taint dependency sequences: A characterization of insecure execution paths based on input-sensitive cause sequences[END_REF], producing explicit dependency chains pondered by risk quantifiers. Demay et al generate security monitors based on fine-grained feedback from the Value plugin [START_REF] Demay | Sidan: A tool dedicated to software instrumentation for detecting attacks on non-control-data[END_REF]. Berthome et al [START_REF] Berthomé | Attack model for verification of interval security properties for smart card c codes[END_REF] propose a source-code model for verifying physical attacks on smart cards, and use Value to verify it. Bouajjani et al [START_REF] Bouajjani | On inter-procedural analysis of programs with lists and data[END_REF] automatically synthesize invariants of sequential programs with singly-linked lists.

The variety of objectives a static analyzer can have, and the variety of design choices for a given objective, make it difficult to benchmark static analyzers. Chatzieleftheriou and Katsaros [START_REF] Chatzieleftheriou | Test-driving static analysis tools in search of C code vulnerabilities[END_REF] have valiantly produced one such comparison, including the Value plug-in. Industrial reception has been extremely positive. Delmas et al verify the compliance to domain-specific coding standards [START_REF] Delmas | Taster, a Frama-C plug-in to encode Coding Standards[END_REF]; their plug-in is undergoing deployment and qualification. At the same company, the value analysis is used to verify the control and data flows of a DAL C, 40-kloc ARINC 653 application [START_REF] Cuoq | Fan-C, a Frama-C plug-in for data flow verification[END_REF]. Pariente and Ledinot [START_REF] Pariente | Formal Verification of Industrial C Code using Frama-C: a Case Study[END_REF] verify flight control system code using a combination of Frama-C plugins, including Value and Slicing. Their contribution includes a favorable evaluation of the cost-effectiveness of their adoption compared to traditional verification techniques. Yakobowski et al use Value in collaboration with WP to check the absence of runtime errors in a 50 kloc instrumentation and control (I&C) nuclear code [START_REF] Cuoq | Formal verification of software important to safety using the Frama-C tool suite[END_REF]. Through these successes and over the past few years, Frama-C has demonstrated its adoptability within numerous industrial environments.

Conclusion

This article attempts to distill a synthetic presentation of the Frama-C platform from a software analysis perspective. Frama-C answers the combined introductory challenges of scalability, interoperability, and soundness with a unique architecture and a robust set of analyzers. Its core set of tools and functionalitiesabout 150 kloc. developed over the span of 7 years -has given rise to a flourishing ecosystem of software analyzers. In addition to industrial achievements and partnerships, a community of users and developers has grown and strived, contributing to the dissemination of the tools. This growth, fostered by a number of active communication channels1 , should be interpreted as a testimony to the health and stability of the platform, and good omens for its future success.

/

  * @ requires \valid(a) && \valid(b); requires \separated(a,b); assigns * a, * b; ensures * a == \at( * b,Pre) && * b == \at( * a,Pre); * / void swap(int * a, int * b);

Fig. 2 .

 2 Fig. 2. Example of ACSL specification

(Fig. 3 .

 3 Fig. 3. The PathCrawler test generation method

Fig. 4 .

 4 Fig. 4. Aoraï's instrumentation

  Create a logical variable for each input. Add constraints for the precondition into the constraint store. Let the initial path prefix π be empty (i.e. the first test case can be any test case satisfying the precondition). Continue to Step (G 2 ). (G 2 ) Symbolically execute the path π: add constraints and update the memory according to the instructions in π. If some constraint fails, continue to Step (G 5 ). Otherwise, continue to Step (G 3 ). (G 3 ) Call the constraint solver to generate a test case t, that is, concrete values for the inputs, satisfying the current constraints. If it fails, go to Step (G 5 ). Otherwise, continue to Step (G 4 ) . (G 4 ) Run a traced execution of the program on the test case t generated in the previous step to obtain the complete execution path. The complete path must start by π. Continue to Step (G 5

  5. Sante first calls Value to analyze p and to generate an alarm for each potentially unsafe statement. Next, Slicing is used to reduce the whole program with respect to one Overview of the Sante method or several alarms. It produces one or several slices p 1 , p 2 , . . . , p n . Then, for each p i , PathCrawler explores program paths and tries to generate test cases confirming the alarms present in p i . If a test case activating an alarm is found, the alarm is confirmed and classified as a bug. If all feasible paths were explored for some slice p i , all unconfirmed alarms in p i are classified safe, i.e. they are in fact false alarms. If PathCrawler was used with a partial criterion (k-path), or stopped by a timeout before finishing the exploration of all paths of p i , Sante cannot conclude and the statuses of unconfirmed alarms in p i remain unknown.

	Program p	Precondition
		Value analysis
		p, Alarms
	Option: all, each,
	min, smart	Program slicing
	p1	p2	. . .	pn	For smart: try smaller slices
					if necessary
		Test generation
		Diagnostic
	Fig. 5.

f in the original C code. This instrumentation is summarized in Fig.4. It consists of two main parts: prototypes whose specification takes care of performing

See http://frama-c.com/support.html.

* This work was partly supported by ANR U3CAT and Veridyc, and FUI9 Hi-Lite projects