Romain Lagrange Den-Service

Philippe Piteau

Xavier Delaune

Jose Antunes

FLUID-ELASTIC COEFFICIENTS IN SINGLE PHASE CROSS FLOW: DIMENSIONAL ANALYSIS, DIRECT AND INDIRECT EXPERIMENTAL METHODS

The importance of fluid-elastic forces in tube bundle vibrations can hardly be over-emphasized, in view of their damaging potential. In the last decades, advanced models for representing fluid-elastic coupling have therefore been developed by the community of the domain. Those models are nowadays embedded in the methodologies that are used on a regular basis by both steam generators providers and operators, in order to prevent the risk of a tube failure with adequate safety margins. From an R&D point of view however, the need still remains for more advanced models of fluid-elastic coupling, in order to fully decipher the physics underlying the observed phenomena. As a consequence, new experimental flowcoupling coefficients are also required to specifically feed and validate those more sophisticated models. Recent experiments performed at CEA-Saclay suggest that the fluid stiffness and damping coefficients depend on further dimensionless parameters beyond the reduced velocity.

In this work, the problem of data reduction is first revisited, in the light of dimensional analysis. For single-phase flows, it is underlined that the flow-coupling coefficients depend at least on two dimensionless parameters, namely the Reynolds number Re and the Stokes number Sk . Therefore, reducing the experimental data in terms of the compound dimensionless quantity r V Re Sk  necessarily leads to impoverish results, hence the data dispersion. In a second step, experimental data are presented using the dimensionless numbers Re and Sk . We report experiments, for a 3x5 square tube bundle subjected to water transverse flow. The bundle is rigid, except for the central tube which is mounted on a flexible suspension allowing for translation motions in the lift direction.

The evolutions of the flow-coupling coefficients with the flow velocity are determined using two different experimental procedures: [START_REF] Tanaka | Fluidelastic Vibration of Tube Arrays in Cross-Flow[END_REF] In the direct method, an harmonic motion of increasing frequency is imposed to the tube. [START_REF] Chen | Instability Mechanisms and Stability Criteria of a Group of Circular Cylinders Subjected to Cross-flow. Part 1: Theory[END_REF] In the indirect method, the coefficients are obtained from the modal response of the tube (frequency, damping). The coefficient identification was performed well beyond the system instability boundary, by using active control, allowing an exploration of a significant range of flow velocity.

For a given Sk , the results show that: (a) at low Re , the flow- coupling coefficients are close to zero; (b) at intermediate Re , the flow stabilizes the tube; (c) at high Re , the flow destabilizes the tube, leading to a damping-controlled instability at a critical Re . Reducing the data in terms of Re and Sk clarifies the various experimental "branches", which are mixed when using r V . The two identification techniques lead to reasonably compatible fluid-elastic coefficients.

NOMENCLATURE  2 f D f g C c DLV    Dimensionless damping coefficient 2 2 f K f g K c LV    Dimensionless stiffness coefficient 0 s f T s f C C c C C    Dimensionless total damping coefficient f C Fluid-added damping 0 f C Fluid-added damping in still fluid  0 f f f C C C   Fluid-

INTRODUCTION

The knowledge of the fluid force acting on a structure subject to a cross flow is a crucial information that must be accounted for when designing heat-exchanger tube bundles. The large vibrations resulting from a fluid-elastic instability may lead to some mechanical degradation of the concerned tube, which may affect the power plant operation and safety. This instability can be described as a self-excited feedback mechanism between the motion of the structure and the fluid forces. Since the pioneering work of Tanaka and Takahara [START_REF] Tanaka | Fluidelastic Vibration of Tube Arrays in Cross-Flow[END_REF], several authors [START_REF] Chen | Instability Mechanisms and Stability Criteria of a Group of Circular Cylinders Subjected to Cross-flow. Part 1: Theory[END_REF][START_REF] Chen | Instability Mechanisms and Stability Criteria of a Group of Circular Cylinders Subjected to Cross-flow. Part 2: Numerical Results and Discussion[END_REF][START_REF] Caillaud | The Measurement of Fluidelastic Effects at Low Reduced Velocities Using Piezoelectric Actuators[END_REF][START_REF] Caillaud | Measurement of Critical Velocities for Fluidelastic Instability Using Vibration Control[END_REF][START_REF] Caillaud | Active Vibration Control for the Measurement of Fluidelastic Effects[END_REF][START_REF] Mureithi | Dynamics of an Inline Tube Array Subjected to Steam-Water Cross Flow. Part 2: Unsteady Fluid Forces[END_REF][START_REF] Inada | An Experimental Study on the Fluidelastic Forces Acting on a Square Tube Bundle in Two-Phase Cross Flow[END_REF][START_REF] Sawadogo | Fluidelastic Instability Study in a Rotated Triangular Tube Array Subject to Two-Phase Cross-Flow. Part 1: Fluid Force Measurements and Time Delay Extraction[END_REF][START_REF] Sawadogo | Fluidelastic Instability Study in a Rotated Triangular Tube Array Subject to Two-Phase Cross-Flow. Part 2: Experimental Tests and Comparison with Theoretical Results[END_REF] measured the fluid-elastic force to feed the stability criterion models developed by Connors [START_REF] Connors | Fluidelastic Vibration of Tube Arrays Excited by Cross-Flow[END_REF], Blevins [START_REF] Blevins | Fluidelastic Whirling of a Tube Row[END_REF], Chen [START_REF] Chen | Instability Mechanisms and Stability Criteria of a Group of Circular Cylinders Subjected to Cross-flow. Part 1: Theory[END_REF][START_REF] Chen | Instability Mechanisms and Stability Criteria of a Group of Circular Cylinders Subjected to Cross-flow. Part 2: Numerical Results and Discussion[END_REF], Lever and Weaver [START_REF] Lever | A Theoretical Model for the Fluidelastic Instability in Heat Exchanger Tube Bundles[END_REF][START_REF] Lever | On the Stability of Heat Exchanger Tube Bundles. Part 1: Modified Theoretical Model[END_REF][START_REF] Lever | On the Stability of Heat Exchanger Tube Bundles. Part 2: Numerical Results and Comparison with Experiments[END_REF], Price and Païdoussis [START_REF] Price | Fluidelastic Instability of an Infinite Double Row of Circular Cylinders Subject to a Uniform Cross-Flow[END_REF][START_REF] Price | Improved Mathematical Model for the Stability of Cylinder Rows Subject to Cross-Flow[END_REF][START_REF] Païdoussis | The Mechanisms Underlying Flow-Induced Instabilities of Cylinder Arrays in Cross-Flow[END_REF], Granger et al. [START_REF] Granger | Motion-Dependent Excitation Mechanisms in a Square In-Line Tube Bundle Subject to Water Cross-Flow: An Experimental Modal Analysis[END_REF][START_REF] Granger | An Improvement to the Quasi-Steady Model with Application to Cross-Flow-Induced Vibration of Tube Arrays[END_REF] and Tanaka et al. [START_REF] Tanaka | Fluidelastic Vibration of Tube Arrays in Cross-Flow[END_REF].

Still, further experimental work is needed to accurately understand the effects of changing the tube mass, damping, frequency or diameter, as well as the bundle configuration, even for a single flexible tube within a rigid bundle subject to singlephase flow.

There are basically two techniques for obtaining the fluidelastic force. The direct method [START_REF] Tanaka | Fluidelastic Vibration of Tube Arrays in Cross-Flow[END_REF][START_REF] Mureithi | Dynamics of an Inline Tube Array Subjected to Steam-Water Cross Flow. Part 2: Unsteady Fluid Forces[END_REF][START_REF] Inada | An Experimental Study on the Fluidelastic Forces Acting on a Square Tube Bundle in Two-Phase Cross Flow[END_REF][START_REF] Sawadogo | Fluidelastic Instability Study in a Rotated Triangular Tube Array Subject to Two-Phase Cross-Flow. Part 1: Fluid Force Measurements and Time Delay Extraction[END_REF][START_REF] Sawadogo | Fluidelastic Instability Study in a Rotated Triangular Tube Array Subject to Two-Phase Cross-Flow. Part 2: Experimental Tests and Comparison with Theoretical Results[END_REF] is based on imposing a controlled oscillatory motion to a given tube within a rigid bundle and measuring the forces exerted by the flow as a function of the flow velocity and of the motion frequency. The fluid-elastic force is obtained from the transfer function between the tube displacement and the measured force from which the structural inertia term is subtracted. The indirect method [START_REF] Caillaud | The Measurement of Fluidelastic Effects at Low Reduced Velocities Using Piezoelectric Actuators[END_REF][START_REF] Caillaud | Measurement of Critical Velocities for Fluidelastic Instability Using Vibration Control[END_REF][START_REF] Caillaud | Active Vibration Control for the Measurement of Fluidelastic Effects[END_REF][START_REF] Piteau | Experiments and Computations of a Loosely Supported Tube in a Rigid Bundle Subjected to Single-Phase Flow[END_REF][START_REF] Piteau | Experimental identification of fluid-elastic coupling forces on a flexible tube within a rigid square bundle subjected to single-phase cross-flow[END_REF] applies to an instrumented tube, for which the fluid-elastic force is extracted from the changes in the modal frequency and damping of the vibrating tube. Both techniques have advantages and drawbacks, concerning the complexity of the setup, the sensitivity to external perturbations and the range of parameters that may be explored. Given all these constraints, an experimental rig allowing both measuring methods, as well as exploring beyond the instability threshold, has been developed at CEA, see [4-6, 22,23]. In this paper, we analyse these experimental results in the light of a dimensional analysis. In particular, we investigate the influence of the tube frequency and the water flow velocity on the fluid-elastic force.

DIMENSIONAL FLUID-ELASTIC COEFFICIENTS

We consider the vibration of a flexibly mounted tube (diameter D , length L ), part of a 3x5 square tube bundle (pitch P ), immerged in a viscous fluid of volume mass density f  , kinematic viscosity  and gap velocity g V , see Fig. 1. We note s M , s C and s K the mass, damping and rigidity coefficients of the first mode of vibration of the flexible tube in air. The modal displacement X in the lift direction is assumed to satisfy the equation:

      s f s f s f a M M X C C X K K X F         , (1)
where a F is an excitation fluid force considered as independent on the tube motion. The modal frequency of the flexible tube in still water is

    0 2 s s f F K M M    . ( 2
)
As most often assumed, we postulate that the added mass f M does not depend on flow velocity, but only on fluid density and bundle geometry. In other words,

f M is an unknown function f M H of   , , , f D L P  :   , , , f f M f M H D L P   . (3)
On the other hand, we consider that the fluid added damping f C and rigidity f K depend, at least, on the flexible tube frequency 0 F , the fluid material properties f  , , the gap velocity g V , as well as the bundle geometry

  0 , , , , , , f f C f g C H F V D L P    , (4)
  0 , , , , , ,

f f K f g K H F V D L P    . ( 5
)
Let 0 f C be the fluid added damping coefficient in still fluid, defined as the value of f C as 0

g V  . Then,  0 f f f C C C   is a
measure of the effect of the fluid velocity g V on the fluid added

damping coefficient f C , and is a function of   0 , , , , , , f g F V D L P   :     0 , , , , , , f f f g C C H F V D L P    . (6)
The relations of dependence ( 3) to ( 6) constitute a minimal model, based on experimental observations, bibliography reporting and physical intuition. More advanced models would also consider the effect of some other parameters, for e.g. the roughness of the tubes.

DIMENSIONLESS FLUID-ELASTIC COEFFICIENTS

The dimensional analysis is based on the Vaschy-Buckingham theorem. The theorem states that an equation involving n physical variables with k fundamental units (usually 3 k  in classical mechanics) can be reduced to an equation involving n k  dimensionless parameters. Thus, introducing a scale of length  , mass M and time  , the equations (3), ( 5) and ( 6) are physically meaningful if 3 , , ,

f f f M M D L P H M M              , ( 7
)
0 2 1 3 2 1 1 , , , , , , f f f g K K V F D L P H M M                        , (8)
  0 1 1 3 2 1 1 , , , , , , f f f g C C V F D L P H M M                        . (9)
As [START_REF] Mureithi | Dynamics of an Inline Tube Array Subjected to Steam-Water Cross Flow. Part 2: Unsteady Fluid Forces[END_REF] involves five dimensional quantities with two fundamental dimensions (length and mass), it can be reduced to a relation between three dimensionless quantities. Similarly, as ( 8) and ( 9) involve eight dimensional quantities with three fundamental dimensions (length, mass and time), they can be reduced to a relation between five dimensionless quantities. These dimensionless quantities are not unique and derive from a specific choice for the characteristic length  , mass M and

time  . Picking D   , 2 f M D L   and g D V   , the dimensionless equations rewrite   2 , f f m f M H l p D L   , ( 10
)   2 2 , , , f f K k r f g K c H V l Re p LV     , ( 11
)    2 , , , f f D c r f g C c H V l Re p DLV     , ( 12
) with l L D  , p P D  , g Re DV   , and 0 r g V V DF
 the tube aspect ratio, the pitch ratio, the Reynolds number and the reduced velocity, respectively. In what follows, we shall also make use of the Stokes number, obtained from the ratio between Re and r V :

2 0 r Sk Re V D F    .
To study the stability of the tube, we introduce the total damping coefficient

    0 T s f s f c C C C C   
, which also rewrites as

    0 0 0 0 1 1 2 s f f f D f g T s f s f s f C C C C c DLV c C C C C C C             . ( 13
) Introducing     2 0 2 s f f Sc M M D L     and      0 0 0 2 2 s f s f C C M M F         
as the Scruton number and the reduced damping parameter for a flexible tube in a still fluid, (13) simplifies to In the indirect method, the motion of the flexible tube is free and the measurement of its displacement and velocity is provided by a laser transducer (Keyence LK-G500). In order to investigate the variations of the fluid-elastic coefficients well beyond the stability threshold, a feedback control loop made of an electro-dynamical shaker was used, see [START_REF] Caillaud | The Measurement of Fluidelastic Effects at Low Reduced Velocities Using Piezoelectric Actuators[END_REF][START_REF] Caillaud | Measurement of Critical Velocities for Fluidelastic Instability Using Vibration Control[END_REF][START_REF] Caillaud | Active Vibration Control for the Measurement of Fluidelastic Effects[END_REF][START_REF] Piteau | Experimental identification of fluid-elastic coupling forces on a flexible tube within a rigid square bundle subjected to single-phase cross-flow[END_REF]. The coefficients D c and K c are obtained from the variations of the tube motion frequency and damping. This technique for extraction of the fluid-elastic coefficients from the modal parameters is simple and proved robust enough. Since three coefficients are to be extracted from only two modal parameters, the assumption of a constant added mass (e.g. a velocity-independent flow inertia coefficient) has to be somewhat arbitrarily enforced, which is a disadvantage of this identification method. Any possible changes in f M due to the fluid velocity will then be reflected in the coupling coefficient f K . The modal frequency of the moving tube is set by changing the thickness of the blades or by adding a suitable mass to the moving fixture. Several configuration tests, denoted L1, L2, L3, L3B, L4B and L4TI, see Table 1, have been performed, covering a large range of frequencies 0 13 ~39 Hz F  . The second modal frequency is much larger than 0 F , such that the tube dynamics is mainly a rigid translation.

  1 , , , , . 4 T D r T c r c V c H V l Re p Sc Sc    (14

EXPERIMENTAL SETUP AND MEASUREMENT METHODS

In the direct method, an harmonic motion of imposed frequency is directly applied to the tube, thanks to a PRODERA shaker, see Fig. 2. The fluid force acting on the flexible tube is only measured in the lift direction with a KISTLER sensor. The coefficients D c and K c are directly extracted from the measure of the fluid force. This method, first introduced by Tanaka in the 80's and mainly followed by Chen in the 90's has yield some interesting results. However, despite its apparent simplicity, it has been progressively abandoned due to its difficult experimental implementation. Three measurement campaigns, Table 1. Tested configurations. Indirect method.

EXPERIMENTAL RESULTS

In this section, we present the experimental measurements of D c , K c and T c , obtained from the indirect and direct methods. We focus on the evolution of these dimensionless coefficients with the Reynolds and the Stokes numbers. On the following figures, curves with identical colours correspond to configuration tests with similar Stokes numbers, see Table 1. In this range of R e , some kinetic energy is conveyed by the vibrating tube to the fluid, which in turn propagates this energy through the far domain. This corresponds to an energy loss for the tube leading to a damping of its vibrations. It is believed that the energy propagation through the far domain is enhanced by some fluid vortices whose existence still needs to be proven. At high Reynolds numbers, D c increases and the tube becomes unstable for some critical c Re , corresponding to 0 T c  . In this range of R e , some elastic and kinetic energy is conveyed by the fluid motion to the tube. This corresponds to an energy gain for the tube whose vibrations are amplified. In this range, we note that D c decreases with the Stokes number Sk , such that a tube with a high frequency is more stable than a tube with a low frequency.

The measurements of

In figures 4, 5 and 6, we compare the measurements of D c obtained from the indirect method and the direct method. For 6), significant differences are observed, especially at low Reynolds numbers. We attribute these differences to a bad signal to noise ratio due to some parasitic frequencies in the experimental setup as the forcing frequency 0 F of the direct method is increased. Also, at high forcing frequencies, the precise determination of the fluid-elastic force is complicated as most of the measured force has an inertia origin.

The measurements of K c obtained from the indirect method are shown in Fig. 7. At low Reynolds numbers, the data are scattered and difficult to analyze in the sense that no special trend is clearly observed. However, if the fluid velocity has a negligible effect on the dynamics of the structure, as it is believed from the analysis of D c , then one would expect that 0 K c  . At intermediate Reynolds numbers, K c becomes positive such that the relative rigidity of the tube diminishes. As already pointed out, this could be related to the existence of some fluid vortices interacting with the structure. At high Reynolds numbers, K c decreases such that the relative rigidity of the tube is enhanced. Still, K c being positive, the total rigidity s f K K  of the tube in a flowing fluid is smaller than its rigidity s K in a fluid at rest. Thus, in this range of R e , the frequency of the tube vibrations is smaller than 0 F .

The evolution of the total damping coefficient T c is shown in Fig. 8. The variations of this coefficient are directly related to those of D c through the equation [START_REF] Lever | On the Stability of Heat Exchanger Tube Bundles. Part 1: Modified Theoretical Model[END_REF]. Consequently, at low Reynolds numbers, 0 Indirect and direct methods. Indirect and direct methods.

CONCLUSIONS

In this work, we have studied the lift vibration of a flexible tube subject to a single phase cross flow. The flexible tube is located in the central position of a square rigid tube bundle. A dimensional analysis shows that the fluid-elastic coefficients are not a function of the reduced velocity only, but also depend on the Reynolds number. This observation is confirmed in our experiments by using two different methods of measurement. In the direct method, an harmonic motion of increasing frequency is imposed to the tube. In the indirect method, the coefficients are obtained from the changes in tube vibration frequency and damping. Both methods suggest the existence of three different dynamics for the flexible tube. At low Reynolds numbers, the fluid velocity has no effect on the stability of the tube. At moderate Reynolds numbers, the tube loses some energy and gets stabilized. At large Reynolds numbers, the tube gains some energy from the fluid and becomes unstable at a critical Reynolds. The experiments show that a tube with a high frequency is more stable than a tube with a low frequency. Finally, it shall be noted that some significant differences are observed in comparing the results of measurement from the two experimental methods as the Stokes number is increased and as the Reynolds numbers is decreased (i.e. low reduced velocities). We attribute these differences to a bad signal to noise ratio due to some parasitic frequencies in the experimental setup as the forcing frequency of the direct method is increased. Still, these first comparisons are very encouraging and should foster further developments of the direct method.

  analysis shows that the fluid-elastic coefficients D c and K c are not a function of r V only, but also depend on the Reynolds number. The reduced velocity has the disadvantage to encapsulate both 0 F and g V in a same dimensionless number, whereas these two parameters are independent. Consequently, to distinguish the effect of the frequency and the fluid flow on the variation of the fluid-elastic coefficient, it is preferable to use the couple of dimensionless numbers p . Similarly, the total damping coefficient T c depends on   , , , , Sc Sk l Re p .

Figure 1 .

 1 Figure 1. Bundle and flexibly mounted tube.

Figure 2 .

 2 Figure 2. Experimental setup.

  , L1TI 2018 and L3TI have been conducted, with imposed frequencies covering the range 0 13 ~39 Hz F  .

  D c obtained from the indirect method are shown in Fig. 3. Whatever the Stokes number, a clear general trend is observed. At low Reynolds numbers, 0 D c  , such that the fluid velocity has a negligible effect on the stability of the tube. At intermediate Reynolds numbers, D c becomes negative such that the tube is getting stabilized.

  In this range of Re the fluid velocity does not affect the stability of the structure. At intermediate Reynolds numbers, T c increases (i.e. D c decreases), meaning that the tube loses some energy and gets stabilized. At high Reynolds numbers, T c decreases (i.e. D c increases) meaning that the tube is gaining some energy from the fluid, leading to amplified vibrations. Eventually, at a critical Reynolds number c R e , the total damping coefficient vanishes and the tube becomes unstable. The figure 8 clearly shows that the stability threshold c R e increases with the Stokes number S k . The determination of the exact relation between c R e and S k is however difficult from the present experimental results as the configuration tests listed in Table 1 do not have exactly the same Scruton numbers. Still, we show in figure 9 that the representation of T c versus the reduced velocity r V R e S k  does not yield a collapse of the experimental data on a master curve, in particular close to the stability threshold 0 T c  . In other words, the flow-coupling coefficients do not only depend on r V , and contrarily to the prediction of Connors [11], c R e is probably not a linear function of S k .Finally, the dependence of the fluid-elastic coefficients with the Stokes number might be explained by several physical reasons. The most likely one is that the flow regime might be triggered by the tube frequency. The tube frequency might also affect the time-lag between the tube motion and the fluid forces. Whatever the physical explanation, the results presented suggest a strong dependence on both the flow Reynolds number and the Stokes number, beyond the classicaly assumed dependence on the reduced velocity.

Figure 3 .

 3 Figure 3. Evolution of the flow added damping coefficient D c with the Reynolds and the Stokes numbers. Indirect method.

Figure 4 .

 4 Figure 4. Evolution of the flow added damping coefficient D c with the Reynolds number.Indirect and direct methods.

Figure 5 .

 5 Figure 5. Evolution of the flow added damping coefficient D c with the Reynolds number.

Figure 6 .

 6 Figure 6. Evolution of the flow added damping coefficient D c with the Reynolds number.

Figure 7 .

 7 Figure 7. Evolution of the flow added rigidity coefficient K c with the Reynolds and the Stokes numbers. Indirect method.

Figure 8 .

 8 Figure 8. Evolution of the total damping coefficient T c

Figure 9 .

 9 Figure 9. Evolution of the total damping coefficient T c with the reduced velocity and the Stokes number.Indirect method

ACKNOWLEDGMENTS

The authors acknowledge financial support for this work, which was performed in the framework of a joint research program co-funded by FRAMATOME, EDF and CEA (France).


R Lagrange

F López Jiménez

D Terwagne

M Brojan

P M Reis
email: preis@mit.edu

From wrinkling to global buckling of a ring on a curved substrate

We present a combined analytical approach and numerical study on the stability of a ring bound to an annular elastic substrate, which contains a circular cavity. The system is loaded by depressurizing the inner cavity. The ring is modeled as an Euler-Bernoulli beam and its equilibrium equations are derived from the mechanical energy which takes into account both stretching and bending contributions. The curvature of the substrate is considered explicitly to model the work done by its reaction force on the ring. We distinguish two different instabilities: periodic wrinkling

of the ring or global buckling of the structure. Our model provides an expression for the critical pressure, as well as a phase diagram that rationalizes the transition between instability modes.

Towards assessing the role of curvature, we compare our results for the critical stress and the wrinkling wavelength to their planar counterparts. We show that the critical stress is insensitive to the curvature of the substrate, while the wavelength is only affected due to the permissible discrete values of the azimuthal wavenumber imposed by the geometry of the problem. Throughout, we contrast our analytical predictions against finite element simulations.

I. INTRODUCTION

Wrinkling is a stress-driven mechanical instability that occurs when a stiff and slender surface layer, bonded to a compliant substrate, is subject to compression. This universal instability phenomenon is found in numerous natural and technological/engineering examples, over a wide range of length scales, including: carbon nanotubes [START_REF] Lourie | [END_REF], pre-stretched elastomers used in flexible electronics applications [2], human skin [3], drying fruit [4], surface morphology of the brain [5] and mountain topographies generated due to tectonic stresses [START_REF] Price | Analysis of Geological Structures[END_REF][START_REF] Huddleston | [END_REF].

Over the past decade, there has been an upsurge of interest in the study of the mechanics of wrinkling, along with a change of paradigm in regarding surface instabilities as an opportunity for functionality, instead of a first step in the route to structural failure [8,9].

The first mechanical studies of wrinkling were motivated by the stability of sandwich panels [START_REF] Allen | Analysis and Design of Structural Sandwich Panels[END_REF], used in lightweight structural applications, in which the core acts as a soft substrate for the much stiffer skin. More recently, [START_REF] Bowden | [END_REF] showed how the wrinkling of a thin film on an elastomeric substrate can be used to produce complex self-organized patterns. Their seminal work has instigated the realization of wrinkling through several different actuation mechanisms, including thermal mismatch [12], tissue growth/atrophy [5,13,14], swelling by a liquid [15] or vapor solvent [16], and pneumatics [START_REF] Terwagne | Advanced Materials[END_REF]. The opportunities in applications opened by such a wide range of external stimuli have enabled the usage of wrinkling in photonics [START_REF] Kim | [END_REF], optics [19], self-assembly [20], microfluidics [21] and morphogenesis [22].

In order to provide a theoretical background to these recent developments, several authors have built on the pioneering work of [START_REF] Allen | Analysis and Design of Structural Sandwich Panels[END_REF], who first provided close form solutions for the critical stress and wavelength obtained when an initially straight beam, adhered to an infinite plane substrate, is placed under a state of uniaxial compression. [23] extended this work to consider the case of a plate adhered to a flat substrate under equi-biaxial compression and performed a nonlinear analysis of the Föppl-von Kármán equations [START_REF] Landau | Theory of Elasticity[END_REF][START_REF] Timoshenko | Theory of Elastic Stability[END_REF]. [START_REF] Huang | [END_REF] further refined these efforts by considering the effect of a finite substrate. Both studies showed the existence of multiple buckling modes associated with the same value of critical stress.

The stability of these modes under different loadings conditions has been addressed by [27][28][29], who produced a stability diagram covering the evolution from low to high values of overstress. However, experiments by [30] found disagreement at low values of overstress, suggesting that a finite intrinsic curvature of their experimental system, even if small, may play an important role in dictating pattern selection.

Early studies of wrinkling on curved substrates, as in the flat configuration, were also motivated by a structural problem; in this case, in the context of the stability of the outer shell of rockets [START_REF] Kachman | Test report on buckling of propellant cylinders under compressive loads[END_REF][32][33]. More recent studies that consider instabilities as a possible source of functionality have led to applications of curved configurations in adhesion [34], microfluidics [35], morphogenesis of microparticles [36], optics [37] and aerodynamic drag reduction [START_REF] Terwagne | Advanced Materials[END_REF].

Curvature also plays a relevant role in the growth of biological systems [38]. Despite these important emerging applications, the mechanics of wrinkling on curved substrates remains poorly understood, when compared to the planar counterpart.

Systematic Finite Element simulations of wrinkling in curved systems have been performed [3,4,39,40] that highlighted a complex pattern formation process. These numerical studies also suggested the possibility for curvature to affect the selected patterns and modify the relevant characteristic length scales, which calls for a robust theoretical backing. Analytical predictions are challenged by the difficulty of modeling the stiffness of the substrate, even in two-dimensional configurations. [START_REF] Cheng | Weight optimization of cylindrical shells with cellular cores[END_REF] and [30] used the stiffness provided by [START_REF] Allen | Analysis and Design of Structural Sandwich Panels[END_REF] for the flat case, such that their model therefore neglects the contribution of curvature on the response of the substrate. [4] used the prediction provided by [START_REF] Brush | Buckling of Bars, Plates and Shells[END_REF], which accounts for curvature but does not consider its influence on the wrinkling wavelength and their prediction does not converge to the classical planar case when the curvature tends to zero. As such, there is a need to quantify the effect of curvature on the stiffness of the substrate and its subsequent influence on wrinkling.

Here, to the best of our knowledge, we provide the first analytical work that accounts for both the curvature of a (2D) shell-substrate system, as well as the finite size of the substrate.

As an initial step, we focus our study on a curved film adhered to a cylindrical substrate, instead of dealing with non-zero Gaussian curvature geometries, which is left for a future study. We assume axial-symmetry to further simplify the system to the 2D problem of a ring on an annular substrate. Mechanical loading is applied by depressurizing a circular cavity inside the substrate, which places the system under a state of compression. This geometry is motivated by recent experiments [START_REF] Terwagne | Advanced Materials[END_REF] that demonstrated the usage of wrinkling on spherical samples for switchable and tunable aerodynamic drag reduction. In our simplified 2D system, we solve the elasticity problem for the substrate and derive a close form expression for its stiffness, which is then used in the stability analysis of the ring to quantify the buckling patterns.

The paper is organized as follows: In § II, we introduce our system along with its material and geometrical parameters. We also describe the possible instability modes, and present a simplified phase diagram, with the aim of providing physical insight on the problem. In § III, we then introduce the kinematics of the ring attached to the substrate and determine the stiffness of the substrate. We proceed by defining a strain energy that includes both bending and stretching of the ring, as well as the effect of the substrate. Energy minimization yields the equilibrium equations of the problem. An asymptotic expansion is then used to calculate the principal solution and the bifurcation at the onset of instability. In § IV, we describe the finite element simulations that we have performed for this same system.

The results of our investigation are presented in § V. Throughout, we directly compare the analytical predictions to the numerical simulations. We start with the fundamental solution and the critical conditions that lead to instability. We then construct a phase diagram which rationalizes the dependence of the instability modes on the governing parameters.

The results for our system are then quantitatively compared to those for wrinkling of a film on a planar substrate, highlighting the effect of curvature. Finally, § VI summarizes our findings and provides perspectives for potential extensions of our work in future studies.

II. DEFINITION OF THE PROBLEM

We study the stability of a thin elastic ring, bound to an equally curved 2D substrate that contains an inner cavity, a schematic diagram of which is presented in Fig. 1(a). The system is initially at equilibrium, with identical pressures inside and outside of the sample.

Motivated by recent experiments on spherical specimens [START_REF] Terwagne | Advanced Materials[END_REF], the system is then loaded by applying a depressurization, P, to the inner cavity. The thickness of the ring is H, its Young's modulus E F and its Poisson's ratio ν F . We refer to

E F = E F /(1 -ν F
2 ) as the reduced Young's modulus of the film. The substrate is made of a linearly elastic material with Young's modulus E S , Poisson's ratio ν S and reduced Young's modulus

E S = E S /(1 -ν S 2 ).
The thickness of the substrate is R -R 0 , where R 0 is the radius of the inner cavity.

For convenience, we now introduce new rescaled quantities to reduce the number of parameters of the problem. As such, we use R and E F to normalize lengths and pressures and define

h = H R , β = R 0 R , ξ = E F E S , p = P E F , (1)
as the dimensionless thickness, cavity size, stiffness ratio, and pressure, respectively.

The principal solution corresponds to an axisymmetric deformation that leads to a decrease of both the inner and outer radii. As the depressurization increases, the onset of instability is reached. In Fig. 1(b-e) we show representative results obtained from Finite Element Modeling (FEM), of the three possible instability configurations of the ring-substrate system, for different values of the dimensionless ring thickness, h, and ratio of stiffness, ξ. For these results, all the other mechanical properties were kept constant: the Poisson's ratios of the film and substrate are ν F = ν S = 0.5, and the dimensionless size of the cavity is β = 0.2. In the simulations, this is achieved by fixing R = 100 units of length and E S = 1 units of pressure, while changing the values of H, R 0 and E F accordingly. More details of our numerical simulations are provided in §IV.

The first mode, a representative example of which is shown in Figure 1(b) for h = 10 -2 and ξ = 10 3 , corresponds to periodic wrinkling of the film with a well defined wavelength.

The displacements are localized in the region close to the film. The second instability mode, for example at h = 10 -2 and ξ = 10 6 in Figure 1(c), corresponds to a global buckling of the structure, where both the ring and the cavity deform into an ellipse such that the wavelength is λ = πR.

In addition to these two instability modes (wrinkling and global buckling), we have also numerically observed an instability on the inner surface of the cavity. However, this third mode does not affect the ring and is only found for low values of cavity size and stiffness ratio (ξ = 10 2 and h = 10 -2 in Fig. 1c-d). This instability was first discussed by [START_REF] Biot | Mechanics of Incremental Deformation[END_REF] and

we therefore refer to it as the Biot mode; it is local in nature and only depends on the compressive strain at the inner surface. This type of instability mode has been recently studied in the case of elastomeric materials with voids by [START_REF] Michel | [END_REF] and [45]. Understanding the specifics of this Biot mode is however outside the scope of our work and we shall not take it into account in our analytical model and systematic numerical investigation.

A schematic phase diagram of our system is provided in Fig. 2. For low values of h and ξ, the ring wrinkles with a short wavelength. As either h or ξ are increased, the wavelength also increases. Once these parameters reach a critical value, represented by the dashed line in Fig. 2, the instability transitions from wrinkling to global buckling. In what follows, we focus on rationalizing how the wavelength of the wrinkling mode, and the threshold value for the transition to global buckling, evolve with the elastic and geometrical parameters of the system.

III. ANALYTICAL MODEL

The ring is treated as an Euler-Bernoulli beam. The effect of the substrate is modeled as a restoring force that acts on the ring and is determined by solving the elasticity problem of the substrate with adequate boundary conditions. Minimization of the potential energy provides the equilibrium equations of the problem, which are solved using an asymptotic expansion that yields the principal and bifurcated solutions.

κ/R = |dT/ds|.
Here, κ is dimensionless and can be written in terms of v and u as

κ = 1 + (-1 + 2u + 2v) v + v 2 -v -u 2 + (-u + v ) u + 1 2 v 2 -u 2 + h.o.t., (2)
where the prime notation represents derivation with respect to θ and high order terms (h.o.t.) are neglected under the assumption of small displacements and moderate rotations.

We now define the elongation of the ring as e = ds/ds 0 to express the stretching deformation, η = (e 2 -1)/2, in terms of v and u as

η = 1 2 u 2 + v 2 + (u -v ) 2 + (1 + v) u + v, (3)
so that the hoop stress in the film is σ 0 = E F η.

Following Euler-Bernoulli beam theory [START_REF] Timoshenko | Theory of Elastic Stability[END_REF], the total energy of deformation E of the ring is the sum of a stretching energy E S and a bending energy E B ,

E = E S + E B = 2πR 0 Eds 0 , (4)
with

E S = 2πR 0 E F H 2 η 2 ds 0 , (5a)
E B = 2πR 0 E F H 3 24R 2 (κ -1) 2 ds 0 , ( 5b
)
and E is the energy of deformation per unit length of the initial configuration of the ring.

Assuming that the reaction force of the substrate derives from a potential 2πR 0 W ds 0 , the equilibrium states of the ring are the solutions of

δE S + δE B - 2πR 0 δW ds 0 = 0, ( 6
)
where δA is the variation of quantity A, for an arbitrary displacement field R (δv, δu), which is 2π periodic. The computation of the variations in Eq. ( 6) leads to the Euler-Lagrange equations for the equilibrium of the ring,

∂E ∂v - ∂E ∂v + ∂E ∂v - ∂δW ∂δv = 0, ( 7a
) ∂E ∂u - ∂E ∂u + ∂E ∂u - ∂δW ∂δu = 0, (7b)
along with static boundary conditions that are naturally satisfied due to the 2π periodicity condition on the displacements v and u. All derivative terms in Eq. (7b) are explicitly reported in Appendix A.

B. Asymptotic expansion and reactive force of the substrate

We seek a solution of Eq. ( 7) as an expansion of the form

v = v 0 + εA sin (mθ) + O ε 2 , (8a) u = εB cos (mθ) + O ε 2 , (8b)
where (v 0 , 0) corresponds to a radial pre-buckling deformation and ε(A sin (mθ) , B cos (mθ))

represents an instability of azimuthal wavenumber m. From the requirement of 2π periodic functions v and u, m has to be an integer. Also, we consider m > 1 since m = 1 corresponds to a solid body translation of the ring. The instability displacement field has amplitudes εA and εB, where ε is a small parameter.

Before substituting Eq. ( 8) into Eq. ( 7) and solving at each order in ε, we first need to determine the work, δW , done by the reaction force, F = -σe ⊥ + τ e , that the substrate exerts on the ring. The normal stress σ and the tangential stress τ at the interface between the ring and the substrate are computed by solving the corresponding two-dimensional elasticity problem using an Airy function, with pressure P at r = R 0 and the displacement field given by Eq. ( 8), at r = R.

It is worth to note that in our computation, we enforce the continuity of the displacement and stress fields at the interface ring/substrate, similarly to [46]. This is a main difference with previous studies of wrinkling surfaces, e.g. in [33], who assumed zero shear stress and continuity of the normal stress. Such an approach yields a stretching energy in the ring much larger than the bending energy, in contradiction to what is expected in the wrinkling of a thin film [27]. Our explicit solution for the boundary value problem of the substrate is reported in B.

In short, we find that the stresses σ and τ at the interface are

σ E F = k 0 v 0 + kεA sin(mθ) -γp, (9a)
τ E F = µεB sin(mθ), (9b)
with the following governing parameters

k 0 = (1 -ν S ) (1 -β 2 ) 1 -2ν S + β 2 1 ξ , (10a) k
= (1 -ν S ) 2 S A + S B B A 1 ξ , (10b) µ
= (1 -ν S ) 2 S B A B + T B 1 ξ , (10c)
γ = 2β 2 (1 -ν S ) 1 -2ν S + β 2 . ( 10d
)
Here, k 0 , k, and µ are the dimensionless pre-wrinkling, wrinkling and shear stiffnesses of the substrate, respectively. The coefficient γ quantifies the effect of pressure, p, on the ring through the term -γp in Eq. (9a). The value of γ is smaller than one, reflecting the fact that the transmitted pressure decreases within the substrate. There are however two exceptions for which γ = 1: the limit where there is no substrate, β → 1, and the limit of a perfectly incompressible substrate, ν S → 0.5, in which the volumetric pressure remains constant throughout the substrate. The quantities S A , S B and T B used in Eq. ( 10) are functions of both the cavity size β and the azimuthal wavenumber m, and their full expressions are reported in B.

From the solution of the linear elasticity problem for the substrate (see B), the reaction force, F, has a constant direction along e r . However, based on physical intuition, one would expect the restoring force to change its direction as the ring deforms. Accounting for this scenario in full would have required solving the elasticity problem for the substrate, with nonlinear kinematics. Here, for simplicity, we assume that the magnitude of the restoring force is given by the solution of the linear elasticity problem, while the direction of the force is given by the vectors e ⊥ and e . Depending on which configuration of the ring is used (initial or deformed), two cases need to be considered to define the normal and tangential vectors e ⊥ and e , respectively. If we define these vectors from the initial configuration, then e ⊥ , e = (e r , e θ ), so that the reaction force F keeps a constant direction while the ring deforms. If we use the deformed configuration, then F is modeled as a follower force whose direction changes with the ring deformation. In this case, e ⊥ , e = (N, T), where N is the inward normal vector orthogonal to the ring center-line, in its deformed state. In order to take both of these options into account in the same model, we define a parameter, χ, which can take the values χ = 1 or χ = 0, when either the undeformed or deformed configurations are used. Thus, we express the reaction force and its elementary work per unit length of the ring as

F = -[(σe r + τ e θ ) χ + (σN + τ T) (χ -1)] , (11a)
δW = F • R (δve r + δue θ ). ( 11b
)
where • is the Euclidean dot product.

1. The principal solution: order 0 in ε

We now proceed to obtain the principal solution and the critical instability modes by substituing Eq. ( 8) into Eq. ( 7) and solving the resulting equations for each order of ε. At order 0, Eq. ( 3) for the stretching deformation η yields η 0 = v 0 = σ 0 /E F and the linear approximation of the Euler-Lagrange Eq. (7a)

v 0 = γp h + h 3 /12 + k 0 = σ 0 E F , (12)
relates the dimensionless pressure, p, and the dimensionless hoop stress, σ 0 /E F , in the ring.

Note that Eq. (7b) is automatically satisfied at order 0 in ε. In the absence of a substrate, the ratio σ 0 /P reads lim β→1

σ 0 P = 1 h + O (h) , (13)
and in the limit of the ring that is much stiffer than the substrate we obtain lim ξ→∞

σ 0 P = γ h + O (h) , (14)
so that σ 0 P < lim ξ→∞ σ 0 P ≤ lim β→1 σ 0 P . For ν S = 0.5, we have γ = 1 through the definition in Eq. (10d) and both limits in Eqs. ( 13) and ( 14) are equal. Thus, the hoop stress in a very stiff ring lying on a soft incompressible substrate is the same as the hoop stress in a ring with no substrate, which serves as a verification of the rationale thus far.

The instability: order 1 in ε

At order 1 in ε, the Euler-Lagrange Eqs. ( 7) write

(a 1 p + a 1 ) A + b 1 p + b 1 B = 0, ( 15a
) (a 2 p + a 2 ) A + b 2 p + b 2 B = 0, ( 15b
)
where a i , a i , b i and b i are functions of the azimuthal wavenumber m, the stiffness ratio ξ, the cavity size β and the substrate Poisson's ratio ν S , reported in full in C. The linear system in Eq. ( 15) has a nontrivial solution when its determinant vanishes, which for a given value of m occurs at the dimensionless pressure

p m = P m E F = -a 1 b 2 + b 1 a 2 a 1 b 2 + a 1 b 2 -b 1 a 2 -b 1 a 2 . (16)
The critical azimuthal wavenumber and the dimensionless critical pressure can now be obtained from Eq. ( 16) by minimizing over all possible values of m

[m c , p c ] = min m=2,3,... (p m ) . (17)
For the general case, this minimization must be performed numerically. However, in the absence of the substrate, we obtain that m c = 2 and Eq. ( 16) simplifies to lim

β→1 (p c ) = lim β→1 (p 2 ) = - 1 4 -χ h 3 - 4 (4 -χ) 2 h 5 + O h 7 . ( 18
)
The term of order h 3 corresponds to the classical dimensionless critical pressure for a ring with no substrate [START_REF] Timoshenko | Theory of Elastic Stability[END_REF][47][48][START_REF] Combescure | Calcul analytique de la déformation d'un anneau avec défauts[END_REF]. The term in h 5 is a correction also reported by [50]. In the limit of a ring much stiffer than the substrate, we also have

m c = 2, with lim ξ→∞ (p c ) = lim ξ→∞ (p 2 ) = - 1 2 (-1 + 2ν S -β 2 ) h 3 β 2 (-1 + ν S ) (4 -χ) - 2 (-1 + 2ν S -β 2 ) h 5 β 2 (-1 + ν S ) (4 -χ) 2 + O h 7 , (19)
so that lim ξ→∞ (p c ) ≤ lim β→1 (p c ) < 0. Again, for ν S = 0.5, we observe that the limits equal, corroborating the physical intuition that as the stiffness ratio ξ increases, the effect of the substrate becomes negligible.

C. Critical stress, wavelength and comparison to wrinkling on a planar substrate

To determine the effect of curvature on the instability, we compare the critical hoop stress, σ c , in the ring, given by Eq. ( 12) with p = p c , as well as the wavelength λ c = 2πR/m c of the wrinkling mode, against their counterparts for an initially planar film on a plane substrate of infinite thickness [23] σ

P lane = E F 4 3 E S * E F 2/3 = E F 4 4(1 -ν S ) 2 3 -4ν S 3 ξ 2/3 , (20a)
λ P lane = 2πH 1 3 E F E S * 1/3 = 2πH 1 3 3 -4ν S 4(1 -ν S ) 2 ξ 1/3 , (20b)
where

E S * = 4E S (1 -ν S ) 2 /(3 -4ν S )
is the effective stiffness of the substrate. Together, Eqs. ( 20) and ( 12) yield

σ c σ P lane = 4γp c h + h 3 /12 + k 0 12(1 -ν S ) 2 3 -4ν S 1 ξ -2/3 , (21a)
λ c λ P lane = 1 hm c 3 -4ν S 12(1 -ν S ) 2 ξ -1/3 , (21b)
and, in the limit of absence of the substrate, Eq. (21a) reduces to lim

β→1 σ c σ P lane = 1 σ P lane lim β→1 (σ c ) = 1 σ P lane lim β→1 (σ 2 ) = 4 4 -χ 12(1 -ν S ) 2 (3 -4ν S ) ξ -2/3 h 2 +O h 4 , (22)
where ν S is the Poisson's ratio of the substrate in the planar case.

In § V D, we shall make use of Eqs. (21) and Eq. ( 22) to further quantify the sensitivity of the stress and wavelength to the curvature of the substrate.

IV. NUMERICAL SIMULATIONS

In § V, we will contrast the predictions from the above analysis with the results of a series of finite element simulations performed using the commercial package Abaqus, with the BUCKLE analysis, which provides the buckling load and the corresponding eigenmodes.

As in the analytical study, the cylindrical structure was modeled as an annulus of a soft substrate with a stiffer thin film adhered to its exterior, under plane strain conditions.

The bonding between the substrate and the film is assumed to be perfect, such that both share nodes. Negative pressure is applied on the interior surface, to model the effect of the pressure differential between the inner surface of the annulus (the cavity) and the exterior of the system. The buckling analysis provides the value of the critical pressure, p c , as well as the corresponding critical mode. Rigid body motions are removed by constraining the displacement of two points on the film.

Both substrate and film are modeled as incompressible linearly elastic materials, ν S = ν F = 0.5. These results were compared to additional simulations using a Neo-Hookean model but no difference was observed given the low values of strain involved. The substrate was modeled using quadrilateral plane strain elements. Due to the incompressibility, the corresponding hybrid element, CPE4H, was used. The film was modeled using B21 beam elements. In order to account for the effect of plane strain, the stiffness of the beam is defined as

E F / (1 -ν 2 F ).
All of the results presented were obtained using 1000 elements in the circumferential direction, and 150 (R -R 0 ) elements in the radial direction. The mesh size was validated with a convergence analysis. By way of example, differences less than 0.5% in critical pressure were obtained when comparing results for a mesh with twice the elements in each direction, even in the cases of wrinkling with the shorter wavelengths. The deviations between the two meshes were, however, larger (∼ 5%) for the case when the critical buckling mode is a Biot instability, due to the infinite number of wavelengths associated to the same buckling mode. However, as stated in § II, this mode shall not be studied in detail as we focus on the wrinkling and global buckling loads.

To test the validity of using beam theory to describe the ring in our problem, we also performed numerical simulations with 2D solid elements (CPE4H, the same used for the substrate) and found excellent agreement with the simulations using beam elements. However, in order to achieve such agreement, the mesh needs to be greatly refined, resulting in a significant increase in computational cost. Attempts to use a mesh size similar to that of our previous simulations showed clear disagreement. Given the excellent agreement between the two versions -either using CPE4H (with a fine mesh) or B21 for the ring -as well as the significantly lower computational cost of the B21 elements, we have decided to use beam elements for the ring in our analysis.

V. RESULTS

Having introduced our analytical and numerical methods, we now present the results of a systematic exploration of the mechanical response of our system for different geometric and material parameters. Throughout, we provide a direct comparison between analytical results and numerical simulations, finding good agreement. A few instances of discrepancy will also be discussed.

For the geometric parameters, we have varied the cavity size, β = R 0 /R, and the dimensionless thickness, h = H/R. Three representative values were chosen for β (depicted in the insets of Fig. 3a-c): a small cavity, β = 0.2; a cavity with size half of the external radius, β = 0.5; and a large cavity, β = 0.8. Moreover, h was varied in the range 10 -3 to 10 -1 .

This parameter has two different physical interpretations. On one hand, for a substrate with given curvature, i.e. fixed R, increasing h is equivalent to increasing the thickness of the ring. On the other hand, for a ring of given thickness H, the value of h decreases with the curvature, 1/R.

For the material properties, we have considered values for the stiffness ratio between the film and the substrate, ξ = ĒF / ĒS , spanning over five orders of magnitude, from 10 2 to 10 7 .

The substrate is taken to be incompressible, ν S = 0.5, since most of the relevant experiments that have motivated our study [START_REF] Bowden | [END_REF]12,[START_REF] Terwagne | Advanced Materials[END_REF]51] use nearly incompressible elastomeric substrates. Moreover, it is important to note that, even if the analytical model has been derived assuming a general value of ν S , it is expected to be less accurate for increasing deviations from incompressibility. This was analyzed by [30], who showed that for the wrinkling of plates deviations when ν S = 0.3 are of just a few percent.

In the presentation of our results, we first consider the effect of h and ξ on the hoop-stress and on the critical pressure. Then, we rationalize the transition from wrinkling to global buckling, shown in Fig. 1. We finally compare the critical stress and wavelength of the wrinkling mode to their planar substrate counterparts and discuss the effect of curvature.

A. Hoop stress prior to wrinkling

In Fig. 3, we plot the hoop stress in the ring, normalized by the pressure, σ 0 /P , as a function of the dimensionless thickness, h = H/R. We find that the hoop stress decreases In short, our model exhibits limitations if the ring and the substrate have comparable stiffness, or when there is a Biot mode (small cavity). Apart from this extreme combination of parameters, rarely observed in experimental configurations, the model performs successfully. (c) β = 0.8. Theory (lines) is given by Eq. ( 17). The Poisson's ratio of the substrate is ν S = 0.5.

The legend (bottom right) is common to all three plots.

C. Phase diagram

We proceed by focusing on the transition from wrinkling to global buckling, towards first constructing a phase diagram in the (ξ, h) parameter space and then quantifying the dependence of the boundary, h * (ξ), between the two modes on the size of the cavity, β.

Given that we do not have a closed form expression for h * , we use a numerical method which tracks any jump from m c = 2 to m * c > 2 in Eq. ( 17), when h is decreased.

In Fig. 5(a), we plot h * as a function of ξ, for β = 0.5. The boundary between modes is consistent with a power-law, h * ∼ ξ -3 , which divides the phase diagram into a wrinkling domain (h < h * ) and a global buckling domain (h > h * , shaded region). Again, there is good agreement between numerical and analytical solutions.

In Fig. 5(b), we extend the phase diagram to cavity sizes from β = 0.2 to β = 0.8. We predict that the power-law with exponent -3, mentioned above for β = 0.5, is still valid as β increases, even if there are some deviations towards the higher values. To quantify the appropriateness of the -3 power-law, we fit a curve of the form ξ = f 2 h f 1 to the analytically calculated boundaries. In Fig. 5(c) we plot the fitting parameters f 1 and f 2 as a function of β. We find that the exponent is f 1 ≈ -3 for β < 0.7, and then decreases rapidly to a value of f 1 ≈ -4.25 for β = 0.8. The decrease of the prefactor f 2 is more pronounced and reflects the fact that h * decreases as β increases.

The evolution of m * c (i.e. azimuthal wavenumber of the wrinkling mode for h = h * ) as a function of β and ξ is presented in Fig. 5(d). We observe that m * c increases with β, while it is nearly insensitive to ξ. For example, for β = 0.5, we find that m * c = 9 for ξ > 10 2 and m * c = 10 for ξ = 10 2 . Predictions for ξ > 10 2 are in relatively good agreement with finite element simulations, which have shown that m * c = 10. For ξ = 10 2 , FEM simulations yield m * c = 15, pointing out, once again, the limitation of our analytical approach in the limit when the stiffness of the ring and the substrate become comparable.

D. Comparison of wrinkling in our curved system with that on infinite planar substrate Thus far, we have shown that the ring may wrinkle or buckle globally, depending on the curvature and stiffness of the substrate. We now focus on the wrinkling mode of the ring, with the aim of comparing the critical stress and wavelength of our curved system to their counterparts for a infinite planar substrate. We shall center our discussion of this comparison for β = 0.5. In D, we report the results for β = 0.2 and β = 0.8, which are qualitatively similar.

In Fig. 6(a), we plot σ c /σ P lane given by Eq. (21a), as a function of h. For h < h * (wrinkling domain), we find that σ c /σ P lane ≈ 1, in agreement with the FEM simulations.

To further quantify how close to unity is this ratio, in Fig. 6(b) we plot the dimensionless critical normal force σ c h/E F as a function of the planar result, σ P lane h/E F . We obtain a line with unit slope, indicating that the substrate curvature has no significative effect on the critical stress for wrinkling. From this observation, and using Eq. ( 12) with σ 0 = σ P lane , we write the following approximation for the dimensionless critical pressure,

|p c | = h + k 0 γ σ P lane E F + O h 3 , (23)
with k 0 and γ given by Eq. [START_REF] Allen | Analysis and Design of Structural Sandwich Panels[END_REF]. Finally, back to Fig. 6(a), for h > h * (global buckling domain), we find that σ c /σ P lane first decreases with h, then reaches a local minimum and eventually increases as a power-law with slope 2, towards the asymptotic limit given by Eq. (22). The evolution of σ c /σ P lane in the global buckling domain is also reproduced well by FEM simulations. We now investigate the effect of curvature on the wavelength of the instability mode.

In Fig. 7(a), we plot λ c /λ P lane , given by Eq. (21b), as a function of h. Focusing on the wrinkling domain, we find that λ c /λ P lane ≈ 1, in agreement with FEM simulations. In

VI. CONCLUSION

We have considered the two-dimensional problem of a ring bound to an elastic substrate which contains a cavity that is depressurized. An energy formulation was used to derive the Euler-Lagrange equations that govern the equilibrium of the ring, and solved them via an asymptotic expansion. As an improvement to previous results in the literature, our analytical approach accounts for the effect of curvature in modeling the reaction force of the substrate. These analytical results were compared with numerical simulations.

We first studied the principal solution, obtaining an expression for the hoop stress in the ring as a function of the applied pressure. We then performed a stability analysis of the problem to determine the critical pressure, P c , and the corresponding instability mode.

Depending on the dimensionless thickness and stiffness ratio (h and ξ) we have identified two different regimes: local wrinkling of the ring, and global buckling of the structure. The boundary between both regions of instability was described via a detailed phase diagram, which quantifies the value of h and ξ at which the transition between instabilities occurs and takes into account the cavity size, β. Our results can be used as a design guideline to target a desired mode. Finally, we have shown that the critical stress for wrinkling and the resulting wavelength do not depend significantly on the curvature of the substrate. However, curvature imposes a discretization of the wrinkling wavelength due to the periodic closing conditions of the ring.

Our study focused on a 2D curved system which exhibits instability modes analogous to the cylindrical pattern found for uniaxial compression of a film on a flat infinite substrate.

Considering more complex loading conditions (e.g. also introducing axial loading) or shells with non-zero Gaussian curvature, should lead to more complex patterns that deserve to be investigated further. Having validated the FEM analysis, as well as carefully considering the elastic response of the curved substrate, extending the study to these other scenarios and addressing the issue of pattern selection is an exciting avenue for future research. Applying the boundary conditions Eq. (B2) to Eq. (B3) yields a linear system for A i and B i , with solution

B 1 = 1 2 (β 2 + 1 -2ν S ) E S v 0 (1 -ν S ) -P β 2 , (B4a)
B 2 = -R 2 β 2 β 2 + 1 -2ν S E S v 0 (1 -ν S ) + P (1 -2ν S ) , (B4b)
k 0 = K 0 R E F = (1 -ν S ) (1 -β 2 ) 1 -2ν S + β 2 1 ξ , (B8c)
k = KR E F = (1 -ν S ) 2 S A + S B B A 1 ξ , (B8d)
µ = MR E F = (1 -ν S ) 2 S B A B + T B 1 ξ , (B8e)
as indicated in Eq. [START_REF] Allen | Analysis and Design of Structural Sandwich Panels[END_REF].

We note that in the case of an inextensible wrinkling mode, A/B = m, the stiffness K simplifies to K given by The cavity size is β = 0.8 and the Poisson's ratio of the substrate is ν S = 0.5. The legend is common to both plots.

K = E S 1 R (

FIG. 1 .

 1 FIG. 1. (a) Schematic diagram of our system: a ring is bound to a curved substrate which contains a circular cavity. The system is loaded by applying a pressure differential between the inside of the cavity and the outside of the ring. (b-d) Representative examples of the three possible instability modes of a ring on a curved substrate which contains a cavity that is depressurized. (b) Wrinkling mode (h = 10 -2 and ξ = 10 3 ), (c) global buckling mode (h = 10 -2 and ξ = 10 6 ) and (d) Biot mode (h = 10 -2 and ξ = 10 2 ). (e) is a zoom in of (d) that exhibits the deformation of the surface of the inner cavity in the Biot mode. The adjacent colorbar applies to pictures b-e) and refers to the maximum principal component of the strain tensor of the mode, which has been normalized by the maximum value of each configuration.

FIG. 2 .

 2 FIG.2. Schematic phase diagram of the instability modes of the system, obtained for a critical value of the pressure differential: i) wrinkling of the ring or ii) global buckling of the structure (shaded regions). The primary parameters that govern this transition are the stiffness and thickness ratios: ξ = E F /E S and h = H/R, respectively. The depicted examples from FEM simulations are for an incompressible film and substrate, ν F = ν S = 0.5. They were obtained for R = 100 units of length and E S = 1 units of pressure, while varying H and E F . The colorbar refers to the maximum principal component of the strain tensor of the mode, which has been normalized by the maximum value of each configuration.

FIG. 3 .

 3 FIG. 3. Pre-instability hoop stress, σ 0 , in the ring, normalized by the pressure, P , as a function of the dimensionless ring thickness, h = H/R. Cavity sizes are: (a) β = 0.2, (b) β = 0.5 and (c) β = 0.8. The Poisson's ratio of the substrate is ν S = 0.5. Analytical predictions are given by Eq. (12) as solid lines and FEM results are plotted as data points. The legend (bottom right) is common to all three plots.

FIG. 4 .

 4 FIG. 4. Dimensionless critical pressure |p c | = |P c |/E F versus the dimensionless ring thickness h = H/R, for several stiffness ratios ξ = E F /E S and cavity sizes β: (a) β = 0.2, (b) β = 0.5 and

FIG. 5 .

 5 FIG. 5. (a) Phase diagram in the (ξ, h) parameter space, showing the transition boundary, h * , from wrinkling to global buckling, for β = 0.5. (b) Phase diagram for β = 0.2 to β = 0.8. Solid lines are analytical predictions and dotted lines correspond to the numerical fit ξ = f 2 h f 1 . (c) Fitting coefficients, f 1 and f 2 , as functions of β. (d) Azimuthal wavenumber of the wrinkling mode for h = h * , as a function of β. The Poisson's ratio of the substrate is ν S = 0.5.

FIG. 6 .

 6 FIG. 6. (a) Critical stress σ c normalized by its planar substrate counterpart, σ P lane , as a function of the dimensionless ring thickness h = H/R. (b) Dimensionless critical normal force versus planar substrate counterpart. Analytical predictions (solid lines) are given by Eq. (21a) and FEM results are shown as data points. Insets are sketches of the wrinkling and global buckling modes. The cavity size is β = 0.5 and the Poisson's ratio of the substrate is ν S = 0.5. The legend is common to both plots.

Fig. 7 (

 7 Fig. 7(b), once again, we quantify how close to unity this ratio is by plotting the critical wavelength λ c as a function of the planar result, λ P lane . To first approximation, noting the large dynamic range (at least two orders of magnitude) in both axes of the plot, we find a line with unit slope that passes through the origin, suggesting that the substrate curvature has no significative effect on the wavelength of the wrinkling mode. The only deviations arise from the discrete nature of the wavenumber, since the geometry of the ring enforces 2π periodic wrinkling modes. To highlight this phenomenon, in Fig.7(c), we plot the critical azimuthal wavenumber m c , given by Eq. (17), as a function of h and we superimpose the planar substrate result (dashed line), m P lane = 2πR/λ P lane . We find that m c is a decreasing stair function of h. The deviation in wavelength between the curved and planar cases scales as (|λ c -λ P lane |)/λ P lane ∼ 1/m, which is maximum for m * c . As shown previously in Fig.

FIG. 7 .

 7 FIG. 7. (a) Critical wavelength, λ c , normalized by its planar substrate counterpart, λ P lane , as a function of the dimensionless ring thickness h = H/R. Analytical predictions (solid lines) are given by Eq. (21b). (b) Zoom of the wrinkling domain by plotting λ c versus λ P lane . (c) Azimuthal wavenumber m c of the wrinkling mode as a function of h. Analytical predictions (solid lines) are given by Eq. (17). Dotted lines show m c for the planar substrate case. Insets on (a) and (b) are sketches of the wrinkling and global buckling modes. Inset on (c) is a zoom in showing the difference in m c for the curved and planar substrates. The cavity size is β = 0.5 and the Poisson's ratio of the substrate is ν S = 0.5. FEM results shown as data points. The legend (bottom right) is common to all three plots.

the 2 A 1 1 -A 1 1 

 21111 Fulbright Program. P.M.R. is grateful to financial support from the National Science Foundation, CMMI-1351449 (CAREER) and Saint-Gobain.σ rθ (r, θ) =   -A 1 m (m + 1) r m + A 2 m (m -1) r -m -A 3 m (m -1) r m-2 + A 4 m (m + 1) r -m-(α -m -1) r m+1 + A 2 (α + m -1) r -m+1 -A 3 mr m-1 + A 4 mr -m-(α + m + 1) r m+1 + A 2 (α -m + 1) r -m+1-A 3 mr m-1 -A 4 mr -m-

A 1 =

 1 (m -1) β 2 m+2 + (α -m + 1) β 2 m + β 2 2 (-m 2 + 1) β 2 m+2 + m 2 β 2 m+4 + β 4 m+2 α + (m 2 -1 + α 2 ) β 2 m + αβ 2 GR -m A + (m -1) β 2 m+2 + (-m -α + 1) β 2 m -β 2 2 (-m 2 + 1) β 2 m+2 + m 2 β 2 m+4 + β 4 m+2 α + (m 2 -1 + α 2 ) β 2 m + αβ 2 GR -m B,(B4c)A 2 = -(m + 1) β 2 m+2 + β 4 m+2 + β 2 m (α + m + 1) 2 (-m 2 + 1) β 2 m+2 + m 2 β 2 m+4 + β 4 m+2 α + (m 2 -1 + α 2 ) β 2 m + αβ 2 GR m A + (m + 1) β 2 m+2 + β 4 m+2 + β 2 m (α -m -1) 2 (-m 2 + 1) β 2 m+2 + m 2 β 2 m+4 + β 4 m+2 α + (m 2 -1 + α 2 ) β 2 m + αβ 2 GR m B,(B4d)A 3 = -β 2 ((α -m + 1) (m + 1) β 2 m + m 2 β 2 m+2 + α + 1 + m) (2 (-m 2 + 1) β 2 m+2 + m 2 β 2 m+4 + β 4 m+2 α + (m 2 -1 + α 2 ) β 2 m + αβ 2 ) m GR -m+2 A β 2 (m 2 β 2 m+2 -1 + (-m -α + 1) (m + 1) β 2 m -m + α) (2 (-m 2 + 1) β 2 m+2 + m 2 β 2 m+4 + β 4 m+2 α + (m 2 -1 + α 2 ) β 2 m + αβ 2 ) m GR -m+2 B,The dimensionless stress at the interface and the dimensionless stiffness parameters are obtained from Eq. (B5), by dividing it with E F , σ E F = k 0 v 0 + kεA sin(mθ) -γp,

1 -ν S ) 2 S A + S B 1 m = E S 1 mR 2 ( 1 -K = E S 1 mR 2 ( 1 -

 122121 ν S ) 2 (β -2 m -β 2 m + 2m (1 -β -2 )) (m 2 -1) (β -1 -β) 2 (m 2 -1) + (β -m -β m ) 2 (3 -4ν S ) + 3-4ν Sβ a substrate with no cavity, leads to the limiting case limβ→0 ν S ) 2 (m 2 -1) 3 -4ν S ,(B10)and, for an infinite plane substrate, yields limR→∞ lim β→0 K = K P lane = E S 4(1 -ν S ) 2 3 -4ν S

FIG. 9 .

 9 FIG. 9. (a) Critical stress σ c normalized by its plane substrate counterpart σ P lane , as a function of the dimensionless ring thickness h = H/R. Analytical prediction is given by Eq. (21a). (b) Critical wavelength λ c normalized by its plane substrate counterpart λ P lane , as a function of h. Analytical prediction is given by Eq. (21b). Insets are sketches of the wrinkling and global buckling modes.

ACKNOWLEDGEMENTS D.T. thanks the Belgian American Education Foundation (B.A.E.F.), the Fulbright Program and the Wallonie-Bruxelles International Excellence Grant WBI World. M.B. thanks

monotonically with the dimensionless thickness, and increases with both the stiffness ratio and the cavity size. When the ring is much stiffer than the substrate (i.e. ξ → ∞), the hoop stress scales as σ 0 /P ∼ (H/R) -1 , with a prefactor given by Eq. (14). For a given thickness of the ring and a given cavity size, the hoop stress in the ring decreases as the curvature of the substrate increases. This observation is consistent with the classic result for the hoop stress, σ 0 = P R/H, for a depressurized thin-walled cylindrical pressure vessel. This result can be recovered from Eq. ( 12) by taking γ = 1, k 0 = 0 and performing a Taylor expansion in h, about 0. In this regime, there is a clear asymptote as ξ → ∞, given by Eq. (19), which corresponds to the buckling of a ring with no substrate.

Depending on the model used for the reaction force of the substrate (constant direction force, χ = 1, or follower force, χ = 0, see § III B), Eq. ( 11) yields two different analytical predictions, shown in Fig. 4 as a solid line for χ = 0, and a dashed line for χ = 1. We note that the agreement between FEM simulations and analytical predictions is superior for χ = 0, in particular for the global mode. In other words, it is better to model the reaction force of the substrate as a pressure field which remains normal to the ring center-line C as it deforms, than as a force with constant direction. Thus, from now on, all analytical results will be only presented for χ = 0.

Despite the overall good agreement between the analytical and numerical results, there are noticeable discrepancies in the global modes for low values of the stiffness ratio. The reason is that, for ξ = 10 2 , the strains in the substrate can become significant, such that the assumptions for our linear theory are no longer valid (see B). A particularly extreme case is the appearance of Biot modes for β = 0.2 and ξ = 10 2 . These modes, although possible, only occur in a small region of our parameter space and, as mentioned above, are beyond The terms in the Euler-Lagrange equation for the equilibrium of the film, Eq. (7b), are

Appendix B: Response of the substrate

In this appendix we consider the boundary value problem of the substrate, subjected to the ring displacement R (v (θ) , u (θ)), v and u given by Eq. ( 8), at the interface r = R between the ring and the substrate and to the pressure P at r = R 0 = βR. The substrate is assumed to be in a state of plane strain. We note α = 3 -4ν S and introduce the shear

The components of the stress in the substrate are represented by σ rr , σ rθ , and the displacement field is (U r , U θ ). The 2D problem of elasticity is solved by finding an Airy function of the form [52],

where A i and B i are unknown constants determined by the boundary conditions

which are assumed to apply at r = R 0 and r = R. The first two equations stand for the continuity of the stress at the boundary of the cavity, whereas the last two stand for the continuity of the displacement at the interface between the substrate and the ring.

The stress and displacement fields resulting from the Airy function Eq. (B1) are [START_REF] Barber | Elasticity: 2nd Edition[END_REF] σ rr (r, θ)

Substituting A i and B i into Eq. (B3) yields the stress at the interface r = R σ rr (R, θ) = K 0 Rv 0 + KRεA sin (mθ) -γP, (B5a)

where

and

Appendix C: Terms of the linear stability analysis

The terms a i , a i , b i and b i that appear in Eq. ( 16) are

where γ, k 0 , S A , S B and T B were given in B.


Claude Guérin

Richard K Albert

Jeremy R Beitler

Luciano Gattinoni

Samir Jaber

John J Marini

Laveena Munshi

Laurent Papazian

Antonio Pesenti

Antoine Vieillard-Baron

Prone position in ARDS patients: why, when, how and for whom

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