Salma Mesmoudi
email: salma.mesmoudi@iscpif.fr

Isabelle Alvarez
email: irsteaisabelle.alvarez@lip6.fr

Sophie Martin Irstea

Mariette Sicard
email: mariette.sicard@grignon.inra.fr

Nathalie Perrot
email: nathalie.perrot@grignon.inra.fr

Coupling geometric analysis and viability theory for system exploration. Application to a living food system

Keywords: model exploration, viability, geometric analysis, food system, optimal control, dynamical system

This paper addresses the issue of studying a food complex system in a reverse engineering manner with the aim of identifying the set of all possible actions that makes it reach a quality target with respect to manufacturing constraints. Once the set of actions is identified, several criteria can be considered to identify interesting trajectories and control policies. A viability approach, coupling the viability theory and a geometric approach of robustness, is proposed to study complex dynamical systems. It can be implemented for several types of systems, from linear to non linear or hybrid systems. The proposed framework was adapted to a living food system: a ripening model of Camembert cheese to identify the set of states and actions (capture basin) from which it is possible to reach a predefined quality target. Within the set of viable trajectories, particular trajectories that improve the Camembert cheese ripening process are identified using the proposed approach. The results are applied at a pilot scale and are discussed in this paper.

Introduction

The cheese ripening process, such as the one used for Camembert, is considered to be a complex system. Numerous interactions take place over time at different scales, from the microscopic to the macroscopic level. To better understand and eventually enhance Camembert ripening control, numerous studies have been carried out in food sciences. But despite the number of experimental databases collected, process knowledge remains incomplete. Experimental trials are very costly (the necessary time is 41 days per trial). However, models have been developed to help us to more effectively understand such complex processes [START_REF] Perrot | Modelling and analysis of complex food systems: state of the art and new trends[END_REF].Cheese processing has been modeled by means of mechanistic models [START_REF] Riahi | Model describing Debaryomyces hansenii growth and substrate consumption during a smear soft cheese deacidification and ripening[END_REF], the partial least square method [START_REF] Cabezas | Application of Partial Least Squares regression to predict sensory attributes of artisanal and industrial Manchego cheeses[END_REF], neuronal methods [START_REF] Jimenez-Marquez | Impact of modeling parameters on the prediction of cheese moisture using neural networks[END_REF],dynamic Bayesian networks [START_REF] Baudrit | A dynamic bayesian network to represent a ripening process of a soft mould cheese[END_REF],genetic algorithms [START_REF] Barriere | Modeling human expertise on a cheese ripening industrial process using gp. Parallel Problem Solving From Nature -Ppsn X[END_REF], stochastic models [START_REF] Aziza | Stochastic, compartmental, and dynamic modeling of cross-contamination during mechanical smearing of cheeses[END_REF], finite element methods [START_REF] Bona | Multicomponent diffusion modeling and simulation in prato cheese salting using brine at rest: The finite element method approach[END_REF] and the fuzzy symbolic approach [START_REF] Perrot | Fuzzy concepts applied to food product quality control: a review[END_REF][START_REF] Ioannou | Experimental analysis of sensory measurement imperfection impact for a cheese ripening fuzzy model[END_REF].Cheese making has been modeled by means of microorganism kinetics, contamination evolution, substrate consumption, mineral diffusion, sensory property prediction, ripening time prediction and expert knowledge. These models may become a key source to integrate the knowledge from experimental databases. Simulations can be performed with these models to investigate and to better understand food processes. The aim of this study is to use a viability approach, developed within a complex system approach, to more effectively explore the Camembert ripening process. It relies on an approach coupling the viability theory and a geometric analysis [START_REF] Alvarez | Describing the Result of a Classifier to the End-User: Geometric-based Sensitivity[END_REF]. For basic information on set valued analysis and viability theory, we refer for instance to [START_REF] Aubin | Set-Valued Analysis[END_REF] or [START_REF] Aubin | Viability theory[END_REF] and [START_REF] Aubin | Viability Theory[END_REF]. The viability theory aims at controlling dynamical complex systems with the goal of maintaining them within a given constrained set. Such problems are frequently encountered in ecology or economics where the systems badly deteriorate when they leave some regions of the state space. This theory has been applied to ecological problems such as prey-predator dynamics studied by [START_REF] Bonneuil | Viable populations in a prey-predator system[END_REF] to determine the necessary conditions to allow prey and predator coexistence. It was also applied to the renewable resource domain, for example, to the viability of trophic interactions in a marine ecosystem by [START_REF] Chapel | Defining yield policies in a viability approach[END_REF] or to the restoration cost of an eutrophic lake by [START_REF] Martin | The cost of restoration as a way of defining resilience: a viability approach applied to a model of lake eutrophication[END_REF] considered as socio-ecological systems. Other applications can also be found in the areas of finance [START_REF] Bonneuil | The hybrid guaranteed capture basin algorithm in economics[END_REF], highway traffic fluxes [START_REF] Aubin | A viability approach to hamilton-jacobi equations:application to concave highway traffic flux functions[END_REF] and sociology [START_REF] Bonneuil | Viability in dynamic social networks[END_REF]. This is the first time that viability theory has been applied to a living food system. Moreover the coupling of such an approach with a geometric study, brings original tools to deal with decision help and reverse engineering problems in food science.

Viability theory proposes a framework to calculate the set of all possible states (the capture basin) from which it is theoretically possible to reach a quality target with respect to manufacturing constraints [START_REF] Aubin | Viability kernels and capture basins of sets under differential inclusions[END_REF]. The ripening process is modeled as a dynamical system. This work is detailed in [START_REF] Sicard | Methods, concepts and tools of system complex to control food processes. Camembert cheese ripening application[END_REF] and [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF]. Contrary to optimal control,the viability approach computes the set of all possible states without considering optimization aims at first. The main objective of this method is to keep new and unexpected trajectories that would be disregarded otherwise. The capture basin can be seen as an explicit representation of the process itself.

Once the capture basin has been computed, it is interesting to study its geometric properties. Outside the capture basin, it is not possible to control the ripening process to make an acceptable cheese. So, if a cheese is in a state near the boundary, small perturbations (for example, dysfunction of sensors or actuators, measurements errors, control errors, etc. ) can throw the processed cheese out of the capture basin, where it cannot recover. Obviously, trajectories that stay far from the capture basin boundary are safer than trajectories that go near the boundary. The distance to the boundary can thus be seen as a measure of the robustness of a state to perturbations or errors.

Geometric analysis of the capture basin allows to define the robustness of states and trajectories. Robustness can then be used to qualify particular control sequences. It can also be used explicitly as a selection criterion (as well as control cost and other classical optimization criteria). The global distance analysis can also be used to identify unsafe areas in the process.

This paper presents the geometric study of the capture basin applied to a food living system: the Camembert cheese ripening process, based on the viability study made in [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF].

Section 2 introduces the viability framework. Presenting in detail the principles and the algorithm of the geometric method, and the definition of robustness. Section 3 presents the results that were obtained with the Camembert cheese ripening application following by conclusions in section 4.

Material and Methods

2.1

Capture basin of the Camembert cheese ripening process.

The domain of viability theory is the control of the dynamical systems [START_REF] Aubin | Viability theory[END_REF].

The evolution of the state variable vector x ∈ X ⊂ R n is described by:

x (t) = f (x(t), u(t)) u(t) ∈ U (x(t)) (1)
where U (x) ⊂ R p is the set of admissible controls when the state of the system is x. Given a constraint set K ⊂ X, the viability theory methods and tools first aim at determining the viability kernel, that is the subset of K gathering all states of the system such that there exists at least one control function that makes it possible to remain in the constraint set indefinitely:

V iab(K) = {x ∈ K | ∃u(.) | x(t) ∈ K ∀t ∈ [0, +∞[}.
The states belonging to the viability kernel are called viable states.

Given a constraint set K and a target set C, another concept is the capture basin define as all the states for which a control function exists and makes it possible to reach the target set while remaining in the constraint set:

Capt(K, C) = {x ∈ K | ∃u(.), ∃T | x(T ) ∈ C and x(t) ∈ K ∀t ∈ [0, T ]}.
In the case of the Camembert cheese ripening application [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF], the state space is 3-dimensional: m, the cheese mass, T S , the cheese temperature and the microorganism respiration rate, r O 2 ,CO 2 . The cheese mass loss is linked to the evaporation phenomena and the carbon loss through respiration of microorganisms. The control dynamics for m and T S are two control variables, the ripening room temperature T ∞ (Kelvin) and the relative humidity rh [START_REF] Hélias | Modeling of Camemberttype cheese mass loss in a ripening chamber: Main biological and physical phenomena[END_REF]:

dm dt = s {w o 2 • r o 2 -w co 2 • r co 2 -k [a w • p sv (T s ) -rh • p sv (T ∞ )]} (2)
dT S dt = s m • C h(T ∞ -T s ) + εσ(T 4 ∞ -T 4 s ) -λk [a w • p sv (T s ) -rh • p sv (T ∞ )] +α r o 2 + r co 2 2 ( 3
)
a w is the cheese surface water activity (dimensionless), p sv is the saturation vapor pressure (P a), k is the average water transfer coefficient (kg.m -2 .P a -1 .s -1 ), C is the cheese specific heat (J.kg -1 .K -1 ), h is the average convective heat transfer coefficient (W.m -2 .K -1 ), ε is the cheese emissivity (dimensionless), σ is the Stefan-Boltzmann value (W.m -2 .K -4 ), α is the respiration heat for 1 mol of carbon dioxide release (J.mol -1 ) and λ is the latent vaporization heat of water (J.kg -1 ). This set of previously established differential equations represents a simple but accurate model to predict cheese mass changes according to available online measurements. The main hypotheses underlying the model are 1) the cheese water activity is constant during ripening, 2) the respiratory activity of the microflora plays a major role by inducing heat production, combined with important water evaporation, 3) the temperature gradient inside the cheese is negligible, and the limiting phenomenon is the convective transfer. The water activity and the specific heat of the cheeses are assessed by offline measurements. The others parameters in the model are obtained from literature.

The controlled dynamics of the microorganism respiration, r O 2 ,CO 2 , are extracted from expert knowledge and experimental results [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF].

We fix the following constraints on the state variables: Cheese mass can vary from 250g to 310g, the cheese temperature from 7 • C to 17 • C and the respiration rate of microorganisms from 0 to 50 g/m 2 /day.

The sets of admissible controls are ripening room temperature from 8 • C to 16 • C and relative humidity from 84% to 98%.

The considered manufacturing and quality objectives regarding the sets of admissible controls are described below:

• The final mass of the cheese, must belong must range from 250to270g

(the lower limit is fixed by administrative rule),

• The final temperature must remain between 8 • -10 • C

• The final respiration must remain between 23 and 50g/m 2 /day

• A specific profile for the gas rate evolution of (CO 2 , O 2 ) is fixed. The respiration rates should begin at level 0 on day 1 (microbial growth latency), reach a maximum between day 3 and day 8 and decrease slowly during the last days of ripening. This temporal constraint is chosen to embed a quality requirement. (see [35][36]).

These four objectives are described by subsets of the state space, also called targets for the viability theory framework [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF].

The discretization of the control variables is set upon a balance between (a) the sensitivity of the microorganisms in the cheese to an increment of control variables (known to be significant for an increment of 1 • C and 1% of relative humidity [START_REF] Leclercq-Perlat | Controlled production of Camembert-type cheeses. part i: Microbiological and physicochemical evolutions[END_REF], [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF] ) and (b) the precision of the regulation. In ripening rooms, as well as cold chambers [START_REF] Baudrit | A Joint treatment of imprecision and variability in food engineering: Application to cheese mass loss during ripening[END_REF], spatial variations of humidity and temperature are always observed due to climate control. The resulting gradients deriving from the shape of the room and air regulation devices. Place from an industrial point of view, since it is impossible to place sensors in every point of the ripening chamber, the precision of the regulation is estimated to be around 2%. For the state variables, 1 g in mass loss and 1m -2 .s -1 in respiration rate are fixed. It corresponds to the level of attempted precision described by the experts of the domain for those variables [START_REF] Leclercq-Perlat | Controlled production of Camembert-type cheeses. part i: Microbiological and physicochemical evolutions[END_REF], [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF].

Given those dynamics, constraints and targets, the viability theory algorithm [START_REF] Saint-Pierre | Approximation of the viability kernel[END_REF] allows to compute the sequence of T -capture basins, considering constraints depending on time, with backward computation from the target followed by a forward computation from initial conditions (see [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF] for more details).

A T -capture basin gathers all states for which there exists a control function that allows to reach the target at time T while remaining in the constraint set between time 0 and T :

Capt(T, K, C) = {x ∈ K | ∃u(.), x(T ) ∈ C and x(t) ∈ K ∀t ∈ [0, T ]}.
Calculation requires a nearly exhaustive search in the control space at each time step. In this study, the computation is performed by the distributing the effort on a cluster of 200 CPU (Central Processing Unit). A dimension 5 (plus time) set is produced, the capture basin with the three state variables (cheese mass, respiration rate and cheese temperature) and the two associated control variables (ripening room temperature and relative humidity of the air) [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF]. Actually an exhaustive search is made to record all the control values (with respect to the discretization scheme), that lead from one state in

Capt(T -t, K, C) to one state in Capt(T -t + dt, K, C), with t ∈ [0; T -dt].
In the 12 days-capture basin, the traditional Camembert ripening process strategy (T s = 12°C, rh = 92%) gathers a subset of states (which are all reached using these control values, with exception of initial states, and from which these controls leads to states that are also in the capture basin). All the others states correspond to different possible strategies which eventually lead to the targeted cheese quality.

Once the capture basin is computed, it is possible to perform a geometric analysis in order to identify risky areas. It is also possible to explore it in order to find new control sequences.

Geometric analysis of a capture basin.

The state space is divided into two areas within the capture bassin: belongs or not to the capture basin. If usual existence assumptions are done on the dynamics, the constraints and the target (compactness, convexity, continuity, for more details see [START_REF] Aubin | Viability Theory[END_REF]), then the corresponding capture basin is a closed set. We will assume in the following that these assumptions are satisfied and that the capture basin is a closed set. We denote Γ the boundary of the T -capture basin.

In this way, a T -capture basin can be seen as a classification system where the decision boundary is the boundary of the T -capture basin. We thus propose to use technics belonging to the field of classification in order to analyse the capture basin. In particular, the distance to the decision boundary and other geometric concepts can be used to obtain relevant information about a particular state [START_REF] Alvarez | Describing the Result of a Classifier to the End-User: Geometric-based Sensitivity[END_REF].

Geometric robustness

We assume in this section that it is possible to define a meaningful distance in the state space (which dimension is finite), see [START_REF] Alvarez | Viability and resilience of complex systems: concepts, methods and case studies from ecology and society[END_REF] for more details. In practice, when performing a viability study on a dynamical system, homogenization and rescaling are part of the state variables discretization of in order to run the viability algorithms with efficiency. In the case of the Camembert cheese ripening, the discretization step of state and control variables is de-termined by domain experts (see [START_REF] Sicard | A viability approach to control food processes: Application to a Camembert cheese ripening process[END_REF]), taking into account sensor accuracy. Discretization steps are considered as distance unit on each variable axis. Euclidean distance or sup-norm distance are then defined accordingly. With the sup-norm, the combination of perturbations of the same size on different variables does not change the size of the resulting perturbation.

In this work, we define perturbations one-time instantaneous "shocks" in the state space. Consequently, the occurrence of a perturbation causes a jump in the state space and can be described by a function D that associates state x with state D(x) the reachable state from x after this perturbation [START_REF] Martin | The cost of restoration as a way of defining resilience: a viability approach applied to a model of lake eutrophication[END_REF]. The size of a perturbation is then defined by the size of the jump, the distance d(D(x) -x). A perturbation can be, for instance, an accidental loss of mass, or a sudden change of respiration rate of microorganism (due to external influences such as change of activity or mortality). Former errors of measurement can also be seen as perturbations (once they are discovered).

Considering such perturbations, the distance to the capture basin boundary Γ is a useful indicator of the robustness of a state. If a perturbation occurs when the state of the system is x, as long as its size is smaller than d(x, Γ), the new state of the system, D(x), will still belong to the capture basin. Consequently, the distance to the boundary can be seen as a measure of the robustness of a state to perturbations (for instance errors of measurement concerning the localization of the state variable) in the state space. This robustness indicator could be used to assist the process operator when the process is sufficiently monitored: If the state of the system has a low robustness value, it means that it is close to the boundary, so monitoring should be intensified.

We define p(x) as a point of the boundary for which the distance d(x, Γ) is reached (for the Euclidean distance it is the orthogonal projection). p(x) gives the direction and size of the smallest perturbation (the sensitive move) at state x that will move the system outside the capture basin.

State robustness can be used to define a robustness of the trajectory followed by the system in the capture basin. Several definitions can be proposed, depending on the manager's perception of risks [START_REF] Plous | The Psychology of Judgment and Decision Making[END_REF]. The risk of adverse decision can be considered in many different ways. We will consider the most risk-adverse indicator which is the minimum of the distance to the boundary over the trajectory: Definition 1 Min-robustness of a Trajectory.

Let T > 0, x ∈ Capt(T, K, C) and x : t -→ x(t) be a trajectory starting at

x such that ∀t ∈ [0; T ], x(t) ∈ Capt(T -t, K, C), we note Γ T -t the boundary of Capt(T -t, K, C
), then the robustness value of Min-robustness of x(.) is defined by:

r(x(t)) := min{d(x(t), Γ T -t )|t ∈ [0; T [}.
With this definition it is possible to assign a robustness value to each cheese ripening trajectory that reaches the target set while satisfying the constraints. Trajectories that stay always far from the T -capture basin boundaries have a higher Min-robustness value, as it can be seen in figure 1. In the following sections we use the term robustness for Min-robustness when dealing with trajectories.

[Figure 1 around here] When a trajectory has a low Min-robustness value, there are times for which x(t) is close to the boundary of Capt(T -t, K, C) which means that a perturbation (or a measurement error in the state space) with a low intensity may cause a jump outside Capt(T -t, K, C) and the target will no more be reached in the allotted time. Hence, trajectory Min-robustness can be seen as a criterion to choose a particular trajectory among the set of trajectories that reach the target while remaining in the constraint set.

Global quality indicator of a capture basin

Geometric indicators give useful information about the feasibility of controlling the system under the constraints set.

Relative volume. The size of the capture basin compared with the size of the constraint set is a first basic indicator. Since the capture basin V is a subset of the constraint set K, it is easy to compare the size of the corresponding hypervolumes: vol(V ) vol(K) . This indicator is generally used in sensitivity analysis. However, this indicator can gives misleading information about the structure of the problem, even when the capture basin is a simply connected set (which means roughly that it has no holes). because the same volume can enclose rather different shapes. For instance, in a dimension d space, a hypercube with side 2r and a hyperrectangle with one side size 2αr and all other sides 2r (α)

1 d-1 with α >> 1,
share the same volume (2r) d . In the first case, the distance to the boundary of the hypercube reaches its maximum at the center with a value of r. In the second case, the points of the hyperrectangle cannot be farther from the boundary than r (α)

1 d-1
.

Maximal maximal ball.

Maximal balls are open balls included in the capture basin, that are maximal. This means that it is not possible to find a larger open ball centered on the same state that will be totally included in the capture basin. Centers of maximal balls form the skeleton. The maximal maximal ball is the maximal ball with the largest radius among maximal balls. Its center is the farthest point of the capture basin from the decision boundary. The larger the radius of the maximal maximal ball is, the more robust to perturbation the system is: among all the closed sets with the same inner volume, the maximum value of the distance to the boundary is reached for the closed ball.

A relatively large maximal maximal ball suggests that it should be easy to find control sequences that guarantee viability constraints verified over time.

Distance and projection algorithm

When the decision boundary Γ is described by an analytical formula, or by a set of constraints, it is possible to compute the exact value of d(x, Γ) for all points of the input space, with efficient algorithms. When it is not the case, then it is necessary to compute an approximation of the distance. Considering a grid G of N points per axis in dimension n, Γ N is the set of points that belong to the discretized boundary. When Γ N is known, efficient algorithms make it possible to compute an approximation of the distance. In practice, Γ N is often not known and what is computed is the distance to the complementary set (points that are not in the capture basin). An algorithm from [START_REF] Meijster | A General Algorithm for Computing Distance Transforms in Linear Time[END_REF] has been adapted in [START_REF] Alvarez | Describing the Result of a Classifier to the End-User: Geometric-based Sensitivity[END_REF] to compute on a grid of N points per axis the exact distance to the complementary set. This algorithm is optimal since its complexity is in O(N n ), that is a linear complexity in the number of points of the grid where the distance is computed.

V = ∅ is the set of points that are not viable. The distance and projection algorithm DistanceMapToSet consists in computing first the distance and a point where the distance is reached along the first axis. The first axis is run in both directions, a function labels points of V if any point is discovered, and a counter records the distance to the last encountered point. This procedure is linear in the number of grid points, actually 2N n .

Then a second procedure AdditionalAxis, is called repeatedly for the n-1 remaining axes. Its complexity is also in O(N n ), so the complete algorithm is in O(n • N n ). The basic idea, coming from [START_REF] Hirata | A Unified Linear-Time Algorithm for Computing Distance Maps[END_REF], and detailed in [START_REF] Meijster | A General Algorithm for Computing Distance Transforms in Linear Time[END_REF] for dimension 2, consists in considering at step k the distance map g k-1 (X) computed at step k -1 on G (the square distance in the Euclidean case). Building functions F X specific to the distance are considered at each coor-dinate x k along the kth axis, based on the map g k-1 (X) computed on the first k -1 axes and the distance along axis k. When the difference between two building functions is a monotonic function, then the lower envelope of these functions gives the distance at step k. This monotonicity property is verified in particular in the case of the Euclidean square distance and of the sup-norm distance (we note it D depending on the distance d: D = d 2 in the Euclidean case and D = d in the sup-norm case).

In the case of the Euclidean distance, g n (X) = d 2 (X, p n (X)). In the case of the sup-norm, g n (X) = d(X, p n (X)). (See Appendix for details).

Theoretical example

We propose here a simple model to illustrate the interest of computing the capture basin and then use the distance to the boundary to select a trajectory robust to perturbations in the state space.

Description of the viability problem

We use in this example a controlled dynamical system, the state space is two-dimensional (x, y) and the control space is one -dimensional u:

x (t) = x 0 x 0 > 0 y (t) = y 0 + u(t) y 0 > 0 u(t) ∈ [-u 0 ; +u 0 ] u 0 > 0 (4)
The constraint set is described by:

K := {(x, y) | y ∈ [y 1,K ; y 2,K ]} (5)
The target is described by:

C := {(x, y) | x = x 1 and y ∈ [y 1,C ; y 2,C ]} (6)
The dynamics, the constraint set and the target are displayed in figure 2a.

[Figure 2 around here]

Calculation of the capture basin

In this example, the capture basin is easy to determine as shown in figure 2b.

Outside the capture basin Capt(K, C), it is not possible to reach the target. For states (x, y) which belong to the capture basin (as a topological set), every control value u ∈ [-u 0 ; +u 0 ] applies.

Finding robust evolutions

We consider perturbations in the state space, such as in 2.2.1: a state (x, y) can be translated to D(x, y) = (x , y ). The state space is associated with the Euclidean norm. The distance to the boundary of the capture basin as then an index of robustness to these perturbations.

Given an initial point A and knowing the capture basin, we can compute the evolution which is viable and reaches the target, and which also maximizes the integral of the distance to the boundary of the T-capture basins. Such an evolution is displayed in figure 2c.

It is worth noting that the computation of the capture basin should be computed before the determination of this robust evolution.

Indeed, maximizing the distance to the boundary of the constraint set leads to choose a much less robust evolution displayed in figure 2d. This evolution spends a long time on the boundary of the capture basin and the smallest perturbation can provoke a jump of the system state outside the capture basin and the target can no longer be reached.

Results: Geometric analysis of the ripening process

The proposed viability algorithm was used to compute the capture basin of the cheese ripening process. We then applied the previously described distance and projection algorithm to this capture basin. The main objectives are to offer a visualization of the ripening process as a set of possible states, to provide global indicators for the capture basin, to study the robustness of theoretical trajectories in order to select the best experimental trajectories.

Complete capture basin robustness

The complete calculated capture basin is organised in 11th days completed by the capture basin of the 12th day which is the target itself. The dimension of the state space is 3 (see section 2.1).

[Figure 3 around here]

Figure 3 shows the twelve sections (one section per day) of the complete capture basin. The Euclidean square distance (d) to the boundary is used to color the points of each section.

The more the color goes towards red, the more robust the state is (far from boundary), on the contrary, the more the color goes towards blue, the less robust the state is (closer to the boundary).

The shape of the capture basin suggests that the process is rather robust to perturbations in the state space: the daily sections have relatively large maximal maximal balls, with a square radius of the maximal maximal ball of 25 steps which corresponds to a distance d of 5 steps. This is also the theoretical maximum of distance to the boundary of the constrained domain. This distance can be expressed mathematically, as half the smallest constraint range of the state variables which corresponds to the cheese temperature range, 10 2 where 10 represents the number of steps between 17 • C and 7 • C . The maximum value ( 5) is reached for an important part of the viable states inside each section from day 4 to day 12, as it can be seen in Figure 3, the red area of these daily sections covers a significant proportion of the capture basin.

The different sections of the capture basin show the evolution of the possible states of the cheese from one day to the next.

The viability set of day 1 is reduced to a two-dimension set, since the breathing dimension is not active, the microorganisms involved in the process having not germinated yet. All the mass values are viable at the starting day. From the 2nd day until the 6th day (see Figure 3), the microorganisms breathing develops and the viability domain spreads over the constraint set. The percentage of the viable states increases rapidly (see Figure 4) to 89 %.

[Figure 4 around here] This shows that the exploration domain is correctly defined, since the capture basin is strictly inside the exploration domain but one section occupies almost all the domain.

During the 5th and 6th day, the proportion of viable states is maximal, but the capture basin drifts progressively to higher values of the breathing and higher surface temperature. This can be seen easily in Figure 3. After Day 8 it is no longer possible to reach the target if the microorganisms haven't reached a minimal development. From the 6th day, the viability domain begins converging towards the target. This convergence implies a decrease in the volume of the viable state area, with the fading of high mass values and a slight decrease of the values of the breathing of the microorganisms. The 11th day is very constrained: the percentage of the viable points reaches its lowest level, 17 %, just before the target itself (6 %).

Even in the smaller capture basins, the maximum distance to the boundary remains almost unchanged. This characteristic of the shape of the capture basin suggests that the ripening process is particularly robust to perturbations.

Contribution to a reverse engineering approach

The capture basin provides the set of all states from which it is possible to reach the target with respect to the constraints. For each state x in the (T -t)-capture basin the viability analysis gives also at least one vector u of control values such that f (x, u) is in the (T -t + dt)-capture basin. Actually, in the continuous case, for states that are not on the boundary, all the control values are available. However, some control values may be used only for a very short time. With a constant time step discretization, all the control values may not be available anymore. In our case, with respect to the discretization scheme, we compute all the possible control values that can be associated to a state x. In this enlarged (T -t)-capture basin (state space x with possible controls), we used the same definition. In practice, this means that we prefer a state far from the capture basin boundary with control values far from the boundary of available controls.

In accordance with experts, we focus our study on a viability tube computed in less than twelve days, since the analysis of the capture basin showed that it was possible to reach the target in less time. A gain in time of 4 days was considered by the experts as extremely desirable. This is the reason why we have worked on a more constrained set than the complete 12 day tube. We have calculated a new 8 day capture basin and all the associated trajectories.

Inside this capture basin, a trajectory that reaches a weighted compromise between reduced initial mass, reduced operational costs and robustness is selected.

A geometric calculus is applied to this trajectory. On this basis, we show, through simulations and experimental validations on a pilot, that it is possible to predict and quantify the disturbances that could lead to disrupted trajectories of the capture basin, or conversely the perturbations that ensure future states to stay within the boundaries.

added value of the geometric robustness and simulations

The result of the geometric calculus on the 8 day trajectories is presented Table 1.

The geometric robustness of the trajectory is computed with the DistanceMapToSet algorithm with the sup norm in an extended state space. With regards to the identical results for T and T s , we only use one of the temperatures, labeled T , for the computation of the distance. This equality can be easily explained by the time to reach T s knowing the control temperatures. It can be considered as instantaneous as regard to the process reaction time manipulated (1 day). The extended state space is then the 4-dimensions space of the capture basin with normalized state and control variables. The sup norm is used instead of the Euclidean norm, since measurement and control errors can combine but not substitute for one another. A distance (d) of 3, like in day 3, means that one or more variables are 3 steps away from the boundary. Any 2 steps perturbation is not enough to move the cheese state outside.

The distance d can be compared to the maximal maximal ball radius of the capture basin(mmbr) for each day. At day 3 for example, a value of 4 means that the thickness of the capture basin at this day (radius of the maximal maximal ball) is 4, that is at least four steps along one variable at least are needed to reach the boundary of the capture basin from its center.

In Table 1 we can see that the 8-days trajectory is very robust, since each day is almost as large as possible (the thickness of the tube is 4).

To assess the validity of our approach, we subsequently test the predicted 8-day trajectory under perturbations.

Since the accuracy of common actuators concerning the control of the relative humidity is rather limited, we decided to apply the disturbance to this variable. Small errors (1%) can induce a great change in the water transfer coefficient and as a consequence in the mass loss [START_REF] Baudrit | A Joint treatment of imprecision and variability in food engineering: Application to cheese mass loss during ripening[END_REF]. In our model, the discretization step on relative humidity is 2%, to take into account the imprecision of sensors and the difficulty of control.

Since the robustness of the trajectory shown in Table 1 is 3 after the two first days, we apply a perturbation to the relative humidity on day 3 and 4, to be sure to be near the boundary without leaving the viability domain. Instead of the normal value of 94%, we fix a set point of 96% for the relative humidity during these two days. The temperature is unchanged. According to the theoritical results, this disturbance should not bring the trajectory out of the capture basin. Table 2 shows the disrupted trajectory simulated by the model. As expected, the robustness of the disrupted trajectory is now only 2 on day 3 and 4 (it is only two steps from the boundary). The change of state values induced by the change of control values has no other effect on the robustness of the trajectory. 5 around here]

The 8-day trajectory also performs very well in the 8-day capture basin. Moreover, the proposed approach also helps us to know to which extent this new way of control is robust, and what is the acceptable uncertainty with regards to the boundary of the capture basin. This is very encouraging for further developments of control strategies at 8-day ripening rather than the standard of 12-day process.

[Figure 6 around here] [Figure 7 around here]

Conclusions

Thanks to the viability theory framework applied to the ripening process, we are able to compute the set of initial points included in a given constrained set from which starts at least one evolution (1) satisfying the manufacturing constraints and (2) reaching the quality target at the allotted time. We evaluate the robustness of these trajectories by using recent original mathematical developments, finally choosing a trajectory with low operational costs among the more robust ones. This trajectory has a 8-day ripening time and an initial mass of 0.284 kg, whereas the standard one was ripened in 12 days from an initial mass of 0.3 kg. Moreover, we show that a coupling between the viability method and an adapted geometrical robustness computation leads to interesting information for reverse engineering purposes. We show that it is possible to predict and quantify the disturbances that can lead to disrupted trajectories that go outside the viable set and conversely that can ensures the state to stay in the viable set. This allows us to effectively control the cheese ripening process in an optimized way. The method developed in this work can be applied to other processes for which a dynamical model is available, in order to explore its state space and propose new ways of controlling it.

Robustness is crucial for the dairy industry to avoid process drift. Other measures of robustness based on the capture basin boundary should be developed to identify the areas where the process is more sensitive to perturbation and, consequently, where it should be more carefully monitored.

In future works, we will investigate in more detail the enlarged state space with viable controls, and the links with the time of exit from the capture basin (since some control values can apply for very short time if a state is near the capture basin boundary).

Further studies will be focus on the exploitation of this approach to consider a multi-objective optimization problem related to the modeling of a food industrial process, that is the Camembert cheese ripening process. One aim is to construct the Pareto boundary [START_REF] Pareto | Manuel d'économie politique[END_REF] of the optimal strategies, i.e. the optimal paths, among each possible viable path to the target. defined by:

F X (i) = F X k x k (i) := g k-1 (X) + (i -x k ) 2 (0 ≤ x k < N ) intersect(F X k x k , F X k y k ) = truncate(x 2 k -y 2 k + g k-1 (X) -g k-1 (Y )) ÷ 2(x k -y k ) (9)
In the case of the sup norm, the building functions ( 8) are truncated V-shaped function defined by:

F X (i) = F X k x k (i) := max(g k-1 (X), |i -x k |), (0 ≤ x k < N )
The intersection of truncated V -shaped functions F X and F Y used by the algorithm is given by the following formula (with x k ≤ y k ), contrary to [START_REF] Meijster | A General Algorithm for Computing Distance Transforms in Linear Time[END_REF]:

intersect(F X , F Y ) =      max((x k + g k-1 (Y )), ((x k + y k ) ÷ 2)) if g k-1 (X) ≤ g k-1 (Y ) min((y k -g k-1 (X)), ((x k + y k ) ÷ 2)) otherwise (10)
The algorithm can be used for other distances, as long as the difference of building functions verifies the monotonicity property (see [START_REF] Hirata | A Unified Linear-Time Algorithm for Computing Distance Maps[END_REF]).

We now consider the kth axis u k , 2 ≤ k ≤ n. F X k l (x k ), with 0 ≤ l < N , gives the distance (or square distance in the Euclidean case) between X = (x 1 , .x k , .x n ) and p X k k-1 (l). We have

F X k l (x k ) = F (X+(l-x k )u k ) (x k ). If V ∩ (X + E k ) = ∅,
then there is at least one l, 0 ≤ l < N such that g X k k-1 (l) and p X k k-1 (l) are defined. Then we have:

   g k (X) = min F (X+(l-x k )u k ) (x k ), 0 ≤ l < N p k (X) = p k-1 (X + (m -x k )u k ) with m = argmin F (X+(l-x k )u k ) (x k ), 0 ≤ l < N (11)
But the computation of all these values of F is suboptimal, so the algorithm considers the set of building functions F X k l , and computes the lower envelope of this set, using the intersection between two building functions F to switch from one function to another.

For all points X k ∈ [0, N -1] n-1 , the procedure AdditionalAxis recruits building functions along axis k. A candidate function f = F X k i is recruited if there is no other recruited function or if it intersects the previous recruited function F X k j with j < i on the grid. From the definition of the building functions, we then have: Let 0 ≤ j < i < N and w = max(0; 1 + intersect(F X k j , F X k i )), (w < N ) ⇒ ∀l, w ≤ l < N, F X k i (l) < F X k j (l) (12)

The candidate function is then recruited after F X k j from the coordinate w if w < N . If there is no previous recruited function it is recruited from w = 0.

Let us suppose that there are already two building functions F X k x 1 and F X k

x 2 with x 1 < x 2 < i. F X k x 2 was recruited from w(x 2 ) with:

0 ≤ w(x 2 ) = 1 + intersect(F X k x 1 , F X k x 2 ) < N
. Then if f (w(x 2 )) = F X k i (w(x 2 )) < F X k x 2 (w(x 2 )), the monotonicity property insures that for all l ≥ w(x 2 ), then f (l) < F X k

x 2 (l) (This can be easily verified from equations ( 9) and ( 10) respectively). So the candidate f replaces the last function F X k

x 2 . This is done until all previously recruited functions are replaced or alternatively until the last recruited function F X k j , recruited from w j with j < i is such that F X k j (w j ) < f (w j ). In that case f is recruited with w as in [START_REF] Baudrit | A dynamic bayesian network to represent a ripening process of a soft mould cheese[END_REF].

Let (F X k x i ) i∈I be the sequence of building functions recruited along axis k from (w i ) i∈I respectively. The functions are labeled such that if i < j then x i < x j . The sequence of recruitment points (w i ) defines the lower envelope of the recruited building function. We then have:

∀x k , 0 ≤ x k < N, let m = max {i ∈ I, w i ≤ x k } then g k (X) = g X k k (x k ) = F X k xm (x k ) p k (X) = p X k k-1 (x m ) (13)
In the case of the Euclidean distance, g n (X) = d 2 (X, p n (X)). In the case of the sup-norm, g n (X) = d(X, p n (X)). p n (X) is a point where the distance is reached (the orthogonal projection in the Euclidean case).

List of figures: Figure 1: Geometric robustness of trajectories in the capture basin. The robustness of x(t) is low whereas the robustness of y(t) is high.

Figure 2: Construction of a capture basin. , (a) dynamics (4), constraint set (5) and target set ( 6), (b) Capture basin, (c) From initial point A, the blue evolution, viable and which reaches the target, also maximises the integral of the distance to the boundary of the capture basin. The green curve describes the skeleton of the capture basin, (d) From initial point A, the red evolution, viable and which reaches the target, also maximises the integral of the distance to the boundary of the constraint set. The black dotted curve describes the skeleton of the constraint set.

Figure 3: Distance square map of the capture basin (12 sections corresponding to the different days of the process). Viable states that are near the boundary are not shown (the yellow bounding box shows the maximal extent of viable states). Axes x(red), y(green) and z(blue) correspond to the mass(g), temperature ( • C) and the respiration rate of microorganisms (g/m 2 /day) respectively. The black box delimits the target: final mass between 251 and 276g, final temperature between 8 and 10 • C, respiration rate between 34 and 55g/m 2 /day. Figure 7: comparative profile of the respiration rate variation reached for the three experiments TVA, DT and SRT led on the pilot.

Figure 4 :

 4 Ratio of the number of viable states.

Figure 5 :

 5 Evolution of the microorganisms during the cheese ripening for TVA (viable optimized trajectory) and DT (disrupted trajectory), (a) K. marxianus, (b) G. candidum, (c) B. aurantiacum.

Figure 6 :

 6 Comparative profile of the mass loss reached for three experiments led on the experimental pilot TVA (viable optimized trajectory), DT (disrupted trajectory) and SRT (standard ripening trajectory).

Table 1 :

 1 The new 8-day trajectory: control values, predicted state values, robustness (distance d to the boundary of the viability tube, dim 4, sup norm), maximal distance to the boundary (mmbr).

		Control value Ripening state		Distance
	Day T • C	rh%	RR	ML	T s	• C d(mmbr)
	1	12	84	0	284	12	1 (-)
	2	13	94	0	280	13	1 (2)
	3	14	94	1.7	279	14	3 (4)
	4	14	94	6.4	277	14	3 (4)
	5	12	94	17.6	276	12	3 (4)
	6	12	94	32.1	271	12	3 (4)
	7	9	94	44	270	9	3 (4)
	8	9	94	43.2	268		-	-

Table 2 :

 2 8-days disrupted trajectory: control values, predicted state values, predicted robustness (d), maximal distance to the boundary (mmbr).

		Control value Ripening state Distance
	Day T •	RH%	RR	ML	d(mmbr)
	1	12	84	0	284	1 (1)
	2	13	94	2	280	1 (4)
	3	14	96	7	279	2 (4)
	4	14	96	19	278	2 (4)
	5	12	94	33	277	3 (4)
	6	12	94	48	275	3 (4)
	7	9	94	41	272	3 (4)
	Validation of the simulations on the experimental pilot ripening
	room
	Three trajectories are tested on the INRA's experimental pilot ripening room
	[36]. The two first are organized in 8 days while the third one is organized
	in 12 days for a ripening trial close to factory practices. Control and state
	variables are recorded: temperature (T • ) and relative humidity (rh%), mass
	loss (ML), respiration rate (RR). The trajectory corresponding to table (1)
	control set points, for a ripening in 8 days, is called TVA (viable optimized
	trajectory). The trajectory corresponding to table (2) control set points,
	for a ripening in 8 days close to the boundary of the capture basin is called
	DT (disrupted trajectory). The 12-day trajectory, corresponding to standard
	controls in a factory for 12 days of ripening is called SRT (standard ripen-
	ing trajectory, for a control temperature set point of 12 • C and a relative

Acknowledgment

We thank Leclercq-Perlat M.N., Lecornue, F., Guillemin H., Savy, M., Picque D. for the experiments, Bourgine, P. for the ideas, Tonda, A. for the rephrasing and English, the French ANR for the grant for the INCALIN project and the funding from the European Community's Seventh Framework Programme (FP7/2007-2013) under the grant agreement n°FP7-222 654-DREAM.

humidity of the air of 92%). This last trajectory is shown to provide a reference for traditional practices. In particular, the robustness of this standard trajectory, as shown in table [START_REF] Alvarez | Describing the Result of a Classifier to the End-User: Geometric-based Sensitivity[END_REF], is almost always at the possible maximum.

The experimental results show that both trajectories reach the designed target. The main impact of the differences of control the slope of mass loss, as it can be seen in figure 6.

The mass loss is smaller for the disturbed trajectory from the 3rd day. Considering the quality of the cheese, it is well known that the range of control of relative humidity, from 96% to 100 especially if this relative humidity is maintained all along the duration of the process, can have a negative impact on the microorganism kinetics and, as a consequence, on the global state of the cheese [START_REF] Leclercq-Perlat | Controlled production of Camembert-type cheeses. part i: Microbiological and physicochemical evolutions[END_REF]. In fact, we can see that the respiration rate is modified at the beginning of the ripening process (see figure 7). The profiles of respiration rate for the disturbed and undisturbed trajectories are rather different during days 3 to 5.The respiration rate reached is maximum 1 day earlier by comparison to the SRT and TVA trajectories. Nevertheless, the consequence on the microorganism kinetics for the complete duration of the process, including the time of ripening after wrapping of the cheese (from day 12 to day 45), is very limited, as it can be seen in figure 5.

Appendix: Geometric Algorithms

Algorithms are adapted from [START_REF] Meijster | A General Algorithm for Computing Distance Transforms in Linear Time[END_REF]. We consider a grid G of N points per axis in dimension d. V is the set of points that are not viable. Since the boundary of V is not known explicitly, this kind of algorithm is appropriate (when the boundary is known, other algorithm are more suitable). We note

The definition of the building functions depends on D and is such that:

The minimum of a building function F X k x k (i) is reached for i = x k . In the case of the square Euclidean distance, the definition of the building functions from equation ( 8) is given by the Pythagorean theorem and the coordinate of their intersection on the kth axis (noted u k ) on the grid is


We call model and note m a triplet composed of a set of instances corresponding to vertices noted V , and two sets of edges noted E a and E r (see figure 2). The first set of edges denotes attributes. It relates instances to scalar values through attribute names. The second one denotes references. It relates instances to each other through reference names. We note m.V , m.E a and m.E r the V , E a and E r components of a given model m. m

V , E a , E r where

   V ⊆ I E a ⊆ V × A × S E r ⊆ V × R × V Figure 2: Model graph
As an illustration, figure 3 shows a model conforming to an Ecore metamodel together with its corresponding graph-based representation. The metadata are willfully out of the scope of the graph-based representation of models.

Migration specification

For a given model m, we call migration specification and note -→ m a quadruplet composed of m and of three sets noted D i , D a and D r . These sets specify the instances, attributes and references of m that are intended to be deleted (see figure 4). We note -→ m .D i , -→ m .D a and -→ m .D r the D i , D a and D r components of a given migration specification -→ m . The migration specification is intended to be used together with the model to produce a migrated model. This migration specification is computed from a Modif specification. A Modif specification makes it possible to state the deletion of a class, an attribute or a reference (corresponding to D i , D a and D r at the model level). See appendix A for more details about the computation of a migration specification from a Modif specification.

   D i ⊆ m.V D a ⊆ m.V × A D r ⊆ m.V × R

Migration

We call migrator and note M the tool aiming at producing a migrated model from a migration specification. It is defined in figure 5. It simply commits the deletion of instances and of features (attributes and references) on the source model to produce a target model. We note m the model of figure 3. It is formally and explicitly defined by its components as follows:

m = M( - → m ) m .V = m.V \ - → m .D i m .E a = m.E a \ {(i, a, s) ∈ I × A × S | i ∈ - → m .D i ∨ (i, a) ∈ - → m .D a } m .E r = m.E r \ {(i, r, i ) ∈ I × R × I | i ∈ - → m .D i ∨ i ∈ - → m .D i ∨ (i, r) ∈ - → m .D r }
m.V = {s 1 , s 2 , s 3 , sa 0 , sa 1 , sa 2 } m.E a = {(s 1 , name, A), (s 2 , name, B), (s 3 , name, C), (s 1 , ini, true), (s 2 , ini, true), (s 3 , ini, true), (sa 0 , name, a 0 ), (sa 1 , name, a 1 ), (sa 2 , name, a 2 )} m.E r = {(s 1 , inside, s 2 ), (s 2 , inside, s 3 ), (s 3 , contexte, s 2 ), (s 2 , context, s 1 ), (s 1 , entry, sa 0 ), (s 2 , entry, sa 2 ), (sa 2 , sin, s 2 ), (sa 0 , sin, s 1 ), (sa 0 , next, sa 1 ), (sa 1 , prev, sa 0 )}

We note -→ m the migration specification of m to a model conforming to the metamodel of figure 6,

i.e. a model where actions have been deleted. This specification is formally and explicitly defined by its components as follows:

- → m .D i = {sa 0 , sa 1 , sa 2 } - → m .D a = ∅ - → m .D r = ∅
Once the migrator M is applied to the migration specification -→ m , we obtain the migrated model m depicted and formally defined by figure 7.

A B C s1 s2 context inside s3 A A name B B name C C name context inside true true ini ini ini m = M( - → m ) m .V = {s 1 , s 2 , s 3 } m .E a = {(s 1 , name, A), (s 2 , name, B), (s 3 , name, C), (s 1 , ini, true), (s 2 , ini, true), (s 3 , ini, true)} m .E r = {(s 1 , inside, s 2 ), (s 2 , inside, s 3 ), (s 3 , contexte, s 2 ), (s 2 , context, s 1 )}
Figure 7: Concrete syntax and semantics of a migrated model without actions

Rewriting tools

In this paper, the focus is put on the reuse of rewriting tools. To begin with, this kind of tools has to be specified within our semantic domain.

We call rewriting tool and note T the tool aiming at producing an output model from an input model and conforming to the same metamodel. T can be taken for an endogenous model transformation. We note T(m) the output model produced by T from an input model m.

So far, we consider T as a black box. As a consequence, the action of T can only be specified by the differences between a given input model and the corresponding output model. Within our semantic domain, the action of T is therefore specified by means of deleted and added graph elements for a given model m. This specification is formally stated in figure 8. As an illustration, we introduce three simple examples of rewriting tools whose scope corresponds to the metamodel of figure 6.

Deleted elements :

   T i D (m) = m.V \ T(m).V (deleted instances) T a D (m) = m.E a \ T(m).E a (deleted attributes) T r D (m) = m.E r \ T(m).E r (deleted references)
Added elements :

   T i A (m) = T(m).V \ m.V (added instances) T a A (m) = T(m).E a \ m.E a (added attributes) T r A (m) = T(m).E r \ m.E r (added references)
Figure 8: Specification of a rewriting tool

Identity

We call identity and note Id the rewriting tool performing no actions. It applies to a statechart conforming to the metamodel of figure 6 and it provides it unchanged as a result. This trivial example will be useful later to underline the benefits of our approach.

We use m the model depicted by figure 7. By definition, Id(m ) = m . The action of Id on m is therefore formally stated as follows:

Deleted elements :

Id i D (m ) = Id a D (m ) = Id r D (m ) = ∅ (no deleted elements) Added elements : Id i A (m ) = Id a A (m ) = Id r A (m ) = ∅ (no added elements)

Copy

We call copy and note Cp the rewriting tool performing a complete copy of an input model. It applies to a statechart conforming to the metamodel of figure 6 and it provides a new model where all the instances are new.

We use m the model depicted by figure 7 and we note m = Cp(m ), the result of the copy. The formal definition of m and it's graphical representation is given by figure 9.

A B C s4 s5 context inside s6 A A name B B name C C name context inside true true ini ini ini m = Cp(m ) m .V = {s 4 , s 5 , s 6 } m .E a = {(
s 4 , name, A), (s 5 , name, B), (s 6 , name, C), (s 4 , ini, true), (s 5 , ini, true), (s 6 , ini, true)} m .E r = {(s 4 , inside, s 5 ), (s 5 , inside, s 6 ), (s 6 , contexte, s 5 ), (s 5 , context, s 4 )} Figure 9: Copy of the model from figure 7 By definition, the action of Cp on m is formally stated as follows:

Deleted elements :

   Cp i D (m ) = m .V Cp a D (m ) = m .E a Cp r D (m ) = m .E r
Added elements :

           Cp i A (m ) = {s 4 , s 5 , s 6 } Cp a
A (m ) = {(s 4 , name, A), (s 5 , name, B), (s 6 , name, C), (s 4 , ini, true), (s 5 , ini, true), (s 6 , ini, true)} Cp r A (m ) = {(s 4 , inside, s 5 ), (s 5 , inside, s 6 ), (s 6 , contexte, s 5 ), (s 5 , context, s 4 )}

Flattener

We call flattener and note Fl the rewriting tool performing a flattening of a hierarchical statechart conforming to the metamodel of figure 6. It provides a new model where the superstates are deleted and potentially replaced by one new state at least. We suppose we have no more information about the algorithm implemented by the flattener.

We use m the model depicted by figure 7 and we note m = Fl(m ), the result of the flattener. The formal definition of m and it's graphical representation is given by figure 10. We can note that only one instance remains in the result. By definition, the action of Fl on m is formally stated as follows:

ABC s3 ABC ABC name true true ini m = Fl(m ) m .V = {s 3 } m .E a = {(s 3 , name, ABC), (s 3 , ini, true)} m .E r = ∅
Deleted elements :

   Fl i D (m ) = {s 1 , s 2 } Fl a D (m ) = m .E a \ {(s 3 , ini, true)} Fl r D (m ) = m .E r
Added elements :

   Fl i A (m ) = ∅ Fl a A (m ) = {(s 3 , name, ABC)} Fl r A (m ) = ∅

Reverse migration

We call reverse migrator and note R the tool aiming at undoing the migration operations that have been applied to a given initial model m. Thus, R applies to a model which is supposed to have been migrated and it produces a new model which corresponds as far as possible to the initial model.

If the migration specification -→ m associated to the initial model m is unknown, then the reverse migration cannot be computed.

If the migration specification is known, and if the migrated model has not been modified, then the reverse migration is already and trivially available:

R(M( - → m )) = m
We focus now on the more challenging problem of an initial model m which has been migrated through a known migration specification -→ m , and which has been later modified by a rewriting tool T. We note m the model produced by the reverse migrator in this case:

R(T(M( - → m ))) = m (reversed migrated model)
We aim at recovering as much as possible the initial model without undoing the action of T.

The first solution to this problem is a well known and rather obvious approach which consists in considering a reverse migration as a normal case of migration. For a given migration operation, if the corresponding reverse operation is available, then it is applied (e.g. renaming of attributes, which is undone by the reverse renaming). If it is not available, then it can leads to a loss of data. This is actually the case on which we focus here, and this is the reason why we only consider deletion. As a result, nothing is supposed to be deleted nor added by this reverse migration:

R(T(M( - → m ))) = M( ------→ T(M( - → m )))
where:

------→ T(M( - → m )) = (T(M( - → m )), ∅, ∅, ∅)
The reverse migration is obtained from this migration specification according to the principles of figure 5. By definition, it leads to identity:

R(T(M( - → m ))) = M( ------→ T(M( - → m ))) = T(M( - → m ))
To illustrate this first solution, we consider again the three examples of rewriting tools that have been introduced in the previous section, namely identity (Id), copy (Cp), and flattener (Fl). They apply to models conforming to the metamodel of figure 6. We aim at performing a reverse migration of the model of figure 7 after its processing by identity, copy and flattener. The result is a model which is supposed to conform to the metamodel of figure 3.

Identity

As mentioned in section 4.1, identity keeps the model M( -→ m ) of figure 7 unchanged. Thus, in this case:

R(Id(M( - → m ))) = Id(M( - → m )) = M( - → m ) (model of figure 7)
As a result, the reverse migration leads to the model of figure 7 which corresponds to the model of figure 3 without actions. Ideally, these actions should have been recovered but actually, they are finally lost.

Copy

According to the definition of copy (cf. section 4.2), Cp(M( -→ m )) leads to the model m of figure 9, and thus:

R(Cp(M( - → m ))) = Cp(M( - → m )) = m (model of figure 9)
The result of the reverse migration is a model which is equivalent to the model of figure 7 except the instances have been recreated. As for identity, this resulting model corresponds to the model of figure 3 without actions. Thus in this case as in the previous one, these actions should have been recovered but actually they are lost.

Flattener

The result of flattener applied to the model of figure 3 is the model noted m in figure 10 (cf. section 4.3), and thus:

R(Fl(M( - → m ))) = Fl(M( - → m )) = m (model of figure 10)
The result of the reverse migration in this case contains a unique state named ABC and without actions. Actually, this state has been modified by the flattener through the initial states named A, B and C. These states had entry actions. These original actions are not taken into account by the reverse migration and then, they cannot be recovered.

Recontextualization by keys

So far, the reverse migration does not enable the recovery of elements that are deleted by the initial migration. We call initial context of a given model m the set of model elements from m that are deleted by a migration specification -→ m . By definition of migration specifications (see figure 4),

this initial context corresponds to - → m .D i ∪ - → m .D a ∪ - → m .D r .
We call recontextualization the tool aiming at putting back the initial context of a given model m in its reverse migrated version after the processing of a rewriting tool.

In a first approach, we use instances as keys to reconnect a reverse migrated model to its initial context. We note C k this tool:

C k (R(T(M( - → m ))) = m (recontextualized reversed migrated model)
Its action can be informally described as follows:

1. deleted instances are recovered

m r = R(T(M( - → m )) (reversed migrated model) m = C k (m r ) (recontextualized reversed migrated model) m .V = m r .V ∪ - → m .D i m .E a = m r .E a ∪ {(i, a, s) ∈ m.E a | i ∈ - → m .D i ∨ ((i, a) ∈ - → m .D a ∧ i ∈ m r .V )} m .E r = m r .E r ∪ {(i, r, i ) ∈ m.E r | (i, i ) ∈ m .V 2 ∧ (i ∈ - → m .D i ∨ i ∈ - → m .D i ∨ (i, r) ∈ - → m .D r )}
Figure 11: Computation of a recontextualized model by keys

Identity

As previously mentioned, identity keeps the input model unchanged and the reverse migration leads to the model of figure 7 which corresponds to the model of figure 3 without actions:

m r = R(Id(M( - → m ))) (model of figure 7) m = C k (m r ) (recontextualized reversed migrated model)
By definition of the recontextualization by keys (see figure 11), m is defined as follows:

m .V = m r .V ∪ {sa 0 , sa 1 , sa 2 } m .E a = m r .E a ∪ {(sa 0 , name, a 0 ), (sa 1 , name, a 1 ), (sa 2 , name, a 2 )} m .E r = m r .E r ∪ {(s 1 , entry, sa 0 ), (s 2 , entry, sa 2 ), (sa 2 , sin, s 2 ), (sa 0 , sin, s 1 ), (sa 0 , next, sa 1 ), (sa 1 , prev, sa 0 )}

The deleted actions are actually recovered and as expected, the recontextualized reversed migrated model corresponds to the initial model of figure 3:

C k (R(Id(M( - → m ))) = m
This situation is depicted by figure 12. In this figure, the initial context appear in red and the links between the reversed migrated model and the initial context appear in blue.

Copy

As previously mentioned, the copy of M( -→ m ) leads to the model of figure 9 which is equivalent to the model of figure 7 except the instances have been recreated :

m r = R(Cp(M( - → m ))) (model of figure 9) m = C k (m r )
(recontextualized reversed migrated model)

As a consequence and by definition of the recontextualization by keys (see figure 11), the initial context can be put back but not reconnected to the reversed migrated model. The resulting model actually contains the recovered instances and the possible links between them (attributes or references), but it does not contains reference links between the recovered instances and the instances newly created by copy:

m .V = m r .V ∪ {sa 0 , sa 1 , sa 2 } m .E a = m r .E a ∪ {(sa 0 , name, a 0 ), (sa 1 , name, a 1 ), (sa 2 , name, a 2 )} m .E r = m r .E r ∪ {(sa 0 , next, sa 1 ), (sa 1 , prev, sa 0 )}
The deleted actions are actually recovered but they are not connected to states. This situation is depicted by figure 13. In this case, the graph is not connected. Each connected component of such graphs is the semantics of a sub-model conforming to the initial metamodel (of figure 3 in this example), and which can be valid or not with regard to multiplicity or containment.

Flattener

The flattening of M( -→ m ) leads to the model of figure 10 which contains a unique state named ABC and without actions:

m r = R(Fl(M( - → m ))) (model of figure 10) m = C k (m r )
(recontextualized reversed migrated model)

The only instance s 3 was in the initial model. Its name has been changed and by definition of the recontextualization by keys (see figure 11), this change is not undone. The initial instances s 1 and s 2 have been deleted by the flattener and they are obviously not recovered by the recontextualization. Actually, as in the case of copy, the resulting model only contains the recovered instances corresponding to actions. However it does not contains reference links between the recovered instances and the remaining instance s 3 :

m .V = m r .V ∪ {sa 0 , sa 1 , sa 2 } m .E a = m r .E a ∪ {(sa 0 ,
name, a 0 ), (sa 1 , name, a 1 ), (sa 2 , name, a 2 )} m .E r = m r .E r ∪ {(sa 0 , next, sa 1 ), (sa 1 , prev, sa 0 )} As in the previous case, it leads to a graph which is not connected. This result is depicted by figure 14. 7 Recontextualization by graphs So far the rewriting tool T is considered as a black box and as a consequence, its action can only be specified by the differences between a given input model and the corresponding output model.

In order to improve the recontextualization of a reverse migrated model, T is considered now as a grey box where the computation functions are not known but where the computation dependencies are known. For that purpose, for a given input model m, we suppose now that the tool T

G(T, m)

T(m).V → P(m.V ) (edges from output instances to input instances)

Figure 15: Dependency graph provides the resulting model T(m) together with a dependency graph noted G(T, m). For a given instance i of T(m), this graph makes it explicit the instances of m that have been used to create or update i. It is formally defined in figure 15 by means of a total function over the output instances (T(m).V ) to parts of the input instances. Thus, G(T, m)(i) provides the (possibly empty) set of input instances that are used to compute i.

The use of dependency graphs enhances the recontextualization. Indeed, the links between the reverse migrated model and it initial context can be created not only over the remaining instances, but also over new instances through the instances that have been used to compute them.

In this second approach, we note C g the tool using instances as vertices of dependency graphs to reconnect a recontextualized reverse migrated model to its initial context. Its actions complete the actions performed by C k (i.e. recontextualization by keys). It applies to a model together with a dependency graph:

C g C k (R(T(M( - → m ))) , G(T, M( - → m )) = m (enhanced recontextualized reversed migrated model)
The actions performed by C g can be informally described as follows:

1. deleted attributes of a non-deleted instance A which has been afterwards deleted by the tool are put on each instance B computed from A (according to the dependency graph, B depends on A)

2. references relating a deleted instance X and a non-deleted instance A which has been afterwards deleted by the tool are put between the recovered instance X and each instance B computed from A (according to the dependency graph B depends on A)

3. deleted references relating a non-deleted instance A1 and a non-deleted instance A2 which has been afterwards deleted by the tool are put between the remaining instance A1 and each instance B computed from A2 (according to the dependency graph B depends on A2)

4. deleted references relating the non-deleted instances A1 and A2, afterwards both deleted by the tool are put between each instance B1 computed from A1 and each instance B2 computed from A2 (according to the dependency graph B1 depends on A1 and B2 depends on A2)

More formally, the recontextualization by a dependency graph is specified by figure 16. To illustrate this approach, we consider again the three examples of rewriting tools namely identity (Id), copy (Cp), and flattener (Fl). We suppose that the dependency graph of identity is defined as follows. It simply states that each instance is computed (without any change) from and only from itself:

m k = C k (R(T(M( - → m ))) (recontextualized reversed migrated model) g = G(T, M( - → m )) (dependency graph) m = C g (m k , g) (enhanced recontextualized reversed migrated model) m .V = m k .V m .E a = m k .E a ∪ {(i, a, s) ∈ T(M( - → m )).V × A × S | ∃i ∈ g(i) \ T(M( - → m )).V, (i , a, s) ∈ m.E a ∧ (i , a) ∈ - → m .D a } m .E r = m k .E r ∪ {(i 1 , r, i 2 ) ∈ m k .V × R × m k .V | ∃i 1 ∈ g(i 1 ) \ m k .V, (i 1 , r, i 2 ) ∈ m.E r ∧ ((i 1 , r) ∈ - → m .D r ∨ i 2 ∈ - → m .D i )} ∨ ∃i 2 ∈ g(i 2 ) \ m k .V, (i 1 , r, i 2 ) ∈ m.E r ∧ ((i 1 , r) ∈ - → m .D r ∨ i 1 ∈ - → m .D i )} ∨ ∃(i 1 , i 2 ) ∈ (g(i 1 ) \ m k .V ) × (g(i 2 ) \ m k .V ), (i 1 , r, i 2 ) ∈ m.E r ∧ (i 1 , r) ∈ - → m .D r }
G(Id, M( - → m )) : Id(M( - → m )).V → P(M( - → m ).V ) s 1 → {s 1 } s 2 → {s 2 } s 3 → {s 3 }
By definition of the recontextualization by dependency graphs (fig. 16), we have:

m = C g C k (R(Id(M( - → m ))) , G(Id, M( - → m )) m .V = C k (R(Id(M( - → m ))).V m .E a = C k (R(Id(M( - → m ))).E a ∪ ∅ m .E r = C k (R(Id(M( - → m ))).E r ∪ ∅ And thus : m = C k (R(Id(M( - → m ))).
In the case of identity, the recontextualized model by dependency graph brings nothing compared to the recontextualized model by keys. This situation is depicted by figure 17.

Copy

We consider again the initial model m of figure 3, and the migrated model M( -→ m ) of figure 7 corresponding to the initial model without actions. The application of copy to this migrated model leads to an equivalent model where instances have been recreated. In this case, the dependency graph is defined as follows. It simply maps new instances to the corresponding initial instances:

G(Cp, M( - → m )) : Cp(M( - → m )).V → P(M( - → m ).V ) s 4 → {s 1 } s 5 → {s 2 } s 6 → {s 3 }
By definition of the recontextualization by dependency graphs (fig. 16), we have:

m = C g C k (R(Cp(M( - → m ))) , G(Cp, M( - → m )) m .V = C k (R(Cp(M( - → m ))).V m .E a = C k (R(Cp(M( - → m ))).E a ∪ ∅ m .E r = C k (R(Cp(M( - → m ))
).E r ∪ {(s 4 , entry, sa 0 ), (s 5 , entry, sa 2 ), (sa 2 , sin, s 5 ), (sa 0 , sin, s 4 )} Thanks to this mapping, entry actions can be connected to new instances as far as they are a copy of a deleted instance which was connected to actions before the migration. Finally, an equivalent model of the initial model is actually recovered. This situation is depicted by figure 18.

Flattener

We consider again the initial model m of figure 3, and the migrated model M( -→ m ) of figure 7 corresponding to the initial model without actions. The application of flatten to this migrated model leads to the model of figure 10 which contains a unique state named ABC and without actions. As an illustration in this case, we suppose that the dependency graph is defined as follows.

The flat state corresponds to a collapsed sequence of nested states. Thus, the remaining state was already a simple state before the flattening and it has been updated by all the super-states to which it belonged:

G(Fl, M( - → m )) : Fl(M( - → m )).V → P(M( - → m ).V ) s 3 → {s 1 , s 2 , s 3 }
By definition of the recontextualization by dependency graphs (fig. 16), we have:

m = C g C k (R(Fl(M( - → m ))) , G(Fl, M( - → m )) m .V = C k (R(Fl(M( - → m ))).V m .E a = C k (R(Fl(M( - → m ))).E a ∪ ∅ m .E r = C k (R(Fl(M( - → m ))
).E r ∪ {(s 3 , entry, sa 0 ), (s 3 , entry, sa 2 ), (sa 2 , sin, s 3 ), (sa 0 , sin, s 3 )} Thanks to this mapping, entry actions can be connected to the remaining instance s 3 . However, the resulting graph does not corresponds to a valid model with regard to the initial metamodel (see figure 3). Indeed, the model we obtain contains two entry actions whereas at most one is expected. This situation is depicted by figure 19. At this point, an update of the recontextualized model is needed in order to make it conform to the initial metamodel. This update requires domain specific knowledge and is therefore out of the scope of automatic recontextualization. In this example, action a2 clearly needs to be added to the first sequence of actions a0 ; a1. By this way, there is only one entry action associated to the unique state ABC. As a result, we obtain a flattened version of the initial statecharts including entry actions whereas the reused flattener does not take them into account. This last transformation rely on the usual operational semantics of statecharts. It is illustrated by figure 20.

Properties

The recontextualization by graphs is intended to extend the recontextualization by keys. The recontextualization by graphs (and therefore especially by keys) must not undo the modifications performed by the rewriting tool. These two properties can be verified within the formal framework of our semantic domain.

Recontextualization by graphs extends recontextualization by keys

Let m be a given model, -→ m a migration specification and T a rewriting tool. Let G(T, M( -→ m )) be a dependency graph over T(M( -→ m )).

An edge of G(T, M( -→ m )) relating an instance i to itself has no added value. Indeed, this edge states that i is an instance of the result which was already in the input of T, and which value depends (obviously) on itself.

We say that G(T, M( -→ m )) is semantically empty if all its edges have no added values. More formally, by definition:

G(T, M( - → m )) is semantically empty ∀i ∈ T(M( - → m )).V, G(T, M( - → m ))(i) ⊆ {i} Property 1 If G(T, M( - → m
)) is semantically empty, then the recontextualization by graphs corresponds to the recontextualization by keys:

C g C k (R(T(M( - → m ))) , G(T, M( - → m )) = C k (R(T(M( - → m ))) Proof Let m k = C k (R(T(M( - → m
m .V = m k .V
By definition of C g , we also have:

m .E a = m k .E a ∪ {(i, a, s) ∈ T(M( - → m )).V × A × S | ∃i ∈ g(i) \ T(M( - → m )).V, (i , a, s) ∈ m.E a ∧ (i , a) ∈ - → m .D a }
By definition of a semantically empty graph, ∀i ∈ T(M( -→ m )).V, g(i) ⊆ {i}. Thus we have:

∀i ∈ T(M( - → m )).V, g(i) \ T(M( - → m )).V = ∅
And then we have:

m .E a = m k .E a ∪ {(i, a, s) ∈ T(M( - → m )).V × A × S | ∃i ∈ g(i) \ T(M( - → m )).V, . . . } = m k .E a ∪ {(i, a, s) ∈ T(M( - → m )).V × A × S | ∃i ∈ ∅, . . . } = m k .E a ∪ ∅ = m k .E a
By definition of C g , we finally have:

m .E r = m k .E r ∪ {(i 1 , r, i 2 ) ∈ m k .V × R × m k .V | ∃i 1 ∈ g(i 1 ) \ m k .V, (i 1 , r, i 2 ) ∈ m.E r ∧ ((i 1 , r) ∈ - → m .D r ∨ i 2 ∈ - → m .D i )} ∨ ∃i 2 ∈ g(i 2 ) \ m k .V, (i 1 , r, i 2 ) ∈ m.E r ∧ ((i 1 , r) ∈ - → m .D r ∨ i 1 ∈ - → m .D i )} ∨ ∃(i 1 , i 2 ) ∈ (g(i 1 ) \ m k .V ) × (g(i 2 ) \ m k .V ), (i 1 , r, i 2 ) ∈ m.E r ∧ (i 1 , r) ∈ - → m .D r }
For the same reason and by definition of semantically empty graphs, if

i 1 ∈ m k .V , then g(i 1 ) ⊆ {i 1 } ⊆ m k .V and if i 2 ∈ m k .V , then g(i 2 ) ⊆ {i 2 } ⊆ m k .V . And thus g(i 1 ) \ m k .V = ∅ and g(i 2 ) \ m k .V = ∅. Hence: m .E r = m k .E r ∪ {(i 1 , r, i 2 ) ∈ m k .V × R × m k .V |∃i 1 ∈ ∅, • • • ∨ ∃i 2 ∈ ∅, • • • ∨ ∃(i 1 , i 2 ) ∈ ∅ 2 , . . . } = m k .E r ∪ ∅ = m k .E r
The three components of m and m k coincide and thus we actually have m = m k .

Recontextualization does not undo the rewriting

As mentioned in section 4, the action of a rewriting tool T is specified by the differences between a given input model and the corresponding output model. This specification is formally stated in figure 8 by deleted elements T i D (m), T a D (m), T r D (m) and added elements T i A (m), T a A (m), T r A (m).

The recontextualization aims at resetting the initial context on the result of T without undoing its specific actions. Thus, the model elements that have been added by the tool must not be deleted by the recontextualization, and the the model elements that have been deleted by the tool must not be recovered by the recontextualization. This property holds under the following condition which means that the instances that are created by the tool (T i A (M( -→ m ))) have new names with regard to instances of the initial context.

T i A (M( - → m )) ∩ m.V = ∅
This condition can be easily satisfied by a renaming of the instances created by the tool (thanks to a kind of post-processor ).

Property 2 Let m be a given model and T a rewriting tool. Let m be the recontextualized reversed migrated model by graphs: m = C g C k (R(T(M( -→ m ))) , G(T, M( -→ m )) . Then we have:

(1) deleted elements are not added:

T i D (M( - → m )) ∩ m .V = ∅ ∧ T a D (M( - → m )) ∩ m .E a = ∅ ∧ T r D (M( - → m )) ∩ m .E r = ∅
(2) added elements are not removed:

T i A (M( - → m )) ⊆ m .V ∧ T a A (M( - → m )) ⊆ m .E a ∧ T r A (M( - → m

A Modif

A.1 Modif operators

A.2 Modif specification

We note C the name space for classes. This name space is only relevant at the metamodel level. We call Modif specification and note modif a possibly empty sequence of operators noted op (see figure 21).

FormalFigure 1 :

 1 Figure 1: ISAR: name spaces for model graphs

Figure 3 :

 3 Figure 3: Metamodel, conforming model, concrete syntax and semantics

Figure 4 :

 4 Figure 4: Migration specification

Figure 5 :

 5 Figure 5: Computation of a migrated model

Figure 6 :

 6 Figure 6: Metamodel of statecharts without actions

Figure 10 :

 10 Figure 10: Flattening of the model from figure 7

  2. attributes of deleted instances are recovered 3. deleted attributes are recovered as far as their source instances still exist 4. references between two deleted instances are recovered 5. references from deleted instances targeting still existing instances are recovered 6. references from still existing instances targeting deleted instances are recovered 7. deleted references are recovered as far as their source and target instances still exist More formally, the recontextualization by beys on graphs is specified by figure 11. To illustrate this approach, we consider again the three examples of rewriting tools namely identity (Id), copy (Cp), and flattener (Fl).

Figure 12 :

 12 Figure 12: Recontextualization by keys of the model from figure 7

Figure 13 :

 13 Figure 13: Recontextualization by keys of the model from figure 9

Figure 14 :

 14 Figure 14: Recontextualization by keys of the model from figure 10

Figure 16 :

 16 Figure 16: Computation of a recontextualized model by a dependency graph

Figure 17 :

 17 Figure 17: Recontextualization by graph of the model from figure 7

Figure 18 :

 18 Figure 18: Recontextualization by graph of the model from figure 9

Figure 19 :

 19 Figure 19: Recontextualization by graph of the model from figure 10

Figure 20 :

 20 Figure 20: Transformation to a valid model with regard to the initial metamodel

  ))) be the recontextualized reversed migrated model by keys. Let g = G(T, M( -→ m )) be a semantically empty dependency graph. Let m = C g (m k , g) be the enhanced recontextualized model by g. By definition of C g (fig. 16), we have:

  )) ⊆ m .E r Proof sketch The second point is trivial since neither C k nor C g remove anything from the model T(M( -→ m )) resulting from T (cf. figures 11 and 16).About the first point, we consider the elements that are added by C k , and then by C g . By definition of C k and T i D , knowing that R is identity and that M( -→ m ).V = m.V \ -→ m .D i we have:m .V = C k (R(T() ∩ m .V = (m.V \ -→ m .D i ) \ T(M( -→ m )).V ∩ T(M( -→ m )).V ∪ -→ m .D i = m.V \ ( -→ m .D i ∪ T(M( -→ m )).V ) ∩ T(M( -→ m )).V ∪ -→ m .D i = ∅The proofs ofT a D (M( -→ m ))∩m .E a = ∅ and of T r D (M( -→ m ))∩m .E r = ∅follow the same principles. By definition of M( -→ m ).E a and of M( -→ m ).E a (cf. figure 5), T does not takes as input any edge (attribute or reference) which is deleted ( -→ m .D a or -→ m .D r ), or whose source and/or target is deleted ( -→ m .D i ). T cannot delete such unavailable edges. Now by definition of C k and C g (cf. figures 11 and 16) the only added edges are taken from -→ m .D a and -→ m .D r , or they are related to a deleted instance from -→ m .D i .

modif op 1 A

 1 Figure 21: Modif specification

Table 1 :

 1 Among the metamodel refactoring operators that are available in Modif, we focus on delete. It applies to classes, attributes and references. Its scope and parameters are defined in table 1. Modif operators: scope and parameters

	Operator Scope	Parameters
	delete	class	class: c
	delete	attribute class: c; name: a
	delete	reference class: c; name: r

• together with the metamodel to produce a refactored metamodel

• together with the metamodel and a conforming model to produce a migration specification

A.3 Generation of migration specifications

We suppose we have a metamodel noted mm, a conforming model noted m and a valid Modif specification noted modif. The scope of modif is mm.

We call migration generator and note G M the tool aiming at producing a valid migration specification from mm, m and modif:

This generator compile the sequence of operator calls in modif to produce the resulting migration specification.


F Leblanc

J.-Y Chaufray

R Modolo

UVSQ/ S Hess

M Yagi

F Gonzalez-Galindo

L Lorenzato

C Mazelle

G M Chanteur

MARS-SOLAR WIND INTERACTION: COUPLING BETWEEN HYBRID, IONOSPHERIC, THERMOSPHERIC AND EXOSPHERIC MODELS

Introduction

The solar wind interaction with the Martian neutral environment is investigated by means of three dimensional hybrid magnetospheric simulations. In such formalism, ions have a kinetic description while electrons are treated as an inertialess fluid, ensuring the neutrality of the plasma and contributing to currents and pressure terms. This model has been successfully used to describe the near ionized environment of Mars [START_REF] Modolo | Influence of Solar EUV flux on the Martian plasma environment[END_REF][START_REF] Modolo | Simulated solar wind plasma interaction with the Martian exosphere: influence of the solar EUV flux on the bow shock and the magnetic pile-up boundary[END_REF][START_REF] Modolo | Dynamic Martian magnetosphere: Transient twist induced by a rotation of the IMF[END_REF]. The main drawback of the hybrid formalism is the coarse spatial resolution that a sequential program with a uniform grid description can afford (about 130-150 km), mainly restricted by computational limitation (memory and CPU). In the frame of the HELIOSARES project (PI F. Leblanc) dedicated to the modeling of Mars environment (neutral and charged species) from the surface up to the solar wind, a modeling effort of parallelization has been conducted. In parallel, developments concerning a three-dimensional exospheric model and an ionospheric extension of a Global Circulation Model [START_REF] Forget | Improved general circulation models of the Martian atmosphere from the surface to above 80 km[END_REF] have been performed. The 3D exospheric model provides a description of the Martian exosphere of CO2, H, H2 and O from the exobase level to few Martian radii [START_REF] Yagi | Mars exospheric thermal and non-thermal components: seasonal and local variation[END_REF], including thermal and non-thermal (as produced by dissociative recombination and sputtering) contributions to the O component. On the other hand, a 3D multifluid dynamical core describing the dynamic of the Martian ionospheric plasma has been included in the LMD Martian General circulation Model [START_REF] Chaufray | 3D Martian Ionosphere model : II Effect of transport processes without magnetic field[END_REF][START_REF] Gonzalez-Galindo | 3D Martian ionosphere model : I. The photochemical ionosphere below 180 km[END_REF] that has been extended up to the exobase.

Methodology

The exospheric 3D model used as inputs the descriptions of the thermosphere and ionosphere provided by the extended Martian General Circulation Model [START_REF] Chaufray | 3D Martian Ionosphere model : II Effect of transport processes without magnetic field[END_REF]. In a similar approach than [START_REF] Chaufray | Mars solar wind interaction: Formation of the Martian corona and atmospheric loss to space[END_REF], simulations outputs of the exospheric model are then used as inputs of the magnetospheric model. A realistic and three dimensional description of the neutral environment is therefore used to constrain the production in the magnetospheric model. This approach is also used for the magnetosphere-ionosphere models coupling, the ionospheric density distribution of the ionospheric model [START_REF] Chaufray | 3D Martian Ionosphere model : II Effect of transport processes without magnetic field[END_REF] being used as inputs of the magnetospheric model.

The exospheric model

This model provides the description of the thermal exospheric components (H, H2, O and CO2) as well as a description of the non-thermal oxygen component. The non-thermal oxygen contribution is estimated by means of a Monte-Carlo approach with a description of the thermosphere determined by the extended general circulation model LMD-MGCM [START_REF] Chaufray | 3D Martian Ionosphere model : II Effect of transport processes without magnetic field[END_REF]. In this model, hot oxygen atoms are either produced by dissociative recombination of the O2 + as modeled by the extended LMD-GCM model or by sputtering using as inputs the flux of pick-up ions reimpacting Mars' exobase as calculated by the magentospheric model. Hot test-particles are then followed through the extended LMD-GCM thermospheric model. The contribution of the thermal part of the exosphere is obtained by using a Chamberlain approach extended to three-dimension and adapted to include planet's rotation and using the extended LMD-GCM model as bottom conditions. Simulation results of this model are presented in [START_REF] Yagi | Mars exospheric thermal and non-thermal components: seasonal and local variation[END_REF].

The ionospheric model

The LMD-GCM [START_REF] Forget | Improved general circulation models of the Martian atmosphere from the surface to above 80 km[END_REF][START_REF] Gonzalez-Galindo | A ground-to-exosphere Martian general circulation model: 1. Seasonal, diurnal, and solar cycle variation of thermospheric temperatures[END_REF] describing the Martian thermosphere as been completed with the implementation of photoionization and secondary ionization by X-rays and photo electron impact. These inputs allow describing the Martian ionosphere at low altitudes where the ion transport can be neglected [START_REF] Gonzalez-Galindo | 3D Martian ionosphere model : I. The photochemical ionosphere below 180 km[END_REF]. [START_REF] Chaufray | 3D Martian Ionosphere model : II Effect of transport processes without magnetic field[END_REF] extended this model by implementing an ambipolar diffusion in the ionosphere. This improvement allows investigating the upper ionosphere and its coupling with the magnetosphere.

The hybrid model

The hybrid model used in this work has been presented in [START_REF] Modolo | Influence of Solar EUV flux on the Martian plasma environment[END_REF], parallelized and completed with a realistic ionospheric description with the implementation of a simplified chemistry model [START_REF] Modolo | Solar Wind interaction with Mars neutral environment from hybrid simulations: improved spatial resolution[END_REF]. An updated version of the model with a crustal field description is presented by [START_REF] Hess | Hybrid simulation of the Solar wind interaction with the martian ionosphere and magnetic field[END_REF]. This simulation model describes the electromagnetic environment and the dynamic of 6 ion species (H + sw , He ++ , H + pl , O + , CO + 2 and O + 2 ). Simulations are performed on a meso-scale calculation server CICLAD ( http://cicladweb.ipsl.jussieu.fr) dedicated to Institut Pierre Simon Laplace modeling effort : The project has access to 3 servers of 32-cores and 3 servers of 64-cores with 4Go memory per CPU. Communications are perfomed with infinite-band cables.

Results

In this presentation, we will describe the results of a first coupling between these three models. Simulations are performed at mean solar activity, at Ls=90 • and for nominal solar wind conditions (with the crustal field placed at midnight) for all exospheric, ionospheric-thermospheric and magnetospheric models. The potential state of Mars' environment from the solar wind up to the surface will be presented. Starting from the magnetosphere, we will discuss the possible signatures of Mars' season on the main structures in the magnetosphere (Bow Shock and magnetopause), in the exosphere (thermal CO2, O, H and H2 exospheric seasonal signatures and non-thermal O densities as induced by both sputtering and dissociative recombination) and in the ionosphere/thermosphere (role of the ionospheric transport).

Conclusion

Future developments are still planned in order to take into account the feedback from the magnetosphere onto the ionosphere/thermosphere and therefore the exosphere. These will be done during the forthcoming years. All these developments are done in the frame of MAVEN science preparation.

Acknowledgements

Authors are indebted to the ANR HELIOSARES (ANR-09-BLAN-0223), The French National Program PNST, the French Agency CNES, the CSA and IPSL for their support.