## Abstract

The equations governing the appearance of flexural static perturbations at the edge of a semi-infinite thin elastic isotropic plate, subjected to a state of homogeneous bi-axial pre-stress, are derived and solved. The plate is incompressible and supported by a Winkler elastic foundation with, possibly, wavenumber dependence. Small perturbations superposed onto the homogeneous state of pre-stress, within the three-dimensional elasticity theory, are considered. A series expansion of the plate kinematics in the plate thickness provides a consistent expression for the second variation of the potential energy, whose minimization gives the plate governing equations. Consistency considerations supplement a constraint on the scaling of the pre-stress so that the classical Kirchhoff–Love linear theory of pre-stretched elastic plates is retrieved. Moreover, a scaling constraint for the foundation stiffness is also introduced. Edge wrinkling is investigated and compared with body wrinkling. We find that the former always precedes the latter in a state of uni-axial pre-stretch, regardless of the foundation stiffness. By contrast, a general bi-axial pre-stretch state may favour body wrinkling for moderate foundation stiffness. Wavenumber dependence significantly alters the predicted behaviour. The results may be especially relevant to modelling soft biological materials, such as skin or tissues, or stretchable organic thin-films, embedded in a compliant elastic matrix.

## 1. Introduction

The edge buckling phenomenon is ubiquitous in nature and it can be observed at the boundary of almost all biological thin structures. Examples include, among many, lettuce leaves, flower petals or the gut tube in animals, see the review by Li *et al.* [1]. There, the driving force behind edge wrinkling is undoubtedly growth. In Kabir *et al.* [2], growth-induced buckling of cell microtubules embedded in a compressed kinesin substrate is shown experimentally and it is modelled as instability of a beam-plate supported by a Winkler elastic foundation. Alongside biological systems, thin-film flexible materials, with special regard to organic films, easily develop stress-induced instability and, with it, a great potential for integrated applications in moving parts and complex geometries. In particular, stretchability (i.e. elasticity under tensile strain) and flexibility have been identified as key material properties required to develop collapsible and portable devices [3], bio-sensors and textile integration [4], energy scavengers [5] and embedded capacitors and batteries. Lipomi *et al.* [6] experimentally investigated pre-stressed stretchable organic photovoltaic cells laid on an elastic substrate as an application of body buckling to increase integration compliance in portable devices.

From a mechanical standpoint, modelling of edge buckling is related to edge wave propagation, whose consideration dates back to 1960 and is now credited to Konenkov [7], despite a long history of discovery and rediscovery, see the overview by Norris *et al.* [8] and the more recent contributions [9–13]. However, to move from edge waves to edge *wrinkles*, we must add the effects of a large enough pre-deformation such that the edge wave speed drops to zero. In this way, a localized static solution might exist in the neighbourhood of this pre-deformation.

Alongside some established mathematical tools, such as Gamma convergence and the asymptotic method [14], results which consistently separate flexural and extensional effects can be obtained through a Taylor expansion of the potential energy in powers of the plate thickness *h*, as in [15–17]. Kaplunov *et al.* [18] used an asymptotic technique to analyse vibrations of thin elastic pre-stressed incompressible plates in the low-frequency limit *η*=*kh*≪1, where *k* is the wavenumber, for the special case of plane deformation. Pichugin & Rogerson [19] provided an extension to the three-dimensional case. Tovstik [20] considered the vibrations of a pre-stressed transversely isotropic infinite thin plate that is supported by an elastic foundation with inertial contribution. Remarkably, no consistent attempt at considering edge wrinkles in pre-stressed elastic plates can be traced in the literature, to the best of the authors’ knowledge.

In this paper, we derive the equations governing edge wrinkling of a homogeneously pre-stressed plate made of incompressible isotropic hyperelastic material, together with the corresponding boundary conditions, when the plate is bilaterally supported by a Winkler elastic foundation [21,22]. As in [17], we adopt a through-the-thickness expansion for the second variation of the plate energy, with the differences that our material is incompressible and the plate is elastically supported. These assumptions are introduced to better model soft solids and thin-films embedded in an elastic matrix. It is worth emphasizing that consideration of different scalings for the pre-stress ** σ** leads to a diverse mechanical response. In this paper, we assume that the pre-stress is small as it scales as

*h*

^{2}and the plate is thin, i.e.

*h*≪1. Besides, the Winkler foundation is soft and its stiffness

*κ*scales as

*h*

^{3}. As a result, a flexural behaviour for the supported plate is considered. By contrast, in [18,19] attention is set on a large pre-stress, for

*h*. Consequently, in-plane deformation (membrane regime) takes over. Indeed, we show that the scaling assumed for the pre-stress determines the leading term in the energy expansion, while consistency considerations suggest the proper scaling for the foundation stiffness.

The paper is structured as follows. Section 2 introduces the problem and presents the variational framework. We carry out the through-thickness energy expansion in §3 and minimize the second variation of the potential energy in §4. The plate governing equation as well as the boundary conditions are given in §5. We draw a comparison with the classical Kirchhoff–Love theory of pre-stressed plates in §6. In §7, we seek solutions in the form of edge wrinkles and derive the corresponding bifurcation curve. Body wrinkling is considered in §8 and its occurrence is compared with that of edge wrinkling as a function of the foundation stiffness and wavenumber dependence. Finally, conclusions are drawn in §9.

## 2. Formulation of the problem

Consider a hyperelastic plate *B* occupying the region *B*, named *equilibrium configuration*, of the three-dimensional Euclidean space *e*_{1}, *e*_{2}, *e*_{3}} denote a fixed orthonormal basis set for *x*_{1},*x*_{2},*x*_{3}}. The plate is incompressible and it has been homogeneously pre-deformed. The equilibrium configuration takes the form
*h*>0 denotes the plate thickness and the region *ω* in the plane *x*_{3}=0 is named the *plate mid-plane*. Here, we assume that *e*_{3} is a principal axis for the homogeneous pre-deformation, while no such provision is taken for *e*_{1} and *e*_{2}. Hence, the plate is pre-deformed by the application of a constant Cauchy stress ** σ** such that

*σ*

_{13}=

*σ*

_{23}=0, whereas the other shear stress components are generally non-zero (see figure 1 for the case of a uni-axial stress). For simplicity, we further assume

*σ*

_{33}=0.

Having been homogeneously pre-stretched, the plate undergoes a small incremental motion. Thus, the deformation reads
*χ*^{(0)} is the homogeneous pre-stretch and *ϵ**χ*^{(1)} the small incremental deformation, where |*ϵ*|≪1. Let ** F**=grad

*χ*^{(0)}be the homogeneous gradient of the pre-deformation.

In this paper, we focus on *flexural edge wrinkles* arising in a thin plate, where *thin* is to be understood in the sense that the plate thickness *h* is small compared with the wrinkle wavelength ℓ=2*πk*^{−1}, i.e. *kh*≪1. The plate is uniformly and bilaterally supported along *e*_{3} by an elastic Winkler foundation with stiffness *κ*>0. Besides, to fix ideas, we assume that the foundation reaction is directly applied to the plate mid-plane, although this restriction will prove unnecessary. Consequently, the plate top and bottom faces are stress-free, i.e.

The gradient of the incremental displacement (deprived of the small parameter *ϵ*) is
*instantaneous* elastic moduli, which is endowed with the major symmetry property [24], eqn (2.10)
*p* is a Lagrange multiplier due to the internal constraint of incompressibility and *Γ*_{lk} and summing over repeated indexes, we get
** χ** is a one-parameter family in

*ϵ*, we can write the

*potential energy*

*E*of the system as a function of

*ϵ*,

*W*=

*W*(

**+**

*F**ϵ*

**) is the elastic energy stored in**

*Γ**B*,

**denotes the traction applied on ∂**

*t**B*, the boundary of

*B*, and the last integral accounts for the contribution of the foundation. Following [15,17], when the potential energy

*E*(

*ϵ*) is expanded as a Taylor series about

*ϵ*=0, the first variation vanishes because the current configuration is an equilibrium one. The second variation of the potential energy is

*E*

_{B}′′(0), is developed in §4 while the Winkler foundation contribution is simply

*κ*>0 is named the

*Winkler modulus*or the foundation stiffness.

## 3. Through-the-thickness expansions

Unless otherwise stated, the summation convention over twice repeated indexes is adopted, with the understanding that all Greek subscripts take on values in the set {1,2}, while Roman subscripts range in the set {1,2,3}. A comma is used to denote partial differentiation with respect to the relevant coordinate, i.e. *w*_{,1}=∂*w*/∂*x*_{1}. We assume that the incremental fields admit the following through-the-thickness expansions in *x*_{3}∈[−*h*/2,*h*/2]
*a*
and
*b*
where ** v**,

**,**

*a***,**

*b***and**

*c**x*

_{1}and

*x*

_{2}(see also [25,26]). Note that, at leading order, the displacement of the mid-plane has been decomposed as the in-plane displacement

*v*

_{α}

*e*_{α}, plus the transverse displacement

*w*

*e*_{3}. Consequently, it may be assumed that

*v*

_{3}=0 without loss of generality.

For the gradient operator, we use the decomposition
** f**=

*f*

_{α,α}, of a tensor, (∇⋅

**)**

*Σ*_{j}=

*Σ*

_{αj,α}, and the two-dimensional gradient of a vector, ∇

**=**

*f**f*

_{i,β}

*e*_{i}⊗

*e*_{β}, may be defined. Then, from (2.4) we obtain

*a*

*b*

*c*or, in terms of components,

*a*

*b*

*c*Upon substituting (3.3) into (2.5), we obtain

*a*and

*b*and so forth.

The stress-free boundary conditions at the top and bottom surfaces of the plate, equations (2.3), extend to the incremental stress and give
*x*_{3}=±*h*/2. Then, adding and subtracting together the two conditions, we get, to leading order [17], eqn (51),
*a*
and
*b*
We observe that equation (3.8) shows that the assumption that the foundation reaction acts directly at the plate mid-plane may be abandoned with no harm provided that, as it will appear later, the foundation is soft and its reaction is

The incremental incompressibility condition, div *χ*^{(1)}=0, applied to the expansion (3.1a), gives
*x*_{3}, the coefficients of this polynomial in *x*_{3} vanish independently, i.e.

We only consider flexural deformations, so we may set *v*_{α}=0. It then follows, from the first of equations (3.11) that
*i*=*α* in (3.9a) and making use of equations (3.5a) and (3.7a), we obtain
*σ*_{33}=0, that
*i*=3 in (3.9a) yields *i*=*α* in equation (3.9b) and using (3.5b), (3.7b) and (3.12), we obtain
*i*=3 gives the leading term in the incremental pressure
** c** rests undetermined.

A comparison of the results with the literature shows that the plate kinematics (3.20) encompasses eqns (3.17,23,29,30) of Kaplunov *et al.* [18]. For instance, the linear-through-the-thickness-*ζ* expression for the axial displacement _{1}, using equation (3.31) and (3.33) corresponds to the one-dimensional version of the second term in (3.20) here, given that *w*. Likewise, equation (3.23)_{2} brings in the quadratic term in the transverse displacement *u*_{2}, corresponding to the third term in (3.20), which is proportional to the curvature. Finally, equation (3.23)_{3} gives the last of equations (3.19) for the leading term in the pressure increment. Conversely, the governing equation for pre-stressed plates (7.4) cannot be directly obtained from the static limit of (3.56) of Kaplunov *et al.* [18], in the light of the fact that the latter equation is obtained under the assumption of pre-stress *h* and in the absence of the foundation. However, once such assumptions are modified, correspondence can be achieved. We observe that in the works of Dai & Song [27] and Wang *et al.* [28] a theory for, respectively, compressible and incompressible thin plates is developed through an expansion of the plate kinematics about the lower surface *x*_{3}=−*h*/2, which is then fed into the governing equations.

## 4. Variational formulation

We begin with the following general expression for the second variation of the total potential energy of a hyperelastic body *B*:
*B* may be the configuration of any pre-stressed body and, in the following, it is identified with the set (2.1). Hereinafter, a general material and a general state of pre-stress *σ*_{αβ} are considered. Besides, in order to restrict the formulation to bending, we assume that the pre-stress scales as *h*^{2}, i.e. *h*^{3}, and then obtain the reduced boundary-value problem by energy minimization [17]. Consequently, in the following derivation, terms of order higher than *h*^{3} are neglected.

With the help of the results established in the previous sections, we can proceed to simplify the second variation (4.1). First, by substituting (3.3) and (3.6) into (4.1) and integrating along the thickness of the plate, we obtain, up to *n*_{α} is the unit vector normal to the mid-plane boundary ∂*ω*. On the account of the incremental equilibrium equation in the absence of incremental body forces, i.e. div ** Σ**=

**o**and with (3.9b), it is

Finally, we have, with the help of equations (3.5b) and (3.18) and up to *w* and ∇∇*w*, and where
*h* about the undeformed state *λ*_{i}=1 only the leading order term may be consistently retained, i.e. *p*=*μ*, whence equation (4.9) gives
*plate flexural rigidity* [30]
*μ*=*E*/(2(1+*ν*)) between the shear modulus, *μ*, and Young’s modulus, *E*. In this context, the only non-zero elements of

## 5. Plate governing equation

Let us define the moments (per unit length) by
*bending moments* are
*twisting moments* are
*shearing force* (per unit length)
*w*, we obtain

Setting to zero the first variation inside *ω* gives the *plate governing equation*
^{2}*w*=*w*_{,1111}+2*w*_{,1122}+*w*_{,2222} is the biharmonic operator applied to *w*. In the same fashion, setting to zero the first variation of the boundary integral yields the *natural boundary conditions*
*V* _{n}=*q*_{α}*n*_{α}+(∂/∂*τ*)(*τ*_{β}*M*_{βα}*n*_{α}) is the well-known *Kirchhoff equivalent shearing force* and *M*_{n}=*n*_{β}*M*_{βα}*n*_{α} is the *bending moment*, both acting on a surface with unit normal *n*_{α}. Setting

## 6. Stress distribution and comparison with the classical theory

Equations (3.15), (4.7) and (4.9) give the normal stress distribution (no sum over *α* in this section)
*Σ*_{3α} consistently at *c*_{α} may be obtained, to leading order, through the incremental equilibrium equation div ** Σ**=

**o**, i.e.

*Q*

_{x}and

*Q*

_{y}, respectively, in eqns (106,107) of [30], corrected to incorporate the pre-stress. The bending moments (5.2) as well as the twisting moment (5.3) correspond to the classical definitions (101,102) of Timoshenko & Woinowsky-Krieger [30], §21. The Kirchhoff shearing force (5.11) amounts to the corresponding definition (g) of [30], provided that we take

*M*

_{12}=

*M*

_{yx}=−

*M*

_{xy}. The governing equation (5.9) coincides with the classical equation for combined bending and compression (or tension) of thin plates, first derived by Saint Venant [32]. Finally, we observe that equation (3.20), up to

## 7. Edge wrinkling solution

Let us consider the hyperelastic plate *B* to occupy the semi-infinite region (figure 1)
*e*_{i} are directed along the principal axes of the underlying homogeneous deformation, in which case the deformation gradient has a diagonal representation, namely
*x*_{2}=0, we have *θ*=*π*/2 and the boundary conditions (5.11) and (5.12) now specialize to
*E*′′(0) to be consistent, we need to assume that *κ*_{1}=*κ*/*h*^{3} is of order

We look for a *wrinkling solution* to equation (7.4), which varies sinusoidally along the edge and decays away from it as *A*_{1}, *A*_{2} are yet-undetermined amplitude constants, *γ*_{1}, *γ*_{2} the attenuation coefficients, such that ℜ(*γ*_{α})>0 and *k*>0 is the wavenumber (see figure 2 for an illustration). Substitution into the plate equation (7.4) shows that *γ*_{1} and *γ*_{2} are the roots of the bi-quadratic equation
*d*_{0}=D/*h*^{3}=*μ*/2, *σ*=(*σ*_{11}+*σ*_{22})/(*kh*)^{2} and

Enforcing the boundary conditions (7.5) gives a homogeneous linear algebraic system of two equations in the two unknowns *A*_{1} and *A*_{2}
*bifurcation condition*,

The bifurcation condition (7.9) is satisfied whenever *γ*_{1}=*γ*_{2} or when the term in square brackets vanishes. It can be shown that the former case is spurious, whereas the latter yields
*γ*_{α})>0. Equation (7.10) expresses the *bifurcation criterion* for the appearance of wrinkles on the edge of a semi-infinite plate compressed by a lateral stress. It can be rationalized by squaring and then solved to yield
*κ*_{1}=0), equation (7.11) is trivially satisfied by *kh*=0, which amounts to a zero pre-stress condition and, consequently, to a critical buckling stretch of *λ*=1. In other words, without the substrate, the plate would buckle as soon as it is laterally compressed.

Substituting the expansion *h*^{2}, we find the bifurcation condition
*characteristic length* *h*/ℓ_{s}. Typically, curves go through a maximum which determines the effective critical stretch of contraction as well as the edge wrinkles wavenumber. A similar pattern is shown in figure 4 where the foundation stiffness is taken to be proportional to the wavenumber *k*, which is the situation of a linear elastic half-space foundation considered in [34,20]. By contrast, when the foundation stiffness scales as the wavenumber squared, bifurcation curves become straight lines and edge wrinkling starts at zero wavenumber. This is indeed the nonlinear correction considered by Brau *et al.* [34] and it may be argued that appearance of edge wrinkles at zero wavenumber is a good test for such an assumption. Furthermore, any dependence of the foundation stiffness on powers of the wavenumber greater than 2 leads to bifurcation at *λ*_{1}=1 and zero wavenumber.

## 8. Body versus edge wrinkling

The body wrinkling solution takes the form
*n*_{2}, thus giving the bifurcation curves for bulk wrinkling. In the special case of uni-axial pre-stress along *x*_{1} (i.e. *n*_{2}=0, *n*_{1}=1), equation (8.2) simplifies to
*n*_{2}=0 at a critical stretch *λ**_{1 body} which is a little smaller than the corresponding threshold for edge wrinkling *n*_{2}=0, equation (8.2) gives
*λ*_{2}>1. Figure 5 shows that body buckling may be preferred to edge buckling in a bi-axially pre-stretched scenario, where transverse extension *λ*_{2}=1.05 favours body wrinkle formation in compression, i.e. *λ*_{1}<1, up to moderate values of foundation compliance. Similarly, figure 6 compares body and edge wrinkling for a thin plate supported by a Winkler foundation whose stiffness is proportional to the wavenumber *k*: in this situation, body wrinkling is preferred to edge wrinkling up to large values of foundation compliance.

## 9. Conclusion

In this paper, a consistent model for flexural edge wrinkling of bi-axially pre-stressed thin incompressible elastic plates, supported by a local elastic foundation, is developed and solved. The governing equations and boundary conditions are derived from minimizing a reduced form for the second variation of the potential energy, which is obtained by expanding the three-dimensional kinematics through the plate thickness. Small deviations from the homogeneously pre-stressed state are investigated. Consistency of the expansion demands that, to obtain purely flexural deformations, the pre-stress scales as the plate thickness squared. Besides, it demands that the foundation stiffness scales as the plate thickness cubed. Within such assumptions, the classical Kirchhoff–Love theory of pre-stressed elastic plates is obtained. Furthermore, the ad hoc assumptions on the stress distribution are also retrieved, while parabolic and cubic terms are introduced to correct the classical linear (along the thickness) plate kinematics. Edge wrinkling is described by means of bifurcation curves and it is compared with body wrinkling. It is found that edge wrinkling always occurs prior to body wrinkling in a uni-axially pre-stressed situation, regardless of the foundation stiffness. The bifurcation landscape is more involved in a bi-axial condition and body wrinkling may precede edge wrinkling for moderate foundation stiffness. The situation where the foundation reaction depends on the wavenumber *k* is also discussed. In particular, it is observed that dependence on *k*^{2} determines buckling at zero wavenumber while dependence on *k*^{β}, *β*>2 produces buckling at *λ*=1. Such results may be employed to infer the mechanical behaviour of the supporting matrix in flexible embedded systems, with special regard to biological tissues or organic thin-films.

## Authors' contributions

M.D. and A.N. developed the model and the incremental approach; Y.F. carried out the variational reduction. All authors gave final approval for publication.

## Competing interests

We declare we have no competing interests.

## Funding

M.D. gratefully acknowledges partial funding from the Università degli Studi di Modena e Reggio Emilia. A.N. is grateful to the National Group of Mathematical Physics (GNFM-INdAM) for partial support through the ‘Progetto Giovani Ricercatori 2015’ scheme, prot. U2015/000125.

- Received May 27, 2016.
- Accepted August 4, 2016.

- © 2016 The Author(s)

Published by the Royal Society. All rights reserved.