Marc-Antoine Campana

Periodic

M Pezerat

Charles Professeur

Président M Michon

Guilhem Professeur

Supaero Rapporteur

M Le Bot

M Scarpa Fabrizio

Periodic inserts in auxetic media for a vibroacoustic control

FIGURES 5.2 Nonlinear forces and moments related to out-of-plane wave propagation (u z , θ x , θ y ) of to the unit cell described 4.31 for ( f req, µ) = (27.5, 0.65) (green case) and ( f req, µ) = (222.3, 1.86) (blue case) . . . . . . . . . . . . . . . . . . .

Inserts en milieu périodique pour le contrôle vibroacoustique

DEDICATION AND ACKNOWLEDGEMENTS

T he ambition of working in the domain of research which led me to conduct this PhD started almost five years ago. At that time I was a student at the engineering school ENSMM and all the studies concerning the field or mathematics, physics and mechanics fields fascinated me at a point that I really wanted to work deepen on these problematics during my early working career. Everything started in February 2015 when my first research internship in acoustics began at the lab associated to my engineering school, FEMTO-ST. I am very thankful to Morvan and Emeline, who understood my interest for the research at this period and let me have this first experience as an early researcher working with them, experience that I enjoyed so much that I decided to reiterate the year after for my master thesis, still working with them to finally having the great opportunity of starting the PhD in this very challenging project. Thank you for the trust you placed on me, making me able to have three years of rich experience of life and work through this PhD and for all the help and support you provided during these three years.

Many thanks to my supervisor Fabrizio, who always has been present to all the help me in a technical, experience or human point of view. It has been an honor to work with a so acknowledge professor as well as a very human being person. Thanks to my thesis reporters Professors Le Bot and Michon for all the precious advice they give me after the two oral presentation in order to improve the quality of the manuscript. Thanks to my colleagues from Bristol, Simone and Rita with who we started to discover the PhD world together, making me feel very confident to start this experience in the best possible conditions. Thanks to my colleagues from Besançon, Dario Safi and Taoufik making me even more enjoyable the few months I passed in France. Thanks to my colleague Sepide for the 6 months she passed in Bristol for her secondment. Thanks to all my other PhD colleagues from the Viper project that I could meet for several occasions during these three years within project meetings or conferences. Thanks to Massimo who received Simone and me for our 3 months placement in GeorgiaTech, allowing me to discover a country I never had the occasion to go before and the chance to have work in such a great university. Thanks to everyone who was involved in the project, including all the professors and Christophe for all the help and support they provided.

Thanks to all the people I met during the time I passed in Bristol, and more specifically to all my French friends I had the luck to meet here, helping me to have so much great time and remembering a piece of homeland, my experience in UK would have been really less enjoyable if I didn't meet you.

And last but not least, thank you for the support of all my family, sisters, my brother and my parents.

Thesis objectives

Periodic structures such as honeycombs have been widely used in aerospace industry [START_REF] Wittenauer | Structural honeycomb materials for advanced aerospace designs[END_REF], having the advantage to be lighter than commonly used panels providing an overall decent stiffness over the mass ratio ensuring a good mechanical strength for a lower cost than usual structures used traditionally. Later, an alternative of honeycomb structures has been investigated called re-entrant structures [START_REF] Spagnoli | On the auxetic behaviour of metamaterials with re-entrant cell structures[END_REF] presenting other interesting features in the structural design. Contrary to honeycombs presenting an anticlastic behaviour when solicited in out-of-plane bending, reentrant structures have the particularity to be auxetic, meaning that the equivalent homogeneous 1 CHAPTER 1. INTRODUCTION structure acts as it has a negative Poisson's ratio. This property provides a synclastic behaviour to the re-entrant structure, allowing the construction of slightly curved panels while keeping the interesting mechanical properties of an honeycomb.

Over the structural properties, such periodic systems give also interesting dynamic properties concerning the wave propagation. One of the way to study wave propagation in a periodic structure is the usage of the Floquet-Bloch theorem applied to one of the unit cells of the structure, giving the relation between the wave number and the frequency of an assumed infinite periodic system, called dispersion relation. This dispersion relation is graphically represented by a dispersion curve, allowing us to identify some frequency areas where waves are strongly attenuated ; those areas are called stop bands or bandgaps and are one of the main feature we are looking to improve the dynamic of our system. Honeycombs-like structures present a very high stiffness over the mass ratio, meaning that their impact on the dynamic behaviour will mainly be identified in high frequencies as well as for stop band location, which is not compatible with the expectations of aerospace industry in which low frequency areas are targeted.

In order to improve the wave propagation in low frequency, one of the solution presented in this work is the use of inclusions such as resonators in the honeycomb structure. This solution has the advantage to not modify the hosting structure in order to keep its mechanical properties, giving more freedom to the design of the resonator to control the stop band in a targeted frequency area. Also, the resonators act as a solution for passive control, meaning that the structure and the physics of the resonator itself is sufficient to have an impact on the wave propagation behaviour, without adding any electronic device for instance. In principle, the resonators are acting as absorbers for the hosting structure, meaning that almost all the energy coming from the waves will be concentrated on them as soon as they become evanescent. This high amount of energy can sometimes cause great deformations and displacements to the resonators, meaning that the hypothesis of geometrical linearity might not be respected any more. In that case, it is necessary to analyse the dynamic behaviour of the system in nonlinear regime, developing tools to identify the impact on wave propagation and therefore on the dispersion relation.

The first part of the thesis is focused on a literature review in periodic structures for linear and nonlinear regime. The second part will show examples of linear periodic structures analyses using the finite element method combined with the Floquet-Bloch Theorem to a re-entrant unit cell including resonators to see their impact to the stop band opening in the dispersion curves for different configurations. The nonlinear results applied to analytical system is exhibited in the third part, as well as first applications to finite element models, using both Lindstet-Poincaré and the inverse method developed to confirm the theoretical approach. Finally, the fifth and sixth chapter will open a discussion about the performance achieved with resonators in complex structures such as the re-entrant configuration, and the limitations of the nonlinear methods by suggesting some ways of improving the, in future works.

C H

A P T E R 2 LITERATURE REVIEW

Periodic Structures

Periodic structures in engineering

A periodic structure is a mechanical system which can be represented by copying a part of the structure, in one, two or three dimensions. At the end of 19 th century, the concept of unit cell has started to be considered to study the mechanical behaviour of periodic systems, starting with the study of continuous structures [START_REF] Rayleigh | Xvii. on the maintenance of vibrations by forces of double frequency, and on the propagation of waves through a medium endowed with a periodic structure[END_REF], in particular for the analysis of elastic and acoustic wave propagation through periodic structures.

In our environment, a lot of engineering structures can be considered as periodic, even if they were not initially designed to use the physical properties of periodic systems. Common example of those structures could be rail trains for transports, fuselage of an airplane for aerospace or bridge, skyscraper for civil engineering: all have in common this particularity to be represented as a certain number of repetition of one unit cell. Thus, the particularity of engineering system to be periodic can be taken as an advantage to implement other periodic system embedded on them, in order to improve the mechanical characteristic of those. A very common example is the usage of honeycombs, initially designed to improve the stiffness of structures while light weighting replacing the traditional solutions like plate structures. In addition to the mechanical properties, the wave propagation specific to periodic structures can also be exploited to improve the dynamic of the initial structure, like vibrations disturbing airplanes or during the rockets take-off.

Two main interests can be identified by studying a unit cell of a periodic structure instead of the full structure itself. The first interest is the development of some techniques applied to the unit cell in order to get dynamic information of the full structure having a different representation of the description of the dynamic behaviour that is usually done, analysing the wave propagation

Auxetic Structures

Auxetic structures are periodic structures with the particularity to have negative Poisson's Ratio, meaning that they will be expanded when they are stretched and contracted when they are compressed, contrary to classic structure having an opposite behaviour [START_REF]Xxiii. on the flexure of a flat elastic spring[END_REF]. Auxetics structures present many interesting characteristics in a structural point of view, leading to many mechanical engineering applications in structural mechanics and dynamics, such as composites laminate using honeycombs [START_REF] Herakovich | Composite laminates with negative through-the-thickness poisson's ratios[END_REF], increased shear modulus and indentation resistance [START_REF] Evans | Auxetic materials: the positive side of being negative[END_REF], dynamics of fan blade in aerospace field [START_REF] Lira | A gradient cellular core for aeroengine fan blades based on auxetic configurations[END_REF] and also use of special foams [START_REF] Alderson | The interpretation of the straindependent poisson's ratio in auxetic polyethylene[END_REF][START_REF] Scarpa | Dynamic properties of high structural integrity auxetic open cell foam[END_REF] or wave propagation directionality [START_REF] Malischewsky | Unusual behaviour of wave propagation in auxetic structures: P-waves on free surface and s-waves in chiral lattices with piezoelectrics[END_REF][START_REF] Ruzzene | Wave propagation in sandwich plates with periodic auxetic core[END_REF]. Another very interesting property coming from auxetics structures is their shape deflection when they are solicited in out-of-plane bending. Contrary to classic structures having a anticlastic behaviour when they are bent, auxetic structures adopt a synclastic behaviour since the curvature is defined by the Poisson's Ratio [START_REF] Ashwell | The anticlastic curvature of rectangular beams and plates[END_REF] making them able to get a dome shape once solicited (Fig. 2. [START_REF] Alderson | The interpretation of the straindependent poisson's ratio in auxetic polyethylene[END_REF]. This specificity open some possibilities to define the geometry of the structure, making able to do slightly curved sandwich panels which can be adapted to cylindrical shapes like airplane fuselage, or can be used for morphing structures [START_REF] Olympio | Zero poisson's ratio cellular honeycombs for flex skins undergoing one-dimensional morphing[END_REF].

Floquet-Bloch theorem and Brillouin Zone

The Floquet-Bloch theorem is a combination of contributions from two different scientific fields, starting from the Floquet theorem [START_REF] Floquet | Sur les équations différentielles linéaires à coefficients périodiques[END_REF] in mathematics field, initially developed to find solutions about ordinary differential equations implying a periodic square matrix. Bloch extended and generalized this theorem 35 years later [START_REF] Bloch | Über die quantenmechanik der elektronen in kristallgittern[END_REF] from one dimensional to three dimensional periodic structures and applied it to physic fields, especially to describe the energy eigenstates for and electron in crystallography. Floquet demonstrated that if an ODE is written under the form:

(2.1)

∂ y ∂x = M(x)y
with M a periodic matrix i.e. M(x + x 0 ) = M(x) ∀x, the solution of the ODE can be expressed such as:

(2.2) Y (x) = A(x)e Rx
with A(x) a periodic matrix and R a constant matrix. From Eq. 2.2, one can write:

(2.3) Y (x + x 0 ) = Y (x)e Rx 0 .
Eq. 2.3 is actually the relation used to express the Floquet-Bloch theorem for periodic structures in analytical studies where examples are given in [START_REF] Collet | Floquet-Bloch decomposition for the computation of dispersion of two-dimensional periodic, damped mechanical systems[END_REF][START_REF] Scarpa | Control of Wave Propagation in Sandwich Beams with Auxetic Core[END_REF]. Indeed, if we fix x = 0, Y (0) and Y (x 0 ) refers respectively to the displacement at the left and right boundary of the unit cell, while Rx 0 is the Floquet multiplicator. With Felix Bloch's contribution, the Floquet-Bloch theorem has been written in a more compact form in one dimensional case:

(2.4) q R = e jkd q L with q L and q R respectively the degrees of freedom at the left and right of the unit cell, k the constant propagation also known as wave number, d the length of the unit cell and j the imaginary number such as j 2 = -1. This theorem gives a relation between the boundary displacements of the unit cell and the wave number, and by extension another relation between the frequency of the system and the propagation constant. Eq. 2.4 can also be extended to 2D and 3D periodic structures. Some of the techniques used to get the relation between frequency and wave number will be explained in chapter 3. This relation is very useful to plot a diagram called dispersion curve, giving essential information mainly about the group velocity of the different type of waves travelling inside the finite periodic structure, just using properties of the unit cell.

For 2D periodic structures, the relation is extended from one parameter of wave number to 2 different orientations in space, usually called k x and k y in case we are working in the Cartesian coordinates, leading to a dispersion surface (k x , k y , f requenc y) instead of a single curve (k, f requenc y).

Examples of a dispersion curve for 1D periodic structure and dispersion surface are shown Fig.

Group velocity and Bragg stop bands

The usage of the dispersion relation can lead to various information about the properties of periodic structures in engineering, such as the group velocity of the different types of waves travelling through the structure by looking the mode shape associated to the corresponding dispersion curve. By definition [START_REF] Lightill | Group Velocity[END_REF], the group velocity c g is defined as:

(2.5)

c g = 1 2π d f dk
where f and k denote respectively the frequency and the wave number. c g represents the velocity of a group of waves travelling through a non-dissipative structure. In a practical way, 2.1. PERIODIC STRUCTURES if the periodic system is harmonically excited by a source at a certain frequency, the waves generated by the excitation will travel along the structure with a velocity c g and the wave peaks will be separated by a certain wave length λ = 2π k , arrive to the other side of the structure and be reflected. Sometimes, a phenomena can occur if the periodic structure is designed in such a way that for specific ranges of frequency, the reflected and transmitted waves become destructive leading to an absence of wave propagation, which can be mathematically translated as c g = 0.

The phenomena is well known on the physics domain under the name of Bragg scattering [START_REF] Bragg | The diffraction of short electromagnetic waves by a crystal[END_REF].

An example of a diatomic mass is shown Fig. 2.4, representing the two branches of the dispersion curve of such a periodic system. We can observe that the group velocity is effectively equal to zero at the boundaries of the band gap, leading to the apparition of evanescent waves. Before finding applications in engineering structures, the Bragg scattering has been investi-CHAPTER 2. LITERATURE REVIEW gated in the field of solid state physics [START_REF] Kittel | Introduction to solid state physics[END_REF], using the first type of metamaterial called photonic crystals for optical wave propagation [START_REF] Yablonovitch | Inhibited spontaneous emission in solid-state physics and electronics[END_REF], inspiring later phononic crystals [START_REF] Kushwaha | Acoustic band structure of periodic elastic composites[END_REF] for acoustic and elastic waves. The area of frequency obtained that way is called a stop band or band gap and has very interesting properties in the design of mechanical system for dynamics since it can be used to create stop band filters to attenuate the amplitude of vibrations for one or several frequency ranges. The stop band can also be found by using more classical way of vibration analysis as frequency response functions, building a periodic structure with a finite number of unit cells and performing an harmonic analysis to excite a specific type of wave. If a stop band is present, the amplitude of displacement is exponentially attenuated and this attenuation is amplified if the number of unit cell used to create the periodic structure is increasing. Fig.

Nonlinearities in periodic media

In the field of dynamics and more specifically in vibration, a system can be defined as nonlinear if the elastic restoring force of the vibrating structure is not proportional to the displacement it occurs, contrary to linear structures. This non proportionality makes a dependence of the eigenfrequencies of the system with the amplitude of the structure [START_REF] Mook | Nonlinear Oscillations[END_REF], possibly leading to instabilities when the frequency is increased or decreased, characterized by a "jump" in the 2.1. PERIODIC STRUCTURES amplitude response. This phenomena is represented by a backbone curve, giving the evolution of the eigenfrequency measured versus the variation of amplitude [START_REF] Londoño | Identification of backbone curves of nonlinear systems from resonance decay responses[END_REF][START_REF] Renson | Robust identification of backbone curves using control-based continuation[END_REF]. The nature of the nonlinearity of a system can be multiple, coming for instance from contact [START_REF] Manini | Friction and nonlinear dynamics[END_REF], material [START_REF] Wang | Effect of material nonlinearity on buckling and postbuckling of fiber composite laminated plates and cylindrical shells[END_REF] or inertia [START_REF] Damaren | Simulation of Flexible-Link Manipulators With Inertial and Geometric Nonlinearities[END_REF]. Another type of nonlinearity existing is related to geometrical nonlinearities, appearing when large amplitude of vibrations are considered and making the hypothesis of small strain not viable any more. For instance, the Green-Lagrange terms for a beam in the plane (x, z)

under geometrical nonlinear hypothesis are written under the form [START_REF] Hui | Geometrically nonlinear analysis of beam structures via hierarchical one-dimensional finite elements[END_REF]:

(2.6)

E xx = ∂u ∂x + 1 2 ( ∂ 2 u ∂x 2 + ∂ 2 w ∂x 2 ) E zz = ∂w ∂z + 1 2 ( ∂ 2 u ∂z 2 + ∂ 2 w ∂z 2 ) E xz = ∂u ∂z + ∂w ∂x + ∂u ∂x ∂u ∂z + ∂w ∂x ∂w ∂z
where u and w represent respectively the displacements in x and z direction. If the displacement is assumed to be small, the crossed terms of Eq. 2.6 are neglected leading to classical linear equations of motion. One of the most known important type of nonlinearity occurring in this kind of problem is the Duffing nonlinearity [START_REF] Kovacic | The Duffing Equation: Nonlinear Oscillators and their Behaviour[END_REF], represented by a cubic term for the displacement, leading to amplitude-dependent frequency shifting and therefore can be exploited to change the resonance of the system. For weak nonlinearities, related to small strain and deflections, several methods have been developed to solve these kind of problem using mainly asymptotic developments. Among the most classical methods used, one can cite the perturbation method [START_REF] Jiang | Nonlinear perturbation solution of duffing equation and the application in control system[END_REF][START_REF] Telban | Hybrid-galerkin perturbation method for forced oscillations of the duffing equation[END_REF], the harmonic balance method [START_REF] Hu | Solution of a duffing-harmonic oscillator by the method of harmonic balance[END_REF] or the method of multiple scales [START_REF] Hill | Comparing the direct normal form method with harmonic balance and the method of multiple scales[END_REF].

More recently, some studies have been developed to study nonlinearities in periodic structures to take advantage of the frequency shifting to improve the vibroacoustics results in term of stop band width and/or location. Many of the research done have been based on the extension of the classical nonlinear techniques cited before, starting with the Harmonic Balance Method (HBM) applied to electromagnetic and acoustic waves [START_REF] Nayfeh | The application of the method of multiple scales to wave propagation in periodic structures[END_REF], first order perturbation using an asymptotic development on the wave number for a periodic spring-mass system [START_REF] Chakraborty | Dynamics of a weakly non-linear periodic chain[END_REF] or the effect of nonlinear damping [START_REF] Marathe | Wave attenuation in nonlinear periodic structures using harmonic balance and multiple scales[END_REF]. Perturbation like techniques such as the Lindstet-Poincaré method [START_REF] Verhulst | The Poincaré-Lindstedt method[END_REF] have been adapted to periodic structures especially for low amplitude nonlinear vibrations of discrete systems and finite element method [START_REF] Narisetti | A perturbation approach for predicting wave propagation in one-dimensional nonlinear periodic structures[END_REF][START_REF] Manktelow | Finite-element based perturbation analysis of wave propagation in nonlinear periodic structures[END_REF].

Nonlinear metamaterials have been investigated to observe any potential impact on the dispersion relation using different concepts as large amplitude of oscillations with pendulums [START_REF] Georgiou | An invariant manifold approach for studying waves in a one-dimensional array of non-linear oscillators[END_REF] or displacement with more conventional resonating systems [START_REF] Hussein | Dispersion characteristics of a nonlinear elastic metamaterial[END_REF]. Other more conventional metamaterials like negative-stiffness structures using bi-stable elements [START_REF] Nadkarni | Dynamics of periodic mechanical structures containing bistable elastic elements: From elastic to solitary wave propagation[END_REF] or negative-mass represented by nonlinear resonators in a linear periodic chain [START_REF] Fang | Wave propagation in one-dimensional nonlinear acoustic metamaterials[END_REF][START_REF] Lazarov | Low-frequency band gaps in chains with attached non-linear oscillators[END_REF][START_REF] Wang | Lower band gaps of longitudinal wave in a one-dimensional periodic rod by exploiting geometrical nonlinearity[END_REF] have been evaluated using the HBM. The HBM has been used in particular to examine the relative displacement between the resonant mass and the principal mass of the periodic chain.

Since we will be interested in the effect of large displacement of resonators in the following

Resonators in periodic systems

Periodic resonating systems can be presented as another class of periodic structures used for wave propagation attenuation that can be described as systems containing two different parts, the hosting structure representing the body of the periodic structure, and the resonator, usually represented as a spring-mass system as a first approximation and whose role is to improve the vibroacoustic performances of the overall system acting as an energy absorber. The basic principle behind the efficiency of the resonator is called the negative mass-effect [START_REF] Huang | On the negative effective mass density in acoustic metamaterials[END_REF][START_REF] Yao | Experimental study on negative effective mass in a 1d mass-spring system[END_REF], consisting on condensing the equation of motion of the resonator with the hosting structure to create a new combined single structure, commonly called metamaterial, for which the overall total mass is not constant anymore but depends on the frequency. Considering the notations shown Fig. 2.6, the effective mass m e f f can be expressed as:

(2.7)

m e f f (ω) = m h + m r ω 2 r ω 2 r -ω 2
where ω r = k r m r is the eigenfrequency of the resonator. The complete study to understand how to get this equation is shown in subsection 4.2.1. As shows Eq. 2.7, the effective mass depends on the frequency and so can become negative when ω > ω r , creating a negative phase shift between the displacement of the hosting structure and the resonator resulting also to the creation of the band gap, called resonant band gap. Due to the displacement phase shift caused by the resonators, the boundaries of the resonant stop band differs from the Bragg stop band creating a discontinuity for the wave number in the dispersion curve. A representation of the difference between Bragg and resonant band gaps and respective mode shapes taking the example of resonators attached to a beam structure is explained in [START_REF] Sharma | Local resonance and bragg bandgaps in sandwich beams containing periodically inserted resonators[END_REF]. Although resonators are present in a unit cell and are therefore also repeated in a periodic way, it is actually not necessary to have a exact periodic repetition in space for the resonators to observe the phenomena of resonant band gap [START_REF] Sheng | Locally resonant sonic materials[END_REF]. However, the choice of using a resonator in each cell is very convenient since it allows to use tools for periodic media, making studies faster and giving access to the dispersion curve for the stop band identification.

Based on that theory, several studies have been performed on periodic structures with resonators systems. Claus et al. studied the evolution of the resonating stop band varying the parameters of a spring-mass resonator attached to a beam and a plate [START_REF] Claeys | On the potential of tuned resonators to obtain low-frequency vibrational stop bands in periodic panels[END_REF] and proposed later a design of a waveguide composed with a square cross-section duct in which are attached several resonators [START_REF] Claeys | Design and validation of metamaterials for multiple structural stop bands in waveguides[END_REF]. Some studies are also focused on beam sandwich structures composed by homogenised continuum media with local resonators [START_REF] Chen | Dynamic behaviour of sandwich structure containing spring-mass resonators[END_REF]. Qian [START_REF] Kundu | Bandgap properties in locally resonant phononic crystal double panel structures with periodically attached spring-mass resonators[END_REF] studied two-dimensional periodic sandwich plate continuum systems using resonators between the two skin panels. The influence of using multi-resonators in series [START_REF] Xiao | Longitudinal wave band gaps in metamaterial-based elastic rods containing multi-degree-of-freedom resonators[END_REF][START_REF] Zhou | Effects of relevant parameters on the bandgaps of acoustic metamaterials with multi-resonators[END_REF] or in parallel [START_REF] Hall | Modelling and experimental validation of complex locally resonant structures[END_REF][START_REF] Xiao | Sound transmission loss of metamaterial-based thin plates with multiple subwavelength arrays of attached resonators[END_REF] tuning resonators in different ways to have close eigenfrequencies and increase the width of the resonating band gap on that way, or the effect of a graded meta material [START_REF] Banerjee | Frequency graded 1d metamaterials: A study on the attenuation bands[END_REF][START_REF] Pai | Acoustic metamaterial structures based on multi-frequency vibration absorbers[END_REF] to generate pseudo-periodic systems. The usage of more complex hosting structures like auxetic and/or chiral unit cell with ring resonators has been investigated [START_REF] Bacigalupo | Optimal design of lowfrequency band gaps in anti-tetrachiral lattice meta-materials[END_REF][START_REF] Bacigalupo | Optimal design of auxetic hexachiral metamaterials with local resonators[END_REF] performing optimization studies to get the widest band gap at the lowest frequency. More recent studies are focused on the active control of the resonator to make appear or disappear a band gap [START_REF] Billon | Design and experimental validation of a temperature-driven adaptive phononic crystal slab[END_REF][START_REF] Candido De | Adaptive locally resonant metamaterials leveraging shape memory alloys[END_REF] creating resonators based on temperature-dependent materials, softening the whole structure as the temperature increases to switch the band gaps to lower frequencies. The effect of damping on resonators has been also investigated to see the effect on the dispersion diagram [START_REF] Belle | On the impact of damping on the dispersion curves of a locally resonant metamaterial: Modelling and experimental validation[END_REF] and on the sound transmission loss (STL) with experimental validation [START_REF] Belle | The impact of damping on the sound transmission loss of locally resonant metamaterial plates[END_REF] and an industrial application of the effect of local resonant material is proposed for noise reduction on an automotive [START_REF] Jung | Realisation of a locally resonant metamaterial on the automobile panel structure to reduce noise radiation[END_REF].

C H A P T E R

LINEAR RESONATORS IN AUXETIC STRUCTURES

T his chapter will focus on the linear structure analysis of periodic structures including resonators for stop band investigation. The first part is dedicated to the principal theories and tools recently found in the literature and adapted specifically to the studies of our type of structures. The second part is dedicated to the study of some resonating systems example, starting with beam element models as a first approach to get fast results, followed by shell elements analysis to get results closer to the practical application and a final part will conclude with volume elements including skins to get a model close to something that could be manufactured and tested.

Wave finite element methods for linear periodic analyses

Wave finite elements methods (WFEM) have been developed for periodic structures using implicitly (Transfer Matrix Method [START_REF] Lin | Dynamics of Beam-Type Periodic Structures[END_REF][START_REF] Zhong | On the direct solution of wave propagation for repetitive structures[END_REF]) or explicitly (direct and inverse WFEM) the Floquet-Bloch theorem to study free wave propagation with eigenvalue formulations in one and two dimensions, or forced vibrations for one dimensional cases. The theory behind WFEM has been initially

proposed by Mead [START_REF] Mead | A general theory of harmonic wave propagation in linear periodic systems with multiple coupling[END_REF] while the first numerical investigations have been done by Orris [START_REF] Orris | A finite element study of harmonic wave propagation in periodic structures[END_REF] with h-version of finite element method and later by Mead with a p-version [START_REF] Mead | Free vibration of a thin cylindrical shell with discrete axial stiffeners[END_REF] before being studied more deeply in the 2000's with the apparition of more powerful computers allowing also system response at mid frequency range [START_REF] Mace | Finite element prediction of wave motion in structural waveguides[END_REF][START_REF] Mencik | Multi-mode propagation and diffusion in structures through finite elements[END_REF].

Inverse WFEM

The inverse WFEM (IWFEM) has been the most used tool to analyse the dispersion curves of a periodic structure in this research. This method gives the advantage to be more computationally stable than the direct WFEM as it will be explained in subsection 3.1.2, but does not give as much information as the latter one since it is based on the imposition of the wave number to get the frequency. That means that only the real part of the wave number is available for the dispersion curve plotting, while the direct method will give also the imaginary part. However, the imaginary part of the wave number might not be necessary here since the majority of simulations are performed assuming 2D wave propagation.

Considering the 2D unit cell presented in Fig. 3.1, one can write the dynamic relation of the system:

(3.1) (K u -ω 2 M u )q = f
K u , M u being respectively the stiffness and mass matrix of the unit cell, q the degrees of freedom and f the forces coming from the neighbour unit cells constituting the periodic structure, written under the form:

(3.2) q = q i q B q T q L q R q LB q RB q LT q RT T ,

f = f i f B f T f L f R f LB f RB f LT f RT T .
Assuming there is no internal forces in the system, one can write that f i = 0. The Floquet-Bloch conditions linking the degrees of freedom at the boundaries of the unit cell are given such as:

(3.3)

q RB = λ x q LB , q LT = λ y q LB , q RT = λ x λ y q LB , q R = λ x q L , q T = λ y q B ,

WAVE FINITE ELEMENT METHODS FOR LINEAR PERIODIC ANALYSES

with λ x = e jµ x and λ y = e jµ y where µ x and µ y are the reduced wave number in the two directions of propagation. The equilibrium of forces at interfaces gives:

(3.4)

f LB + λ -1 x f RB + λ -1 y f LT + λ -1 x λ -1 y f RT = 0, f L + λ -1 x f R = 0, f B + λ -1 y f T = 0.
Eq. 3.3 can be written under a matrix form verifying q = Λ a q (r) such as:

(3.5) Λ a =                     I 0 0 0 0 I 0 0 0 Iλ y 0 0 0 0 I 0 0 0 Iλ x 0 0 0 0 I 0 0 0 Iλ x 0 0 0 Iλ y 0 0 0 Iλ x λ y                     , q (r) =        q I q B q L q LB        .
Eq. 3.4 can also be written under matrix form verifying the relation Λ b f = 0:

(3.6) Λ b =        I 0 0 0 0 0 0 0 0 0 I Iλ -1 y 0 0 0 0 0 0 0 0 0 I Iλ -1 x 0 0 0 0 0 0 0 0 0 I Iλ -1 x Iλ -1 y Iλ -1 x λ -1 y        . Pre multiplying Eq. 3.1 by Λ b and writing K r = Λ b K u Λ a (M r = Λ b M u Λ a for the mass matrix)
we get the following equation:

(3.7) K r (λ x , λ y ) -ω 2 M r (λ x , λ y ) q (r) = 0
Eq. 3.7 represents an eigenvalue problem in ω 2 , establishing the relation between the wave numbers (µ x , µ y ) in input and the frequency in output. Imposing values of µ x and µ y in the IFBZ gives access to the contour of the dispersion surface.

Direct WFEM

The direct form can be obtained by inverting Eq. 3.7 in order to get an eigenvalue problem in λ x or λ x as a frequency in output. Let consider matrices

K λ x , K λ y , K 1 λx , K 1 λy , K λ x λ y , K λx λy , K λy λx , K 1 λxλy ,
K 1 be partitions of the K r matrix, such as:

(3.8) K r = λ x K λ x + λ y K λ y + 1 λ x K 1 λx + 1 λ y K 1 λy + λ x λ y K λx λy + λ y λ x K λy λx + λ x λ y K λ x λ y + 1 λ x λ y K 1 λxλy + K 1 .
Processing exactly the same way partitioning the M r , Eq. 3.7 can be rewritten in the following form:

(3.9)

λ x D λ x +λ y D λ y + 1 λ x D 1 λx + 1 λ y D 1 λy + λ x λ y D λx λy + λ y λ x D λy λx +λ x λ y D λ x λ y + 1 λ x λ y D 1 λxλy +D 1 q (r) = 0 with D X = K X -ω 2 M X , (X being λ x ,λ y ,...)
. By imposing an arbitrary value for λ y = c (same reasoning for λ x ) and multiplying Eq. 3.9 by λ x , it comes:

(3.10) αλ 2 x + βλ x + γ q (r) = 0 with α = D λ x + cD λ x λ y + 1 c D λ y , β = D 1 + cD λ y + 1 c D 1 λy and γ = D 1 λx + cD λx λy + 1 c D 1 λxλy . Eq. 3.
10 is a quadratic eigenvalue problem in λ x , meaning that the dispersion relation is found by imposing a value of frequency and get the value of λ x as an output.

Transfer Matrix Method

Inverse and direct WFEM are both designed for infinite structure methods, meaning that the analysis of the unit cell applying Floquet-Bloch will give a dispersion relation as a result. The Transfer Matrix Method (TMM) can be useful to both finite and infinite structure analysis but is only usable for 1D Periodic system and presents some ill-conditioning problems. However, the TMM has the advantage to do fast computing simulation to get the FRF of a structure, compared to classical direct methods used in general vibration theory.

TMM theory

Let consider a unit cell belonging to a 1D periodic structure, where q L , q R and q I denote respectively the degrees of freedom on the left, right and internal degrees of freedom. Considering the dynamic matrix of the unit cell as D u = K u -ω 2 M u and assuming that no internal forces are present, we can decompose D u such that :

(3.11)     D LL D LI D LR D IL D II D IR D RL D RI D RR         q L q I q R     =     f L 0 f R    
By isolating q I from the second line of the matrix in Eq. 3.11, one can reduce the matrix doing a condensation of internal degrees of freedom and create a new dynamic matrix such as:

q L q R = f L f R
where :

(3.13)

DLL = D LL -D LI × D -1 II × D IL , DLR = D LR -D LI × D -1 II × D IR , D RL = D RL -D RI × D -1 II × D IL , DRR = D RR -D RI × D -1 II × D IR .
We can then construct a Transfer Matrix by rewriting the matrix 3.12 in order to have a relation between displacement and force on the right side with displacement and force on the left side:

(3.14) q R -f R (n) = T 11 T 12 T 21 T 22 q L f L (n)
.

Where n represents the n th unit cell of the structure and the Transfer Matrix coefficients are given by:

(3.15)

T 11 = -D-1 LR × DLL , T 12 = D-1 LR , T 21 = -DRL + DRR × D-1 LR × DLL , T 22 = -DRR × D-1 LR .
. The Transfer Matrix is symplectic, i.e. it can be written under the for T T JT = J with

(3.16) J = 0 I -I 0 ,
and I being the identity matrix. Thereby, the eigenvalues of the Transfer Matrix are linked with the wave numbers of the system [START_REF] Zhong | On the direct solution of wave propagation for repetitive structures[END_REF], implying that all the eigenvalues are going in pair under the form λ and 1 λ where λ = e jµ . The continuity relation and the force equilibrium between each cell of the periodic structure can be respectively written q (n+1)

L = q (n) R and f (n+1) L = -f (n) R , or in vector form: (3.17) q L f L (n+1) = q R -f R (n)
If we consider a finite structure with N unit cells, combining Eq. 3.14 and 3.17, we get the following relation between the extremities left and right of the finite structure:

(3.18) q L f L (N) = T N-1 q R -f R (1)
.

Inverting the T N-1 matrix and applying boundary conditions, one can obtain the value of displacement for each point of the structure and get the FRF. This method is theoretically faster than usual direct methods to calculate FRFs because it involves to invert only p × p systems (where p is the number of degrees of freedom of the unit cell) while direct methods require the inversion of N p×N p systems. However, it appears that the matrix T N-1 is full and ill-conditioned as the number N increases. This comes from the fact that the eigenvalues of T are written under the form λ = e jµ and 1 λ = e -jµ . Due to the symplectic properties of the Transfer Matrix, one can diagonalize the matrix T and write it under the form:

(3.19) T = ΦLΦ -1
where L is the diagonal matrix of the eigenvalues of T and Φ the matrix of eigenvectors.

Multiplying T by itself N -1 times in relation 3.19, we get T N = ΦL N Φ -1 where:

(3.20) L N =                e N jµ 1 0 • • • • • • • • • 0 0 . . . . . . 0 . . . . . . e N jµ k 2 . . . . . . . . . . . . e -N jµ 1 . . . . . . . . . . . . . . . 0 0 • • • • • • 0 e -N jµ k 2               
and k is the size of the L matrix. Inverting T N is then the same operation of inverting L N , which can lead to numerical problems when N is high due to the exponential terms. Numerical ill-conditioning can also be caused by the truncation of inertia terms in the dynamic stiffness matrix [START_REF] Waki | Numerical issues concerning the wave and finite element method for free and forced vibrations of waveguides[END_REF].

A TMM ill-conditioning solution

A very similar issue has been found in acoustics for multi-layering purpose, and a solution has been proposed by Dazel et al. [START_REF] Dazel | A stable method to model the acoustic response of multilayered structures[END_REF]. However, the problem can be adapted to elastic wave problems, and here will be explained how to get the FRF of the full periodic structure using this method.

20 3.1. WAVE FINITE ELEMENT METHODS FOR LINEAR PERIODIC ANALYSES Let Q R = q R -f R -respectively Q L = q L
f L -being state vectors decomposed in a translation vector noted Ω and an information vector noted X such as Q L = Ω L X L for each interface between unit cells of the periodic system. We can write the iterative equation:

(3.21) Ω n+1 L X n+1 L = TΩ n L X n L
equivalent to the Transfer Matrix relation found in Eq. 3.14. The first step consists in initializing the operators by applying boundary conditions. Let consider a free-free boundary condition for the full periodic structure. At the left extremity of the periodic system, we get:

(3.22) q 1 0 Q 1 L = I 0 Ω 1 L q 1 X 1 L
where q 1 represents the degrees of freedom vector at the left extremity of the periodic structure. The idea to fix the ill-conditioning issue is to separate exponentially growing eigenvalue terms and their inverse in order to hide the non-stable terms. Eigenvalues of T noted λ i = e jµ i are computed and using a L,N formulation [START_REF] Nobrega | Vibration band gaps for elastic metamaterial rods using wave finite element method[END_REF] and are sorted such as ℜ(µ 1 ) > ℜ(µ 2 )... > ℜ(µ 2m ), 2m being the dimension of T and also considering ℜ(µ p ) > 0 for p ∈ [1, m] and ℜ(µ p ) < 0 for

p ∈ [m + 1, 2m].
Then the eigenvector matrix of T noted Φ is built considering the same order than the sorted eigenvalue previously done. Let Ψ = Φ -1 , Φ k the k th column of Φ and Ψ k the k th row of Ψ. We can decompose T such as:

(3.23) T = 2m i=1 e jµ i Φ i Ψ i .
The eigenvector matrix is then decomposed in the following way:

(3.24) Φ = [Φ q |Φ r ], Ψ = [Ψ q |Ψ r ]
where Φ q and Ψ q represent the matrices of the m first columns (respectively rows) of Φ and Ψ while Φ r and Ψ r represent the m last ones.

T can now be written:

(3.25) T = m i=1 e jµ i Φ i Ψ i + e jµ m α
where :

(3.26) α = Φ r        1 e j(µ m+1 -µ m ) . . . e j(µ 2m -µ m )        Ψ r .
We can notice doing this separation that the argument of each term in the α matrix is negative, implying that they will not generate exponentially growing term. Combining Eq. 3.21 and 3.25 leads to:

(3.27) Q n+1 L = Ω n+1 L X n+1 L = TΩ n L X n L = m-1 i=1 e jµ i Φ i Ψ i Ω n L + e jµ m αΩ n L X n L .
The next step is to create operators X n+1

L

and Ω n+1 L in order to respect Eq. 3.27. Exponentially growing terms will be included in X n+1 L since it will not be taken in consideration for the numerical computation. We consider the following relation:

(3.28) X n+1 L =        e jµ 1 .
. .

e jµ m-1 e jµ m               Ψ 1 [Ω] n L . . . Ψ m-1 [Ω] n L Ψ m [Ω] n L        [Ξ] n X n L .
[Ξ] n can be interpreted as the projection of the translation matrix on the m first eigenvectors of the Transfert Matrix, and got the advantage to be numerically stable contrary to the latter one.

Also, the operator Ω n+1 L allowing the relation in Eq. 3.27 has to be written in the following way:

(3.29) Ω n+1 L = Φ r + αΩ n L Ξ -1 n        e j(µ m -µ 1 )
. . . (3.31)

e j(µ m -µ m-1 ) e j(µ m -µ m )        . The relation Q n+1 L = Ω n+1 L X n+1
X N+1 L = (Ω N+1 L b ) -1 F. Now that X N+1
L is known, the information about the displacement value at the extremity right of the structure is also known. The final step is to use the iterative relation between X N+1 L and X N L using Eq. 3.28 and inverting it:

(3.32) X n L = Ξ -1 n        e -jµ 1 . . . e -jµ m-1 e -jµ m        D m X n+1 L .
The matrix D m is again stable because it only contains negative argument terms, which does not give any computation problem. By applying this relation until n = 1, we get the final equation giving the value of displacement at the left extremity of the structure:

(3.33)

X 1 L = N+1 i=1 Ξ -1 i D m X N+1 L .

Beam elements

Unit cell models

Prior to analyse periodic structures including resonators, a choice has to be done for the shape of the hosting unit cell. So far, two main patterns have been used in literature to represent sandwich cores: the first one represents classical hexagonal honeycombs, and the second one the re-entrant models, having the particularity to be auxetic. For a 2D periodic structure, there is infinite manners to chose the unit cell in order to perform Floquet analysis. Usually, the choice is done taking in account how to write the periodic conditions as easily as possible, i.e. less nodes at the interface for less conditions giving more convenience and facility to apply the periodicity conditions. The way the unit cell has been chosen is represented in Fig. 3.2. Also, Fig. under Matlab to do the calculations and post-processing, mainly concerning the dispersion curves and mode shapes of the periodic system. The element type used here is BEAM4, corresponding to beams under Euler-Bernoulli hypothesis [START_REF] Timoshenko | History of strength of materials[END_REF]. The beams are all connected on a plane at a fixed height (z = 0) meaning that the equations of motion for in-plane (u x , u y , θ z ) and out-of-plane (u z , θ x , θ y ) are uncoupled [START_REF] Tian | Numerical simulation of elastic wave propagation in honeycomb core sandwich plates[END_REF], making able the separation of studies for in-plane and out-of-plane wave propagation blocking accordingly the degrees of freedom to have an easiest reading of dispersion curves.

Simulation and Results

Hosting unit cell selection

After applying Floquet boundary conditions accordingly to the two unit cells, the dispersion relation is computed using Inverse WFEM technique and the dispersion diagrams are plotted.

Segment [O -A] [A -B] [B -C] [C -O] µ x [0, π] π [π, 0] 0 µ y 0 [0, π] π [π, 0]

Results including resonators

Results observed in previous subsection shows that the re-entrant configuration exhibits some Bragg band gaps for out-of-plane movement. However, one can notice that these band gaps are not located in low frequency areas since they are not related to the first modes of the periodic structure. For that reason, resonators will be included in the hosting unit cell in order to improve the vibroacoustic behaviour of the structure by trying to open band gaps in lower frequencies.

For this model, the resonator is represented by a simple beam attached to the corners of the unit cell. The second study focuses on the efficiency related to the number of resonators present in the unit cell, in order to see if it is necessary to chose between 1,2 or 4 resonators to get the This observation can be explained since the first case represents a configuration where

ρ 1 = ρ 2 ,
meaning 2 resonators with same density, while the second case represents a configuration where ρ 2 0 for the 4 resonators configuration, which is equivalent to say that the 2 extra resonators do not play any role in the dynamic behaviour of the system and overall behaves like a system with 2 resonators. With the same reasoning, the configuration n • 1 of 3.11a describes a case where Based on the previous results, a parametric analysis can be done based on the density variation of the resonators for the unit cell including 2 resonators in order to analyse the evolution where m res and m aux denote respectively the total mass of the resonators and the mass of the hosting structure. As we could expect, increasing the mass of the resonators makes decrease the band gaps to lower frequencies. However, after around 30% of added mass, the effect seem to be less impacting. As a consequence, a compromise has to be found around this value trying to not outreach too much added mass while having interesting vibroacoustic properties. Frequency (Hz) FIGURE 3.12. Parametric study on the mass variation of the resonator

Conclusion

Beam element analyses have been performed here for infinite and finite periodic structure hypotheses. Although beams do not represent the reality of sandwich structures, the aim was to CHAPTER 3. LINEAR RESONATORS IN AUXETIC STRUCTURES introduce a simple model showing example of studies that can be done to analyse the impact of the resonators on the dispersion curves. It has been shown that re-entrant structures present band gaps for out-of-plane type of waves contrary to hexagonal pattern. Also, the optimal number of resonators in each unit cell for band gap opening for a constant mass is two, corresponding to one resonator in each butterfly pattern of the periodic structure.

Shell elements

Unit cell models

Models involving shell elements give a better representation of the reality, since we can have access to real thickness/height for the hosting structure, and give the possibility to add skins to the core. The element type used to model the unit cell is under ANSYS is called SHELL181,based on the Midlin-Reissner shell theory [START_REF] Mindlin | Influence of Rotary Inertia and Shear on Flexural Motions of Isotropic, Elastic Plates[END_REF][START_REF] Reissner | he Effect of Transverse Shear Deformation on the Bending of Elastic Plates[END_REF] taking in account through-the-thickness shear deformations. The geometry of the unit cell designed under ANSYS APDL is shown Fig. 3.13 and the geometrical and material properties associated in table 3.4. For such a model, in-plane and out-of-plane equations are not uncoupled as they were for the beam element models, meaning that the dispersion curves will include all the type of waves (compression, bending, shear).

Results without resonators

As done for the beam element analysis, the analysis of the core of the re-entrant unit cell without including the resonator is done to get the dispersion curve using Floquet conditions on the same way as previously. Each eigenfrequency we compute using the Inverse WFEM for a fixed couple (µ x , µ y ) belongs to a certain type of wave associated to a wave mode φ. However, when several type of waves are present, it is not possible to separate the out-of-plane and in-plane studies in the pre-processing, making complicated the following of each branch of the dispersion curve, especially when they are crossing at a particular points (ω 0 , (µ x 0 , µ y 0 )). The method to sort and follow the different branches is called branch tracking, and a method based on the computation of the group velocity for one dimensional case has already been developed [START_REF] Billon | Design and experimental validation of a temperature-driven adaptive phononic crystal slab[END_REF]. Since the equations of motion in-plane and out-of-plane are coupled, a method based on the Modal Assurance Criterion (MAC), initially introduced for the collinearity / orthogonality of modes in classical vibration study [START_REF] Allemang | The modal assurance criterion (mac): Twenty years of use and abuse[END_REF][START_REF] Collet | Floquet-Bloch decomposition for the computation of dispersion of two-dimensional periodic, damped mechanical systems[END_REF], is used in order to track the branches of the dispersion curve making easier the reading of the graph. The principle of the MAC is to compare two consecutive set of wave modes φ p (ω i ) and φ p+1 (ω j ) for the p th couple (µ x p , µ y p ) to see if they are close to be collinear. The equation to have the value of the criterion is:

(3.35) M AC = |φ H p (ω i )φ p+1 (ω j )| ||φ p (ω i )||.||φ p+1 (ω j )|| .
If we consider the computation of the N first eigenfrequencies, the evaluation between φ p (ω i ) for a fixed i and φ p+1 (ω j ) for j ∈ [1, N] is done. Then, M AC represent a vector of size (1, N) whose all the values are included between 0 and 1 thanks to the normalization. Taking the maximum of M AC gives the index where φ p (ω i ) and φ p+1 (ω j ) are collinear which allows to conclude that the eigenfrequencies ω i (µ x p , µ y p ) and ω j (µ x p+1 , µ y p+1 ) belong to the same branch. The dispersion curve of the re-entrant shell structure before and after the MAC algorithm is given Fig.

Results with resonators

Since shell elements are used, the resonator will be represented as a single long beam connecting the sides of the hosting structure for simplification. The shape of the unit cell is represented in Fig. Another possible solution is to modify the stiffness of the resonator by keeping the mass constant. If the thickness (h r ) is decreased, the out-of-plane thickness of the resonator will also decrease and so the branches of the dispersion curve will switch to lower frequency to get rid of the obstructing branch. Since the resonators are also part of the structure contrary to theoretical or purely homogeneous models, this modification will also have consequences on the branches of the hosting structure but with a limited impact since the modifications are done locally. Results both for the dispersion curves and the FRF using the TMM are shown in Fig. 3.18 for ρ r = 14ρ and h r = h 1. [START_REF] Ashwell | The anticlastic curvature of rectangular beams and plates[END_REF] . One can see that the branches related to the band gap decreased in frequency to not be disturbed by the structural wave mode. As a consequence, the width of the band gap has Hz, but as also decreased in term of frequency passing from 2750 Hz to 1580 Hz for the middle frequency of the band gap. To understand quantitatively what result is the most interesting, a performance index f i for the stop band i can be calculated doing the ratio between the width of the band gap and its position [START_REF] Scarpa | Directional and band-gap behavior of periodic auxetic lattices[END_REF]:

(3.36) f i = ∆ω i ω m,i
where ∆ω i and ω m,i represent respectively the width and the middle location of the band gap.

By comparing the two results, it appears that f 1 = 0.148 for the first configuration and f 2 = 0.165 meaning that the results Fig. 3.18 provide better results considering this criterion.

Results including skins and resonators

According to the previous study, it has been shown that using resonators in a re-entrant core improves the wave propagation attenuation by creating resonant stop bands for out-of-plane waves. Here skins will be added to the model to get something closer to the reality and observe how they are influencing the wave propagation in such type of cores. In theory, adding skins to the core structure should create some modifications in the dynamic behaviour, since the global stiffness of the structure is drastically increased and therefore a shift of dispersion curves to higher frequencies should appear.

VOLUME ELEMENTS

The skins are represented by plates with same material properties as the core structure with a thickness T s = 0.6 mm, and are attached to the top and bottom of the core. Additional Floquet boundary conditions are taken in consideration, adding also the corners of the unit cell and the dispersion analysis is done with a resonator height h r = 0.3 mm and resonator density such as ρ r = 10ρ. The dispersion contour and the mode shapes are shown Fig. 3.19 and exhibits that, contrary to the model containing only the core of the periodic structure, the separation between in-plane wave modes in red and out-of-plane modes in blue is not always done using the same algorithm of displacement sorting used before. This phenomena can be explained by observing the mode shape associated to the green dot in the figure, showing a combined wave behaviour associated x,y and z directions. As a consequence, no clear band gaps related to a certain direction is created hypothetically meaning that the design of the used resonators does not have enough impact to create band gaps for the unit cell containing the core and the skins. This result shows the limitation of designing resonator using only shell element, implying that the model has to be improved in order to get relevant vibroacoustic performances, using for instance volume elements to have more flexibility on the shape of the resonators.

Conclusion

Compared to beam elements, shell elements makes able the representation of more realistic structures using 3D shapes and also including the skins to see the impact on the wave propagation.

Results without skins are promising since a clear out-of-plane stop band is appearing using flat resonators even though their added mass stay relatively high. Results with skins does not show interesting results for stop band analyses, principally because coupling waves are travelling through the structure and the resonators, making the wave propagation analysis much more complicated due to the unusual shape of the dispersion contour and showing the limitations of the models for the resonator design.

Volume elements

Preliminary studies

Since shell elements are not sufficient to get relevant results including resonators, the model of the unit cell including skins is done using volume element. Linear 8-nodes elements solid185 containing only 3 degrees of freedom (u x ,u y and u z ) are used. The geometry of the hosting structure is the same as used for shell element in previous section since the only parameters changed are the resonators, which are made with a combination of a simple piece of plate in which is attached a parallelepipedic mass at its end. The geometry of the resonators and their location in the unit cell is shown Fig. 3.20 for 2 different configurations. The blue configuration is equivalent to the one used for shell elements i.e. using a resonator in each butterfly pattern, instead the red configuration correspond to 4 resonators put in one over two butterfly pattern. The 3.5 as well as the material properties.

For the following studies, the resonators are made of steel and their material properties remain unchanged. Only geometrical parameters are susceptible to change since it is the most convenient way to change their mass and stiffness in a practical point of view since using volume elements makes more flexibility on the shape model of the resonators, contrary to beams or shell elements with which the topological possibilities were much more limited. In order to simplify the finite element model design, four resonators are included in each unit cell. Indeed, the final shape of the unit cell has been created just doing one quarter of the unit cell including one of the resonator, Keeping the same performance considerations as previously, the main constraint we are looking to respect in this study is to keep the ratio between the width of the band gap over its location as low as possible in frequency keeping a constant mass. A first study using the red resonator configuration is performed using standard parameters to get the overall shape of the dispersion curve and the mode shapes of the resonator at the boundary of the band gap.

E [GPa] 210 L s [mm] 2.5 h m [mm] L m + h s ρ [k g.m -3 ] 7800 L m [mm] 1 b s [mm] 1 ν 0.33 h s [mm] 0.2
The results are presented in Fig. 3.21 and one can observe that contrary to the model used in subsection 3.3.4, this time the resonator mode shapes does not present coupled directionalities but are clearly defined as out-of-plane bending. This difference leads to a more clear and familiar resonating dispersion curve a distinguishable band gap between 7600 and 8000 Hz. Moreover, one can observe that the boundaries of the resonant stop band looks to be defined at the point O (µ x = µ y = 0) since the resonating waves are not obstructed by a structural wave as it was the case in the results shown in the shell element models Fig. 3.17, most likely due to the usage of the skins and meaning that the studies can be largely reduced to the interval [O -A].

A second analysis is done comparing the blue configuration and the red configuration for the same total mass of resonators, keeping the same mass m R we used for the previous study with the red configuration. One have to determine the mass of the blue resonator m B such as m B = 2m R since the red configuration contains 2 times more mass than the blue one. Generally speaking, the mass of the resonator is defined as:

(3.37)

m res = ρ(V s + V m ) = ρb s (h s L s + h s L m + L 2 m ),
with V s the volume of the resonator's beam and V m the volume of the resonator's mass. If we decide to change the value of L m to get a new value of mass, we obtain after calculation L b m = 1.6 mm for the blue configuration to get the same total mass as the red configuration considering L r m = 1 mm, L b m and L r m being the parameter L m for the red and blue configuration respectively. Under these consideration, we get a percentage of added mass due to the resonator m res m aux = 13.7% (see Eq. 3.34 for notations) for both configurations. The [O -A] part of the dispersion curve is plot on both cases. The results are shown Fig. 3.22 and we can observe that the band gap is located at lower range of frequency for the blue resonator case compared to the red case for a very close bandwidth. In a more precise way, with the criteria Eq. 3.36 used in previous section, we got f r = 0.045 for the red configuration and f b = 0.057 for the blue configuration meaning that the blue configuration provides better results.

Parametric study of the resonators

A parametric study is done on the geometric variables defining the resonators in order to get larger band gaps located in lower frequencies. Theoretically speaking, the shift of frequency could be achieved by increasing the concentrated mass of the resonator while decreasing its stiffness.

In other words, the ideal solution would be to maximize the value of L m (and so h m ) to get the biggest mass possible while minimizing the thickness h s to get the lowest stiffness, keeping an overall constant mass for the resonator. However, compared to the shell model, the adding of skin and the more complex shape of the resonator could be less predictable as before, encouraging the use of a parametric study. Some geometrical considerations also have to be respected in order to create the model, mainly considering the fact that the total length L s + L m has to be less than the height H s to stay inside the unit cell space design as we can see in Fig. 3.23.

As a first approach, the study can be done only considering the parameters h s and b s as design variables. Under that consideration, one does not have to take in account the limit length mentioned before, and the only constraint is given by the mass of the resonator. Taking as a starting point the standard parameters considered for the blue configuration study in previous subsection, the function to take in consideration to have a constant mass will be under the form b s = f (h s ) where the expression is obtained using Eq. 3.37: and L m to keep the mass constant: However, the width of the stop band could not be predicted that way which make the usage of several parametric analysis essential to get optimised results.

(3.39) L s = 1 h s (h s L s 1 + h s L m 1 + L 2 m 1 -h s L m -L 2 m ).

Conclusion

With a large variability in the model design, volume element give a better representation of the reality compared with shell of beam elements, giving also more possibility for the design of the resonators. Results obtained with skins are this time more promising than shell elements, and a parametric analysis on the geometrical properties of the resonators has been performed to show which is the best configuration in term of width and low frequency stop band under some geometrical constraints. This parametric has been useful since the resonator's shape obtained for the best results was not intuitively predictable as it could have been using flat resonators (section 4.1).

C H A P T E R

NONLINEAR RESONATORS IN PERIODIC STRUCTURES

This chapter is dedicated to the analysis of periodic structures in nonlinear regime, assuming geometric nonlinearities due to the large displacements that generate resonators close to their eigenfrequencies. The first part focuses on the methods found on the literature, such as the Lindstedt-Poincaré method which has been modified and adapted to perform analyses with a commercial software, and a method developed during the thesis called the inverse method in order to identify the apparent band gap of a finite structure under nonlinear hypotheses. The second part will treat some academic example to validate the tools presented, and a case of a nonlinear periodic system is studied to see the impact on the dispersion curve in the case of a resonating periodic structure.

Theories and tools for nonlinear periodic structures

Lindstedt-Poincaré Method adapted

The Lindstedt-Poincaré Method has originally been introduced to solve nonlinear differential equations in different scientific fields, such as quantum physics, electromagnetic or mechanics.

The principle of the method is to do an asymptotic development to the degrees of freedom and the frequency of the system in order to have a corrected frequency term, generally depending on the amplitude of displacement. Narisetti has extended this nonlinear method to periodic discrete structures [START_REF] Narisetti | A perturbation approach for predicting wave propagation in one-dimensional nonlinear periodic structures[END_REF] while Manktelow [START_REF] Manktelow | Finite-element based perturbation analysis of wave propagation in nonlinear periodic structures[END_REF] adapted it for finite element method applications.

The same idea will be explained here, using a slightly different approach for the periodicity hypothesis in the literature ; previously the unit cell was selected with its neighbours to apply the Floquet-Bloch conditions and the Lindstedt-Poincaré method, instead here the inverse WFEM presented section 3.1.1 is used, but both ways are leading to the same final result. Rewriting the problem that way was more a question of convenience to extend programs developed for linear structure for nonlinear analysis.

Let consider the dynamic equation of the unit cell including nonlinear hypothesis:

(4.1) M u ẍ + K u x + ε f nl = f l
where x denotes the displacement vector, ẍ the acceleration, f l the linear forces applied to the unit cell, f nl the nonlinear forces assumed as proportional to cubic displacements, which is the most encountered in geometrical nonlinearities, and ε a small parameter that tune the level of nonlinearity. Considering the dimensionless time τ = ωt, the Lindstedt-Poincaré asymptotic expansions at first order are written as:

ω = ω 0 + εω 1 + o(ε), x = x (0) + εx (1) + o(ε). (4.2)
Including Eq. 4.2 in 4.1 leads to the following expression :

(4.3) (ω 0 + εω 1 ) 2 M u d 2 (x (0) + εx (1) ) dτ 2 + K u (x (0) + εx (1) ) = f l -ε f nl .
Developing Eq. 4.3 by neglecting all the terms associated to order 2 or higher leads to two equations:

(4.4) o(ε 0 ) : ω 2 0 M u d 2 x (0) dτ 2 + K u x (0) = f l , (4.5) o(ε 1 ) : ω 2 0 M u d 2 x (1) dτ 2 + K u x (1) = -2ω 0 ω 1 M u d 2 x (0) dτ 2 -f nl .
Eq. 4.4 represents the linear part of the dynamic equation for the unit cell. Assuming the harmonic regime, forces and displacements can be written under the form x (0) = X 0 e jτ and f l = F l e jτ , X 0 and F l being the amplitude of displacement and force respectively. Eq. 4.4 becomes:

(4.6) (K u -ω 2 0 M u )X 0 = F l .
Eq. 4.6 is equivalent to 3.1, so we can apply the Inverse WFEM seen in subsection 3.1.1 to transform the equation:

(4.7) (K r -ω 2 0 M r )φ i = 0
where φ i is the i th eigenvector of the system. The left-hand-side of Eq. 4.5 represents the linear kernel of the equation and is similar to the one in Eq. 4.4, meaning that the wave modes associated are similar. If we note y (0) and y (1) the reduced displacement vectors such that x (0) = Λ a y (0) and x (1) = Λ a y (1) , we premultiply Eq. 4.5 by the matrix Λ b already seen in subsection 3.1.1 leading to the new equation:

(4.8) ω 2 0 M r d 2 y (1) dτ 2 + K r y (1) = -2ω 0 ω 1 M r d 2 y (0) dτ 2 -Λ b f nl .
Introducing the modal coordinates y (0) = Φz (0) and y (1) = Φz (1) into Eq. 4.8 where Φ =

[φ 1 | ... | φ p ], p being the dimension of K r and M r , and pre multiplying 4.8 by Φ H being the conjugate transpose of Φ, it comes :

(4.9) ω 2 0 Φ H M r Φ d 2 z (1) dτ 2 + Φ H K r Φz (1) = -2ω 0 ω 1 Φ H M r Φ d 2 z (0) dτ 2 -Φ H Λ b f nl .
Φ matrix is orthogonal in respect to K r and M r giving the following relations:

(4.10) Φ H M r Φ = m i j δ i j Φ H K r Φ = k i j δ i j
where δ i j represents the kronecker symbol, m i j and k i j represent respectively the modal mass and stiffness associated to M r and K r . Eq. 4.9 is a system of p decoupled scalar equation where the i th term is:

(4.11) ω 2 0 m ii d 2 z (1) i dτ 2 + k ii z (1) i = -2ω 0 ω 1 m ii d 2 z (0) i dτ 2 -φ H i Λ b f nl .
Assuming the harmonic regime as been established, we can write z (0) i = A i 2 e jτ + Āi 2 e -jτ where A i is an arbitrary amplitude and expand the nonlinear force into a Fourier series such as: (4.12)

f nl = +∞ n=-∞ c n e -jnτ
where c n represent the complex Fourier coefficients. Here only terms c 1 and c -1 are taken in consideration since all the other coefficients would provide secular terms leading to nonconvergent solutions, meaning that the nonlinear force considered here is written under the form f nl = c 1 e -jτ + c -1 e jτ . Assuming z (1) i is also written under the form z (1) i = B i 2 e jτ + Bi 2 e -jτ with B i an arbitrary amplitude, Eq. 4.11 becomes:

(4.13) (-ω 2 0 m ii + k ii )z (1) i = (ω 0 ω 1 A i m ii -φ H i Λ b c 1 )e -jτ .
The equation related to the complex conjugate terms of Eq. 4.13 is not taken in consideration since it would lead to exactly the same relation. Left-hand-side of Eq. 4.13 has the same linear kernel than Eq. 3.1, which means that right-hand-side of 4.13 has to be equal to 0 in order to avoid secular terms for the value of z (1) i [START_REF] Narisetti | A perturbation approach for predicting wave propagation in one-dimensional nonlinear periodic structures[END_REF]. This condition leads to:

(4.14) ω 0 ω 1 A i m ii -φ H i Λ b c 1 = 0
which can be rearranged to get the value of the corrected frequency ω 1 :

(4.15)

ω 1 = (Λ H b φ i ) H c 1 ω 0 A i m ii .
One can notice that Eq. 4.15 only contains terms that are already known using linear Floquet-Bloch analysis, excepted for c 1 coming from the nonlinear forces. This term can be either calculated analytically if the value of f nl is explicitly known or can be found using transient analysis, typically for cases involving finite element methods. The method developed to find this coefficient using a commercial software (ANSYS APDL) will be explained in chapter 4. Finally, the expression of the dispersion relation for nonlinear system using Lindstedt-Poincaré approach is:

(4. [START_REF] Campana | Impact of non-linear resonators in periodic structures using a perturbation approach[END_REF])

ω = ω 0 1 + ε (Λ H b φ i ) H c 1 ω 2 0 A i m ii +o(ε).

Inverse method for dispersion curve identification on finite structure

This subsection will present a new method introduced in this work to identify the dispersion curve of an infinite structure using a model of finite structure. Usually, identification methods are based on harmonic analysis, the structure is excited and the wave number is identified performing a spatial FFT to the deformed shape at a fixed time [START_REF] Manktelow | Analysis and experimental estimation of nonlinear dispersion in a periodic string[END_REF]. An inverse approach is used here since the wave number will be imposed and the frequency of the system is measured. The main advantage of this approach is to avoid space aliasing due to the discretization of the system, and allows to find directly the boundaries of the bandgap of a periodic structure just by imposing the wave number equal to the contour of the IFBZ (µ = 0 or µ = π) while the harmonic analysis can not predict at which frequency the band gap is appearing leading to a less precise estimation of the stop band boundary.

However, one have to notice that contrary to harmonic analysis, the inverse approach can not be used easily for experimental tests. So far, imposing a wave number to a structure is relatively though, so this technique can be mainly used as a numerical pre-test to identify at which frequencies the boundaries of the band gap should be localized or to validate infinite Let consider an arbitrary unit cell belonging to a periodic structure. Imposing an arbitrary value for µ 0 and using Eq. 3.7, one can find the eigenvalues ω 0,k of the system (K r (µ 0 ), M r (µ 0 ))

where k is the branch number of the dispersion curve observed. The eigenvector φ(µ 0 , ω 0,k ) contains the informations describing the motion of the unit cell through time. Considering a periodic structure containing N unit cells, the vector of global mode shape of the whole structure can be constructed concatenating the eigenvector multiplied by λ n = e n jµ 0 for n ∈ [0, N] :

(4.17) Θ = [φ | λφ | ... | λ N-1 φ | λ N q N bd ] T
where q bd is the vector of degrees of freedom at the boundaries of the unit cell which have not been reduced by the Floquet theorem (q L for 1D wave propagation and q L ,q B ,q LB for 2D wave propagation). Fig. 4.1 illustrates an example of initial conditions application in a 1D case with q bd = q 1 . The displacement value of all the degrees of freedom of the finite structure for a fixed value of wave number and frequency can be then written under the form:

(4.18) U(ω 0 , µ 0 , t) = A 2 Θe jω 0 t + Ā 2 Θe -jω 0 t .
Considering a real amplitude A, Eq. 4.18 can be simplified to:

(4.19) U(ω 0 , µ 0 , t) = A(ℜ(Θ) cos(ω 0 t) -ℑ(Θ) sin(ω 0 t))
ℜ and ℑ being respectively the real part and imaginary part function. For t = 0, the initial conditions that we have to impose to the full periodic structure is :

(4.20) U(ω 0 , µ 0 , t = 0) = Aℜ(Θ).
Once the initial conditions 4.20 have been applied to the finite structure for a value of amplitude A, the system is released for a certain number of periods which can be different depending on the studied system. After waiting a certain time, one node belonging to the periodic structure (ideally around the middle of the structure to avoid boundary effects) is chosen and a time analysis is performed in order to get the value of the frequency. For linear systems, a Fast

Fourier Transform [START_REF] Oberst | The fast fourier transform[END_REF] (FFT) can be used to get the value of the frequency. However, release the structure can present drops of amplitude as the system lose energy through time, meaning that the global amplitude might not remain constant through time ; this can cause issues for the study of nonlinear systems since the frequency depends to the amplitude. For that reason, it can be preferable to perform an Hilbert Transform [START_REF] Feldman | Hilbert transform in vibration analysis[END_REF] (noted H ) to the time signal to get the value of the average amplitude and frequency :

(4.21) Û(t) = H (U(t)) = α(t)e jβ(t)
where α(t) represent the instantaneous amplitude or envelope of U(t) and ψ(t) the instantaneous phase. The instantaneous frequency f i can be derived from the the instantaneous phase such as :

(4.22) f i (t) = 1 2π dβ dt (t).
The average amplitude and frequency of U(t) can then be found by doing an average of the instantaneous amplitude and frequency during the period T of signal analysis :

(4.23) α avg = 1 T T 0 α(t) dt , f avg = 1 T T 0 f i (t).
Repeating the operation for all the values of µ 0 and ω 0,k will give a panel of values f avg (µ, ω 0 ) which can be plotted to represent the apparent dispersion curve for a finite structure on the nonlinear case, and be eventually be compared with the Lindstedt-Poincaré method to confirm theoretical analyses.

Periodic discrete systems including resonators

Here is presented the case of a simple periodic resonating system using spring and masses. An introduction presenting the linear system is done to understand how behave the system and what type of dispersion curve we obtain on that case, then a non-linear duffing spring is added and a Lindstedt-Poincaré analysis is performed to observe what is the impact of the geometric nonlinearity on the wave propagation analysis. Finally, the results are compared with the inverse method applied to the finite structure developed in the thesis to validate the Lindstedt-Poincaré method.

Linear spring-mass resonator

Let us consider a one-dimensional principal periodic chain made of spring-mass units (mass m and stiffness k). A resonator represented by another spring-mass system (mass m R and stiffness k R ) is attached to each mass. This chain is infinite, and an appropriate unit cell is defined to apply the Floquet-Bloch Theorem (Fig. 4.2).

By writing the equation of motion of the unit cell linked to its neighbours we obtain the following expression:

(4.24)    m d 2 u j dt 2 + k(2u j -u j-1 -u j+1 ) -k R (q j -u j ) = 0 m R d 2 q j dt 2 + k R (q j -u j ) = 0
where u j and q j respectively denote the displacements of the j th principal mass and the j th resonating mass. In harmonic regime the displacements can be expressed as x j = X j e iωt , with i being a complex number such as i 2 = -1. Eq. 4.24 becomes:

(4.25) -ω 2 mu j + k(2u j -u j-1 -u j+1 ) -k R (q j -u j ) = 0 -ω 2 m R q j + k R (q j -u j ) = 0
The term q j depends on u j from the second line of Eq. 4.25. By applying the Floquet-Bloch

Theorem to the first line (u j+1 = e jµ u j , µ being the reduced wave number) one can obtain the following dispersion equation:

(4.26) -mω 2 + 2k(1 -cos(µ)) + (k R - k 2 R k R -m R ω 2 ) = 0
The roots of Eq. 4.26 are ω 1 and ω 2 , representing the two branches of the dispersion curve: ) The FRF exhibits the presence of a frequency region in which the amplitude is considerably attenuated. This area represents the band gap (resonant in this case). A resonant band gap can be differentiated from a Bragg one by observing the boundaries of the propagating part [START_REF] Sharma | Local resonance and bragg bandgaps in sandwich beams containing periodically inserted resonators[END_REF]. In the resonant case, a wave will stop propagating for a value of µ = π and then propagates again for µ = 0. The boundaries of the band gap can be expressed by calculating the values of ω 1 (µ = π) and ω 2 (µ = 0), leading to the following expression:

ω 1 = (m+m R )Ω 2 R +4k sin 2 (µ/2)-((m+m R )Ω 2 R +4k sin 2 (µ/2)) 2 -16km sin 2 (µ/2)Ω 2 R 2m ω 2 = (m+m R )Ω 2 R +4k sin 2 (µ/2)+ ((m+m R )Ω 2 R +4k sin 2 (µ/2)) 2 -16km sin 2 (µ/2)Ω 2 R 2m with Ω R = k R m R
ω 2 g = (m+m R )Ω 2 R +4k-((m+m R )Ω 2 R +4k) 2 -16kmΩ 2 R 2m ω 2 d = Ω 2 R (1 + m R m )
with ω g and ω d being the pulsation of the left and the right boundary, respectively. It is also important to notice that the value of Ω R is bounded by the values of ω g and ω d , and this has consequences on the value of the amplitude of the resonant mass. The rewriting of the second line of Eq. 4.25 leads to

q j = Ω 2 R Ω 2
R -ω 2 u j , and indicates that the value of the amplitude of q j strongly depends on the value of Ω R for a fixed value of the u j amplitude. The term q j is bounded in the propagative zone, since ω -< ω g < Ω R < ω d < ω + , with ω -and ω + representing respectively the value of the pulsation before ω g and after ω d . This property of a resonating system is different from spring-mass systems with Bragg bandgaps, like the diatomic mass system of Hussein et al. [START_REF] Hussein | Dynamics of phononic materials and structures: Historical origins, recent progress, and future outlook[END_REF]. This particular periodic configuration shows that the amplitudes of the two masses depend directly to the boundaries of the bandgap due to the Bragg effect, which is not the case for the resonating system where the boundaries depend on the value of the eigenfrequency of the resonator, changing the way of evaluating the maximum admissible amplitude to stay in low-level vibration domain.

Nonlinear spring-mass resonator system 4.2.2.1 Perturbation method for identification of nonlinear dispersion diagram

A Duffing spring mass resonator is added between the masses m and m R . In this case the nonlinear force is expressed as:

f nl = k R (q j -u j ) + εΓ(q j -u j ) 3
where Γ is the nonlinear stiffness constant and ε a small perturbation parameter. Replacing the linear restoring force -k R (q j -u j ) by f nl in Eq. 4.24 leads to:

(4.30)    m d 2 u j dt 2 + k(2u j -u j-1 -u j+1 ) -k R (q j -u j ) -εΓ(q j -u j ) 3 = 0 m R d 2 q j dt 2 + k R (q j -u j ) + εΓ(q j -u j ) 3 = 0
The objective here is to identify the correction term in the dispersion relation using the Lindstedt-Poincaré adapted method developed in [START_REF] Narisetti | A perturbation approach for predicting wave propagation in one-dimensional nonlinear periodic structures[END_REF]. The main goal is to generate a first order asymptotic development on the displacements of the masses and on the frequency, corresponding to the following expressions: (4.31)

u j = u (0) j + εu (1) j + O(ε 2 ), q j = q (0) j + εq (1) j + O(ε 2 ), ω = ω 0 + εω 1 + O(ε 2 ).
Replacing expressions 4.31 in Eq. 4.30 gives, after development and separation of the equations in the ε 0 and ε 1 orders: (4.32) ε 0 :

     ω 2 0 d 2 u (0) j dτ 2 + (2u (0) j -u (0) j-1 -u (0) j+1 ) -β(q (0) j -u (0) j ) = 0, κ 2 ω 2 0 d 2 q (0) j dτ 2 + (q (0) j -u (0) j ) = 0, (4.33) ε 1 :      ω 2 0 d 2 u (1) j dτ 2 + (2u (1) j -u (1) j-1 -u (1) j+1 ) -β(q (1) j -u (1) j ) = -2ω 0 ω 1 d 2 u (0) j dτ 2 + Γ(q (0) j -u (0) j ) 3 κ 2 ω 2 0 d 2 q (1) j dτ 2 + (q (1) j -u (1) j ) = -2κ 2 ω 0 ω 1 d 2 q (0) j dτ 2 -Γ β (q (0) j -u (0) j ) 3 (4.40) ω 2 0 d 2 u (1) j dτ 2 + (2u (1) j -u (1) j-1 -u (1) j+1 ) -β κ 2 ω 2 0 1 -κ 2 ω 2 0 u (1) j = F(τ)
with:

(4.41)

F(τ) = 2ω 0 ω 1 α + (1 -κ 2 ω 2 0 ) 2 (1 -κ 2 ω 2 0 ) 2 u (0) j -Γ κ 2 ω 2 0 1 -κ 2 ω 2 0 4 (u (0) j ) 3
Replacing u (0) j by its expression 4.34, the equation becomes:

(4.42) F(τ) = ω 0 ω 1 α + (1 -κ 2 ω 2 0 ) 2 (1 -κ 2 ω 2 0 ) 2 A u - 3 Γ 8 κ 2 ω 2 0 1 -κ 2 ω 2 0 4
A 2 u Āu e iτ e i jµ + d 1 e 3iτ e 3i jµ

The terms d 1 are associated to the 3 rd order of the nonlinearity. The linear kernel of Eq. 4.40

is similar to the one in Eq. 4.26; this implies that one needs to have all the coefficients in e i jµ be equal to 0 not to obtain a secular term in the temporal expression of u (1) j . By imposing this condition and rearranging the equation, we obtain the following expression for ω 1 : (4. [START_REF] Kittel | Introduction to solid state physics[END_REF])

ω 1 = 3 Γ|A u | 2 8 κ 8 ω 7 0 (1 -κ 2 ω 2 0 ) 2 (α + (1 -κ 2 ω 2 0 ) 2 )
Equation 4.44 leads to the establishment of the final expression describing the corrected dispersion relation for the periodic structure:

(4.44) ω = ω 0 + ε 3 Γ|A u | 2 8 κ 8 ω 7 0 (1 -κ 2 ω 2 0 ) 2 (α + (1 -κ 2 ω 2 0 ) 2 )
Eq. 4.44 illustrates the importance of the pulsation of the resonator to determine the value of ω 1 . If κ 2 ω 2 0 approaches 1 (i.e., ω approaches the resonance Ω R ), the value of ω 1 will increase and could become larger than ω 0 . This would however contradict the hypothesis underlying the perturbation method used here.

For the following analyses, the parameters α,κ, Γ are considered to be equal to 1 (unless Γ = 0, and this to obtain the linear case). A graphic representation of the corrected dispersion curve considering A u = 1 is shown Fig. 4.6.

Domain of amplitude validity

The perturbation analysis assumes an asymptotic development of the term ω to obtain the nonlinear corrected result. However, the asymptotic development of the term ω 2 = (ω 0 + εω 1 ) 2 = Dispersion curve using the perturbation approach in the nonlinear case for A u = 1 and ε = 0.05

ω 2 0 + 2εω 0 ω 1 leading to Eq. 4.33 is true only if 2εω 0 ω 1 ω 2 0 i.e. 2ε ω 1 ω 0 1
. ω 1 is replaced by its expression 4.43 leading to the following condition:

(4.45) 3|A u | 2 ε 4 ω 6 0 (1 -ω 2 0 ) 2 (1 + (1 -ω 2 0 ) 2 )
1

By assuming that the order of magnitude of is the same as ε, the following condition for the amplitude can be written:

(4.46) |A u |h(ω 0 ) < 1 where h(ω 0 ) = 3 4 ω 6 0 (1-ω 2 0 ) 2 (1+(1-ω 2 0 ) 2 )
. Plotting the function 1 h will give the maximum admissible value that |A u | can assume to satisfy the equation 4.46. Fig. 4.7 shows that for values of µ up to 1.0, almost any value of |A u | would not affect the dispersion in the nonlinear regime. This means that the hypothesis of linearisation in this area remains true even for very high values of amplitude. However, the admissible amplitude is much lower when the wave number increases and approaches a maximum theoretical limit of 0.48. The opposite behaviour however happens in the right branch (Fig. 4.8); in this case it is possible to observe a limit for |A u |, with a minimum value around 0.6 for µ = 0 and reaching 1.8 for µ = π.

Inverse method

In general, techniques used to identify the apparent wave number of a periodic structure are based on applying an harmonic excitation to one of the masses of the periodic structure and then evaluating the response by measuring the wave number after a certain time [START_REF] Lazarov | Low-frequency band gaps in chains with attached non-linear oscillators[END_REF][START_REF] Manktelow | Finite-element based perturbation analysis of wave propagation in nonlinear periodic structures[END_REF]. This paper introduces here an alternative numerical method to compare the results obtained from the theory. The idea underpinning this methodology is to impose the wave number of the whole structure as an initial condition, releasing the system and then observing how it evolves in time to obtain the frequency of the corresponding imposed wave number. One of the advantages of this novel methodology is that a time-domain analysis is required, rather than the typical spacedomain one to obtain the wave number by using direct methods [START_REF] Berthaut | K-space identification of apparent structural behaviour[END_REF]. The use of a time-domain analysis helps to reduce the error made due to the spatial discretisation of the periodic structure.

In particular for spring-mass systems, the distance between two consecutive masses does not need to be represented in the model, and this allows to fix the element length to unity. According to the Shannon theorem, a minimum of 2 elements per wavelength is necessary to avoid space aliasing. However, to represent a sinusoidal excitation in a correct way, it is common to use at least 6 elements per wavelength, and that condition can not be reached when the wave number exceeds 1. Fig. 4.9 illustrates this problem by showing the space representation for two different wave numbers. It is however important to note that this numerical method is representative of a theoretical experimental setup, since imposing a fixed shape to the whole structure as an initial condition may be impractical.

The finite structure here is represented by an assembly of 300 unit cells in a linear chain.

The equation of motion of the finite structure with a free-free boundary condition is written and solved using the ODE45 function in Matlab. A perfectly matched layer (PML) is also used at the boundaries of the structure to avoid wave reflection [START_REF] Berenger | A perfectly matched layer for the absorption of electromagnetic waves[END_REF]. The PML consists in a viscous damping force c(x p ). ẋp applied to all the masses of the principal chain under the form:

(4.47) c(x p ) = Ce 1-1 (1-xp ) n
with C a constant positive value, x p ∈ [0; 2] a value interpolating the mass at the p th position such as c(0) = c(2) = C and c(1) = 0, and n an even number representing the order of the PML.

The higher the value of n is, the less the masses close to the boundaries will be affected by the damping. Fig. 4.10 shows the trend of the PML for two different values of n and C = 2.

The resolution of the linear system (i.e. with Γ = 0) has to be performed in order to obtain the values of ω 0 (µ) corresponding to the dispersion curve of the linear case first. This step is (4.48)

A q = A u 1 -κ 2 ω 2 0 .
Then, the initial conditions applied to the whole structure are written as:

(4.49)                  U j (t = 0) = A u cos(µ(x j -x 0 )) U j (t = 0) = 0 Q j (t = 0) = A q cos(µ(x j -x 0 )) Q j (t = 0) = 0 j ∈ [[1; 300]]
The terms U j and Q j represent the values of the displacements of the j th principal and resonating mass of the system, respectively. The location of the mass considered to obtain the time-displacement response is x 0 . The value of x 0 has to be chosen carefully to obtain a coherent representation of the results: if the chosen mass is too close to the boundaries, the PML and potential residual wave reflections will not provide correct values of amplitude. To avoid that issue, the observed mass should be located in the middle of the structure. The writing of the initial spatial displacement as in Eq. 4.49 ensures that the amplitude A u will be observable for the mass located in x 0 for any value of µ, when x j = x 0 . This condition means that the imposed initial displacement has a phase shift of φ = µx 0 , and that will not affect the results.

After the masses have been released for an imposed wave number µ 0 and a waiting time of t = 25T 0 such as T 0 = 2π ω 0 (µ 0 ) , the frequency of the time domain signal from the mass located in x 0 4.2. PERIODIC DISCRETE SYSTEMS INCLUDING RESONATORS is measured. One example of input and output signal obtained for a µ = 0.1 on the right branch of the dispersion curve is shown Fig. 4.12. Fig. 4.13 shows that the amplitude does not remain constant and equal to A u with time. Consequently, the Fast Fourier Transform commonly used to obtain the frequency of a periodic signal would not give accurate results here, since the signals possess variations of amplitude and therefore variations of the frequency due to the nonlinear behaviour of the resonators. An average of the periods of the signal could also be performed to obtain the global value of the frequency. Depending on the imposed wave number, the value of the output average amplitude will however be different. It will be therefore not possible to compare effectively those results with the ones from the perturbation method, which assumes a fixed amplitude for every wave number. To fix this problem, an algorithm is implemented, which consists in iterating the simulations starting with an imposed amplitude A u 0 . The goal of the algorithm is to reach the targeted value A u by comparing it with the output amplitude of the system (noted A out ). After an iteration, the value of A out is estimated by performing a Hilbert transform of the signal:

(4.50) A out = |H (U 0 (t))| If A out < A u ,

Numerical results

Before obtaining the numerical nonlinear dispersion curves, a verification of the condition estimated in subsection 4.2.2.2 has to be performed to verify that the domain of validity of the amplitudes is consistent with the one of the theoretical model. This verification is done by launching several numerical simulations with a variation of the amplitude and the wave number in the following domains:

A u ∈ [0.1 , 0.6], µ ∈ [1 , π] for the left branch and A u ∈ [0.1 , 1], µ ∈ [0 , π]
for the right branch. The values of wave number are not considered in the domain [0; 1] for the left branch since the behaviour of the structure is almost linear in this area, even for high values of amplitude (see Fig. to observe the difference between the numerical method and the perturbation approach. One can observe that for this particular value of amplitude, the numerical and theoretical results provide a close match, and the dispersion curve is shifted to higher frequencies compared to the linear case. This result is consistent with the fact that a nonlinear cubic spring has been added to the model, hence stiffening the structure and increasing value of eigenfrequencies. A result imposing an amplitude for an high amplitude (A u = 1.5) and the time signal associated to the 4.23 demonstrates that the results will be inaccurate up to a certain value of wave number (around 0.8 in that case), but they will still remain correct after this value. Fig. 4.24 confirms that the perturbation approach will no longer be accurate after a certain value of amplitude.

Comparison with HBM and harmonic analysis

Although the theoretical method presented above provides quite accurate results for relatively low values of amplitude, a comparative study is performed with previous techniques shown in the literature for further benchmark. Another way to obtain the dispersion relation for nonlinear structures is by using the Harmonic Balance Method on the relative displacement between the resonant mass and the principal one noted v j = u j -q j [START_REF] Lazarov | Low-frequency band gaps in chains with attached non-linear oscillators[END_REF][START_REF] Wang | Lower band gaps of longitudinal wave in a one-dimensional periodic rod by exploiting geometrical nonlinearity[END_REF]. One can assume that v j can be written under the following form:

(4.52) v j (τ) = A v e iτ + Āv e -iτ
With A v being the amplitude of the relative displacement. By rewriting Eq. 4.30 and using the HBM form Eq. 4.52 one obtains a new expression of the dispersion relation of the periodic resonant system as:

A u = 0.3, ε = 0.05) (4.53) cosh(µ) = 1 - ω 2 0 2 - ω 2 0 β(κ 2 + 3κ 2 Γ|A v | 2 ) 2(κ 2 + 3κ 2 Γ|A v | 2 -ω 2 0 )
The dispersion relation is here written in the direct form (i.e., frequency in input and wave number in output). Using Eq. 4.48 and knowing that the amplitude of the relative displacement can be written as A v = A q -A u one obtains:

(4.54) Eq. 4.54 indicates that for a constant imposed value of A u the amplitude A v does not remain constant but varies with the frequency. Thanks to Eq. 4.54, 4.53 can be used to compare the dispersion curves obtained with the perturbation method and the HBM by imposing a value of amplitude A u . In addition, a classical harmonic analysis has been performed on the numerical The perturbation method appears to provide a slightly better accuracy in this case, since both the harmonic analysis and the inverse approach are closer to that method compared to the HBM. The relative error in frequency between the two methods is however relatively small here (about 1%).

A v = κ 2 ω 2 0 1 -κ 2 ω 2 0 A u ,
The increase of the amplitude to |A u | = 1.2 shows however that the gap between the perturbation method and the HBM becomes larger (figure 4.26) and reaches a 6% error that shows the higher accuracy of the perturbation method according to the inverse approach. It appears that as the amplitude increase, the perturbation method provides more accurate results than the harmonic balance mainly due to the fact that in HBM only the relative displacement is expanded in Fourier series. Instead, the asymptotic development in the perturbation method is carried out both on the displacements and the frequency providing an overall improved approximation.

Conclusion

The results relative to the nonlinear dispersion curve exhibits a frequency shift on higher frequencies compared to the linear case, showing on that case an hardening effect. This effect could be predictable since the nonlinearity imposed for this model was relative to a nonlinear say that this hardening effect is not very satisfying since it is generally better for aerospace applications to lower the stop band to lower frequencies, therefore having a softening nonlinear effect. However, even if the nonlinearity considered here is weak, the study shows that this nonlinearity has to be taken in consideration when a periodic structure has to be designed, in order to have a correspondence between the stop band predicted by the model and the one we could observe experimentally, especially for resonant structures showing often high amplitude o displacements close to the peak of resonance and therefore geometrical nonlinearities. The Lindstedt-Poincaré method can be theoretically used to get nonlinear dispersion curves generalizing the example seen in previous section for models designed with the finite element method. This section will explain a numerical method developed during the thesis in order to get this dispersion curve using an interaction between the commercial software ANSYS APDL in combination with Matlab. All the relations written are done for a couple (ω 0 , µ 0 ).

where φ N is the normalized wave mode. By doing so, we can ensure that the maximum of φ N is equal to unity for each couple (ω 0 , µ 0 ) and therefore the choice of the maximum amplitude applied to the unit cell is constant for the transient analysis. The computation is performed taking in consideration 10 periods of T = 2π ω 0 to ensure that the harmonic regime is established. The temporal values of forces are then collected for each degrees of freedom of the structure.

Since ANSYS does not make the difference between linear and nonlinear efforts during the post-processing, the output obtained represent the sum of the two type of efforts. We note this quantity F tot = F l + F nl written under the following form: (4.58)

F tot = [F 1 y (t) M 1 θ (t) ... F n y (t) M n θ (t)] T .
A similar analysis is done taking in consideration the hypothesis the geometrical nonlinearities on ANSYS. The nonlinear forces are then found by subtracting the total forces previously calculated with the linear forces such as F nl = F tot -F l . An example of this procedure to get the nonlinear force related to the first degree of freedom F 1 y is given Fig. 4.28. One can notice in this example that the total forces and the linear force seem to be similar at first look, but after subtraction we actually observe that nonlinear force is simply present with much lower intensity than the linear one. This intensity does only depend on the amplitude imposed to build the vector U(t) and an arbitrary value A = 4.5 mm has been chosen here. Also, one can notice that the nonlinear force curve does have the shape of a function under the form: F nl (t)e -jω 0 t dt.

Finally, the value of corrected frequency is calculated using the relation 4.15, and the process is repeated for each value of ω 0 . The final corrected dispersion curve is shown Fig. 4.29 for the four first branches of bending wave modes. For this example, a total of 200 corrected frequencies An example of a 1D periodic structure embedding a resonator is studied in order to see the effect of the geometric nonlinear behaviour on the bandgap of the system. The periodic structure chosen corresponds sequence of periodic truss structures with a rectangular shape and modeled with beam elements, and is made as simple as possible to get a good visualisation of the resonator effect on the whole periodic system. A picture of the full periodic system and the unit cell chosen with the resonator is shown in Fig. 4.31 representing the real shape of the structure on the software and the geometric/material properties are shown in table 4.1. ρ R , E R and ν R appoint the material parameters for the resonator. The study is done considering only in plane waves in this case, meaning that only the degrees of freedom u x ,u y and θ z are conserved here. geometric nonlinearities are applied to the whole studied system. For that reason, it is important to design a resonator which has to have a stronger geometric nonlinear behaviour than the hosting structure in order to consider only a nonlinear resonant system. To this end, a much thin and long structure for the resonator is chosen in comparison with the host structure, explaining the choice of the values for h 2 and b 2 compared to h 1 and b 1 . Also, one has to consider that the amplitude of the resonator increases considerably close to the resonance, making the hypothesis of small perturbations not valid anymore. Therefore, the choice of the input amplitude for the analysis has to be chosen carefully to get coherent results even using Lindstedt-Poincaré method. The

[k g.m -3 ] 7800 ρ R [k g.m -3 ] 5ρ L 1 [mm] 200 h 1 [mm] 5 L R [mm] 100.5 E [GPa] 210 E R [GPa] 5E L 2 [mm] 20 b 1 [mm] 5 h 2 [mm] 1 ν 0.33 ν R ν L 3 [mm] 10 b 2 [mm] 1

Amplitude limit analysis

The value of amplitude limit has to be taken into account when the hypothesis of weak nonlinearities is considered in the analyses. If the amplitude is too high, it can induce large strain or deflection and may overpass the Lindstedt-Poincaré hypothesis, leading to wrong result for the comparison between the theoretical and the practical approach. A method is proposed here to determine the admissible limits of amplitude using the Lindstedt-Poincaré technique. The method is an extension of what has been observed for the discrete system in previous subsection.

Following the same considerations, the hypothesis are respected only if:

(4.61) 2ω 0 ω 1 < ω 2 0 .
Replacing ω 1 by its expression 4.15 leads to the following condition in amplitude:

(4.62) A > 2Λ b φ H c 1 m ii ω 2 0 .
Contrary to Eq. 4.46 seen in the discrete resonator analysis, there is no analytical expression for the amplitude limit since c 1 is known thanks to force extracted from the transient analysis seen in previous subsection, therefore c 1 depends also on the amplitude A. However, since only geometrical nonlinear analyses are performed here, the relation between c 1 and A should be such as c 1 ∝ A 3 . This verification can be done in a simple case performing a sequence of nonlinear analysis doing a variation of amplitude and checking the value of c 1 (A). Using the same model used in subsection 4.3.1, transient analyses are performed for µ = π 2 on the second branch of the dispersion curve, and the c 1 related to the force F 31 y is estimated by doing a variation of amplitude for A ∈ [0.1 , 5] mm. Fig. 4.34 represents the functions ℜ(c 1 ) = f (A) and ℑ(c 1 ) = f (A) superposed with y = αA 3 and y = βA 3 where α = -5.55 × 10 -5 N.mm -3 and α = 5.75 × 10 -8 N.mm -3 , and illustrates that the nonlinear force is indeed proportional the cube of the amplitude, meaning that we can write c 1 = γA 3 with γ = α + jβ. This observation may not be valid for any type of structure, and so such an analysis has to be done for every new case to treat. Moreover, this analysis has been done for a single value of force for a specific couple (ω 0 , µ 0 ), but we can suppose here for simplification that the relation remains true for every cases with eventually a different value for α. Under these considerations, Eq. 4.62 can be simplified leading to:

(4.63) A 2 < m ii ω 2 0 2Λ b φ H γ v
where γ v = [γ 1 , ..., γ p ] with p the number of degrees of freedom of the system. In theory one just need to compute all the values belonging to γ v in order to get the amplitude limit A l im :

(4.64)

A l im = m ii ω 2 0 2Λ b φ H γ v .
A practical way to find A l im is to use a test value with a random amplitude and combine the relations ω l im 1 = ω 0 2 and ω 1 = K A One can also calculate the value ω l im 1 associated to the amplitude limit A l im :

(4.66)

ω l im 1 = K A 2 l im .
Combining Eq. 4.65 and 4.66 with the condition ω l im 1 = ω 0 2 leads to the value of amplitude limit : (4.67)

A l im = A t ω 0 2ω t 1 .

The amplitude limit depends on the wave number and branch analyzed, meaning that it has to be computed for all the values of couple (ω i 0 , µ i 0 ) in the range of interest, and so the value of the global admissible value that A g l im has to be chose such as:

(4.68) A g l im = max i A i l im (ω i 0 , µ i 0 ) 0 1 2 3 4 5
Amplitude (mm)

Inverse Method

Like the Lindstedt-Poincaré method, the inverse method seen in previous section can be also theoretically implemented using finite element method in order to find the apparent dispersion curve of a structure. The principle is the same as seen for discrete systems, instead here the wave The computation of 3 branches of the dispersion curves is done using all the methods (Inverse WFEM, Lindstedt-Poincaré, inverse methods using FFT and Hilbert Transform) in order to be compared, considering an amplitude A = 3 mm which verifies the amplitude limit conditions presented in subsection 4.3.2 and ε = 1 for the Lindstedt-Poincaré method. Fig. 4.37 shows the results of all the dispersion curves obtained. One can observe that there is no much difference between using the FFT and the Hilbert Transform since the results are pretty close one to each other. However, the FFT seem to show more accuracy since it matches more with the Lindstet-Poincaré method, assuming that FFT is a better choice on that case. It is also noticeable that the inverse methods does not give always correct results at that stage, specifically in the first branch corresponding to lowest frequency. This is mainly due to the difficulty to analyse the amplitude time signal on some cases caused by the nonlinearity and the imposed mode shape as an initial condition, as well as the node observed. As an example, 2 different nodes are observed for the analysis done at the point ( f 0 , µ 0 )=(2811 Hz, 0.6491) corresponding to the result not following the trend of the other dots the second branch observed Fig. 4.37. The first node observed (n 1 ) is the one used as a reference for all the computation (middle of the structure) while the second node n 2 is at location n 1 + 2 from this one (after two elements) therefore very close to the first one. Time amplitude plots of the two nodes for u y displacement are shown Fig. 4.38, and show that the signal obtained is very different even for two close nodes. If a FFT is performed for the node n 2 , it appears that the frequency calculated is 3373 Hz which exactly corresponds to the result obtained with the Lindstedt-Poincaré method on the dispersion curve.

Conclusion

An idea of how to perform nonlinear analyses on periodic structures have been presented in this chapter. The theoretical approach is based on the Lindstet-Poincaré method and an study has been performed on a discrete resonant system and then extended with first investigations for finite element models. The theoretical approach used has been validated with an new approach of analysing the wave propagation in nonlinear structures applied to finite structures with the inverse approach, both applied for discrete systems and finite element systems. A way to estimate the amplitude limit of the methods has been described and the coherence of the method has been observed in the case of discrete resonant systems. Even if the first results obtained show good correspondence, improvements have to be done to both methods to be able to study more complex periodic structures. This chapter is treating about the limitations encountered for all those studies in a physical and numerical point of view, exposing some ideas that could be follow in future works related to these topics in order to obtain improved results.

Linear structure analyses

Methods used

Concerning the part of linear periodic structures, all the methods used have been already known and developed during the past years, mainly for the Inverse WFEM or the branch tracking. Only the modified version of the Transfer Matrix Method could be presented as an improved method since no really clear applications have been done yet using this approach in the case of vibrations of periodic structure. Even if the ill-conditioning problem linked to the exponential values has been solved with the new formulation (see subsection 3.1.3.2), some other ill-conditioning problems still persist making impossible the study of some specific periodic structures. The persisting problem was mainly due to the number of degrees of freedom at the boundaries of the structure. Indeed, internal degrees of freedom have been reduced using a wave finite element method CMS using Craig-Bampton [START_REF] Hussein | Bloch mode synthesis: Ultrafast methodology for elastic band-structure calculations[END_REF][START_REF] Zhou | Wave finite element method based on reduced model for one-dimensional periodic structures[END_REF] leading to a very few number of degrees of freedom, but still remained the problem of external degrees of freedom which still let some troubles in the finite structure response calculation. As an example, the study of the unit cell modeled with shell elements proposed in section 3.3 has been done using two different type of meshing for the finite structure and infinite structure cases, reducing the number of degrees of freedom for the TMM (see Fig. 5.1). We can observe after remeshing that the number of nodes has been reduced to 2 (so 12 dofs) at each side of the single unit cell, multiplied by the number of repetition in the y direction (leading to a maximum of 96 dofs when the cell is repeated 8 times). In comparison, if the same meshing was done using Fig. 5.1b, the TMM computation would have lead to wrong results because of the too numerous number of external nodes. These remarks and adjustment made us realize than even before doing the wave propagation study, the way of how to select the unit cell and the type of meshing, as well as the placement of the resonator has to be done very carefully in order to get converging results.

However, this technique is no longer efficient when skins are added to the core or when volume elements are used, since too many nodes have to be taken in consideration making the TMM not interesting anymore. Since the goal of the work was not focused on model reduction or ways to improve the TMM that way, it has been judged not necessary to improve the developed tool and simply use is for simple cases involving few external degrees of freedom. Instead, classical frequency response on finite structure can be performed directly using a commercial software to validate the results obtained with the WFEM.

Results obtained with resonators

Some results have been obtained for wave propagation using resonators in different type of models, and the objective was double concerning the studies. The first objective was to have an overview of different techniques we can apply to the analysis of resonating system, using different type of shapes and model for the resonators going from simple 1D system, representing the theoretical model to complex 3D structures much closer to the reality. Doing so, one could have an idea of how the dispersion curve is evolving taking in consideration certain types of element implying certain different waves, or the effect of the skins on the wave propagation behaviour and therefore on the stop band analysis. The second objective was the study of the parameter of the resonators to get improved results in term of stop bands located in low frequency, by showing an example of how to make a trade off between the added mass of the resonator versus vibroacoustic performances using beam elements, or how to improve these results keeping a constant mass by changing the geometrical parameters of the resonators through the volume element analyses.

The design of the resonators can be expanded to try to improve the results. Indeed, depending on the approach and the constraints imposed, several ways of studies can be carried out, using for example multiple resonators in one single unit cell with different materials of geometric properties to get larger stop bands. Another idea would be to create super cells of resonating system [START_REF] Henneberg | Periodically arranged acoustic metamaterial in industrial applications: The need for uncertainty quantification[END_REF] to have different material and/or geometric properties for each resonator, keeping the same number of resonator in each unit cell compared to the studies done here. The main problems of those approaches is that their require to build models with a very large number of degrees of freedom,

Applications and limits of nonlinear models proposed

Lindstedt-Poincaré Method

The nonlinear analyses done here were a combination between a first approach of studying nonlinear resonators using the Lindstedt-Poincaré method in an analytical way and using a finite element software (Ansys) for more general situations, which required an algorithm development to get an automatic method allowing us to get the dispersion curve, even though a very close approach, but not identical, has been explored before. The idea here was to expand the method combining the Inverse WFEM and transient analyses in one unit cell only instead of the unit cell and 2 neighbours (for 1D case) as done in the literature. However, several issues occurred during the algorithm development which could not be anticipated easily before starting to work on it, leading to sometimes not straightforward analyses of structures and might be improved in the future.

One of the most important issue encountered when testing the program was the apparition of secular nonlinear forces taken from ANSYS for some degrees of freedom belonging to the structure at certain values of (µ, ω), leading to a wrong calculation to the c 1 value for some points of the dispersion curve. This problem occurred for example when trying to get the out-of-plane dispersion curve of the resonant unit cell analysed in Fig. 4.31, getting some infinitely increasing good results meaning that the method is theoretically correct but some cases can not be treated (mostly because of the design of the model) or that the algorithm has to be adjusted to cover a larger panel of those cases.

Inverse approach

The inverse approach is another way to find the dispersion curve using a finite structure, and present some advantages but also drawbacks compared to the direct approach consisting on performing an harmonic or sweep sine analysis. As we could see in the finite element method example, it provides quite good results since it matches almost perfectly with the Lindstedt-Poincaré method for the case of the simple beam, and some of the imperfections observed come from the fact that the analysis is done on a finite structure, and so what we get is an apparent 92 5.2. APPLICATIONS AND LIMITS OF NONLINEAR MODELS PROPOSED dispersion curve and not a real dispersion curve.

However, the method is currently poorly working for more complex cases, even by taking a simple beam with different sections as we did for the Lindstedt-Poincaré analysis (see Fig. 4.29).

The issues are most likely due to the method used to export and treat the results. As we could see in the analysis, most of the outputs calculated were matching with the dispersion curve, but some of them were inaccurate, mainly because a wrong node was observed for a certain wave number

value. An idea to improve the result would be to import the temporal displacement of a group of nodes located around the middle of the structure and compute the mean of the results to get a more precise frequency in output. Even if the FFT or the Hilbert Transform gave quite precise results, the using of signal processing method such as the Wigner-Ville [START_REF] Ville | On the quantum correction for thermodynamic equilibrium[END_REF][START_REF] Wigner | On the quantum correction for thermodynamic equilibrium[END_REF] distribution can be an alternative more adapted for the studies since it would give more details about the distribution of the frequency value over the time instead of a deterministic value like the Hilbert Transform.

CHAPTER 6. CONCLUSIONS AND FUTURE WORKS the type of resonators selected for the study caused by a coupling wave effect between in-plane and out-of-plane propagation. Volume elements have been used to have a better control design of the resonators concerning their shapes, in order to have more flexibility on the geometrical parameters for parametric analyses. A study comparing the location and number of resonators has been made showing that the most efficient way to get the widest band gaps at low frequency for the same mass is most likely using a resonator in each re-entrant pattern as it was for the beam element model. A succession of parametric studies changing the geometrical features of the resonators in a specific range of values to compute the best result, and shows that the width of the stop bands can not be found as intuitively as it was done for the mid frequency location which directly depends on the stiffness over the mass ratio of the resonator.

The second part of the thesis was about the development and study of nonlinear periodic systems with and without resonators. Studies have been done starting with the simple one dimensional resonant spring mass resonator which has been subject of a journal publication [START_REF] Campana | Impact of non-linear resonators in periodic structures using a perturbation approach[END_REF],

where the objectives were multiple. One of the objective was to show that the choice of amplitude to perform the study has to be done very carefully for such a locally nonlinear system since the boundaries of the band gap are dependant on the parameters of the resonator. According to the position of the eigenfrequency of the resonator compared to the stop band boundaries, the choice of the amplitude has to be capped to a maximum value, otherwise it could lead to a non respect of the Lindstedt-Poincaré hypothesis by going into strong nonlinearities assumptions.

The second objective was the introduction of the inverse method, showing close results and confirming the amplitude limit analyses studies imposed by the Lindstedt-Poincaré method. Even if this technique is not applicable by performing experimental tests, it is very suitable to quickly measure the value of the stop band in the nonlinear dispersion curve since the wave number can be imposed in that case. The methodologies have been then expanded to the finite element method, creating an algorithm to compute nonlinear dispersion curves using ANSYS APDL as Finite Element software for both Lindstedt-Poincaré method and the inverse approach, showing convincing results for the studied structures.

Future works

The study of resonators in auxetic structures can be widely investigated for both linear and nonlinear purposes, and some main ideas can be explicited as future works and perspectives to complete the work done in this thesis.

• Different shapes of resonators to get a various panel of results, also attaching them on the core of the periodic structure instead of attaching them in the shells, even it would most likely impact the result for higher range of frequencies since the skins modes are 96 6.1. FUTURE WORKS not considered in the lowest frequencies. Damping can also be considered locally to the resonators to see the influence on the dispersion curve.

• Vibroacoustics tools such as the Transfer Matrix Method could be enhanced in order to be applicable to complex structures including skins, allowing efficient optimization results on a very wide range of parameter choice for the resonators which could not have been done in the thesis since the goal was not focused on this kind of topic.

• Manufacturing the samples using the Kirigami technique [START_REF] Neville | Shape morphing kirigami mechanical metamaterials[END_REF] and creating the resonators using steel or aluminium materials for some applications and performing vibration experimental tests using a shaker to excite the out-of-plane waves of the structure for different shape of resonators in order to compare the results with the dispersion curves obtained with volume element.

• Nonlinear methods applied to finite element methods can be improved in order to work for more complex cases by improving the algorithm in more specific cases, even if the automation the method seems very challenging since the computations can take a lot of time depending on the size of the models, and sometimes does not seem to converge for specific geometries. One can observe on those examples that results related to the resonator A and B are potentially presented stop band located for high frequencies for the first case and relatively low frequencies for the second case, while no clear wave propagation improvement seem to be made with the resonator C, meaning that the effect of the resonator looks to be more impacting if they are attached to the core. This phenomena can be explained since the first frequency modes of the sandwich panel are related to the core of the structure only, meaning that the resonator will not be efficient for low frequency areas because they will not interact with the traveling waves. Also, the dispersion contour Abstract : In the domain of mechanical engineering, control of noise and vibration is one of the most important problem that engineers have to deal to the protection of structures. During the few past years, solutions have been developed to reduce structural vibrations using periodic structures. The main interest of using such structures is linked to the concept of band gaps, areas of frequency where the amplitude of oscillations of structures is strongly attenuated with the creation of evanescent waves caused by the Bragg effect, making them good candidates to realise stop band filters. The work presented in this thesis focuses on the study of inserts in honeycombs structures presenting a negative Poisson's ratio or auxetic. Inserts used that way are called resonators and represent an alternative to create band gaps in low frequency range.

The first part of the thesis is consecrated to the study of different type of resonators using different finite element models to analyse the evolution of wave propagation through the periodic structure, changing the material properties or the number of resonators used. Since resonators are used as dynamic absorbers and therefore are subject to high vibration amplitudes, the small displacement hypothesis used in linear structural analysis might not be available any more in majority of cases, creating geometric nonlinearities. For that reason, the second part of the thesis is devoted to the study of non-linear periodic structures, having principally a Duffing behaviour in order to get a representation of a dispersion diagram depending on the resonator's amplitude of oscillation.

Université Bourgogne Franche-Comté 32, avenue de l'Observatoire 25000 Besançon

  THESE DE DOCTORAT DE L'ETABLISSEMENT UNIVERSITE BOURGOGNE FRANCHE-COMTE ET DE L'UNIVERSITE DE BRISTOL PREPAREE EN CO-TUTELLE AVEC l'UNIVERSITE DE BRISTOL

3 . 1 41 3. 6 46 4. 1 FIGURE Page 2 . 1 4 2. 2 7 2. 3 8 2. 4 9 2. 5 10 2. 6 13 3. 1 16 3. 2 24 3. 3 25 3. 4 25 3. 5 26 3. 6 27 3. 7 Forced 27 3. 8

 314164612142738495106131162243254255266277278 FIGURE Page

3. 22 1 . 1

 2211 Dispersion corresponding to the red configuration (a) and the blue configuration (b) for an equivalent total mass . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.23 Length limit L s + L m < H s tolerated to stay inside the unit cell contour . . . . . . . . . 3.24 Constant mass stop band evolution performing parametric analysis of the resonator for each couple of value (h s , b s ) presented table 3.6 . . . . . . . . . . . . . . . . . . . . . 3.25 Values of the criteria f for different geometrical values of the resonator. The highest value is obtained for the configuration (5,7) with f = 0.0778. . . . . . . . . . . . . . . . 3.26 Dispersion curve for the best value corresponding to the best result obtained in the parametric study and the associated resonator shape . . . . . . . . . . . . . . . . . . . . 4.1 Initial conditions application for inverse method in a finite structure containing N unit cells . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2 1D resonant spring mass system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3 Dispersion curve for the infinite structure. Propagative part in red (ℜ(µ)) and evanescent part in blue (ℑ(µ)) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.4 Forced Response Frequency at one of the extremity of a finite structure of 20 unit cells 4.5 1D nonlinear resonant spring mass system . . . . . . . . . . . . . . . . . . . . . . . . . . 4.6 Dispersion curve using the perturbation approach in the nonlinear case for A u = 1 and ε = 0.05 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.7 Value of 1 h(ω 0 ) for the left branch of the dispersion curve (ω (1) 0 ) . . . . . . . . . . . . . . . 4.8 Value of 1 h(ω 0 ) for the right branch of the dispersion curve (ω (2) 0 ) . . . . . . . . . . . . . . 4.9 Spatial representation of masses for two different values of µ . . . . . . . . . . . . . . . 4.10 Representation of the PML . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.11 Flowchart of the algorithm to get the right amplitude A out . . . . . . . . . . . . . . . . x LIST OF FIGURES 4.12 Initial conditions of the system . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 4.13 Time signal obtained for the mass x 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 4.14 Initial conditions applied to the mass of the principal chain after algorithm . . . . . . 67 4.15 Time signal obtained for the mass x 0 after algorithm . . . . . . . . . . . . . . . . . . . . 68 4.16 Value of 2 ω 0 ω 1 for the left branch of dispersion curve (perturbation method error) . . . . 68 4.17 Model error estimation for the left branch of the dispersion curve . . . . . . . . . . . . 69 4.18 Value of 2 ω 0 ω 1 for the right branch of dispersion curve (perturbation method error) . . . 69 4.19 Relative frequency difference between numerical and perturbation method . . . . . . 70 4.20 Nonlinear dispersion curve (A u = 0.3, ε = 0.05) . . . . . . . . . . . . . . . . . . . . . . . . 70 4.21 Zoom on the left branch of the dispersion curve (A u = 0.3, ε = 0.05) . . . . . . . . . . . 71 4.22 Zoom on the right branch of the dispersion curve (A u = 0.3, ε = 0.05) . . . . . . . . . . 71 4.23 Right branch of the dispersion curve for A u = 1.5 . . . . . . . . . . . . . . . . . . . . . . 72 4.24 Time signal for µ = 0.3 and a targeted amplitude A u = 1.5 . . . . . . . . . . . . . . . . . 72 4.25 Comparison between Perturbation, HBM, inverse approach and harmonic analysis for |A u | = 0.5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 4.26 Comparison between Perturbation, HBM and inverse approach for |A u | = 1.2 . . . . . 74 4.27 Beam unit cell used for the nonlinear analysis . . . . . . . . . . . . . . . . . . . . . . . . 75 4.28 Steps to get the nonlinear force related to F 1 y . . . . . . . . . . . . . . . . . . . . . . . . . 77 4.29 Corrected dispersion curve obtained after Lindstedt-Poincaré Method application (blue) and linear dispersion curve (red) for A = 4.5 mm and ε = 1 . . . . . . . . . . . . . 78 4.31 Resonant structure studied (a) finite structure (b) unit cell considered with the resonator in red . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 4.30 Algorithm used to get the corrected dispersion curve using Lindstedt-Poincaré method 79 4.32 Linear dispersion curve of the resonant system (a) Mode shape for the first branch (b) and the second branch (c) close to the resonance . . . . . . . . . . . . . . . . . . . . . . . 80 4.33 Comparison between linear and nonlinear dispersion curve (a) values of |c 1 | for µ = 1.8634 at the first branch of dispersion curve (b) . . . . . . . . . . . . . . . . . . . . 81 4.34 Comparison between |c 1 | = f (A) and the function y = αA 3 with al pha = 5.55 × 10 -5 N.mm -3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 4.35 Initial conditions applied to the structure using ANSYS (a) and response obtained for the middle node (b) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.36 Results obtained to get the nonlinear frequency using the FFT (a) and the Hilbert Transform (b) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 4.37 Dispersion curves using all the methods explained in the research applied to a simple beam . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 4.38 Time signal obtained after releasing at 2 different nodes for the second branch of dispersion curve computed at µ = 0.6491 . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 5.1 Different type of meshing for the TMM unit cell (a) and the inverse WFEM (b) . . . . 91 xi INTRODUCTION VIPER Project This research is part of the Marie Skłodowska Curie Project under Horizon 2020 European Funding Program VIPER (Vibroacoustic In PERiodic media) and focuses on development of solutions for broadband noise and vibration control using periodic structures. The VIPER Project involves five European universities: Università degli studi di Napoli (IT) "Federico II", Ecole Centrale de Lyon (FR), Université de Franche-Comté (FR), Katholieke Universiteit Leuven (BE) and University of Bristol (UK), where twelve PhDs are conducting their research covering a specific topic of the VIPER objectives. Many different type of research are conducting by the fellows, such as development of numerical techniques to improve the time computing of simulations, fluid structure interferences, effect of uncertainty and quasi-periodicity in the frequency response or topology design. The work presented here focus on the effect of inserts in honeycomb-like periodic structures for vibration control, represented by resonators added to the host structure in order to improve the dynamic performances in targeted frequency areas.

FIGURE 2 . 1 .

 21 FIGURE 2.1. Bending deformed shape for a positive (anticlastic) (a), zero (cylindrical) (b) and negative (synclastic) (c) Poisson's ration. (a) defines an hexagonal honeycomb and (c) a re-entrant honeycomb [47]

2. 2 ,

 2 respectively for a simple beam and a simple plate unit cell. As we can see, the curve and the dispersion surface present a repeated pattern along the axis of the wave number, coming from the definition of the Floquet-Bloch theorem. By doing the substitution µ = kd (respectively µ x = k x d x and µ y = k y d y ), where µ is called the reduced wave number, it has been demonstrated that the study of the wave propagation can be reduced to the interval [-π, π]. Leon Brillouin[START_REF] Brillouin | Wave Propagation in Periodic Structures[END_REF] did this observation some years after the Bloch's contribution, and the reduced interval considered has been called the First Brillouin Zone (FBZ). Brillouin also demonstrated that this zone can be reduced again depending on the symmetric aspect of the reciprocal space associated to the unit cell and the new interval created is known as the Irreducible First Brillouin Zone (IFBZ). Finally, the IFBZ can be reduced a last time for the study of band gaps just by observing the contour of the zone, since it has been claimed that the maximum admissible solutions are observed for value of µ x and µ y equal to 0 or π, even if this result has been contested with few counter examples in a recent study comparing the results with the choice of reciprocal spaces[START_REF] Maurin | Probability that a band-gap extremum is located on the irreducible brillouin-zone contour for the 17 different plane crystallographic lattices[END_REF]. In the case of the squared plate unit cell shown Fig.2.2b, the reciprocal space is also a square and presents 2 axes of symmetry on the x and y direction, making possible to restrict the study of the periodic structure on a triangle contour [O -A -B -O]. An example of step of decomposition from the dispersion surface to the contour of the IFBZ in the case of the simple plate is shown Fig. 2.3.

FIGURE 2 . 2 .FIGURE 2 . 3 .

 2223 FIGURE 2.2. Example of unit cells associated to dispersion curve for a 1D periodic structure (a) and dispersion surface for 2D periodic structure (b)

FIGURE 2 . 4 .

 24 FIGURE 2.4. Periodic system composed with the diatomic mass (m 1 ,m 1 ) (a) and dispersion curve associated (b)

FIGURE 2 . 5 .

 25 FIGURE 2.5. Representation of some stop bands using FRF obtained for 10,15 and 20 unit cell (a) and the dispersion curve with Floquet-Bloch theorem (b)

FIGURE 2 . 6 .

 26 FIGURE 2.6. Periodic structure containing the hosting structure (m h , k h ) and the resonator (m r , k r ) (a) condensed in effective mass m e f f (b)

FIGURE 3 . 1 .

 31 FIGURE 3.1. General representation of a unit cell for 2D wave propagation

3. 1 . 3

 13 WAVE FINITE ELEMENT METHODS FOR LINEAR PERIODIC ANALYSES (

L

  is then respected, and has the advantage to be numerically stable because exponentially growing terms are not involved in the computation. Now the terms [Ξ] n and Q n+1 l have to been successively calculated to link the end and the beginning of the periodic system. If we consider the finite structure with N unit cells and an imposed force at the extremity right of the structure, we get the following boundary condition: with q N+1 the degrees of freedom at the extremity right of the structure.

FIGURE 3 . 2 .

 32 FIGURE 3.2. Unit cell selection for hexagonal (a) and re-entrant (b) configuration

Fig. 3 .

 3 Fig. 3.4 shows dispersion curves for hexagonal unit cell for in-plane and out-of-plane modes while Fig. 3.5 gives the results for re-entrant configuration. The segments [O -A -B -C -O] represent the contour of the IFBZ for that specific configuration whose the values of µ x and µ y are described

FIGURE 3 . 4 .FIGURE 3 . 5 .

 3435 FIGURE 3.4. Dispersion curves for hexagonal unit cell : in-plane wave modes (a) and out-of-plane wave modes (b)

FIGURE 3 . 6 .FIGURE 3 . 7 .

 3637 FIGURE 3.6. Finite structure representation applying the Transfer Matrix Method

Fig. 3 .

 3 FIGURE 3.8. Unit cell including resonators : 2 resonators configuration (a) and 4 resonators configuration (b)

Fig. 3 .

 3 Fig. 3.9 gives the dispersion curves of both cases, zooming on the first modes to have access to the low frequency area. One can observe that the resonators open two band gaps around 2000 and 3000 Hz at the 2 first modes of the unit cell. By extracting the eigenvectors φ of the dispersion analysis, the mode shape of the unit cell is obtained by observing the real value of the eigenmode ℜ(φ) at the boundaries of the band gap. The mode shape at 2421 Hz and (µ x , µ y ) = (1.87, 0) is represented Fig. 3.10 and shows that the global displacement of the structure is done by the resonators, confirming that the band gap created is resonant.

FIGURE 3 . 9 .

 39 FIGURE 3.9. Dispersion curve comparison between no resonator and 2 resonators configuration and ρ 1 = ρ 2 = ρ aux for the re-entrant structure : no resonator (a) 2 resonators (b)

ρ 1 0

 1 and so only one resonator is actually present in the unit cell. Generally speaking, the more interesting result for a constant mass (i.e. larger band gaps at lowest frequencies) is obtained with the configuration n • 10 for 2 resonators, which corresponds exactly to 2 resonators with the same density. On the other hand, the less relevant result is related to the configuration n • 10 for 4 resonators, corresponding to having 4 resonators with the same density.

FIGURE 3 . 11 .

 311 FIGURE 3.11. Variation of density of the resonator analysis according to table 3.3 : 2 resonators configuration (a) and 4 resonators configuration (b)

FIGURE 3 .

 3 FIGURE 3.13. Geometry of the hosting unit cell

3 . 14 .FIGURE 3 . 14 .

 314314 FIGURE 3.14. Dispersion curve of the re-entrant structure without resonator using the MAC algorithm : Dispersion curve before (a) and after (b) MAC algorithm

  3.16, with L r = L 3 , b r = 1.33 mm and h r = 0.3 mm. As done for beam elements, the density of the resonator ρ r is variable while the Young Modulus and Poisson's Ratio have the same values as the hosting structure. One resonator in each butterfly pattern of the periodic structure (two resonators for each unit cell) is placed since it has been demonstrated in previous section that it was the most optimal way to put resonators for a constant mass. Although the resonator representation is simple here, one still need to pay attention of its orientation in the hosting cell. In that configuration, the stiffness of the resonators is less important in the (x, y) direction, encouraging the opening of out-of-plane band gaps in low frequencies. If the resonators were placed such as the sides related to (L r , b r ) are in the (y, z) plane instead of (x, y), that would have encourage more the opening of in-plane band gaps in low frequency, which is not the aim of the study here.

FIGURE 3 . 15 .FIGURE 3 . 16 .

 315316 FIGURE 3.15. Mode shape associated to some curves shown in 3.14b

CHAPTER 3 .

 3 LINEAR RESONATORS IN AUXETIC STRUCTURES also decreased by reducing the stiffness of the resonating, passing from a value of 407 Hz to 260

FIGURE 3 . 18 .

 318 FIGURE 3.18. Results obtained for ρ r = 14ρ and h r = h 1.4 h. Apparent band gap not obstructed visible on the dispersion curve (a) and the FRF using TMM (b).

FIGURE 3 . 19 .

 319 FIGURE 3.19. Dispersion curve for the unit cell including skins and resonators with resonant mode shapes associated

FIGURE 3 . 20 .

 320 FIGURE 3.20. Shape and geometric parameters of the resonators used in the volume element unit cell. Two types of configuration are used for the studies corresponding to either the red resonators only or the blue resonators only

FIGURE 3 . 21 .FIGURE 3 . 22 . 2 m

 3213222 FIGURE 3.21. Dispersion curve with resonant band gap in orange (a) and Resonating Mode shape (b) associated to the red resonator configuration with standard parameter values

FIGURE 3 . 24 .

 324 FIGURE 3.24. Constant mass stop band evolution performing parametric analysis of the resonator for each couple of value (h s , b s ) presented table 3.6

FIGURE 3 . 26 .

 326 FIGURE 3.26. Dispersion curve for the best value corresponding to the best result obtained in the parametric study and the associated resonator shape

FIGURE 4 . 1 .

 41 FIGURE 4.1. Initial conditions application for inverse method in a finite structure containing N unit cells

FIGURE 4 . 2 .

 42 FIGURE 4.2. 1D resonant spring mass system

FIGURE 4 . 3 .

 43 FIGURE 4.3. Dispersion curve for the infinite structure. Propagative part in red (ℜ(µ))and evanescent part in blue (ℑ(µ))

  being the eigenfrequency of the resonator. The dispersion curve is shown Fig. 4.3, with the propagative part in red and the evanescent part in blue for m R = m = 1 k g and k R = k = 1 N.m -1 . To validate numerically the result, a finite spring-mass chain of 20 unit cells has been simulated to obtain the Frequency Response Function (FRF) of the system. The boundary conditions in the numerical model are representative of a free-free state at the two extremities of the chain, while the mass located in one of the ends is excited by imposing a longitudinal propagating wave generated by the force F(t) = F 0 cos(ωt) with F 0 = 1N. An harmonic analysis is performed by varying the value of ω, and the amplitude of the displacement is captured on the mass at the other end of the structure. The results are shown in Fig. 4.4 and 4.3.

FIGURE 4 . 5 .

 45 FIGURE 4.5. 1D nonlinear resonant spring mass system

  FIGURE 4.6. Dispersion curve using the perturbation approach in the nonlinear case for A u = 1 and ε = 0.05

FIGURE 4 . 7 .

 47 FIGURE 4.7. Value of 1 h(ω 0 ) for the left branch of the dispersion curve (ω (1) 0 )

FIGURE 4 . 8 .

 48 FIGURE 4.8. Value of 1 h(ω 0 ) for the right branch of the dispersion curve (ω (2) 0 )

FIGURE 4 . 9 .

 49 FIGURE 4.9. Spatial representation of masses for two different values of µ

FIGURE 4 . 10 .

 410 FIGURE 4.10. Representation of the PML

Fig. 4 .

 4 Fig. 4.11. Fig. 4.15 shows the results using initial conditions shown in Fig. 4.14, with a value of A u 0 = 0.53 to get an average value A out = 0.5. The dashed line represents the value of amplitude equal 0.5.

  4.7). A relative frequency error is estimated as:(4.51) f err = | f num -f per | f perwith f num being the frequency obtained with the numerical simulations and f per the one calculated with the perturbation method. The value of 2 ω 0 ω 1 that represents the condition leading

FIGURE 4 . 11 .

 411 FIGURE 4.11. Flowchart of the algorithm to get the right amplitude A out

FIGURE 4 . 12 .FIGURE 4 . 13 .FIGURE 4 . 14 .

 412413414 FIGURE 4.12. Initial conditions of the system

FIGURE 4 . 15 .

 415 FIGURE 4.15. Time signal obtained for the mass x 0 after algorithm

FIGURE 4 . 16 .

 416 FIGURE 4.16. Value of 2 ω 0 ω 1 for the left branch of dispersion curve (perturbation method error)

FIGURE 4 . 17 . 3 FIGURE 4 . 18 .

 4173418 FIGURE 4.17. Model error estimation for the left branch of the dispersion curve

FIGURE 4 . 19 .

 419 FIGURE 4.19. Relative frequency difference between numerical and perturbation method

FIGURE 4 . 20 .

 420 FIGURE 4.20. Nonlinear dispersion curve (A u = 0.3, ε = 0.05)

FIGURE 4 . 22 .

 422 FIGURE 4.22. Zoom on the right branch of the dispersion curve (A u = 0.3, ε = 0.05)

FIGURE 4 . 23 .FIGURE 4 . 24 .

 423424 FIGURE 4.23. Right branch of the dispersion curve for A u = 1.5

FIGURE 4 . 25 .

 425 FIGURE 4.25. Comparison between Perturbation, HBM, inverse approach and harmonic analysis for |A u | = 0.5

FIGURE 4 . 26 .

 426 FIGURE 4.26. Comparison between Perturbation, HBM and inverse approach for |A u | = 1.2

4. 3 Application to finite element models 4 . 3 . 1 Lindstedt-Poincaré adapted 4 . 3 . 1 . 1

 34314311 Method explanation on a beam with 2 different sections

1 y

 1 (t) = a cos(ωt) + b cos(3ωt) with a and b two coefficients such as b < a, which is coherent with the equation of nonlinear geometrical forces after doing an asymptotic development on the displacement for low amplitude of vibrations, assuring that the physical sense of the model is respected. This procedure is repeated to all the forces and moments for each degrees of freedom, and the c 1 Fourier coefficient is calculated by integrating the nonlinear force for the last period calculated:

ω 1

 1 have been calculated (50 for each branch) and the model contains a total of 61 nodes (122 dofs), meaning that a total of 122 × 200 × 2 efforts have been calculated. The algorithm summarizing the procedure to get the corrected dispersion curve is given Fig. 4.30.

FIGURE 4 . 28 . 4 . 3 . 1 . 2

 4284312 FIGURE 4.28. Steps to get the nonlinear force related to F 1 y . Linear forces (b) are subtracted to linear forces (a) to get the nonlinear forces (c)

FIGURE 4 . 29 .FIGURE 4 . 31 .

 429431 FIGURE 4.29. Corrected dispersion curve obtained after Lindstedt-Poincaré Method application (blue) and linear dispersion curve (red) for A = 4.5 mm and ε = 1

FIGURE 4 . 30 .

 430 FIGURE 4.30. Algorithm used to get the corrected dispersion curve using Lindstedt-Poincaré method

  linear dispersion curve for in-plane wave propagation and the mode shapes close to resonances are shown in Fig.4.32. Looking at the mode shapes of the structure for some values of µ, one can observe at which value of wave number the resonator is the most likely to be impacted by the nonlinear analysis. The nonlinear analysis is then performed using an amplitude A = 0.5 mm in order to stay in the Lindstedt-Poincaré hypothesis and ε = 0.5. The corrected dispersion curve Fig.4.33a is showing an important impact at specific values of wave numbers, and more specifically for µ ∈ [1.8 , π] for the first branch and µ ∈ [0 , 2] for the second branch, corresponding to areas where the resonator got the higher value of amplitude compared to the rest of the structure. Moreover, values of |c 1 | representing the contribution of the nonlinear forces from each node are plotted Fig.4.32, showing that this contribution is much more important for the nodes related to the resonator compared to the rest of the structure. This observation is fundamental to effectively show that, overall, only the resonator has an effective nonlinear contribution. In that case, we can see that only the part of the branches related to the resonating modes are affected, leading to a partial dispersion curve shifting contrary to the example shown in Fig.4.29 where all the branches of the dispersion curves are shifted, eventually leading to less control for the location and width of the band gap.

FIGURE 4 . 32 .

 432 FIGURE 4.32. Linear dispersion curve of the resonant system (a) Mode shape for the first branch (b) and the second branch (c) close to the resonance

Table 4 . 1 :

 41 Material and geometrical properties of the hosting unit cell and the resonator according to Fig.4

FIGURE 4 . 33 .

 433 FIGURE 4.33. Comparison between linear and nonlinear dispersion curve (a) values of |c 1 | for µ = 1.8634 at the first branch of dispersion curve (b)

FIGURE 4 . 34 .

 434 FIGURE 4.34. Comparison between c 1 = f (A) and the function y = γA 3 , real part in (a) and imaginary part in (b)

FIGURE 4 . 35 .FIGURE 4 . 36 .

 435436 FIGURE 4.35. Initial conditions applied to the structure using ANSYS (a) and response obtained for the middle node (b)

FIGURE 4 . 37 .

 437 FIGURE 4.37. Dispersion curves using all the methods explained in the research applied to a simple beam. A = 3 mm for all the models and ε = 1 for the Lindstedt-Poincaré method

FIGURE 4 . 38 .

 438 FIGURE 4.38. Time signal obtained after releasing at 2 different nodes for the second branch of dispersion curve computed at µ = 0.6491

FIGURE 5 . 1 .

 51 FIGURE 5.1. Different type of meshing for the TMM unit cell (a) and the inverse WFEM (b)

FIGURE 5 . 2 .

 52 FIGURE 5.2. Nonlinear forces and moments related to out-of-plane wave propagation (u z , θ x , θ y ) of to the unit cell described 4.31 for ( f req, µ) = (27.5, 0.65) (green case) and ( f req, µ) = (222.3, 1.86) (blue case)

  Other ideas of resonators shape and location for shell element models are presented here but could not have been studied more in detail. The main purpose is to have an overview of different dispersion contour and mode shapes associated attaching resonators in different ways such as fixed-fixed, fixed-free on the core or on the skins. Material used for skins and cores of the models are the same as used for the model studied in chapter 3 and the resonators are made with steel.

  Fig. A.3 looks more familiar to the ones studied in the literature related to resonant structures, contrary to the contour shown in Fig. A.2, most likely coming from the boundary conditions of the resonator on the cores of the structure which are fixed-free on that case, like usual studied structures.

FIGURE A. 1 .

 1 FIGURE A.1. Dispersion contour and mode shapes associated for the sandwich structure (core + skins)

  2.5 gives an example of comparison between the dispersion curve and the FRF using 10, 15 and 20 unit cells for the case of a simple beam with 2 different sections. One can notice in the dispersion diagram that the Bragg band gaps appear when µ is a multiple of π, which is actually corresponding to the Bragg conditions d = n λ 2 , n ∈ N i.e. µ = nπ if we consider k = 2π λ and µ = kd with d the length of the unit cell and λ the wave length.

Table 3 .

 3

3.3 gives the geometrical features of the re-entrant unit cell and table 3.1 gives the values of material and geometrical properties associated. The material properties correspond to a composite made with flax and polypropylene which were initially planned to be used to manufacture the samples, while the dimensions (length L and angle θ) correspond to the dimensions of the mould. Dimensions and CHAPTER 3. LINEAR RESONATORS IN AUXETIC STRUCTURES Density, ρ [k g.m -3 ] 1040 L [mm] 6 Young Modulus, E [GPa] 8.1 h [mm] 0.3 Poisson Ratio, ν 0.2 b [mm] 0.3 θ [ • ] 30 1: Material and geometrical properties of the re-entrant core for the beam element model

Table 3 . 2

 32

: Values of µ x and µ y for the dispersion contour in table 3.2.

Table 3 .

 3 3: Values of ρ 1 and ρ 2 for each configuration of the band gap with mass variation. Fig.3.12 shows the results of the study considering the

	3.2. BEAM ELEMENTS

Table 3 .

 3 4: Material and geometrical properties of the re-entrant core for the shell element model

	3.3. SHELL ELEMENTS

Table 3 .

 3 5: Material and standard geometrical properties of the resonators for the volume element model and symmetries in (O,y,z) and (O,x,y) planes are done to get the entire unit cell.

  Zoom on the left branch of the dispersion curve (A u = 0.3, ε = 0.05)

	4.2. PERIODIC DISCRETE SYSTEMS INCLUDING RESONATORS
	3
	2.5
	2
	1.5
	1						Linear perturbation
							Linear inverse
							Non-linear perturbation
							Non-linear inverse
	0.07	0.08	0.09	0.1	0.11	0.12	0.13	0.14	0.15
				Frequency (Hz)
	0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 FIGURE 4.21. 0.23	0.24	0.25		0.26	0.27	0.28 Linear perturbation 0.29 Linear inverse Non-linear perturbation Non-linear inverse	0.3
				Frequency (Hz)

CHAPTER 3. LINEAR RESONATORS IN AUXETIC STRUCTURES

The dispersion curve including the resonators is computed taking in consideration ρ r = 10ρ and Fig. 3.17 summarizes the results obtained. The results have been separated in two different colors for global in-plane waves in red and out-of-plane waves in blue. To make this separation, an average between the value of φ i in-plane (u x and u y ) is compared with the out-of-plane values (u z ). Since some modes are coupled because of the complexity of the structure, this technique will work only in specific cases in particular in low frequency where the modes are still distinguishable.

Comparing with Fig. 3.14, one can see that a resonant band gap is opened (branches a and b in Fig. 3.17). However, the band gap is obstructed by the branch c representing a structural out-of-plane mode shape. This problem could be improved by adding mass to the resonator increasing the density which would have as an effect to switch the branches of the dispersion curve on lower frequency and therefore the branches related to the resonant stop band. However, this technique can add a lot of weight to the structure, which is not a very satisfying solution in practice if we try to limit the mass addition to the system. with

Assuming an harmonic regime is reached, the displacements u (0) j and q (0) j can be expressed in the following form:

2 e jiµ e iτ + Āu 2 e -jiµ e -iτ q (0) j =

A q 2 e jiµ e iτ + Āq 2 e -jiµ e -iτ

The second order derivatives are then:

(4.35) When adopting the same approach used during the linear analysis, the purpose is to isolate the corrected displacement q (1) j in the second line of Eq. 4.33 and inject it in the first line. In this way one can obtain all the nonlinear terms in the same equation. To do so, it is assumed that the displacement q (1) j can be expressed under the form:

Where B is the amplitude of q (1) j , whose the value does not have importance for the following calculations. Replacing Eq. 4.38 in the second line of Eq. 4.33 we obtain:

Injecting expression of q (1) j in the first line of Eq. 4.39 leads to:

APPLICATION TO FINITE ELEMENT MODELS

The main goal of the method is to find the nonlinear coefficients c 1 present in Eq. 4.15 since their analytical expressions can be very difficult to explicit as we could do for the discrete analysis as seen in the study of 1D nonlinear resonator. The explanation of the method will be done using the example of a beam presenting a variation of section. Only the bending waves will be taken geometries option in the software. To do so, a temporal displacement is imposed to all the nodes of the structure taking in consideration the mode shapes as the amplitude of displacement at a certain frequency. The vector of displacement in our specific case is written under the form:

where A is an arbitrary amplitude, n is the number of nodes of the unit cell,

In order to standardize the values of amplitude, one also have to normalize the mode shape before constructing the displacements u doing the operation :

mode calculated using the Floquet theorem is used as an initial condition for the full structure. As the method developed is recent, it still presents some parts that need to be improved, especially in the signal processing part which is far more harder than what we have seen for the discrete part when it is about structures modeled with FEM. For that reason, the example that will be given to compare with the Lindstedt-Poincaré method is a simple beam in order to confirm that the method theoretically works for a simple FEM case. However, a discussion about possible improvements for the method will be presented the sixth chapter.

The unit cell is identical to the one seen in Fig. 4.27 (except here h 1 = h 2 = 10 mm) and the degrees of freedom analysed remain the same (u y , θ). The first step consists again in calculating the eigenfrequencies and the wave vectors of the unit cell under linear hypothesis using the Inverse WFEM to get the vectors φ. Once this is done, a finite structure is built by repeating the unit cell a certain amount of time. In this example, we consider a finite structure built with 20 unit cells. Under these considerations and according to the method described in 4.1.2, the global mode shape of the finite structure is:

T the degree of freedom vector at the left of the structure and λ 0 = e jµ 0 . The vector of initial condition applied to all the nodes of the finite structure is then written such that:

Once the initial conditions are set up to the finite structure, a transient analysis is done during a time equivalent to a certain number of periods (5 in this example), and the results are obtained by looking at a specific node of the structure, preferably around the middle of the structure to avoid any border effects. An example of the initial condition imposed to the structure and the response in time obtained for the middle node for the second branch of the dispersion curve at µ = π 2 and amplitude A = 4.5 mm is given in Fig. 4.35. The time-amplitude signal is then analysed using FFT and Hilbert Transform to get the frequency. The two different approaches are used to see which one is likely to be used in that type of case. Unlike the example of discrete resonant system studied in previous section, the signal does not seem to lose so much amplitude in time, making justifiable the usage of the FFT in here. The FFT and the instantaneous frequency obtained using the Hilbert Transform are plotted in Fig. 4.36 and one can notice that the frequency obtained using the two methods are quite close (1592 Hz vs 1595 Hz). In comparison, the frequency obtained with the linear model is 1363 Hz exhibiting a significant change compared to the nonlinear model with a frequency increasing of 15 %. As usual, the operations has to be repeated for every couple (ω 0 , µ 0 ) to get the new values of the associated nonlinear frequency.

CONCLUSIONS AND FUTURE WORKS

The main objective of the thesis consisted in analysing the effect of resonators in periodic systems, and more specifically in auxetic structures providing interesting mechanic properties concerning the stiffness over the mass ratio, and also providing flexibility for the design of the structure due to its synclastic behaviour. Also, the improvement and development of tools for the case of non linear periodic structure has been investigated, and some example including resonant and non resonant systems have been used as first applications. All the studies performed have required a great phase of comprehension of wave propagation in periodic structures in order to create house made algorithms and tools covering a large number of methods used around the Floquet-Bloch analysis.

The first part of the thesis was focused on the linear resonant system analysis. The study of linear periodic structures including resonators has been divided in three categories, going from the simplest using beam elements to more realistic one using volume elements. Beam elements have been used as an introduction for the physical understanding of the resonator effect on the hosting structure and having first results with parametric analyses changing the material properties with fast computation. The study of the number of resonators to use has been performed, showing this specific case that using exactly two resonators per unit cell (corresponding to one resonator in each re-entrant pattern) with same material properties were leading to the lowest and widest stop bands in frequency. Shell elements were used to have results closest to manufacturable samples and having the opportunity to include the skins. One of the study carried out without the skins demonstrated the negative effect that can have a structural wave obstructing a resonant band gap with the observation of the dispersion curve, and an example of strategy used to get the full efficiency of the stop band without any mass addition have been shown. Also, the usage of skins shown the limitations of using shell element design, showing no clear band gaps with