Panagiotis Kotronis

Simplified modeling strategies for soil-structure interaction problems: The multifiber beam concept

Simplified modeling strategies for soil-structure interaction problems: The multifiber beam concept

Introduction

The objective of this course is to introduce the differences between the Euler Bernoulli and Timoshenko theories, to present various finite element beam formulations and to show the equations of displacement based Timoshenko multifiber beam elements. The paper follows mainly the work and ideas exposed in [START_REF] Pegon | A timoshenko simple beam element in castem[END_REF], [START_REF] Guedes | A fibre timoshenko beam element in castem[END_REF], [START_REF] Kotronis | Simplified modelling strategies to simulate the dynamic behaviour of r/c walls[END_REF], [START_REF] Mazars | Using multifiber beams to account for shear and torsion: Applications to concrete struc-tural elements[END_REF], [START_REF] Bitar | Un nouvel élément fini poutre pour la modélisation des structures sous sollicitations complexes[END_REF], [CKCed]. The important subject of force based beam finite elements is not discussed hereafter. The reader can find information on this subject in the following references [START_REF] Spacone | Fiber beam-column model for nonlinear analysis of r/c frames. i[END_REF], [START_REF] Spacone | Fiber beam-column model for nonlinear analysis of r/c frames[END_REF].

1 2 Classical beam theories

Kinematics

We consider hereafter a beam of length L and section S(x) (figure 1). G(x, y, z) is the center of gravity of the section S(x) and P (x, y, z) a point in the section. We define as neutral axis the line that links the center of gravities of all the sections. We also suppose that G x , G y , G z are principal axes. For the 2D case studied hereafter (loadings are in the xy plane), the displacements u T = {u x , u y } of the point P (x, y, z) can be expressed as a function of the displacements U x (x), U y (x) and the rotation Θ z (x) of the section S(x) (often defined in the literature as generalised displacements, see also equations ( 8) and ( 16)) considering the following two kinematic hypotheses: • Euler Bernoulli theory: In this theory, the section remains plane and perpendicular to the neutral axis (figure 2).

Figure 2: Euler Bernoulli theory [START_REF] Bitar | Un nouvel élément fini poutre pour la modélisation des structures sous sollicitations complexes[END_REF] A consequence of the previous kinematic assumption is that the rotation of the section Θ z (x) equals U y (x) (the symbol defines hereafter the first derivative with respect to x and the symbol the second derivative with respect to x)

(figure 2). The displacements of the point P (x, y, z) take thus the following form [START_REF] Guedes | A fibre timoshenko beam element in castem[END_REF]:

u x (x, y) = U x (x) -yΘ z (x) = U x (x) -yU y (x) u y (x, y) = U y (x) (1)
and the strains are calculated as (infinitesimal strain theory assumption):

ε x = ∂u x ∂x = U x (x) -yΘ z (x) = U x (x) -yU y (x) γ xy = ∂u x ∂y + ∂u y ∂x = U y (x) -Θ z (x) = 0 (2)
One can notice that due to the adopted kinematic hypothesis shear strains are found equal to zero.

• Timoshenko theory: In this theory, the kinematic assumption is that the section remains plane but not necessarily perpendicular to the neutral axis. In other words, Θ z = U y (x), see figure 3. Displacements and strains (infinitesimal strain theory assumption) are now calculated as [START_REF] Guedes | A fibre timoshenko beam element in castem[END_REF]:

u x (x, y) = U x (x) -yΘ z (x) u y (x, y) = U y (x) (3)
ε x = ∂u x ∂x = U x (x) -yΘ z (x) γ xy = ∂u x ∂y + ∂u y ∂x = U y (x) -Θ z (x) (4)
Figure 3: Timoshenko theory (Wikipedia)

We define hereafter

β y = U y (x) -Θ z (x) (5)
The variable β y represents the rotation of the section due to shear. From equation (4) one can notice that shear strains are now constant in the section (and not as in the Euler Bernoulli theory necessarily equal to zero).

Euler Bernoulli theory: Principal work principle and generalised forces

The virtual work principle is written as (body forces and inertial forces are hereafter neglected. The symbol σ defines stresses):

L 0 S δε x σ x dSdx -w exter = 0 (6)
with w exter the work of the external forces.

Replacing (2) in (6) we get:

L 0 S δ(U x (x) -yΘ z (x))σ x dSdx -w exter = 0 (7)
We define hereafter the generalised forces in the section as:

Normal force:

F x = S σ x dS
Bending moment :

M z = - S yσ x dS (8)
and equation (7) becomes:

L 0 (F x δU x + M z δΘ z )dx -w exter = 0 (9)
Within a beam theory σ y = σ z = σ yz = 0. Furthermore and because of equation (2) σ xz = σ xy = 0. Hooke's law thus becomes (with E the Young's modulus and ν the Poisson's coefficient):

σ x = Eε x ε y = ε z = -νε x (10)
Finally, using equations (2), ( 8) and (10) the generalised forces become:

F x = S Eε x dS = S E(U x (x) -yΘ z (x))dS = S EU x (x)dS = ESU x (x) M z = - S yEε x dS = - S yE(U x (x) -yΘ z (x))dS = S y 2 EΘ z (x)dS = EI z Θ z (x) (11)
Remarks:

• The axis z is a principal axis and therefore S ydS = 0.

• An homogeneous section is considered.

Introducing equation (11) in the virtual work principle (9) we have:

L 0 (δU x (x)ESU x (x) + δΘ z (x)EIΘ z (x))dx -w exter = 0 (12)
If F s T = {F x , M z } the generalised force vector and D s T = {U x , Θ z } the generalised displacement vector we define the stiffness matrix of the section K s as:

F s = K s D s = F x M z = ES 0 0 EI z U x Θ z = ES 0 0 EI z U x U y (13)

Timoshenko theory: Principal work principle and generalised forces

The principle work principle now becomes:

L 0 S (δε x σ x + 2δε xy σ xy )dSdx -w exter = 0 (14)
Using equation (4) we get:

L 0 S (δ(U x (x) -yΘ z (x))σ x + (δU y (x) -δΘ z (x))σ xy dSdx -w exter = 0 (15)
where the generalised forces are:

Normal force:

F x = S σ x dS
F x = S σ x dS = S Eε x dS = S E( dU x dx -y dΘ z dx )dS = S EU x ds = ESU x (18)
In a similar way we have:

F y = GSβ y M z = EI z Θ z (19)
If F s T = {F x , F y , M z } the generalised force vector and D s T = {U x , β y , Θ z } the generalised strain vector we define the stiffness matrix of the section K s as:

F s = K s D s = ⎧ ⎨ ⎩ F x F y M z ⎫ ⎬ ⎭ = ⎡ ⎣ ES 0 0 0 GS 0 0 0 EI z ⎤ ⎦ ⎧ ⎨ ⎩ U x β y Θ z ⎫ ⎬ ⎭ (20)
Remark: The Timoshenko beam theory provides constant shear strains and stresses in the section (see equation (4) and Hooke's law). This result violates the boundary conditions of the beam theory (σ y = σ z = 0) and does not agree with the theoretical distribution of stresses which is parabolic for a rectangular cross section. A simplified way to deal with this inconsistency is to change the definition of the shear force by adding a shear corrector factor (or Reissner corrector factor) k that depends on the cross section geometry and the material characteristics [START_REF] Cowper | The shear coefficient in timoshenko's beam theory[END_REF]. The modified expressions take thus the following form:

Shear force:

F y = S kσ xy dS (21)
F s = K s D s = ⎧ ⎨ ⎩ F x F y M z ⎫ ⎬ ⎭ = ⎡ ⎣ ES 0 0 0 kGS 0 0 0 EI z ⎤ ⎦ ⎧ ⎨ ⎩ U x β y Θ z ⎫ ⎬ ⎭ (22)

A 2 node beam finite element formulation

Consider a 2D beam finite element with two nodes and three degrees of freedom per node (figure 4). Displacements U x (x), U y (x) and rotations Θ z (x) along the beam are discretized using the nodal displacements as follows:

U = N Φ (23) or ⎧ ⎨ ⎩ U x U y Θ z ⎫ ⎬ ⎭ = ⎡ ⎣ N 1 N 2 N 3 N 4 N 5 N 6 N 7 N 8 N 9 N 10 N 11 N 12 N 13 N 14 N 15 N 16 N 17 N 18 ⎤ ⎦ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U x1 U y1 Θ z1 U x2 U y2 Θ z2 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ (24)
where N i (x), i = 1, 18 the shape functions and Φ the nodal displacements. The equation providing the generalised strain vector D s becomes (a for axial, s for shear and b for bending):

D s = ⎧ ⎨ ⎩ U x β y Θ z ⎫ ⎬ ⎭ = ⎡ ⎣ B1 B2 B3 B4 B5 B6 B7 B8 B9 B10 B11 B12 B13 B14 B15 B16 B17 B18 ⎤ ⎦ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U x1 U y1 Θ z1 U x2 U y2 Θ z2 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ = ⎡ ⎣ B a B s B b ⎤ ⎦ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U x1 U y1 Θ z1 U x2 U y2 Θ z2 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ (25)
Using the previous equation in the virtual work principle we get a system of linear equations that take the following form (with F the vector of nodal forces at the element level and K e the element stiffness matrix):

F = K e Φ (26)
The element stiffness matrix K e is given by (where the symbol + means matrix assembly, adding in an adequate way the different degrees of freedom) :

K e = K a + K s + K b (27) K a = L 0 B T a ESB a dx K s = L 0 B T s kGSB s dx K b = L 0 B T b EIB b dx (28)

An Euler Bernoulli 2 node finite element beam

Consider for example the following classical shape functions for the horizontal, vertical and rotational degrees of freedom [START_REF] Frey | Analyse des structures et milieux continus: mécanique des structures[END_REF]:

N 1 = 1 - x l and N 4 = x l N 8 = 1 -3( x l ) 2 + 2( x l ) 3 and N 9 = x -2 x 2 l + x 3 l 2 N 11 = 3( x l ) 2 -2( x l ) 3 and N 12 = - x 2 l + x 3 l 2 N 14 = N 8 and N 15 = N 9 N 17 = N 11 and N 18 = N 12 (29)
In this formulation the rotational and vertical displacements are made interdependent. The other shape functions are considered equal to zero and so equation (24) becomes:

⎧ ⎨ ⎩ U x U y θ z ⎫ ⎬ ⎭ = ⎡ ⎣ N 1 0 0 N 4 0 0 0 N 8 N 9 0 N 11 N 12 0 N 8 N 9 0 N 11 N 12 ⎤ ⎦ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U x 1 U y 1 Θ z 1 U x 2 U y 2 Θ z 2 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ (30)
The different stiffness matrices, for a constant homogeneous section S, take the following form (the matrix due to shear is equal to the zero matrix):

K a = ES l 1 -1 -1 1 K b = EI l 3 ⎡ ⎢ ⎢ ⎣ 12 6l -12 6l 6l 4l 2 -6l 2l 2 -12 -6l 12 -6l 6l 2l 2 -6l 4l 2 ⎤ ⎥ ⎥ ⎦ (31)
K e = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ES/l 0 0 -ES/l 0 0 0 12EI/l 3 6EI/l 2 0 -12EI/l 3 6EI/l 2 0 6EI/l 2 4EI/l 0 -6EI/l 2 2EI/l -ES/l 0 0 ES/l 0 0 0 -12EI/l 3 -6EI/l 2 0 12EI/l 3 -6EI/l 2 0 6EI/l 2 2EI/l 0 -6EI/l 2 4EI/l ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (32)

A Timoshenko 2 node finite element beam (1st formulation)

For the case of a Timoshenko beam finite element, the choice of the shape functions is crucial because of the possible shear locking numerical problem [START_REF] Hughes | A simple and efficient finite element for plate bending[END_REF], [START_REF] Donea | A modified representation of transverse shear in c 0 quadrilateral plate elements[END_REF], [START_REF] Zienkiewicz | The finite element method for solid and structural mechanics[END_REF]. If these functions are not appropriately chosen, the finite element beam presents a spurious stiffness for the case of elongated beams. A simple adequate set of shape functions is presented hereafter [START_REF] Pegon | A timoshenko simple beam element in castem[END_REF], [START_REF] Guedes | A fibre timoshenko beam element in castem[END_REF]:

N 1 = N 8 = N 15 = x 2 -x L N 4 = N 11 = N 18 = x -x 1 L (33) ⎧ ⎨ ⎩ U x U y θ z ⎫ ⎬ ⎭ = ⎡ ⎣ N 1 0 0 N 4 0 0 0 N 8 0 0 N 11 0 0 0 N 15 0 0 N 18 ⎤ ⎦ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U x 1 U y 1 Θ z 1 U x 2 U y 2 Θ z 2 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ (34)
The generalised strains become:

D s = ⎡ ⎣ B a B s B b ⎤ ⎦ Φ = ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ -1 l 0 0 1 l 0 0 0 -1 l -x 2 -x l 0 1 l -x-x 1 l 0 0 -1 l 0 0 1 l ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U x 1 U y 1 Θ z 1 U x 2 U y 2 Θ z 2 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ (35)
In order to avoid shear locking problems, [START_REF] Donea | A modified representation of transverse shear in c 0 quadrilateral plate elements[END_REF] propose to eliminate the linear term in B s or [START_REF] Hughes | A simple and efficient finite element for plate bending[END_REF] to sub-integrate the K s using only 1 Gauss integration point. According to the solution proposed by [START_REF] Donea | A modified representation of transverse shear in c 0 quadrilateral plate elements[END_REF], [START_REF] Guedes | A fibre timoshenko beam element in castem[END_REF], β y is now corrected as:

β y = - 1 l U y 1 - 1 2 Θ z 1 + 1 l U y 2 - 1 2 Θ z 2 = B s Φ (36)
Finally, the element stiffness matrix, for a constant homogeneous section, takes the following form:

K e = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ES L 0 0 -ES L 0 0 0 kSG L kSG 2 0 -kSG L kSG 2 0 kSG 2 EI L + kSGL 4 0 -kSG 2 -EI L + kSGL 4 -ES L 0 0 ES L 0 0 0 -kSG L -kSG 2 0 kSG L -kSG 2 0 kSG 2 -EI L + kSGL 4 0 -kSG 2 EI L + kSGL 4 ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (37)

A Timoshenko 2 node finite element beam (2nd formulation)

Another way to avoid shear locking problems is to use higher order functions and to integrate exactly the stiffness matrix. For a 2 node finite element this leads to shape functions that depend on the material properties [START_REF] Ville De Goyet | L'analyse statique non linéaire par la méthode des éléments finis des structures spatiales formées de poutres à section non symétrique[END_REF], [START_REF] Friedman | An improved two-node timoshenko beam finite element[END_REF], [START_REF] Kotronis | Simplified modelling strategies to simulate the dynamic behaviour of r/c walls[END_REF], [START_REF] Mazars | Using multifiber beams to account for shear and torsion: Applications to concrete struc-tural elements[END_REF]. For example, according to [START_REF] Friedman | An improved two-node timoshenko beam finite element[END_REF] the shape functions become:

N = ⎡ ⎣ N 1 0 0 N 4 0 0 0 N 8 N 9 0 N 11 N 12 0 N 14 N 15 0 N 17 N 18 ⎤ ⎦ (38)
with:

N 1 =1 - x L N 4 = x L N 8 = 1 1 + φ [2( x L ) 3 -3( x L ) 2 -φ( x L ) + 1 + φ] N 9 = L 1 + φ [( x L ) 3 -(2 + φ 2 )( x L ) 2 + (1 + φ 2 )( x L )] N 11 = -1 1 + φ [2( x L ) 3 -3( x L ) 2 -φ( x L )] N 12 = L 1 + φ [( x L ) 3 -(1 - φ 2 )( x L ) 2 - φ 2 ( x L )] N 14 = 6 (1 + φ)L [( x L ) 2 -( x L )] N 15 = 1 1 + φ [3( x L ) 2 -(4 + φ)( x L ) + (1 + φ)] N 17 = -6 (1 + φ)L [( x L ) 2 -( x L )] N 18 = L 1 + φ [3( x L ) 2 -(2 -φ)( x L )] (39)
φ is the ratio between bending and shear stiffnesses. For an homogeneous section we get:

φ = 12 L 2 S Ey 2 dS S kGdS = 12 L 2 EI kGS (40)
The generalised strains now become:

D s = ⎡ ⎣ N 1 0 0 N 4 0 0 0 N 8 -N 14 N 9 -N 15 0 N 11 -N 17 N 12 -N 18 0 N 14 N 15 0 N 17 N 18 ⎤ ⎦ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U x 1 U y 1 Θ z 1 U x 2 U y 2 Θ z 2 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ (41)
Finally, the stiffness matrix of the element takes the following form:

K e = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ES L 0 0 -ES L 0 0 0 12EI (1+φ)L 3 6EI (1+φ)L 2 0 -12EI (1+φ)L 3 6EI (1+φ)L 2 0 6EI (1+φ)L 2 (4+φ)EI (1+φ)L 0 -6EI (1+φ)L 2 (2-φ)EI (1+φ)L -ES L 0 0 ES L 0 0 0 -12EI (1+φ)L 3 -6EI (1+φ)L 2 0 12EI (1+φ)L 3 -6EI (1+φ)L 2 0 6EI (1+φ)L 2 (2-φ)EI (1+φ)L 0 -6EI (1+φ)L 2 (4+φ)EI (1+φ)L ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (42)
For elongated beams, φ is almost zero and the matrix reduces to the classical Euler Bernoulli stiffness matrix, see equation (32). The main problem of this finite element formulation is that it is not appropriate for non linear calculations as the shape functions depend on the material properties that evolve with the loading. Nevertheless, some good results where obtained by keeping the initial shape functions (calculated at the first step -elasticity) unchanged, [START_REF] Kotronis | Simplified modelling strategies to simulate the dynamic behaviour of r/c walls[END_REF], [START_REF] Mazars | Using multifiber beams to account for shear and torsion: Applications to concrete struc-tural elements[END_REF].

Timoshenko finite elements with internal degrees of freedom

Another way to avoid shear locking problems is to enrich the displacements field [START_REF] Ibrahimbegovic | A modified method of incompatible modes[END_REF] and to add internal nodes in the element [START_REF] Ibrahimbegović | Geometrically non-linear method of incompatible modes in application to finite elasticity with independent rotations[END_REF], [CKCed]. In that way higher order shape functions are obtained that do not depend on the material properties. The numerical integration is exactly performed and the elements are suitable for non linear calculations.

The finite element presented in [CKCed] (named "FCQ" Timoshenko beam for "Full Cubic Quadratic") has additional internal degrees of freedom, cubic shape functions for the vertical displacements and quadratic for the rotations. The element is free of shear locking and one element is able to predict the exact tip displacements for any complex distributed loadings and any suitable boundary conditions. One element gives the exact solution for the case of bending of a Timoshenko beam free of distributed loadings. It is also proven that the element presented in [START_REF] Friedman | An improved two-node timoshenko beam finite element[END_REF] is a particular case of the more general FCQ Timoshenko beam element. For more information the reader is invited to read the relevant reference.

Multifiber beam formulation

The section of the finite element beam (Euler Bernoulli or Timoshenko) is divided in different "fibers" [START_REF] Owen | Finite elements in plasticity: Theory and practice[END_REF]. In each fiber a constitutive law is introduced (e.g. concrete, steel...). Depending on the mesh discretization of the section (e.g. with triangular or rectangular finite elements) one or more Gauss points are associated in each fiber, see figure 5.

Figure 5: Multifiber beam modelling [START_REF] Guedes | A fibre timoshenko beam element in castem[END_REF], [START_REF] Mazars | Using multifiber beams to account for shear and torsion: Applications to concrete struc-tural elements[END_REF] The different formulations of the finite element beams presented in section 2 (expressed now in 3D) can be used for the multifiber element. Consider equations (1) and (3) expressing the kinematic assumptions: the y (and z in 3D) are replaced with y f (and z f ), the coordinates of the fiber f in the section. The generalised forces take the following forms (where E f and G f the Young's and the shear moduli of the fiber respectively):

F x = S E f ε x dS = S E f ( dU x dx -y f dΘ z dx )dS = S E f dSU x - S E f y f dSΘ z (43)
and in a similar way:

⎧ ⎨ ⎩ F x F y M z ⎫ ⎬ ⎭ = ⎡ ⎣ S E f dS 0 -S E f y f dS 0 S kG f dS 0 -S E f y f dS S E f y 2 f dS ⎤ ⎦ ⎧ ⎨ ⎩ U x β y Θ z ⎫ ⎬ ⎭ (44)
By introducing the previous equation of the virtual work principle we obtain:

L 0 δD s T K s D s dx -w exter = 0 (45) Considering that D s = ⎡ ⎣ B a B s B b ⎤ ⎦ Φ = B Φ (46)
the virtual work principle becomes:

L 0 δΦ T B T K s B Φdx -w exter = 0 (47)
The section stiffness matrix presented hereafter is valid for homogeneous and non homogeneous sections even if the chosen axes are not the principal ones [START_REF] Guedes | A fibre timoshenko beam element in castem[END_REF]:

K e = L 0 B T K s Bdx (48)
The numerical implementation of a multifiber beam is similar to a classical beam with the main difference that further loops are needed in the section level (scanning all the fibers) in order to construct the section stiffness matrix

K s [GPP94], [KM05], [ MKRC06].
In the following chapter, a case study is presented on a reinforced concrete viaduct considering soil-structure interaction [START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF]. For this, a multifiber Timosheko beam is coupled with macro-elements [START_REF] Grange | A macro-element to simulate 3d soil-structure interaction considering plasticity and uplift[END_REF], [START_REF] Grange | A macro-element to simulate dynamic soil-structure interaction[END_REF], [START_REF] Grange | Simplified modeling strategies for soil-structure interaction problems: The macro-element concept[END_REF]. The reader can find other applications of the multifiber beam concept in the recent literature: non linear shear [START_REF] Combescure | A fibre model accounting for transverse shear in castem[END_REF], non linear torsion [START_REF] Mazars | Using multifiber beams to account for shear and torsion: Applications to concrete struc-tural elements[END_REF], shaking Table 1: Viaduct: characteristics of deck cross section.

A(m 2 ) I x (m 4 ) I y (m 4 ) J(m 4 ) 1.11 0.13 2.26 2.39

Finite element mesh

A finite element model using multifiber beams and concentrated masses is chosen to reproduce the structure (figure 8). The piers are at first considered fixed at the base. The mass and rotational inertia details are given in Table 2 [START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF].

Figure 8: Viaduct: finite element mesh. Non linear Timoshenko multifiber beam elements are used to reproduce the behaviour of the piers [START_REF] Kotronis | Simplified modelling strategies to simulate the dynamic behaviour of r/c walls[END_REF], [START_REF] Mazars | Using multifiber beams to account for shear and torsion: Applications to concrete struc-tural elements[END_REF]. Six (6) elements are used for the piers P1 and P3 and nine (9) elements for the pier P2. Forty (40) concrete fibers and eighty (80) steel fibers are assumed in the sections, (figure 9). The deck is simulated using elastic linear beam elements. Calculations are made with FEDEASLab, a finite element MATLAB toolbox [START_REF] Filippou | Fedeaslab, getting started guide and simulations examples[END_REF].

Material parameters

A damage model with two scalar variables, one in compression and one in tension is adopted for concrete [START_REF] Borderie | Phénomènes unilatéraux dans un matériau endommageable: modélisation et application l'analyse des structures en béton[END_REF]. The model is able to reproduce the unilateral effect, the permanent strains and the stiffness recovery. A modified version of the classical Menegotto-Pinto model [START_REF] Menegotto | Method of analysis of cyclically loaded reinforced concrete plane frames including changes in geometry and nonelastic behaviour of elements under combined normal force and bending[END_REF] with an isotropic hardening is used for steel. It is worth noting that as the tests are pseudo-dynamic, the damping coefficient adopted in the numerical simulations has to be small [START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF].

Loading sequence

The accelerations imposed at the base of the structure derive from a synthetic accelerogram consistent with a 5% damping response specturm selected according to Eurocode 8 for a soil of class B. The peak of accelerations is situated at 0.35g ("weak" earthquake). A second similar accelerogram (dilated) is also imposed at the base of the structure. Its peak of acceleration is equal to 0.7g ("strong" earthquake) [START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF]. One can clearly see that despite the small number of degrees of freedom of the finite element model the non linear behaviour of the viaduct is reproduced quite satisfactorily. Not only the peaks in both directions are well reproduced but the frequency content of the response is correctly matched.

Experimental versus numerical results: dynamic analysis

Dynamic analysis considering soil-structure interaction

Two modelling strategies are studied hereafter to take into account soil-structure interaction. The first uses the macro-element approach [GKM09a], [START_REF] Grange | A macro-element to simulate dynamic soil-structure interaction[END_REF], [START_REF] Grange | Simplified modeling strategies for soil-structure interaction problems: The macro-element concept[END_REF] and the second linear elastic springs applied at the base of each pier. The elastic stiffness of the springs is calibrated such as that they accumulate the same energy as the non-linear macro-element [START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF]. The three types of boundary conditions are denominated hereafter as follows: linear springs (EL), macro-element (ME) and fixed (Fixed).

The results for the weak earthquake are presented in figure 12 for a class C soil. The predicted numerical behaviour of the viaduct differs depending on the assumed boundary conditions. The displacements are strongly amplified for the case of the structure resting on the macro-element and on the linear elastic springs. The results are more pronounced for the internal forces at the base of the piers (moments and shear forces). Loads on the structure are significantly reduced for the case of the macro-element.

Results obtained with the elastic linear springs are similar to the ones found for the fixed piers [START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF].

The limits of the classical engineering approach based on elastic linear springs are thus evident. For the case of the reinforced concrete viaduct internal forces and displacements are higher than the ones obtained using the macro-element, which allows a more appropriate description of the non linear behaviour of the foundation soil system.

Figure 12: Viaduct -soil-structure interaction: comparison of the displacements, moments and shear forces for the weak motion and for a class C soil [START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF].

Conclusions

In this course, the formulation of a multifiber beam element was presented in detail. Adopting an Euler Bernoulli or a Timoshenko kinematic assumption, the element is able to reproduce the non linear behaviour of composite structures. A case study on a reinforced concrete viaduct subjected to earthquake loadings showed the performance of the approach. Combined with macro-elements [START_REF] Grange | Simplified modeling strategies for soil-structure interaction problems: The macro-element concept[END_REF], it can take into account in a efficient, fast and robust way soil-structure interaction phenomena.

Figure 1 :

 1 Figure 1: Beam [Bit13]

Figure 4 :

 4 Figure4: A 2 node finite element beam[START_REF] Bitar | Un nouvel élément fini poutre pour la modélisation des structures sous sollicitations complexes[END_REF]

Figure 6 :Figure 7 :

 67 Figure 6: Viaduct: plan view (scale 1:2.5) [PVP + 96].

Figure 9 :

 9 Figure 9: Viaduct: details of the multifiber beam element mesh (piers P1, P2 et P3).

FiguresFigure 10 :

 10 Figures10 and 11show the comparison between the experimental and the numerical results of the dynamic analysis considering the piers fixed at the base. The two earthquakes (weak and strong) are imposed. The figures show the evolution with time of the shear forces at the base and the lateral displacements at the top of the piers P1, P2 and P3[START_REF] Grange | The effects of soil-structure interaction on a reinforced concrete viaduct[END_REF].

Figure 11 :

 11 Figure 11: Viaduct -fixed base: comparison between experimental and numerical displacements and shear forces for the strong earthquake [GBKT11].

  Eε x , σ xy = Gγ xy , with G the shear coefficient) we obtain:

	Equation (15) becomes:
	0	L	(F x	d dx	δU x + F y	d dx	δβ y + M z	d dx	δΘ z )dx -w exter = 0	(17)
	Using Hooke's law (σ x =
				Shear force: F y =			(16)

S σ xy dS

Bending moment :

M z = -S yσ x dS

table tests

 tests PVP + 96]). Inertial forces were calculated numerically and imposed to the model piers through actuators by applying the adequate displacements. Details of the deck and piers are given (scaled) in figures 7(a) and 7(b). Piers are made of reinforced concrete and present hollow rectangular section shapes. The characteristics of the section of the deck are given in Table1.

	[KRM05],
	[INK + 08], [GKM09c], retrofitting with fiber reinforced polymers [DMKP13]. . .
	4 Case study: A reinforced concrete viaduct
	5 Description of the structure
	A 1:2.5 scaled viaduct was tested pseudo-dynamically in ELSA laboratory (JRC Ispra)
	(figure 6, [

Table 2 :

 2 Viaduct: masses and rotational inertia.Mass M (kg) Rotational inertia I x (kg.m 2 ) Rotat. inertia I z (kg.m 2 )

	M A	27.5	285	234
	M B	32	287	271
	M C	34	288	322
	M D	13.75	143	117


Hai Hoang

G Galliero

Swelling/Shrinkage Induced by Shear in Narrow Pores

published or not. The documents may come

Swelling/shrinkage induced by shear in narrow pores

I. Introduction

Understanding the mechanical response of a tight porous medium in which pores are filled partially or completely by a fluid and are subject to an external phenomena, e.g. external stress, plays an important role in many engineering problems: low permeability reservoirs (Shale, Coal seams) in petroleum engineering, activated carbons/zeolites in chemical engineering processes, cement paste in civil engineering, etc. (Coussy 2010[START_REF] Israelachvili | Intermolecular and Surface Forces[END_REF][START_REF] Karniadakis | Rheological Evidence for a Dynamical Crossover in Polymer Melts via Nonequilibrium Molecular Dynamics[END_REF][START_REF] Anderson | Clay swelling -A challenge in the oilfield[END_REF]. In such situations where the pore sizes are of the order of one to ten molecular sizes, the main problem is how to take into account the solid-fluid interaction in the mechanical response of the whole porous medium composed of the fluid and the solid matrix [START_REF] Anderson | Clay swelling -A challenge in the oilfield[END_REF][START_REF] Hensen | Why Clays Swell[END_REF][START_REF] Vandamme | Adsorption and strain: The CO2-induced swelling of coal[END_REF][START_REF] Pijaudier-Cabot | Revisiting poromechanics in the context of microporous materials[END_REF][START_REF] Malikova | Dynamics in Clays -Combining Neutron Scattering and Microscopic Simulation[END_REF][START_REF] Brochard | Poromechanics of microporous media[END_REF].

In a tight porous medium, the surface effects due to the fluid adsorbed on the solid surface and the molecular packing due to confinement can be dominant over the volume effects [START_REF] Israelachvili | Intermolecular and Surface Forces[END_REF][START_REF] Karniadakis | Rheological Evidence for a Dynamical Crossover in Polymer Melts via Nonequilibrium Molecular Dynamics[END_REF][START_REF] Schoen | Computer simulation of condensed phases in complex geometries[END_REF][START_REF] Hansen | Theory of simple liquid[END_REF]. Among others, this can lead to the swelling of the porous medium [START_REF] Anderson | Clay swelling -A challenge in the oilfield[END_REF][START_REF] Hensen | Why Clays Swell[END_REF][START_REF] Vandamme | Adsorption and strain: The CO2-induced swelling of coal[END_REF][START_REF] Pijaudier-Cabot | Revisiting poromechanics in the context of microporous materials[END_REF][START_REF] Malikova | Dynamics in Clays -Combining Neutron Scattering and Microscopic Simulation[END_REF][START_REF] Brochard | Poromechanics of microporous media[END_REF]. To take into account these effects poromechanics approaches were proposed in the literature. [START_REF] Pijaudier-Cabot | Revisiting poromechanics in the context of microporous materials[END_REF] revisited classical poromechanics in the context of microporous materials with a continuous pore size distribution. To account for the surface effects, the constitutive equations were re-formulated by introducing an apparent porosity and an interaction free energy that are related to the Gibbs adsorption isotherm. In a recent study, [START_REF] Brochard | Poromechanics of microporous media[END_REF] derived poroelastic constitutive equations for a generic porous medium, in which the surface effects are characterized by the amount of fluid adsorbed that was shown to depend on both the fluid bulk pressure and the strain of the medium.

It should be noticed that these constitutive equations for microporous medium always assume that the adsorption stress, i.e. the pressure (normal to the solid phase) of the fluid in the pore, sometimes called solvation pressure, and the density of the fluid in the pore are independent on the deviatoric strain. Thus, the volume of a pore is assumed to be unchanged when the porous medium is sheared. As we will show later, this can be inadequate when dealing with some microporous medium in which the solid phase is highly structured and the solid-fluid interactions are strong. This is simply due to the fact that in a slit nanopore configuration, the relative ordering between the solid surfaces has a significant effect on the properties of the fluid in the pore (Schoen 1993, Hoang and[START_REF] Hoang | Shear behavior of a confined thin film: Influence of the molecular dynamics scheme employed[END_REF].

Thus, in this work using extensive molecular dynamics (MD) simulations we explore the volumetric deformation of a simple saturated (by a liquid) slit pore induced by the solid displacement in the direction parallel with the solid-fluid interface. To do so, we analyze the swelling/shrinkage (i.e. the gap between the two flat walls) due to a relative displacement of the solid walls, see Fig. 1.

II. MOLECULAR DYNAMICS SIMULATIONS

All fluid and solid molecules have been modeled as spherical ones. All interactions are described by a classical truncated Lennard-Jones (LJ) 12-6 potential between pairs of particle:

( )      > ≤               -       = c c LJ r r r r r r r U if 0 if 4 6 12 σ σ ε (1)
where r is the distance between the two particles, ε is the potential depth, σ is the particle diameter, and r c is the cut-off diameter (=3.5σ in this work). In the following the variables have been expressed in usual LJ dimensionless units, noted with a star as superscript. The solid walls are formed by atoms arranged on a faced centered cubic (FCC) lattice with a lattice size 6 . 1 * = a . It is assumed that their structures are unchanged during the swelling/shrinkage, i.e. the Biot coefficient in the conventional constitutive equations is equal to one [START_REF] Biot | General Theory of Three Dimensional Consolidation[END_REF][START_REF] Coussy | Structure and solvation forces in confined films: Linear and branched alkanes[END_REF]). The walls are separated by a distance * Gap L .

To perform the MD simulations, we have employed the constant density Grand Canonical like Molecular Dynamics method proposed by [START_REF] Hoang | Grand canonical-like molecular dynamics simulations: Application to anisotropic mass diffusion in a nanoporous medium[END_REF] which is similar to the one proposed by Gao et al. (1997). It consists in simulating explicitly the fluid in the pore in direct contact with reservoirs maintained at the desired density. For that study we have employed 7 . 0 * = ρ and temperature 1 * = T which corresponds to a liquid state in bulk conditions for the cutoff employed in this work [START_REF] Galliero | Molecular dynamics study of the repulsive form influence of the interaction potential on structural, thermodynamic, interfacial, and transport properties[END_REF]. Figure 1 shows a sketch of the simulation box that contains both fluid and solid particles. Periodic boundary conditions (PBC) have been applied on all three directions. The configuration has been chosen so that the fluid out of the pore space is inappreciably influenced by a small change in the pore studied.

To study the swelling/shrinkage of the pore, the solid walls forming the pore are not allowed to move in the y direction but only in the x direction during the simulations. To model the solid displacement in the direction parallel to the solidfluid interface, the walls are prescriptively moved in the z direction, see Fig. 1. The equation of motion of the centers of mass of the solid walls in the x direction is governed by the classical equation of motion as:

( ) j j 1 1 j j LJ S S S F S i i N i N j i i r x r r U x m N ∑∑ = = ∂ ∂ - = ′ ′ (2)
where S N , S m and S x are the number of solid molecules, the mass of a solid molecule and the coordinate of the center of mass in the x direction for each solid wall respectively. To avoid as much as possible the effects due to the finite size of the two solid walls in the y direction, Eq. ( 2) is computed using only the solid molecules in the sampling region as described in [START_REF] Hoang | Grand canonical-like molecular dynamics simulations: Application to anisotropic mass diffusion in a nanoporous medium[END_REF].

III. RESULTS and DISCUSSIONS 3.1 Swelling/Shrinkage induced by walls displacement

As mentioned previously, the use of the conventional and recently modified constitutive equations leads to an unchanged pore width when the solid walls are moved in the direction parallel to the solid-fluid interfaces. Thus the shear is not assumed to influence the pressure (normal to the walls) in the fluid. However, in a narrow pore for which its size is comparable to the typical distance of interaction between molecules, the relative structural ordering between the solid walls has a significant effect on the fluid in the pore (Schoen 1993, Hoang and[START_REF] Hoang | Shear behavior of a confined thin film: Influence of the molecular dynamics scheme employed[END_REF], in particular on the normal pressure of the confined fluid. This means that a relative displacement of the two solid walls can modify the normal pressure of the fluid. Hence, a pore can swell or shrink if the solid walls are translated of a given distance in the direction parallel to the solid-fluid interface. That effect is named "static" swelling/shrinkage in the following.

To explore such static swelling/shrinkage of a pore, we have computed the variation of the pore width with the amplitude of the solid walls displacement * z ∆ starting from an out-of-registry initial configuration of the solid walls. To do so, we have performed MD simulations in which, once the system has reached the equilibrium, the solid walls are moved in opposite direction of a distance 2

Dynamic Swelling/Shrinkage

To extend the results presented in the previous section, it seems interesting to explore the swelling/shrinkage of a slit pore when the solid walls are moved at a constant velocity, i.e. to induce a dynamic swelling/shrinkage. To do so, we have performed MD simulations in which, once the system has reached the equilibrium, the solid walls are moved at a constant velocity in the direction parallel to the solidsurface interface, see Fig. When the velocity is sufficiently small, it is consistent to assume that the state of the system is completely re-equilibrated at each moment. In other words, the pore size at any time should depend only on the instantaneous relative structural ordering between the solid surfaces. The instantaneous pore size should so vary between the values corresponding to the in-registry and out-of-registry configurations. Thus, we have carried MD simulations for different wall velocities to check this assumption and to explore the behavior of the dynamic swelling/shrinkage for large shearing rate. γ . The latter trend is mainly due to the fact that the pressure of a fluid usually increases when the wall velocity increases (Hoang andGalliero 2013, Heyes 1986). It should be noticed that * c γ is one order of magnitude lower than the shear rate associated to the shear thinning threshold of a bulk LJ liquid (Krögerand andHess 2000, Galliero and[START_REF] Galliero | Shear viscosity of the Lennard-Jones chain fluid in its gaseous, supercritical, and liquid states[END_REF].

Figure 5 displays the variation of the pore size with time for different walls velocities. As expected, all cases exhibit an oscillatory behavior with a period

* Wall * * 2V a = θ
. This behavior is simply due to the periodic structure of the solid surfaces. It is interesting to notice that the values of peaks and valleys are nearly equal to the equilibrium ones corresponding to the in-registry and out-of-registry configurations for low walls velocities, whereas it is not the case for high velocities, see Figs. 3, 4 and5. In other words, the walls velocity has a noticeable effect on the instantaneous size of a pore, which is not so obvious when looking only to the average size, see Figs. 4 and5. In particular, the amplitude of the oscillations is reduced when the walls velocity increases. Such behaviors can be understood in terms of re-equilibration of the fluid in the pore undergoing shear. For low enough shear rate the pore size is governed only by the effect of the relative structural ordering between the solid surfaces, not by the wall velocity. When the shear rate is sufficiently large, the fluid in the pore does not have enough time to re-equilibrate and so the effect of the relative structural ordering is reduced.

IV. CONCLUSIONS

The poro-mechanics approaches have shown to be valuable to describe the deformation of a porous media partially or fully saturated by a fluid. In its present form, such theory assumes that the volume of a pore is unchanged when the porous medium is sheared. However, this is incorrect at the pore scale for systems in which the fluid is highly structured due to a strong solid-fluid interaction. Thus, in this work, using molecular dynamics simulations, we have explored the volumetric deformation of a simple slit nano-pore induced by shear.

First, the swelling/shrinkage of a pore induced by a displacement of the solid walls (immerged in a fluid reservoir) of a given distance in the direction parallel to the solid-fluid interface has been investigated. Results have shown that, in narrow pores, the average pore size is changed by such solid displacement, i.e. there is a shear-induced swelling/shrinkage. This effect is simply due to the fact that the normal pressure of the simulated fluid appreciably varies with the relative structural ordering between solid surfaces (going from out of registry to in registry configurations).

Then, we have explored the situation in which solid walls are moved at a constant velocity. It has been found that the shear rate employed has only a slight effect on the average pore size, whereas it has a noticeable effect on the instantaneous pore size. More precisely, while all cases exhibit an oscillatory behavior of the instantaneous pore size, which is simply due to the periodic structure of the solid surfaces, the amplitude of the oscillations is reduced when the walls velocity increases. Such behaviors can be understood in terms of re-equilibration of the fluid in the pore undergoing shear.

Figure 1 :

 1 Figure 1: (a): Side view of the simulation cell used in the Grand Canonical like Molecular Dynamics. (b): A schematic representation of the solid displacement in the direction parallel to the solid-fluid interface.

Figure 2 :Figure 2

 22 Figure 2: Variation of the pore size for a given solid displacement * z ∆ .

Figure 3 :

 3 Figure 3: (a) Variation in the normal pressure with gap size. Circles correspond to the in-registry configuration. Diamonds correspond to the out-of-registry configuration. Solid line represents the pressure in the reservoirs. (b) Variation in the normal pressure with the relative structural ordering between the solid surfaces Latt L z ∆ for

Figure 4 :

 4 Figure 4: Variation in the average size of the pore with the shear rate * * * 2 Equi Gap Wall L V -

Figure.

  Figure 4 shows the dependence of the average pore size for different shear rates * * * 2 Equi Gap Wall L V -

Figure 5 :

 5 Figure 5: Variation in the pore size with time when the solid walls are moved at a constant velocity. Red curves correspond to 001 . 0 2 * *

Acknowledgments

This work has been supported by the "Failflow" project which is funded by an advanced grant from the European Research Council. We gratefully acknowledge computational facilities (in link with a "Conseil Regional d'Aquitaine" project) provided by Pau University, by the I2M and the MCIA in Bordeaux University.


Omer Nguena Timo
email: omer.nguena-timo@crim.ca

Dimitri Prestat
email: prestat.dimitri@courrier.uqam.ca

Antoine Rollet
email: antoine.rollet@labri.fr

Multiple Mutation Testing for Timed Finite State Machine With Timed Guards and Timeouts

The problem of generating tests detecting all logical and timing faults which can occur in real-time systems is challenging; this is because the number of (timing) faults is potentially too big or innite. As a result, it might be time consuming to generate an important number of adequate tests. The traditional model based testing approach considers a fault domain as the universe of all machines with a given number of states and input-output alphabet while mutation based approaches dene a list of mutants to kill with a test suite. In this paper, we combine the two approaches by developing a mutation testing technique for realtime systems represented with deterministic timed nite state machines with timed guards and timeouts(TFSM-TG). In this approach, fault domains consisting of fault-seeded versions of the specication (mutants) are represented with non-deterministic TFSM-TG. The test generation avoids the one-by-one enumeration of the mutants and is based on constraint solving. We present the results of an empirical proof-of-concept implementation of the proposed approach.

Introduction

This paper deals with mutation testing of timed systems. Traditional modelbased testing approaches consist of using the model of the specication as a base for the generation of a test suite, then applying it on the implementation of the system under test. This is a black-box testing, meaning that the implementation is considered as unknown, but with the assumption that it may be described by a formal model. It is the so-called test hypothesis. Since it is generally impossible to generate an exhaustive test suite, several approaches have been proposed to limit the size of test cases while providing enough condence.

This work was partially supported by MEI (Ministère de l'Économie et Innovation) of Gouvernement du Québec.

Concerning FSM model-based testing, one may consider coverage criteria of the specication: the test may be satisfying if for instance all the transitions or all the states are covered by the test suite. Another approach consists of obtaining a sucient fault coverage: the test suite is guaranteed to detect all faults of a specic nature, e.g. output faults (the output provided by a transition is not correct) or transfer faults (the arrival state of a transition is not correct). As an example, the so called W-method [START_REF] Chow | Testing software design modeled by nite-state machines[END_REF] is known to detect both transfer and output faults whereas TT-method [START_REF] Naito | Fault detection for sequential machines by transitiontours[END_REF] 4 focuses only on output faults. Gener- ally, fault coverage based methods consider the traditional fault domain as the universe of all machines with a specic alphabet and a given number of states.

On the other side, mutation testing is a technique originally developed to verify if a given test suite has a satisfying detecting power. In model-based testing, the principle consists of applying a certain number of small variations, called mutations, on the specication, then to check if these mutations are detected by the test suite. It is generally called killing the mutants. Note that mutation testing techniques are also used in code based testing, but in this case, the mutations are applied directly in the source code. Dealing with real-time systems increases the diculty of testing. For instance, faults related to timing errors should also be detected by a test method. For timed model based testing, many dierent models and conformance relations have already been proposed in the past. Some of them use a discrete representation of the time in the specication model, such as synchronous models [START_REF] Bousquet | Lutess: A specication-driven testing environment for synchronous software[END_REF][START_REF] Raymond | Automatic testing of reactive systems[END_REF][START_REF] Marre | Test sequences generation from LUSTRE descriptions: GATEL[END_REF], some other methods extend the so-called ioco theory [START_REF] Tretmans | Test generation with inputs, outputs, and repetitive quiescence[END_REF] by dening a new conformance relation adapted to timed automata [START_REF] Alur | A theory of timed automata[END_REF] with inputs and outputs [START_REF] Larsen | Testing real-time embedded software using Uppaal-Tron: An industrial case study[END_REF][START_REF] Krichen | Conformance testing for real-time systems[END_REF][START_REF] Mikucionis | T-Uppaal: Online model-based testing of real-time systems[END_REF][START_REF] Bertrand | O-line test selection with test purposes for non-deterministic timed automata[END_REF][START_REF] Nguena Timo | Conformance testing of variable driven automata[END_REF][START_REF] Bohnenkamp | Timed testing with TorX[END_REF]. Another category extends the FSM based testing theory using timed extensions of the FSM model and an adapted fault model [START_REF] Lallali | Timed modeling of web services composition for automatic testing[END_REF][START_REF] Merayo | Formal testing from timed nite state machines[END_REF][START_REF] Derderian | Aiding test case generation in temporally constrained state based systems using genetic algorithms[END_REF][START_REF] Zhigulin | FSM-Based Test Derivation Strategies for Systems with Time-Outs[END_REF].

The work in [START_REF] El-Fakih | Testing Timed Finite State Machines with Guaranteed Fault Coverage[END_REF] proposes a method to generate a test suite with guarantee fault coverage for partial deterministic TFSM extended with timed guards only. TFSM with timed guards and timeouts (TFSM-TG) can express timed behaviors which cannot be expressed with TFSM with timeouts only [START_REF] Bresolin | Deterministic Timed Finite State Machines: Equivalence Checking and Expressive Power[END_REF]. More recently [START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF] proposes a method to generate tests detecting all the nonconforming implementations of an initialized TFSM with timeouts and TFSM-TG. The implementations belong to the traditional fault domain and they can have more states than the specication, up to a given number. The tests are generated by applying the W-method on an (possible huge) abstraction of TFSM-TG with classical FSM and transforming the resulting abstract tests into timed tests.

In this paper, we combine fault coverage based approaches and mutation testing applied to real-time systems. We propose a mutation testing technique for real-time systems represented with complete and deterministic TFSM-TG. We derive complete test suites from mutation machines. A mutation machine can be designed to represent the traditional fault model; but it can also represent customized fault models as well. Customized fault models could be parts of the traditional fault model and they can be covered with a reduced number of tests. The tests generated for the traditional fault model remain valid to test customized fault models. However they could suer from redundancy and could include useless tests. Running implementations with tests is time consuming.

Then, reducing the number of tests to apply is really useful. This is the main motivation of our work. Our test generation approach is based on constraint resolution. It uses the distinguishing automaton of the specication and the mutation machine to determine revealing combs. They characterized the undetected mutants and serve to encode them with Boolean formula. We use a solver to check the satisability of the formula and conclude about the completeness of test and generate new tests if needed. We propose an abstraction of timed states of TFSM-TG and we use it to build the distinguishing automaton. This paper extends the approach proposed in [START_REF] Petrenko | Multiple Mutation Testing from FSM[END_REF][START_REF] Nguena Timo | Multiple Mutation Testing from Finite State Machines with Symbolic Inputs[END_REF] in a timed context. This is the main contribution of our paper. It is adapted to be applied on the TFSM-TG model. This contribution includes the abstraction of timed states, the denition of the distinguishing automaton for TFSM-TG, and the representation of the transitions passed during executions of TFSM-TG with the so-called combs.

The paper is organized as follows. The next section introduces the general theoretical background, a fault model for TFSMs-TG and the coverage of fault models with complete test suites. Section 3 characterizes detected mutants with revealing combs. In Section 4 we introduce an abstraction for executions of TFSM-TG and we dene the distinguishing automaton. Section 5 presents an encoding of undetected mutants with Boolean formulas. We present a method for the test analysis and a method for the complete test suite generation in Section 6; we also present the results of an empirical proof-of-concept tool. We conclude the paper in Section 7.

Preliminaries

A timed guard π is an interval π = l π , u π where l π is a non-negative integer, while u π is either a non-negative integer or ∞, l π ≤ u π , and the symbols and represent ] or [. For example, the interval [0, 5[ contains all the non-negative real numbers smaller than 5. We let Π denote the set of timed guards. Given a timed guard π = l π , u π and a real number x, we dene π -x = max(0, l π -x), max(0, u π -x) .

2.1 Timed FSM with timed guards and timeouts Denition 1. A timed nite state machine with timeouts and timed guards [START_REF] Bresolin | Deterministic Timed Finite State Machines: Equivalence Checking and Expressive Power[END_REF] (TFSM-TG) is a 6-tuple S = (S, s 0 , I, O, λ S , ∆ S ) where S, I and O arenite non-empty set of states, inputs and outputs, respectively, s 0 is the initial state, λ S ⊆ S × I × Π × O × S is an input/output transition relation and ∆ S ⊆ S × N ≥1 ∪ {∞} × S is a timeout transition relation.

Remark that we allow dening multiple timeout transitions in states of TFSM-TG, in the opposite of [START_REF] Bresolin | Deterministic Timed Finite State Machines: Equivalence Checking and Expressive Power[END_REF]. Later this will be used to compactly represent sets of implementations of a TFSM-TG specication. S in state s res an input/output transition (s, i, π, o, s ) to reach the transition's target state s and produces output o if input i is applied in transition's starting state s when timed guard π is respected; π expresses the minimal and the maximal delays in s to re the transition. It res a timeout transition (s, δ, s ) ∈ ∆ S dening timeout δ in s and reaches s if no input is applied in s before δ expires. Let δ max (s) denote the maximal (possibly innite) timeout dened in state s. We require that the maximal nite constant in guards of the transitions dened in every state s should be smaller than δ max (s). Because multiple timeouts can be dened in the same state, S necessarily res a timeout transition dening δ max (s) if no input is applied at s before δ max (s) expires. A clock measuring the amount of the time elapsed in every state is implicitly reset when transitions are red.

A timed state of TFSM-TG S is a pair (s, x) ∈ S × R ≥0 where s ∈ S is a state of S and x ∈ R ≥0 is the current value of the clock and x < δ for some δ ∈ N ≥1 ∪ {∞} such that (s, δ, s ) ∈ ∆ S . The initial timed state of S is (s 0 , 0).

An execution step of S in timed state (s, x) corresponds either to the time elapsing or the ring of an input/output or timeout transition; it is permitted by a transition t of S. Formally, a tuple ((s, x), a, t, (s , x ))

∈ (S × R ≥0 ) × (((I × O) ∪ R ≥0 ) × (λ S ∪ ∆ S )) × (S × R ≥0
) is an execution step if it satises one of the following conditions:

(timeout) t = (s, δ, s ) ∈ ∆ S , a ∈ R ≥0 , x + a = δ and x = 0 (time-elapsing) t = (s, δ, s ) ∈ ∆ S , a ∈ R ≥0 , x + a < δ, x = x + a and s = s (input/output) t = (s, i, l, u , o, s ) ∈ λ S , x ∈ l, u , a = (i, o) with (i, o) ∈ I × O and x = 0
In the time-elapsing step, the target state s of t and s can be dierent. This is because the timeout δ is not expired and t is not red; δ max (s) is never exceeded.

An execution of S in timed state (s 0 , x 0 ) is a nite sequence of steps e = stp 1 stp 2 . . . stp n with stp k = ((s k-1 , x k-1 ), a k , t k , (s k , x k )), k ∈ [1, n] such that stp 1 is not an input/output step and stp k is an input/output step implies that stp k-1 is a time-elapsing step for every k ∈ [1..n]. If needed, the elapsing of zero time unit can be inserted before input/output steps which are not immediately preceded with a time-elapsing step. The timed input/output sequence of execution e is the sequence

(i 1 , o 1 )d 1 (i 2 , o 2 )d 2 . . . (i l , o l )d l in ((I × O) × R ≥0 ) * such that (i 1 , o 1 )(i 2 , o 2 ) . . . (i l , o l )
is the sequence of input/output pairs occurring in the execution and d 1 d 2 . . . d l ∈ R l ≥0 is a sequence of non-negative real numbers and t k is the amount of the time elapsed since the occurrence of i k-1 or the beginning of the execution if no input has occurred. l is smaller than the number of steps in the execution. The timed input sequence and the timed output sequence of the execution e are inp(e

) = i 1 d 1 i 2 d 2 . . . i l d l and out(e) = o 1 d 1 o 2 d 2 .
. . o l d l , respectively; they are said to be applicable and produced in (s, x), respectively. We denote by inp(s, x) the set of timed input sequences applicable in (s, x). Given a timed input sequence α, let out S ((s, x), α) denote the set of all timed output sequences which can be produced by S when α is applied in s, i.e., out((s, x), α) = {out(e) | e is an execution of S in (s, x) and inp(e) = α}. We let Exec S and Exec S (α) denote the set of executions of S and the set of executions with the timed input sequence α.

Transitions starting in the same state are called compatible if they are timeout transitions or have the same input and the intersection of their timed guards is non-empty. A TFSM-TG S is deterministic (DTFSM-TG) if it has no compatible transition; otherwise, it is non-deterministic. S is initially connected if every state of S is part of a reachable timed state. Let λ S (s, i) denote the set of input/output transitions dened in state s with input i. S is complete if the union of the timed guards of the transitions in λ S(s,i) equals [0, ∞[ for every (s, i) ∈ S × I; it implies that every i ∈ I is the input of at least one input/output step from every reachable timed state of S. Note that inp(s, x) = (I × R ≥0 ) * for every timed state (s, x) of a complete machine S.

We dene distinguishability and equivalence relations between timed states of complete TFSMs-TG. Then, we extend these relations to TFSM-TG. Similar notions were introduced in [START_REF] Zhigulin | FSM-Based Test Derivation Strategies for Systems with Time-Outs[END_REF]. Intuitively, timed states producing different timed output sequences in response to the same timed input sequence are distinguishable. Formally, let (s, x s ) and (m, x m ) be the timed states of two complete TFSMs-TG dened on the same input and output sets. Given a timed input sequence α, (s, x s ) and (m, x m ) are distinguishable with α denoted (s, x s ) α (m, x m ), if the sets of timed output sequences in out((s, x s ), α) and out((m, x m ), α) dier; otherwise they are equivalent and we write (s, x s ) (m, x m ), i.e., if the sets of timed output sequences coincide for every timed input sequence α. Two TFSM-TG are equivalent if their initial timed states are indistinguishable; otherwise they are distinguishable with a timed input sequence.

Henceforth the TFSMs-TG are complete and initially connected.

Mutants and fault model

Let S = (S, s 0 , I, O, λ S , ∆ S ) be a complete DTFSM-TG, called the specication machine. A mutant is a variant of the specication and represents a possibly faulty implementation that should be detected by tests; it is also a deterministic and complete TFSM-TG. A mutant can have fewer or more states than the specication, up to a specied bound. The faults can be introduced by performing combinations of the following atomic mutation operations: changing the target state of a transition (transfer fault), changing the output of a transition (output fault), changing a timeout, merging the timed guards of transitions, splitting the timed guard of a transition, adding an extra-state. Similar operations were dened in [START_REF] Vega | Model-based mutation operators for timed systems: A taxonomy and research agenda[END_REF][START_REF] Nilsson | Test case generation for mutation-based testing of timeliness[END_REF] for testing other types of timed machines. Our operations are adapted to TFSM-TG. We compactly represent mutants with a mutation machine. Intuitively, mutants and the specication are all sub-machines of a global mutation machine describing the fault model.

TFSM-TG S = (S, s 0 , I, O, λ S , ∆ S ) is a sub-machine of TFSM-TG M = (M, m 0 , I, O, λ M , ∆ M ) if S ⊆ M , s 0 = m 0 , λ S ⊆ λ M and ∆ S ⊆ ∆ M . Denition 2. A non-deterministic TFSM-TG M = (M, m 0 , I, O, λ M , ∆ M ) is a mutation machine of S if S is a sub-machine of M. Transitions in λ M but not in λ S or in ∆ M but not in ∆ S are called mutated.
M. Let ∆ M (m) denote the set of timeout transitions dened in m. For a state m, every mutant denes exactly one timeout transition of ∆ M and a subset z mi of transitions of λ M (m, i), for every (m, i) ∈ M × I. The subset z mi satises the following three cluster conditions: [START_REF] Alur | A theory of timed automata[END_REF] it is non-empty, (2) the transitions in z mi are not compatible with each other and (3) the union of their timed guards equals [0, ∞[. We let Z mi = {z 1 mi , z 2 mi , . . . , z l mi } be the maximal set of subsets of λ M (m, i) such that z k mi satises the three cluster conditions, for each k = 1...l.

The number of mutants is |Mut(M)| = (m,i)∈M ×I |Z mi |× m∈M |∆ M (m)|-1.
Faults in mutants are represented with mutated transitions which can be viewed as alternatives for transitions of the specication. Note that changes of delays to re input/output transitions cannot be expressed if the specication is a TFSM with timeouts only; this is because the timed guard for every input in a TFSM with timeouts is always [0, ∞[. Some executions of non-deterministic sub-machines of M are not executions of any mutants or the specication. We can show that

P ∈Mut(M)∪{S} Exec P (α) ⊆ Exec M (α).
A transition t is suspicious in M if M denes another transition t compatible with t. In other words, t and t specify possible behaviors of dierent mutants. A transition of the specication is called untrusted if it is suspicious in the mutation machine; otherwise, it is trusted. Every trusted transition is dened in each mutant. The set of suspicious transitions of M is partitioned into a set of untrusted transitions all dened in the specication and the set of mutated transitions undened in the specication. We let Susp O denote the set of suspicious transitions in an artifact O containing transitions.

Figure 1b presents a mutation machine M 1a is deterministic; it denes all the trusted transitions and the un-

trusted transitions t 4 , t 8 , t 12 . ∆ M1 (s 1 ) = {t 2 }, Z s1a = {{t 1 }}, Z s1b = {{t 3 }}, ∆ M1 (s 2 ) = {t 8 , t 11 }, Z s2a = {{t 7 , t 9 , t 10 }}, Z s2b = {{t 4 }, {t 5 , t 6 }}, ∆ M3 (s 1 ) = {t 15 }, Z s3a = {{t 12 },
{t 13 , t 16 , t 17 }} and Z s3b = {{t 14 }}. Mut(M 1 ) contains seven mutants; two of them appear in Figure 2a and Figure 2b. The mutants are obtained from the specication by changing the behavior of suspicious transitions, which is done by replacing the untrusted transitions in the specication by other suspicious transitions they are compatible with. The mutant in Figure 2b has a dierent behavior in state s 3 for input a. The mutants and the specication have the same behavior in s 1 for input a; this is because t 1 is trusted.

Let P be a mutant with an initial state p 0 of the mutation machine M of S.

We use the equivalence relation to dene conforming mutants.

Denition 3. Mutant P conforms to S, if (p 0 , 0) (s 0 , 0); otherwise, it is nonconforming and a timed input sequence α such that (p 0 , 0) α (s 0 , 0) is said to detect P. P survives α if α does not detect P.

s1 start s2 s3 [0, ∞[a/x [t1] ∞ [t2] [0, ∞[b/x [t3] [0, ∞[b/x [t4] [5 , ∞ [a / x [t 7 ] 8 [t8] [2 , 5 [a / x [t 9 ] [0 , 2 [a / y [t 1 0 ] [0, ∞[a/y [t12] [0, ∞[b/x [t14] 5 [t15] (a) The specication machine S1 s1 start s2 s3 [0, ∞[a/x [t1] ∞ [t2] [0, ∞[b/x [t3] [0, ∞[b/x [t4] [0, 3]b/y [t5] ]3, ∞[b/x [t6] [ 5 , ∞ [ a / x [t 7 ] 8 [t8] [ 2 , 5 [ a / x [t 9 ] [ 0 , 2 [ a / y [t 1 0 ] 7 [t 1 1 ] [0, ∞[a/y [t12] [0, 3[a/y [t13] [0, ∞[b/x [t14] 5 [t15] [3, 5]a/y [t16] ] 5 , ∞ [ a / y [t 1 7
] (b) The mutation machine M1 Fig. 1: A mutation machine for a specication machine; s 1 is the initial state The set Mut(M) of all mutants in mutation machine M is called a fault domain for S. If M is deterministic and complete then M includes only the specication and Mut(M) is empty. A general fault model is the tuple S, , Mut(M) following [START_REF] Petrenko | Test Suite Generation from a FSM with a Given Type of Implementation Errors[END_REF][START_REF] Nguena Timo | Checking Sequence Generation for Symbolic Input/Output FSMs by Constraint Solving[END_REF]. The conformance relation partitions the set Mut(M) into conforming mutants and nonconforming ones which we need to detect. Denition 4. A test for S, , Mut(M) is a timed input sequence. A complete test suite for S, , Mut(M) is a set of tests detecting all nonconforming mutants in Mut(M).

We address the problem of checking the completeness of a test suite and the problem of generating a complete test suite for a fault model S, , Mut(M) , where the specication machine S is deterministic and complete and the mutation machine can have more states than S. Our approach consists in eliminating, from the fault model, the mutants detected by tests meanwhile avoiding their one-by-one enumeration. Undetected mutants will serve to generate new tests.

Revealing combs for characterizing detected mutants

We introduce combs to characterize mutants having common executions of the mutation machine since mutation machine includes all the mutants. The mutants dening all the transitions in a comb have common executions and can be detected with the same test. The comb for an execution e is the sequence of transitions permitting the steps in e; it is denoted by π e . The transitions in combs do not necessarily form paths in the state-transition diagram of M. This is because they contain non-red time-out transitions permitting time-elapsing steps when

s1 start s2 s3 [0, ∞[a/x [t1] ∞ [t2] [0, ∞[b/x [t3] [0, ∞[b/x [t4] [5 , ∞ [a / x [t 7 ] [ 2 , 5 [ a / x [t 9 ] [ 0 , 2 [ a / y [t 1 0 ] 7 [t 1 1 ] [0, 3[a/y [t13] [0, ∞[b/x [t14] 5 [t15] [3, 5]a/y [t16] ] 5 , ∞ [ a / y [t 1 7 ] (a) Nonconforming mutant P1 s1 start s2 s3 [0, ∞[a/x [t1] ∞ [t2] [0, ∞[b/x [t3] [0, ∞[b/x [t4] [5 , ∞ [a / x [t 7 ] 8 [t8] [ 2 , 5 [ a / x [t 9 ] [ 0 , 2 [ a / y [t 1 0 ] [0, 3[a/y [t13] [0, ∞[b/x [t14] 5 [t15] [3, 5]a/y [t16] ] 5 , ∞ [ a / y [t 1 7
] (b) Nonconforming mutant P2 1b; it is also an execution of S 1 . The transitions permitting the steps appear in brackets and also serve to identify the comb. The rst, third and fourth steps in e

e 1 = (s 1 , 0)0.5 [t2] (s 1 , 0.5)b/x [t3] (s 2 , 0)7 [t8] (s 2 , 7)0.5 [t8] (s 2 , 7.5)a/x [t7] (s 1 , 0) is an execution of M 1 in Figure

Distinguishing Automaton and Revealing Combs

The distinguishing automaton for a specication S and a mutation machine M is aimed to represent synchronous executions of S and the mutants in M. It will be used to identify unexpected executions of the mutants and the corresponding combs without enumerating the mutants one-by-one. Each of its states is composed of a timed state of S and a timed state of M. Each transition between two states of the automaton represents a synchronization between a transition of S and a transition of M. Because S and M could have an innite number of reachable timed states, a naive denition of the distinguishing automaton could include an innite number of the automaton's states. For this reason our denition is based on a new notion of abstract timed states for TFSM-TG.

Abstract timed state for TFSM-TG

We introduce a new notion of abstract execution which we show to be equivalent to the notion of execution described in Section 2. In abstract executions the range for clock values is nite while it is innite in executions. In particular, by increasing a clock value by one in a state which denes the innite timeout, the value of the clock will continuously increase up to a value which depends on the number of time-elapsing steps. It happens that after a certain threshold value, the exact value of the clock is not necessary to determine transitions which can be red. In state m, the threshold value is the maximal nite integer used in input-output and timeout transitions starting at state m, denoted by max m . We let max + m be a special number representing any real number greater than max m , the maximal nite value in transitions starting from m; it has the following arithmetic properties:

max + m ∈]max m , ∞[, max + m + k = max + m for
every nite real number k ≥ 0 and max + m + ∞ = ∞. The abstract time domain for the clock when the TFSM-TG is in state m is Atd m = [0, max m ]∪{max + m } where [0, max m ] is an integer interval. For machine with states in M , the abstract time domain is Atd S = m∈S Atd m . An abstract execution of a TFSM-TG is dened with abstract execution steps which only consider abstract time domains for clocks and arithmetic properties introduced for the domain. Every execution corresponds to an abstract execution. The latter can be obtained from the former by replacing every timed state (s, x) with (s, max + m ) wherever x is greater than max m . We can show that there is an abstract execution of a TFSM-TG with timed input/output sequence α/β if and only if there is an execution of the TFSM-TG with α/β. This indicates that abstract executions can be used instead of executions.

4.2 Distinguishing automaton Denition 6. Given a specication machine S = (S, s 0 , I, O, λ S , ∆ S ) and a mutation machine

M = (M, m 0 , I, O, λ M , ∆ M ), a nite automaton D = (C ∪ {∇}, c 0 , I, λ D , ∆ D , ∇), where C ⊆ S × S × Atd S × Atd M , λ D ⊆ C × I × Π × C is the input transition relation, ∆ D ⊆ C ×(N ≥1 ∪{∞})×C
is the timeout transition relation and ∇ is the accepting (sink) state, is the distinguishing automaton with timeouts and timed guards for S and M, if it holds that: c 0 = (s 0 , m 0 , 0, 0) For each (s, m, x s , x m ) ∈ C and i ∈ I

(R 1 ) : ((s, m, x s , x m ), i, π -xs s ∩ π -xm m , (s , m , 0, 0)) ∈ λ D if there exists (s, i, π s , o s , s ) ∈ λ S , (m, i, π m , o m , m ) ∈ λ M s.t. π -xs s ∩ π -xm m = ∅ and o s = o m (R 2 ) : ((s, m, x s , x m ), i, π -xs s ∩ π -xm m , ∇) ∈ λ D if there exists (s, i, π s , o s , s ) ∈ λ S , (m, i, π m , o m , m ) ∈ λ M s.t. π -xs s ∩ π -xm m = ∅ and o s = o m
For each (s, m, x s , x m ) ∈ C and the only timeout transition (s, δ s , s ) ∈ ∆ S dened in the state of the deterministic specication (R 3 ) : ((s, m, x s , x m ), δ m -x m , (s , m , 0, 0))

∈ ∆ D if there exists (m, δ m , m ) ∈ ∆ M s.t. δ s -x s = δ m -x m and δ m -x m > 0 (R 4 ) : ((s, m, x s , x m ), δ m -x m , (s, m , max + s , 0)) ∈ ∆ D if there exists (m, δ m , m ) ∈ ∆ M such that δ m -x m > 0, δ s -x s > δ m -x m ,

and

x s + δ m -x m > max s (R 5 ) : ((s, m, x s , x m ), δ m -x m , (s, m , x s + δ m -x m , 0)) ∈ ∆ D if there exists (m, δ m , m ) ∈ ∆ M such that δ m -x m > 0, δ s -x s > δ m -x m , and
x s + δ m -x m ≤ max s (R 6 ) : ((s, m, x s , x m ), δ s -x s , (s , m, 0, max + m )) ∈ ∆ D if there exists (m, δ m , m ) ∈ ∆ M such that δ s -x s > 0, δ s -x s < δ m -x m ,

and

x m + δ s -x s > max m (R 7 ) : ((s, m, x s , x m ), δ s -x s , (s , m, 0, x m + δ s -x s )) ∈ ∆ D if there exists (m, δ m , m ) ∈ ∆ M such that δ s -x s > 0, δ s -x s < δ m -x m ,

and

x m + δ s -x s ≤ max m (∇, i, [0, ∞[, ∇) ∈ λ D for all i ∈ I and (∇, ∞, ∇) ∈ ∆ D
In a state (s, m, x s , s m ) of D, (s, x s ) and (m, s m ) represents timed state of the specication and mutation machines. The two machines synchronize on input actions (see R 1 , R 2 ) and time delays. They do not synchronize on timeout, i.e., the expiration of a timeout in a machine is not necessarily synchronized with the expiration of a timeout in the other machine (case of R 4 to R 7 ); A synchronization of timeout transitions may occur depending on the respective clocks (see R 4 ). The automaton also uses a sink state ∇ to identify synchronizations on an input but with dierent outputs. An execution of D from a timed state (c, x) is a sequence of steps between timed states of D; it can be dened similarly to that for a TFSM-TG. An execution starting from (c 0 , 0) and ending at ∇ is called accepted. D is complete because S and M are complete. Every execution of D corresponds to an execution of the specication S and an execution of the mutation machine M. Indeed, transitions in D are directly dened by execution steps of S and M and indirectly by transitions of S and M permitting the steps. We let e = (e 1 , e 2 ) represent an execution of D, where e 1 and e 2 are the corresponding executions of S and M. Clearly inp(e) = inp(e 1 ) = inp(e 2 ). For every execution e 1 of S there is an execution e 2 of M such that e = (e 1 , e 2 ) is an execution of D. For every execution e 2 of M there is an execution e 1 of S such that e = (e 1 , e 2 ) is an execution of D. Moreover e 1 and e 1 have the same timed input/output sequence whenever (e 1 , e 2 ) and (e 1 , e 2 ) represent executions of D.

It means that D represents the comparisons of timed input/output sequences of M with an expected timed input/output sequence of S.

Lemma 3. An execution e 2 of M is unexpected if and only if there is an execution e 1 of S such that e = (e 1 , e 2 ) is an accepted execution of D.

Only unexpected executions of M having deterministic combs can be executions of mutants; the combs for such executions are revealing; they could be identied in checking whether Mut(M) contains nonconforming mutants. Lemma 4. Rev M (α) = ∅ if and only if there exists an accepted execution e = (e 1 , e 2 ) of D such that π e2 is a deterministic comb and α = inp(e 2 ).

The deterministic revealing combs for a test α can be computed from the distinguishing automaton for the specication and mutation machines. Their computation works as follows. First we compute the accepted executions of the distinguishing automaton D with α; this can be done by dening a product an automaton for α and D, which we do not formalize for the sake of simplicity.

Each accepted execution corresponds to an execution of the specication and an unexpected execution of the mutation machine; the deterministic combs for the unexpected execution of the mutation machine belong to Rev M (α).

The following corollary is a consequence of Lemma 4 and Lemma 1.

Corollary 2. Mut(M) contains nonconforming mutants detectable with α if and only if there exists an accepted execution e = (e 1 , e 2 ) of D such that π e2 is a deterministic comb and α = inp(e 2 ).

5 Boolean Formulas Encoding (Un)Detected Mutants

We encode the test-surviving mutants with Boolean formulas over Boolean variables; each variable corresponds to a transition in M. Henceforth we let variable t represent both a transition in M and the corresponding variable. A solution of a Boolean formula assigns True or False to the variables, which we can use to choose a (mutant) sub-machine in M. Intuitively, we can replace in a mutant a suspicious transition by another (compatible) suspicious transition to obtain a dierent mutant. However, every mutant denes all the trusted transitions. We say that a sub-machine (possibly a mutant) of M is determined by a Boolean formula ϕ dened over the transition variables of M if: (1) the sub-machine includes all the trusted transitions in M and, (2) there exists a solution of ϕ which assigns True to a transition variable if and only if the corresponding transition is in the sub-machine. The encoding Boolean formula is the conjunction of two sub-formulas. The rst sub-formula encodes the complementary of the detected mutants, which can contain mutants and non-deterministic sub-machines of the mutation machine as well. The second sub-formula encodes all the mutants only.

The encoding of the detected mutants uses the suspicious transitions in revealing combs corresponding to the tests in a test suites. This encoding is inspired by Lemma 2. Given a set of revealing combs Rev M (α) for test α, we let ϕ α = π∈Rev M (α) t∈Suspπ t be a Boolean formula over the variables for the suspicious transitions in revealing combs of Rev M (α). ϕ α encodes all the mutants and also non-deterministic sub-machines involved in revealing combs; this means that every solution of ϕ α determines a sub-machine of M that denes all the transitions in a comb π ∈ Rev M (α). So, ϕ α determines all the sub-machines with unexpected outputs for input α. Let ¬ϕ denote the negation of ϕ. Every solution of ¬ϕ α sets variables for some suspicious transitions in every revealing comb to False; this indicates that every mutant determined by ¬ϕ α does not dene some suspicious transitions from each comb in Rev M (α) . Lemma 5. For every P determined by ¬ϕ α and every π ∈ Rev M (α), there exists t ∈ Susp π which does not belong to λ P ∪ ∆ P .

The following lemma is a consequence of Lemma 5 and Corollary 1. Lemma 6. Every mutant determined by ¬ϕ α survives α.

Let T S = {α 1 , α 2 , . . . , α n } be a test suite. We dene ϕ T S = αi∈T S ϕ αi . The formula ϕ T S determines all the sub-machines which produce unexpected outputs for an input α i ∈ T S, i.e., the sub-machines detected by at least one test in T S. A determined sub-machine is not necessarily a mutant; so we need to encode all the mutants in M. We can proof the following lemma by using Lemma 6. Lemma 7. Every mutant determined by ¬ϕ T S survives the test suite T S.

Note that ¬ϕ T S determines not only mutants but also non-deterministic sub-machines of the mutation machine M = (M, m 0 , I, O, λ M , ∆ M ). In other to exclude these non-deterministic sub-machines of M, we encode the set of mutants with a Boolean formula ϕ M . Each of its solutions determines exactly one timeout transition in every state (which is expressed by (Eq 2)) and a subset z mi ∈ Z mi , where Z mi is the maximal set of subsets of λ M (m, i) satisfying the cluster conditions introduced in Section 2.2. Let us dene:

ϕM = m∈M ϕ∆m ∧ i∈I ϕmi ∧ t∈λ S ∪∆ S ¬t ( Eq 1)
where We can prove the following theorem by using Lemma 8 and Lemma 7.

ϕ∆m = t∈∆ M (m)   t ∧ t ∈∆ M (m),t =t ¬t   ( Eq
Theorem 1. Formula ϕ M ∧ ¬ϕ T S determines exactly the mutants undetected by the test suite T S.

For M 1 in Figure 1b, ϕ ∆s2 = (t 8 ∧ ¬t 11 ) ∨ (t 11 ∧ ¬t 8 ), ϕ s2a = t 7 ∧ t 9 ∧ t 10 and ϕ s2b = (t 4 ∧¬t 5 ∧¬t 6 )∨(t 5 ∧t 6 ∧¬t 4 ). ϕ ∆s2 ∧ϕ s2a ∧ϕ s2b determines the transitions starting at s 2 in a mutant. We can make similar formulas for s 1 and s 3 , according to the sets Z s1a , Z s1b , Z s2a , Z s2b , Z s3a and Z s3b presented before Denition 3.

Finally, ϕ M1 = ϕ ∆s1 ∧ ϕ s1a ∧ ϕ s1b ∧ ϕ ∆s2 ∧ ϕ s2a ∧ ϕ s2b ∧ ϕ ∆s3 ∧ ϕ s3a ∧ ϕ s3b .
The test α 1 = b0.5a7.5 triggers executions in the distinguishing automaton for M 1 . These executions correspond to the expected execution e 1 , e 2 and : encodes the mutants surviving α 1 , e,g., the mutant P 2 in Figure 2b. Consider a test suite T S and a fault model M, , M ut(M) . The test analysis problem consists in verifying whether the test suite is complete for the fault model. The test generation problem aims at generating a complete test suite, by adding new tests to the given one. Solving the test analysis problem, we encode the mutants undetected by the tests with the Boolean formula ϕ M ∧ ¬ϕ T S ; then we use a solver to determine a mutant undetected by T S, according to Theorem 1 T S is complete if no mutant can be determined or only conforming mutants can be determined. We exclude encountered conforming mutants from the fault domain; this is done by making the conjunction of ϕ M ∧ ¬ϕ T S with the formula t∈λ P ∪∆ P ¬t encoding the mutants dierent from P.

e 3 = (s 1 , 0)0.5 [t2] (s 1 , 0.5)b/x [t3] (s 2 , 0)7 [t11] (s 3 , 0)0.5 [t15] (s 3 , 0.5)a/y [t12] (s 1 , 0) and e 4 = (s 1 , 0)0.5 [t2] (s 1 , 0.5)b/x [t3] (s 2 , 0)7 [t11] (s 3 , 0)0.5 [t15] (s 3 , 0.5)a/y [t13] (s 1 , 0).
Theorem 2. T S is complete for M, , Mut(M) if and only if ϕ M ∧ ¬ϕ T S is not satisable or all the mutants it determines are conforming.

In order to generate a complete test suite, our procedure analyzes the current test suite T S. If it is complete the procedure stops and return T S; otherwise we generate a test α detecting nonconforming mutants from a mutant surviving T S. α can be the timed input sequence of an execution of the distinguishing automaton for the specication and a surviving mutant (which is a part of the mutation machine), according to Corollary 4. Then we compute ϕ α and we analyze the new test suite T S = T S ∪ {α} by solving the formula ϕ M ∧ ¬ϕ T S ∧ ¬ϕ α ; this may trigger the generation of a new test. We generate iteratively new tests until all the procedure stops and a complete test suite is returned.

For example, the test suite T S = {b0.5a7.5} is not complete for M 2b, where ϕ T S = ϕ α1 and α 1 = b0.5a7.5. An implementation of our approach generates the complete test suite {α 1 = b0.5a7.5, α 2 = b0.5a0.5a3.5a0.5, α 3 = b0.5b0.5} detecting all the seven nonconforming mutants. α 2 detects P 2 ; it is the timed input sequence of an accepted execution of the distinguishing automaton of S 1 and P 2 . A nonconforming mutant undetected by α 1 and α 2 was used to generate α 3 ; it denes transition t 5 .

Let us compare our complete test suite generation approach with the one proposed in [START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF]. We focus on detecting specic faults represented in the mutation machine which can have more states than the specication, so do the mutants.

The specication machine is not necessarily "minimal". We believe that for the traditional fault model, the method in [START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF] could be faster than ours (especially if the time required to minimize the specication is not considered), but the method in [START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF] could generate redundant tests. The reason is that the method in [START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF] is an application of the W-method over a classical TFSM which represents an abstraction of TFSM-TG and the W-method was developed for minimal classical FSM. The size of the abstract TFSM-TG [START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF] could contribute to decrease the eciency of the test generation [START_REF] Dorofeeva | Fsm-based conformance testing methods: A survey annotated with experimental evaluation[END_REF][START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF]. The W-method generates tests by combining input sequences (namely, state cover sets, all input sequences of a certain length and state identication input sequences or characterization sets); (15, 4.45) it does not check whether each generated test actually detects faults and several generated tests can detects only one mutant. Thus the number of generated tests can be bigger and less adequate to focus on faults specied in customized fault models represented with mutation machines. Where [START_REF] Tvardovskii | Deriving tests with guaranteed fault coverage for nite state machines with timeouts[END_REF] will generate a huge test suite to detect all possible faults, our method will generate test suites of reduced sizes for customized fault models.

We performed an empirical evaluation of our approach with a proof-ofconcept tool developed in C++ and randomly generated TFSM-TG with two inputs and outputs, with a maximal timeout of 5 (resp. 10) for the specication (resp. the mutations machines). The tool uses cryptoSAT [START_REF] Soos | Extending SAT Solvers to Cryptographic Problems[END_REF]. Preliminary results appear in Table 1 which presents sizes and generating times for test suites from mutation machines dening multiple of mutants with few states. For each number of states and number of mutants, we generated tests for several mutation machines. The tool can take a long time in encoding the fault domain; this may happen when there are too many compatible (not necessarily equal) timed guards and timeouts dened in the same state. Since this operation is done for every state, it could be distributed. In most of the situations, the tool can rapidly encode the fault domain and the test generation becomes faster.

Conclusion

We have proposed a multiple mutation testing theory to testing real-time systems represented with nite state machine extended with timed guards and timeouts (TFSM-TG). We developed an approach to generate complete test suites for fault models represented with mutation machines. The approach relies on the denition of distinguishing automaton for mutation machine and the construction of Boolean formulas encoding the (faulty) implementations undetected by tests. We implemented the approach in a proof-of-concept tool which we used to evaluate the eciency of the approach. The experimental results show that the approach can be used to derive tests for non-trivial TFSM-TG fault models representing an important number of faults.

Ongoing work includes developing open access benchmarks and use them to compare the existing test generation methods for TFSM-TG. We also plan to generate symbolic tests, i.e., timed input sequences with delay intervals instead of simple delays, and complete test suites consisting of a single test.

Fig. 2 :

 2 Fig.2: Two mutants included in the mutation machine M 1 in Figure1b

  Each solution of ϕ ∆m sets to True the variable of exactly one timeout transition dened in m. All the variables in exactly one set Z mi has all its variables assigned to True by each solution of ϕ mi . The variable for each trusted transition is set to T rue in every solution of ϕ mδ ; this is because they are not compatible with any other transitions and thus dened by every mutant. We can conclude that each solution of ϕ M determines a mutant in Mut(M) and every mutant is determined by a solution of ϕ M , i.e., ϕ M determines all the mutants in Mut(M). According to its denition, every solution of ϕ M never assigns True to the variables for compatible transitions and determines a mutant. Lemma 8. ϕ M determines the mutants in Mut(M).

  1 . Transitions represented with dashed lines are mutated. Transition identiers appear in brackets. M 1 is nondeterministic because, e.g., t 4 and t 6 are compatible or the two timeout transitions t 8 and t 11 start from s 2 . The suspicious transitions in Susp M1 are t 4 , t 5 , t 6 , t 8 , t 11 , t 12 , t 13 , t 16 , t 17 ; the other transitions are trusted. The specication S 1 in

	Figure

  The executions of the mutants correspond to deterministic combs. Only executions of a non-deterministic sub-machine of M can correspond to a nondeterministic comb. Mutants included in M can have common executions. Because we want to avoid the enumeration of the mutants one-by-one, we would like to identify the mutants which have an execution in common, i.e., they are involved in the execution. A mutant P is involved in an execution e of M if P denes all the suspicious transitions in the comb π e of e, i.e., Susp πe ⊆ λ P ∪ ∆ P . This is because every mutant denes all the trusted transitions in M and only suspicious transitions vary between the mutants. We let Rev M (α) denote the set of revealing combs for executions with timed input sequence α. It holds that Rev M (α) = P∈Mut(M) Rev P (α). The following lemma indicates that detected mutants are involved in revealing combs of the mutation machine. Lemma 2. Test α detects mutant P if and only if there exists e ∈ Exec M (α) such that π e ∈ Rev M (α) and Susp πe ⊆ λ P ∪ ∆ P .Corollary 1. Mutant P survives α if and only if for every π ∈ Rev M (α), some transitions in Susp π do not belong to λ P ∪ ∆ P .

	Lemma 1. Mut(M) contains nonconforming mutants if and only if Rev M (α) =
	∅, for some timed input sequence α.

  1 are time-elapsing. t 8 permits the fourth step. A distinct execution, let us call it e 2 , could involve t 11 at the fourth step. Another distinct execution is possible by making a timeout step at the third step, i.e., by ring t 11 . The comb π e1 = t 2 t 3 t 8 t 8 t 7 is for e 1 does not form a path in M 1 ; the symbol " " is inserted between timeout transitions permitting time-elapsing steps which just precede input/output steps. π e1 is deterministic and a mutant can dene all the suspicious transitions in it. π e1 is not revealing because the timed output sequence of e 1 , x0.5x7.5 is expected. The mutants dening only the suspicious transitions occurring in π e1 and all the trusted transitions cannot be detected by b0.5a7.5. π e2 = t 2 t 3 t 8 t 11 t 7 is not deterministic because t 8 and t 11 are compatible and suspicious; a mutant cannot dene these two transitions. π e2 is not revealing since it is not deterministic.

  The comb for e 3 and e 4 are π e3 = t 2 t 3 t 11 t 15 t 12 and π e4 = t 2 t 3 t 11 t 15 t 13 ; they are revealing because the produced timed output sequence x0.5y7.5 is unexpected; so they characterize the mutants detected by test α 1 . We recall that π e1 and π e2 are not revealing. The formula ¬ϕ α1 = ¬((t 11 ∧ t 12 ) ∨ (t 11 ∧ t 13 ))

  1 in Figure 1b because ϕ M1 ∧ ¬ϕ T S determines the nonconforming mutant P 2 in Figure

Table 1 :

 1 Size of the generated complete test suites and generating time ; for an entry (x, y), x is the size of the test suite and y is the generating time in seconds

		#mutants in the fault domain
	#states	10 4	10 8	10 12
	8 states (6, 0.22) (7, 0.32) (12, 45.57)
	10 states (3, 0.16) (8, 2.34) (21, 46.94)
	12 states (5, 0.72) (5, 2.77)