Alessandro Della Corte
email: alessandro.dellacorte@uniroma1.it

Antonio Battista
email: antonio.battista@univ-lr.fr

Ivan Giorgio
email: ivan.giorgio@uniroma1.it

Francesco Dell'isola
email: francesco.dellisola@uniroma1.it

Ivan Giorgio Memocs

Francesco Dell

Modeling deformable bodies using discrete systems with centroid-based propagating interaction: fracture and crack evolution

In this paper we use a simple discrete system in order to model deformation and fracture within the same theoretical and numerical framework. The model displays a rich behavior, accounting for different fracture phenomena, and in particular for crack formation and growth. A comparison with standard Finite Element simulations and with the basic Griffith theory of fracture is provided. Moreover, an 'almost steady' state, i.e. a long apparent equilibrium followed by an abrupt crack growth, is obtained by suitably parameterizing the system. The model can be easily generalized to higher order interactions corresponding, in the homogenized limit, to higher gradient continuum theories.

Introduction

Motivation and basic ideas

Modeling fracture in an effective way has always been a major challenge for solid mechanics. Many sophisticated theoretical and numerical tools have been devel-oped, and considerable progresses have been obtained in recent years in the framework of both classical and generalized continuum models (see [1,2,3,4,5]). However, the richness of the phenomenology is leading to theoretical formulations that are becoming increasingly complex and may be computationally expensive and/or involve a stochastic approach in order to capture the peculiar characteristics displayed by real fractures. The idea behind the present work is to explore a very simple discrete system characterized by a centroid-based propagating interaction, and evolving through actual configurations connected by virtual (i.e., not visible) ones; the system, in its basic version and without including fracture, was introduced in [6], while some preliminary results on fracture have been shown in [7]. The model presents some similarities with Molecular Dynamics (MD) but, as it will be shown, its specific features cannot be recovered within standard MD models. The potential advantages of the proposed approach are mainly linked to the simplicity of the model. First of all this results in a significantly low computational cost; moreover, deformation and fracture are covered here within the same simple model; finally as we will see, more general continuum theories can be easily numerically investigated by means of slight modifications of the algorithm. Since generalized continua are one of the most promising and rapidly evolving areas in modern mechanics, this last feature of the proposed model may be particularly appealing (for a theoretical coverage on higher gradient theories, see e.g.: [8,9,10,11], and specifically for an approach combining higher gradient theories with lattice models, see [12]). A closely related topic is that of micromorphic/microstructured continua (see [START_REF] Eringen | Microcontinuum field theories: I. Foundations and Solids[END_REF][START_REF] Germain | [END_REF]15,16] for classical references and [17,18,19,20,21,22] for interesting applications), which can be viewed as a generalization of higher gradient theories and may benefit as well from the development of new discrete approaches. What is making these subjects critical, in the opinion of the authors, is the advancements in manufacturing possibilities in the last years, as 3D printing and other computeraided manufacturing techniques, which are resulting in new metamaterials requiring a suitable theoretical description (and related numerical techniques) for objects whose richness at the micro-scale cannot be captured by classical continuum models (see e.g.: [23] for a review of recent results). Discrete systems are, of course, very frequently studied in order to address the aforementioned problems. In particular, Molecular Dynamics (MD) is by now a very large research area with specific methods and very sound results (see e.g.: [24,25,26,27]). The model presented herein, while sharing certain basic features with MD ones, is characterized by some relevant differences in the approach. MD is indeed based on the numerical study of systems constituted by a very large number N of elements. The numerical computation of the trajectories of particles in the ordinary 6N-dimensional phase space of positions and momenta employs the classical mechanics laws of motion (for a sample of the numerical method the reader can see: [28,29,30]). The system investigated in the present paper, on the other hand, does not consider any explicit equation of motion. Instead, in order to simplify the model to the maximum, the elements of our discrete system move according to an interaction law which:

• is purely geometrically based;

• propagates along square frames centered in special elements (the 'leaders' of the system); • drives the elements through virtual configurations which are needed for the computation of the evolution, but are not visible in the final output.

An important class of discrete systems for continuous deformable bodies, and in fact a particular case of MD models, is given by masses-springs or beads-springs models, which can be characterized in different ways by changing the geometry of the springs and their type (extensional, torsional). The mass-spring systems are common tools in Computer Graphics for the simulation of soft bodies [31,32].

There are several relevant differences between these systems and the proposed model. The two probably most important ones are the following:

i. In case of mass-spring system, a uniform motion of the leaders results in a global motion r(t) = r * (t) + r * * (t), with r * * (t) being a periodic function and r * (t) a transient term such that ||r * (t)|| → 0 as t → ∞ [START_REF] Bishop | The mechanics of vibration[END_REF]. As we will see, no such decomposition makes sense in case of the proposed system. For instance, it can be proven that a uniform motion imposed to a leader results, asymptotically, in a rigid motion of the system. ii. It can also be proven that uniformly accelerated motion of a set of leaders results in our case in a disaggregation of the system (i.e. there is no R ∈ R such that for every

t ∈ [0, ∞[ the system is contained in B R (x, y) for some x, y ∈ R 2 ).
The characteristics of our model have, as we will see, several advantages, but at the same time make harder a standard variational formulation. As we will show, energy-based investigations are possible in the discrete context here considered, but the identification of an explicit Lagrangian whose minimization leads to the exact dynamics displayed by the discrete system is far from trivial, and will be one of the main objectives of future investigations. One of the main checks we will perform will concern the systematyc comparison with Finite Element simulations. Nowadays FE analysis is indeed probably the most reliable tool for the numerical simulation of the behavior of deformable bodies, also thanks to the possibility of adapting it to the features of the problem studied. Isogeometric analysis (see for instance [START_REF] Hughes | [END_REF]35,36,37]), in particular, can be especially convenient for shape optimization problems that easily arise in the study of multi-agent systems moving in unbounded domains and starting form arbitrary configurations, as is the case in our context.

A summary of the algorithm and of the formalism employed

In this section, we will briefly describe the model studied in the present paper. The algorithm and the formalism will be summarized in their most relevant features; the reader can find a fully detailed description in [6].

Let us consider a discrete system S constituted by a finite number of points which are characterized by their position in the real plane; we will call them the 'elements' of S. The elements are set, in the initial configuration C 0 , in the nodes of a squared grid sized L × L. We consider a set of discrete time steps T m = {0,t 1 , ...,t m , ...} and an orthonormal reference system with axes parallel to the sides of the system, the unit length being equal to the cell side in C 0 . Each Lagrangian element (i, j) of S, placed in (x 1 i j , x 2 i j )(t m ), has a set of neighbors

N n ( ī, j) := (i, j) ∈ C 0 : ρ[( ī, j), (i, j)] = n
where ρ is the R 2 -Chebyshev distance, i.e. the distance in R 2 given by ρ((

x 1 , x 2 ), (y 1 , y 2 )) = max |x 1 -y 1 |, |x 2 -y 2 |
. With this definition, selecting n = 1 and n = 2, we have respectively the first and second neighbors as shown in Fig. 1. Let us select a leader element L whose actual position is defined by a prescribed motion M : t m ∈ T m -→ (x 1 ī j, x 2 ī j)(t m ) ∈ R 2 . We are now ready to describe the interactions between the elements. We will consider virtual steps in between two actual states, i.e. states that are invisible in the real displacements of the system but are necessary for computing its evolution. Let us consider a configuration V t 1 0 such that the leader L is positioned in M (t 1 ): this is defined as the first virtual configuration. Then, recalling that by centroid of a set of points P 1 (x 1 1 , ...,

x 1 n ), ..., P m (x m 1 , ..., x m n ) ∈ R n one means the point P(∑ m i=1 x i 1 m , ∑ n i=1 x i 2 m , ... ∑ n i=1 x i n m ), the next virtual configuration V t 1
1 will be defined as the on in which: 1. the leader L remains in V t 1 0 ; 2. every one of its first neighbors N 1 ( ī, j) moves to the centroid of its own first neighbors in V t 1 0 ; 3. all others elements remain in the position they had in V t 1 0 . Iterating and generalizing to the n-th virtual step, V t 1 n will be the virtual configuration in which:

1. the leader L together with its first n -1-th neighbors are in the same position they had in V t 1 (n-1) ; 2. every one of the n-th neighbors of L has moved to the centroid of its own first neighbors in V t 1 (n-1) ; 3. all others elements are in the same position they had in V t 1

(n-1) . We will get the actual configuration C t 1 when n equals the maximum Lagrangian Chebyshev distance of the elements from the leader (see Fig. 2). It is easy to see that in the model so described edge effects will arise, since spontaneous shrinking will concern boundary elements because of the non-symmetric placement of their neighbors. In order to overcome this problem several standard possibilities can be considered. Choosing probably the most simple one, we introduce a 'fictitious' boundary constituted by an external frame of elements which simply follow, at a fixed distance and always in the same direction, their closest 'true' element. The fictitious elements move only in a specific virtual time step which follows the other ones. In this way, every true element will have a complete set of 8 neighbors (see Fig. 3 for a graphical representation. An explicit description of the motion of the fictitious elements is in [6]). The distance at which the fictitious elements are positioned represents a key characteristic length, determining the lattice step of the system in the equilibrium configuration. In the following numerical simulations, unless differently specified, this length is equal to 1 length unit for all the fictitious elements. The previous algorithm can be easily generalized to the case of second and, in general, n-th neighbors interaction, by simply computing the centroid of sets of points having maximum Lagrangian Chebyshev distance of 2 (in general, n). In this case, additional sets of fictitious elements have also to be introduced in order to have, for every element of the system, a full set of neighbors. In [6] the system is also gen-Fig. 3 Graphical representation of fictitious elements. Every fictitious element (empty red dots) moves in rigid translation with the closest true one (blue dots); every vertex true element carries in rigid translation the vertex fictitious one and also its two fictitious neighbors. eralized to the case of multiple leaders, which will be used throughout the present paper. The geometric centroid of a given set of points P 1 , ..., P n minimizes the sum of the squared distances, i.e. the function:

f (P) = n ∑ i=1 ||P i -P|| 2 (1)
It is therefore possible to see that for a first-(second-)neighbor interaction, a natural deformation energy density can be written as a function of first-(second-)order finite differences of the placement χ for a given point, and therefore, as the step length goes to zero, a second gradient homogenized energy E [∇χ(P), ∇∇χ(P)] can be conjectured for second neighbors interaction systems (for more details see again [6]).

The algorithm described before has much in common with theoretical models conceived for describing collective behaviors (e.g. the evolution of swarms, flocks etc.). Such models have already been used for mechanical investigations, in cases in which the potential energy of a mechanical system is minimal in correspondence with stable equilibrium positions of a swarm system [38,[START_REF] Temür | Recent Advances in Swarm Intelligence and Evolutionary Computation[END_REF]. However, in the present work we prefer not to use the name 'swarm' system, since one usually associates that name, and in particular the expression Particle Swarm Optimization (PSO), to cases in which each particle of the swarm represents a potential solution of a given (a priori formulated) optimization problem [START_REF] Kaveh | Harmony search algorithms for structural design optimization[END_REF][START_REF] Clerc | Particle swarm optimization[END_REF][START_REF] Vaz | [END_REF]. Our aim, in fact, is not the study of numerical tools for the solution of well defined optimization problems, but rather the development of a new model directly accounting for the description of the phenomenology; other features distinguishing the proposed discrete model from PSO ones is the presence of leaders (whereas PSO models are usually anarchic) and the absence of any randomness.

It should also be pointed out that some basic features of our model (in particular the Lagrangian character of the neighborhood) are quite different from standard swarm robotic modeling. In modeling swarms, flocks and schools dynamics, a topological concept of neighborhood (rather than a metric one) is emerging as one of the most promising in order to account for the observed phenomenology (see [43,44]).

The model proposed in the present paper mixes these two ideas, since it is based on a concept of neighborhood which depends on the topological distance between the elements, while the way in which the interaction works depends on the metric (Euclidean) distance in the actual configuration. Finally, another feature of the proposed model which is often met in swarm modeling is that the elements, due to the presence of the fictitious boundary, do not behave all in the same way; this has been proposed as one of the possible discriminating factors between crowds and swarms/flocks/schools models (see e.g. [45,46]), but also in this case it can be pointed out that the property of being (or not) a boundary element, in the model considered herein, is Lagrangian, i.e. it is preserved during the time evolution of the system.

Short summary of previous numerical results

We summarize in this section the numerical results obtained in [6]. With the aim of comparing the discrete models characterized by (respectively) first and second neighbor interaction (FNI and SNI) with first and second gradient continuum theories, some simple cases were investigated. Two squared systems FNI and SNI were considered, in which the leaders L, situated in a vertex of the square, were pulled outside (or pushed inside) in the diagonal direction at 45 degrees with respect to the sides of the system in C * with a uniform motion. The numerical results obtained were compared with Finite Element simulations of a 2D continuous squared body, at the vertex of which a prescribed displacement was imposed. Two cases were considered: a standard energy (for the first gradient case) and the Mindlin general form [15] for the second gradient case.

In Fig. 4 we can appreciate one of the similarities between the continuum models and the discrete system here described (for other features see [6]). As we can see, in the FNI/first gradient continuum cases a loss of boundary convexity around the vertex is observed, whereas this behavior is not present in SNI/second gradient continuum cases; this is highlighted in Fig. 5. Releasing the vertex at a certain time step, i.e. letting the system evolve while the leader is stopped, one observes the return to the original configuration, in agreement with the behavior expected in elastically deformed bodies. As can be seen in Fig. 6, the return to the initial configuration is not instantaneous, which implies the presence of inertia effects even if no explicit variable accounting for the elements mass was introduced. It is possible to prove that actually we have an asymptotic convergence, which means that some viscous effect has to be considered if one want to identify the Lagrangian system exactly corresponding to the evolution of the proposed system ( [6]). Fig. 6 Imposing a stop to the leader (at t = 150) the system tends to return to the original configuration.

These results (among others) suggest that the evolution of the presented discrete system resembles that of elastic deformable bodies, and that varying the order of the interaction, specific characteristics of higher gradient theories are also recovered. In the following we will see a more direct comparison with 2D continuum simulations performed with COMSOL. Finally, we want to devote a few words to the nonlinear character of the evolution of the system. In the following sections, a fracture algorithm (obviously entailing nonlinearities) will be introduced, but even the basic form of the model as described up to this point exhibits a nonlinear behavior. In Fig. 7, we compare two identical systems in which different actions are imposed to one single leader. In the left panel, a speed of 1 length units per unit time is imposed for a total of 10 time steps; in the right panel, we imposed a speed of 0.001 length units per unit time for a total of 10000 time steps. As one can see, the resulting configurations are clearly different, thus implying the nonlinearity with respect to the imposed external action.

Further comparison with Finite Element simulations

In this section we want to directly compare the results of our discrete model with the ones obtained by standard Finite Element simulations.

In Fig. 8 we show the superposition of:

• a simulation with an imposed external action consisting in pulling in opposite directions two opposed vertexes of our system; • a standard FE simulation (performed with COMSOL) in which a similar action (imposed displacement) is applied on a squared 2D continuum.

For the continuous simulation, we considered both a classical Cauchy continuum (left) and a second gradient one (right). In color map the modulus of the displacement is shown. The corresponding discrete simulation involve first-(left) and second-neighbors (right) interaction. As it can be seen, in both cases the shape of the sample is accurately approximated by the discrete system, in the specific case employing only a very limited number of elements (12 × 12 = 144 elements in total). This means that our system produces reliable results with a very limited computational cost, at least with this kind of external action.

Spontaneous evolution and Discrete Energy Equivalent

The previously described algorithm can be adapted to the case of a spontaneously evolving system, i.e. systems in which we do not impose a particular motion to any of the elements. To do so, we start from an arbitrary initial configuration C 0 , in which in general the elements do not lie in the centroid of their neighbors. We apply the algorithm starting from a selected element P (which we will call the 'pseudoleader') and then proceeding through concentric square frames. The only difference with the case of a 'true' leader is that the pseudo-leader P does not obey to an imposed motion, but simply goes, in the first virtual step, in the centroid of its first neighbors as all other elements do.

The observations made at the end of section 1.2 suggest the introduction of a quantity E D (t), which we will call Discrete Energy Equivalent (DEE), defined on a geometrical basis. The DEE will represent a measure of the deformation energy stored in the actual configuration of the system in case of first neighbors interaction. Due to the minimum property of the centroid above recalled (see Eq.( 1)), and considering that an (Euclidean) distance equal to 1 or to √ 2 corresponds to the pairs of first neighbors elements in the equilibrium configuration C * , it is natural to define E D as follows:

E D = ∑ S 1 (d t -1) 2 + ∑ S 2 (d t - √ 2) 2 2 (2)
Here the sums are extended over all the pairs of first neighbors 'true' elements. The sets S 1 and S 2 contain the pairs that in C * lie respectively orthogonally or diagonally with respect to the sides. A unit length is subtracted from d t ( √ 2 in case of diagonal pairs) because the step length is unitary; more precisely, because the fictitious elements are at a fixed unitary distance from the true boundary. Since the sums have the dimension of an area, a coefficient α having the dimension of a force unit over a length unit has to be introduced to define the energy of the system: E = αE D . The first set of numerical simulations will concern the study of the DEE decay in systems evolving without leaders, and since in our simulations E D (C 0 ) = 0, C 0 can be seen as a pre-stressed configuration. A basic case, in which the center element is selected as the pseudo-leader, is represented in Fig. 9. A sharper case in which the initial configuration is more stressed is shown in Fig. 10. In both cases, the energy decay is well interpolated by an exponential behavior; in Fig. 11 the interpolation of the form e ax b , with a ≈ -0.59 and b ≈ 0.52 is shown. This behavior is well known for several elastic systems in which viscous dissipation occurs in both linear and nonlinear cases (see e.g. [47,48,49,50]). One may wonder whether the system always asymptotically converges to the same limit configuration (i.e., the one with E D = 0) independently on the choice of the pseudo-leader. This seems a rather natural request if the system is intended to model the behavior of a deformable body. In the Appendix (see 4) we prove that in the 1D case the answer is yes. In 2D, the conjecture that the asymptotic configuration is independent on the pseudo-leader choice remains more than reasonable, but the proof is more difficult, since in this case it is not true anymore that the total discrete energy decreases in every virtual step. Numerical evidence of this somewhat surprising statement is shown in Fig. 12, where one can observe that, after a first approximately exponential phase, the energy increases slightly before its eventual decay. We want now to perform some quantitative analysis in order to evaluate some magnitudes concerning the evolution of our system on the basis of the DEE above defined. Let us first consider a very simple case in which uniaxial motion is imposed to a set of leaders constituted by the elements belonging to two opposite sides of the system. In particular, the fictitious elements relative to the left side are treated as leaders, as they are motionless by definition (from a mechanical point of view, the side has an imposed displacement equal to zero). Moreover, we impose a motion to the elements belonging to the right side; the motion is uniform and directed along the x axis up to a certain time t s , after which they stop. The system is then left evolving until equilibrium is achieved. The result, with relative energy versus time plotting, is shown in Fig. 13 (here and in all the following simulations, the leaders are represented by red dots). Since our external action is an imposed motion, it is not immediate to derive from it a discrete version of the applied uniaxial force σ D . We can however reasonably define it in various equivalent ways; for instance, we can use the identity

(σ D ) 2 A 2E = E = αE D
where A is the area of the system at the equilibrium, E is Young's modulus and α the previously introduced constant. We will not lose generality by selecting the unit for σ D so as to have 2E = 1, and therefore we get (σ D ) 2 = αE D A . The numerical data for A and E D give:

(σ D ) 2 = αE D A ≈ 12.5α 382 ≈ 0.03α (3)
Using the reasonable assumption that a certain imposed motion will approximately correspond to the same σ D if only small and local changes are considered in C * , we will employ the value now obtained in the simulations on the crack formation and growth considered below.

3 Fracture and crack formation and evolution

Introduction of the fracture

In the proposed model, the fracture is intended as a loss of interaction between neighboring elements. Specifically, when the Euclidean distance (evaluated in actual configurations) between two interacting elements overcomes a certain threshold s f , the two elements do not interact anymore. Obviously, when this happens, the computation of the centroid relative to the considered elements is ill defined, in the same sense as intended for what concerns boundary elements, as seen in section 1.2. We solve the problem in a similar way as done before. Indeed, we introduce for every bond break a new pair of fictitious elements evolving in the same way as the fictitious boundary elements introduced above, i.e. each of them following one of the two elements whose bond has been broken. In this way, the two elements among which the fracture occurred become boundary elements in the same sense as the preexisting ones. The fracture algorithm works as follows: we introduce a prescribed order in the n-th neighbors of a given non-fictitious element starting, for example, from the neighbors in the left upper corner of the frame square centered in the considered element and proceeding clockwise. The labeling of neighbor elements allows the introduction of an adjacency matrix J, that in this case is a tensor of boolean variables J i jk , the first two indexes determining the particle (i, j) of the system, and the third one k, ranging from 1 to 8 in FNI, identifying one of its neighbors. When the distance between two elements is larger than s f , the value of the corresponding element of J is set to FALSE= 0 (otherwise it is TRUE= 1). When the centroid is computed, each neighbor coordinate is multiplied by the relative element of J. Moreover, another term is considered in the centroid computation, consisting in the coordinate of fictitious elements, each multiplied by J (where the bar indicates logical complement). It is clear that in this way, when the centroid is evaluated, a true or a fictitious element enters the calculation according to the fact that the interaction is present or broken. The irreversibility of the fracture, which is a desirable feature when modeling solid bodies, is verified since once an element of J is set to 0, it can not become 1 anymore. Another degree of freedom of the model, concerning the fictitious elements appearing after the fracture, naturally emerges form the described algorithm. Indeed, one can choose to place the new fictitious elements at a distance ρ f that does not correspond to the lattice step in C * (more precisely, we should say that it would not correspond to the fixed distance at which the fictitious boundary is from the true elements. It is indeed this last distance that determines the lattice step at which the system is in equilibrium). In this way, as we will see in the following, relevant features of fracture phenomena can be modeled.

Basic fracture and crack evolution

A well known variational approach to the fracture was developed by Griffith in the '20. From Griffith's model (in 2D case), we have at the equilibrium:

E T OT = 2γa + σ 2 2E A - σ 2 2E β a 2 (4)
where γ represents the energy per unit line required to break atomic bonds, E the Young Modulus, A the area of the sample, a is the crack length, σ the stress and β an non-dimensional parameter accounting for the measure of the part of the surface relaxing as a consequence of the crack opening. The most relevant qualitative aspect of Griffith's theory is that the dependence of the energy on the crack length forecasts that below a critical value of a the system is stable, and a crack growth is possible only providing additional energy to the system. If this critical value is reached, the system becomes unstable and the crack evolves spontaneously. In the next groups of simulations we selected two leaders belonging to the bottom side (red elements) and imposed to them a uniform velocity, with components (0.03, 0.03) and (-003, 0.03) respectively. The leaders were stopped at t = 80. In Fig. 14 we set the fracture threshold s f = 1.17 length units; we can see that the crack keeps growing even after the leaders have stopped, that means that the crack length is reaching the equilibrium value for the given imposed action. The crack length equilibrium value is reached around t = 160, after which the cracks remain stable while the system keeps relaxing and releasing deformation energy. In Fig. 15 the same simulation was performed with the fracture threshold set at s f = 1.12. It can be seen that the crack growth does not stop and indeed goes up to the complete split of the system, which is reached around t = 250. This means that stopping the motion at t = 80 is sufficient for furnishing to the system the energy needed to reach the length threshold after which the crack grows spontaneously.

Uniaxial external action

The previous results indicate that qualitative features of Griffith's theory are recovered through our discrete model. However, in order to be able to compare quantitatively our numerical results with Griffith's theory, we need to apply a simpler

Fig. 15 The same simulation of Fig. 15 is performed with s f = 1.12; in this case, the critical length is overcome and the crack growth keeps up to the complete split of the system.

external action, i.e. a uniaxial one, as done in the simulation shown in Fig. 13. Let us first formulate Eq.( 4) by means of the discrete variables considered here. In our model, the energy can be written as:

E T OT = 3α(s e f f ) 2 l + (σ D ) 2 A -(σ D ) 2 β l 2 (5)
where s e f f = s f -1 is the DEE lost when a single bond is broken, A represents the area of the system (in squared unit lengths) and l is the number of broken pairs of elements in the crack.

In the simulation shown in Fig. 16 we imposed to the leaders the same external action used in the simulation of Fig. 13,i.e. we imposed an uniaxial displacement of the leaders in the x direction up to 180 time steps followed by a stop, and a zero displacement to the fictitious elements relative to the left vertical side. In order to see a crack formation and evolution we created a 'defect' by removing 6 bonds close to the middle points of the bottom side. We can observe that a crack indeed opens in correspondence with the defect, and that it reaches a stable length already at t = 1000.

The DEE E D that we obtained as an output (≈ 5) times the dimensional constant α is the discrete estimate of the deformation energy of the system, and therefore represents the quantity given by the last two terms in the right hand side of Eq.( 5). This allows us to provide an estimate for β . Indeed, substituting the value for σ D taken from the simulation shown in Fig. 13, we have:

β = (σ D ) 2 µn 2 -αE D (σ D ) 2 l 2 ≈ 25.7 (6)
According to the standard interpretation in Griffith's fracture theory, this means that an area measuring approximately β × l 2 = 25.7 × 9 is relaxed as a consequence of the crack formation; this area corresponds to a fraction of approximately 0.58 of the total area. Let us now consider the first term in the right hand side of Eq.( 5). The numerical output tells us that a total of 9 bonds, in Fig. 16, are stretched beyond the fracture threshold (see Fig. 17 for a graphical representation). As we recall, we introduced a 'defect' consisting in removing six interactions close to the edge; this means that a total of three elementary bonds were broken by applying the external action. Since the threshold was set at s f = 2.62, the corresponding DEE, i.e. the energy E b required to break the bonds can be estimated by:

E b = α(s f -1) 2 × 6 = 7.68α (7)
The sum E b +αE D ≈ 12.68α is an estimate of the total energy in the system. As one can see, it is very close to the total deformation energy (≈ 12.5α, see Eq.( 3)) with the same external action in case no fracture threshold is considered. This is consistent since, in the proposed model, no additional dissipated energy is associated to the break of the interactions, and therefore the total energy should coincide with the deformation energy. The small difference has to be related with the fact that we (e) Fig. 16 Uniaxial external action imposed to a set of boundary leaders with a defect consisting in the removal of six bounds around the middle elements of the bottom side (s f = 2.62). The leadres stop at t = 180; a crack opens in correspondence with the defect and reaches a stable length (as before, we use empty red dots to indicate that the fictitious elements relative to the left side have an imposed null velocity, i.e. the left side has an imposed displacement equal to zero). introduced a local defect in the simulation shown in Fig. 16; since we imposed an external action consisting in imposed displacements, this means that, assuming the same value for σ D as the one measured in the simulation of Fig. 13, we are slightly overestimating it. Actually, the ratio 12.68-12.5 12.68 ≈ 0.014 is very close to the ratio between the number of removed bonds over the total bonds present in the system ( 6 520 ≈ 0.012).

Fig. 17 Graphical representation of the bonds broken in the simulation shown in Fig. 16.

Almost-steady state

In the next simulation we want to underline a peculiar behavior shown by the model which can be called as an 'almost steady' state followed by a catastrophic evolution of the crack. We considered the same system of the previous one, making only a small change in the fracture threshold s f . Tuning very finely the threshold (i.e. setting s f = 2.5741), one can obtain that a significantly long almost steady phase is followed by a nearly spontaneous crack opening around 1000 time steps (Fig. 18). We recall that the leaders were stopped at t = 180, which means that the abrupt opening of the crack around t = 900 is a consequence of an internal evolution of the system, though it is hardly visible since only very tiny changes in the positions of the elements (virtually invisible at the 'macro' level) occur for a long time. In the opinion of the authors, this behavior, which is experimentally well known [START_REF] Sih | Mechanics of fracture initiation and propagation: surface and volume energy density applied as failure criterion[END_REF], is quite interesting, especially considering the simplicity of the discrete model employed.

The sensitivity of the model with respect to very fine changes in the parameters suggests that non trivial instability issues may arise from the study of its homogenized form. This can be a very interesting aspect to be investigated, useful theoretical tools being available in recent literature (the interested reader can see, for instance, [START_REF] Luongo | [END_REF]53,54,55,56])

Conclusion and further objectives

In this work a new simple discrete system (introduced in [6]) with a centroid-based propagating interaction has been used in order to account for fracture phenomena in deformable bodies. Some results on pure deformation including a direct comparison (j) Fig. 18 'Almost steady' state followed by a quick and spontaneous evolution of the crack; s f = 2.5741, leaders are stopped at t = 180. Between t ≈ 500 and t ≈ 1000 the energy plot is almost flat, then an abrupt crack opening occurs (as before, we use empty red dots to indicate that the fictitious elements relative to the left side have an imposed null velocity, i.e. the left side has an imposed displacement equal to zero).

with standard Finite Element simulations are presented. A discrete equivalent of the deformation energy is also defined and used to prove the independence of the asymptotic configuration on the choice of the element from which the algorithm starts in case of spontaneous evolution. A fracture algorithm is then introduced, with the definition of a suitable adjacency matrix specifying the pairs of elements among which no interaction exists anymore because the two elements overcame a prescribed threshold in the actual distance. The discrete equivalent of an external stress is also defined in the uniaxial case, and used to compare the numerical results with Griffith's theory of fracture. Different examples of crack formation and growth are considered.

The results indicate that the proposed model is promising for developing new and computationally advantageous tools for the study of fracture, and that the model is rich enough to produce interesting behaviors such as a spontaneous crack evolution after a long almost steady state. Future investigations will concern how, by suitably weighting the contribution of the neighbors in the computation of the centroid, one can obtain different material models, and in particular different kinds of anisotropies. Moreover, one may try to generalize the model by suitably coupling the coordinates in order to account for Poisson's ratios different from 0 (while using directly the centroid-based algorithm we described this is the only possible value). Also interesting will be the study of more cases in which higher order interactions, reproducing behaviors typical of higher order continua, are considered. Finally theoretical results, among which the investigation of the homogenized limit of the considered system and the relative variational formulation, are needed in order to make further and sound progresses in the subject.

Proof of the independence of the asymptotic behavior on the choice of the pseudo-leader for 1D systems

In Section 2 we considered the problem of the uniqueness of the asymptotic configuration of the system evolving spontaneously (i.e., without leaders), independently on the selected pseudo-leader. As this represents an important objectivity requirement for the model, in the present appendix we deserve a more formal treatment to the 1D case.

Let us consider a 1-dimensional system S constituted by N elements mathematically represented by points in the real line, and an infinite set of discrete time steps T = {t 0 ,t 1 , ...,t m , ...}. The system will evolve exactly as explained in the 2D case, with the set of first neighbors being constituted, for every element, by the two closest ones. We call configuration of S a strictly increasing function C mapping i into x i ∈ R, with the integer i : 0 < i ≤ N indicating the Lagrangian label of the elements of the system. We indicate by C t j the configuration at time t j ∈ T, and by S ī the system evolving with the chosen pseudo-leader ( ī). We say that S ī asymptotically converges to the configuration C t ∞ , and write S ī → C t ∞ , if for every ε > 0 there ex-ists a positive integer M such that |x t m ix t ∞ i | < ε for every integer i : 0 < i ≤ N and every t m ∈ T : m > M (we can notice that the given definition, which is natural when considering finite systems, would correspond in the homogenized limit to a uniform functional convergence). Let C * be the set of configurations in which |x i+1x i | = 1 for every integer i such that 0 < i < N; obviously one has E D (C * ) = 0 for every C * ∈ C * (see Eq. ( 2)), and the representative C * is unique up to a translation. We have the following Proposition: A system S i evolving spontaneously tends to an asymptotic configuration C t ∞ , and it is C t ∞ ∈ C * ∀i ∈ N : 1 ≤ i ≤ N.

Proof: Let us first prove that the limit configuration exists. Then we will prove the Lemma by showing that, if C t ∞ is the asymptotic configuration, one has E D (C t ∞ ) = 0 (the thesis follows then immediately). We start by noticing that in every virtual time step one element or two elements which are not neighbors move to the centroid of their first neighbors, while all other elements do not move. Since the centroid locally minimizes the discrete energy, this means that E D (V h ,t j ) ≥ E D (V h+1 ,t j ) and E D (V h ,t j ) ≥ E D (V k ,t j+1 ) for all integers j, and for all h and k for which virtual configurations are defined. Since actual configurations of the system are a subset (preserving the order relation) of the virtual configurations, one has that E D (C t j ) (t j ∈ T) is a monotonically decreasing function of j. Since E D is non negative, this means that there exists lim m→∞ E D (C t m ) = E D (C t ∞ ). Let us now enumerate all virtual configurations (including the ones identified by definition with actual ones): V 1 ,V 2 , ...V n , .... Since E D converges, one has:

∆ E D (V h ) := E D (V h ) -E D (V h + 1) → 0 for m → ∞ (8)
Let ∆ s h be the sum of the moduli of the displacements of the elements that move in the virtual time step V h . Since ∆ E D (V h ) is a continuous function of ∆ s h vanishing if and only if ∆ s h = 0, Eq. ( 8) implies that lim h→∞ ∆ s h = 0, i.e. all virtual displacements of the elements (and thus the actual ones) tend to zero as time goes to infinity. Therefore, there exists a limit configuration C t ∞ . Let us now suppose by contradiction that E D (C t ∞ ) = Ẽ > 0 for some element ī chosen as the pseudo-leader. Let now d max and d min be respectively the maximum and minimum distance between (non fictitious) elements which are first neighbors in C t ∞ . Since Ẽ > 0, one has that at least one between the two inequalities d max > 1 and d min < 1 holds; without losing generality, we select the first possibility. Let ĩ and ĩ + 1 be the two Lagrangian elements such that lim m→∞

(|x t m ĩ+1 -x t m ĩ |) = d max (9)
Recalling that the centroid algorithm moves the element ĩ to the centroid of its neighbors, it easily follows from Eq.( 9) that: lim m→∞ (|x t m ĩx t m ĩ-1 |) = d max (10) Let x t m F be the coordinate of the left fictitious element of the system at time t m . Iterating the previous reasoning, one gets eventually lim m→∞ |x t m 0x t m F | = d max > 1 (11) and since by definition |x 0x F | = 1 at every time, the contradiction completes the proof.

Remark The result is easily generalized to system with second or higher order interaction. As mentioned in the paper, however, the generalization to 2D systems is not straightforward as in that case it is not true that the energy is decreasing.

Fig. 1

 1 Fig. 1 Graphical representation of neighbors: the first and second neighbors of the red element are the ones respectively in yellow and green.

Fig. 2

 2 Fig. 2 Graphical representation of the virtual configurations. In red the elements which are already in the position they will have in C t 1 .

Fig. 4 AFig. 5

 45 Fig. 4 A comparison between FE simulation of first (a) or second (b) gradient continua and first (c) or second (d) neighbors interaction with a similar external action.

Fig. 7

 7 Fig. 7 Nonlinearity of the system: in (a), a speed s = 1 length units per time step is imposed to the leader for an interval I of 10 time steps; in (b), the s = 0.001 and I = 10000. The two resulting configurations are different.

Fig. 8

 8 Fig.8Superposition of the discrete system with a standard FE continuum simulation involving first gradient/first-neighbors (left) and second gradient/second neighbors models (right).

(a) t = 1 (b) t = 30 (Fig. 9 ( 30 (Fig. 10 (Fig. 11

 1309301011 Fig.9 (a-d): spontaneuos evolution of the system from a prestressed configuration; the pseudoleader is the central element. (e): time history of the Discrete Energy Equivalent

(a) t = 1 (b) t = 30 (Fig. 12 Fig. 13

 1301213 Fig.12Spontaneous evolution starting from a severely prestressed configuration: the energy is not a monotonically decreasing function.

Fig. 14

 14 Fig.14Crack growth as a consequence of two leaders pushing from inside at 45 degrees with uniform velocity; with s f = 1.17 the crack length stabilizes around t = 160.

(a) t = 1 (

 1


Franck Techer
email: franck.techer@univ-nantes.fr

Christophe Jallais
email: christophe.jallais@ifsttar.fr

Yves Corson
email: yves.corson@univ-nantes.fr

Fabien Moreau
email: fabien.moreau@ifsttar.fr

Daniel Ndiaye
email: daniel.ndiaye@ifsttar.fr

Bruno Piechnick
email: bruno.piechnick@ifsttar.fr

Alexandra Fort
email: alexandra.fort@ifsttar.fr

Attention and driving performance modulations due to anger state : Contribution of electroencephalographic data

Keywords: Anger, Event-Related Potentials, Visual N1, Lateral variations, Car simulator, Attention

Introduction

Anger, which is a negative and highly arousing emotion, is commonly experienced while driving [START_REF] Mesken | Frequency, determinants, and consequences of different drivers' emotions: An on-the-road study using self-reports, (observed) behaviour, and physiology[END_REF]. A consequence of experiencing this emotional state can be outward expression of aggression through verbal or behavioral means [START_REF] Sullman | The expression of anger on the road[END_REF]. Such expression can include the use of one's vehicle to show own frustration to other road users or to frustrate another driver [START_REF] Deffenbacher | Anger, aggression and risky behavior: A comparison of high and low anger drivers[END_REF].

Furthermore, studies have linked driving anger to traffic rules infringements [START_REF] Nesbit | A quantitative review of the relationship between anger and aggressive driving[END_REF][START_REF] Stephens | Anger-congruent behaviour transfers across driving situations[END_REF], a reduction of lateral control [START_REF] Roidl | Emotional states of drivers and the impact on speed, acceleration and traffic violations-A simulator study[END_REF], and a reduction of following distances [START_REF] Garrity | Relations Among Personality Traits, Mood States, and Driving Behaviors[END_REF]. One frequently reported effect of driving anger is its ability to modulate driving style, leading to higher driving speeds and stronger accelerations [START_REF] Stephens | Anger-congruent behaviour transfers across driving situations[END_REF][START_REF] Roidl | Emotional states of drivers and the impact on speed, acceleration and traffic violations-A simulator study[END_REF][START_REF] Jeon | The effects of social interactions with in-vehicle agents on a driver's anger level, driving performance, situation awareness, and perceived workload[END_REF].

The impact of anger on driving performance involves more than just behavioral modifications.

Several negative effects also occur at a cognitive level, including the tendency to use a heuristic processing style, making drivers rely on superficial cues rather than on the significance of stimuli and leaving drivers unlikely to carefully analyze their environments [START_REF] Lerner | Portrait of the angry decision maker: how appraisal tendencies shape anger's influence on cognition[END_REF]. For example, increased driving speeds caused by anger seem to be predominantly mediated by a situational awareness deficit [START_REF] Jeon | The effects of social interactions with in-vehicle agents on a driver's anger level, driving performance, situation awareness, and perceived workload[END_REF]. Angry drivers are less likely to be aware of critical information or potential hazards on the road.

Although a number of studies have investigated the impacts of anger state on driving behavior, its impact on attentional processing while driving has rarely been studied. Basic research studies assessing the efficiency of attentional processing have revealed that emotions could improve the most basic perceptual abilities [START_REF] Phelps | Emotion facilitates perception and potentiates the perceptual benefits of attention[END_REF]. Concerning anger, Techer et al. [START_REF] Techer | Assessing the impact of anger state on the three Attentional Networks with the ANT-I[END_REF], used the Attention Network Test -Interactions (ANT-I) [START_REF] Callejas | Modulations among the alerting, orienting and executive control networks[END_REF] to investigate the impact of anger on attention. The ANT-I is a modified version of the Attention Network Test (ANT) [START_REF] Fan | Testing the Efficiency and Independence of Attentional Networks[END_REF] that can be used to study the influence of several contexts on attentional processing of neutral stimuli. This test, based on the model of attentional networks [START_REF] Posner | Orienting of attention[END_REF], allows the evaluation of the three attention sub-systems: the alerting, the orienting and the executive control networks. Techer et al. [START_REF] Techer | Assessing the impact of anger state on the three Attentional Networks with the ANT-I[END_REF] found a positive impact of anger on the alerting network, probably attributable to its high level of arousal. Such effect on attention may also be observed using physiological measures such as Event Related Potentials (ERP). This electrophysiological technique is based on the observation of averaged brain electrical signals after a repeated stimulus presentation so as to infer on the underlying cognitive processes [START_REF] Luck | An Introduction to the Event-Related Potential Technique[END_REF]. According to ERP literature, the first negative electrical peak after stimulus onset (N1) mainly reflects the perceptual processing stage of a target [START_REF] Luck | An Introduction to the Event-Related Potential Technique[END_REF][START_REF] Vogel | The visual N1 component as an index of a discrimination process[END_REF]. The amplitude of this component seems partly linked to the quantity of attentional resources allocated to sensorial processing. Additionally, this component is also sensible to the task. Its amplitude and latency are higher during a discrimination task than during a simple detection task [START_REF] Vogel | The visual N1 component as an index of a discrimination process[END_REF][START_REF] Fort | Task-dependent activation latency in human visual extrastriate cortex[END_REF]. For its part, the positive P3 peak is thought to inform on the cognitive processing of information and the inhibition of non-pertinent stimuli or responses [START_REF] Luck | An Introduction to the Event-Related Potential Technique[END_REF][START_REF] Polich | Updating P300: an integrative theory of P3a and P3b[END_REF].

During the ANT, both alerting and orienting signals could impact N1 amplitude [START_REF] Neuhaus | Event-related potentials associated with Attention Network Test[END_REF], the largest amplitude being elicited during trials with alerting and orienting cues. As for P3, changes in amplitude were observed when target stimulus was flanked by incongruent stimuli, suggesting an effect of the executive control network [START_REF] Neuhaus | Event-related potentials associated with Attention Network Test[END_REF]. Thus, changes in the efficiency of one attentional network may be observed with ERP measures.

Modulations of ERP amplitudes can also be observed during emotional information processing.

Literature suggests that the valence dimension of a stimulus seems to affect preferentially early components [START_REF] Rozenkrants | Affective ERP processing in a visual oddball task: Arousal, valence, and gender[END_REF][START_REF] Van Strien | Valence interacts with the early ERP old/new effect and arousal with the sustained ERP old/new effect for affective pictures[END_REF]. For example, unpleasant picture processing can lead to larger P1 amplitude when compared with pleasant picture processing. As for late components, their amplitude can increase when processing highly arousing stimuli [START_REF] Rozenkrants | Affective ERP processing in a visual oddball task: Arousal, valence, and gender[END_REF]. The effects on ERPs provoked by stimulus arousal are also observed during neutral information processing according to individual's arousal level [START_REF] Polich | Updating P300: an integrative theory of P3a and P3b[END_REF]. It has been interpreted as an arousal-related modulation in the quantity of available attentional resources. Therefore, in a driving context, physiological arousal and negative valence evoked by an angry mood may impact the electrophysiological response of drivers.

To our knowledge, no study has ever measured the impact of an anger state on ERPs while driving. Additionally, the relationship between anger and the alerting network [START_REF] Techer | Assessing the impact of anger state on the three Attentional Networks with the ANT-I[END_REF] is of particular interest due to the relative importance of this network during driving. Alerting signals are common and represent critical information for driving safety. Thus, a more efficient alerting network would allow the driver to take a greater advantage of the numerous alerting signals.

However, ERP technique is particularly sensitive to motor activation, which may be a limitation to its usability in a driving context. For that reason, Bueno et al. [START_REF] Bueno | An electrophysiological study of the impact of a Forward Collision Warning System in a simulator driving task[END_REF] developed a simulated driving task (consisting in a motorcycle following paradigm, on a straight rural road) compatible with ERPs collection. It successfully revealed several effects of a forward collision alerting system as well as cognitive distraction on ERPs. [START_REF] Bueno | An electrophysiological study of the impact of a Forward Collision Warning System in a simulator driving task[END_REF]. Thus, the use of this ADAS reduced P3 latency evoked by the motorcycle's braking lights.

The aim of this present study was to investigate the impact of an anger state on attention while driving, using the ERP technique, and its impact on driving. Anger was expected to influence ERPs following auditory alert and braking light of the leading motorcycle due to increased arousal level and the greater efficiency of the alerting network. Such effect would imply an increase of auditory N1, and visual N1 and P3 amplitudes. Additionally, anger was expected to disrupt driving performance as measured by reaction times, control of speed and lateral position.

Participants

Thirty-three participants (19 females) aged between 25 and 40 (M = 32.3; SD = 5.5) were involved in this study and received a financial compensation. They reported a normal or corrected to normal vision, no neurologic disease and no medical treatment. Every participant was right-handed and had more than three years of driving experience. The research protocol was carried out in accordance with The Code of Ethics of the World Medical Association.

Material

Mood induction and measurement

Two experimental sessions were used. In the Anger session, participants were induced using the autobiographical recall procedure [START_REF] Jallais | Inducing changes in arousal and valence: Comparison of two mood induction procedures[END_REF]. They had ten minutes to recall and write down a personal event during which they experienced anger, and were encouraged to provide as many details as they could. In the Control session, participants were not induced. In order to keep their natural mood, they had ten minutes to complete questionnaires about their driving habits.

Mood states were measured using a modified version of the Brief Mood Introspection Scale (BMIS) [START_REF] Corson | Emotions and False Memories Valence or Arousal?[END_REF]. The BMIS is a 16-item self-report questionnaire in which each adjective refers either to anger, happiness, sadness or calmness. It is rated on a 7-point scale, ranging from "not at all" to "absolutely", providing a score for the valence and the arousal dimensions of the emotional state.

Apparatus

The experimental scenario was presented to participants in a driving simulator composed of a 24" screen, an adjustable car seat, a steering wheel and three pedals.

Electroencephalographic data was collected using the Biosemi ActiveTwo system® sampled at 1024Hz. Electrodes were placed on an electrode cap which was organized according to the international 10-20 system. Two mastoids electrodes were also placed on Ma1 and Ma2, and one EOG electrode was placed near the right eye. The reference electrode was placed on the nose. Event related potentials were extracted offline using the ELAN software [START_REF] Aguera | ELAN: A Software Package for Analysis and Visualization of MEG, EEG, and LFP Signals[END_REF].

Procedure

Each participant completed an Anger and a Control experimental session with a weekly interval.

The session order was counterbalanced between the participants. Each session followed the same structure. After completing the informed consent form, the EEG recording apparatus was set on participants. They were seated in the simplified driving simulator at 90cm eye-screen distance and completed a 13.5-minute training scenario based on the procedure of Bueno et al. [START_REF] Bueno | An electrophysiological study of the impact of a Forward Collision Warning System in a simulator driving task[END_REF], in which they were instructed to follow a motorcycle. The driving task was designed to require a low level of motor activation so as to improve the quality of ERP collection.

Participants drove on a straight rural road with no crosswinds and no regular steering was needed. A speed limiter ensured that participants drove at a consistent speed. They were instructed to maintain the accelerator pedal pushed as long as the motorcycle did not brake.

Every time they saw the brake lights of the motorcycle, they had to release the accelerator pedal as quickly as they could. They were also asked to adjust their pace in order to keep a safe distance to the lead vehicle. When participants reached the 70km/h speed limit, the first of the 90 trials started. Each trial started by a random duration comprised between 400 and 1600ms.

In 80 of the 90 trials, a 500ms alerting auditory signal indicated to the participant that the lead vehicle would brake imminently. In the 10 remaining trials, no alerting signal was displayed.

After 500ms, the motorcycle braked for a duration of 2 seconds, decelerating from 70km/h to 35km/h. In 15 of those trials, the deceleration of the lead vehicle was more important, shifting from 70 km/h to 10 km/h, so as to break the monotony of the task and prevent the driver from losing concentration. The brake lights of the motorcycle were lit throughout the 2000ms corresponding to the braking maneuver. At the end of the trial, the motorcycle accelerated to get back to the initial inter-vehicular distance and participants had to push the accelerator pedal at its maximum again. Between trials, the motorcycle was programmed to keep a distance corresponding to 2 seconds at a 70km/h speed so as to reproduce exactly the same driving conditions for each trial. The next trial started immediately. This trial structure ensured that each participant was placed in a comparable situation with one motorcycle braking on average every nine seconds. To counter the high braking predictability, the random time implemented at the beginning of each trial made the interval between two targets irregular, ranging from 7.8 seconds to 10.2 seconds. After the training scenario and before the Mood Induction Procedure (MIP), participants had to fill in the first BMIS questionnaire (Moment 1), then followed the MIP corresponding to the current session. Finally, after the MIP, they had to fill in the second BMIS (Moment 2) followed by the experimental scenario which was identical to the training scenario.

Analyses Event-related potentials

Auditory and visual ERPs were computed for trials with auditory alert. They were filtered offline at 1Hz -30Hz, computed on a time window from -200ms to 800ms and baselined from 200ms to 0ms pre-stimulus window. Each ERP component was analyzed at the traditional locations observed in the literature [START_REF] Luck | An Introduction to the Event-Related Potential Technique[END_REF].Auditory N1 analyzed windows ranged from 80ms to 160ms after auditory alert onset at Cz, FC1, FC2 and Fz. Visual N1 window was set from 160ms to 230ms after the braking lights onset and was analyzed on IMA, IMB, P7, P8, PO3, PO4, O1 and O2 sites. Finally, P3 was recorded during the 250ms to 400ms window for Pz, P3, P4, PO3, PO4, O1 and O2 electrode sites.

All ERP analysis were conducted using ELAN software and ERP waves visualization was conducted using ERPA software [START_REF] Aguera | ELAN: A Software Package for Analysis and Visualization of MEG, EEG, and LFP Signals[END_REF]. For each component, according to its polarity, the most negative or positive value on each electrode site was recorded. Those values were then averaged between electrodes to obtain the maximum or minimum ERP component amplitude. The latencies of negative or positive peaks were recorded and averaged between electrodes to obtain the latency of the ERP component. Considering the little number of trials without auditory alert, these trials were not included in the analyses and were not compared to cued trials. Thereby the analysis was not on the impact of anger on the alerting effect, but rather on its impact on an alerting cue processing during driving.

Driving performance

Several driving indicators were recorded to analyze driving performance. Response times were calculated as the duration between braking lights onset and participants' reaction. Number of anticipations was also recorded as the number of trials for which participants reacted before the braking lights. Finally, intervehicular distance and standard deviation of lateral position (SDLP) in meters were recorded to give an indication of longitudinal and lateral control of participants.

Results

Mood measurement

The emotional state was assessed before and after the MIP of each session. A repeated measure ANOVA was carried out with the Session (Anger/Control) and the Moment (Moment 1/Moment_2) as within-subject factors, on the mean anger ratings or the mean arousal ratings as dependent variables.

The mean anger ratings analysis revealed a main effect of the Session (F(1, 32) = 11.08, p < .001, η² = .258) and a main effect of the Moment (F(2, 64) = 25.08, p < .001, η² = .439). A significant Session × Moment interaction was also found (F(2, 64) = 20.07, p < 001, η² = .385).

Planned comparisons revealed no difference before induction (F(1, 32) < 1, n.s., η² = .002). At Moment 2, a significant difference was found (F(1, 32) = 14.75, p < .001, η² = 0.357) with a higher anger rating after induction during the Anger session (M = 2.11, SD = 1.22) than during the Control session (M = 1.18, SD = 0.39). Similar effects were observed for the mean arousal ratings. The analysis revealed a main effect of the Session (F(1, 32) = 10.81, p < .01, η² = .252) and a main effect of the Moment (F(2, 64) = 17.22, p < .001, η² = .350). A significant Session × Moment interaction was also found (F(2, 64) = 6.28, p < .01, η² = .164). Planned comparisons revealed no difference in Moment 1 (F(1, 32) = 3.62, n.s., η² = 0.102). After induction, a significant difference was found (F(1, 32) = 16.28, p < .001, η² = .384) with a higher arousal rating in Moment 2 during the Anger session (M = 4.10, SD = 0.39) than during the Control session (M = 3.80, SD = 0.32).

Table 1 Anger and Arousal mean ratings (SD) before and after induction according to the session

Anger ratings

Arousal ratings

Discussion

This study assessed the influence of an anger state on attentional processing while driving, as measured by Event Related Potentials (ERPs) and its impact on driving performance. Anger was expected to influence attentional processes as a result of high levels of arousal and the ERP mean amplitudes and latencies for each measured ERP component according to the session (*=p<.05) alerting network efficiency associated therewith [START_REF] Techer | Assessing the impact of anger state on the three Attentional Networks with the ANT-I[END_REF]. Additionally, anger was expected to alter driving performance by disrupting reaction times as well as control of speed and lateral position.

First of all, it can be noted that N1s latencies are quite late as in the Bueno's studies [START_REF] Bueno | An electrophysiological study of the impact of a Forward Collision Warning System in a simulator driving task[END_REF] compared to basic studies. This may be due to our dynamic environment associated to low salience of stimuli because of the auditory and visual context within the car simulator. Further research will be necessary to test this hypothesis.

Concerning the induction, results confirm that the MIP was efficient to induce participants in anger since after the MIP, participants reported higher anger and arousal ratings.

However, ERP data did not support our initial hypothesis concerning an increase of the auditory N1, and visual N1 and P3 amplitudes. On the contrary, during the anger session, the amplitude of visual N1 was smaller than in the control session.

The generation of mind wandering (MW) could explain this unexpected result. MW is a reorientation of the attention to internal thoughts not linked to the driving task [START_REF] Smallwood | The science of mind wandering: empirically navigating the stream of consciousness[END_REF]. According to the literature, a N1 amplitude reduction could be linked with a perceptual decoupling during MW episodes [START_REF] Smallwood | The science of mind wandering: empirically navigating the stream of consciousness[END_REF][START_REF] Chaparro | Behavioral and neural correlates of spontaneous attentional decoupling : towards an understanding of mind wandering[END_REF][START_REF] Kam | Slow fluctuations in attentional control of sensory cortex[END_REF]. A strong bidirectional relationship exists between negative emotions and MW [START_REF] Smallwood | The science of mind wandering: empirically navigating the stream of consciousness[END_REF][START_REF] Smallwood | Shifting moods, wandering minds: Negative moods lead the mind to wander[END_REF]. On the one hand, a negative mood can lead to MW [START_REF] Smallwood | Shifting moods, wandering minds: Negative moods lead the mind to wander[END_REF] and thinking about negative life events can be used to induce a negative emotional state [START_REF] Jallais | Inducing changes in arousal and valence: Comparison of two mood induction procedures[END_REF]. On the other hand, the emergence of MW can also depend on the task. Normally, it appears foremost when a task is highly repetitive and monotonous [START_REF] Berthié | The restless mind while driving: drivers' thoughts behind the wheel[END_REF], which corresponds to the protocol of the current study. Therefore, it is possible that our monotonous task associated to the autobiographical induction procedure had generated MW. A reduction of the visual N1 amplitude may indicate a reduction of the attention allocated to the sensorial processing of braking lights through an attenuation of sensorial sensitivity [START_REF] Chaparro | Behavioral and neural correlates of spontaneous attentional decoupling : towards an understanding of mind wandering[END_REF]. Moreover, the procedure used here was previously developed to show an impact of cognitive distraction on ERPs [START_REF] Bueno | An electrophysiological study of the impact of a Forward Collision Warning System in a simulator driving task[END_REF]. These authors observed that a cognitive distraction reduced the visual N1 and P3 amplitudes. However, our participants

were not asked to be engaged in a cognitive distraction task and as a result of the fog used in the driving environment, only stimuli used for the following task were displayed on the simulator screen. Therefore, our results seem to be consistent with the idea that drivers were distracted due to their attention being allocated towards internal thoughts.

The presence of MW could also give another point of view to the fact that induced anger in the present experiment did not show the expected effect on the alerting network as observed by Techer et al. [START_REF] Techer | Assessing the impact of anger state on the three Attentional Networks with the ANT-I[END_REF] using the ANT-I. This task requires a continuous attentional focus which could not be suitable for the generation of MW. On the contrary, as mentioned above, the driving task used in the present study was monotonous and favorable to MW which could be responsible for the reduced amplitude of visual N1 and for minimizing or suppressing the alerting effect.

The N1 effect may reflect a lack of sensitivity to stimuli present in the environment [START_REF] Chaparro | Behavioral and neural correlates of spontaneous attentional decoupling : towards an understanding of mind wandering[END_REF].

Therefore, angry drivers may process critical information in a more superficial way which is consistent with a reactivity reduction in case of unexpected hazards [START_REF] Stephens | Drivers Display Anger-Congruent Attention to Potential Traffic Hazards[END_REF].

Concerning the expected effects on driving performance, we observed no effect of anger on response times and speed control. It may be a direct consequence of the highly controlled situation in which participants were placed, reducing the possible sources of variability.

However, analyses revealed that participants had a degraded lateral control of their vehicle during the anger session, as indicated by the increased standard deviation of lateral position (SDLP). Such increase in variations of lane positioning was previously observed for drivers induced in anger by the driving situation [START_REF] Stephens | Anger-congruent behaviour transfers across driving situations[END_REF]. It seems also consistent with the literature stating that angry driving could lead to higher lateral accelerations [START_REF] Roidl | Emotional states of drivers and the impact on speed, acceleration and traffic violations-A simulator study[END_REF]. In both of these studies, it was interpreted to be a result of drivers engaging in aggressive behaviors to avoid being impeded.

However, in our study, the lead motorcycle was not impeding and drivers were instructed to follow it. Moreover, no active steering was required to complete the driving task. Therefore, the increased variations in lateral position may be a reflection of changes occurring at an attentional level. As previously mentioned, ERP data suggested that drivers' attention was impacted by anger in a comparable way to mind-wandering. However, previous studies using comparable driving conditions seem to indicate that drivers involved in distractive tasks or in MW episodes tend to have a reduced SDLP [START_REF] He | Mind Wandering Behind the Wheel Performance and Oculomotor Correlates[END_REF][START_REF] Lemercier | Inattention behind the wheel: How factual internal thoughts impact attentional control while driving[END_REF]. SDLP is also observed with an increase of cognitive workload [START_REF] Cooper | The Impact of Eye Movements and Cognitive Workload on Lateral Position Variability in Driving[END_REF]. It may indicate that drivers prioritized other driving sub-tasks more demanding than lateral control [START_REF] Cooper | The Impact of Eye Movements and Cognitive Workload on Lateral Position Variability in Driving[END_REF]. Thus, the increased SDLP observed here does not seem to be induced by MW. Then, the other possible explanation could be linked to the high levels of arousal which defines anger. As stated by Logan and Crump [START_REF] Logan | The Left Hand Doesn't Know What the Right Hand Is Doing The Disruptive Effects of Attention to the Hands in Skilled Typewriting[END_REF], when a conscious monitoring of a low-demanding task is possible, the performance obtained for this task can be disrupted.

Thus, if arousal can rise the amount of available attentional resources [START_REF] Kahneman | Attention and Effort[END_REF], the participants in the anger session may have allocated attention to their internal thoughts for one part and to a conscious monitoring of the vehicle lateral control for another part. This last can explain the increased variation of lateral position provoked by anger. However, Ünal et al. [START_REF] Ünal | Driving with music: Effects on arousal and performance[END_REF] asked to their participants to listen to music of their choice while performing the same car following task. They observed that a high level of arousal could reduce the standard deviation of lateral position [START_REF] Ünal | Driving with music: Effects on arousal and performance[END_REF]. This effect was found with a musical playlist set up by the participants.

Unfortunately Ünal et al. [START_REF] Ünal | Driving with music: Effects on arousal and performance[END_REF] did not asses the valence dimension in their study. This issue raises the question about the role of valence. Thus, future research should take into account this dimension so as to disentangle the effects of moods defined by a similar level of arousal on lateral variations.

Considering the electrophysiological and behavioral results, the effects of anger observed here may not be explained only by the presence of MW episodes or high levels of arousal, but from a combination of these two factors: one having an effect on the sensory processing and the other one on the attentional resources management revealed by the lateral control of the vehicle.

Consequently, to identify the effects of a such emotion, several indicators have to be taken into account. Previous studies revealed that perceptual decoupling from visual stimuli is not a continuous phenomenon [START_REF] Chaparro | Behavioral and neural correlates of spontaneous attentional decoupling : towards an understanding of mind wandering[END_REF]. That is why drivers may constantly switch between a physiological activity corresponding to aroused states and a perceptual decoupling corresponding to MW. Further research on this topic is needed to better understand the influence of anger which may have a more complex impact on ERPs than expected.

In any case, the results obtained from the present study open new investigation opportunities regarding the impact of emotional states on attention during a driving task. Future studies may focus on the effects of other emotional states in different driving tasks. As an ultimate goal, those results may be used in the design of future driving assistance systems in order to provide support corresponding to drivers' current internal state.

Figure 1 :

 1 Figure 1: Grand average of the visual ERP on electrode P8 according to the session

Acknowledgements: Authors want to thank the IFSTTAR-COSYS-LEPSiS team for their support with the driving simulator, and Kristy-Anne Leith for her suggestions to improve English language.