A Antonelli

B Desmorat

B Kolev

R Desmorat

DISTANCE TO PLANE ELASTICITY ORTHOTROPY BY EULER-LAGRANGE METHOD

Keywords: Mathematics Subject Classification. 74B05, 74E10. Key words and phrases. anisotropy, distance to a symmetry class, orthotropy, plane elasticity. 1

Constitutive tensors are of common use in mechanics of materials. To determine the relevant symmetry class of an experimental tensor is still a tedious problem. For instance, it requires numerical methods in three-dimensional elasticity. We address here the more affordable case of plane (2D) elasticity, which has not been fully solved yet. We recall first Vianello's orthogonal projection method, valid for both the isotropic and the square symmetric (tetragonal) symmetry classes. We then solve in a closed-form the problem of the distance to plane elasticity orthotropy, thanks to the Euler-Lagrange method.

Introduction

Experimental techniques have been developed in order to measure all the components of an elasticity tensor, such as the ultrasonic measurements [START_REF] Krishnan | Elastic constants of triclinic copper sulphate pentahydrate crystals[END_REF][START_REF] Arts | General anisotropic elastic tensors in rocks: Approximation, invariants, and particular directions[END_REF][START_REF] Arts | A study of general anisotropic elasticity in rocks by wave propagation[END_REF][START_REF] François | Une nouvelle analyse des symétries d'un matériau élastique anisotrope. exemple d'utilisation à partir de mesures ultrasonores[END_REF][START_REF] Dellinger | Computing the optimal transversely isotropic approximation of a general elastic tensor[END_REF][START_REF] Brown | Triclinic elastic constants for low albite[END_REF]. Once its components C raw ijkl are identified (in a working basis), one may wonder what is its relevant material symmetry (i.e. symmetry group). For single crystals or for composite materials, a material symmetry is expected (for example cubic symmetry or orthotropy) so that the question becomes to find a tensor C, with a given material symmetry, the closest to the measured tensor C raw (usually with no material symmetry at all, i.e. triclinic in 3D, biclinic in 2D).

The precise identification of the symmetry class -or the symmetry group -of a constitutive tensor is not an easy task, mainly for two reasons [START_REF] François | Une nouvelle analyse des symétries d'un matériau élastique anisotrope. exemple d'utilisation à partir de mesures ultrasonores[END_REF]:

• the tensor is measured in a given orientation, possibly experimenter dependent, which may not coincide with an expected symmetry group, not allowing then for the recognition by eye (on the Kelvin matrix representation) of the orthotropic, tetragonal, cubic, . . . well-known expressions for C raw (the well-known normal forms of elasticity tensors [START_REF] Fedorov | Theory of Elastic Waves in Crystals[END_REF]), • the experimental discrepancy/errors makes the material symmetries approximate. Sufficient conditions (in [START_REF] Bóna | Coordinate-free characterization of the symmetry classes of elasticity tensors[END_REF]) and necessary and sufficient conditions (in [START_REF] Olive | Characterization of the symmetry class of an elasticity tensor using polynomial covariants[END_REF]) have been formulated to characterize, in an arbitrarily oriented coordinate system, the symmetry class of a three-dimensional elasticity tensor. For plane (2D) elasticity, the symmetry classes characterization is much simpler than in 3D [START_REF] Blinowski | Two-dimensional Hooke's tensors-isotropic decomposition, effective symmetry criteria[END_REF][START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF][START_REF] Vannucci | Stiffness design of laminates using the polar method[END_REF]. Methods furthermore exist to bring back a rotated elasticity in its normal form [START_REF] Baerheim | Classification of symmetry by means of maxwell multipoles[END_REF][START_REF] Zou | Identification of symmetry type of linear elastic stiffness tensor in an arbitrarily orientated coordinate system[END_REF][START_REF] Abramian | Recovering the normal form and symmetry class of an elasticity tensor[END_REF]. Unfortunately, they are a priori useless in the common case of a triclinic/biclinic measured elasticity tensor. When experimental discrepancy has to be dealt with, the literature approaches are based on the concept of distance of a tensor to a considered symmetry class [START_REF] Gazis | The elastic tensor of given symmetry nearest to an anisotropic elastic tensor[END_REF][START_REF] Fedorov | Theory of Elastic Waves in Crystals[END_REF][START_REF] Francois | Détermination des symétries matérielles de matériaux anisotropes[END_REF][START_REF] François | Une nouvelle analyse des symétries d'un matériau élastique anisotrope. exemple d'utilisation à partir de mesures ultrasonores[END_REF][START_REF] Helbig | 3. representation and approximation of elastic tensors[END_REF][START_REF] François | Determination of the symmetries of an experimentally determined stiffness tensor: Application to acoustic measurements[END_REF][START_REF] Moakher | The closest elastic tensor of arbitrary symmetry to an elasticity tensor of lower symmetry[END_REF][START_REF] Kochetov | On obtaining effective orthotropic elasticity tensors[END_REF][START_REF] Diner | Identifying symmetry classes of elasticity tensors using monoclinic distance function[END_REF], starting from a given (usually measured) elasticity tensor C raw with no material symmetry, and sometimes from the additional quantification of the measurement errors [START_REF] Bóna | Symmetry characterization and measurement errors of elasticity tensors[END_REF][START_REF] Danek | Effective elasticity tensors in context of random errors[END_REF][START_REF] Guilleminot | A stochastic model for elasticity tensors with uncertain material symmetries[END_REF].

The 3D case has, by far, been the most studied. It remains the most challenging case, and determining the distance to a 3D symmetry class is usually done numerically, with the risk of reaching a local minimum instead of the global minimum [START_REF] François | Determination of the symmetries of an experimentally determined stiffness tensor: Application to acoustic measurements[END_REF][START_REF] Kochetov | On obtaining effective orthotropic elasticity tensors[END_REF][START_REF] Danek | Effective elasticity tensors in context of random errors[END_REF][START_REF] Moakher | The closest elastic tensor of arbitrary symmetry to an elasticity tensor of lower symmetry[END_REF]. Recently, in [START_REF] Stahn | Distances of stiffnesses to symmetry classes[END_REF], following [START_REF] Baerheim | Harmonic decomposition of the anisotropic elasticity tensor[END_REF], the authors have used the basis carried by a second-order deviatoric tensor, a linear covariant t(C asym ), built from the Backus asymmetric part C asym = C raw -(C raw ) s of the raw elasticity tensor. They obtain, this way, closed-form expressions for upper bounds estimates of the distances to the elasticity symmetry classes. In 2D, the asymmetric part C asym is always isotropic (see remark 2.2) so that, in order to apply this procedure, another (non isotropic) second-order covariant of C raw has to be used.

The 2D case has been shown to be more affordable than the 3D case [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF]. Indeed:

• There are only four symmetry classes, easily characterized (instead of eight in threedimensional elasticity [START_REF] Forte | Symmetry classes for elasticity tensors[END_REF]); • Both the isotropic and the tetragonal (square symmetry) strata1 are linear subspaces of the vector space of plane elasticity tensors Ela [START_REF] Abud | The geometry of spontaneous symmetry breaking[END_REF]. Closed form-expressions for the distance to these symmetry strata are obtained using orthogonal projections on these subspaces.

The problem of calculating the distance to plane elasticity orthotropy has also been posed in [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF], but not solved explicitly. An upper bound estimate, expressed in terms of invariants, of this distance has been obtained in [START_REF] Oliver-Leblond | Continuous anisotropic damage as a twin modelling of discrete bi-dimensional fracture[END_REF]. It is the aim of the present work to derive a closed form solution for the exact 2D distance problem. For the sake of self-completeness, we will recall Vianello's orthogonal projection method for the isotropic and the square symmetric cases (section 3). We will determine, by the Euler-Lagrange method, the closest orthotropic plane elasticity tensor to a given raw tensor C raw ∈ Ela(2). The method is detailed in section 4 in a sufficiently general manner so that it also applies to other situations. For instance, the case of the distance to a transversely isotropic second order (3D) constitutive tensor (such as a thermal expansion or a thermal conductivity tensor) is treated in Appendix B. A discussion on different equivalent definitions of the distance to a symmetry class is provided in section 5. In Appendix C, we give a proof that, concerning the distance to orthotropic plane elasticity, the upper bounds estimates proposed in [START_REF] Stahn | Distances of stiffnesses to symmetry classes[END_REF] and in [START_REF] Oliver-Leblond | Continuous anisotropic damage as a twin modelling of discrete bi-dimensional fracture[END_REF] are equal.

Using the Euclidean structure of space, no distinction will be made between covariant, contravariant or mixed tensors. All tensor components will be expressed with respect to an orthonormal basis.

Plane elasticity tensors

In this paper, we consider the natural action of the orthogonal group O(2) on the vector space of elasticity tensors in two dimensions (sometimes, called plane elasticity tensors [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF][START_REF] Vannucci | Stiffness design of laminates using the polar method[END_REF])

Ela(2) = {C; C ijkl = C jikl = C ijlk = C klij } ,
which describe constitutive equations for two-dimensional linearly elastic bodies. The orthogonal group O(2) consists in two components. The subgroup, noted SO(2), of rotations r θ ,

r θ := cos θ -sin θ sin θ cos θ , (2.1)
and the set (which is not a subgroup), noted SO(2)s, consisting of axial symmetries with respect to the axis with angle θ/2 with Ox. These reflections are written as the composition

r θ s = cos θ sin θ sin θ -cos θ , where s := 1 0 0 -1 , (2.2)
(s is the symmetry with respect to the x-axis) and the identity element is denoted by e := 1 0 0 1 .

The action of the full orthogonal group O(2) on the vector space Ela(2) is then written as

(ρ(g)C) ijkl = g ip g jq g kr g ls C pqrs , g ∈ O(2), C ∈ Ela(2).

2.1.

Plane elasticity symmetry classes. The vector space of plane elasticity tensors splits into four symmetry strata Σ [H] [START_REF] Verchery | Une méthode d'aide graphique à la conception des séquences d'empilement dans les stratifiés[END_REF][START_REF] Blinowski | Two-dimensional Hooke's tensors-isotropic decomposition, effective symmetry criteria[END_REF][START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF][START_REF] Vannucci | Plane anisotropy by the polar method[END_REF], a symmetry stratum being the set of all tensors having the same symmetry class. They are indexed by the conjugacy class [H] of some symmetry group H (see Appendix A) and are naturally ordered as follows biclinic -→ orthotropic -→ square symmetric -→ isotropic (2.3)

In the following we will refer to a symmetry stratum by its mechanical name, setting Σ ortho , Σ square and Σ iso , and generically Σ, instead of the notation Σ [H] (which uses the definitions of subgroups of the orthogonal group O(2)).

Raw plane elasticity tensors are in general biclinic but one may expect, for several mechanical/material science reasons, a particular symmetry, the orthotropic, the square symmetry (also called tetragonal in [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF]) or the isotropic symmetry. A natural question is then: how far is a given raw elasticity tensor C raw ∈ Ela(2) from an orthotropic, a square symmetric or an isotropic tensor? In the following, non biclinic plane elasticity tensors, which are either orthotropic, square symmetric or isotropic, will be said to be at least orthotropic and the set of all tensors which are at least orthotropic will be called the closed orthotropic stratum and noted by Σ ortho . Hence

Σ ortho = Σ iso ∪ Σ square ∪ Σ ortho .
Similarly, at least square symmetric elasticity tensors are either square symmetric or isotropic, and belong to the closed symmetry stratum

Σ square = Σ iso ∪ Σ square .
The arrows in (2.3) represent this "at least of a given symmetry" ordering.

Normal forms.

It is well-known that any symmetric second-order tensor a can be rotated onto a diagonal tensor A. This diagonal representation is called a normal form for a under the action of the rotation group. Similarly, any plane elasticity tensor C (which is at least orthotropic) can be written as C = ρ(g)A where A ∈ Ela(2) has the following Kelvin matrix representation 2.3. Harmonic decomposition. Following [START_REF] Blinowski | Two-dimensional Hooke's tensors-isotropic decomposition, effective symmetry criteria[END_REF][START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF], we denote by C = (λ, µ, h, H) the following explicit harmonic decomposition of a plane elasticity tensor C, into two scalars λ, µ, a secondorder deviatoric (harmonic) tensor h and a fourth order harmonic tensor H:

A =   A 1111 A 1122 0 A 1122 A 2222 0 0 0 2A 1212   . ( 2
C := 2µ I + λ 1 ⊗ 1 + 1 2 (1 ⊗ h + h ⊗ 1) + H, (2.5)
where I is the fourth order identity tensor with components I ijkl = 1 2 (δ ik δ jl + δ il δ jk ). Explicitly, one has

λ := 1 8 (3C iijj -2C ijij ) = 1 8 (C 1111 + 6C 1122 -4C 1212 + C 2222 ), (2.6
)

µ := 1 8 (2C ijij -C iijj ) = 1 8 (C 1111 -2C 1122 + 4C 1212 + C 2222 ), (2.7
)

h := (tr 12 C) ′ i.e. h ij = C kkij - 1 2 C kkll δ ij , (2.8)
H := C - 1 2 (1 ⊗ h + h ⊗ 1) -2µ I -λ1 ⊗ 1, (2.9)
where (•) ′ means the deviatoric part. Recall that an harmonic tensor is a totally symmetric and traceless tensor.

The two scalars λ and µ and

I 2 (C) := h 2 = 1 2 (C 1111 -C 2222 ) 2 + 2(C 1112 + C 2212 ) 2 , J 2 (C) := H 2 = 1 8 (C 1111 -2C 1122 -4C 1212 + C 2222 ) 2 + 2 (C 1112 -C 1222 ) 2 ,
are invariants of the elasticity tensor. Let (•) s be the total symmetrization and ⊙ be the (totally) symmetric tensor product. For symmetric second-order tensors a, b, we have

a ⊙ b = (a ⊗ b) s , (a ⊙ b) ijkl = 1 6 (a ij b kl + a kl b ij + a ik b jl + b ik a jl + a il b jk + b il a jk ),
The explicit harmonic decomposition (2.5) can be rewritten as

C = 2 3 (µ -λ)(I -1 ⊗ 1) + (2µ + λ)1 ⊙ 1 + 1 2 (1 ⊙ h + h ⊙ 1) + H, (2.10)
thanks to the two following properties,

2 I + 1 ⊗ 1 = 3 1 ⊙ 1, 1 ⊗ h + h ⊗ 1 = 1 ⊙ h + h ⊙ 1, (2.11)
the second being specific to the 2D case when h = h ′ is harmonic [START_REF] Desmorat | Second order tensorial framework for 2D medium with open and closed cracks[END_REF]. The totally symmetric part S and the asymmetric part C asym -in the sense of Backus [START_REF] Backus | A geometrical picture of anisotropic elastic tensors[END_REF]-of a plane elasticity tensor C are then

S = C s = (2µ + λ)1 ⊙ 1 + 1 2 (1 ⊙ h + h ⊙ 1) + H, C asym = C -S = 2 3 (µ -λ) (I -1 ⊗ 1) , Remark 2.2.
We deduce from the last equality that, in 2D, the asymmetric part C asym of a plane elasticity tensor is always isotropic.

Closest tensor in a given plane elasticity symmetry class

Given an experimental raw tensor C raw , the aim is to find a plane elasticity tensor C Σ in a given (exact) symmetry stratum Σ which is the closest to C raw , i.e., to calculate

C Σ = argmin C∈Σ C raw -C , (3.1)
a problem which has already been extensively studied [START_REF] Gazis | The elastic tensor of given symmetry nearest to an anisotropic elastic tensor[END_REF][START_REF] Francois | Détermination des symétries matérielles de matériaux anisotropes[END_REF][START_REF] François | Une nouvelle analyse des symétries d'un matériau élastique anisotrope. exemple d'utilisation à partir de mesures ultrasonores[END_REF][START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF][START_REF] François | Determination of the symmetries of an experimentally determined stiffness tensor: Application to acoustic measurements[END_REF]. Here, we choose to work with the Euclidean norm C = C, C (see [START_REF] Moakher | The closest elastic tensor of arbitrary symmetry to an elasticity tensor of lower symmetry[END_REF][START_REF] Morin | Generalized euclidean distances for elasticity tensors[END_REF] for other norms), derived from the O(2)-invariant scalar product

C 1 , C 2 = C 1 :: C 2 = C 1 ijkl C 2 ijkl . Once a tensor C Σ solution of (3.1) is found, the distance to the considered symmetry class Σ is just d(C raw , Σ) = C raw -C Σ .
In general, several methods allow to solve such a minimization problem [START_REF] François | Une nouvelle analyse des symétries d'un matériau élastique anisotrope. exemple d'utilisation à partir de mesures ultrasonores[END_REF][START_REF] Kochetov | On obtaining effective orthotropic elasticity tensors[END_REF][START_REF] Danek | Effective elasticity tensors in context of random errors[END_REF][START_REF] Moakher | The closest elastic tensor of arbitrary symmetry to an elasticity tensor of lower symmetry[END_REF], most of them appealing to numerical minimization schemes. However, concerning the problem for the isotropic and the square symmetry plane elasticity symmetry classes, observe that both of the corresponding closed strata are linear subspaces of Ela(2), which are described respectively by the following linear equations

• h = 0 and H = 0 for Σ iso = Σ iso , • h = 0 for Σ square = Σ iso ∪ Σ square ,
where h and H are the second and fourth order components of the harmonic decomposition C = (λ, µ, h, H). The distances to these symmetry strata are then obtained using orthogonal projections on these subspaces and no further calculations are required. The fact that the isotropic stratum is a linear subspace is very general but the fact that square symmetry stratum is a linear subspace is exceptional (symmetry strata are not in general linear spaces [START_REF] Abud | The geometry of spontaneous symmetry breaking[END_REF][START_REF] Auffray | On anisotropic polynomial relations for the elasticity tensor[END_REF]) and seems to have been noticed first by Vianello [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF]. For the sake of self-completeness, we recall with full details these analytical solutions below, which are based on the harmonic decomposition of elasticity tensors.

3.1. Closest isotropic tensor. The isotropic elasticity tensor C iso the closest to C raw is the orthogonal projection on the isotropic stratum

C iso = 2µ I + λ1 ⊗ 1, (3.2)
where the Lamé constants of C iso are determined from C raw , as

λ := 1 8 (3C raw iijj -2C raw ijij ) = 1 8 (C raw 1111 + 6C raw 1122 -4C raw 1212 + C raw 2222 ), µ := 1 8 (2C raw ijij -C raw iijj ) = 1 8 (C raw 1111 -2C raw 1122 + 4C raw 1212 + C raw 2222 ).
This has allowed Vianello to obtain a nice expression for the distance of C raw to isotropy

d(C raw , Σ iso ) = I 2 (C raw ) + J 2 (C raw ) 3.2.
Closest square symmetric tensor. In the same manner, the at least square symmetric (tetragonal) elasticity tensor C square the closest to C raw is the orthogonal projection of C raw on the square symmetry stratum

C square = C iso + H, (3.3)
with both C iso (by (3.2)) and H (by (2.6)-(2.9)) determined from C raw . If H = 0, C square is square symmetric. The Kelvin matrix representation of the square symmetric elasticity tensor which is the closest to C raw is then

C square =    1 2 (C raw 1111 + C raw 2222 ) C raw 1122 √ 2 2 (C raw 1112 -C raw 2212 ) C raw 1122 1 2 (C raw 1111 + C raw 2222 ) √ 2 2 (C raw 2212 -C raw 1112 ) √ 2 2 (C raw 1112 -C raw 2212 ) √ 2 2 (C raw 2212 -C raw 1112 ) 2C raw 1212    .
and the distance of C raw to square (tetragonal) symmetry, obtained by Vianello, is the nice formula d(C, Σ square ) = I 2 (C raw ).

Closest orthotropic tensor using Euler-Lagrange method

Starting from a raw (usually experimental) plane elasticity tensor C raw ∈ Ela(2), assumed biclinic, it remains to determine the closest orthotopic tensor C ortho to C raw . Since the orthotropic stratum is not a linear subspace, a different approach is required. In what follows, we will apply the Euler-Lagrange method, using the parametrization C = ρ(g)A [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF][START_REF] Dellinger | Computing the optimal transversely isotropic approximation of a general elastic tensor[END_REF], where g ∈ O(2) and A is the normal form (2.4).

Remark 4.1. The parametrization C = ρ(g)A is not one to one, because several orthogonal transformations g may correspond to the same tensor C = ρ(g)A, depending on the symmetry group of A. In particular, redundant solutions are obtained as C = ρ(gh)A, when h is in the symmetry group of A.

Since the reflection s belongs to the symmetry group of A, namely the dihedral group D 2 , we deduce that for each solution (g, A) of the optimization problem, it corresponds a redundant solution (gs, A). Therefore, we can reduce the minimization problem to the set SO(2)× Fix(D 2 ), rather than O(2) × Fix(D 2 ), where Fix(D 2 ) denotes the subspace of normal forms for at least orthotropic elasticity tensors. We shall thus consider the optimization problem

d(C raw , Σ ortho ) = min g,A C raw -ρ(g)A (4.1) over g = r θ ∈ SO(2) and A ∈ Fix(D 2 ) (with N = 4 independent components A ijkl ).
Solving problem (4.1) can be done numerically with the risk to reach a local minimum instead of the global minimum [START_REF] François | Determination of the symmetries of an experimentally determined stiffness tensor: Application to acoustic measurements[END_REF][START_REF] Kochetov | On obtaining effective orthotropic elasticity tensors[END_REF][START_REF] Danek | Effective elasticity tensors in context of random errors[END_REF][START_REF] Moakher | The closest elastic tensor of arbitrary symmetry to an elasticity tensor of lower symmetry[END_REF]. To avoid this difficulty in 3D, François and coworkers did propose pole figures for the given elasticity tensor C raw [START_REF] Francois | Détermination des symétries matérielles de matériaux anisotropes[END_REF][START_REF] François | Une nouvelle analyse des symétries d'un matériau élastique anisotrope. exemple d'utilisation à partir de mesures ultrasonores[END_REF] (renamed plots of the monoclinic distance in [START_REF] Kochetov | On obtaining effective orthotropic elasticity tensors[END_REF][START_REF] Moakher | The closest elastic tensor of arbitrary symmetry to an elasticity tensor of lower symmetry[END_REF]). Accordingly, they got an initial value for C = ρ(g)A, not too far from C raw , which was then optimized by a standard numerical (iterative) scheme. In the following, one will avoid numerical schemes and a closed form solution will be sought.

4.1.

Euler-Lagrange first order equations. Since the method is very general, we explain it below for a linear action ρ(g)T of the orthogonal group SO(d) (where d is 2 or 3) on a space T of tensors T. We assume that we know a normal form A for a given symmetry class (with N independent parameters). Moreover, we use a scalar product on tensors T, for which the action is isometric, meaning that

ρ(g)T 1 , ρ(g)T 2 = T 1 , T 2 , (4.2
) for all g ∈ SO(d) and all tensors T 1 , T 2 . The problem is now, given a raw tensor T raw , to calculate the critical points of the functional f (g, A) = T raw -ρ(g)A 2 = T raw -ρ(g)A, T raw -ρ(g)A .

(4.3)

To do so, we consider a path g(t) ∈ SO(d) with g(0) = g and ġ(0) = δg and a path A(t) of normal forms with A(0) = A and Ȧ(0) = δA, and we get

df (g, A).(δg, δA) = d dt t=0 f (g(t), A(t)) = -2 T raw -ρ(g)A, (δρ(g))A + ρ(g)δA . (4.4)
To calculate δρ(g), we use the infinitesimal action ρ ′ induced by ρ (see Appendix A) and defined by

ρ ′ (u)T = d dt t=0 ρ(h(t))T, (4.5)
where h(t) is a path in SO(d) with h(0) = I and ḣ(0) = u ∈ so(d) is an infinitesimal rotation.

Example 4.2. For an elasticity tensor T = C ∈ Ela(d), we get

(ρ ′ (u)C) ijkl = u ip C pjkl + u jp C ipkl + u kp C ijpl + u lp C ijkp . (4.6)
Remark 4.3. Due to the fact that ρ is isometric (4.2), we have

T 1 , ρ ′ (u)T 2 = -ρ ′ (u)T 1 , T 2 , (4.7)
for every infinitesimal rotation u ∈ so(d) and every tensors T 1 , T 2 . Now, introducing the infinitesimal rotations v = g -1 δg and u = δg g -1 , and, thanks to the fact that ρ is a linear representation, we have (see Appendix A)

δρ(g) = ρ(g)ρ ′ (v) = ρ ′ (u)ρ(g), where v = g -1 u g.
Therefore, the Euler-Lagrange first-order equations (deduced from (4.4)) df (g, A).(δg, δA) = 0, ∀δg, δA, recast as T raw -ρ(g)A, ρ(g)δA = 0, ∀δA, and T raw -ρ(g)A, ρ ′ (u)ρ(g)A = 0 ∀u.

Moreover, by (4.2), the first equation rewrites as

ρ(g t )T raw -A, δA = 0, ∀δA, (4.8)
and by (4.7), we have ρ(g)A, ρ ′ (u)ρ(g)A = -ρ ′ (u)ρ(g)A, ρ(g)A = 0 so that the second equation simplifies as

T raw , ρ ′ (u)ρ(g)A = 0, ∀u. (4.9)
Consider now a basis (u α ) of the space of skew symmetric matrices (of dimension d(d -1)/2) and an orthonormal basis (E I ) of the space of tensors, such that the first N vectors E 1 , . . . , E N span the subspace of normal forms A and set

A = N I=1 Λ I E I . (4.10)
Then, the first order equations (4.8)-(4.9) can be written directly in terms of the unknown variables, the rotation g = (g ij ) and the components Λ I of A, as the two sets of equations

Λ I =< ρ(g t )T raw , E I > 1 ≤ I ≤ N, (4.11)
< ρ(g t )T raw , ρ ′ (u α ) N I=1 Λ I E I >= 0 1 ≤ α ≤ d(d -1)/2. ( 4.12)
These equations are polynomial in the variables g ij and Λ I . Using the genericity of T raw and some arguments of algebraic geometry, one can prove that they have only a finite number of solutions (g k , A k ) when the symmetry group corresponding to the symmetry stratum under consideration is finite. These solutions are the critical points of the function f introduced in (4.3). Therefore, only a finite number of critical points (g k , A k ) needs to be compared. The global minimum of T raw -T with T of a given symmetry class Σ, will simply be

d(T raw , Σ) = min k T raw -T k , T k = ρ(g k )A k ,
and the closest tensor(s) to the considered symmetry class will be

T Σ = argmin k T raw -T k .
4.2. Distance to plane elasticity orthotropy. We now apply the preceding method to calculate the distance to the orthotropic stratum of a given plane elasticity tensor C raw . Three upper bounds estimates of this distance [START_REF] Gazis | The elastic tensor of given symmetry nearest to an anisotropic elastic tensor[END_REF][START_REF] Stahn | Distances of stiffnesses to symmetry classes[END_REF][START_REF] Stahn | Distances of stiffnesses to symmetry classes[END_REF] are compared in Appendix C. The Lie algebra so(2) (the space of infinitesimal 2D rotations) of the rotation group SO(2) is the one-dimensional vector space spanned by

u 1 = 1 √ 2 0 -1 1 0 .
It was shown in [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF] that C and C raw have the same isotropic part C iso defined by (3.2). Therefore, according to (2.4), the parametrization of C ortho can be written as

C ortho = C iso + T, T = ρ(g)(Λ 1 E 1 + Λ 2 E 2 ), g = cos θ -sin θ sin θ cos θ , (4.13)
where the orthonormal basis of tensors (E 1 , E 2 ) have Kelvin matrix representations

E 1 =    1 2 √ 2 -1 2 √ 2 0 -1 2 √ 2 1 2 √ 2 0 0 0 -1 √ 2    , E 2 =    1 √ 2 0 0 0 -1 √ 2 0 0 0 0    , with C iso , E 1 = C iso , E 2 = 0.
The later equalities can be checked using the Kelvin matrix representation of C iso ,

C iso =   2µ + λ λ 0 λ 2µ + λ 0 0 0 2µ   .
We are then looking for the critical points of the functional (4.3) with T raw = C raw -C iso given and where

A = Λ 1 E 1 + Λ 2 E 2 is bi-dimensional. Equation (4.11
) is used to determine the components Λ 1 and Λ 2 as functions of the angle of the rotation (2.1), giving

Λ 1 (θ) = X 1 cos 4θ + Y 1 sin 4θ, Λ 2 (θ) = X 2 cos 2θ + Y 2 sin 2θ, (4.14)
with

X 1 = 1 2 √ 2 (C raw 1111 -2C raw 1122 -4C raw 1212 + C raw 2222 ), Y 1 = √ 2(C raw 1112 -C raw 2212 ), X 2 = 1 √ 2 (C raw 1111 -C raw 2222 ), Y 2 = √ 2(C raw 1112 + C raw 2212
). Then, these results are injected in Equation (4.12), using (4.6), in order to determine θ. This gives A cos 8θ + B sin 8θ + C cos 4θ + D sin 4θ = 0, (4.15) with the following closed form expressions for the constants

A = 4X 1 Y 1 , B = 2(Y 2 1 -X 2 1 ), C = 2X 2 Y 2 , D = Y 2 2 -X 2 2 .
The sought tensors C k are finally given by (4.13) with θ k solution of (4.15). Remark 4.4. Equation (4.15) has several solutions θ (some corresponding to maxima, others to saddle points, and other to minima, possibly local). By remark 4.1, if (θ, A) is a solution, then (θ + π, A) is also a solution. Therefore, to obtain all the solutions of the Euler-Lagrange equations, it is sufficient to seek for θ in an interval of length π. Furthermore, if θ is a solution of (4.15), then θ -π 2 is also a solution of (4.15) (with identical Λ 1 but with opposite Λ 2 ). This allows us to seek for θ in an interval of length π/2, but to take account, for each solution

θ k ∈] -π 4 , π 4 [, of a twin solution θ -k = θ k -π 2 ∈] -3π 4 , -π 4 [. Check first if θ 1 = π/4
and θ -1 = -π/4 are solutions (an improbable event in case of an experimental, therefore noisy, elasticity tensor C raw ). Seek for the other solutions as θ = 1 2 arctan t, where t is a real root of the degree 4 polynomial (generically,

A -C = 0) (A -C)t 4 + (2D -4B)t 3 -6At 2 + (2D + 4B)t + A + C = 0,
deduced from (4.15) when θ = π 4 + n π 2 . Finally, choosing among the solutions C k , corresponding to solutions θ k and their twins θ -k = θ k -π 2 , the plane elasticity tensor

C ortho = C ℓ = C iso + ρ(g(θ ℓ )) Λ 1 (θ ℓ )E 1 + Λ 2 (θ ℓ )E 2 ,
which realizes the global minimum

d(C raw , Σ ortho ) = min k C raw -C k = C raw -C ortho ,
completes the solving.

Different equivalent definitions of the distance to a symmetry class

We have considered (and solved in a closed form for plane elasticity) the distance to a symmetry class problem formulated thanks to the normal form A / rotation g parameterization (Eq. (4.1)), as [START_REF] Vianello | An integrity basis for plane elasticity tensors[END_REF][START_REF] Dellinger | Kelvin notation for stabilizing elastic-constant inversion[END_REF][START_REF] Dellinger | Computing the optimal transversely isotropic approximation of a general elastic tensor[END_REF]

d(C raw , Σ [H] ) := min g∈SO(d), A∈Fix H C raw -ρ(g)A , C raw , A ∈ Ela(d), (5.1
R H (C) := 1 |H| h∈H ρ(h)C, (5.2)
where |H| is the order of the group H. The Reynolds operator is an orthogonal projector, meaning that, for all tensors C,

R H (R H (C)) = R H (C) , and C -R H (C) , R H (C) = 0. ( 5.3)
For the conjugate subgroup gHg -1 of H, we get

R gHg -1 (C) = 1 |H| h∈H ρ(ghg -1 )C = ρ(g)R H ρ(g -1 ) (C) . (5.4)
In particular, for plane elasticity and H = D 2 (of order |D 2 | = 4), we have

R D 2 (C) = 1 4 (C + ρ(r π )C + ρ(s)C + ρ(r π s)C) = 1 2 (C + ρ(s)C) , R gD 2 g -1 (C) = 1 2 C + ρ(gsg -1 )C ,
where the elements e, r π , s and r π s of the dihedral group D 2 have been defined in section 2.

Another definition of the distance to a symmetry class is due to François and coworkers [START_REF] Francois | Détermination des symétries matérielles de matériaux anisotropes[END_REF][START_REF] François | Determination of the symmetries of an experimentally determined stiffness tensor: Application to acoustic measurements[END_REF],

d(C raw , Σ [H] ) := min g∈SO(d) C raw -R gHg -1 (C raw ) . (5.5)
and a slightly modified definition of it has recently been used by Weber et al [START_REF] Weber | Distance of a stiffness tetrad to the symmetry classes of linear elasticity[END_REF],

d(C raw , Σ [H] ) = min g∈SO(d) ρ(g)C raw -R H ((ρ(g)C raw ) . (5.6)
Remark 5.2. The three distances to a symmetry class defined by (5.1), (5.5) and (5.6) are equal when the norm • is SO(d)-invariant (meaning that ρ(g)C = C for all rotations g and tensors C), which is the case for the Euclidean norm.

We provide now a proof of this fact. Consider first the definition (5.1),

d(C raw , Σ [H] ) = min g∈SO(d) min A∈Fix H C raw -ρ(g)A ,
because the variables g and A are independent. We have then min

A∈Fix H C raw -ρ(g)A = min A∈Fix H ρ(g) -1 C raw -A , 2
For infinite compact groups H such as O(2) or SO(2), the finite sum is replaced by an integral and the Haar measure, a (bi-invariant) probability measure on H [START_REF] Sternberg | Group theory and physics[END_REF], is used.

and hence, by (5.3), the minimizer A ∈ Fix H is the orthogonal projection

A = R H ρ(g) -1 C raw of ρ(g) -1
C raw on Fix H. Thus, using (5.4), we get min

A∈Fix H C raw -ρ(g)A = ρ(g) -1 C raw -R H ρ(g) -1 C raw = C raw -R gHg -1 (C raw ) ,
which proves the equivalence of definitions (5.1) and (5.5). Now, since ρ(g) preserves the norm of elasticity tensors, we have moreover

d(C raw , Σ [H] ) = min g∈SO(d) ρ(g -1 )C raw -R H ρ(g -1 )C raw = min g∈SO(d) ρ(g)C raw -R H (ρ(g)C raw ) ,
which shows the equivalence with definition (5.6).

Conclusion

The problem of calculating the distance of a raw tensor to a symmetry class has been posed. It has been fully solved for plane (bidimensional) elasticity tensors. First, Vianello's orthogonal projection method, valid for both the isotropic and the square symmetric (tetragonal) symmetry classes, has been recalled. Then, the remaining case of the distance to plane elasticity orthotropy has been solved, thanks to Euler-Lagrange method. The solution proposed is analytical, it requires only to find the roots of a degree four polynomial and to compare the at most eight closed-form solutions of the first-order Euler-Lagrange equations. The method is general and relies on the use of the infinitesimal action of Lie algebra of the rotation group SO(d) to solve the first-order equations. This use seems to be new in the present context. Another important feature of this method is that, dealing with polynomial functions, there is generically only a finite number of critical points, provided that the symmetry group associated to the symmetry class under consideration is finite. In that case, finding the global minimum is immediate. To illustrate further this method, an application to constitutive (3D) symmetric second-order tensors is provided in Appendix B.

Appendix A. Basic concepts in representation theory

A linear action of a group G on a vector space V (also usually called a linear representation of G on V ) is a mapping ρ : G → GL(V ), g → ρ(g). where ρ(g) is an invertible linear transformation of V and GL(V ) is the group of invertible, linear mappings of V into itself, and such that ρ(e) = Id, and ρ(g

1 g 2 ) = ρ(g 1 )ρ(g 2 ),
where g 1 , g 2 ∈ G and e is the unit element of G. The orbit of a vector v ∈ V is the set

Orb(v) := {ρ(g)v; g ∈ G} .
The symmetry group of a vector v ∈ V is the subgroup of G defined by

G v := {g ∈ G; ρ(g)v = v} ,
and its symmetry class, noted [G v ], is defined as the conjugacy class of G v in G, i.e., [G v ] := gG v g -1 ; g ∈ G .
Observe that all vectors in a same orbit Orb(v) have the same symmetry class, since

G ρ(g)v = gG v g -1 .
Finally, to each symmetry class [H], corresponds a symmetry stratum Σ [H] (which, in general, is not a linear subspace of V ) which is the set of all vectors v which have exactly the symmetry

[H] Σ [H] := {v ∈ V ; G v ∈ [H]} . Remark A.1. To each symmetry stratum Σ [H] corresponds a closed symmetry stratum Σ [H] ,
which is the set of vectors v ∈ V which have at least the symmetry [H].

Given a symmetry class [H] and choosing a representative H in this conjugacy class, we can build the fix point set

Fix(H) :

= {v ∈ V ; ρ(h)v = v; ∀h ∈ H} .
This set is a linear subspace of V and intersects each orbit Orb(v) when H ⊂ G v . The subspace Fix(H) defines a normal form for Σ [H] , if the representative H in [H] has been well chosen.

We will suppose now that G is a Lie group (which means that G is not only a group but also a differentiable manifold and that the group operations are smooth). In practice, all groups G considered are closed subgroup of the general linear group (and hence G can be considered as a matrix group). Then, we can assume that the representation ρ of G on V is smooth. In that case, it induces a linear mapping, called the infinitesimal action

ρ ′ : g → gl(V ), u → T e ρ.u,
where g = T e G, the tangent space at the identity element is called the Lie algebra of G. This infinitesimal action is sufficient to describe the linear tangent mapping at every point g since the relation ρ(g 1 g 2 ) = ρ(g 1 )ρ(g 2 ) leads to

T g ρ.δg = ρ(g)ρ ′ (g -1 δg) = ρ ′ (δg g -1 )ρ(g).

Appendix B. 3D symmetric second-order tensors

In this section, Euler-Lagrange method, described in section 4, is illustrated on the space T = S 2 (R 3 ) of three-dimensional symmetric second-order tensors. We consider a given (for example experimental) constitutive symmetric second-order tensor T raw = a raw in 3D (for example an anisotropic thermal expansion tensor or an anisotropic conductivity tensor), which is furthermore assumed to have three distinct eigenvalues (i.e. to be orthotropic).

The action of the rotation group SO(3) on

T = S 2 (R 3 ) = {a; a ji = a ij } ,
is written ρ(g)a = g a g t (g ∈ SO(3)), or in components (ρ(g)a) ij = g ip g jq a pq .

Note that, on even order tensors, the action of the full orthogonal group O(3) cannot be distinguished from the action of the 3D rotation group SO(3). The space of three-dimensional symmetric second-order tensors splits into three symmetry classes: orthotropy (three distinct eigenvalues), transverse isotropy (two distinct eigenvalues) and isotropy (one single eigenvalue). Since a generic symmetric second-order tensor T raw = a raw is orthotropic, a natural question then arise of how far it is from a transversely isotropic or an isotropic tensor ?

The distance to isotropy is obtained straightforwardly using the orthogonal projection of a raw onto the space of spherical tensors

a iso = 1 3 (tr a raw )1
and the distance of a raw to isotropy is

d(a raw , isotropy) = a raw -a iso = √ a raw ′ : a raw ′ .
We will thus focus on the distance of a raw to transverse isotropy, and seek by the Euler-Lagrange method for the closest transversely isotropic symmetric second-order tensor a to a raw (assumed orthotropic). The parametrization of (at least) transversely isotropic tensors is given by a = ρ(g)A = gAg t , where g ∈ SO(3) and

A = Λ 1 E 1 + Λ 2 E 2 , with E 1 = 1 √ 6   1 0 0 0 1 0 0 0 -2   , E 2 = 1 √ 3   1 0 0 0 1 0 0 0 1   .
This choice corresponds to an axis n = e 3 of transverse isotropy for A.

Remark B.1. Note that, contrary to the case of the distance to orthotropy for plane elasticity tensors, the symmetry group involved here for transverse isotropy is not finite. We expect thus to find an infinity of solutions (g, A) but they will lead, anyway, to a finite number of tensors a = gAg t , candidate to be a global minimum for the distance to the transversely isotropic stratum.

The problem is to determine the minimum of the functional

f (g, A) = a raw -gAg t 2 , A = Λ 1 E 1 + Λ 2 E 2 , g ∈ SO(3).
Without loss of generality, we can assume (after a diagonalization of a raw ) that the proper basis of a raw is the canonical basis of R 3 , (e i ), which will be done in the sequel.

The Lie algebra so(3) of the rotation group SO(3) is the vector space of 3 × 3 skew symmetric matrices (infinitesimal rotations), with as basis

u 1 = 1 √ 2   0 -1 0 1 0 0 0 0 0   , u 2 = 1 √ 2   0 0 -1 0 0 0 1 0 0   , u 3 = 1 √ 2   0 0 0 0 0 -1 0 1 0   .
To solve the problem, we calculate first the infinitesimal action ρ ′ on S 2 (R 3 ), which is written as (ρ ′ (u)a) ij = u ip a pj + u jp a ip , or, in a more intrinsic form, as

ρ ′ (u)a = ua -au = [u, a]. (B.1)
Then, we recast (4.9) using this expression, which leads to 0 = a raw , ρ ′ (u)ρ(g)A = tr a raw u(ρ(g)A) -(ρ(g)A)u = tr u ρ(g)A a rawa raw ρ(g)A

= u : [ρ(g)A, a raw ] = u : [a, a raw ] .
Since this last equality is true for every skew symmetric matrix u, it implies that the commutator [a, a raw ] vanishes and thus that the symmetric second-order tensors a and a raw commute. The given tensor a raw being orthotropic, it has three distinct eigenvalues and [a, a raw ] = 0 means that a is diagonal in the basis (e i ). Since a = gAg t , either both a and A are isotropic, in which case, all rotations g ∈ SO(3) are solutions or both a and A are transversely isotropic, in which case, each solution g of a = gAg t must send the transverse isotropy axis of A onto the one of a and writes thus g = g k h, where

g 1 = R(e 2 , π/2), g 2 = R(e 1 , -π/2), g 3 = 1,
R(n, θ) denotes the 3D rotation of angle θ around the vector n and h belongs to the subgroup of SO(3) of rotations which do not change the Oz axis. This subgroup of 3D rotations is isomorphic to the orthogonal group O(2). We have moreover

A = Λ 1 E 1 + Λ 2 E 2 ,
where, by (4.11), for each rotation g solution of a = gAg t , Λ 1 = (g t a raw g) :

E 1 , Λ 2 = 1 √ 3 tr a raw .
This achieves the determination of the critical points of f . Remark B.2. One recognizes Λ 2 E 2 = a iso as the isotropic part of a raw .

Isotropic solutions have already been calculated and correspond to the unique solution a iso = 1 3 (tr a raw )1, but which may not be a global minimum in the present situation (and, indeed, will not be one generically). Transversely isotropic solutions (g, A) correspond finally to the three following candidates for the global minimum a trans-iso

  λ 1 0 0 0 λ 2 +λ 3 2 0 0 0 λ 2 +λ 3 2   or   λ 1 +λ 3 2 0 0 0 λ 2 0 0 0 λ 1 +λ 3 2   or   λ 1 +λ 2 2 0 0 0 λ 1 +λ 2 2 0 0 0 λ 3  
where λ 1 < λ 2 < λ 3 are the distinct eigenvalues of a raw . The global minimum is obtained by comparing the three candidate distances a rawa associated with these critical points. We will summarize these results with the following conclusion. Let a raw = diag(λ 1 , λ 2 , λ 3 ) ∈ S 2 (R3 ) be a given symmetric second-order tensor, with λ 1 < λ 2 < λ 3 . Then, the transversely isotropic second-order tensor a trans-iso ∈ S 2 (R 3 ) closest to a raw (for the Frobenius norm) commutes with a raw and writes

a trans-iso =   λ 1 +λ 2 2 0 0 0 λ 1 +λ 2 2 0 0 0 λ 3   if λ 2 < 1 2 (λ 1 + λ 3 ), a trans-iso =   λ 1 0 0 0 λ 2 +λ 3 2 0 0 0 λ 2 +λ 3 2   if λ 2 > 1 2 (λ 1 + λ 3 ).
in the proper basis of a raw . The distance of a raw to the transversely isotropic class is then

d(a raw , transverse isotropy) = 1 √ 2 min i =j |λ i -λ j | .
Appendix C. Upper bounds estimates of the distance to a symmetry class Upper bounds estimates of the distance d(C raw , Σ [H] ) to a symmetry class have been proposed in the literature, for instance in [START_REF] Gazis | The elastic tensor of given symmetry nearest to an anisotropic elastic tensor[END_REF] and [START_REF] Stahn | Distances of stiffnesses to symmetry classes[END_REF] for 3D elasticity symmetry classes, and in [START_REF] Oliver-Leblond | Continuous anisotropic damage as a twin modelling of discrete bi-dimensional fracture[END_REF] for 2D orthotropic symmetry class. They are easier to calculate than the distance itself.

First, Gazis and coworkers did simply suggest the estimate

M g 0 (C raw , Σ [H] ) := ρ(g 0 )C raw -R H (ρ(g 0 )C raw ) , (C.1)
with g 0 a "well-chosen" rotation 3 , for example related to the manufacturing process or to the observed microstructure. Obviously,

M g 0 (C raw , Σ [H] ) ≥ d(C raw , Σ [H] ), (C.2)
since the minimization in (5.6) is not performed when calculating

M g 0 (C raw , Σ [H] ).
Studying the 3D elasticity problem, Stahn and coworkers astutely defined in [START_REF] Stahn | Distances of stiffnesses to symmetry classes[END_REF] a smaller upper bound M (C raw , Σ [H] ) of the distance to a 3D elasticity symmetry class, in a few steps,

(1) by computing

t = C asym : 1 - 1 4 (1 : C asym : 1) 1,
which is a second-order covariant of the asymmetric part C asym = C raw -(C raw ) s of the given elasticity tensor C raw , (2) following [START_REF] Baerheim | Harmonic decomposition of the anisotropic elasticity tensor[END_REF], by computing an eigenbasis of t and a rotation g 0 that brings it into its diagonal form (i.e., such that ρ(g 0 )t = g 0 t g t 0 is diagonal),

(3) by computing ρ(g 0 )C raw (and assuming in fact that ρ(g 0 )C is close to the normal form of the sought tensor with higher symmetry), (4) finally, by setting

M (C raw , Σ [H] ) := min g∈O ρ(g)ρ(g 0 )C raw -R H (ρ(g)ρ(g 0 )C raw ) (C.3)
with R H the Reynolds operator (5.2). The procedure introduces the proper cubic (octahedral) group O which is defined by

O = {g ∈ SO(3); ge i = ±e j } ,
where (e i ) is the canonical orthonormal basis R 3 . This group is of order |O| = 24. Generically, the second-order tensor t is orthotropic, and, in 3D, there are exactly 24 rotations which bring t into a diagonal form. Given such a rotation, say g 0 , the other ones are written as gg 0 with g ∈ O, so that definition (C.3) does not depend on a particular choice for g 0 .

Remark C.1. As the min in (C.3) is over the finite group O (and not over the rotation group SO(3)), one has

M (C raw , Σ [H] ) = min g∈O ρ(gg 0 )C raw -R H (ρ(gg 0 )C raw ) ≥ min g∈SO(3) ρ(g)C raw -R H (ρ(g)C raw ) = d(C raw , Σ [H] ),
as g = gg 0 is a rotation, and where the min over the subgroup O ⊂ SO(3) in the first line is generally strictly larger than the min over SO(3) in the last line. Hence,

M (C raw , Σ [H] ) is an upper bound of the distance d(C raw , Σ [H] ).
Remark C.2. In 2D, the asymmetric part C asym of an elasticity tensor is isotropic (see remark 2.2), so that the Stahn and coworkers' procedure applies only if the definition of the second-order tensor t is changed, for example into the second-order harmonic component h of C raw defined in section 2. where H is furthermore the fourth-order harmonic part of C raw = (λ, µ, h, H) (determined by (2.9)), and Remark C.3. In the particular case of orthotropic plane elasticity, using definition (C.2) for t, we have the equality of the upper bounds estimates proposed in [START_REF] Oliver-Leblond | Continuous anisotropic damage as a twin modelling of discrete bi-dimensional fracture[END_REF] and in [START_REF] Stahn | Distances of stiffnesses to symmetry classes[END_REF],

C upper := 2µ I + λ1 ⊗ 1 + 1 2 (1 ⊗ h + h ⊗ 1) + 2K 3 I 2 2 h ⊙ h - 1 4 h 2 1 ⊙ 1 ,
∆(C raw , Σ ortho ) = C raw -C upper = J 2 I 2 2 -2K 2 3 I 2 = M (C raw , Σ [H] ).
A proof of this result is as follows. The harmonic decomposition of the raw plane elasticity tensor is C raw = (λ, µ, h, H) (by formulas (2.6)-(2.9)). Let g 0 be a rotation such that h 0 := ρ(g 0 )h = g 0 h g t 0 is diagonal. In 2D, there are only two such rotations, g 0 and g 1 = r π 2 g 0 , the second one giving ρ(g 1 )h = -h 0 (as h is deviatoric). Since the harmonic decomposition is equivariant, we get ρ(g 0 )C raw = (λ, µ, h 0 , H 0 ).

ρ(r π 2 ) -1 R D 2 ρ(r π 2 ) = R r π 2 D 2 r -1 π 2 = R D 2 .
Using the harmonic decomposition of ρ(g 0 )C raw , the upper bound (C.4) is

∆(C raw , Σ ortho ) = ∆(ρ(g 0 )C raw , Σ ortho ) = H 0 -2 h 0 : H 0 : h 0 h 0 4 h 0 ⊙ h 0 - 1 4 h 0 2 1 ⊙ 1
and we are left to show that 2 h 0 :

H 0 : h 0 h 0 4 h 0 ⊙ h 0 - 1 4 h 0 2 1 ⊙ 1 = R D 2 (H 0 ) ,
for all h 0 = 0 deviatoric and diagonal, and all fourth-order harmonic tensor H 0 , which can be checked by comparing these expressions in components. Indeed, a general harmonic tensor K, has for matrix components

K =   K 1111 -K 1111 √ 2K 1112 -K 1111 K 1111 - √ 2K 1112 √ 2K 1112 - √ 2K 1112 -2K 1111   . Now set T = 2 h 0 : K : h 0 h 0 4 h 0 ⊙ h 0 - 1 4 h 0 2 1 ⊙ 1 ,
where h 0 = λ 0 0 -λ , with λ = 0 (otherwise, the elasticity tensor C raw has already the square symmetry). We get then

T =   K 1111 -K 1111 0 -K 1111 K 1111 0 0 0 -2K 1111   ,
which is indeed the Kelvin matrix form of R D 2 (K).

We have therefore shown that for a biclinic elasticity tensor C raw ∈ Ela(2)

∆(C raw , Σ ortho ) = M (C raw , Σ ortho ) = J 2 I 2 2 -2K 2 3 I 2 ≥ d(C raw , Σ ortho ).
Funding. Two of the authors, R. Desmorat and B. Kolev, were partially supported by CNRS Projet 80-Prime GAMM (Géométrie algébrique complexe/réelle et mécanique des matériaux).

  3,t = h = (tr 12 C raw ) ′ = (1 : C raw ) ′ . M (C raw , Σ [H] ) built from t = h is then, by remark C.1, an upper bound estimate of the distance d(C raw , Σ [H] ).Finally, an invariant formula has been obtained in 2D for an upper bound estimate of the distance to orthotropic elasticity[START_REF] Oliver-Leblond | Continuous anisotropic damage as a twin modelling of discrete bi-dimensional fracture[END_REF] Corollary 3.2]. Indeed, settingh := (tr 12 C raw ) ′ , I 2 := h 2 = h : h, J 2 := H 2 = H :: H, K 3 := h : H : h,

2

 2 then the positive invariant∆(C raw , Σ ortho ) := C raw -C upper =J bound of the distance d(C raw , Σ ortho ) of the plane elasticity tensor C raw to orthotropy.

where H 0 2 =

 02 := ρ(g 0 )H. The Stahn and coworkers upper bound is then (by definition (C.3))M (C raw , Σ ortho ) = min ρ(g 0 )C raw -R D 2 (ρ(g 0 )C raw ) , ρ(r π 2 )ρ(g 0 )C raw -R D 2 ρ(r π 2 )ρ(g 0 )C raw = min ρ(g 0 )C raw -R D 2 (ρ(g 0 )C raw ) , ρ(g 0 )C raw -ρ(r π 2 ) -1 R D 2 ρ(r π 2 ) (ρ(g 0 )C raw ) = ρ(g 0 )C raw -R D 2 (ρ(g 0 )C raw ) = H 0 -R D 2 (H 0 ) , since r π 2 D 2 r -1π D 2 in 2D and thus

  .[START_REF] Arts | General anisotropic elastic tensors in rocks: Approximation, invariants, and particular directions[END_REF] This normal form depends on four independent parameters if it is orthotropic, three parameters if it has square symmetry (A 2222 = A 1111 ), and two parameters if it is isotropic (A 222 = A 1111 and 2A 1212 = A 1111 -A 1122 ). We denote by N this number of independent parameters.

	Remark 2.1. Observe that in (2.4), the tensor A is fixed by the dihedral group D 2 := {e, r π , s, r π s}. Hence, the space of tensors A having Kelvin matrix representation (2.4), is noted Fix(D 2 ). A
	systematic way to calculate such a normal form for a tensor in a given symmetry class is described
	in Appendix A.

  )with Σ [H] = Σ [D 2 ] = Σ ortho the orthotropic stratum.Remark 5.1. Here, the notation [H] stands for the conjugacy class of H, that is the set of all subgroups gHg -1 , where g ∈ SO(d) (see Appendix A). Indeed, a symmetry class is described by all the symmetry groups gHg -1 rather than the subgroup H itself. For instance, the orthotropy is described by the conjugacy class[D 2 ], where the dihedral group H = D 2 defined in remark 2.1 is just a representative in this class. The set Σ [H] correspond to all tensors having the same symmetry class [H] and for this reason Σ [gHg -1 ] = Σ [H] . Alternative definitions of d(C raw , Σ [H] ) introduce the Reynolds operator R H : Ela(d) → Fix(H) associated here with the finite group 2 H [34, Chapter 2],

A symmetry stratum is the set of all tensors having the same symmetry class.

in fact these authors did take g0 = e the identity in their example.

(Adrien Antonelli) Université Paris-Saclay, ENS Paris-Saclay, CNRS, LMT -Laboratoire de Mécanique et Technologie, 91190, Gif-sur-Yvette, France.

Email address: adrien.antonelli@ens-paris-saclay.fr Email address: rodrigue.desmorat@ens-paris-saclay.fr


Gerhard Münz

Ali Fessi

Georg Carle

Oliver Paul
email: olivier.paul@int-evry.fr

Dušan Gabrijelčič

Jozef Stefan

Yannick Carlinet
email: yannick.carlinet@francetelecom.com>francetelecom

Sherif Yusuf

Vrizlynn Morris Sloman

Jan Van Lunteren

Patricia Sagmeister

Gero Dittmann

DIADEM Firewall: Web Server Overload Attack Detection and Response

High-profile web servers often become the victim of web server overload Distributed Denial-of-Service (DDoS) attacks. Motivations of such attacks range from technical challenge (e.g. script kiddies) to financial profit (e.g. blackmailing the web server's owner).

This paper presents the DIADEM Firewall architecture that allows an ISP to protect its customers from being the target of a DDoS attack. Additionally, it provides protection against usage of customer hosts for attacks. Furthermore, the use-case of the web server overload attack detection and response mechanism will be explained in more details. Finally, we outline the integration an FPGA based highspeed classifier engine integrated into the Linux Netfilter firewall as well as its deployment during a response action against the DDoS attack.

Introduction

Today's networking environments are becoming more and more complex. The number and diversity of deployed network devices and configurations is rising and the diversity of services offered and used is increasing. On the other hand, broadband Internet access is becoming a commodity for more and more Internet users.

However, networking environments have not only increased the utility provided for the users, they have also become prone to various kinds of malicious activities and security threats. Attacks and the tools used to perform them are getting more and more sophisticated and harmful by exploiting the increasing number of permanently connected and often badly secured customer hosts. Such hosts are susceptible to worm and virus infections, and once infected, they can be easily misused as attack sources in Distributed Denial-of-Service (DDoS) attacks. Firewalls as well as Intrusion Detection and Prevention Systems (IDS/IPS) are used to detect such malicious activities and to protect the networking environment and other potential victims, such as web servers, corporate networks etc. However, these systems are usually operated as standalone devices protecting an isolated subsystem or network only. Possibilities for coordinated operation are very limited and mostly realized on a manual basis. Hence, the existing attack detection and prevention systems are not sufficiently effective as countermeasures against sophisticated attacks and upcoming threats.

In the scope of the DIADEM 1 Firewall project [START_REF]DIADEM Firewall Homepage[END_REF], we are developing an architecture that enables an ISP (Internet Service Provider) to better protect its own networking environment as well as the connected hosts and servers of its customers. The main goal of the project is to establish a tight loop of detection, decision and response across multiple, distributed firewall and networking devices.

In this paper, we present our approach on the basis of one selected use-case, the web server overload DDoS attack. During such an attack, multiple distributed attack sources attempt to overload a web server by simultaneously requesting a given object located on the server. Depending on the number of attackers and the performance of the web server, the attack causes service degradation or even Denial-of-Service for legitimate users. This attack is common today and has various properties we can expect from new and forthcoming attacks.

The rest of this paper is structured as follow. An overview of the overall DIADEM Firewall architecture is presented in section 2. The general attack detection architecture and particularly the web server overload attack detection mechanism are explained in section 3, followed by a description of the response architecture and the web server use-case related response actions in section 4. Section 5 outlines how a high-speed packet classifier engine implemented on an FPGA is deployed for efficient mitigation of the attack. Finally, conclusions and directions for future work are given in section 6.

Overall Architecture of DIADEM Firewall

This section gives an overview of the architecture of DIADEM Firewall and presents its high-level components. As shown in Figure 1, we separate the different components of the architecture according to three levels: data level, element level and administrative level.

The Data Level

The data level contains the equipment that deals with the actual traffic in the operator's network. Considered devices include routers, network monitors, and firewalls that are used for network monitoring and response purposes. Project results include implementations of innovative components with special monitoring and response functionality that are evaluated in scenarios with commercial and open-source (Linux) routers and firewalls.

The Element Level

The element level is an abstraction layer for the data level equipment, which hides device specific details by providing abstract interfaces (APIs) to control the equipment in the data level independently of its specific control interface. Two types of components are considered: Monitoring Elements and Firewall Elements.

Monitoring Element (ME)

The Monitoring Element is an abstraction of a router or a metering device that monitors traffic at a specific point in the network and exports both flow and packet-based monitoring data to a remote system called Violation Detection (c.f. section 2.3) for further analysis. Flow-based data consists of statistics such as the number of packets and the number of octets observed during the reporting interval for a specific source-destination flow. Packet-based information includes packets that are selected with filtering and packet sampling methods. The Monitoring Element can also provide statistics about packets with specific header fields, for instance, the number of observed TCP SYN packets or ICMP echo request packets within a given time interval. The Monitoring Element can also be configured to focus on application specific data such as HTTP traffic to detect web server attacks. Section 3 explains the corresponding monitoring and detection process in more detail.

The Monitoring Element provides a Monitoring API to the administrative level for dynamically setting and changing the monitoring configuration. Thus, the network monitoring process and exported monitoring data can be adapted to the current needs of the Violation Detection system.

Monitored data is exported from the Monitoring Element to the Violation Detection system using the IPFIX (IP Flow Information Export) protocol [START_REF] Claise | IPFIX Protocol Specifications[END_REF] which is currently standardized in the IETF. The strength of the IPFIX protocol resides in its capability to carry any kind of data records in a bandwidth-efficient way.

The Monitoring API is realized using the Netconf protocol [START_REF] Enns | NETCONF Configuration Protocol[END_REF] which is also currently standardized in the IETF. The Netconf protocol is very flexible and can be easily adapted to different configuration tasks.

The Firewall Element (FE)

The Firewall Element is an abstraction for a router or a firewall device that allows enforcing response policies in order to stop or mitigate an attack. Response policies may require that one or more specific data flows are blocked, redirected or rate-limited. The Firewall Element offers a Firewall API and a Response API to the administrative level. Hiding the device-specific details of the data level components, the Firewall API can be used to trigger immediate actions. On the other hand, the Response API delivers response policies to the Firewall Element that allow for making fast local decisions.

The Administrative Level

The administrative level allows the system administrator of to enter new policies, new attack counter-measure logic, and take decisions when needed for a counter-measure. It comprises two components: Violation Detection and System Manager.

Violation Detection (VD)

The Violation Detection system receives the monitoring data from Monitoring Elements which are distributed in the network. It then performs a variety of attack detection methods deploying both signature-based and anomaly detection techniques. If an attack is detected, information about the type of attack, the target, recognized attack sources etc. is encoded into IDEMF (Intrusion Detection Exchange Message Format) [START_REF] Debar | The Intrusion Detection Message Exchange Format[END_REF] messages and sent as an alert notification to the System Manager.

The Violation Detection is described in more details in Section 3.

System Manager (SM)

The System Manager integrates a policy service for specifying, storing and distributing policies to Policy Management Agents (PMAs). PMAs can be deployed in the Violation Detection, the Firewall Elements, and in the System Manager itself.

PMAs receive events via a publish-subscribe system. Upon reception of a new event, a PMA searches for a matching policy and enforces the corresponding action(s) if necessary. For example, the PMA of the System Manager might receive an attack notification event from the Violation Detection. A policy specifies how to react, e.g. by calling the Firewall API in order to perform a specific response action on the Firewall Element, or by triggering a reconfiguration of the Violation Detect to adapt its detection strategy.

The functionality of the System Manager is explained in more details in Section 4.

Deployment

Figure 2 sketches the deployment of the DIADEM Firewall architecture in an operator's network. The network is equipped with several highly distributed Monitoring Elements that deliver monitoring data to the Violation Detection which analyzes the data in order to detect ongoing attacks. If an attack is detected or suspected by the Violation Detection, a notification message is sent to the System Manager which uses this information to determine and trigger the appropriate response to the attack. A typical response consists of blocking, redirecting, or rate-limiting the suspicious traffic. This requires the reconfiguration of the appropriate Firewall Elements.

Ideally, if the source(s) of the attack can be correctly identified, the attack traffic can be blocked as close as possible to the source generating the attack. This prevents attack traffic from reaching the victim and reduces unwanted traffic within the network.

Like the Monitoring and Firewall Elements in the element level, the components of the administrative level (Violation Detection and System Manager) can also be distributed for performance reasons or to avoid single points of failures.

Inter-Domain Deployment

Cooperation between different administrative domains is advantageous for network monitoring as well as for attack detection and response. However, a strong trust relationship between cooperating network providers is required. Handing out monitoring data would also reveal the network topology and information about network performance and utilization etc. Furthermore, customer's privacy comes into question here. Less difficult is the dissemination of alert notifications between distinct administrative domains. Upon receiving an alert from another domain, the System Manager can decide to enforce appropriate countermeasures immediately, or to ask the local Violation Detection at first if the attack can also be observed in the own domain. Depending on the trust relationship with the foreign domain, the System Manager might also decide to ignore a foreign alert completely. Note that a System Manager will not usually be permitted to enforce response actions or policies on Firewall Elements belonging to a different domain.

Violation Detection

This section describes the Violation Detection architecture within DIADEM Firewall. In particular, a detection module for web server overloading is presented.

Architecture

As represented in Figure 4, the Violation Detection system is organized in several modules: -The IPFIX collector receives monitoring data from Monitoring Elements. -Detection modules apply detection methods on the monitoring data. Detection results are then passed to the PMA as IDEMF events. Existing detection tools (e.g. Snort [START_REF]Snort Homepage[END_REF]) may be integrated in the violation detection function as detection modules. -The Policy Management Agent (PMA) stores and interprets detection policies which may be triggered by IDMEF events or can be used to filter events. Possible actions are notify, log, store, defer the event to an investigation module, or discard the event. Other actions might be to reconfigure the Monitoring Elements to get more information about current traffic. -The event database is used to store interesting IDMEF events for later investigation and correlation with other events. -The investigation modules examine events that cannot be clearly classified by detection modules. Possible investigation methods include: o compare/correlate an event with earlier events stored in the event database o reconfigure detection modules and/or Monitoring Elements If the investigation module concludes that an event is an attack, a new IDMEF event is passed to the PMA.

Web Server Overload Detection Module

The overall behavior of the Web Server Overload Detection Module is to build a view of application level actions from the network-level measurement parameters provided by the IPFIX collector. These inferred application level actions are used to detect overloading attacks by matching detected results with a detection model. When an attack is detected, an IDMEF event is sent to the PMA within the Violation Detection system. Each Web Overload Detection Module is expected to monitor traffic directed to and originating from a specific web server. As a result, models and training data used by each module are server specific.

The Web Server Overload Detection Module consists of four components as shown in Figure 3.

-The parameters extraction component receives IPFIX records from the IPFIX collector, selects those that are relevant for the monitored targets, and extracts the useful measurement parameters. Monitored targets include potential victims (web server address and port number) as well as potential attackers (e.g. IP addresses). The definition of monitored targets is included in the target database. -The parameter estimation components transform the extracted measurement parameters into input parameters that are adapted to the inference engine. For example, IP addresses may be transformed into country codes in order to limit the size of the inference model. Each parameter may need a different form of correction or adaptation. For each parameter, the estimation function looks up the corresponding model from the estimation database which defines the necessary transformation depending on the type of parameter as well as the target web server. In the case of the IP address to country code transformation, this model would be a data structure representing relationships between these two parameters.

-The inference engine infers information about the invoked HTTP operation (HTTP method, targeted URI, operation result, etc.) from the parameters received by the estimation components. In contrast to classical application level processing schemes, the inference engine retrieves application level information based on probabilistic models applied to the measurement results from lower levels. Thus, the inference approach does not require expensive packet processing including reassembly and pattern matching and allows for a faster, yet less precise application level parameters estimation [START_REF] Paul | Tradeoffs for Web Communications Fast Analysis[END_REF]. For each estimated element, the inference engine selects the corresponding inference model from the inference database. This model derives a tuple including the identity of the user (IP address), the performed action (HTTP header values), and the target web server from the estimated input parameters.

Although inference models are target specific, we separate the definition of targets from the definition of models for efficiency reasons. The definition of targets may change frequently (several times an hour) while the modification of models should happen less frequently (for example, once a month or more often depending on the frequency of the server's content update).

The detection engine checks if the tuple received from the inference engine conforms to normal target behavior. Therefore, profiles of normal behavior are stored in the detection database. The detection engine retrieves the appropriate profile from the database and compares it to the received tuple. A profile represents the normal HTTP server behavior as well as the acceptable deviation. Profiles are built incrementally using received tuples and a change detection algorithm. The profile is used, along with the identity of the source, to decide whether an overload attack is going on. If an attack is detected, an alert is sent to the notification component which sends a corresponding IDMEF event to the PMA. The IDMEF event includes the identity of the victim (address and port), the identity of the attack sources (IP addresses and port numbers), the type of attack (used method, impacted URI, operation result etc.), the level of deviation between the attacker's behavior and a normal behavior, and the number of packets and bytes related to the communication.

Attack Response Mechanisms

This section describes the System Manager architecture and how it communicates with other modules in the DIADEM Firewall architecture in Figure 1.

The System Manager

The System Manager mainly communicates with the Violation Detection system and the Firewall Elements, as indicated in Figure 5.

The System Manager receives IDMEF notifications from the Violation Detection system indicating the type of attack, the target and the presumed source(s) of the attack, the signature of the attack etc. The System Manager uses this information to determine the appropriate countermeasures against the attack. Notifications may result in the reconfiguration of the Violation Detection or the Firewall Elements in order to respond to the detected attack. The appropriate response actions are specified by response policies enabled by the System Manager. When notifications related to a given attack are no longer received over a given period of time, the System Manager stops the enforced response actions and restores the prior state of the system.

The System Manager provides reports to the DIADEM Firewall architecture administrator as well as to the impacted customer regarding the detection and the treatment of the attack.

Response Policies

Response policies specify the actions needed to be undertaken in order to stop or mitigate the attack, making use of any of the parameters provided by the triggering alert notification. In case of an automated response, they equate to one of the following:

• perform an action on an Firewall Element, such as redirection, rate limiting etc., • perform an action on the Violation Detection component; e.g. reconfigure the Violation Detection for further investigation of the attack, or • generate and send a new event to the Violation Detection system or to a Firewall Element. In addition, human network administrators can be informed if manual intervention is required.

Policy Management Agents

Automated policy enforcement is entrusted to Policy Management Agents (PMAs) that act as local policy interpreters in the Violation Detection system and the Firewall Elements. Once policies are compiled, they are loaded and enabled in the PMAs. The Elvin event system [START_REF]Content Based Messaging[END_REF] is used for communication with and between the PMAs. PMAs only subscribe to those events which are relevant to their stored policies. On receiving an event, the PMA queries a domain service to determine the target objects on which to perform the policy action(s). A policy can optionally specify constraints which must evaluate to true before the action can be performed.

Web Server Overload Response

When the Violation Detection system detects a web server attack, it notifies the System Manager, which inspects the attack signature and may decide to reconfigure the Violation Detection to improve the detection result. For example, the Violation Detection system can be asked to look for a more specific attack signature or to increase the confidence in the detection result in order to reduce the probability of false positives. As a consequence, it might be necessary that the Violation Detection system also reconfigures the associated Monitoring Elements to provide more detailed monitoring data.

If the alert notification is severe enough, the corresponding response policy may trigger a response mechanism on the Firewall Elements. The web server overload attack requires established TCP connections between the attackers and the victim. Therefore, the sources of the attack can be identified easily. An effective response consists in: -setting new firewall rules that drop all packets originated by the attacking sources and directed to the web server and closing the TCP connection by inserting a TCP FIN or RST packet towards the server. Note here, that the inserted packet must carry the right sequence number.

While the firewall rules prevent the attackers from sending more HTTP requests, closing the existing TCP connections immediately releases the web server by stopping the currently served requests. Both functionalities are provided by a Linux based firewall implementation. In order to speed up the packet filtering, we make use of a packet classifier engine as presented in section 5. Figure 6 shows an example of the policy that triggers the closing of a TCP connection under the condition that the confidence in the detection result exceeds a predefined minimum confidence level.

Whether a response to an attack is successful or not, does not only depend on the functionality of the firewall devices but also on the points in the network where the Finally, response actions should be stopped when the attack has finished. Though, it is often not possible to determine if the attack has finished as long as the response actions are enforced. In case of the web server overload attack, the firewall rules are removed after a certain period of time. Incoming traffic from the attacker will be analyzed further for a faster detection and response in case the attack persists.

High-Speed Packet Classification

As described in section 2.2, the Firewall Element is an abstraction for a router or a firewall device that allows enforcing response policies in order to stop or mitigate an attack. Such a firewall device can be a software firewall such as Linux Netfilter or a commercial hardware firewall.

In the scope of the DIADEM Firewall project, a highspeed packet Classifier Engine (CE) is developed and integrated into a Netfilter-based software firewall in order to offload the processing-intensive packet classification function from the firewall system into hardware for very high-speed links. The CE is based on the P2C algorithm [START_REF] Van Lunteren | Fast and Scalable Packet Classification[END_REF] which is a multi-field packet classification scheme that searches the fields independently in a first phase and determines the final classification result using a TCAM in a second phase. The P2C concept is shown in Figure 7 for 5tuple classification.

A key element of the P2C scheme is novel encoding of the intermediate search results, which allows a significant reduction of the storage requirements and minimizes the dependencies within the field search structures, thus enabling fast incremental updates. The encoding involves several styles that can be applied simultaneously and allow several performance aspects to be balanced at the granularity of individual rules. The scheme can be used in combination with various search algorithms and memory technologies. The P2C implementation for the CE performs the field searches using the Balanced Routing Table search (BaRT) algorithm [START_REF] Van Lunteren | Searching Very Large Routing Tables in Wide Embedded Memory[END_REF]. BaRT is a scheme for exact-and prefix-match searches that achieves wire-speed search performance in combination with high storage-efficiency and fast dynamic updates, through application of a novel hash function. BaRT can efficiently process the search key in segments of about 8 bits, requiring in the worst-case only four memory accesses to search a 32-bit IP address and two accesses for a 16-bit port number, all of which can be performed efficiently in a pipelined fashion. A detailed description of BaRT is presented in [START_REF] Van Lunteren | Searching Very Large Routing Tables in Wide Embedded Memory[END_REF].

FPGA prototypes based on the combination of P2C and BaRT have demonstrated the feasibility of achieving wirespeed classification performance for OC-192 and OC-768 with off-the-shelf technology. The prototype supported 1733 rules of a commercial firewall using only 2 KB of SRAM in combination with 5 KB of TCAM.

The CE will be implemented on an FPGA board and connected to a Linux Netfilter firewall through a PCI bus. It will be configured and initialized through the upper layer software of the Firewall Element. The lookup tables, which are a set of packet filter rules, are embedded in SRAM. The multi-field search, currently performed through a TCAM can also be implemented using TCAM emulation in SRAM. The communication between the Firewall Element and the CE is realized through an API located in user space and a loadable kernel module embedded in the Linux kernel. The basic architecture can be seen in Figure 8.

Conclusions

In this paper, we discussed DDoS attack detection and response, in particular the web server overload DDoS attack. To date, no efficient and cost-effective solution for the detection and mitigation of DDoS attacks is available. This is mostly for the reason that DDoS attacks can be detected at the victim, but effective countermeasures need to be enforced closer to the attack sources. Thus, protection systems deployed at a single point in the network cannot cope with DDoS attacks efficiently.

The DIADEM Firewall architecture is a significant improvement as it provides distributed detection and mitigation of DDoS attacks within the network. We did present how to deploy the DIADEM Firewall architecture as a protection against web server overload attacks, which present a common threat to high-profile web servers. Finally, we showed how response actions are accelerated using a high-speed packet classifier engine.

Next steps of the project include experiments in a test environment that resembles an operational broadband network. Furthermore, the DIADEM Firewall architecture will be used to implement efficient detection and response mechanisms against other highly distributed DoS attacks, such as TCP SYN flood attacks, profiting from the distributed nature of the DIADEM Firewall architecture.

Figure 1 :

 1 Figure 1: DIADEM Firewall Components

Figure 2 :Figure 3 :

 23 Figure 2: Deployment of the DIADEM Firewall Architecture

Figure 4 :

 4 Figure 4: Web Server Overloading Detection Module

Figure 7 :

 7 Figure 8: Linux Netfilter and Classifier Engine

short for DIstributed Approach for a high-speed programmable firewall and attack DEtection systeM.

Acknowledgments

We gratefully acknowledge support from EU FP6 for DIADEM Firewall Project FP6 IST-2002-002154.


Sacha Bensoussan

Kévin Premat

Eimad Shotar

Évelyne Cormier

Amira Al Raasi

Jean-Philippe Spano

Laetitia Morardet

Raphaël Bonaccorsi

Véronique Morel

Jugurtha Mathout

Jacques Chiras

Frédéric Clarençon
email: frederic.clarencon@aphp.fr

Percutaneous Reinforced Cementoplasty Using Spindles as a Palliative Option for Malignant Fractures of the Humerus

Keywords: Percutaneous cementoplasty, Augmented osteoplasty, Pathological fracture, Humerus, Polymethyl Methacrylate CT: Computed tomography IV: Intravenous PMMA: Polymethyl methacrylate RC: Reinforced cementoplasty VAS: Visual analogue scale 2

The humerus is the second most frequent long bone site of metastatic disease [START_REF] Frassica | Evaluation and Treatment of Metastases to the Humerus[END_REF]. Percutaneous cementoplasty is a minimally invasive alternative for patients with malignant humeral fracture unsuitable for surgery, but polymethyl methacrylate (PMMA) resists poorly to torsional forces and shear stress and is considered less effective for treating long bone tumors, with an increased risk of secondary fractures [START_REF] Garnon | Percutaneous extra-spinal cementoplasty in patients with cancer: A systematic review of procedural details and clinical outcomes[END_REF][START_REF] Cornelis | Augmented osteoplasty for proximal femur consolidation in cancer patients: Biomechanical considerations and techniques[END_REF]. Various strategies have been proposed to increase the resistance of PMMA in humeral fractures, but these studies are limited to case-reports [START_REF] Sun | Percutaneous cementoplasty for painful osteolytic humeral metastases: initial experience with an innovative technique[END_REF][START_REF] Kawai | Percutaneous Osteoplasty with Use of a Cement-filled Catheter for a Pathologic Fracture of the Humerus[END_REF][START_REF] Kelekis | Percutaneous Augmented Osteoplasty of the Humeral Bone Using a Combination of MicroNeedles Mesh and Cement[END_REF][START_REF] Anselmetti | Painful Pathologic Fracture of the Humerus: Percutaneous Osteoplasty With Bone Marrow Nails Under Hybrid Computed Tomography and Fluoroscopic Guidance[END_REF]. The present study describes a procedure of percutaneous reinforced cementoplasty (RC) for humeral pathological fractures using intralesional spindles based on a technique described in the cervico-trochanteric region [START_REF] Premat | Reinforced cementoplasty using dedicated spindles in the management of unstable malignant lesions of the cervicotrochanteric region[END_REF].

Procedure

This is a monocentric observational retrospective study. All the spindle-reinforced cementoplasty procedures of the humerus performed at our Institution between January 1 st , 2017 and January 1 st , 2020 were reviewed. All the patients were discussed during a weekly multidisciplinary meeting dedicated to cancerous bone disease involving interventional radiologists, orthopedic surgeons, medical and radiation oncologists to validate the indication of RC. To be considered eligible for RC, patients had to meet the following criteria: (1) pathological fracture of the humerus causing pain and functional impairment (2) inability to undergo surgical treatment because of poor general condition, need for a systemic chemotherapy, or refused surgery.

All procedures were performed under general anesthesia, in supine position with both arms placed along the body, in a hybrid angiosuite that combines a flat panel C-arm and a dedicated CT scan (16row) (Miyabi Emotion 16, Siemens, Erlangen, Germany), by interventional radiologists with expertise in musculoskeletal interventional radiology ranging from 2 years to 30 years. At the beginning of the procedure, 1 g intravenous (IV) cefazolin (or other broad-spectrum antibiotics) was administered to minimize any risk of infection. An 8-gauge 150-mm bone needle (Thiebaud, Margencel, France) was inserted in a cranio-caudal fashion into the humeral head, laterally to the acromion, in the long axis of the humerus and advanced distally to the fracture site in the medullary cavity, under fluoroscopic guidance. A biopsy at the fracture site was systematically performed to confirm the malignant nature of the disease. Additional 11-gauge bone needles were inserted, under fluoroscopic guidance, transversally to the shaft of the humerus, into the lytic lesion, at the discretion of the operator for additional PMMA bone cement injection. 2.5-mm diameter spindles with length of 5-cm, 7-cm or 8cm cm (material: stainless steel in accordance with the international standard AISI 316 L [Thiebaud]) were inserted through the 8-gauge needle and pushed by means of the mandrel to be placed in the medullary cavity of the humerus, across the fracture site (Figure 1). The choice in size of the spindles was left at the discretion of the operator depending on the size of the lesion. For patients with other lytic lesions adjacent to the fracture site, additional spindles were inserted to form a steel frame in the medullary cavity across all the lesions. After control X-ray, PMMA bone cement (Biomet UK Ltd, Bridgend, UK) was injected into the needles to ensure a satisfactory filling (defined as superior to 50%) of the lesion and the fracture site. A post-operative CT-scan covering the treated humerus was systematically acquired to assess cement distribution, look for potential leakages, and to serve as a reference for the follow-up.

For retrospective analysis, medical records were screened until May 15 th , 2021. The following data was sought: secondary fractures, pre-and post-procedure pain estimation at 6 weeks using a visual analogue scale (VAS), technical success of the procedure, number of needles used during the procedure, number of spindles inserted at the fracture site, volume of PMMA injected, additional radiation therapy, and procedure-related complications. The occurrence of secondary fractures was defined as any one of the following events during the follow-up, attested by CT-imaging or plain radiographs: (1) any new fracture of the treated humerus (2) fracture of the cement plug (3) secondary displacement at the fracture site. Procedure-related complications were divided into major (death, permanent sequela, or symptomatic leak) and minor (transient pain, asymptomatic leak, or puncture site hematoma) complications and were graded according to the Cardiovascular and Interventional Radiological Society of Europe (CIRSE) classification [START_REF] Filippiadis | Cirse Quality Assurance Document and Standards for Classification of Complications: The Cirse Classification System[END_REF].

Baseline, per and post-procedural data are summarized in Table 1.

Discussion

Different strategies have been proposed to increase the resistance of PMMA cement to torsional forces and shear stress in long bones, but only a few case-reports emphasize the humerus. The use of cement-filled catheter has been reported [START_REF] Sun | Percutaneous cementoplasty for painful osteolytic humeral metastases: initial experience with an innovative technique[END_REF][START_REF] Kawai | Percutaneous Osteoplasty with Use of a Cement-filled Catheter for a Pathologic Fracture of the Humerus[END_REF]. Anselmetti et al. [START_REF] Anselmetti | Painful Pathologic Fracture of the Humerus: Percutaneous Osteoplasty With Bone Marrow Nails Under Hybrid Computed Tomography and Fluoroscopic Guidance[END_REF] described a percutaneous osteoplasty technique using metallic tunnelers designed for central venous catheter placement. Kelekis et al. [START_REF] Kelekis | Percutaneous Augmented Osteoplasty of the Humeral Bone Using a Combination of MicroNeedles Mesh and Cement[END_REF] performed augmented cementoplasty of the humerus using a combination of microneedles and bone cement. Percutaneous screw fixation has been described in the shoulder girdle for fractures involving the scapula [START_REF] Garnon | Percutaneous CT and Fluoroscopy-Guided Screw Fixation of Pathological Fractures in the Shoulder Girdle: Technical Report of 3 Cases[END_REF], but no case-report was found in the literature for humeral pathological fractures. We used the same approach previously described by Premat et al. [START_REF] Premat | Reinforced cementoplasty using dedicated spindles in the management of unstable malignant lesions of the cervicotrochanteric region[END_REF] for the treatment of malignant lesions of the cervicotrochanteric region, with the use of dedicated long bone spindles in addition to PMMA cement. The main advantage of this technique of RC is its simplicity: the 8-gauge bone needle can be placed accurately using only fluoroscopic guidance, one bone needle was needed in most procedures (4/6), and intralesional spindling and cement injection can be performed through the same needle. All complications observed in this series were considered minor complications, as they were asymptomatic. Cement leakages around the fracture site represented most of the complications (3/6). In one patient, two spindles migrated through the fracture site immediately after their insertion with their tip protruding in soft tissues. This is a potential limitation of this technique as it is almost impossible to reposition a misplaced spindle. This patient didn't report any residual pain during follow-up. No secondary fracture was observed, which suggests the efficacy of this procedure on bone consolidation. However, advanced imaging techniques and post-processing tools will help in the future to better evaluate the absence of secondary fracture, especially cement plug fragmentation [START_REF] Blum | 3D reconstructions, 4D imaging and postprocessing with CT in musculoskeletal disorders: Past, present and future[END_REF][START_REF] Garnon | Cement Plug Fragmentation Following Percutaneous Cementoplasty of the Bony Pelvis: Is it a Frequent Finding in Clinical Practice?[END_REF].. To the best of our knowledge, this is the largest series of RC for malignant fractures of the humerus.

With experimental studies questioning the efficacy of augmented osteoplasty in long bones [START_REF] Cazzato | Biomechanical effects of osteoplasty with or without Kirschner wire augmentation on long bone diaphyses undergoing bending stress: implications for percutaneous imaging-guided consolidation in cancer patients[END_REF], larger comparative trials would be needed to assess the efficacy of these procedures.

Disclosure

All the authors declare no conflict of interest regarding this study.

Figure 1 :

 1 Figure 1: Fifty-one-year-old woman with multiple myeloma. a. multiple osteolytic formations of the

Table 1 :

 1 Baseline, procedural and follow-up data

	Patient			1	2	3	4	5	6
	Age (years)		46	70	61	51	70	58
	Sex			M	F	M	F	M	M
	Primary tumor Colon cancer Multiple	Multiple	Multiple	Multiple	Multiple
					myeloma	myeloma	myeloma	myeloma	myeloma
	Side			Left	Left	Right	Right	Left	Left
	Fracture		Metaphysis	Metaphysis	Metaphysis Diaphysis	Metaphysis	Metaphysis
	location
	Displacement No	Yes	No	Yes	Yes	Yes
	Number	of	1	1	1	1	3	3
	bone needles
	inserted
	Number	of	NA	NA	7	NA	30	12
	spindles
	inserted
	Volume		of	NA	NA	25	NA	20	21
	PMMA (ml)
	Vertebroplasty	No	No	Yes (T11,	No	Yes (T2, T5, T9) No
	in the same			T12, L4)
	intervention
	Total			112	80	130	58	122	140
	procedure
	time (minutes)
	Total			17	11.8	39	15	NA	NA
	fluoroscopy
	time (minutes)
	Procedure-		Minor:	Minor:	No	No	Minor:	Two	Minor:
	related			Leakages	Leakages			spindles	Leakages
	complications	(subacromial-	(subacromial-			protruding	(soft
				subdeltoid bursa)	subdeltoid bursa)			fracture from	site the	tissue)
								(asymptomatic)
	Radiation		Yes	No	Yes	No	Yes	No
	therapy
	Duration	of	35.8	4.5 (death)	1.7	46.5	14 (death)	15.6
	follow-up
	(months)
	Baseline VAS	10	10	10	9	8	10
	VAS	at	6	0
	weeks