Dario Magliacano
email: dario.magliacano@univ-fcomte.fr

Sepide Ashani

Morvan Ouisse

Elke Deckers

Giuseppe Petrone

Wim Desmet

Sergio De Rosa

Formulation and validation of the shift cell technique for acoustic applications of poro-elastic materials described by the Biot theory

Keywords: vibroacoustics, Biot theory of poro-elasticity, shift cell, FEM, dispersion diagram

The inclusion of vibroacoustic treatments at early stage of product development through the use of poro-elastic media with periodic inclusions, which exhibit proper dynamic filtering effects, is a powerful strategy for the achievement of lightweight sound packages and represents a convenient solution for manufacturing aspects. This can have different applications in transportation (aerospace, automotive, railway), energy and civil engineering fields, where weight, space and vibroacoustic comfort are still critical challenges. This paper develops the shift cell operator approach as a numerical tool to investigate the dispersion characteristics of periodic poro-elastic media. It belongs to the class of the k(ω) (wave number as a function of the angular frequency) methods and leads to a quadratic eigenvalue problem, even when considering frequency-dependent materials, contrarily to the ω(k) approach that would lead to a non-linear eigenvalue problem for frequency-dependent materials.

Introduction

Fast urbanization and transport development cause serious noise-induced health risks, such as annoyance, sleep disturbance, or even ischemic heart disease [START_REF] Cao | Porous materials for sound absorption[END_REF]. Therefore, nowadays, environment noise control is becoming a subject of great interest. Generally, common sound absorbing materials could be divided into two categories: resonant [START_REF] Zhao | Improving low-frequency sound absorption of micro-perforated panel absorbers by using mechanical impedance plate combined with helmholtz resonators[END_REF] and poro-elastic materials. Resonant materials for sound absorption mainly involve Helmholtz resonators [START_REF] Cai | Noise attenuation capacity of a helmholtz resonator[END_REF] and/or perforated panels [START_REF] Yildiz | Properties of sound panels made from recycled footwear treads[END_REF]. These materials show good performances at low frequencies, but they often have the disadvantage of narrow frequency stop-bands [START_REF] Lv | Effect of micro-slit plate structure on the sound absorption properties of discarded corn cob husk fiber[END_REF]. Poro-elastic materials for acoustic applications are composed of channels, cracks or cavities that allow the sound waves entering the materials. Sound energy is dissipated by thermal and viscous losses; these energy consumption principles assure sound absorption over broader frequency ranges [START_REF] Berardi | Acoustic characterization of natural fibers for sound absorption applications[END_REF][START_REF] Xinzhao | Electrically conductive graphene-coated polyurethane foam and its epoxy composites[END_REF]. Poro-elastic materials suffer from a lack of performance at low frequencies compared to their efficiency at higher ones [START_REF] Groby | Using simple shape three-dimensional inclusions to enhance porous layer absorption[END_REF]. This difficulty is usually overcome by multi-layering [START_REF] Yang | Wave and finite element method for predicting sound transmission through finite multi-layered structures with fluid layers[END_REF]; however, the efficiency of such devices relies on the allowable thickness [START_REF] Weisser | Acoustic behavior of a rigidly backed poroelastic layer with periodic resonant inclusions by a multiple scattering approach[END_REF][START_REF] Gaborit | A simplified model for thin acoustic screens[END_REF]. An efficient way to enhance the low frequency performances of sound packages consists in embedding periodic inclusions in a poro-elastic layer [START_REF] Groby | Acoustic response of a rigid-frame porous medium plate with a periodic set of inclusions[END_REF][START_REF] Xiong | Sound attenuation optimization using metaporous materials tuned on exceptional points[END_REF], in order to create wave interferences or resonance effects that may be advantageous for the dynamics of the system. In this context of increasingly complex material systems, numerical tools to properly design sound packages are more and more useful. Several theoretical models are available to estimate the physical behavior of poro-elastic materials, and the most complex of them require the definition of more than ten parameters. For example, one of the most accurate models is the Biot theory of poro-elasticity [START_REF] Biot | Mechanics of deformation and acoustic propagation in porous media[END_REF], which takes into account both the mechanical and the acoustical behaviors of the material [START_REF] Detournay | Fundamentals of poroelasticity[END_REF]. Furthermore, the measurement of all the necessary parameters, which usually constitutes the first step in the definition of a model, is already a specific issue in the case of poro-elastic.

In addition, numerical simulations, usually carried out through the Finite Element Method (FEM), are often problematic, in terms of computational times and convergence. On the other hand, analytical models constitute a powerful instrument to quickly catch physics and general trends of the problem, but they are partially limited by restrictive approximating hypotheses and come short considering non-trivial geometries. In this context, the present work investigates the application of the shift cell approach to poroelastic media; this allows to obtain dispersion characteristics of frequencydependent damped materials through the resolution of a quadratic eigenvalue problem, whose accuracy only depends on the FEM meshing. This technique has already been successfully applied to describe the mechanical behavior of periodic structures embedding visco-elastic materials [START_REF] Collet | Floquet-bloch decomposition for the computation of dispersion of two-dimensional periodic, damped mechanical systems[END_REF][START_REF] Billon | Design and experimental validation of a temperature-driven adaptive phononic crystal slab[END_REF], piezoelectric materials [START_REF] Collet | Adaptive metacomposites for vibroacoustic control applications[END_REF] and foams modeled as equivalent fluids [START_REF] Magliacano | Computation of dispersion diagrams for periodic porous materials modeled as equivalent fluids[END_REF]. The main novelty of the present work consists in the formulation and application of the shift cell technique to Biot-modeled poro-elastic media. Materials modeled in this way account for wave propagation and interaction in both fluid and solid phases, thus leading to the fact that diphasic models are the most comprehensive ones in order to describe the vibroacoustics of porous media. However, compared to equivaled fluid models, they require more parameters to be used (a set for each of the two phases), and therefore the process of extension of the shift cell technique is definetly not trivial and requires a specific dissertation, which is herein provided for the first time in literature. This paper is organized as it follows. Section 2 recalls the fundamentals of Biot theory and introduces the shift cell operator formulation for Biotmodeled foams. Section 2.2.2 defines a weak formulation of the problem, and Section Appendix B describes its FE implementation. In Section 3 two validations of the method are shown. At last, Section 4 provides conclusions and future perspectives.

Shift cell operator technique for Biot-modeled foams

Biot theory

Although for many porous materials the frame can be considered almost rigid for a wide range of acoustical frequencies, thus allowing the use of models with motionless skeleton assumption [START_REF] Deckers | Modelling techniques for vibro-acoustic dynamics of poroelastic materials[END_REF][START_REF] Allard | Propagation of sound in porous media: Modelling sound absorbing materials[END_REF], this is not generally true: for example, for a poro-elastic material attached to a vibrating structure and for many other similar situations, frame vibrations are induced by those of the elastic structure.

The wave propagation through a poro-elastic media can be analyzed only considering a solid-fluid coupled behavior; such description is provided by the Biot theory of sound propagation in poro-elastic media [START_REF] Biot | Mechanics of deformation and acoustic propagation in porous media[END_REF]. In this context, two compressional waves and a shear wave propagate. The parameters that characterize a poro-elastic material are: φ is the open porosity; σ is the static flow resistivity; α ∞ is the tortuosity; Λ is the viscous characteristic length; η visc is the viscosity; q 0 = η visc σ is the static viscous permeability; ν visc = η visc ρ 0 ; ν therm = ν visc P r ; P r is the Prandtl number. Furthermore, additional quantities are defined in Appendix A [START_REF] Atalla | A mixed displacement-pressure formulation for poroelastic materials[END_REF]. Zienkiewicz et al. proposed a simplified u -p formulation [START_REF] Zienkiewicz | Dynamic behaviour of saturated porous media; the generalized biot formulation and its numerical solution[END_REF], where u is the solid phase displacement and p is the pressure of the fluid phase. In particular, by neglecting the second time derivatives of the relative fluid displacement from the original Biot u -U formulation [START_REF] Detournay | Fundamentals of poroelasticity[END_REF], the u -p formulation [START_REF] Atalla | A mixed displacement-pressure formulation for poroelastic materials[END_REF][START_REF] Atalla | Enhanced weak integral formulation for the mixed (u,p) poroelastic equations[END_REF] is deduced in order to reduce the primary variables in the context of finite element analysis; indeed, if one considers a 3D model, instead of the 3+3 nodal variables that are in the u -U formulation, in the case of the u -p one there are only 3+1 nodals variables. In addition, the solid displacement u and the pore fluid pressure p are always the most interesting quantities. In an infinite homogeneous isotropic poro-elastic media, three waves propagate (two compressional waves and one shear wave):

k shear = ω ρ11 ρ22 -ρ12 2 N ρ22 , (1)
k f ast,slow = A 1 2 ± A 2 1 4 -A 2 , with (2)
A 1 = ω 2 ρ11 R -2 ρ12 Q + ρ22 P RP -Q 2 , A 2 = ω 4 ρ11 ρ22 -ρ12 2 RP -Q 2 . ( 3
)
The symbols introduced in Eq. 1-3 are defined in Appendix A. The two phases present in a poro-elastic material behave in a different manner, respect to the pure elastic case (where the onyl compressional wave is fluid-born): the main difference is the existence of a second (solid-born) compressional wave, which is highly attenuated in the low frequency range. Each of the waves propagates both in the solid and in the fluid phases of the poro-elastic medium [START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF].

Shift cell operator technique 2.2.1. Introduction

Herein, the shift cell operator technique applied to Biot-modeled foams is presented, providing details on its implementation [START_REF] Billon | Composites périodiques fonctionnels pour l'absorption vibroacoustique large bande[END_REF]. The shift cell approach provides a reformulation of classical Floquet-Bloch periodic conditions [START_REF] Floquet | Sur les équations différentielles linéaires à coefficients périodiques [On the linear differential equations with periodic coefficients[END_REF], and its major advantage is that it allows the introduction of a generic frequency dependence of visco-elastic material behavior [START_REF] Collet | Floquet-bloch decomposition for the computation of dispersion of two-dimensional periodic, damped mechanical systems[END_REF]; this is fundamental, if one looks for the computation of the dispersion curves of a porous material, modeled as an equivalent fluid or with the Biot theory. Indeed, even if the usage of Floquet-Bloch (F-B) periodic conditions actually allows it, a very powerful non-linear solver is required in that case. The shift cell operator [START_REF] Collet | Floquet-bloch decomposition for the computation of dispersion of two-dimensional periodic, damped mechanical systems[END_REF][START_REF] Magliacano | Computation of dispersion diagrams for periodic porous materials modeled as equivalent fluids[END_REF], instead, leads to a quadratic eigenvalue problem even in the presence of frequency-dependences and/or damping. The main mathematical difference with respect to the classical F-B approach is that, in the case of the shift cell operator, the phase shift of the boundary conditions and the exponential amplitude decrease, related to wave propagation, are integrated into the partial derivative operator. As a consequence, the periodicity is included in the overall behavior of the structure, while simple continuity conditions are imposed at the edges of the unit cell. Considering a poro-elastic layer modeled through Biot's theory [START_REF] Biot | Mechanics of deformation and acoustic propagation in porous media[END_REF], the coupled starting system is constituted by the equation of motion of the solid part and the classical Helmholtz equation, respectively:

     ∇ • σ(u) + ω 2 ρu + γ∇p = 0 ∆p + ω 2 ρ 22 R p -ω 2 ρ 22 φ 2 γ∇ • u = 0 , (4)
where u = (u, v, w) is the solid phase displacement vector and p = p(x, ω) is the acoustic pressure [START_REF] Bensoussan | Asymptotic Analysis of Periodic Structures[END_REF]. The following quantities are defined [START_REF] Atalla | A mixed displacement-pressure formulation for poroelastic materials[END_REF]: ω is the angular frequency; σ(u) = Cε(u) is the stress tensor of the frame in vacuum, whose generic element can be written as

σ ij = (µ 1 -Q 2 R )δ ij ε kk + 2µ 2 ε ij , where δ ij is the Kronecker delta and ε kk = tr(ε) = ε ux + ε vy + ε wz ; C is the Hooke elasticity tensor with C 11 = µ 1 -Q 2 R + 2µ 2 and C 12 = µ 1 -Q 2 R ; ε(u) = 1 2 (∇u + ∇u T )
is the symmetric strain tensor; µ 1 = 2ν 1-2ν N and µ 2 = N are respectively the first and second Lamé parameters.

For each physical property of the system, the periodicity is described by α(x -rn) -α(x) = 0, where α is a generic physical property, n is a vector of integers normal to the face considered, r = (r 1 ; r 2 ; r 3 ) is a matrix containing the three vectors defining the cell periodicity directions and lengths, and Ω is the domain of interest. This applies everywhere except on the discontinuity surfaces, where appropriate boundary conditions apply [START_REF] Magliacano | Computation of dispersion diagrams for periodic porous materials modeled as equivalent fluids[END_REF]. By further developing the latter equation and applying the Bloch theorem [START_REF] Bloch | Uber die quantenmechanik der elektronen in kristallgittern [On the quantum mechanics of the electrons in crystal lattices[END_REF], which extends Floquet's theory to 3D systems, one obtains:

           (∇ + jk) • C 1 2 ((∇ + jk)u + (∇ + jk)u T )+ + ω 2 ρu + γ(∇ + jk)p = 0 (∇ + jk) T • (∇ + jk)p + ω 2 ρ22 R p -ω 2 ρ22 φ 2 γ(∇ + jk) • u = 0 , (5)
with the wave vector k defined as:

k = kθ = k   θ x θ y θ z   = k   cos θ cos φ cos θ sin φ sin θ   (6)
and k = -jλ, where λ is an eigenvalue of the problem.

Weak formulation

The solution approach follows a common weak formulation of a differential problem in a discrete coordinate scheme. A (u, p) formulation, in its classical form, can be found in literature [START_REF] Atalla | Enhanced weak integral formulation for the mixed (u,p) poroelastic equations[END_REF]:

                                         Ω σ(u) : ε(δu)dΩ -ω 2 Ω ρu • δudΩ+ - Ω (γ + φ(1 + Q R ))∇p • δudΩ - Ω φ(1 + Q R )p∇ • δudΩ+ - Γ (σ T (u, p) • n) • δudΓ = 0 Ω φ 2 ω 2 ρ22 ∇p • ∇δpdΩ - Ω φ 2 R pδpdΩ+ - Ω (γ + φ(1 + Q R ))∇δp • udΩ - Ω φ(1 + Q R )δp∇ • udΩ+ - Γ φ(U n -u n )δpdΓ = 0 , (7)
where δu and δp are admissible variations of the solid phase displacement vector and the interstitial fluid pressure of the poro-elastic medium, respectively. Considering that σ(u) = Cε(u) = C 1 2 (∇u + ∇u T ), and introducing the shift cell operator as explained above, one obtains:

145                                          Ω (C 1 2 ((∇ + jk)u + (∇ + jk)u T )) : ((∇ -jk)δu + (∇ -jk)δu T )dΩ+ -ω 2 Ω ρu • δudΩ - Ω (γ + φ(1 + Q R ))(∇ + jk)p • δudΩ+ - Ω φ(1 + Q R )p(∇ -jk) • δudΩ = 0 Ω φ 2 ω 2 ρ22 (∇ + jk)p • (∇ -jk)δpdΩ+ - Ω φ 2 R pδpdΩ - Ω (γ + φ(1 + Q R ))(∇ -jk)δp • udΩ+ - Ω φ(1 + Q R )δp(∇ + jk) • udΩ = 0 , (8)
where the boundary condition caused the integral on the boundary to vanish. Therefore, one can define the following quantities:

• σ θ (u) = Cε θ (u), whose generic term is σ θij = (µ 1 -Q 2 R )δ ij ε θkk +2µ 2 ε θij ; • ε θ (u) = 1 2 (θu + θu T ).
Therefore:

150                                                          Ω σ(u) : ε(δu)dΩ + jk Ω σ θ (u) : ε(δu)dΩ+ -jk Ω σ(u) : ε θ (δu)dΩ + k 2 Ω σ θ (u) : ε θ (δu)dΩ+ -ω 2 Ω ρu • δudΩ - Ω (γ + φ(1 + Q R ))(∇ + jk)p • δudΩ+ - Ω φ(1 + Q R )p(∇ -jk) • δudΩ = 0 Ω φ 2 ω 2 ρ22 ∇p • ∇δpdΩ + jk Ω φ 2 ω 2 ρ22 θ • p∇δpdΩ+ -jk Ω φ 2 ω 2 ρ22 θ • ∇pδpdΩ + k 2 Ω φ 2 ω 2 ρ22 pδpdΩ - Ω φ 2 R pδpdΩ+ - Ω (γ + φ(1 + Q R ))∇δp • udΩ + jk Ω (γ + φ(1 + Q R ))θ • δpudΩ+ - Ω φ(1 + Q R )δp∇ • udΩ -jk Ω φ(1 + Q R )θ • δpudΩ = 0 , (9)
which can be written in a more structured form, as:

                                         Ω σ(u) : ε(δu)dΩ + jk Ω ( σ θ (u) : ε(δu) -σ(u) : ε θ (δu))dΩ+ + k 2 Ω σ θ (u) : ε θ (δu)dΩ -ω 2 Ω ρu • δudΩ - Ω γ∇p • δudΩ+ -jk Ω γθ • pδudΩ - Ω φ(1 + Q R )(∇p • δu + p∇ • δu)dΩ = 0 Ω φ 2 ω 2 ρ22 ∇p • ∇δpdΩ + jk Ω φ 2 ω 2 ρ22 (θ • p∇δp -θ • ∇pδp)dΩ+ + k 2 Ω φ 2 ω 2 ρ22 pδpdΩ - Ω φ 2 R pδpdΩ - Ω γ∇δp • udΩ+ + jk Ω γθ • δpudΩ - Ω φ(1 + Q R )(∇δp • u + δp∇ • u)dΩ = 0 . (10)
Finally, one can discretize the weak formulation through the FE Method: considering that ϕ s and ϕ f are the eigenvectors of the solid and fluid parts respectively, the system of equations can be written in its matrix form:

(K s + jkL s + k 2 H s -ω 2 M s )ϕ s -(N s + jkO s + T s )ϕ f = 0 ((K f + jkL f + k 2 H f -ω 2 M f )ϕ f -ω 2 (N f -jkO f + T f )ϕ s = 0 , (11)
with the following matrices (∝ means "proportional to"):

• K s ∝ Ω σ(u) : ε(δu)dΩ; • L s ∝ Ω ( σ θ (u) : ε(δu) -σ(u) : ε θ (δu))dΩ; • H s ∝ Ω σ θ (u) : ε θ (δu)dΩ; • M s ∝ Ω ρu • δudΩ; • N s ∝ Ω γ∇p • δudΩ; • O s ∝ Ω γθ • pδudΩ; • T s ∝ Ω φ(1 + Q R )(∇p • δu + p∇ • δu)dΩ; • K f ∝ Ω φ 2 ρ 22 ∇p • ∇δpdΩ; • L f ∝ Ω φ 2 ρ 22 (θ • p∇δp -θ • ∇pδp)dΩ; • H f ∝ Ω φ 2 ρ 22 pδpdΩ; • M f ∝ Ω φ 2
R pδpdΩ;

• N f ∝ Ω γ∇δp • udΩ; • O f ∝ Ω γθ • δpudΩ; • T f ∝ Ω φ(1 + Q R )(∇δp • u + δp∇ • u)dΩ.
Here, M s,f and K s,f are respectively the symmetric mass and symmetric stiffness matrices, L s,f are skew-symmetric matrices, H s,f are symmetric matrices and N s = N T f , O s = O T f and T s = T T f are the matrices that couple the solid and fluid behaviors; all of them are complex and frequency-dependent. Therefore, the coupled system can be written as it follows:

(K s + jkL s + k 2 H s -ω 2 M s ) -(N f -jkO f + T f ) -(N s + jkO s + T s ) 1 ω 2 (K f + jkL f + k 2 H f -ω 2 M f ) ϕ s ϕ f = = 0 0 . ( 12
)
The details of the FE implementation are given in Appendix B.

Validation of the method

In order to validate the shift cell technique implementation for Biotmodeled foams and for waves propagating along the x-axis, two different comparisons are provided: one with an application of the shift cell approach to an equivalent fluid [START_REF] Magliacano | Computation of dispersion diagrams for periodic porous materials modeled as equivalent fluids[END_REF], and another one with a WFEM analysis performed on a Biot-modeled foam [START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF].

Biot model with shift cell vs. JCA model with shift cell

The first considered system is a homogeneous foam with material properties shown in Table 1, represented by a cubic unit cell having a volume of 8 cm 3 , with periodicity in three directions and mesh composed by 10 thetraedral elements along each side of the cube. The second case is constructed by introducing a rigid cylindrical inclusion with radius equal to 0.5 cm at the center of the previous unit cell, as shown in Figure 1. In Figure 2 and Figure 3, dispersion curves of two different systems with an artificially high value of frame Young modulus (E = 10 15 Pa) and nullified loss factor, such that the rigid frame assumption would be valid, are calculated using the shift cell approach and compare the results obtained through the Biot model with those calculated using a Johnson-Champoux-Allard (JCA) [START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF][START_REF] Champoux | Dynamic tortuosity and bulk modulus in air-saturated porous media[END_REF] equivalent fluid [START_REF] Magliacano | Computation of dispersion diagrams for periodic porous materials modeled as equivalent fluids[END_REF]. Therefore, the elasticity of the skeleton is neglected and the Biot model essentially describes the behavior of the equivalent fluid one. The distinction between propagative and evanescent waves is obtained, in a first approximation, through the application of the 1 st classifying criterion described by Magliacano et al. [START_REF] Magliacano | Computation of dispersion diagrams for periodic porous materials modeled as equivalent fluids[END_REF] for equivalent fluids. Looking at Figure 2 and Figure 3, it can be noticed that the comparison shows an almost perfect agreement between the results of the shift cell technique applied on the two different foam models. The advantage of using Biot model, for which the shift cell approach is developed herein, relies on the fact that, as already introduced in Section 1, in some cases (for example: lowfrequency acoustic loads, or mechanical excitations) waves can propagate in both fluid and solid phases. In those contexts, motionless skeleton models cannot be used and a more general diphasic model (like Biot's one) is required in order to describe the poro-elastic behavior of the foam [START_REF] Stinson | The propagation of plane sound waves in narrow and wide circular tubes, and generalization to uniform tubes of arbitrary cross-sectional shape[END_REF]. Moreover, if the frequency range of the study is under the decoupling frequency, which is located at high frequencies for foams with high value of flow resistivity, the equivalent fluid model prediction deviates significantly from the Biot theory; therefore, also in these cases it is necessary to use the latter, in order to have an accurate overview of the wave propagation in the medium. In addition, as it is more clear in Figure 4 and Figure 5, the shift cell approach is capable to catch the behavior of the three types of waves propagating in a porous material with elastic frame.

Biot model with shift cell vs. Biot model with WFEM

In this validation case, shift cell results are compared to those obtained by Serra et al. [START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF] using the Wave Finite Element Method [START_REF] Zhong | On the direct solution of wave propagation for repetitive structures[END_REF] (labeled as "reference" in Figure 4 and Figure 5). Parameters of foam and air used in this validation case can be found in Appendix B of Serra et al. [START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF], and are reported here in Table 2 andTable 3. In the case of poro-elastic media the rigidity of the material is very low, leading to very small wavelengths, and a high dissipation rate occurs within the pores; despite these difficulties, in the paper by Serra et al. [START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF] it is shown that the WFEM provides an efficient tool to compute the waves propagating through poro-elastic media. This validation is also performed with curves computed through the analytical model described in Section 2.1. As it is clear in Figure 4 and Figure 5, wavenumbers calculated using the shift cell approach applied to a Biotmodeled foam completely agree with those calculated through the analytical model; moreover, it can be seen that the slow compressional wave is highly attenuated. The shift cell approach has several advantages, in terms of linearity and convergence, compared to the WFEM. Indeed, as described by Serra et al. [START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF], the WFEM applied to Biot-modeled foams leads to a transcendental eigenvalue problem that can be solved only by using a nonlinear solver. However, there are still a lot of numerical difficulties, and robust solutions have not yet been developed [START_REF] Dauchez | Convergence of poroelastic finite elements based on biot displacement formulation[END_REF]. In the case of WFEM, the use of 10 elements per wavelength in the three directions is recommended as a rule of the thumb [START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF]. Under the hypotheses of plane wave, the use of the shift cell approach leads directly to a quadratic eigenvalue problem, with no assumption on the nature of the waves, whose accuracy only depends on the mesh chosen to discretize the system.

Conclusions

An efficient way to enhance the low frequency performances of sound packages consists in embedding periodic inclusions in a poro-elastic layer, in order to create wave interferences or resonance effects that may be advantageous for the dynamics of the system. This work develops the shift cell technique as a numerical tool to investigate the dispersion characteristics of periodic Biot-modeled poro-elastic media, providing details on its FEM implementation too; this approach allows to obtain dispersion characteristics of frequency-dependent damped materials through the resolution of a quadratic eigenvalue problem, whose accuracy only depends on the FEM meshing. A first validation of the shift cell approach for Biot-modeled poro-elastic materials has been obtained through a comparison with the results obtained on a JCA-modeled 3D unit cell, both in a homogeneous configuration and with a perfectly rigid cylindrical inclusion. For this purpose, the elasticity of the foam skeleton has been neglected and therefore the Biot model essentially described the behavior of an equivalent fluid, thus allowing the comparison between dispersion curves obtained through the application of the shift cell approach to Biot-modeled foams and equivalent fluids.

An additional validation has then been carried out through a comparison of the shift cell results with those obtained using the Wave Finite Element Method, and those computed through an analytical model that is valid for infinite homogeneous isotropic poro-elastic media; in this context, compared to the WFEM, the shift cell technique shows significant computational advantages. The outcome of this research is very promising, since the methodological basis and its validations are given in order to trace future characterizations and applications of periodic poro-elastic media in acoustics.

• b = σφ 2 G(ω) is the viscous drag;

• G(ω) = 1 + 4jα 2 ∞ η visc ρ 0 ω (σΛφ) 2
is the relaxation function, as predicted by JCA model [START_REF] Johnson | Theory of dynamic permeability and tortuosity in fluid-saturated porous media[END_REF][START_REF] Champoux | Dynamic tortuosity and bulk modulus in air-saturated porous media[END_REF];

• ρ = ( ρ11 -ρ 12 2 ρ 22 ); • P, Q, R are elasticity coefficients to be determined by "gedanken experiments" [START_REF] Allard | Propagation of sound in porous media: Modelling sound absorbing materials[END_REF]; in particular [START_REF] Deckers | Modelling techniques for vibro-acoustic dynamics of poroelastic materials[END_REF]:

P = (1-φ)(1-φ-KB KS )KS+φ KBKS KF 1-φ-KB KS +φ KS KF - 2 3 N ∼ = (1 + ν 1-2ν )2N + 1-φ 2 φ KF , Q = (1-φ-KB KS )φKS 1-φ-KB KS +φ KS KF ∼ = (1 -φ)KF , R = φ 2 KS 1-φ-KB KS +φ KS KF ∼ = φKF ; • N = |N |(1 + jη) = Y 2(1+ν)
is the complex shear modulus of the frame;

• Y = |Y |(1 + jη) is the complex Young modulus of the frame;

• η is the loss factor of the frame;

• ν is the Poisson's ratio of the frame;

• KB = 2N (ν+1) 3(1-2ν) is the bulk modulus of the the solid phase in vacuum;

• KS = KB 1-φ is the bulk modulus of the solid phase; • KF is the bulk modulus of the fluid phase, computed starting from the equivalent one (e.g.: KF = φK JCA );

• γ = φ( ρ 12 ρ 22 -Q R ); • µ i = P δ 2 i -ω 2 ρ 11 ω 2 ρ 12 -Qδ 2 i , i = 1,
ε(u) =   ∂u ∂x 1
ε θ (u) =   θ x u 1 2 (θ y u + θ x v) 1 2 (θ z u + θ x w) 1 2 (θ y u + θ x v) θ y v 1 2 (θ z v + θ y w) 1 2 (θ z u + θ x w) 1 2 (θ z v + θ y w) θ z w   , (B.3) σ(u) =   C 11 ∂u ∂x + C 12 ( ∂v ∂y + ∂w ∂z ) (C 11 -C 12 ) 1 2 ( ∂u ∂y + ∂v ∂x ) (C 11 -C 12 ) 1 2 ( ∂u ∂z + ∂w ∂x ) (C 11 -C 12 ) 1 2 ( ∂u ∂y + ∂v ∂x ) C 11 ∂v ∂y + C 12 ( ∂u ∂x + ∂w ∂z ) (C 11 -C 12 ) 1 2 ( ∂v ∂z + ∂w ∂y ) (C 11 -C 12 ) 1 2 ( ∂u ∂z + ∂w ∂x ) (C 11 -C 12 ) 1 2 ( ∂v ∂z + ∂w ∂y ) C 11 ∂w ∂z + C 12 ( ∂u ∂x + ∂v ∂y )   (B.4) σ θ (u) =   C 11 θ x u + C 12 (θ y v + θ z w) (C 11 -C 12 ) 1 2 (θ y u + θ x v) (C 11 -C 12 ) 1 2 (θ z u + θ x w) (C 11 -C 12 ) 1 2 (θ y u + θ x v) C 11 θ y v + C 12 (θ x u + θ z w) (C 11 -C 12 ) 1 2 (θ z v + θ y w) (C 11 -C 12 ) 1 2 (θ z u + θ x w) (C 11 -C 12 ) 1 2 (θ z v + θ y w) C 11 θ z w + C 12 (θ x u + θ y v)   (B.5)
The numerical model is based on the following matrix weak formulation, proposed to provide an expression optimized for the FE implementation:

)( ∂δu ∂z + ∂δw ∂x )))dΩ; • L s,u ∝ Ω ((C 11 θ x u+C 12 (θ y v +θ z w)) ∂δu ∂x +(C 11 -C 12 ) 1 4 ((θ y u+θ x v)( ∂δu ∂y + ∂δv ∂x ) + (θ z u + θ x w)( ∂δu ∂z + ∂δw ∂x )) -(C 11 ∂u ∂x + C 12 ( ∂v ∂y + ∂w ∂z ))θ x δu -(C 11 - C 12 ) 1 4 (( ∂u ∂y + ∂v ∂x )(θ y δu + θ x δv) + ( ∂u ∂z + ∂w ∂x )(θ z δu + θ x δw)))dΩ; • H s,u ∝ Ω ((C 11 θ x u+C 12 (θ y v+θ z w))θ x δu+(C 11 -C 12 ) 1 4 ((θ y u+θ x v)(θ y δu+ θ x δv) + (θ z u + θ x w)(θ z δu + θ x δw)))dΩ; • M s,u ∝ Ω ρuδudΩ; • N s,u ∝ Ω γ ∂p ∂x δudΩ; • O s,u ∝ Ω γθ 1 pδudΩ; • T s,u ∝ Ω φ(1 + Q R )( ∂p ∂x δu + p ∂δu ∂x )dΩ; • K s,v ∝ Ω ((C
• M f ∝ Ω φ 2
R pδpdΩ;

• N f ∝ Ω γ(u ∂δp ∂x + v ∂δp ∂y + w ∂δp ∂z )dΩ;

• O f ∝ Ω γ(θ 1 u + θ 2 v + θ 3 w)δpdΩ;

• T f ∝ Ω φ(1 + Q R )((

  Young modulus [kPa] 70+j9 Resistivity [Pa*s/m 2 ] 3750 Shear modulus [kPa] 25+j7 Viscous char. length [mm] 0.11 Loss factor 0.265 Thermal char. length [mm] 0.742 Poisson ratio 0.39

Figure 1 :

 1 Figure 1: 3D unit cell constituted by a 2 cm cube, homogeneous (on the left) and with a 5 mm radius cylindrical hole (on the right).

Figure 4 :

 4 Figure 4: Dispersion curve comparison with the reference (WFEM by Serra et al. [25]), and analytical model; real part of the wavenumber.

Figure 5 :

 5 Figure 5: Dispersion curve comparison with the reference (WFEM by Serra et al. [25]), and analytical model; imaginary part of the wavenumber.

2

 2 is the ratio of the velocity of the air over the velocity of the frame for the two compressional waves and indicates in what medium the waves propagate preferentially.Appendix B. Finite element implementationIn order to numerically implement the shift cell technique for Biot-modeled foams, the vector equation related to the motion of the solid part is split into three scalar equations. The following matrices are defined accordingly: θ x v θ x w θ y u θ y v θ y w θ z u θ z v θ z w

Table 1 :

 1 Properties of a PU 60 foam.

Table 2 :

 2 Properties of the foam used in the validation with the work by Serra et al.[START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF].

	Frequency [Hz]	Prop. Re(k x ) JCA Prop. Im(k x ) JCA Evan. Re(k x ) JCA Evan. Im(k x ) JCA
		Prop. Re(k x ) Biot
		Prop. Im(k x ) Biot
		Evan. Re(k x ) Biot
		Evan. Im(k x ) Biot
	-1	-0.8	-0.6	-0.4	-0.2	0	0.2	0.4	0.6	0.8	1
						k x r/

Figure 2: Dispersion curves validation with JCA plots; here, the Biot curves are computed for a homogeneous PU 60 foam, with E = 10 15 Pa and structural loss factor = 0.

Table 3 :

 3 Properties of the air used in the validation with the work by Serra et al.[START_REF] Serra | Wave properties in poroelastic media using a wave finite element method[END_REF].

  • L s,v ∝ Ω ((C 11 θ y v +C 12 (θ x u+θ z w)) ∂δv ∂y +(C 11 -C 12 ) 1 4 ((θ y u+θ x v)( ∂δu ∂y + ∂δv ∂x ) + (θ z v + θ y w)( ∂δv ∂z + ∂δw ∂y )) -(C 11 ∂v ∂y + C 12 ( ∂u ∂x + ∂w ∂z ))θ y δv -(C 11 -C 12 ) 1 4 (( ∂u ∂y + ∂v ∂x )(θ y δu + θ x δv) + ( ∂v ∂z + ∂w ∂y )(θ z δv + θ y δw)))dΩ;• H s,v ∝ Ω ((C 11 θ y v+C 12 (θ x u+θ z w))θ y δv+(C 11 -C 12 ) 1 4 ((θ y u+θ x v)(θ y δu+ θ x δv) + (θ z v + θ y w)(θ z δv + θ y δw)))dΩ; • M s,v ∝ Ω ρvδvdΩ; • N s,v ∝ Ω γ ∂p ∂y δvdΩ; • O s,v ∝ Ω γθ 2 pδvdΩ; • L s,w ∝ Ω ((C 11 θ z w+C 12 (θ x u+θ y v)) ∂δw ∂z +(C 11 -C 12 ) 1 4 ((θ y w+θ z v)( ∂δw ∂y + ∂δv ∂z ) + (θ z u + θ x w)( ∂δu ∂z + ∂δw ∂x )) -(C 11 ∂w ∂z + C 12 ( ∂u ∂x + ∂v ∂y ))θ z δw -(C 11 -C 12 ) 1 4 (( ∂w ∂y + ∂v ∂z )(θ y δw + θ z δv) + ( ∂u ∂z + ∂w ∂x )(θ z δu + θ x δw)))dΩ; • H s,w ∝ Ω ((C 11 θ z w+C 12 (θ x u+θ y v))θ z δw+(C 11 -C 12 ) 1 4 ((θ y w+θ z v)(θ y δw+ θ z δv) + (θ z u + θ x w)(θ z δu + θ x δw)))dΩ;• M s,w ∝ Ω ρwδwdΩ;• N s,w ∝ Ω γ ∂p ∂z δwdΩ; • O s,w ∝ Ω γθ 3 pδwdΩ; 22 p( ∂δp ∂x θ 1 + ∂δp ∂y θ 2 + ∂δp ∂z θ 3 ) -( ∂p ∂x θ 1 + ∂p ∂y θ 2 + ∂p ∂z θ 3 )δpdΩ;

	• T s,w ∝ Ω φ(1 + Q R )( ∂p ∂z δw + p ∂δw ∂z )dΩ;
	• K f ∝ Ω	φ 2 ρ 22 ( ∂p ∂x	∂δp ∂x + ∂p ∂y	∂δp ∂y + ∂p ∂z	∂δp ∂z )dΩ;
	11 ∂z + ∂w ∂x ) + ( ∂v ∂δv ∂y )( ∂δv ∂v ∂y + C 12 ( ∂u ∂x + ∂w ∂z )) ∂δv ∂y + (C 11 -C 12 ) 1 4 (( ∂u ∂y + ∂v ∂x )( ∂δu ∂y + φ 2 • L f ∝ Ω ∂z + ∂δw ∂y )))dΩ; • H f ∝ Ω φ 2 ρ 22 pδpdΩ;
	• T s,v ∝ Ω φ(1 + Q R )( ∂p ∂y δv + p ∂δv ∂y )dΩ;
	• K s,w ∝ Ω ((C 11	∂w ∂z + C 12 ( ∂u ∂x + ∂v ∂y )) ∂δw ∂z + (C 11 -C 12 ) 1 4 (( ∂v ∂z + ∂w ∂y )( ∂δv ∂z +
	∂δw ∂y ) + ( ∂u ∂z + ∂w ∂x )( ∂δu ∂z + ∂δw ∂x )))dΩ;

ρ

  ∂δp ∂x u + ∂δp ∂y v + ∂δp ∂z w) + δp( ∂u ∂x + ∂v ∂y + ∂w ∂z ))dΩ. curves are computed for a PU 60 foam with a perfectly rigid cylindrical inclusion, with E = 10 15 Pa and structural loss factor = 0. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 Dispersion curve comparison with the reference (WFEM by Serra et al. [25]), and analytical model; real part of the wavenum-465 ber. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 Dispersion curve comparison with the reference (WFEM by Serra et al. [25]), and analytical model; imaginary part of the wavenumber. . . . . . . . . . . . . . . . . . . . . . . . . . . .

Acknowledgments

This project has received funding from the European Union's Horizon 2020 research and innovation program under the Marie Skłodowska-Curie grant agreement No. 675441.

Appendix A. Quantities defined in Biot model of poro-elasticity

RP -Q 2 ; • ρ11 , ρ12 and ρ22 are parameters depending on the nature and the geometry of the poro-elastic medium and the density of the fluid; in particular:

the bulk density of the fluid phase;

• ρ 1 is the bulk density of the solid phase;

• ρ a = φρ 0 (α ∞ -1) is an inertial coupling term;

List of Tables


L R Cunha

M Ouisse

D A Rade
email: rade@ita.br

Robust smart periodic truss

Smart and periodic structures have received the attention of researchers by virtue of their great potential. These structures have powerful properties like adaptiveness and the ability to operate as mechanical filters. Although, the presence of uncertainties must be taken into account to guarantee robustness. Thus, a finite element model is proposed to elucidate the importance of stochastic aspects and to present the concept of robust frequency bandgap. The smart part consists of piezoelectric actuators connected to resonant circuits in a tridimensional truss unit cell. The periodic part is the replication of this cell to assemble the final structure. Floquet/Bloch conditions are used to model the infinite representation. Then, a Monte Carlo Simulation is carried out and the bandgaps' bounds are analyzed considering frequency responses and dispersion diagrams. The goal being to evaluate the influence of uncertainties affecting the prediction of the attenuation zones. Likewise, the consequences of increasing the uncertainty level are evaluated.

Introduction

Wave propagation [START_REF] Graff | Wave Motion in Elastic Solids[END_REF] is one of the basis of vibration analysis. If this type of propagating energy arrives in a boundary condition or regardless of the impedance mismatch nature, there is a part that will be reflected and another transmitted. The presence of this kind of discontinuity, located periodically in a structure, can creates destructive interferences and, consequently, this movement is blocked and this energy is trapped inside this periodicity. Repetitive or periodic structures have being recently in focus by researchers of vibration control scientific community. It is explained by this unusual mechanical filter behavior for wave dispersion. A detailed review about these kinds of structures can be found in [START_REF] Mead | Wave Propagation in continuous periodic structures: Research contributions from Southampton, 1964-1995[END_REF] and some overviews about their future in [START_REF] Hussein | Dynamics of Phononic Materials and Structures: Historical Origins, Recent Progress, and Future Outlook[END_REF]. One of these perspectives is related to adaptiveness and tunable frequency bandgaps.

Smart material and structures capacity of self-sensing and self-changing to adapt to new conditions according to design requirements is their most important characteristic. Some fundamental aspects of these intriguing structures can be found in [START_REF] Leo | Engineering Analysis of Smart Material Systems[END_REF]. Therefore, periodic and smart structures seem to be a good combination to help the resolution of complex vibration problems. More information about periodic and smart structures is shown in [START_REF] Thorp | Attenuation and localization of wave propagation in rods with periodic shunted piezoelectric patches[END_REF] and [START_REF] Spadoni | Vibration and Wave Propagation Control of Plates with Periodic Arrays of Shunted Piezoelectric Patches[END_REF]. Piezoelectric (PZT) actuators [START_REF] Leo | Engineering Analysis of Smart Material Systems[END_REF][START_REF] Hagood | Damping of Structural Vibrations with Piezoelectric Materials and Passive Electrical Networks[END_REF][START_REF] Preumont | Vibration Control of Active Structures[END_REF] are used in the numerical example of this study. If this actuator is linked in series with a shunt circuit, it will dynamically behaves as adding a Dynamic Vibration Absorber (DVA). This kind of mechanism creates resonant bandgaps if placed periodically in a structure. More information about DVAs can be found in [START_REF] Hartog | Mechanical Vibrations[END_REF][START_REF] Korenev | Dynamic Vibration Absorbers, Theory and Technical Application[END_REF][START_REF] Rade | Optimisation of dynamic vibration absorbers over a frequency band[END_REF]. The principle behind this strategy is that the vibrational energy is transformed into electric energy through the direct piezoelectric effect and is transferred to the circuit where it is partially dissipated and/or dispersed. Among the types of used electric circuit, RL, known as resonant circuits, are considered as some of the most efficient [START_REF] Hagood | Damping of Structural Vibrations with Piezoelectric Materials and Passive Electrical Networks[END_REF][START_REF] Sales | Passive vibration control of flexible spacecraft using shunted piezoelectric transducers[END_REF]. Such circuits comprise an inductor and a resistor that are connected to the piezoelectric transducer that is assimilated to a capacitor, thus forming an RLC circuit. When coupled to a dynamic system, this device operates similarly to a dynamic vibration absorber (DVA). Distributing these devices may lead to multimodal control [START_REF] Dell'isola | Passive damping of beam vibrations through distributed electric networks and piezoelectric transducers: prototype design and experimental validation[END_REF].

Claeys describes the difference between placing a localized mass and a DVA in a periodic structure in [START_REF] Claeys | On the potential of tuned resonators to obtain lowfrequency vibrational stop bands in periodic panels[END_REF]. Normally, increasing the mass density in one specific degree of freedom (dof) in a unit cell may create a large band gap but with weak attenuation. Oppositely, by adding a spring-mass dof, it is possible to create resonant band gaps whose attenuation zones are narrower but with strong dispersion grace of lower frequency response amplitudes. Although, as it is not yet possible to add just an inertia effect on an electrical circuit with a PZT actuator, because only the actuator is already a spring equivalent, the resonant circuit is the best choice. The main advantage of using piezoelectric actuators, rather than DVAs, is the characteristic of no addition of significant mass to the main structure and the convenience of electronically tuning without changing geometric properties. Moreover, these circuits can be redesigned and unusual behaviors can be included. One example of this is the negative capacitance shunting [START_REF] Park | Vibration Control of Beams with Negative Capacitive Shunting of Interdigital Electrode Piezoceramics[END_REF], which aims at removing the intrinsic capacitive effect of the piezoelectric patch [START_REF] De Marneffe | Vibration damping with negative capacitance shunts: theory and experiment[END_REF]. This may be combined with resonant circuits [START_REF] Casadei | Vibration control of plates through hybrid configurations of periodic piezoelectric shunts[END_REF], opening the way to new strategies with wideband efficiency [START_REF] Lossouarn | Multimodal vibration damping of a beam with Periodic array of piezoelectric patches connected to a passive electrical network[END_REF][START_REF] Tateo | Experimental characterization of a bi-dimensional array of negative capacitance piezo-patches for vibroacoustic control[END_REF][START_REF] Tateo | Design variables for optimizing adaptive metacomposites made of shunted piezoelectric patches distribution[END_REF].

Signorelli in [START_REF] Signorelli | Wave Propagation, Power flow, and Resonance in a Truss Beam[END_REF] and [START_REF] Signorelli | Wave Propagation in Periodic Truss Structures[END_REF] shows the wave propagation behavior in truss structures by using beam finite element and the transfer matrix method. This sort of structure also has the weightlessness as its major characteristic. In virtue of this, the use of piezoelectric actuator in these lattices like structures seems to be a good compromise because it favors the lightness design requirements. Nevertheless, uncertainty and robustness analysis are rare in literature. Near-periodic structures, defects, impurities on periodic structures and the localization phenomena are well detailed by [START_REF] Wu | Vibration properties of piezoelectric square lattice structures[END_REF][START_REF] Yong | Propagation of decaying waves in periodic and piecewise periodic structures of finite length[END_REF][START_REF] Yong | Dynamics of Complex Truss-Type Structures[END_REF][START_REF] Mester | Periodic and Near-Periodic Structures[END_REF], but robustness analysis are scarce.

In this paper, one addresses the robustness of repetitive truss structures whose attenuation zones are created and passively controlled by using piezoelectric stack actuators associated with electrical shunt circuits [START_REF] Leo | Engineering Analysis of Smart Material Systems[END_REF][START_REF] Hagood | Damping of Structural Vibrations with Piezoelectric Materials and Passive Electrical Networks[END_REF]. The finite element method is used to find the movement equation for unit cells [START_REF] Orris | A finite element study of harmonic wave propagation in periodic structures[END_REF][START_REF] Mace | Modelling wave propagation in two-dimensional structures using finite element analysis[END_REF] whose three members contain a piezoelectric stack actuator connected to this circuit. Just like DVAs, shunt circuits must be tuned, which means that the values of their electric parameters must be precisely chosen for vibration attenuation in a narrow frequency band. However, the characteristic values of electronic components are prone to variability, due to manufacturing process and temperature, which can lead to mistuning and, consequently, decreasing of the control performance.

In this scenario, it becomes essential to evaluate the probability that the system will comply with the design requirements, given the probability density functions ascribed to the uncertain variables considered. For that, three different cases are analyzed based on typical operation and performance requirements. Uncertainties in the values of mass joints and inductance of the shunt circuit are modeled as Gaussian random variables. Following, a Monte Carlo Simulation [START_REF] Melchers | Structural Reliability Analysis and Prediction[END_REF][START_REF] Lemaire | Structural Reliability[END_REF][START_REF] Haldar | Probability, Reliability and Statistical Methods in Engineering Design[END_REF] is performed considering an infinite and a finite model by using the dispersion constants and frequency response function bandgaps data as output variables. The uncertainty level is varied by increasing the standard deviation of input variables as a percentage of their means. Subsequently, the results are presented and discussed.

Methodology

A unit cell is the smallest part of a periodic structure that is repeated side by side to assemble the complete structure. As long as these cells are identical, appropriate periodic boundary conditions can be taken into account and an infinite model can be created. The equation (1) presents the Floquet [START_REF] Floquet | Sur les équations différentielles linéaires à coefficients périodiques[END_REF] theorem for a 1D periodic structure:

( ) ( ) x X l x X λ = + , (1)
where, in an equivalent mechanical structure, X can be the displacement, x the position, l the cell length, or spatial period, and λ is the Floquet multiplier. For example, using a discrete model, the periodic boundary conditions are described in equations ( 2) and (3):

) ( ) 1 ( ) ( i L i L i R u e u u µ = = + , ( 2
) ) ( ) 1 ( ) ( i L i L i R f e f f µ - = - = + , (3)
for continuity and equilibrium to the left ( L ) and right ( R ) of the unit cell, respectively, where µ is the dispersion constant with δ and ε its real and imaginary components. These variables are the attenuation and phase constants, respectively [START_REF] Mead | Wave Propagation in continuous periodic structures: Research contributions from Southampton, 1964-1995[END_REF][START_REF] Hussein | Dynamics of Phononic Materials and Structures: Historical Origins, Recent Progress, and Future Outlook[END_REF]. For a purely complex µ , the waves propagate and the frequency zone is a pass band. If, in a given frequency band, all waves correspond to purely real values of µ , the waves are attenuated and this zone is a stop band or well kwon in periodic structure scientific community as bandgap. In this article, a simple truss unit cell, composed of 12 longitudinal rods, is shown in Figure 1. The embedded smart devices will be described in the next section. This structure has 3 boundary nodes on each side, the rods' properties are the (common) cross section area A , the Young's modulus E , density ρ and we denote i l the length of the cell i while l is the total length of a structure with n cells. For a sake of simplicity and low computational cost, the repetitive truss structure in Figure 1 is composed only by bar finite elements with 2 nodes and 3 degrees of freedom (dof) per node.

Smart Cell

One category of smart material are the piezoelectric material [START_REF] Leo | Engineering Analysis of Smart Material Systems[END_REF]. These materials produce an electrical output when a mechanical strain is imposed. This effect is called direct. Oppositely, a mechanical strain can also be induced by applying an electrical signal to this material and it is reciprocally called the inverse effect. These properties are expressed mathematically by the following equation:

            =       T E s d d S D E T 33 33 33 ε , (4) b) a)
where D and S are the electric displacement and the mechanical strain, E and T are the electrical field and mechanical stress. The material properties T ε (F/m), d33 (C/N or V/m) and E s (m²/N) are, respectively, the dielectric permittivity coefficient, piezoelectric strain coefficient and the mechanical compliance of the piezoelectric material of the active members. The indexes E and T represent the values in short-circuit (constant electric field) and open-circuit (varying electric field) with stress-free condition, while 33 is the direction along the axis of piezoelectric material polarization and it is one of the most common operating modes of a piezoelectric device. These characteristics are used to model the semi-active finite element present on this repetitive structure.

Periodic Cell with a resonant shunt circuit

In this section, the representation of the finite element model of the three-dimensional truss cell is depicted in Figure 2. This unit cell is composed by 6 nodes, 9 passive elements and 3 active elements. The nodes displacement are connected and the joint has a spherical function. Each of the active members is considered as being composed of a stack-type actuator placed between two passive segments. The actuator is assumed to be composed of piezoelectric discs poled in the axial direction (direction 3), wired in such a way that those discs are electrically connected in parallel. Application of voltages to the actuators will generate an internal longitudinal force enabling to attenuate transversal and longitudinal waves in the truss. As detailed in reference [START_REF] Leo | Engineering Analysis of Smart Material Systems[END_REF][START_REF] Hagood | Damping of Structural Vibrations with Piezoelectric Materials and Passive Electrical Networks[END_REF] the development of the finite element-based equations of motion of the electromechanical system can be made and it is equivalently expressed as:

( ) ( ) ( ) ( ) t t V t t F K KU U M = - + & & , (5) ( ) ( ) ( ) t Q t V t = + Γ U K ~, ( 6
)
where M is the mass matrix, K is the stiffness matrix, K ~ is the electromechanical coupling matrix, F is the vector of external loads, Γ is the matrix of dielectric permittivity, U is the vector of mechanical degrees of freedom,

) (t Q is the electric charge and ) (t V
is the voltage across the electrodes of the piezoelectric patches. The main structures is considered as conservative and the damping matrix C is zero. A RLC shunt circuit can be connect in series by adding the voltage condition as described in equation [START_REF] Hartog | Mechanical Vibrations[END_REF].

( ) ( ) ( ) t Q L t Q R t V & & & + = . (7)
Associating Eqs. ( 5), ( 6) and ( 7), the electromechanical equations of motion are found under the form:

( ) ( ) ( ) ( ) t t t t F Z K Z C Z M = + + & & & , (8)
where ( )

t Z , M , C and K are equal to ( ) ( ) { } T t Q t U , [ ] T L L Γ K 0 M ; - , [ ] T R R Γ K 0 0 ; - and [ ] T T I 0 K K - ;
~, respectively. An equivalent movement equation of a truss model can be found in [START_REF] Hagood | Damping of Structural Vibrations with Piezoelectric Materials and Passive Electrical Networks[END_REF].

In this work, Preumont shows a method by using the admittance of the shunt circuit. Similarly to what happens with a dynamic vibration absorber, the resonant shunt circuits must be tuned, which means that the values of their electrical resistance and inductance parameters must be accurately determined for the attenuation of vibrations of the host structure in a given range of frequencies. According to Hagood and Von Flotow [START_REF] Hagood | Damping of Structural Vibrations with Piezoelectric Materials and Passive Electrical Networks[END_REF], the electromechanical coupling coefficient plays the same role as the mass ratio in tuning a DVA. Indeed, this coefficient can be approximated as follows:

Tuning a resonant RLC circuit

( ) ( ) ( ) ( ) 2 2 2 2 E n E n D n ij K ω ω ω - =
( ) ( ) 2 1 2 ij ij opt K K r + = , ( ) E n S opt opt C r R ω = , 2 1 ij opt K + = δ and ( ) ( ) 1 2 - = S E n opt C L ω
. In this article, the electrical resistance is considered as 0 and no optimal values are used to tune. Thus, equation [START_REF] Rade | Optimisation of dynamic vibration absorbers over a frequency band[END_REF] shows how the inductance value is obtained for a specific tune frequency:

S n C L 2 1 ω = . (9)
However, in practical conditions, the values of the electric characteristics of the shunt circuit are inevitably affected by uncertainties resulting from material composition, manufacturing process and temperature variations. Such uncertainties can lead to mistuning of the shunt circuit and, as a result, deterioration of the performance of the damping/resonating device. In the case of periodic structures, it changes the bandgap ij K USD -UNCERTAINTY IN DESIGN AND MANUFACTURING bounds and, for as strong impedance mismatch in mechanical and electrical domains can creates natural frequencies inside that increases the amplitude and, consequently, decrease the attenuation zone efficacy.

Frequency Response Function and Dispersion Diagrams

To give a sense on the influence of the joint mass and shunt inductance, 8 scenarios were compared by observing the dispersion diagrams and frequency responses and they are illustrated in Figure 3. The equation 8 can be solved considering a harmonic solution to obtain the frequency response function (FRF) in equation [START_REF] Claeys | On the potential of tuned resonators to obtain lowfrequency vibrational stop bands in periodic panels[END_REF]. Moreover, the equation ( 11) is achieved by solving equation ( 10) considering periodic boundary conditions of equations ( 2) and (3).

( ) ( ) 1 2 - + + - = K C M H ω ω ω j , ( 10
) [ ] 0 ) ( ) ( ) ( ) ( ) ( 2 ) ( ( = - + µ µ ω ω µ r r r r j U M C K . (11)
) c) d) e) f) g) h)
bounds. It is important to mention that as the model of one cell has 18 and 21 d.o.f., considering reduced matrices (r), the same number of mode branches can be found. As long as it does not represent all the dof that one finite structure can have, the area after the last branch is not considered as an attenuation zone and they are not exposed in these diagrams. The Figure 4 shows their equivalent respective frequency response by considering a transversal (x or y) and longitudinal (z) excitation on the first node cell and observing the response of equivalent node of last cell. The number of cells considered to model the finite structure is equal to 10. This choice was done concerning the lowest number of cells that the bandgap is easily depicted in frequency response.

The added mass ( a m ) is equal to 0.3 kg and the value of the inductance ( L ) to tune to a 7000 Hz, considering equation ( 9) is 0.1084 H.

It is possible to observe that infinite and finite models have the same band gaps and it is verified by comparing Figures 3 and4. Some important aspects can be observed. Firstly, these models without added joint mass and piezoelectric actuator have no band gaps (Figure 3.a). The same can be observed by adding mass in all joint or adding actuators as shown in Figures 3 b) and3 e). The fact of adding mass in joints, or considering the actuators in half cell, creates the attenuation zones that can be analyzed in Figures 3 c) and 3 f). By comparing these diagrams, a larger one is found in the case of added mass and a deeper one with added mass and PZT in shunt circuit. It is similar to the added mass and mass-spring effects in [START_REF] Claeys | On the potential of tuned resonators to obtain lowfrequency vibrational stop bands in periodic panels[END_REF].

Monte Carlo Simulation and Localization Phenomena

The Monte Carlo Simulation (MCS) [33, 34 and 35] is a method that consists in sampling input random variables, defined as a probability density function, to obtain predefined outputs. As mentioned in reference [START_REF] Melchers | Structural Reliability Analysis and Prediction[END_REF], convergence can be reached with approximately 95% confidence, by using from 10,000 to 20,000 samples, depending on the function being evaluated. This convergence is characterized by the stabilization of mean and standard deviation values. For problems featuring low standard deviations, very large numbers of samples are necessary. Optimized sampling methods can be used to reduce the computation effort necessary for MCS convergence. In this work, the sampling technique referred as Latin Hypercube Sampling (LHS) was used. More details about this technique can be found in reference [START_REF] Choi | Reliability-Based Structural Design[END_REF].

A specific probability density function should be chosen to represent the stochastic variable considering the input variable characteristics. In this article, for the sake of simplicity and easy comprehension, the normal distribution is used. Considering a Gaussian distribution, if the standard deviation has high values, the gap between the properties of neighbor cells can creates the localization phenomena. The localization phenomena can be briefly explained as normal mode that is confined in a specific place of the structure. As the number of cells tends to infinity, the influence of this mode tends to zero. Consequently, it does not appear in the infinite model. It is explained due to the stronger impedance mismatch compared to other cells and, grace of it, the localization phenomenon appears. These phenomena can be considered as valuables as long as they reduce the variation of bandgap borders.

In other words, this kind of unpredictable event in a Monte Carlo Sampling adds robustness to the bandgap behavior in meanings of borders variation but it becomes hard to observe the amplitude on FRF inside bandgaps' zones and to identify the effective attenuation zones. Consequently, a robust interval confidence envelope will not take into account the localization phenomena since the worst configuration may takes part in any cell in a random analysis.

Numerical Example

In practice, it is not easy to design locally or to modulate periodically the inertia and the stiffness of a mechanical structure without weakening it. Normally, the properties as Young modulus and the density do not vary considerably if the material comes from same pack. Therefore, the principal incoming of uncertainties are the imperfections or defects caused by the manufacturing process. Accordingly, Figures 3 (c), (f) and (h) illustrate the 3 situations of interest in which the added joint mass and inductance values are stochastically analyzed. The Table 1 shows the properties of the truss unit cell for the numerical example.

Probabilistic Analysis

A periodic structure is composed by quasi-identical cells. This "quasi" is for the reason that there is no perfectly identical structures even they were designed to be. Thus, a robustness analysis is employed for this periodic structure with the goal of finding a robust bandgap. For this case, to get closer to reality, the variables of each cell should be independents. For the sake of comparison, infinite and finite models are used in this probabilistic analysis. To realize this investigation, for each standard deviation of random input variable (RV), a convergence of Monte Carlo Simulation is achieved and the outputs are obtained. In this article, the band gap edges are the outputs. Figures 3 show an illustration of these variables where the red dots are the bounds for infinite structure and a specified magnitude, different for each case, is used as a threshold for obtaining the bandgaps limits for a finite structure. The stochastic properties of the input random are shown in Table 2.

Variables Distribution Mean ) (µ Standard Deviation ) (σ 1 X ] [ kg m a Normal (Gaussian) 3 . 0 1 X µ γ × 2 X ] [H L Normal (Gaussian) 1084 . 0 2 X µ γ ×
Table 2: Probabilistic variables and their distributions.

Two different cases are considered for each of the 3 scenarios of figures 3 c), f) and g):

• the uncertainty variable (1 RV) in the unit cell is repeated for each cell and, consequently, the structure still periodic (infinite model);

• all cells are independents (10 RV) from the uncertainty point of view (the structure is not perfectly periodic anymore, meaning quasiperiodicity, and is strictly finite).

For both cases, FE models were considered. The equation [START_REF] Signorelli | Wave Propagation in Periodic Truss Structures[END_REF] presents the stochastic coefficient γ . It permits the variation of uncertainty level, given by the standard deviation, which is obtained from a percentage of probability distribution mean, according to:

n n X X = µ γ σ × . ( 12
)
This parameter varies from 2 % to 20 % with a step of 2 % for this numerical example. Then, for each considered value, a Monte Carlo Simulation with Latin Hypercube is performed until its convergence (number of samples around 5000). Considering a Normal approach for 95% of confidence level envelope, the output standard deviation value must multiplied by 1.96. It means that the probability of the output mean of the output samples has 95% of chance to be in this interval.

Numerical Results

The Figure 6 shows the envelope of the frequency band gaps for three considered scenarios. Considering a normal approach, this picture shows the means and the 95% confidence intervals. It is possible to observe that the envelopes widths increase according to the increasing of uncertainty level. In general, the envelopes are narrower for the models with finite structures (10 RV), approximately less than a half of the bandgaps obtained with infinite structures (1 RV).

If an uncertainty level of 10 % is considered for the random variables, the dark gray areas represents the robust bandgaps (bandgap almost certain considering these uncertainties).

The localization phenomena perturbed the result in Figure 6 e) and equivalent results were not found as related to Figures 6 d) and 6 f). This proves that this phenomenon is strongly sensible to low variations of inductance values in a shunt circuit and an erroneous result is obtained.

c) a) b) d) e) f)
The Figure 7 elucidates the concept of robust band gap by situating the classical band gaps and robust band gaps inside the frequency response and dispersion diagrams, for infinite and finite models, respectively. The light red areas (envelopes) in Figure 6 indicate that the attenuation zone calculated with a single random variable is narrower than the dark gray, which was obtained with ten random variables for gamma, considering γ equals to 10%. Thus, when one uncertainty is considered for each cell of a finite structure, which corresponds to the quasi-real case, the effects on the responses are compensated and the structure are intrinsically more robust. This envelope does not take into account the localization phenomena. However, it can appears just in the case of finite model simulation with 10 RV. For high gamma values, the chances of appearing this situation is higher.

It is also possible to observe in figures 7 (b) and 7(d) that there is a branch correspondent to a longitudinal wave that inhibits the appearance of two others bandgaps which can be verified by the appearance of on figures 7 (a) and 7 (c).

Conclusion

A low computational cost three-dimensional periodic structure truss was used as a numerical example. Firstly, solutions for finite and infinite structures were found likewise their frequency response and dispersion diagrams, respectively. This structure was discretized by a bar finite element with 3 degrees of freedom per node. A number of cells equals to ten has been selected for the finite structure. This choice was made in order to obtain credible location of the attenuation frequency zones without substantially increasing the computational cost.

The localization phenomenon was identified related to a place with strong impedance mismatch. This important remark was clarified and it explains the behavior of uncertainties in this kind of repetitive structure. A relevant impedance mismatch, which highly differs from the other cells, localized in some place of the repetitive structure can reveal interesting properties.

An uncertainty analysis was performed using the finite and infinite models of this structure. The added joint mass and inductance values were purposely chosen as random variables with a Gaussian probability density function owning a mean and standard deviation as stochastic properties.

The choice of this probability distribution is explained due to its facility of understanding. Although the mass and inductance values may achieve negative values for this model, a strongly high value of standard deviation is necessary.

A coefficient of variation, namely gamma, was specified to vary the standard deviation as a percentage of the mean value. Numerous Monte Carlo Simulations with Latin Hypercube were conducted until their convergences for each value of the coefficient mentioned before.

The results for infinite and finite models are different. The results for infinite structures are not completely reliable since the same uncertainty is repeated indefinitely. The results for finite structures are more reliable since the uncertainty is considered ten times for the finite structure with ten cells.

The confidence envelope for the complete structure with ten stochastic variables is narrower than those found with one random variable. Besides, they are more reliable due to the model being closer to reality for the same case. The fact of considering ten random variables, instead of one, proves that the uncertainties self-compensate and the periodic structures are intrinsically tough.

Finally, a more accurate and reliable robust smart bandgap was obtained.

Figure 1 :

 1 Figure 1: a) 3D truss unit cell and b) its finite model.

Figure 2 :×Table 1 :

 21 Figure 2: a) Smart unit cell and b) its piezoelectric stack actuator described as a finite element model. As indicated in Figure 2, E (N/m 2 ), A (m 2 ) and ρ (kg/m) are the Young's modulus, cross-section area and linear mass density of the passive members, respectively. The others properties, E Y 33 (N/m 2 ), Ap (m 2 ), tp (m) and ρp (kg/m) are the Young's modulus, cross-section area, thickness of the piezoelectric discs and linear mass density of the piezoelectric material of the active members. The physical and geometrical properties of the main structure and the actuator are provided in Table 1, in which subscripts b and p indicate the properties related to the metallic and piezoelectric material.

Figure 2 b

 2 Figure 2 b) illustrates a piezoelectric transducer connected to a resonant (RL) shunt circuit and bonded to a host vibrating structure.Similarly to what happens with a dynamic vibration absorber, the resonant shunt circuits must be tuned, which means that the values of their electrical resistance and inductance parameters must be accurately determined for the attenuation of vibrations of the host structure in a given range of

ω

  are the n -th natural frequencies in open and closed circuit, respectively. Knowing the value of ij K , one can calculate the optimum values of resistance opt R

Figure 3 :

 3 Figure 3: Considered unit cells a), b), c), d), e), f), g) and h) with their respective dispersion diagrams.

Figures 3 a

 3 Figures 3 a), b), c), d), e), f), g) and h) shows the added mass and actuator position with their correspondent imaginary part of dispersion constants. Transversal and longitudinal wave's modes are not sorted and the bandgaps highlighted in gray comprehend these two kind of waves. The red dots indicate the bandgap

Figure 4 :

 4 Figure 4: a), b) and c) frequency responses of these periodic finite structures with 10 cells by exciting the first node on direction x, y and z and observing the correspondent dofs on the last unit cell.

Figures 3 d

 3 Figures 3 d), g) and h) are the combination of the others. The bandgap width is increased by combining added mass and actuator in scenarios of figures 3 g) and 3 h). Uncertainties can drastically change these bandgaps behavior and it is elucidate in the following.

Figure 5 :

 5 Figure 5: a) FRFs with localization phenomena L1, L2, L3 and L4 and b) their correspondent mode shape.

  Figures 6 a), b) and c) for infinite model obtained from imaginary part of propagation constant considering with 1 RV (strictly periodic structure) and Figures d), e) and f) for finite model obtained from frequency responses with 10 RV (quasi-periodic structure).

Figure 6 :

 6 Figure 6: Bandgap bounds envelopes relatives to (c), (f) and (h) on figure 3, with: a) mass as RV, b) and c) inductance as RV considering an infinite model (1 RV) with their correspondent d), e) and f) to finite model results with 10 RV.

Figure 7 :

 7 Figure 7: a), b), c) and d) Robust band gaps for a) finite and b) infinite structures.

Acknowledgements

The authors are grateful to the Brazilian National Research Council -CNPq and to the National Institute of Science and Technology of Smart Structures in Engineering -INCT-EIE for the financial support to their research work. This work was realized in the context of ANR-12-JS09-008-Covia project, in cooperation with the Labex ACTION (ANR-11-LabX-0001-01).


An ancient age of open-canopy landscapes in northern Madagascar? Evidence from the population genetic structure of a forest tree

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

An ancient age of open-canopy landscapes in northern Madagascar? Evidence from the population genetic structure of a forest tree

Miguel de Navascués based on reviews by Yurena Arjona and Katharina Budde

A recommendation of:

How ancient forest fragmentation and riparian connectivity generate high levels of genetic diversity in a micro-endemic Malagasy tree Jordi Salmona, Axel Dresen, Anicet E. Ranaivoson, Sophie Manzi, Barbara Le Pors, Cynthia Hong-Wa, Jacqueline Razanatsoa, Nicole V. Andriaholinirina, Solofonirina Rasoloharijaona, Marie-Elodie Vavitsara, Guillaume Besnard(2021), bioRxiv, 2020.11.25.394544, ver. 4

Recommendation

We currently live in the Anthropocene, the geological age characterized by a profound impact of human populations in the ecosystems and the environment. While there is little doubt about the action of humans in the shaping of present landscapes, it can be difficult to determine what the state of those landscapes was before humans started to modify them. This is the case of the Madagascar grasslands, whose origins have been debated with arguments proposing them either as anthropogenic, created with the arrival of humans around 2000BP, or as ancient features of the natural landscape with a forest fragmentation process due to environmental changes pre-dating human arrival [e.g. [START_REF] Quéméré | Genetic data suggest a natural prehuman origin of open habitats in northern Madagascar and question the deforestation narrative in this region[END_REF][START_REF] Joseph | Madagascan highlands: originally woodland and forest containing endemic grasses, not grazing-adapted grassland[END_REF]. One way to clarify this question is through the genetic study of native species. Population continuity and fragmentation along time shape the structure of the genetic diversity in space. Species living in a uniform continuous habitat are expected to show genetic structuring determined only by geographical distance. Recent changes of the habitat can take many generations to reshape that genetic structure [START_REF] Landguth | Quantifying the lag time to detect barriers in landscape genetics[END_REF]. Thus, we expect genetic structure to reflect ancient features of the landscape.

The work by Jordi Salmona and collaborators [START_REF] Salmona | How ancient forest fragmentation and riparian connectivity generate high levels of genetic diversity in a microendemic Malagasy tree[END_REF] studies the factors determining the population genetic structure of the Malagasy spiny olive (Noronhia spinifolia). This narrow endemic species is distributed in the

Open Access

Published: 12 November 2021 discontinuous forest patches of the Loky-Manambato region (northern Madagascar). Jordi Salmona and collaborators genotyped 72 individuals distributed across the species distribution with restriction associated DNA sequencing and organelle microsatellite markers. Then, they studied the population genetic structure of the species. Using isolation-by-resistance models [START_REF] Mcrae | Isolation by resistance[END_REF], they tested the influence of several landscape features (forest cover, roads, rivers, slope, etc.) on the connectivity between populations. Maternally inherited loci (chloroplast and mitochondria) and bi-parentally inherited loci (nuclear), were analysed separately in an attempt to identify the role of pollen and seed dispersal in the connectivity of populations.

Despite the small distribution of the species, Jordi Salmona and collaborators [START_REF] Salmona | How ancient forest fragmentation and riparian connectivity generate high levels of genetic diversity in a microendemic Malagasy tree[END_REF] found remarkable levels of genetic diversity. The spatial structure of this diversity was found to be mainly explained by the forest cover of the landscape, suggesting that the landscape has been composed by patches of forests and grasslands for a long time. The main role of forest cover for the connectivity among populations also highlights the importance of riparian forest as dispersal corridors. Finally, differences between organelle and nuclear markers were not enough to establish any strong conclusion about the differences between pollen and seed dispersal.

The results presented by Jordi Salmona and collaborators [START_REF] Salmona | How ancient forest fragmentation and riparian connectivity generate high levels of genetic diversity in a microendemic Malagasy tree[END_REF] contribute to the understanding of the history and ecology of understudied Madagascar ecosystems. Previous population genetic studies in some forest-dwelling mammals have been interpreted as supporting an old age for the fragmented landscapes in northern Madagascar [e.g. [START_REF] Quéméré | Genetic data suggest a natural prehuman origin of open habitats in northern Madagascar and question the deforestation narrative in this region[END_REF][START_REF] Rakotoarisoa | Late Quaternary climatic vegetational shifts in an ecological transition zone of northern Madagascar: insights from[END_REF]. To my knowledge, this is the first study on a tree species. While this work might not completely settle the debate, it emphasizes the importance of studying a diversity of species to understand the biogeographic dynamics of a region.

Reviews

Author's Reply

Download author's replyDownload tracked changes file

Decision by Miguel de Navascués

Dear Jordi Salmona, Guillaume Besnard and coauthors, Thank you for the thorough revision of your manuscript. You have satisfactorily addressed the concerns raised by the reviewers.

However, in the evaluation of your reply it has come to my attention an important problem that was overlooked in the previous round of review. At different parts of the manuscript you refer to "genetic diversity", "heterozygosity", "Ho", "He", and "% of polymorphic sites" regarding nuclear RADseq data. After carefully reading all the relevant text and supplementary materials I conclude that in all of these instances you are referring to the same measurement: the individual heterozygosity per base pair, that is, for each individual the number of heterozygous sites over the total of sequenced data (after application of relevant filters). First, the text needs to state this more clearly and use a consistent name for that measurement throughout the entire work. Second, individual heterozygosity is a poor measure of genetic diversity. For instance, it can be affected by the reproductive system (i.e. clonal or selfing reproduction). Expected heterozygosity (per bp, i.e. π) calculated among individuals of the same forest patch (for comparisons among populations) or over all individuals (for comparisons with other species) is a more appropriate measure and what most readers expect as a measure for genetic diversity, please add this or an equivalent measure to your work. Since the levels of genetic diversity are one of the main results from the work this problem needs to be addressed before the work can be recommended.

In addition to this main point there are some other minors points that I would like you to consider:

In the discussion you suggest that current landscape structure must have been in place "at least 40 generations" corresponding to "at least ~800 years". I think there is little justification for those numbers even in the light of the results by [START_REF] Landguth | Quantifying the lag time to detect barriers in landscape genetics[END_REF] and Mona et al. (2014). The dynamics of how pattern change through time in non-equilibrium models depends on many unknown parameters (effective population sizes, migration rate…), giving such precise numbers seems misleading. While I agree with the general message and I understand that this sentence was intended to address a point Yurena Arjona raised in her review, I think it would be more appropriate to discuss that the number of generations for the period with data on forest cover (1953-present) represent less that 5 generations, which seems too low to erase the signal of previous population structure (or lack of).

Please verify that the description of header for Table S1 is complete and accurate.

In line 617, I suggest to write "geographical distance" instead of "Euclidean distance".

Last, I would like also to suggest you considering an alternative wording for fragmentation/forest fragments throughout the text. Since one of the conclusion of the work is that the landscape structure with open areas is a relatively old feature, it is unclear to me whether there was ever an ancient continuous forest in the region subsequently fragmented (presumably by humans). Thus forest patches might be just that, "patches" and not "fragments".

Evaluation round #1

Author's Reply

Download author's replyDownload tracked changes file

Decision by Miguel de Navascués

Dear Jordi Salmona, Guillaume Besnard and coauthors, I have now received the reviews for your manuscript. While they are in general positive, both point out few concerns regarding the methods, conclusions and writing that require further work on the manuscript. Please consider all their comments carefully before submitting a revised version for evaluation. I believe most of these points can be addressed by working on the text and no further analysis will be necessary. However, pay special attention to properly justify the absence of new analysis regarding the concerns on the methodology.

Both reviewers highlight the lack of explanation for the reduction of samples from sampling to genotyping. Please revise your text so it is completely clear for the reader why there are such differences. In addition, clarify how the survey of Noronhia was conducted and if N. spinifolia was identified in the field. Were all Noronhia spp. sampled? In the main text it is unclear if the 220 samples are from Noronhia spp. or from N. spinifolia only (in the supplementary text it is clearly stated, but I think the reader should not need to go there for this basic information).

Another question that seems to require a more clear explanation is the discussion of the age of forest fragmentation in the area. I believe that some of the confusion stems from the better fit of the genetic distances with data on recent vegetation cover (2000's) but a conclusion that forest fragmentation structuring the genetic diversity of N. spinifolia must be ancient. It might not be clear for the reader why data on older vegetation cover does not explain better the genetic structure of N. spinifolia. I also have some points I would like you to address:

1) Regarding the discussion on the role forest cover for seed dispersal (lines 497-509), it is not obvious to me why patterns of mtDNA and cpDNA haplotype networks would be incompatible with the signal found by IBR analysis on organelle data. I would appreciate if you could develop this point.

2) Supplementary data availability: In addition to making available the missing Table S1, please consider making supplementary tables S6, S7 and S8 available through reputed repositories (such as institutional repositories, Zenodo, DataDryad, supplementary material in bioRxiv...) that can assure long term accessibility to those files (instead of Google Drive/Docs service). Also, I could not find reference PRJNA632767 in the Sequence Read Archive (SRA) NCBI database, make sure the reference is correct and it is accessible.

3) line 59: "tropical biodiversity hotspots" instead of "tropical hotspots" S5) should be included in the main text when and in the supplementary materials (Method S6) when they are first mentioned. 6) line 466: instead of "prehistoric", which strictly refers to the time period before the invention of writing, use "non recent" or some other wording. S2: What is the difference between columns "Allele number" and "Number of alleles"?

7) Table

Best regards, Miguel Navascués

Additional requirements of the managing board:

We would like to receive your revision within 2 months. If you need more time, just tell us. As indicated in the 'How does it work?' section and in the code of conduct, please make sure that: -Data are available to readers, either in the text or through an open data repository such as Zenodo (free), Dryad or some other institutional repository. Data must be reusable, thus metadata or accompanying text must carefully describe the data.

-Details on quantitative analyses (e.g., data treatment and statistical scripts in R, bioinformatic pipeline scripts, etc.) and details concerning simulations (scripts, codes) are available to readers in the text, as appendices, or through an open data repository, such as Zenodo, Dryad or some other institutional repository. The scripts or codes must be carefully described so that they can be reused. -Details on experimental procedures are available to readers in the text or as appendices.

-Authors have no financial conflict of interest relating to the article. The article must contain a "Conflict of interest disclosure" paragraph before the reference section containing this sentence: "The authors of this preprint declare that they have no financial conflict of interest with the content of this article." If appropriate, this disclosure may be completed by a sentence indicating that some of the authors are PCI recommenders: "XXX is one of the PCI XXX recommenders."

Reviewed by Katharina Budde, 06 Jan 2021

The Manuscript "How ancient forest fragmentation and riparian connectivity generate high levels of genetic diversity in a micro-endemic Malagasy tree" written by Salmona et al. deals with an interesting topic in a very interesting and understudied ecosystem. The authors used chloroplast SSRs and RADseq loci from the chloroplast, mitochondrial and nuclear genomes to study the genetic diversity and genetic structure of a narrow endemic tree species, Noronhia spinifolia in Madagascar. The study presents interesting data and results and is well written, nevertheless I have some concerns that I will lay out in detail below and that I hope will be helpful for the authors to improve the manuscript.

Study species: In this section I am missing some more information e.g. about the flowering and fruit phenology of the study species. Is this a monoecious or dioecious species? Are the flowers hermaphroditic or are female and male flower parts separated on the same individual? Please also add information about the possible pollen dispersal and seed dispersal mechanisms. If no literature is available for this particular species the flower and fruit morphology can already give an idea.

Sample numbers: 220 leaf samples collected, DNA extracted from 137 samples and after quality control only 72 samples were used for genotyping and sequencing. Why could only so few samples be used at the end?

The methods appear mostly sound but sometimes I am missing details:

A minimum coverage of 4x for the RADseq seems very low to me and I am wondering why the authors chose such a low threshold. Anyway, it is good that the authors used both called genotypes and genotype likelihoods! Could a higher threshold for the minimum coverage have improved the genotype calling in STACKs? What about linkage disequilibrium at nuclear loci? I could only find information for the organellar loci. And which markers were used for NgsAdmix and ADMIXTURE? Did the authors take care to remove loci showing high LD? I am missing Table S1 (if there is a link in blue, it did not work for me).

The authors mention that they focus on the mature trees. Why did they make this choice? For the future, I think a comparison between mature and young trees might be interesting to evaluate the effect of more recent forest fragmentation.

l. 289 I suggest to substitute "diversity" with "structure" as the authors do not test the effect of environmental factors on genetic diversity estimates but rather the effect on the genetic distance or relationship matrices. Is this correct? It might be helpful to actually mention which measures were tested for correlation. In the following paragraph the authors only mention "genetic matrices". l. 447-450 Is it meaningful to directly compare the genetic diversity levels between species with such different life history traits as a tree species and mammal species (here mouse lemurs and plains zebras)? Also, why did the authors pick these species and not others? Tree species are characterized by harboring huge amounts of standing genetic variation (e.g. Petit & Hampe 2006, https://doi.org/10.1146/annurev.ecolsys.37.091305.110215). The authors have already illustrated well the big difference in diversity levels in the sentences before by comparing with other tree species. Therefore, I would suggest to remove these sentences.

The authors explain their results by a long-standing forest fragmentation in the LM region (l. 521) which is in line with other studies suggesting that forests were already fragmented before intensified human use in some regions of Madagascar. Do they have any idea if the areas surrounding the current forest patches might have been historically covered by a different vegetation type where N. spinifolia was not present? Are there other environmental factors that are correlated with forest cover, e.g. soil type that might be less suitable for N. spinifolia? Another or additional possible explanation for the high genetic diversity and structuring, as also mentioned by the authors, might be hybridization with closely related species and I very much agree with this. Could it also be a case of cryptic radiation as it has been reported for several groups on Madagascar (but see also Pillon et al. 2014; doi: 10.1111/nph.12677)? Noronhia habors many micro-endemic species in Madagascar and I agree that a more exhaustive sampling of all possible closely related species and taking into account morphological and ecological features in future analyses might shed further light on hybridization, ecotype differentiation and speciation.

  of the preprint: https://www.biorxiv.org/content/10.1101/2020.11.25.394544v1 Version of the preprint:

4 )

 4 line 72: Specify what 40-50% are referring to (area?) 5) lines 324-5: It is unclear what the values "m=4, M=5, N=8" are. Please clarify: a brief description/definition (such as in figure

  l. 391What is the null IBD layer? Do you mean a null model only taking into account IBD? Please clarify.

Figure

  Figure 4C Is "Percent of tree cover geographic distance" referring just to a pair-wise distance or is it a conductance measure? Please clarify! I find it sometimes difficult to follow the authors to which kind of estimate (pair-wise distance, resistance etc.) they are referring. Please revise the manuscript accordingly. l. 441" …polymorphic cpSSRs (35)…" what exactly is the 35 referring to?

  l. 231 …site frequency spectrum/ spectra?.. l. 244 Do you mean pair-wise estimates? l. 252 remove "the" before Nei´s and Yang´s l. 267 remove s at the end of dispersals l. 319 and 320 I would suggest using no LD instead of null LD l. 377 from and to? l. 383 are shown l. 429 I suggest to change genetic diversity to genetic structure l. 443 Add an s to level l. 473 …accelerated evolutionary rate… l. 531 over a larger timeframe l. 846-851 Prunier et al 2015 is listed twice in the references Reviewed by Yurena Arjona, 04 Jan 2021 Download the review

Toggle reviews Evaluation round #2 26 May 2021 DOI or URL of the preprint: https://www.biorxiv.org/content/10.1101/2020.11.25.394544v2 Version of the preprint: https://www.biorxiv.org/content/10.1101/2020.11.25.394544v2