Designing necks and wrinkles in inflated auxetic membranes

This article presents the potentiality of inflatable, functionally-graded auxetic membranes to produce wrinkles and necks. We obtain elastic instabilities at desired locations in axisymmetric membranes and with prescribed patterns in square membranes. First, we use an analytical approach to obtain a series of universal results providing insights into the formation of wrinkles and necks in inflated, axisymmetric membranes. For example, we prove analytically that necks and wrinkles may never overlap in pressurized, axially symmetric membranes. Second, we implement the relaxed strain energy of tension field theory into a Finite Element solver (COMSOL). By tuning spatial inhomogeneities of the material moduli, we corroborate our universal results, describe the onset of wrinkling in an averaged way, and also generate non-trivial instabilities at desired locations. This study on membranes with morphing or corrugation on demand has potential applications in Braille reading and haptics.


Introduction
Artificially designed materials with negative Poisson's ratio thicken when stretched, in contrast to classical elastic materials, which become thinner in the directions lateral to an applied loading direction.These so-called auxetic materials have been studied in several fields, including ferromagnetics (Popereka and Balagurov, 1970), crystal elasticity (Milstein and Huang, 1979), foam structures (Lakes, 1987), microporous materials (Evans and Caddock, 1989), and composites (Milton, 1992).
Auxetic materials can be fabricated using additive manufacturing techniques, such as powder bed fusion (King et al., 2015;Sun et al., 2017) or subtractive manufacturing techniques, such as laser cutting for thin structures (Mueller et al., 2013;Bhullar et al., 2017).With the rapid advancement in additive manufacturing techniques, even complex auxetic structures can be produced at large scales in a short time.With a careful design of voids or holes at the micro-scale, auxetic properties can be obtained at the continuum level (Lakes, 1993;Bertoldi et al., 2017).
These structures have a wide range of applications, from soft robotics, with compliant actuators (Lazarus and Reis, 2015) or compliant grippers (Kaur and Kim, 2019), to biomedical applications, with auxetic stents (Dolla et al., 2006), Figure 1: Examples of applications and behaviors for auxetic materials.(a) Hip implants can be designed to include regions with a negative Poisson ratio to minimize retraction from the bone under biomechanical loading (Kolken et al., 2018).(b) Auxetic tubular lattice stents exhibit increased ductility compared to conventional diamond tubular lattices (Jiang et al., 2022).(c) Auxetic thin membranes under tension wrinkle in the neighborhood of clamps, in contrast to conventional membranes, which wrinkle in their central region (Bonfanti and Bhaskar, 2019).(d) Auxetic patches glued onto balloon membranes can undergo larger deformations than conventional patches, a property which can help with the healing of puncture wounds to the bladder for example (Chansoria et al., 2022).
Due to their thin-walled, lightweight, and impressive tensile properties, nonlinear elastic membranes play a prominent role in many fields such as automobile, aerospace, civil and biomedical engineering (Beatty, 1987;Evans, 2009;Fu et al., 2016).However, thin membranes lose their mechanical stability under in-plane compressive stresses due to negligible bending rigidity, leading to interesting behaviors (Timoshenko, 1983;Roddeman et al., 1987;Cerda et al., 2002).Inflatable membranes under large elastic deformations experience various kinds of bifurcation phenomena such as limit-point (snap-through), wrinkling, and necking instabilities.Some of the most common and significant applications of inflatable auxetic materials can be seen in stent deployment for angioplasty (Amin et al., 2015), in morphing structures (Sun et al., 2014), and for smart biomaterials (Bhullar et al., 2015).
Over the past several decades, there has been extensive research on finite inflation in incompressible circular membranes.Early experimental works by Flint and Naunton (1937), Treloar (1944), and Rivlin and Saunders (1951) present a detailed analysis of the deformation shape and strain distributions over the surface of the inflated incompressible isotropic balloons.Adkins and Rivlin (1952) studied theoretically pressurized circular and spherical thin shells using incompressible neo-Hookean and Mooney-Rivlin strain energy functions, to obtain the relationship between inflating pressure, extension ratio, and radius of curvature near the pole.The book by Green and Adkins (1960) covers the mathematical theory of finite strain deformations in nonlinear membranes extensively for both isotropic and anisotropic materials.More recent experimental works include that by Machado et al. (2012), who proposed a method to determine curvatures and membrane stresses using an axisymmetric bulge test for isotropic circular membranes and a three-dimensional digital image correlation technique.
Finite element analysis (FEA) is very useful in understanding the behavior of complex geometries and loading conditions.Many researchers have developed numerical methods based on FEA to study finite inflation in nonlinear membranes using different material models and for different geometries.For example, see the works on axisymmetric membranes by Oden and Sato (1967), Wriggers and Taylor (1990), Gruttmann and Taylor (1992), Jiang andHaddow (1995), Rumpel et al. (2005), Eriksson andNordmark (2012), andSelvadurai (2022), on rectangular membranes using finite difference iterative scheme by Yang and Lu (1973), and on square membranes using FEA by Adler and Mikulas (2001), Lee and Youn (2006), Barsotti andLigarò (2014), andChen et al. (2014).
Analytical solutions associated with the finite inflation of nonlinear membranes are scarce in the literature due to strong material and geometrical nonlinearities.However, analytical solutions play an important role in providing simplified and direct solutions to predict the mechanical behavior of inflated membranes.With an assumption of linear elastic constitutive behavior and spherical deformation shape for a pressurized incompressible isotropic circular membrane, Fichter (1997) provided closed-form analytical solutions for small deformations, later extended by Coelho et al. (2014) for finite strains.Relaxing the constraint of linear elastic material behavior but still assuming the spherical deformation shape of the membrane, Yuan et al. (2021), Foster (1967), and Yang et al. (2021) derived analytical solutions for pre-stretched circular membranes under inflation using different incompressible hyperelastic material models.Dropping the hypothesis of spherical deformation shape and considering compressibility, pressure-deflection formulae for inflated isotropic circular membranes with compressible Mooney-Rivlin material model are provided by Pelliciari et al. (2022) without pre-stretch, by Sirotti et al. (2023) with pre-stretch, and by Pelliciari and Tarantino (2022) for anisotropic pressurized graphene membranes without prestretch.
Many of the above-mentioned works dealt with limit-point instability in inflatable circular membranes, but little attention has been drawn in the past to the study of wrinkling and necking instabilities for functionally-graded inflatable auxetic circular membranes under mechanical loads.
Several theories exist in the literature for studying wrinkles in nonlinear elastic membranes.For example, the theory of incremental deformations (Haughton and Ogden, 1978a,b), the Föppl-von Kármán theory of plates (Dym et al., 1973;Puntel et al., 2011), the reduced-order finite element membrane theory (Damil and Potier-Ferry, 2010), and the numerical bifurcation-continuation analysis (Healey et al., 2013).Although these advanced and refined models provide comprehensive details on the wavelength and amplitude of wrinkles, they are computationally expensive to implement.
In this study, we focus on determining the average deformation in the wrinkled region along with finding the location and orientation of wrinkles, but not their amplitude and wavelength.Therefore, we employ tension field theory, originally proposed by Wagner (1931) and Reissner (1938).Tension field theory has the advantage of being computationally more viable and mathematically elegant.According to tension field theory, membranes are assumed to have zero out-of-plane bending stiffness and cannot sustain in-plane compressive stresses.In the 1980s, in order to account for compressive stresses, Pipkin (1986b) extended the theory by introducing the concept of "relaxed strain energy function", see also Pipkin (1986a), Steigmann and Pipkin (1989b), Steigmann (1990), andPipkin (1994).Relaxed models have been used to study wrinkles in anisotropic membranes (Pipkin, 1995), electroactive elastomeric membranes (De Tommasi et al., 2011;De Tommasi et al., 2012;Greaney et al., 2019;Khurana et al., 2022a,b), pressurized magnetoelastic circular membranes (Saxena et al., 2019), and inflatable isotropic membranes under uniform pressure load (Dhavale et al., 2014;Barsotti, 2015;Pamulaparthi Venkata and Saxena, 2019).
Alongside wrinkling, necking is another interesting bifurcation phenomenon observed in nonlinear membranes, although less explored in hyperelastic isotropic circular membranes.Chaudhuri and DasGupta (2014) found negative Gaussian curvature and circumferential wrinkling at the fixed rim of the hyperelastic isotropic circular membrane under inflation.Necking in pressurized elasto-plastic spherical membranes has been investigated by Needleman (1976) and in pressurized elasto-plastic circular membranes by Chater and Neale (1983).However, a study on necking and multilayered bubbling phenomenon in pressurized functionally-graded hyperelastic isotropic auxetic circular membranes is still missing in the literature.
In this work, as a proof-of-a-concept, we study the effect of varying material properties such as the Young modulus and the Poisson ratio on limit-point instabilities (snap-through), necking, and wrinkling in pressurized isotropic auxetic membranes (circular and square geometries).We use the Blatz-Ko strain energy function to model the membrane's hyperelastic mechanical behavior.We study the effect of pre-stretches on wrinkling instabilities in auxetic circular membranes.For the onset of limit-point instabilities, we investigate the effect of pre-stretches and material parameters for circular geometries.Finally, we show how to obtain wrinkling patterns in specific areas of the membrane by tuning the spatial distribution of Young's modulus and Poisson's ratio, both for square and circular membranes.
The paper is organized as follows.In Section 2, we introduce the problems of interest and the constitutive equations, and we derive the kinematics of the deformation.We also briefly summarize the main features of the relaxed strain energy functional based on tension field theory to derive the membrane stresses.Finally, we write down the equilibrium conditions along with the applied boundary conditions.In Section 3, with the help of membrane theory, we derive the equations linking curvatures and principal stresses, which are necessary for wrinkling and necking in the circular membrane.We establish several universal insights, valid for all hyperelastic isotropic membranes.These include the results that regions of necking and wrinkling cannot overlap, and that necking and wrinkling cannot occur in the center of a pre-stretched circular membrane.
In Section 4, we compare the results of our finite element simulations in COMSOL (COMSOL Multiphysics® Version 6.0, 2021) with the universal predictions of Section 3. Conclusions and limitations of the current work, along with possible directions for future works are detailed in Section 5.These include some preliminary results on square membranes which can be inflated to exhibit a desired wrinkling pattern.These results could be used in applications involving haptics.

Membrane deformations, energy and stress
In this work, we describe the behavior of inflated elastic membranes that are rotationally symmetric about an axis.In the so-called "membrane approximation", the membrane thickness is small in comparison to its diameter and bending effects are neglected.In this approximation, the 3D deformation is deduced from the 2D deformation of the membrane mid-surface.
We consider a circular membrane with radius  in and we use a cylindrical coordinate system to represent the kinematics of the deformation.We identify the position of a point on the mid-surface of the membrane in its undeformed configuration with P 1 (, Φ, 0).As the membrane is radially stretched axisymmetrically, the point P 1 moves to position P 2 (  0 , Φ, 0 ) in the membrane with radius  fin .Upon axisymmetric inflation of the pre-stretched membrane (with a fixed circumference) under a uniform pressure  from the side  = 0 − , the membrane bulges out of the plane towards  = 0 + , and the point P 2 is displaced to position P 3 ((), Φ, ()).Here, () and () represent the radial and transverse deflections of the point P 1 in the final configuration, respectively, as shown in Figure 2. Deformed profile of the circular membrane under radial pre-stretch with the radius  fin .Note that the thickness of the membrane is increased with the in-plane stretching due to auxeticity.(c): Deformed profile of the inflated circular membrane subjected to uniform pressure  and fixed on its circumferential boundary.Point P 1 in the undeformed configuration is displaced to the positions P 2 and P 3 with pre-stretch and inflation, respectively.A point in the final configuration has coordinates (, Φ, ), where  is radial position and  is the vertical deflection.
The in-plane stretches of the 2D membrane are calculated as where the index 1 refers to meridians and 2 to parallels in the current configuration.
In the membrane approximation, it is assumed the deformation is normal preserving.Therefore the gradient of deformation admits the diagonal representation  = diag( 1 ,  2 ,  3 ) in a local basis, where the eigenvalue  3 is relative to the membrane normal direction.By introducing the volumetric variation coefficient  = det , we may therefore write With F written this way, we may now compute the 2D membrane energy from any 3D energy.Specifically here, to model the compressible behavior of the membrane, we use the three-dimensional Blatz-Ko strain energy density (Blatz and Ko, 1962;Brockman, 1986), where The infinitesimal shear modulus  and infinitesimal Poisson ratio  are related through  =  ∕2(1+), where  is the infinitesimal Young modulus.
With the help of the plane-stress state condition (Haughton andMcKay, 1995, 1996),  3D 33 = 0, where  3D =  3D ∕ is the first Piola-Kirchhoff stress, we find the out-of-plane principal stretch ratio  3 as Now, by substituting Equation (4) in Equation (3), we obtain  , i.e. the membrane strain energy function, in terms of the two in-plane principal stretch ratios  1 and  2 , as follows Using Equation ( 5), we then compute the non-zero components of the first Piola-Kirchhoff stress  =  ∕ in the membrane as )) and the components of the principal Cauchy stress  =  −1   associated with the membrane energy, as Note that in axially symmetric membranes, the Cauchy stress has only in-plane components ( 1 ,  2 ) representing the principal stresses along meridians and parallels, respectively.With stretches given by Equation (1), these are functions of  only.
It is worth noting that the Cauchy stress components (7) satisfy the Baker-Ericksen inequality ( 2 − 1 )( 2 − 1 ) > 0 for all choices of the material parameters.This condition is fundamental to ensure that an inflated, elastically homogeneous spherical membrane in its reference configuration remains spherical for all values of applied pressure, see De Tommasi et al. (2012) andDe Tommasi et al. (2013).In this paper we demonstrate that an initially flat membrane with elastic inhomogeneities can attain non-trivial geometries (other than the spherical configuration) under a uniform pressure load.Furthermore, wrinkles and necks may be achieved at desired locations.

Tension Field Theory: Relaxed strain energy functional
Ideally, membranes with negligible bending stiffness may only achieve non-negative stress states as they cannot sustain in-plane compression.Real membranes, however, possess a small (albeit non-negligible) bending stiffness which effectively slightly delays the onset of wrinkling when compressive stresses arise.
Lack of resistance to compression may be seen as a unilateral constitutive constraint.This is embedded in tension field theory by constructing a "relaxed strain energy density"  ⋆ ( 1 ,  2 ) from the parent energy  ( 1 ,  2 ), which sets an in-plane stress component to zero whenever it would be negative in the parent energy.Clearly, in taut regions  ⋆ ≡  ( 1 ,  2 ), whereas in completely slack regions  ⋆ = 0.
Following A.C. Pipkin (Pipkin, 1986b), we may formalize the above by taking The function  ⋆  (  , ) is called "natural width in tension" and is the main player of the relaxed energy construction.Now consider a strip of membrane with energy Equation ( 5), oriented along the principal directions.If the membrane is pulled by so that  1 ≥ 1 along the direction 1, while being free in the direction 2, the membrane contracts laterally with stretch ratio  2 =  ⋆ 2 ( 1 , ).This value is found by solving  2 ( 1 ,  2 ) = 0 for  2 .Therefore if  2 ≤  ⋆ 2 ( 1 , ), then  2 ≤ 0, and the membrane would be compressed in the non-relaxed energy.In tension field theory, this problem is avoided by assuming that if the membrane is compressed further than  ⋆ 2 , its energy does not change once  * 2 is attained: namely,  ⋆ ( 1 ,  2 ) =  ( 1 ,  ⋆ 2 ( 1 , )).Physically, if  1 is kept fixed, there would be no energetic expenditure in shortening the membrane in the direction 2, below the natural width  ⋆ 2 .The same considerations apply to the perpendicular direction, whereas in the fully slack region, the energy is directly set to zero.For auxetic materials ( < 0), "lateral contraction" is changed to "lateral expansion", while all the remaining considerations are unchanged.
Remarkably, for membranes with energy defined by Equation ( 5), the natural width depends on the Poisson ratio  only, as follows, Note that for  = 1 ∕2, we recover the expression for the natural width obtained in Steigmann and Pipkin (1989a) and Steigmann and Pipkin (1989b) for incompressible and isotropic membranes.
The dependence of the natural width in tension on the Poisson coefficient  has not been explored much so far.However, this feature can be powerfully exploited to achieve non-trivial wrinkling patterns on demand.Indeed, as already illustrated by Venkata et al. (2023) for stretched membranes, by carefully tuning the spatial distribution of  and , one can achieve unusual wrinkling patterns.This has a great potential in technological applications (Naebe and Shirvanimoghaddam, 2016;Ren et al., 2018).
Auxetic materials have interesting properties: when  1 > 1, then  ⋆ 2 > 1, indicating that the membrane expands in all directions.This is a very strong difference between auxetic and classical membranes.Classical membranes always contract laterally when pulled in the perpendicular direction.In the next section, we discuss the implications of these features for both classical and auxetic functionally graded membranes.

Equilibrium of a pressurized membrane
Equilibrium equations of axially symmetric pressurized membranes are written along the meridian and normal directions (Green and Adkins, 1960;Gurtin and Ian Murdoch, 1975;Libai and Simmonds, 1998).If the reference membrane is flat, all involved fields depend on the reference radius  only, and the equilibrium equations take the form where (•) ′ = ∕ and where the membrane curvatures  1 (curvature along a meridian line) and  2 (curvature along a parallel line) in axial symmetry may be calculated as Note that for  2 we have also used the alternative expression based on the outward normal  to the current surface, and on the radial unit vector   pointing outwards radially and perpendicularly to the membrane axis .This expression exemplifies that the curvature of a parallel does not coincide, in general, with the curvature of the membrane along the parallel unless of course  ≡   .
Because the membrane is flat in its reference configuration, the fields ,  must satisfy the following boundary conditions, Also, because point loads are not applied to the membrane, the suitable boundary conditions to be imposed on the functions ,  to avoid stretch and curvature singularities at the origin are Finally, to account for the presence of wrinkling, it is sufficient to write the equilibrium equations in terms of the relaxed counterparts ( ⋆ 1 ,  ⋆ 2 ) as obtained from the relaxed energy through  ⋆ =  −1 ( ⋆ ∕)  .This will automatically ensure that no compressive states can be achieved on the inflated membrane.

General insights into necks and wrinkles
Due to constitutive and geometric non-linearities, it would appear that little can be said in general on the placement of wrinkles and necks in an inflated membrane undergoing large deformations.Surprisingly, a careful analysis of the equilibrium equations of an inflated membrane, together with the implementation of tension field theory, provides a deep and fully general characterization of what type of instability patterns may or may not be expected in an inflated membrane.
Even more remarkably, these characterizations are universal, in the sense that they are independent on the choice of the (isotropic and elastic) constitutive behavior of the material.One such remarkable universal result is that necks and wrinkles can never overlap: this insight will be used in the sequel to produce alternating patterns of regions with wrinkling and necking.
In this section we collect these universal characterizations, calling them "insights", that apply to all isotropic, axisymmetric and inflated membranes.Their relevance to the present study is that one can use them as general guidelines to design the spatial distribution of elastic moduli, in order to obtain desired patterns of wrinkles and necks in the inflated membrane.The insights also provide interesting characterizations of the membrane shape when wrinkles or necks occur, therefore helping to understand what type of shapes may, or may not, be obtained in an inflated membrane undergoing such instabilities.
As we deal with thin membranes that offer no resistance to compression, we use tension field theory and base our analysis on the equilibrium equations Equation ( 10), expressed in terms of the relaxed Cauchy principal stress components  ⋆ 1 ,  ⋆ 2 .
First note that the equilibrium of a cap above a parallel can be written, in global form, by balancing pressure  > 0 and membrane tension, to give Because, by definition,  ⋆ 1 ≥ 0, equilibrium imposes  ⋆ 1 > 0 and  2 < 0 everywhere, which concludes the proof.
• Insight 2: Wrinkles can never be aligned with parallels.
If such wrinkles were to occur it would result in  ⋆ 1 = 0, but from Insight 1 we know this is not possible.Therefore, wrinkles can only be aligned with meridians.
• Insight 4: When wrinkles occur, the membrane curvatures are linked through: This result follows from the previous point with  ⋆ 2 = 0.
• Insight 6: For necking to happen, the principal stresses must satisfy: By using Equations ( 14) to ( 15), we obtain the principal curvatures in terms of stresses For necking we need   < 0, and the result follows2 .
This insight may be also used to prove that wrinkles can never overlap necks.Indeed, for necks 2 ⋆ 1 ≤  ⋆ 2 , for wrinkles  ⋆ 2 = 0, and combining both we get 2 ⋆ 1 ≤  ⋆ 2 = 0, which cannot occur because  ⋆ 1 must be strictly positive.
• Insight 7: Neither wrinkling nor necking can take place at the apex of the membrane.
According to Equation ( 13), both principal stretches are equal at the center of the membrane, i.e.  1 =  2 at  = 0. Therefore, we have  ⋆ 1 =  ⋆ 2 and  1 =  2 at  = 0.Then, from Insights 4 and 6, we conclude that neither wrinkling nor necking can occur at the apex of the membrane.
Based on Insight 2, we therefore conclude that slack regions and regions of wrinkles aligned with parallels can never occur in the inflated membranes.These observations help us to simplify the relaxed energy function in Equation ( 8) as follows, We emphasize here that not all insights are directly helpful in providing guidance on how to design the functionally graded membrane to get desired instabilities.From Insights 1, 2, 3, 4, 5, 7 we deduce important information on the shape of an inflated membrane, with remarkable bounds on the membrane curvatures in regions with and without instabilities.These insights dictate what shapes may be expected (or not expected) in an inflated membrane, regardless of its constitutive response.
From Insight 6, together with the consideration that  ⋆ 2 = 0 leads to the appearance of wrinkles, we conclude that it is possible to create necks and wrinkles by tuning the ratio between  ⋆ 2 and  ⋆ 1 .This method is used in the following sections to create instability patterns at desired locations.

Designing instabilities
In principle, it is logical to assume that the spatial distribution of the elastic moduli (), () may be used to achieve wrinkles and necks.However, in practice, the inverse problem is extremely difficult due to geometric and constitutive nonlinearities.To circumvent the analytical complexity of this inverse problem, we use numerical simulations to devise regions of spatial inhomogeneity with sharp variations of the elastic moduli to get desired instabilities.
Across a boundary of radius  between elastically inhomogeneous regions, the following jump conditions have to be satisfied, where [ ] ∶= lim →0 ( ( + ) −  ( − )).We see that across such boundary, we may sharply act on the jump of the tension [ ⋆ 2 ] in the direction of parallels, while at the same time, the value of  ⋆ 1 along meridians shall remain continuous across the phase boundary.By tuning the ratio  stiff ∕ soft between the Young's modulus of two neighboring regions (stiff and soft), we can act on the ratio  ⋆ 2 ∕ ⋆ 1 which, as the insights above reveal, regulates the onset of wrinkles (when  ⋆ 2 ∕ ⋆ 1 is decreased toward 0) or the onset of necks (when  ⋆ 2 ∕ ⋆ 1 is increased above 2).The Poisson coefficient  is also a key player.Indeed, by taking this coefficient to be negative, the membrane will expand laterally and create compressive stresses ( ⋆ 2 = 0).The complexity of the problem requires a trial and error approach to determine appropriate elastic inhomogeneity patterns to achieve the sought instabilities.By using the commercial software COMSOL, we implement the Blatz-Ko relaxed strain energy function, Equation ( 8), for functionally-graded, inflated membranes.
In the simulations below, we inflate pre-stretched membranes by keeping the external rim fixed.Then we solve Equation ( 10) to obtain the Pressure-Volume (P-V) curves and the deformation profiles of the membranes.The numerical implementation is explained in Appendix A. Our simulations are consistent with the universal insights developed above: notably, we see that necking cannot occur in wrinkling regions and that only wrinkles aligned with meridians can occur.By playing on the ratio  stiff ∕ soft and on the location of elastically inhomogeneous regions, we show that a desired number of wrinkling and necking regions can be obtained at prescribed locations.
For the sake of brevity, the effects of pre-stretch and of the Blatz-Ko material parameters are not discussed in detail.Here, we focus on the most relevant design parameters:  and  stiff ∕ soft .

Necks and no wrinkles
To ensure that the inflated membrane develops necks without wrinkles, we have to create regions where the ratio  ⋆ 2 ∕ ⋆ 1 exceeds 2, while maintaining  ⋆ 2 > 0 everywhere in the membrane.To achieve this goal, in Section 4.1, we use a one-step variation of the Young modulus.The reference membrane has a central region with low Young's Modulus ( soft ) surrounded by a region of high Young's Modulus ( stiff ) which extends up to the boundary, as shown in Figure 3.Moreover, to avoid that  ⋆ 2 approaches zero, we set 0 <  < −0.6 which ensures a low-to-moderate auxeticity.Our simulations show that if  stiff ∕ soft is small (< 1.5), the membrane displays a smooth necking during inflation.If  stiff ∕ soft is moderate (1.5 <  stiff ∕ soft < 3), the necking region is more prominent.
On the left panel of Figure 4, we plot the P-V profiles of pre-stretched circular membranes with different sizes of the  soft region.From Case A to Case D, the radius or size of  soft region increases (see Equation (B.1b)).We observe that, as the size of the soft region increases, the value of the limit-point pressure  lim decreases.This clearly occurs because the bigger the soft region, the smaller the effective Young modulus.
The deformation profile at the end of inflation process is shown for each case on the right panel of Figure 4.It can be clearly noticed that, by adjusting the size of the low Young modulus region, we can shift the boundary of the necking region.Therefore necks may be obtained at any desired location.It has to be noted that the simulations confirm that necking cannot occur at the tip of the membrane, in agreement with Insight 7. Membranes with a high stiffness ratio ( stiff ∕ soft >3) undergo a transition into a double-bubbled shape, due to a sudden inflation of the softer region after the limit-point pressure, as shown in Figure 5.In this figure, we plot the P-V curves for two different cases of one-step variation of the Young modulus.

Special case: Double bubbling
It is important to note that the behavior of the P-V curve is strongly dependent on the size of the softer core.In Figure 5(a), the volume decreases immediately after the limit point, whereas in Figure 5(b) the volume keeps increasing while the pressure decreases after the limit point.Note that case (a) is also studied by Selvadurai (2022) for pressurized incompressible elliptical membranes with the Ogden strain energy model.
We also show the deformed shapes at four different stages of inflation.Note that the formation of the second bubble corresponds to the necking point on the P-V curve.
Highly auxetic membranes (with  < −0.8) develop wrinkles in the region between the two bubbles.Indeed, high auxeticity results into high compressive stresses.Note that in line with Insight 5, wrinkles (black color regions) are never observed in the necking region, but only slightly above the neck.

Necks and wrinkles
In this section, we extend the findings from Section 4.1 to the general case where the Young modulus has an -step variation across the membrane, while  = −0.9 is uniform throughout the membrane.Again, we use alternating stiff and soft regions to create multiple necks.Depending on the size of these regions we may or may not achieve lateral compression, and therefore wrinkling, in the membrane.Interestingly, the complexity of the problem leads to wrinkling regions appearing and disappearing during the inflation process.The simulations show that membranes with -step variations of the Young modulus develop  − 1 wrinkled regions during the pre-stretching process, consistent with the findings of Venkata et al. (2023).We recall that during prestretching, the membrane remains flat.These simulations are reported in Figure 6, showing the onset of wrinkling patterns in laterally pre-stretched membranes with two-and three-step variations of the elastic moduli.
Following the pre-stretch, we simulate the inflation.The deformation profiles for the inflation process are shown in Figure 7(b).Wrinkling and necking regions are represented by black (right half of the figure) and blue (left half) colors, respectively.Our simulations show that stiff regions inflate relatively less than soft regions.Moreover, this difference leads to the formation of a neck at the boundary between the two regions.Interestingly, the wrinkles generated during the pre-stretch phase tend to be preserved during the inflation process.Also, the wrinkles never overlap with the necking regions, which is consistent with the findings of Insight 5.
Furthermore, depending on the stiffness ratio between soft and stiff regions, the spacing between the necks and the wrinkles can be modulated.In particular, the higher the ratio  stiff ∕ soft , the closer the two patterns are, and vice-versa.The spatial heterogeneity of the elastic moduli affects also the P-V curve, as shown in Figure 7(a).Interestingly, for membranes with two-step variation of the Young modulus, the pressure increases until the limit-point instability and later drops, whereas it remains nearly constant in the three-step variation case.In all simulations, no further increase of pressure past the limit point is observed, irrespective of the material constant 0 <  < 1.
Finally, we highlight two important features.First, configurations containing multiple regions of necks and wrinkles can be obtained by increasing the number of stiff-soft regions in the original membrane.Second, the pre-stretch   may be used as a parameter to unwrinkle, during inflation, the regions of the membrane that are close to the boundary.

Wrinkles and no necks
Our results from the previous sections suggest that a high contrast between stiff and soft regions (i.e.high  stiff ∕ soft ) leads to necking.Therefore, to avoid necks, we now consider a homogeneous (or nearly homogeneous) distribution of the Young modulus across the membrane.Although different types of distributions (constant, linear, step, and Gaussian) of the Poisson ratio are explored for this case, they all yield similar deformation behaviors.
For the sake of brevity, we consider the case where the material properties ,  are spatially uniform throughout the membrane.We also compare our results with existing solutions in the literature.
In Figure 8, we see that for an auxetic membrane with homogeneous Young's modulus, the pressure increases with the volume up to the limit point, where the membrane loses stability and wrinkles appear near the base of the membrane.Beyond the limit point, the pressure decreases with increasing volume.The results in Figure 8 also show that the membrane attains a spherical shape, and indeed no necks are developed.These findings are consistent with previous experiments on pressurized incompressible circular membranes (Treloar, 1944;Machado et al., 2012;Zhou  et al., 2018) with several analytical and numerical studies (Yang and Feng, 1970;Patil and DasGupta, 2013;Selvadurai, 2022).
For auxetic membranes, we predict the formation of wrinkles aligned with meridians.We note that Chaudhuri and DasGupta (2014) predict wrinkles aligned with parallels in inflated incompressible membranes with homogeneous material properties.However, according to our Insight 2, such wrinkles are not possible in any axisymmetrically inflated isotropic circular membranes.
Additionally, we find that for a Blatz-Ko material with homogeneous material properties, tensile stresses exist everywhere in the membrane, except near the fixed circumferential edge.Hence, wrinkles are obtained only near the fixed boundary.
Finally, we mention that wrinkles, if present, may also be suppressed by tuning the Poisson ratio of the membrane and the pre-stretch applied to the membrane before inflating.We also observed separately that the limit-point pressure  lim increases linearly with the Blatz-Ko material constant  and nonlinearly with the pre-stretch, respectively.

Conclusion and perspectives
With this work we investigated the possibility of inducing various forms of instabilities in inflated hyperelastic membranes.Our results are valid for axially symmetric configurations.In particular, we investigated the effects of spatial material inhomogeneities on membranes with auxetic properties ( < 0).
We showed that wrinkling and necking instabilities can be induced in the membrane during inflation, by tuning the material inhomogeneities.We avoided solving the complex inverse problem of identifying the elastic heterogeneity patterns that deliver the desired instabilities at desired locations.Instead, we developed the following novel methods: 1. We implemented tension-field theory based on the relaxed strain energy function into COMSOL to study the stability of an inflated compressible, hyperelastic membrane of the Blatz-Ko type for an arbitrary geometry; 2. We obtained a series of universal results for the formation of wrinkles and necks in thin inflated membranes; 3. Building upon the universal results and the numerical simulations, we identified spatial inhomogeneity distributions across the undeformed membrane that result in wrinkles alone, necks alone, and simultaneous (but not overlapping) wrinkles and necks in the inflated membrane.
Our study is limited by various factors.For example, typical auxetic membranes are anisotropic in nature and their material properties are deformation-dependent: these features might greatly affect the results presented in this work.This aspect is not addressed in the current work but could be studied with the methods we have developed here.
Another avenue of interest is to investigate the inflation of non-axisymmetric membranes.For example, in Figure 9, we show 2D distributions of material properties in square membranes.The membrane is fixed on its edges and inflated under a uniform pressure load.We show that material properties can be tuned in the membranes to produce desired wrinkling patterns, as shown in the Figure 9.This concept could potentially be used in Braille reading and haptics.
Here, the left bottom vertex of the quarter square membrane is located at ( = 0,  = 0).The side length of the quarter square membrane in the undeformed configuration is ∕2.The smoothed Heaviside function with a continuous second derivative is represented by F , it is an inbuilt function in the COMSOL software.The function, |(•)|, returns the absolute value of any variable, (•).Similarly, for the full square membrane (FG-SQ-2) with side length  in the undeformed configuration, the Young modulus can be mathematically expressed as

Figure 2 :
Figure 2: Deformation profiles of the axisymmetric auxetic circular membrane.(a): Circular membrane with initial radius of  in and thickness .(b):Deformed profile of the circular membrane under radial pre-stretch with the radius  fin .Note that the thickness of the membrane is increased with the in-plane stretching due to auxeticity.(c): Deformed profile of the inflated circular membrane subjected to uniform pressure  and fixed on its circumferential boundary.Point P 1 in the undeformed configuration is displaced to the positions P 2 and P 3 with pre-stretch and inflation, respectively.A point in the final configuration has coordinates (, Φ, ), where  is radial position and  is the vertical deflection.

Figure 3 :
Figure 3: Schematic representation of the one-step variation of the Young modulus distribution  in the reference configuration.The circular membrane has a softer region in the center surrounded by a stiffer region that extends to the disc circumference, see Equation (B.1a).

Figure 4 :
Figure 4: (a) Inflation process: Pressure-Volume curves of a pre-stretched circular membranes with four different cases of one-step variation of the Young modulus.From Case A to Case D, the radius or size of  soft region is assumed to be increasing (see Equation (B.1b)),(b) Deformation profiles of the membrane at the end of inflation process are shown here for each case.Inflation process of each membrane is continued until the apex of the membrane reaches a height of five times the value of radius  in .For each profile, the region in blue highlights the necking area ( G < 0).Blatz-Ko parameter  = 0.4, pre-stretch   = 2, Poisson's ratio  = −0.4,and referential radius  in = 1 cm.

Figure 5 :
Figure 5: Pressure-Volume curves of pre-stretched circular membranes for two different cases of one-step variation of the Young modulus (see, Equation (B.1c)) are shown: (a) One-step variation-I, (b) One-step variation-II.The deformation profiles of the membrane are shown at four different points along the inflation path.The black region in the deformation profiles (at the end of inflation process) denotes wrinkling regions.Parameters:  in = 1 cm,   = 2,  = 0.4 and  = −0.9.

Figure 7 :
Figure 7: (a) Inflation process: Pressure-Volume curves of a pre-stretched circular membrane with two-step (solid curve) and three-step (dash-dotted curve) variations of the Young modulus, (b) Deformation profiles of the membrane at three different locations are shown for each variation case.For each profile, the left half highlights the necking area in blue ( G < 0) and the right half highlights the wrinkling area in black ( 2 = 0).Parameters:   = 2,  = 0.4 and  = −0.9.

Figure 8 :
Figure 8: Pressure-Volume curve of pre-stretched circular membranes with constant material properties.The deformation profiles of the membrane are shown at four different points along the inflation path.The black region in the deformation profiles denotes wrinkling regions.Parameters:  in = 1 cm,   = 2,  = 0.4,  = 70 kPa, and  = −0.9.

Figure 9 :
Figure 9: (a) Square membranes with two different Young's modulus profiles (top) and wrinkling patterns (bottom); (b) Associated Pressure-Volume curves.FG-SQ-1: The circular region in the center and the rest of the membrane have  = 0.3 MPa and  = 0.03 MPa, respectively.FG-SQ-2: The triangular and rectangular strips have  = 0.63 MPa and  = 4.03 MPa, respectively, while the rest of the membrane has  = 0.03 MPa.The Poisson ratio is  = −0.1,pre-stretch   = 1, and the Blatz-Ko coefficient is  = 0.4.Wrinkling patterns at two points Q 1 and Q 2 on the inflation curves corresponding to FG-SQ-1 and FG-SQ-2, respectively, are shown at the bottom of (a).