Thibaut Weller
email: thibaut.weller@univ-montp2.fr

Christian Licht
email: clicht@univ-montp2.fr

Asymptotic modeling of thin linearly quasicrystalline plates Analyse asymptotique de plaques minces quasicristallines linéaires

published or not. The documents may come

Asymptotic modeling of thin linearly quasicrystalline plates

Introduction

As a new phase of solid matter, quasicrystals -henceforth abbreviated as QCs -were first discovered and studied by D. Shechtman. At the time of their discovery, in 1982, they presented a fundamental challenge to theoretical physicists to re-examine traditional concepts and devise new analytical methods for determining their physical properties. Indeed, it had been generally assumed that pentagonal symmetry is impossible in solids, that Bragg diffraction peaks imply periodic translational order, that all pure solids are either crystalline or glassy, and that the (classically) stable state of a solid is crystalline. The discovery of icosahedral phases in metallic alloys by Shechtman et al. [START_REF] Shechtman | Metallic phase with long-range orientational order and no translational symmetry[END_REF] has compelled theorists to reconsider these assumptions and to confront many new problems (see [START_REF]The Physics of Quasicrystals[END_REF] for a comprehensive introduction to the subject). Therefore, the physical properties such as the structural, electronic, magnetic, optical and thermal characteristics of QCs have been investigated intensively. Most of these properties combine effectively to give interesting technological applications which have been protected recently by several patents [START_REF] Blaaderen | Quasicrystals from nanocrystals[END_REF][START_REF] Dubois | Useful Quasicrystals[END_REF][START_REF] Dubois | New prospects from potential applications of quasicrystalline materials[END_REF]. As an example among others, the combination of such kind of properties as high hardness, low friction and corrosive resistance with bio-compatibility is very promising for introducing QCs in surgical applications. It is thus of major technological interest to provide efficient modelings of quasicrystalline structures. Here we perform a mathematical modeling of linear quasicrystalline plates by regarding the thickness as a small parameter denoted by ε. Then we study the behavior of the solution of the physical problem as ε tends to 0. We show that depending on the type of boundary conditions, 26 different models indexed by a triplet p = (p 1 , p 2 , p 3 ) ∈ {0, 1, 2} 3 appear at the limit. Comparing to our previous studies devoted to the mathematical modeling of thin plates in the framework of multi-physical couplings [START_REF] Weller | Asymptotic modeling of thin piezoelectric plates[END_REF][START_REF] Weller | Mathematical modeling of piezomagnetoelectric thin plates[END_REF][START_REF] Weller | Asymptotic modeling of piezoelectric plates with electric field gradient[END_REF], we are here in front of a real 'explosion' of the number of limit models. This multiplication of models, however, has its roots in the very structure of QCs. As shown in [START_REF] Ding | Generalized elasticity theory of quasicrystals[END_REF], the constitutive law of quasicrystalline media present a coupling between two different kinds of displacement fields: the phonon field, denoted by u, and the phason field, denoted by w (see [START_REF]The Physics of Quasicrystals[END_REF]). Roughly speaking, the existence of the phason field w is due to the fact that the dimension of the reciprocal space of a three-dimensional QC is greater than three. Whereas phonons are related to translation of atoms and therefore to classical elastic displacements, phasons are associated with atomic rearrangements and appear in the constitutive equations only through their gradients and not their symmetrized gradients. As it will be explained infra, this is the main cause of the huge amount of models we obtain for a single quasicrystalline thin plate. For all the models we get, we are able to express the limit constitutive law as a Schur complement of the tensor accounting for the coupling between phonon and phason fields. Finally, we study the influence of the icosahedral order on all our models and show that a decoupling between phasons and phonons may appear in the structure.

Setting the problem

We will denote phonon fields by the letters u, v and v while the phason fields will be denoted by w, ψ and ψ .

Depending on the nature of our formulation, these letters and some other symbols may be indexed by ε which stands for the thickness of the plate, regarded as a parameter. Classically, the linearized strain tensor is written e(u) ∈ S 3 where S N indicates the set of all N × N real and symmetric matrices. Used as indices, letters i and j take their values in {1, 2, 3} while α and β take their values in {1, 2}. When there is no ambiguity, the convention of summation over repeated indices is understood. We recall that 2e ij (u) = ∂ i u j + ∂ j u i where the symbol ∂ i refers to the partial derivative with respect to the i-th coordinate. The gradient of a phason field w will be denoted by ∇ w ∈ M 3 where M N is the set of all N × N real matrices. The components of ∇ w are given by (∇ w) ij = ∂ j w i . Given a physical state (u, w) living in a quasicrystalline medium, we therefore have (e(u), ∇ w) ∈ H, with:

H = S 3 × M 3 (1)
In the sequel, any element of H will be denoted by h = (e, g). For the sake of simplicity, the classical symbol • will stand for the inner product in H, S 3 , M 3 and R 3 . The set of all linear mappings from a space X into a space Y is denoted by L(X, Y ) and if X = Y we simply write L(X). In the sequel, for all domain D of R N , H 1 Γ (D) refers to the subset of the Sobolev space H 1 (D) whose elements vanish on a part Γ of the boundary ∂ D of D.

The reference configuration of a linearly quasicrystalline thin plate is the closure in R 3 of the set Ω ε = ω × (-ε, ε) whose outward unit normal is n ε . Here ε is a small positive number and ω is a bounded domain of R 2 with a Lipschitz boundary ∂ω. The lateral part of the plate ∂ω ×

(-ε, ε) is denoted Γ ε lat , while Γ ε ± = ω × {±ε} refers to the upper or lower face, respectively. Let (Γ ε uD , Γ ε uN ), (Γ ε w D,i , Γ ε w N,i ) i=1,2,
i = w ε 0 i on Γ ε w D,i , where w ε 0 = (w ε 0 1 , w ε 0 2 , w ε 0 3
) is a smooth enough given field defined in Ω ε . Moreover, the plate is subjected to body forces f ε and generalized body forces g ε , but also to surface forces F ε and generalized tractions G ε of components G ε i on Γ ε uN and Γ ε w N,i , respectively. Finally, we assume that Γ ε uD = γ 0 × (-ε, ε), with γ 0 ⊂ ∂ω. The equations determining the physical state s ε = (u ε , w ε ) at equilibrium are (see [START_REF] Ding | Generalized elasticity theory of quasicrystals[END_REF][START_REF] Fan | Mathematical Theory of Elasticity of Quasicrystals and Its Applications[END_REF] for example):

P Ω ε : ⎧ ⎪ ⎨ ⎪ ⎩ div σ ε + f ε = 0 in Ω ε , σ ε n ε = F ε on Γ ε uN , u ε = 0 on Γ ε uD div τ ε + g ε = 0 in Ω ε , τ ε ij n ε j = G ε i on Γ ε w N,i , w ε i = w ε 0 i on Γ ε w D,i , i = 1, 2, 3 σ ε , τ ε = Q ε (x) e u ε , ∇ w ε in Ω ε
where u ε , w ε , σ ε and τ ε respectively stand for the phonon field, the phason field, the phonon stress tensor and the phason stress tensor. The operator Q ε is an element of L(H) such that:

σ ε ij = C ε ijkl e kl u ε + R ε ijkl ∇ w ε kl τ ε ij = R ε kli j e kl u ε + K ε ijkl ∇ w ε kl (2)
In these constitutive equations, C ε ijkl , K ε ijkl and R ε ijkl stand for the phonon, phason and phonon-phason coupling coefficients, respectively. As an example, for icosahedral quasicrystals, the stored energy function takes the form (see [START_REF] Levine | Elasticity and dislocations in pentagonal and icosahedral quasicrystals[END_REF]):

Q ε e u ε , ∇ w ε • e u ε , ∇ w ε = 1 2 λ ε e ii u ε e ii u ε + μ ε e ij u ε e ij u ε + 1 2 K ε 1 ∂ i w ε j ∂ i w ε j + 1 2 K ε 2 ∂ ε 3 w ε 3 2 -∂ 1 w ε 3 2 -∂ 2 w ε 1 2 -∂ 3 w ε 1 2 -∂ 3 w ε 2 2 + 2 ∂ 1 w ε 3 + ∂ 3 w ε 1 ∂ 1 w ε 1 -2 ∂ 1 w ε 3 -∂ 3 w ε 1 ∂ 2 w ε 2 -2 ∂ 2 w ε 3 + ∂ 3 w ε 2 ∂ 1 w ε 2 -2 ∂ 2 w ε 3 -∂ 3 w ε 2 ∂ 2 w ε 1 + r ε e 11 u ε + e 22 u ε -2 e 33 u ε ∂ 3 w ε 3 + e 11 u ε -e 22 u ε ∂ 1 w ε 1 + ∂ 2 w ε 2 + ∂ 3 w ε 1 + 2 e 12 u ε ∂ 1 w ε 2 -∂ 2 w ε 1 -∂ 3 w ε 2 + 2 e 13 u ε ∂ 1 w ε 1 + ∂ 1 w ε 2 -∂ 2 w ε 2 + 2 e 23 u ε ∂ 2 w ε 3 -∂ 1 w ε 2 -∂ 2 w ε 1 ( 3
)
where λ ε , μ ε , r ε , K ε 1 and K ε 2 are five physical coefficients which may depend on x and see [START_REF] Ding | Generalized elasticity theory of quasicrystals[END_REF] for the corresponding form of the tensors C ε , R ε and K ε .

We are looking for the state (u ε , w ε ) living in the quasicrystalline plate at equilibrium. For this purpose, we first make the following regularity hypothesis on the exterior loading:

(H 1 ) : f ε , g ε , F ε , G ε i ∈ L 2 Ω ε 3 × L 2 Ω ε 3 × L 2 Γ ε uN 3 × L 2 Γ ε w N,i , i = 1, 2, 3 w ε 0 ∈ H 1 Ω ε
Now, on the space of physical states:

V ε = H 1 Γ ε uD Ω ε 3 × H 1 Γ ε w D,1 Ω ε × H 1 Γ ε w D,2 Ω ε × H 1 Γ ε w D,3 Ω ε (4)
we define a bilinear form m ε :

m ε (ξ, ζ ) = m ε (v, ψ), v , ψ = Ω ε Q ε e(v), ∇ψ • e v , ∇ψ dx ε ( 5
)
and a linear form L ε :

L ε (ξ ) = L ε (v, ψ) = Ω ε f ε • v + g ε • ψ dx ε + Γ ε uN F ε • v ds ε + 3 i=1 Γ ε wN,i G ε i ψ i ds ε (6)
The coupled physical problem then takes the form

P Ω ε : Find s ε = u ε , w ε ∈ 0, w ε 0 + V ε such that m ε s ε , ξ = L ε (ξ ), ∀ξ ∈ V ε
Thus, with the additional and realistic assumptions of boundedness of C ε , R ε and K ε and of uniform ellipticity of C ε and K ε :

(H 2 ) : Q ε ∈ L ∞ Ω ε , L(H) , ∃κ ε > 0; Q ε x ε h • h κ ε |h| 2 H , ∀h ∈ H, a.e. x ε ∈ Ω ε
the Stampacchia theorem (cf. [START_REF] Brézis | Functional Analysis, Sobolev Spaces and Partial Differential Equations[END_REF]) implies the Theorem 2.1. Under assumptions (H 1 )-(H 2 ), the problem P(Ω ε ) has a unique solution.

To derive simplified and accurate models, the true question is to study the behavior of s ε when ε, regarded as a parameter, tends to zero.

The different models

We will show that 26 different limit behaviors appear according to both the type and magnitude of boundary conditions in P(Ω ε ). These limit behaviors will be indexed by a triplet p = (p 1 , p 2 , p 3 ) of {1, 2, 3} 3 \ {0, 0, 0}.

Classically (see [START_REF] Ciarlet | [END_REF]) we come down to a fixed open set Ω = ω × (-1, 1) through the mapping π ε :

x = (x 1 , x 2 , x 3 ) ∈ Ω → π ε x = (x 1 , x 2 , εx 3 ) ∈ Ω ε (7)
Also, we drop the index ε for the images by (π ε ) -1 of the geometric sets defined in the preceding section. To get physically meaningful results, we have to make various kinds of assumptions. They deal with (i) the electromechanical coefficients:

(H 3 ) : Q ε π ε x = Q (x) with Q ∈ L ∞ Ω, L(H) , ∃κ > 0; Q (x)h • h κ|h| 2
H , ∀h ∈ H, a.e. x ∈ Ω

(ii) the electromechanical loading:

(H 4 ) : ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ∃( f , g, F , G i ) ∈ L 2 (Ω) 3 × L 2 (Ω) 3 × L 2 (Γ uN ) 3 × L 2 Γ i w N , i = 1, 2, 3 f ε α π ε x = ε f α (x), f ε 3 π ε x = ε 2 f 3 (x), ∀x ∈ Ω F ε α π ε x = ε 2 F α (x), F ε 3 π ε x = ε 3 F 3 (x), ∀x ∈ Γ uN ∩ Γ ± F ε α π ε x = ε F α (x), F ε 3 π ε x = ε 2 F 3 (x), ∀x ∈ Γ uN ∩ Γ lat g ε i π ε x = ε 2-p i g i (x), ∀x ∈ Ω G ε i π ε x = ε 3-p i G i (x), ∀x ∈ Γ w N,i ∩ Γ ± G ε i π ε x = ε 2-p i G i (x), ∀x ∈ Γ w N,i ∩ Γ lat w ε 0 i π ε x = ε p i w 0 i (x), ∀x ∈ Ω (iii)
the boundedness of the "work of the exterior loading":

(H 5 ) : ⎧ ⎪ ⎨ ⎪ ⎩ p i = 0 : w 0 i = 0 p i = 1 : w 0 i does not depend on x 3 p i = 2 :
-the closure δ of the projection of Γ w D,i on ω coincides with ω

-either G i = 0 on Γ w N,i ∩ Γ lat or Γ w N,i ∩ Γ lat = ∅
Also, with the true physical state s ε = (u ε , w ε ) defined on Ω ε , we associate a scaled physical state s p (ε) = (u(ε), w p (ε)) defined on Ω by:

u ε α x ε = ε u(ε) α (x), u ε 3 x ε = u(ε) 3 (x), w ε i x ε = ε p i w p (ε) i (x), ∀x ε = π ε x ∈ Ω ε (8)
so that s p (ε) is the unique solution of the following mathematical problem:

P(ε, Ω) p : Find s p (ε) ∈ (0, w 0 ) + V such that: m p (ε) s p (ε), ξ = L(ξ ), ∀ξ ∈ V = H 1 Γ uD (Ω) 3 × H 1 Γ w D,1 (Ω) × H 1 Γ w D,2 (Ω) × H 1 Γ w D,3 (Ω)
equivalent to the genuine physical one, where:

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ m p (ε)(s, ξ) = Ω Q (x) k p (ε, s) • k p (ε, ξ) dx k p (ε, ξ) = k p ε, (v, ψ) = e(ε, v), ∇ p (ε, ψ) e αβ (ε, v) = e αβ (v), e α3 (ε, v) = ε -1 e α3 (v), e 33 (ε, v) = ε -2 e 33 (v) ∇ p (ε, ψ) iα = ε p i -1 ∂ α ψ i , ∇ p (ε, ψ) i3 = ε p i -2 ∂ 3 ψ i L(ξ ) = L(v, ψ) = Ω ( f • v + g • ψ) dx + Γ uN F • v ds + 3 i=1 Γ wN,i G i ψ i ds (9)
As explained in our previous papers [START_REF] Weller | Asymptotic modeling of thin piezoelectric plates[END_REF][START_REF] Weller | Mathematical modeling of piezomagnetoelectric thin plates[END_REF][START_REF] Weller | Asymptotic modeling of piezoelectric plates with electric field gradient[END_REF][START_REF] Weller | Asymptotic modeling of linearly piezoelectric slender rods[END_REF], the signs of the various powers of ε in the components of k p (ε, ξ) induce an orthogonal decomposition of H in subspaces H p with ∈ {-, 0, +}, which allows a comprehensive description of our plate models in any theoretically admissible quasicrystallographic classes. We denote by h p the projection on H p of any element h of H. As an example, for p = (0, 1, 2), we have:

H - p = (e, g) ∈ H; e αβ = 0, g 2α = 0, g 3i = 0 H 0 p = (e, g) ∈ H; e i3 = 0, g 1i = 0, g 23 = 0, g 3α = 0 H + p = (e, g)
∈ H; e ij = 0, g 1i = 0, g 2i = 0, g 33 = 0 Then, for a given triplet p, the operator Q can be decomposed in nine elements Q p ∈ L(H p , H p ) with , ∈ {-, 0, +}.

Hypothesis (H 3 ) on the coefficients implies that Q 00 p and Q -- p are positive operators on H 0 p and H - p , respectively. Therefore, the Schur complement:

Q p = Q 00 p -Q 0- p Q -- p -1 Q -0 p ( 10
)
is an element of L(H 0 p ). Note that neither Q 00 p nor Q p are necessarily symmetric, but nevertheless:

κ h 0 p 2 H Q p (x) h 0 p • h 0 p , ∀h 0 p ∈ H 0 p , a.e. x ∈ Ω ( 11
)
This is implied by the coercivity of Q (see (H 3 )) and by the fundamental relation:

(Q h) - p = h + p = 0 ⇒ Q p h 0 p = (Q h) 0 p and Q p h 0 p • h 0 p = Q h • h (12)
The key point of the asymptotic study is to show that if k p is the limit (in a suitable topology) of k p (ε, s p (ε)), then (Q k p ) - p = (k p ) + p = 0, which exhibits Q p as the operator governing the 26 limit problems for any triplet p ∈ {0, 1, 2} 3 \ {0, 0, 0}.

The limit space of phonon fields will be the space of Kirchhoff-Love displacements defined by V K L = {v ∈ H 1 Γ uD (Ω) 3 ; e i3 (v) = 0} while the limit phonon fields need the introduction of nine spaces Ψ 0

,i = {0}, Ψ 1,i = {φ ∈ H 1 Γ w D,i (Ω); ∂ 3 φ = 0} and Ψ 2,i = {φ ∈ H 1 ∂ 3 (Ω); φ = 0 on Γ w D,i ∩ Γ ± } where H 1 ∂ 3 (Ω) = {φ ∈ L 2 (Ω); ∂ 3 φ ∈ L 2 (Ω)}.
The space of limit physical states is defined as

V p = V K L × Ψ p 1 ,1 × Ψ p 2 ,2 × Ψ p 3 ,3 .
Under assumptions (H 1 )-(H 5 ), we have the following convergence result:

Theorem 3.1. Let p ∈ {1, 2, 3} 3 \ {0, 0, 0}, X 0 = X 1 = H 1 (Ω) and X 2 = H 1 ∂ 3 (Ω). When ε tends to 0, the family (s p (ε)) ε>0 of the unique solutions of P(ε, Ω) p strongly converges in X = H 1 Γ uD (Ω) 3 × X p 1 × X p 2 × X p 3 to the unique solution sp of : P(Ω) p : Find s p ∈ (0, w 0 ) + V p such that Ω Q p (x) k(s p ) 0 p • k(ξ ) 0 p dx = L(ξ ), ∀ξ ∈ V p
where k(ξ ) stands for (e(v), ∇ψ), while k(ξ ) 0 p denotes the projection of k(ξ ) on the space H 0 p , for all ξ = (v, ψ) ∈ V p .

Proof. Actually, the problem can be considered as a "vectorial variant" of the asymptotic behavior of thin piezoelectric plates [START_REF] Weller | Asymptotic modeling of thin piezoelectric plates[END_REF] where each component of the phason plays the role of an electrical potential. Hence all the arguments of [START_REF] Weller | Asymptotic modeling of thin piezoelectric plates[END_REF] can be used; the additional case p i = 0 being obviously handled through the assumptions (H

P Ω ε p : Find s ε p ∈ 0, w ε 0 + V ε p such that Ω ε Q ε p (x) k s ε p 0 p • k(ξ ) 0 p dx = L ε (ξ ), ∀ξ ∈ V ε p
The superscript ε in V ε p and L ε means that Ω is replaced by Ω ε in the definition of V p and L, respectively. This transported problem is our proposal to model thin linearly quasicrystalline plates of thickness 2ε. In these models, constitutive laws are given by Q ε

p such that Q ε p (π ε x) = Q p (x)
, for all x in Ω. We enlighten some remarkable properties of these constitutive laws in the next section.

Some properties of Q ε p for icosahedral quasicrystals

The convergence result stated in Theorem 3.1 is of course valid for any quasicrystalline class. However, in all about 200 individual quasicrystals observed to date, a wide majority presents an icosahedral order. This suggests that this symmetry plays a central role in the study of the physical behavior of quasicrystalline media. For such symmetry, there are two elastic coefficients λ ε and μ ε in tensor C ε , one coupling coefficient r ε in tensor R ε and two coefficients K ε 1 and K ε 2 in tensor K ε (see relation [START_REF] Blaaderen | Quasicrystals from nanocrystals[END_REF]). We now focus our attention on the limit operators Q ε p which supply the constitutive relations of the quasicrystalline plate. These limit operators can be computed directly from the explicit expression [START_REF] Ding | Generalized elasticity theory of quasicrystals[END_REF]. In spite of the dimension reduction process, it is possible to show that they keep the same structure than Q ε (see [START_REF]The Physics of Quasicrystals[END_REF]) in the sense that

Q ε p = C ε p R ε p R ε T p K ε p where C ε p , R ε
p and K ε p respectively denote the limit tensors of phonon, phonon-phason coupling and phason coefficients.

These coefficients involve a mixture of the five coefficients λ ε , μ ε , r, K ε 1 and K ε 2 . More precisely, it is possible to show that we have the following properties:

-C ε p always involves a mixture of the later five physical coefficients, except for p = (2, 2, 2) for which K ε 2 is not present, -for p ∈ {(1, 1, 2), (1, 2, 2), (2, 1, 2), (2, 2, 2), (2, 0, 2), (0, 1, 2), (0, 2, 2)}, K ε p involves all the five physical coefficients, -for p ∈ {(1, 1, 1), (1, 1, 0), (1, 2, 1), (1, 2, 0), (1, 0, 1), (1, 0, 2), (1, 0, 0), (2, 1, 1), (2, 1, 0), (2, 2, 1), (2, 0, 1), (0, 1, 1), (0, 1, 0), (0, 2, 1), (0, 0, 1), (0, 0, 2)}, K ε

p involves μ ε , r ε , K ε 1 and K ε 2 but not λ ε , -for p ∈ {(2, 2, 0), (2, 0, 0), (0, 2, 0)}, K ε p involves only K ε 1 and K ε 2 ,
-for p ∈ {(1, 1, 2), (1, 2, 2), (1, 2, 0), (1, 0, 2), (2, 1, 2), (2, 2, 2), (2, 0, 2), (0, 1, 2), (0, 2, 2)}, R ε p involves all the five physical coefficients, -for p ∈ {(1, 2, 1), (1, 0, 1), (1, 0, 0), (2, 1, 1), (2, 1, 0), (0, 1, 1), (0, 1, 0)}, R ε p involves μ ε , r ε , K ε 1 and K ε 2 but not λ ε , -for p ∈ {(1, 1, 1), (1, 1, 0), (2, 2, 1), (2, 2, 0), (2, 0, 1), (2, 0, 0), (0, 2, 1), (0, 2, 0)}, R ε p involves r ε , K ε 1 and K ε 2 but not λ ε neither μ ε , -for p ∈ {(0, 0, 1), (0, 0, 2)}, R ε p vanishes: for quasicrystalline plates submitted to such type and magnitude of boundary conditions, the physical coupling between phonons and phasons disappears at the structural level! Moreover, the tensor K ε p is diagonal.

As an example, for p = (2, 2, 2) we have

C ε (2,2,2) = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ 2(2μ ε K ε 1 -r ε 2 )(λ ε +2μ ε ) (λ ε +2μ ε )K ε 1 2μ ε λ ε K ε 1 +2r ε 2 (λ ε +2μ ε ) (λ ε +2μ ε )K ε 1 0 0 2μ ε λ ε K ε 1 +2r ε 2 (λ ε +2μ ε ) (λ ε +2μ ε )K ε 1 2(2μ ε K ε 1 -r ε 2 )(λ ε +2μ ε ) (λ ε +2μ ε )K ε 1 0 0 0 0 2 μ ε K ε 1 -r ε 2 K ε 1 -2 r ε 2 K ε 1 0 0 -2 r ε 2 K ε 1 2 μ ε K ε 1 -r ε 2 K ε 1 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ R ε (2,2,2) = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ (K ε 1 -2K ε 2 )r ε K ε 1 0 r ε (3λ ε +2μ ε ) λ ε +2μ ε - (K ε 1 -2K ε 2 )r ε K ε 1 0 r ε (3λ ε +2μ ε ) λ ε +2μ ε 0 - (K ε 1 -2K ε 2 )r ε K ε 1 0 0 - (K ε 1 -2K ε 2 )r ε K ε 1 0 ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ K ε (2,2,2) = ⎛ ⎜ ⎜ ⎜ ⎜ ⎝ K ε 2 1 -K ε 2 K ε 1 -2K ε 2 2 K ε 1 0 0 0 K ε 2 1 -K ε 2 K ε 1 -2K ε 2 2 K 1 0 0 0 (λ ε +2μ ε )(K ε 1 +K ε 2 )-4r ε 2 λ ε +2μ ε ⎞ ⎟ ⎟ ⎟ ⎟ ⎠

3

 3 four suitable partitions of ∂Ω ε with Γ ε uD and Γ ε w D,i of strictly positive surface measures. The plate is clamped along Γ ε uD and the phason field w ε satisfies w ε

  1 )-(H 5 ) which, together with Poincaré inequality, imply the boundedness of |(∇ p (ε, w p(ε))) ij | L 2 (Ω) , 1 j 3.The case p = (0, 0, 0) is not considered here because it leads to a Kirchhoff-Love plate model without contribution of phason fields. 2To get physically meaningful results it is necessary to define a physical state s ε (π ε x) = s p (x), for all x in Ω. This physical state is the unique solution of a problem posed over the real plate,

p over the genuine plate Ω ε through the descaling s ε p namely:


Ruijing Hu
email: ruijing.hu@lip6.fr

Efficient Probabilistic Information Broadcast Algorithm over Random Geometric Topologies

This paper studies reliability of probabilistic gossip algorithms over the random geometric topologies which model ad hoc networks. We propose an efficient algorithm that ensures higher reliability at lower message complexity than the three families of gossip algorithms. Such an improvement is reasonably estimated by our reliability model. The results obtained by OMNET++ simulator confirm the prediction that our algorithm is the best choice for random geometric networks.

I. INTRODUCTION

Information broadcast is a basic service for many distributed applications in ad hoc networks. A message generated by a source is transferred to all other sites by successive retransmissions, since sites have no direct connection to each other in the network. In particular, even though global information about the topology is hardly available in large-scale system, an efficient broadcast protocol should still be able to provide high reliability and low latency with tolerable message exchanges.

The simplest solution is using a pure flooding protocol [START_REF] Jetcheva | A simple protocol for multicast and broadcast in mobile ad hoc networks[END_REF]. Upon the first reception of a message, every site delivers it to the application and also forwards it to all sites in its retransmission range once. Evidently, this form of dissemination is very wasteful and may give arise to a broadcast storm [START_REF] Ni | The broadcast storm problem in a mobile ad hoc network[END_REF], entailing for instance packet collision and loss.

A common alternative to the flooding is to require the sites of the system to relay the message only to some of their neighbors. The decision of the neighbors is either carried out in a deterministic or probabilistic way. The problem of the former ensuring high reliability is that its implementation is commonly very hard (i.e., in some schemes it is proven to be NP-hard to reach the optimal performance) and requires more information even if the mobility is not considered [START_REF] Lipman | Neighbor aware adaptive power flooding in mobile ad hoc networks[END_REF]. Hence, in our study we are interested in probabilistic schemes that highlight their simplicity and scalability [START_REF] Friedman | On reliable dissemination in wireless ad-hoc networks[END_REF], [START_REF] Ustebay | Efficient decentralized approximation via selective gossip[END_REF]. More precisely, three basic probabilistic gossip algorithms are widely-deployed in ad hoc networks without overlay construction for reliable dissemination: (1) Probabilistic Edge Gossip (GossipPE) [START_REF] Shen | Directional broadcast for mobile ad hoc networks with percolation theory[END_REF], (2) Probabilistic Broadcast Gossip (GossipPB) [START_REF] Haas | Gossip-based ad hoc routing[END_REF], and (3) Probabilistic Inverse Self-degree Broadcast Gossip (GossipPISB) [START_REF] Drabkin | Rapid: Reliable probabilistic dissemination in wireless ad-hoc networks[END_REF]. Furthermore, we propose an efficient algorithm: Probabilistic Inverse Neighbor-degree Edge Gossip (GossipPINE), which exploits local information about sites in a one-hop neighborhood. It performs better than the basic gossip algorithms over two kinds of random geo-metric topologies with and without border effects, which are commonly used to analyze the behavior of ad hoc networks.

On the other hand, we introduce a model to explain the best reliability given by our algorithm GossipPINE amongst all gossip algorithms. Interestingly, we find that the reliability of GossipPINE is accurately estimated for both topologies, and not affected by complex feature of the topology such as edge dependency and border effect. In contrary, the reliability of the other algorithms shows a significant variation from one topology to another. We therefore argue that GossipPINE is robust even in networks whose global characteristics are unknown beforehand.

Extensive simulations on top of OMNET++ [START_REF]Omnet++: Discrete event simulation system[END_REF] almost conform to our theoretical results of reliability. The latency, another important metric for broadcasting service, is also evaluated.

The road-map of this paper is organized as follows. An overview of the random networks and the gossip algorithms is given in Sections II and III. Section IV presents our gossip algorithm. The performance metrics are defined in Section V, while we explain our model in Section VI. Section VII evaluates our simulation and related work is discussed in Section VIII. Section IX concludes this paper.

II. SYSTEM TOPOLOGIES

A network underlying a large-scale dissemination system can be viewed as a bidirectional or undirected graph. It is comprised of N sites {s 1 , s 2 , • • • , s N }. The set of sites s j , connected to s i (i.e., s i ∼ s j ), is called s i 's neighbors, represented by Λ i whose cardinal is the degree of s i , denoted V i . P (k) represents the degree distribution of sites with k neighbors (i.e., the fraction of sites with degree k) in the graph and V is the mean degree V = N -1 k=0 P (k) • k . Moreover, another important property of the graph is defined as follows.

Clustering Coefficient, denoted C of a given random graph, for distinct sites s i ,s j ,s k , is the conditional probability that, given the existence of edges s i ∼ s k and s j ∼ s k , an edge s i ∼ s j also exists (i.e., P connect (s i ∼ s j |s i ∼ s k , s j ∼ s k )).

Our study focuses on the random geometric graphs [START_REF] Penrose | Random Geometric Graphs[END_REF], which model ad hoc networks in [START_REF] Drabkin | Rapid: Reliable probabilistic dissemination in wireless ad-hoc networks[END_REF], [START_REF] Haas | Gossip-based ad hoc routing[END_REF], [START_REF] Shen | Directional broadcast for mobile ad hoc networks with percolation theory[END_REF].

More precisely, the sites are distributed uniformly at random in a bounded region, which composes a random geometric graph. In this article, such a region is a rectangular plane with length a and width b. Furthermore, two sites are connected, as long as the distance between them is at most ρ. Based on [START_REF] Philips | Connectivity properties of a packet radio network model[END_REF], we can fine-tune ρ >

(1+ε)•ln(N )•a•b N •π
with a positive constant ε in order to ensure that the graph is connected.

As observed in [START_REF] Hu | Fair comparison of gossip algorithms over large-scale random topologies[END_REF], this random graph can have Border Effect, where the degree of the sites close to boundaries is much smaller than the mean degree. If we ignore it, the topology, denoted T (N, ρ), follows a Poisson-distribution [START_REF] Kong | Characterization of the critical density for percolation in random geometric graphs[END_REF]:

P t (k) = exp -Vt Vt k k! with mean degree Vt = N •π•ρ 2 a•b .
In [START_REF] Avin | Distance graphs: From random geometric graphs to bernoulli graphs and between[END_REF], the edge dependency is calculated as C = 0.5865.

Otherwise, the sites with low degree near the boundary, which in fact have an important impact on reliability of gossip algorithm, should be taken into account. As shown in [START_REF] Hu | Epidemic dissemination algorithms in large-scale networks: comparison and adaption to topologies[END_REF], the topology with the border effect, denoted G (N, ρ), has the degree distribution:

P g (k) = P t (k) (a-2ρ)(b-2ρ) ab + 2aρ+(b-2ρ)ρ ab ψ(k) where ψ(k) = 1 0 exp -Vt (F (x) -1) • F (x) k dx with F (x) := 1 π x √ 1 -x 2 -arccos(x) +1.
Trivially, the mean degree can be calculated as Vg =

N -1 k=0 P g (k) • k.

III. GOSSIP ALGORITHMS

Broadcast in large-scale network is commonly studied on the basis of Algorithm 1. Initially, the source sends a message to all of its neighbors (Lines 2 and 3). A site delivers and retransmits a received message provided it has not previously received it; otherwise, the message is discarded. Sites that have received the message at least once are called infected sites. Basically, probabilistic gossip algorithms, which disseminate information in ad hoc networks without overlay construction can be classified into three families, to implement Gossip() procedure: (1) Probabilistic Edge Gossip (GossipPE) [START_REF] Shen | Directional broadcast for mobile ad hoc networks with percolation theory[END_REF], (2) Probabilistic Broadcast Gossip (GossipPB) [START_REF] Haas | Gossip-based ad hoc routing[END_REF], and (3) Probabilistic Inverse Self-degree Broadcast Gossip (GossipPISB) [START_REF] Drabkin | Rapid: Reliable probabilistic dissemination in wireless ad-hoc networks[END_REF]. Besides the received message, the Gossip() procedure takes one parameter which is identical for all sites. In the sequel, Random() generates a random number in the interval [0, 1]. Send( msg , s j )

13

In GossipPE (Algorithm 2), a site s i chooses to send msg over an edge independently from the other edges with regard to a fixed probability p e (see Line 12). Note that when p e = 1 for all sites, we obtain the flooding algorithm.

Algorithm 3: Probabilistic Broadcast Gossip (at s i ) / * p v : probability to broadcast * / 14 GossipPB ( msg ,p v ) 15 if Random() p v then 16 foreach s j ∈ Λ i do 17
Send( msg , s j )

18
Unlike Algorithm 2, in GossipPB (Algorithm 3), each site, except the source, diffuses msg to all its neighbors with probability p v (see line 16). In particular, when p v = 1 this protocol becomes the flooding algorithm.

Algorithm 4: Probabilistic Inverse Self-degree Broadcast Gossip (at s i ) / * c v : constant to control broadcast * / 19 GossipPISB ( msg ,c v ) 20 if Random() min{ cv Vi , 1} then 21 foreach s j ∈ Λ i do 22
Send( msg , s j )

23

In GossipPISB (Algorithm 4), every site except the source rebroadcasts msg to all its neighbors independently with probability min{ cv Vi , 1} (see Line 21). When

c v max {V 1 , V 2 , • • • , V N } this protocol is the flooding algorithm.

IV. OUR ALGORITHM

The idea behind our algorithm is to ensure that sites with very low degree are infected, if one of their neighbors is infected. On the other hand, fewer redundant message copies are received at the sites with high degree. This will be analyzed in Sections VI and VII.

Algorithm Algorithm 5 shows our Probabilistic Inverse Neighbordegree Edge Gossip, denoted GossipPINE. We assume that every site knows the degree of its neighbors. Like GossipPE, site s i randomly chooses its edges over which msg should be transmitted. However, the probability to send on one edge depends on the degree of the connected neighbor (see Line 27).

V. PERFORMANCE METRICS

For performance evaluation of the gossip algorithms, the following metrics are commonly used [START_REF] Eugster | δ-reliability: A probabilistic measure of broadcast reliability[END_REF], [START_REF] Leitao | Epidemic broadcast trees[END_REF], [START_REF] Luo | Route driven gossip: Probabilistic reliable multicast in ad-hoc networks[END_REF]:

Message Complexity, denoted M: measures the mean number of messages received (or sent, since no message loss is taken into consideration) by each site:

M = Ω N -1 ,
where Ω is the total number of messages exchanged during the dissemination.

Reliability, denoted R: is defined as the percentage of messages generated by a source that are delivered by all sites. A reliability value of 100% is indicative that the algorithm was successful in delivering any given message to all sites (i.e., every site is infected for any given message) ensuring thus atomicity similarly to pure flooding algorithms [START_REF] Kermarrec | Probabilistic reliable dissemination in large-scale systems[END_REF].

Latency, denoted L: measures the number of hops required to deliver a message to all recipients, i.e., the number of hops of the longest path among all the shortest paths from the source to all other sites that received the message.

VI. RELIABILITY MODELING

In this section, we compare the probabilistic gossip algorithms by our reliability model, which explains the best performance of our algorithm GossipPINE in random geometric graphs.

Since we are only interested in whether a site has been infected by the end of the broadcast, a dissemination graph is exploited for modeling the reliability.

Given a graph G = (S; E), and a source site s 0 , the probabilistic gossip algorithms follow the same principle: when a new message arrives at a site, this site chooses a subset (possibly empty) of its neighbors in G and forwards the message to these neighbors.

In our analysis, we assume that every site selects a set of neighbors to receive the message in an initial phase. We thus obtain a directed graph

- → G = (S; - → E ) with sites identical to G. For an edge (s i ∼ s j ) ∈ E, the arc --→ s i s j is a part of - → G ,
if s i has selected s j in the initial phase. Then, a site s i is infected during dissemination, if there is a path from s 0 to s i in -→ G . Furthermore, if every site receives the message, the number of arcs in -→ G expresses the message complexity. We call -→ G Dissemination Graph.

We call a site isolated in the dissemination graph, if it has no incoming arcs. Evidently, an isolated site is not infected by a given gossip algorithm. We call the probability that an edge in E becomes an arc towards one site in S the forwarding probability, denoted p forw .

We assume as a precondition that the random network G = (S; E) is connected. Thus, any message will reach all sites by pure flooding. We model the reliability as the probability that no site is isolated in the dissemination graph. Hence, for a random topology with N sites and degree distribution P (k) we model the probability R as

R = N -1 i=1 1 -(1 -p forw (i)) i P (i)•N (1)
where p forw (i) denotes the probability that a site with degree i has an incoming arc from a neighbor. Our evaluation shows that this equation provides a good estimate for the reliability while assuming that every site is isolated independently as done in [START_REF] Friedman | On reliable dissemination in wireless ad-hoc networks[END_REF]. However, due to high clustering coefficient in random geometric graphs, one isolated site can even make adjacent sites never receive the message, particularly in sparse density area. Our model shown in Equation (1) will consequently overestimate the reliability, especially in G (N, ρ) with border effect. Yet, our simulations show that our model still gives a valuable estimate and correctly indicates which algorithm is more or less reliable. GossipPE enables every site to be selected at random by neighbors with probability p e . Therefore, its forwarding probability: p P E forw = p e . In the same way, p P B forw = p v and p P IN E forw (i) = min{ ce i , 1} for GossipPB and GossipPINE respectively. The message reception by GossipPISB in a site s i depends on the degree of its neighbors, while the latter is highly associated with the degree of s i . Lemma 1 studies the average degree of an arbitrary neighbor of s i whose degree is given, ignoring the border effect of random geometric graphs.

Lemma 1: In T (N, ρ), given a site s i with degree V i and one adjacent site s j with V j neighbors, we can get the average value of V j as:

N eighD = 1 + (V i -1) • C + ( Vt -1) • (1 -C) (2)
where

Vt = N •π•ρ 2 a•b
is the mean degree of the graph. Proof: Figure 1 shows s i , a neighbor s j , and the retransmission zones (i.e., Disc(s i ), Disc(s j )). We call the overlapped part between two discs lune (i.e., Lune(s i , s j )). Let disc be the area of Disc(s j ). Then, the average area of the lune is C • disc. Clearly, the density of sites inside Disc(s i ), and also inside Lune(s i , s j ) is Vi disc . If disc is small, we can assume that the density outside Disc(s i ) is Vt disc . Therefore, From Lemma 1, we can derive p P ISB forw (i) = min{ cv 1+(i-1)•C+( Vt-1)•(1-C) , 1} for T (N, ρ), while we approximate it for G (N, ρ) in replacing Vt by Vg that is the mean degree of the random topology with the border effect (see Section II).

s j has Vi disc • C • disc neighbors in Lune(s i , s j ), V disc • (1 -C) • disc neighbors outside Lune(s i , s j ), plus s i .
Theorem 2: If all algorithms generate the same message complexity, then for k < V

p P IN E forw (k) > p P ISB forw (k) > p P E forw = p P B forw . (3)
holds.

Proof: According to the study in [START_REF] Hu | Fair comparison of gossip algorithms over large-scale random topologies[END_REF], we infer input parameters of gossip algorithms to reach a given message complexity. Therefore, when identical message complexity holds for all gossip algorithms, c e = c v , p e = p v , and forw . Furthermore, we can deduce that the sites with degree smaller than V are the majority in T (N, ρ) and G (N, ρ), which dominate the site isolations. In regard to our reliability model described in Equation ( 1), GossipPINE, which ensures the message reception in sites with lower degree, is thus more reliable in random geometric graphs than the other gossip algorithms.

c v p e • V . Trivially, p P E forw = p P B forw . Since k < V implies k < N eighD < V ,

VII. PERFORMANCE EVALUATION

In this section, we evaluate performance of the four probabilistic gossip algorithms in Sections III and IV: GossipPE, GossipPB, GossipPISB, and our algorithm GossipPINE. These algorithms are executed over two kinds of random geometric graphs described in Section II.

Our simulations are conducted on top of OMNET++. We consider networks composed of N = 1000 sites and, in order to ensure connectivity, ε = 1 for every random geometric graph, where a = 7500, b = 3000, ρ = 330, such that the mean degree of G (N, ρ) (resp., T (N, ρ)) is about 14 (resp., 15). The graphs without border effect are generated by bending the rectangular region into a Torus T (N, ρ), while the Euclidean distance between every site is measured.

The results, obtained by our reliability modeling, are also compared with our simulation. We can observe that in terms of reliability, our algorithm outperforms the other gossip algorithms, which conforms to Theorem 2.

For each gossip algorithm, we fixed the parameter values to reach a given message complexity in using the method introduced in [START_REF] Hu | Fair comparison of gossip algorithms over large-scale random topologies[END_REF], and then evaluated the reliability and the latency of the algorithm. 200 different messages are generated by 200 different sources that are chosen uniformly amongst 1000 sites over 50 different graphs related to each of the topologies. Then, the results for each message complexity are averaged by the 200 × 50 = 10000 message disseminations.

As a matter of fact, we can configure the topology parameters with other values, whereas the best reliability of our algorithm GossipPINE is always obtained, as estimated by our model for the geometric graphs. Due to lack of space in our paper, we do not give exhaustive presentation.

Reliability: We now discuss the reliability values for the gossip algorithms in T (N, ρ) and G (N, ρ), obtained by the simulations as well as our model.

The results are plotted in Figure 2, where the models for GossipPINE and GossipPISB are called ModelPINE and ModelPISB respectively. As explained in Section VI, our model gives the same result for GossipPE and GossipPB. This is depicted as ModelPE.

GossipPINE exhibits the best performance in both Figures 2(a) and 2(b), while GossipPB and GossipPE show the worst reliability. In regard to T (N, ρ), the reliability is more or less degraded in G (N, ρ) for every gossip algorithm. Especially, GossipPB is far less reliable than the other algorithms in G (N, ρ). Compared with T (N, ρ), on the boundaries of the rectangular field G (N, ρ) lie a large number of sites with low degree (see Section II). These sites are easily isolated in GossipPB, which results from the clustering effect explained in [START_REF] Hu | Fair comparison of gossip algorithms over large-scale random topologies[END_REF].

Due to the high clustering coefficient, as shown in Lemma 1, sites with low degree are very likely to be clustered together in random geometric graphs (i.e., small groups of interconnected sites with relatively few edges to other groups). As a result, a single isolated site can render many adjacent sites isolated. This effect is further reinforced by the border of the rectangular region comprising sites with low degree. GossipPISB leverage this clustering, so that the forwarding probability p P ISB forw (i) for small i is larger than p P E forw . It thus mitigates isolating sites in a cluster with low degree and achieve higher reliability than GossipPE.

Since the site isolation dependency is not included in our model, the results from our modeling overestimate the performance obtained by simulation in Figures 2(a) and 2(b). We see that such dependency in G (N, ρ) has a greater impact on GossipPB and GossipPISB, for which the reliability falls behind the predicted model, in regard with GossipPINE and GossipPE. The reason is that the message retransmitted by GossipPB or GossipPISB stays within a cluster few hops away from the source, when the retransmission probability is not high enough. In accordance with the percolation theory [START_REF] Grimmett | Percolation[END_REF], GossipPE is never worse than GossipPB on any topology.

It is worth pointing out that in GossipPINE the probability for a message to be spread out of a cluster with high degree sites and into a cluster with low degree sites is larger than in GossipPISB and GossipPE (see Theorem 2). As expected, GossipPINE is accurately described by our model over both T (N, ρ) and G (N, ρ).

Latency: An efficient broadcasting algorithm aims at providing high reliability, while minimizing both message complexity and latency. We analyze the latency performance when the reliability reaches at least 85%. Figures 3(a Over all random topologies, after a given message complexity, latency does not decrease anymore, but converges towards the pure flooding approach (i.e., the shortest routes between the source and the other sites), and therefore, towards the minimum latency. The latency of our algorithm GossipPINE reaches a lower value in T (N, ρ) with the smallest message complexity amongst all gossip algorithms. A similar behavior of the latency curves is observed in G (N, ρ), except that the minimum latency value lies around 22 hops which is larger than in G (N, ρ), since the diameter of G (N, ρ) is greater than that of T (N, ρ).

VIII. RELATED WORK

The initial study on the probabilistic information broadcast in ad hoc networks was conducted in [START_REF] Haas | Gossip-based ad hoc routing[END_REF] and later implemented on DES-Testbed [START_REF] Blywis | Gossip routing in wireless mesh networks[END_REF], which evaluates GossipPB and its variants. A bimodal behavior of the reliability has been heuristically analyzed with regard to the retransmission probability p v . More precisely, from a threshold value of p v , all the sites can receive almost all the messages. Otherwise, the message can rarely reach all sites in the system. In [START_REF] Sasson | Probabilistic broadcast for flooding in wireless mobile ad hoc networks[END_REF], such a threshold is given on the basis of percolation theory. The percolation driven flood routing algorithm proposed in [START_REF] Vakulya | Energy efficient percolation-driven flood routing for large-scale sensor network[END_REF] can be seen as an advanced version of GossipPB in largescale ad hoc networks. It gives an answer how to choose p v to ensure high reliability. In Section VI, we model the reliability as a function of the p v , which at the same time takes into account the properties of underlying networks, such as degree distribution. While assuming the site isolation independency in the end of the dissemination, our model accurately establishes such a relation in infinite random geometric graphs. Even in the graphs with border effect, reasonable estimation can be adaptively given by the corrective degree distribution. Furthermore, it explains the higher reliability obtained by our algorithm than the other gossip algorithms.

GossipPISB proposed in [START_REF] Drabkin | Rapid: Reliable probabilistic dissemination in wireless ad-hoc networks[END_REF] tailors the retransmission probability to the local topology, where sites in sparse density area retransmit with high probability. We have shown that the sites around boundaries of the region are hardly infected when small message complexity is allowed. In [START_REF] Friedman | On reliable dissemination in wireless ad-hoc networks[END_REF], the authors model the reliability, and improves the protocol by adding a pull phase. But pull algorithms are difficult to carry out [START_REF] Kermarrec | Probabilistic reliable dissemination in large-scale systems[END_REF]. Similarly, the authors in [START_REF] Abdulai | Route discovery in wireless mobile ad hoc networks with adjusted probabilistic flooding. adaptive science technology[END_REF] distinguish all sites by four levels. The probability in each level is reversely proportional to the number of levels predefined in ad hoc network. Nonetheless, no theoretical analysis are proposed. The Smart Gossip protocol in [START_REF] Kyasanur | Smart gossip: An adaptive gossip-based broadcasting service for sensor networks[END_REF] constructs local relationship trees in two-hop neighborhoods, which consist of a parent, children and siblings. Every parent decides on an optimal forwarding probability that is computed by its children under reliability and topology constraint. Maintenance of the relationship trees is required. Moreover, in [START_REF] Cartigny | Border node retransmission based probalistic broadcast protocols in ad-hoc networks[END_REF], the authors specify forwarding probability of a site in ad-hoc networks as function of both its degree and the mean degree of its neighbors, while introducing a mechanism where if a site suspects that some of its neighbors have not received the message, it rebroadcasts the message regardless of the initial decision. The reliability is heuristically studied. In contrary, our algorithm is efficient and very simple to apply in random geometric networks, which is confirmed by both our model and simulation.

As explained in [START_REF] Hu | Fair comparison of gossip algorithms over large-scale random topologies[END_REF] that GossipPE is never worse than GossipPB in any topology, the directional antenna is applied to implement GossipPE in ad hoc network in [START_REF] Shen | Directional broadcast for mobile ad hoc networks with percolation theory[END_REF], which shows a better performance than GossipPB. In comparison, our algorithm GossipPINE is robust against border effect that is a common characteristic in the realistic environment with long boundary.

Other probabilistic approaches include counter-based [START_REF] Ni | The broadcast storm problem in a mobile ad hoc network[END_REF], distance-based [START_REF] Huang | Efficient lightweight broadcasting protocols for multi-hop ad-hoc networks[END_REF], and location-based mechanisms [START_REF] Tseng | Adaptive approaches to relieving broadcast storms in a wireless multihop mobile ad hoc network[END_REF]. In regard to GossipPB, variant protocols arm a time counter as shown in [START_REF] Ni | The broadcast storm problem in a mobile ad hoc network[END_REF] for copy accumulations. Hop Counter Aided Broadcast (HCAB) protocol in [START_REF] Huang | Efficient lightweight broadcasting protocols for multi-hop ad-hoc networks[END_REF] makes the sites, upon the first reception of a message, start a random timer and record the value of hop counter. The message will be relayed by the site if, when the timer expires, no message with a hop counter higher than the first was received. Implicitly, every site in HCAB attempts to know whether each relay is done by some other sites, hopefully covering (i.e., infecting) additional regions of the network. Another distance-based protocol is also proposed in [START_REF] Huang | Efficient lightweight broadcasting protocols for multi-hop ad-hoc networks[END_REF], which is called Self-Adaptive Probability Broadcasting (SAPB) protocol. It stipulates that the probability of a site to relay a message depends on coverage (i.e., infection) contribution brought by this site. However, those works contain little theoretical evaluation of the schemes.

IX. CONCLUSION

In this article, we have proposed an efficient probabilistic broadcast algorithm and reliability model to compare probabilistic gossip algorithms in random geometric networks. Our model gives a reasonable estimation, and explains the best reliability obtained by our gossip algorithm GossipPINE, which is confirmed by our simulations. GossipPINE also shows desirable latency performance for high reliability. In particular, our algorithm has the best resilience to the border effect degrading the reliability.

As a near future work, the random mobility in the network will be evaluated and studied.

Figure 1 :

 1 Figure 1: Two kinds of s j 's neighbors

  ) and 3(b) show the simulation results in T (N, ρ) and G (N, ρ) respectively.

Figure 2 :Figure 3 :

 23 Figure 2: Reliability comparison of gossip algorithms


Yahia Salhi
email: yahia.salhi@univ-lyon1.fr

Pierre-E Thérond

Julien Tomas
email: julien.tomas@univ-lyon1.fr

A Credibility Approach of the Makeham Mortality Law *

Keywords: Credibility, Makeham law, Mortality, Life insurance, Graduation, Extrapolation

) and Hardy and Panjer (1998), we propose a credibility approach which consists on reviewing, as new observations arrive, the assumption on the mortality curve. Unlike the methodology considered in Hardy and Panjer (1998) that consists on updating the aggregate deaths we have chosen to add an age structure on these deaths. Formally, we use a Makeham graduation model. Such an adjustment allows to add a structure in the mortality pattern which is useful when portfolios are of limited size so as to ensure a good representation over the entire age bands considered. We investigate the divergences in the mortality forecasts generated by the classical credibility approaches of mortality including Hardy and Panjer (1998) and the Poisson-Gamma model on portfolios originating from various French insurance companies.

Introduction

Recently, interest from life insurers to assess their experienced mortality risk has considerably increased. The new regulation and norms, Solvency II, shed light on the need of life tables that best reect the experience of insured portfolios so as to reliably quantify the underlying mortality risk. Insurers, in France for example, are used to rely on regulatory life tables for pricing purposes, which are sometimes too conservative. In general, ill-suited mortality assumptions and life tables, especially being too conservative, lead to two eects:

(i) Increase of Best Estimate technical provisions (and thus decrease Basic Own-Funds );

(ii) Increase of the base gure used for calculating the capital charge for mortality risk (15 % increase scenario of the conditional death rates under the Solvency II framework).

Therefore, the question of which mortality table can be considered for pricing and reserving purposes is of substantial importance. A rst attempt, to handle this issue, is to use the available data at the portfolio level and build a specic mortality table. However, practitioners may face technical diculties related to the size of the portfolio and the heterogeneity of the guarantees (for the same underlying risk). For instance, an insurer may detain a fairly big portfolio but with insureds holding dierent policies: pure endowment contracts, unit-linked contracts with minimum death guarantees, loan insurance and so on. In such a case, it is dicult to build mortality tables only based on the sole experience of each policy. Especially since it may induce signicant impacts on the technical reserves if the table has to be updated more frequently over time. In this paper, we consider an insurer with exposures to dierent policies and aiming at establishing an experience-based mortality table for each policy.

In the academic literature, various methodologies have been proposed to built and graduate mortality rates at the insured portfolio level. They are usually divided into non-parametric and parametric techniques. The latter are very useful in practice especially when there is sucient data, see [START_REF] Forfar | On graduation by mathematical formula[END_REF] for a comprehensive introduction to the use of parametric models for graduation. These approaches t the parametric structure to the mortality of interest over a given period. The graduated mortality is then used to project future liabilities related to the underlying population. By doing so, the evolution of the ow of data related to latest available information is not taken into account. This should be, for example, used to update the graduated mortality. However, if one decides to re-calibrate the parametric model each year, the forecasts are likely to be unstable. This is mainly due to the instability of parameters estimation due to the lack of sucient data.

In this context and following the work of [START_REF] Bühlmann | A course in credibility theory and its applications[END_REF] and [START_REF] Hardy | A credibility approach to mortality risk[END_REF], we propose a credibility approach which consists on reviewing, as new observations arrive, the parameters of a Makeham t. The framework considered in [START_REF] Hardy | A credibility approach to mortality risk[END_REF] focuses on the update of the aggregate deaths recorded over the whole portfolio. However, such an approach may be not eective in situations where the insurer liability is highly dependent on the age structure of the underlying portfolio. Thus, using an adjustment makes possible to add a structure in the mortality pattern which is useful when portfolios are of limited size so as to ensure a good representation over the entire age-band considered. Note that, adding an age structure is also benecial given the heterogeneity observed in the cost of the guarantees according to the age. To recap, as we can see in Section 5, the proposed adjustment approach is intended to enhance the predictive ability of the credibility-based revisions at the age-level and not on the aggregate portfolio level.

The remainder of the paper is organized as follows. Section 2 has still an introductory purpose. It species the notation, assumptions and the Makeham settings used in the following. Section 3 introduces the Makeham credibility approach and assess the estimation of the credibility model. Section 4 describes the classical credibility approaches of mortality including [START_REF] Hardy | A credibility approach to mortality risk[END_REF] and the Poisson-Gamma model. Section 5 presents an application with experience data originating from French insurance companies. Finally, some remarks in Section 6 conclude the paper.

2 A Credibility Model for Makeham's Law 2.1 Data Structure and Notation. We suppose that we have at our disposal age-specic mortality statistics originating from n portfolios. For each portfolio i ∈ {1, • • • , n}, we observe the deaths of exposures over a period T i . Denote the number of individuals at attained age x during calendar year t = 1, • • • , T i by L i x,t and D i

x,t represents the number of deaths recorded.

We also introduce the following notation,

D i x,• = T i t=1 D i x,t , L i x,• = T i t=1 L i x,t , and D i •,t = x x=x D i x,t , L i •,t = x x=x L i x,t ,
which refer respectively to the aggregate deaths and individuals over the age-band {x, x + 1, . . . , x} and calendar years 1 to T i for each portfolio i. Henceforth, the • indexation refers to the summation over the index of interest. For example, D •

x,• refers to the aggregate deaths over the period [1, T i ] and over the n portfolios, i.e. D

• x,• = n i=1 T i t=1 D i x,t .
2.2 Mortality Law. We consider the (rst) Makeham law of mortality, which generalizes the Gompertz law. Omitting the time dependency, [START_REF] Makeham | On the law of mortality[END_REF] assumes that the force of mortality ϕ i x at attained age x during calendar year t has the following form:

ϕ x = A + B × C x , (2.1)
with A, B and C are some constants. These parameters capture the essential properties of the progression of mortality. For instance, the dominant eect, i.e. the aging eect, over the age is captured by the multiplicative component factor B × C x . The non-age dependent parameter A can be interpreted as the non-senescent mortality, for instance, due to accidents. Both of these capture the exponential increase in the forces of mortality observed for adult mortality, see [START_REF] Bongaarts | Long-range trends in adult mortality: Models and projection methods[END_REF] for more details.

Various modication of the above law have been proposed, especially, to encounter for the time dependency of the mortality, see for e.g. [START_REF] Keytz | The limits of population forecasting[END_REF] among others. Indeed, as soon as age-specic mortality patterns over time are concerned, the time series records of the latter show a discernible downward trend with minor uctuations around. In order to correct this deciency in the model (2.1), we suppose that the time trend is incorporated in the parameter B i denoted henceforth B i t . Therefore, the force of mortality ϕ i x,t for portfolio i writes now as the following expression

ϕ i x,t = A i + B i t (C i ) x . (2.2)
This model should capture the behavior of the probability of death over years through the time-dependent parameter B i t . This also make possible the prediction of future mortality, i.e.

for t = T + 1, T + 2, . . ., through the study of the time series B i t for t = 1, . . . , T . When it comes to small portfolios, the model in (2.2) is not easy to implement. Indeed, as discussed later in this paper, the temporal behavior the factor B i cannot be accurately extracted. Nevertheless, one can use the estimated values of B i t even over the few periods to predict the future behavior of B i .

Note that in order to estimate the parameters of the model (2.2), given the growth of the forces of mortality with the age, we must have a constant C greater than 1 and a positive B. Then,

q i x,t = 1 -exp - x+1 x ϕ y,t dy = 1 -exp - x+1 x A i + B i t × (C i ) y dy = 1 -exp(-A i ) exp - B i t ln C i (C i ) x (C i -1) , (2.3)
where q i x,t denotes the one-year probability of death at attained age x during calendar year t for portfolio i. Consistent estimates A i , B i t and C i of the parameters are obtained by minimizing the following weighted distance:

x x=x L i x,t q i x,t (1 -q i x,t ) (q i x,t -q i x,t ) 2 , with q i x,t = D i x,t /L i
x,t is the crude mortality rates.

2.3 Dierential Mortality Law. It is common in modeling specic portfolio's mortality to consider an adjustment with regard to a baseline mortality. Generally, this implicitly assumes that both populations share common features up to a random eect. Relational models stipulate a deterministic relationship in the form q i x = f (q b x ) links the two mortalities, where q b

x refers to the baseline mortality. The function f : [0, 1] → [0, 1] is a known and deterministic function, see [START_REF] Delwarde | Modèles linéaires et additifs géneralisés, maximum de vraisemblance local et méthodes relationelles en assurance sur la vie[END_REF] for more details. A simple example would suggest that the death rate is common for all companies. Specically, q i x = q b x for any i ∈ {1, • • • , n}. However, such an assumption does not appreciate the specic characteristic of each portfolio's portfolio.

In other words, portfolios having lives in poorer or better conditions than the baseline mortality do not behave in a similar fashion than the baseline mortality. This implies that one should encounter for dierential mortality that arises due to portfolio specic features, e.g. particular socioeconomic groups involved, average income level, etc. However, when it comes to the study of the mortality at a single portfolio level, some specic issues arise:

(i) Size of populations: Insured population are generally of small size, so none or very few deaths are observable at some ages. This may not only bias the estimation of the force of mortality but also lead to a mis-estimation of the parameters in (2.3). This may cause high uctuations for q i x,t and consequently for A i , B i t and C i .

(ii) Length of historical data: Available age-specic mortality statistics lacks of deepness. This makes dicult to isolate a possible time trend as it may be captured by B i t . The latter may be uctuating due to the small size of the dataset as noted before.

(iii) Scale of available data: Insured portfolios show a typical behavior compared to a national mortality. The mortality of insured population is signicantly lower than the national popula-tion from which it is drawn. This could make the use of a baseline mortality based on national demographic statistics as a substitute useless as it may not have the same characteristics of the initial population.

All these characteristics make forecasting of future mortality evolution problematic. In order to overcome these issues when implementing and tting the model (2.3) for each portfolio i ∈ {1, . . . , n} we will make the following assumptions: (i) The baseline mortality q b

x,t is described by the Makeham model in (2.3).

(ii) The age eect is similar on the n portfolios and companies specic model is assumed to share the same parameters A i and C i . Those are set equal to the baseline ones, i.e. A i = A b and C i = C b for any i ∈ {1, • • • , n}.

(iii) The time-dependent parameter B i t is tted at each period. This is given by the following formula:

B i t = D i •,t -A b L i •,t x x=x (C b ) x L i x,t
.

The assumptions (i) and (ii) allows to overcome potential estimation bias of the parameters A i and C i . Indeed, basing the estimation on a large population allows to avoid erroneous inferences of the parameters. Also, if the portfolio i is a subset of the baseline population composed of the aggregated portfolios, we may think that both the non-senescent factor A i and the slope C i are equivalent and thus normalized with the baseline mortality. Empirical evidence of a normalized slope can be found in [START_REF] Thatcher | The long-term pattern of adult mortality and the highest attained age[END_REF]. It is shown that relative rate of increase is the same at all ages and is a shared feature with over subset populations, see also the empirical study of [START_REF] Zhu | Logistic regression for insured mortality experience studies[END_REF]. As for the non-senescent parameter, the assumption is relevant to the extent that this eect is generally of small impact and sometimes ignored (especially for industrialized countries), see [START_REF] Gavrilova | Ageing and longevity: Mortality laws and mortality forecasts for ageing populations[END_REF]. The unique parameter that captures the specic mortality at the portfolio level is B i t , which would a priori not be the same over companies due to the heterogeneity of the underlying populations as explained above. This can be regarded as an unobservable random factor and similar to the so-called frailty factor. Such a methodology is widely understood in the literature as well as in life insurance practice. Assumption (iii) gives an estimate of the time-dependent parameter. By time-dependent we only track the uctuation of B i t over time that might be caused by the small size and length of data. Thus our main aim is to sequentially adjust the estimation of B i t over time in view of the ow of information at our disposal.

3 Credibility of the Makeham Mortality 3.1 Next Period Prediction. In the following, we are interested in the behavior, over time, of the random variable

X i t = B i t B b t , (3.1)
and specically on its next period prediction X i T +1 merging information from other portfolios j = 1, • • • , n with j = i. Specically, suppose that we are at the end of the year T , i.e. at time T + 1, and we want predict the next period deaths D i

x,T +1 in the portfolio (equivalently the probability of death q i

x,T +1 ). Naturally, we can assume that this ratio is constant over time and thus invoke a widespread practice that applies a single factor of reduction/increase to the baseline mortality. On the other hand, one could propose a dynamic model on the same line as [START_REF] Plat | Stochastic portfolio specic mortality and the quantication of mortality basis risk[END_REF]. The latter proposes a modeling framework of the relative ratio of an experienced mortality (death rates) to a baseline and consider that this can be diused using either an autoregressive model or a decomposition similar to the one introduced by [START_REF] Lee | Modeling and forecasting us mortality[END_REF].

Other methodologies have been also proposed, see [START_REF] Ngai | Longevity risk management for life and variable annuities: The eectiveness of static hedging using longevity bonds and derivatives[END_REF] and Hyndman et al.

(2013) among others. However, random eects that constitute the decomposition of the experienced mortality have to be projected using their temporal and statistical features. In our case, we are not only interested in handling populations of small size but also with potentially limited historic period of observation. Therefore, such a methodology would typically not be useful in our setting as it requires a long experience.

Note that the behavior of X i t is broadly related to the so-called basis risk. This refers to the fact that the evolution of the policyholders mortality is usually dierent from that of the national population (baseline), due to some selection eects. This selection eect has dierent impacts on dierent insurance companies portfolios, as mortality improvements and accelerations are very heterogeneous in the insurance industry, see [START_REF] Barrieu | Understanding, modelling and managing longevity risk: key issues and main challenges[END_REF].

3.2 Heterogeneity and Makeham's Law Adjustment. As noted above, we are interested in the accurate adjustment of the portfolio-dependent parameter in (2.3), i.e. B i t . Given the specic parameterization of the problem, one may think of the n portfolios as a subset of the reference population and thus each population is characterized by a risk prole Θ i . In addition,

it is benecial to borrowing information across the dierent portfolios to enhance the knowledge and estimation of the mortality at the single portfolio level. Furthermore, these subpopulations may, for example, share a common mortality feature, while showing some specicity in their mortality prole. This can be seen as a random variable eect or heterogeneity characterizing the specic prole of each portfolio, for i = 1, • • • , n. Therefore, we implicitly assume that each portfolio is endowed by a risk prole θ i which is a realization of a random variable Θ.

In view of the various stylized facts presented above and in order to predict D i x,T +1 , for each age x, we focus on the projection of X i T +1 . Therefore, we suppose that this relative trend level of portfolio i with respect to the baseline mortality (trend) is characterized by the risk prole θ i which is a realization of Θ i . In other words, X i t is viewed as a function of a random element Θ i representing the unobserved characteristics of the portfolio mortality trend (with respect to the baseline). By doing so, we implicitly take into account the heterogeneity of the portfolio i's portfolio mortality prole. It thus remains to predict X i T +1 taking into account this random het- erogeneity. By doing so, we naturally invoke the use of a credibility approach to estimate X i T +1 .

3.3 Credibility Based Adjustment. As noted above, the objective is to estimate the next period projection of the relative ratio for each portfolio i. More precisely, in view the available data up to time T

i , i.e. X i t , for t = 1, • • • , T i , one aims to nd the best estimate of E[X i T i +1 |Θ i ] = µ(Θ i )
, which is unknown. Let µ(Θ i ) be this estimation. For this purpose and using the usual credibility setting, we shall make the following hypotheses:

(H1) Conditionally on Θ i , the random variables X i t , for t ∈ {1, . . . , T i }, are independent with mean and variance given as follows

E X i t | Θ i = µ Θ i and Var X i t | Θ i = σ 2 Θ i ω i t , ω i t = x x=x (C b ) x L i x,t n i=1 x x=x (C b ) x L i x,t
, measures the weight given to the period t experience from the portfolio i.

(H2) The pairs Θ i , X i t , Θ k , X k l , k = i are independent and identically distributed.

The rst assumption (H1) implies that for each risk prole i (portfolio), the true relative ratio µ Θ i (conditionally on the knowledge of the risk prole Θ i ) does not change over time and its variance given Θ i , Var X i t | Θ i changes in proportion to the relative size of the portfolio ω i t .

The latter expresses dierent concerns outlined earlier. Specically, it links the variability of the estimation of the parameter B i t to the size of the underlying population: very small portfolios are subject to larger variability on the estimation of B i t and vice versa.

The second assumption (H2) means that the risk proles are independent. The successive realizations of the relative ratio X i t for any portfolio are independent of each other except through the risk parameter Θ i . Moreover, using the random variable X i t instead of B i t permits to avoid data adjustment for Intuitively, assumption (H2) implicitly suggests that portfolios are comparable as they random sub-groups of a reference (national) population, but not entirely similar which induces the conditional independence.

In view of these assumptions, the following results are straightforward:

(i) The expected prediction of X i T +1 unconditionally on the risk prole Θ i is given by E X i T +1 = E X T +1 (Θ i ) = 1.
In other words, in the absence of any information on the heterogeneity level on the parameter B i t , the best next-period prediction of the latter is the reference one, i.e.

E B i t = B b t .
(ii) Using the law of total variance, the dependence structure of portfolio i associated risk factor over time, is, for l, t ∈ {1, . . . , T },

Cov X i l , X i t = Cov E X i l |Θ i , E X i t |Θ i + E Cov X i l , X i t |Θ i = Var µ Θ i + E Cov X i l , X i t |Θ i =    τ 2 if l = t τ 2 + σ 2 ω i t if l = t, (3.2)
where Var µ

Θ i = Var Θ i := τ 2 , while E σ 2 Θ i = E Θ i := σ 2 .
3.4 Credibility Estimator. Following the Bühlmann-Straub credibility approach, the aim is to nd the best estimate of the actual to expected mortality ratio E X i T i +1 | Θ i = µ Θ i which is linear in the observations. For each portfolio, due to the assumption (H2), µ(Θ i ) depends only on the observations and the linear credibility estimator is of the form

µ Θ i = a i 0 + T i t=1 a i t X i t , (3.3)
where the coecients a i t , for t = 0, • • • , T i , are those minimizing the mean squared errors crite-rion

a i t t=0,••• ,T i = argmin (a i t ) (t=0,••• ,T i ) E µ Θ i -a i 0 - T i t=1 a i t X i t 2
.

In view of Equation 3.2, taking the derivatives of the above criterion with respect to the a i,t 's and equating to zero gives,

a i 0 = 1 - τ 2 ω i • σ 2 + τ 2 ω i • and a i t = τ 2 ω i t σ 2 + τ 2 ω i • , with ω i • = T i t=1 ω i t .
(3.4)

Then, substituting Equation 3.4 into (3.3), leads to the following the Bühlmann-Straub credibility estimator of X i

T i +1 X i T i +1 (Θ i ) = α i X i • + (1 -α i ), with α i = ω i • τ 2 /(ω i • τ 2 + σ 2 ), (3.5)
where

X i • = ( T i t=1 ω i t X i t )/ω i • .
Note that the ratio σ 2 /τ 2 represents the credibility coecient. The parameter α i is called the credibility factor or credibility weight for portfolio i and takes values in [0, 1]. For each portfolio i, note that the larger the volume of historical data, the larger α i will be, see Equation 3.5.

3.5 Estimators of the Structure Parameters. As the risk parameters, Θ i , for i ∈ {1, . . . , n}, are assumed to be identically distributed, their moments are identical. Therefore τ 2 and σ 2 are the same for all portfolios and measure the residual heterogeneity of the risk proles and the pure randomness respectively. These parameters are the key determinants of the credibility estimator, i.e. Equation 3.5. In the following, special attention is addressed to the estimation of these quantities. Recall the denition of the structure parameters,

σ 2 = E σ 2 Θ i = ω i t E Var X i t |Θ i , and τ 2 = Var E[X i t |Θ i .
Then, it is reasonable to propose the estimators σ 2 and τ 2 in the same vein as [START_REF] Bühlmann | A course in credibility theory and its applications[END_REF] based on the observations X i t :

σ 2 = 1 n n i=1 s 2 i , with s 2 i = 1 T i -1 T i t=1 ω i t X i t -X i • 2 , (3.6) and τ 2 = ω • • (ω • • ) 2 -n i=1 (ω i • ) 2 n i=1 ω i • X i • -X • • 2 -(n -1) σ 2 , with X • • = 1 ω • • n i=1 ω i • X i • and ω • • = n i=1 ω i • .
These estimators are unbiased and consistent, see [START_REF] Bühlmann | A course in credibility theory and its applications[END_REF] for more details.

Note that τ 2 can be negative. This would mean that there would be no dierence between the risks. In this case, τ 2 is set to 0. Hence we use as estimator τ 2 = max τ 2 , 0 .

3.6 Empirical Credibility Estimator. The empirical credibility estimator is obtained from the credibility formula (3.3) by replacing the structural parameters σ 2 and τ 2 by their estimators derived in Subsection 3.5. Hence, we have

     X i T i +1 = α i X i • + (1 -α i ), α i = τ 2 ω i • σ 2 + τ 2 ω i • . (3.7)
It follows from Equation 3.5, that the mortality time varying coecient is successively updated as follows

B i T i +1 = B b T +1 1 + α i (X i • -1) , (3.8)
and similarly, the forces of mortality and the probabilities of death are given respectively by

ϕ i x,T i +1 = α i 1 -X i • A b + α i X i • -1 + 1 ϕ b x,T +1 ,
and

q i x,T i +1 = q b x,T +1 1 -q b x,T +1 exp(-A b ) 1 α i (X i
• -1) .

(3.9)

4 Classical Credibility Approaches to Mortality

Next, we wish to compare our model to the [START_REF] Hardy | A credibility approach to mortality risk[END_REF] and Poisson-Gamma credibility analysis to mortality. The actual to expected mortality ratio is the key observation that is the focus of the two following approaches. Specically, the a priori expected number of deaths for portfolio i in calendar year t in the age-band x , x is denoted by

ω i t = E D i •,t = x x=x q b x,t L i x,t .
The actual to expected mortality ratios denoted by X i t are computed for each calendar year t in aggregate for each portfolio i,

X i t = D i •,t E D i •,t = D i •,t ω i t .
Both estimators for τ 2 and σ 2 . The Hardy and Panjer (1998) credibility approach to mortality estimates the structure parameters from the aggregated data using the estimators derived in

σ 2 0 = 1 n i=1 C i n i=1 C i s 2 i , where C i = 1 1 + 2 T i -1 φ with φ = E σ 4 Θ i V σ 2 Θ i .
Again, as the risk parameters, Θ i n i=1

, are assumed to be identically distributed, the factors E σ 4 Θ i and V σ 2 Θ i are independent of the portfolio. Hence, the only portfolio dependent variable in C i is T i , the number of years data available for the portfolio.

Both methods give the same result. In addition, the latter approach, derived from [START_REF] Centeno | The bühlmann-straub model with the premium calculated according to the variance principle[END_REF], allows to obtain a credibility estimator for the variance part of the credibility premium which has the form:

σ i 2 = C i s 2 i + (1 -C i ) σ 2 0 .
The estimate of V µ Θ i denoted by τ 2 is

τ 2 = ω • • W -σ 2 ω • • Ω
, where Ω = 1

n i=1 T i -1 n i=1 ω i • ω • • 1 - ω i • ω • • , and W = 1 n i=1 T i -1 n i=1 T i t=1 ω i t ω • • X i t -X • • 2 , with ω • • = n i=1 ω i • , and X • • = 1 ω • • n i=1 ω i • X i • .
Then, the estimate of E σ 4 Θ i denoted by σ 4 is

σ 4 = 1 n i=1 T i + 1 n i=1 T i -1 s 2 i 2 ,
and the estimate of V σ 2 Θ i denoted by υ σ 2 is

υ σ 2 = 1 R n i=1 T i -1 s 2 i -β 2 2 -2 σ 4 (n -1) , where R = n i=1 (T i -1) - n i=1 (T i -1) 2 n i=1 (T i -1)
, and

β 2 = 1 n n i=1 s 2 i .
4.2 The Poisson-Gamma Approach. A priori, we could assume that E Θ i = 1 so that the baseline mortality produces the a priori expected number of deaths,

E D i •,t = E ω i t Θ i = ω i t .
We suppose here that the parametric distribution for the number of deaths D i

•,t is Poisson conditional to the relative risk level Θ i , so that E D i •,t |Θ i = V D i •,t |Θ i = ω i t Θ i .
Then, under assumption H1, Subsection 3.3, the conditional mean and variance of the actual to expected mortality ratios become:

E X i t | Θ i = µ Θ i = Θ i and V X i t | Θ i = σ 2 Θ i ω i t = Θ i ω i t ,
and the p.d.e with respect to a i,0 and a i,t , Equation 3.4 are:

a i,0 = 1 - τ 2 ω i • 1 + τ 2 ω i • and a i,t = τ 2 ω i t 1 + τ 2 ω i • , since σ 2 = E Θ i = 1.
Then the linear credibility estimator is given by

µ Θ i = X i T i +1 = 1 1 + τ 2 ω i • + τ 2 ω i • 1 + τ 2 ω i • 1 ω i • T i t=1 ω i,t X i,t . (4.1)
And, the expected number of deaths for portfolio i for next year

T i + 1 is ω i T i +1 X i T i +1 = ω i T i +1 1 + τ 2 D i •,• 1 + τ 2 ω i • .
Then, we need to obtain the structure parameter τ 2 = V Θ i . As the distribution of the total number of deaths in portfolio i is D i •,• ∼ MP ω i • Θ i and using the variance decomposition principle,

V D i •,• = V E D i •,• | Θ i + E V D i •,• | Θ i = V ω i • Θ i + E ω i • Θ i = τ 2 (ω i • ) 2 + ω i • .
And,

n i=1 V D i •,• = ω 2 n i=1 (ω i • ) 2 + n i=1 ω i • , leads to τ 2 = n i=1 V D i • -ω i • / n i=1 (δ i • ) 2 . Thus, the estimator of τ 2 writes τ 2 = n i=1 D i •,• -ω i • 2 -D i •,• n i=1 (ω i • ) 2
.

5 Numerical Analysis 5.1 Data Quantitative Analysis. The data come from studies conducted by Institut des

Actuaires. These studies include in total 14 portfolio covering the period 2007-2011 with each companies contributing data for at least 4 of a possible 5 years. Table 1 presents the observed characteristics of the male population of the portfolios. For this dataset, we are considering respectively T i = 3 and T i = 4 for all companies. The remaining years serve to test the predictive feature of the model through an in-sample analysis. The age band for all companies ranges from 30 to 95 years old. Figure 1 shows the age distribution of two portfolios. It graphically depicts the heterogeneity observed between the portfolios with insureds holding dierent policies.

5.2 The Baselines Mortality. We consider two prospective tables as baselines for our credibility models. One is the national demographic projections for the French population over

A b t = A b T and C b t = C b
T remain xed. Table 3 presents the estimated parameters for each year and baselines considered.

5.4 Proximity Between the Observations and the Model. We assess the overall deviation with the observed mortality by comparing criteria measuring the distance between the observations and the models with the χ 2 applied by [START_REF] Forfar | On graduation by mathematical formula[END_REF], the mean average (2002) as well as the standardized mortality ratio (SMR) and the number of standardized residuals larger then 2 and 3, see [START_REF] Tomas | Constructing entity specic projected mortality table: adjustment to a reference[END_REF]. In addition, we nd useful to use the SMR test proposed by [START_REF] Liddell | Simple exact analysis of the standardised mortality ratio[END_REF] and the likelihood ratio test. The tests and quantities summarizing the proximity between the observations and the model are described in the following. The χ 2 allows to measure the quality of the t of the model. It writes,

χ 2 = (x,t) D x,t -L x,t q x (t) 2 L x,t q x (t) 1 -q x (t) .
The MAPE is the average of the absolute values of the deviations from the observations,

MAPE = (x,t) D x,t /L x,t -q x (t) / D x,t /L x,t (x,t) D x,t
× 100.

We can also determine if the t corresponds to the underlying mortality law (null hypothesis H 0 ) with the likelihood ratio test. The statistic, ξ LR , writes

ξ LR = (x,t) D x,t ln D x,t L x,t q x (t) + L x,t -D x,t ln L x,t -D x,t L x,t -L x,t q x (t)
. .

If H 0 is true, this statistic follows a χ 2 law with a number of degrees of freedom equal to the number of observations n: ξ LR ∼ χ 2 (n). Hence, the null hypothesis H 0 is rejected if ξ LR > χ 2 1-α (n), where χ 2 1-α (n) is the (1 -α) quantile of the χ 2 distribution with n degrees of freedom. The p-value is the lowest value of the type I error (α) for which we reject the test. We will privilege the model having the p-value

= P χ 2 1-α (n) > ξ LR = 1 -F χ 2 (n) (ξ LR ) closest to 1.
The SMR is computed as the ratio between the observed and tted number of deaths:

SMR = (x,t) D x,t (x,t) L x,t q x (t)
.

Hence, if SMR > 1, the tted deaths are under-estimated and vice-versa if SMR < 1. Note that we can consider the SMR as a global criterion which does not take the age structure into account, compared to the chi2 and MAPE for instance. We can also apply a test to determine if the SMR is signicatively dierent from 1. [START_REF] Liddell | Simple exact analysis of the standardised mortality ratio[END_REF] proposes to compute the statistic,

ξ SMR = 3 × D 1 2 1 -(9D) -1 -(D/E) 1 3 If SMR > 1, 3 × D * 1 2 (D * /E) 1 3 + (9D * ) -1 -1 If SMR < 1, where D = (x,t) D x,t , D * = (x,t) D x,t + 1 and E = (x,t) L x,t q x (t).
If the SMR is not signicatively dierent from 1 (null hypothesis H 0 ), this statistic follows a standard Normal law, ξ SMR ∼ N(0, 1). Thus, the null hypothesis H 0 is rejected if ξ SMR > N 1-α (0, 1), where N 1-α (0, 1) is the (1 -α) quantile of the standard Normal distribution. The p-value is given by p-value = 1 -F N(0,1) (ξ SMR ).

5.5

In-Sample Numerical Analysis. We tted the approaches over a history covering 3 and 4 years (2007-2009 and 2007-2010 respectively) and compared the overall deviation between the observations and the models (for the year 2010 and 2011 respectively). Table 4 displays the estimates of the structure parameters for the three approaches.

Table 5 presents the tests and quantities summarizing the overall deviation between the observations and the credibility analysis for the male population of portfolio 1 obtained by the The Hardy-Panjer and Poisson-Gamma approaches produce relatively similar graduations.

However, we notice some dierences with the Makeham credibility model which displays more favorable results whatever the baseline mortality considered for the two periods tted.

It is also apparent that using the market baseline mortality IA2013 produces better results than the national demographic projections originating from INSEE, see Subsection 5.2. It illustrates the importance of using an adequate baseline mortality when adjusting the models.

When looking at criteria and quantities which take the age structure of the error into account, the Makeham credibility approach is a benet. The quality of the t increases, sometimes For clarity, the graphical comparisons only consider the market baseline mortality IA2013 as it leads to better results than using the national demographic projections.

(a) Fitted probabilities of death in the log scale. Figure 2a displays the the tted probabilities of death in the log scale for portfolio 1 for the year 2010. Figure 3 and 4 in Appendix C and D display the comparisons for all the portfolios and for the years 2010 and 2011 respectively. It gives us the opportunity to visualize the similarities and dierences between the ts obtained by the approaches. It is again apparent that the Hardy-Panjer and Poisson-Gamma models lead to similar results. In addition, we observe that these approaches have a tendency to strongly overestimate the probabilities of death for the age band 30, 60 and reciprocally underestimate them for the age band 60, 95 . This is explained by the fact that the age structure is not taken in account by the Hardy-Panjer and Poisson-Gamma approaches, conversely to the Makeham credibility model. We can visualize this lack of t in the plots of the tted number of deaths, Figure 2b for portfolio 1 and Figure 5 and 6 for all portfolios in Appendix E and F.

In conjunction with looking to the plots of the ts, we should study the residuals plots.

Such residual plots provide a powerful diagnostic that nicely complements the analysis. The diagnostic plots can show lack of t locally and we have the opportunity to judge the lack of t based on our knowledge on the data and of the performance of the models. We superimposed a smooth curve on the standardized residuals. This smooth helps search for clusters of residuals that may indicate a lack of t. The plots of the standardized residuals, for the male population, are display in Figure 2c for portfolio 1 and

Concluding Remarks

We considered the periodic adjustment of a mortality graduated curve using a Makeham parametric model. This relies on the revision of a single parameter the two remaining been xed.

The framework considered here is closely related to the one introduced in Hardy and Panjer (1998). The main dierence is the age-structure included through the parametric Makeham model. By doing so, we showed that adding an age structure enhances the predictive ability of the death forecast especially when we consider age-sensitive proxies. If one is only interested in predicting deaths at the aggregate portfolio level our methodology yields to the same forecast as in the [START_REF] Hardy | A credibility approach to mortality risk[END_REF] framework. Moreover, we should note that in our methodology especially using the ratio of the considered Makeham parameters allows to overcome the de-trending step recommended in [START_REF] Hardy | A credibility approach to mortality risk[END_REF].

In order to assess the predictive power of our methodology, various other measures of risk and goodness-of-t should be taken into account. Especially, we should consider the age-structure's impact on the prices and reserves and potential benet of our model compared to the current market practice. There are also several piratical we do not address here which we openly acknowledge and leave for future research.

A Tests and quantities summarizing the deviation between the observations and the models for the year 2010

Tables 6 and7 present the tests and quantities summarizing the overall deviation between the observations and the credibility analysis for the male population obtained by the Hardy-Panjer,

Poisson-Gamma and the Makeham credibility approaches with the two baselines mortality considered for the year 2010.

B Tests and quantities summarizing the deviation between the observations and the models for the year 2011

Tables 8 and9 present the tests and quantities summarizing the overall deviation between the observations and the credibility analysis for the male population obtained by the Hardy-Panjer,

Poisson-Gamma and the Makeham credibility approaches with the two baselines mortality considered for the year 2011.

C tted probabilities of death in the log scale for the year 2010 E Fitted number of deaths for the year 2010

Figure 5 displays the tted number of deaths for the male population for the year 2010.

F Fitted number of deaths for the year 2010

Figure 6 displays the tted number of deaths for the male population for the year 2010.

G Standardized residuals for the year 2010

Figure 7 displays the standardized residuals for the male population for the year 2010.

H Standardized residuals for the year 2011

Figure 8 displays the standardized residuals for the male population for the year 2011.

INSEE IA2103

  the Hardy and Panjer (1998) and the Poisson-Gamma credibility approaches are using the Bühlmann-Straub set-up, see Subsection 3.3. Again, the key determinants of the credibility estimator (3.5) are the structure parameters, i.e. the variance part of the credibility premium E σ 2 Θ i denoted by σ 2 and the uctuation part V µ Θ i denoted by τ 2 .4.1 The Hardy-Panjer Approach. As in general we have no knowledge or, at least, no exact knowledge of the parametric distributions for the number of deaths or of the structure distribution, we need estimators for the two components of the credibility estimator (3.5), i.e.

Figure 1 :

 1 Figure 1: Distribution of age groups in portfolios 3 (left panel) and 8 (right panel), male population.

  drastically, compared to the Hardy-Panjer and Poisson-Gamma model in terms of having the minimum χ 2 and MAPE values. The Makeham credibility model leads to the lowest number of standardized residuals lower than 2 and 3. It exhibits as well the highest p-value for the likelihood ratio test.Even when we considering a global indicator of the quality of the t such as the SMR which does not take the age structure into account, the Makeham credibility model seems to perform better than the Hardy-Panjer and Poisson-Gamma approaches. The statistic ξ SMR of the SMR test is the smaller 8 times over 14 for the year 2010, see Tables6 and 7in Appendix A, and 6 times over 12 for the year 2011, see Tables8 and 9in Appendix B. number of deaths, having a SMR lower than 1 for 9 portfolios over 14 in 2010 and for 8 portfolios over 12 in 2011.In the following, these quantitative diagnostics are supplemented by a range of visual comparisons. Besides the tests and quantities, the comparison involves graphical analysis. It consists of representing graphically the tted values against the observations for the years 2010 and 2011.

  (b) Fitted number of deaths. (c) Standardized residuals.

Figure 2 :

 2 Figure 2: Fitted values against the observations for portfolio 1 for the year 2010, male population.

  Figure 7 and 8 in Appendix G and H for all the portfolios and for the years 2010 and 2011 respectively. The standardized residuals, obtained by the Hardy-Panjer and Poisson-Gamma models, present a high curvature for most of the portfolios in Figure 7 and 8. It indicates a clear lack of t. These models overestimate the number of deaths for the age band 30, 60 et underestimate them for the age band 60, 95 , as observed in the plots of the ts previously. Conversely, no strong patterns appear in the standardized residuals retrieved for the Makeham credibility model. The smooth curves over the standardized residuals is meanly at, meaning that no systematic reproducible lack of t has been detected and that the Makeham credibility model captures adequately the variability of the data.

Figure 3

 3 Figure3displays the tted probabilities of death in the log scale for the male population for the year 2010.

Figure 3 :Figure 6 :

 36 Figure 3: Fitted probability of death, log scale, for the year 2010, male population

Table 1 :

 1 the period 2007-2060, provided by the French National Oce for Statistics, INSEE, see[START_REF] Blanpain | Projections de populations 2007-2060 pour la France métropolitaine: méthode et principaux résultats[END_REF]. These projections are based on assumptions concerning fertility, mortality and migrations. We choose the baseline scenario among a total of 27 scenarios. The baseline Observed characteristics of portfolios population.

		Period of observation	Mean age Average	Mean age
		Beginning	End	In	Out	exposure	at death
		1/1/07	12/31/11	36.96	39.74	2.77	68.78
		1/1/07	12/31/11	69.3	73.35	4.05	80.34
		1/1/07	12/31/10	40.16	43.1	2.94	71.77
		1/1/07	12/31/11	37.5	41.13	3.63	54.08
		1/1/07	12/31/11	36.9	39.1	2.2	59.31
		1/1/07	12/31/10	48.5	52.11	3.62	82.34
		1/1/07	12/31/11	66.65	71.29	4.64	73.68
		1/1/07	4/13/11	67.51	71.38	3.86	80.72
		1/1/07	6/30/11	45.97	49.6	3.62	73.17
	10	1/1/07	12/31/11	62.97	67.64	4.67	79.77
	11	1/1/07	12/31/11	38.89	42	3.11	56.44
	12	1/1/07	12/31/11	37.05	39.2	2.15	57.41
	13	1/1/07	12/31/11	43.01	46.89	3.88	71.03
	14	1/1/07	12/31/11	50.12	54.16	4.04	72.37

scenario is based on the assumption that until 2060, the total fertility rate is remaining at a very high level (1.95). The decrease in sex and age-specic mortality rates is greater for men over 85 years old. The baseline assumption on migration consists in projecting a constant annual net-migration balance of 100, 000 inhabitants. The second external reference table, denoted IA2013, is a market table constructed for the French insurance market provided by Institute des Actuaires, see

[START_REF] Tomas | Construction et validation des références de mortalité de place[END_REF]

. It is worth to mention that this table is derived on mortality trends originating from the INSEE table and covers the period 2007-2060.

Following, assumption (i) in Subsection 2.3, the baseline mortality q b

x,t is described by the Makeham model in (2.3). Table

2

presents the estimated parameters for each of the baselines considered.

Table 2 :

 2 Estimated parameters of the Makeham model (2.3) for the baselines of mortality

	considered, male population.
		INSEE		IA2013
		20072009 20072010 20072009 20072010
	A b T	4.2835e -03 4.2787e -03 2.1577e -04 2.4355e -04
	B b T	7.9564e -07 7.7199e -07 4.0863e -06 3.9935e -06
	C b T	1.1484	1.1487	1.1211	1.1213
	5.3 Adjustment of the Makeham model. Following assumptions (ii) in Subsection 2.3, we
	t the Makeham model (2.3) for the baselines of mortality considered so as to estimate B b t for
	each calendar year while the parameters

Table 3 :

 3 Estimated

			INSEE		IA2013
		20072009 20072010 20072009 20072010
	A b T	4.2835e -03 4.2787e -03 2.1577e -04 2.4355e -04
	B b 2007	8.0826e -07 7.9035e -07 4.1740e -06 4.1204e -04
	B b 2008	7.9554e -07 7.7790e -07 4.0843e -06 4.0319e -06
	B b 2009	7.8318e -07 7.658e -07 4.0009e -06 3.9496e -06
	B b 2010	-	7.5406e -07	-	3.8729e -06
	C b T	1.1484	1.1487	1.1211	1.1213
	percentage error (MAPE) applied by Felipe et al.

parameters of the Makeham model (2.3) for each year and baselines of mortality considered, male population.

Table 4 :

 4 Estimates of the structure parameters, male population. Gamma and the Makeham credibility approaches with the two baselines mortality considered for the year 2010. Tables 6, 7 and 8, 9 in Appendix A and B display the results for all the portfolios and for the years 2010 and 2011 respectively.

		Hardy-Panjer	Poisson-Gamma Makeham-Credibility
		INSEE	IA2013	INSEE	IA2013	INSEE	IA2103
	2007-09	µ 0 3.5521 σ 2 44.4032 92.1668 16.3290 τ 2 6.8368 44.0092	1 1 10.7485	1 1 367.5029	1 4.0552e-04 0.1935	1 2.3198e-03 3.5960e-02
	2007-10	µ 0 3.6495 σ 2 65.9649 116.0159 15.7865 τ 2 7.0772 43.4966	1 1 10.9684	1 1 338.4440	1 5.1034e-04 0.2217	1 2.6285e-03 5.0281e-02
	Hardy-Panjer, Poisson-

Table 5 :

 5 Tests and quantities summarizing the deviation between the observations and the models for portfolio 1, calendar year 2010, male population.

				INSEE			IA2103
			Hardy-Panjer	Poisson-Gamma	Makeham-Credibility	Hardy-Panjer	Poisson-Gamma	Makeham-Credibility
	Standardized	> 2	60	60	35	46	46	15
	residuals	> 3	48	48	28	32	32	5
	χ 2		5481.86	5542.82	3569.97	1705.25	1747.25	208.81
	MAPE (%)	233.22	230.94	373.89	117.01	115.42	42.35
	Likelihood	ξ LR	946.98	947.72	443.16	463.48	468.46	88.03
	ratio test	p-value	0	0	0	0	0	0.0364
		SMR	1.1792	1.1919	0.5265	1.7629	1.7957	1.0532
	SMR test	ξ SMR	4.0379	4.2939	12.1893	13.0352	13.4202	1.2845
		p-value	0	0	0	0	0	0.0995

Figure 4: Fitted probability of death, log scale, for the year 2011, male population

Figure 8: Standardized residuals, calendar year 2011, male population

* This work beneted from the nancial support of the ANR project LoLitA (ANR-13-BS01-0011) and the BNP Paribas Cardif Chair "Management de la Modélisation".