Tunable wrinkling of thin nematic liquid crystal elastomer sheets

TL;DR

This paper develops a theoretical framework for understanding and controlling wrinkle formation in thin nematic liquid crystal elastomer sheets, revealing how nematic order influences wrinkle wavelength and critical stress.

Contribution

It introduces a plate theory for nematic elastomers and analyzes how nematic order affects wrinkle wavelength and critical stress under compression.

Findings

Wrinkle wavelength scales with material parameters and compression.

Nematic order causes non-monotonic wrinkle wavelength behavior.

Critical stress for wrinkling can be significantly higher due to soft modes.

Abstract

Instabilities in thin elastic sheets, such as wrinkles, are of broad interest both from a fundamental viewpoint and also because of their potential for engineering applications. Nematic liquid crystal elastomers offer a new form of control of these instabilities through direct coupling between microscopic degrees of freedom, resulting from orientational ordering of rod-like molecules, and macroscopic strain. By a standard method of dimensional reduction, we construct a plate theory for thin sheets of nematic elastomer. We then apply this theory to the study of the formation of wrinkles due to compression of a thin sheet of nematic liquid crystal elastomer atop an elastic or fluid substrate. We find the scaling of the wrinkle wavelength in terms of material parameters and the applied compression. The wavelength of the wrinkles is found to be non-monotonic in the compressive strain owing to the presence of the nematic. Finally, due to soft modes, the critical stress for the appearance of wrinkles can be much higher than in an isotropic elastomer and depends nontrivially on the manner in which the elastomer was prepared.

Equations

$$\mathbf{\epsilon }$$

$${\mathbf{\epsilon }}^{\ast }$$

$$E$$

$$\mathbf{Y}(p)=\mathbf{S}(p^{\text{mbox}(0)})+x_3{\mathsf{F}}^{\text{mbox}(0)}{\mathbf{e}}_3+\cdots ,$$

$$\mathbf{\epsilon }={\mathbf{\epsilon }}^{\text{mbox}(0)}-x_3{\mathsf{F}}^{\text{mbox}(0)\mathrm{T}}\mathsf{b}{\mathsf{F}}^{\text{mbox}(0)}+\frac{1}{2}{x}_3^2{\mathsf{F}}^{\text{mbox}(0)\mathrm{T}}\mathsf{c}{\mathsf{F}}^{\text{mbox}(0)},$$

$${\mathbf{\epsilon }}^{\text{mbox}(0)}=\frac{1}{2}\left[{\mathsf{F}}^{\text{mbox}(0)\mathrm{T}}{\mathsf{F}}^{\text{mbox}(0)}-1\right]=\frac{1}{2}\left(a_{\alpha \beta }-{\delta }_{\alpha \beta }\right){\mathbf{e}}^{\alpha }\otimes {\mathbf{e}}^{\beta },$$

$${\mathsf{F}}^{\text{mbox}(0)\mathrm{T}}\mathsf{b}{\mathsf{F}}^{\text{mbox}(0)}=b_{\alpha \beta }{\mathbf{e}}^{\alpha }\otimes {\mathbf{e}}^{\beta },$$

$${\mathsf{F}}^{\text{mbox}(0)\mathrm{T}}\mathsf{c}{\mathsf{F}}^{\text{mbox}(0)}=c_{\alpha \beta }{\mathbf{e}}^{\alpha }\otimes {\mathbf{e}}^{\beta },$$

$${\mathbf{\epsilon }}^{\ast }$$

$${\mathbf{\epsilon }}^{\text{mbox}(0)\ast }=\frac{1}{2}\left[1^{\ast }-{\mathsf{F}}^{\text{mbox}(0)-\mathrm{T}}{\mathsf{F}}^{\text{mbox}(0)-1}\right]=\frac{1}{2}\left(a_{\alpha \beta }-{\delta }_{\alpha \beta }\right){\mathbf{a}}^{\alpha }\otimes {\mathbf{a}}^{\beta }$$

$$\mathsf{b}=b_{\alpha \beta }{\mathbf{a}}^{\alpha }\otimes {\mathbf{a}}^{\beta }$$

$$\mathsf{c}=c_{\alpha \beta }{\mathbf{a}}^{\alpha }\otimes {\mathbf{a}}^{\beta }.$$

$${\sigma }_{ij}=\delta E/\delta {\epsilon }_{ij}=\lambda {\delta }_{ij}{\epsilon }_{kk}+2\mu {\epsilon }_{ij}-\alpha \Delta Q_{ij},$$

$$\mathbf{n}=n_{\alpha }{\mathbf{g}}^{\alpha }=\left[n_{\alpha }{\left({\mathsf{T}}^{-1}\right)}_{\beta }^{\alpha }\right]{\mathbf{a}}^{\beta }\equiv n_{\beta }^{\text{mbox}(0)}{\mathbf{a}}^{\beta },$$

$$n_{\alpha }{n}_{\beta }$$

$${\mathcal{E}}_{\mathrm{s}}$$

$${\mathcal{E}}_{\mathrm{b}}$$

$$+h\left(B_{\alpha \beta }^2+C_{\alpha \beta }\Delta Q_{\alpha \beta }^{\text{mbox}(0)}\right),$$

$${\sigma }_{xx}^{\text{mbox}(0)}$$

$${\sigma }_{yy}^{\text{mbox}(0)}$$

$${\sigma }_{xy}^{\text{mbox}(0)}$$

$$\tan 2\varphi =\frac{2{\epsilon }_{xy}^{\text{mbox}(0)}}{{\epsilon }_{xx}^{\text{mbox}(0)}-{\epsilon }_{yy}^{\text{mbox}(0)}+h/\alpha }.$$

$${\epsilon }_{xy}^{\text{mbox}(0)}$$

$${\epsilon }_{yy}^{\text{mbox}(0)}$$

$${\varphi }_{\mathrm{c}}=\frac{1}{2}{\cos }^{-1}\left[\frac{8h\mu -{\alpha }^2-12\alpha \gamma \mu }{3{\alpha }^2}\right].$$

$${\partial }^2{\mathcal{E}}_{\mathrm{s}}/\partial {\varphi }^2=-3{\alpha }^2\cos (4\varphi )-({\alpha }^2-8h\mu +12\alpha \gamma \mu )\cos (2\varphi )\ge 0,$$

$${\gamma }_1$$

$${\gamma }_2$$

$${\sigma }_{xx}^{\text{mbox}(0)}({\varphi }_c)=\alpha -2\frac{h}{\alpha }\mu ,$$

$${\mathcal{E}}_{\mathrm{b}}=\frac{t^3}{24}\frac{q}{2\pi }{\int }_0^{2\pi /q}\mathrm{d}x\mathcal{B}(\varphi ){\left({\partial }_x^2\zeta \right)}^2=\frac{t^3\mathcal{B}(\varphi )q^4{A}^2}{48},$$

$$\mathcal{B}(\varphi )=4\mu -\frac{\alpha }{2}\left[3+5\cos (2\varphi )\right]+\frac{h}{2}\left[5+3\cos (2\varphi )\right]{\cos }^2\varphi .$$

$${\mathcal{E}}_{\mathrm{s}}$$

$${\mathcal{E}}_{\mathrm{e}}=\frac{q}{2\pi }{\int }_0^{2\pi /q}\mathrm{d}x\frac{1}{2}Kq{\zeta }^2=\frac{1}{4}Kq{A}^2,$$

$$4\mathcal{E}=t\mu q^4{A}^4+q^2{A}^2\left(\frac{t^3\mathcal{B}(\varphi )q^2}{12}+t{\sigma }_{xx}^{\text{mbox}(0)}+\frac{K}{q}\right).$$

$${\xi }_{\mathrm{e}}=\frac{(tq{)}^2}{12}+\frac{K}{\mathcal{B}(\varphi )qt}.$$

Full text

I Introduction

The wrinkling of a thin film is a common occurrence, but elucidation of its underlying physics has far-reaching implications Li et al. (2012). These include suggesting better ways to fabricate functional surfaces that change shape in predictable ways when subject to appropriate stimuli Tokareva and Minko (2009), finding new and more accurate methods for measuring material properties Huang et al. (2007); Chung et al. (2011), determining growth and form in biological tissues Thompson (1992); Liang and Mahadevan (2011), and bringing insights into the fundamental mechanisms of pattern formation Amar and Goriely (2005); Davidovitch et al. (2011); Vella et al. (2011); King et al. (2012); Paulsen et al. (2015). Many recent studies have focused on how wrinkling arises in isotropic materials as they respond to external effects, for example, when tension is applied to a freely suspended thin sheet Cerda et al. (2002); Cerda and Mahadevan (2003); Géminard et al. (2004), from compressing a film atop a soft substrate Chen and Hutchinson (2004); Huang et al. (2005); Audoly and Boudaoud (2008) or a fluid interface Huang et al. (2010); Brau et al. (2013); Oshri et al. (2015), or through swelling mechanisms Guvendiren et al. (2009). However, less attention has been given to wrinkling instabilities in anisotropic materials, such as Nematic Liquid Crystal Elastomers (NLCE) Küpfer and Finkelmann (1991); Warner and Terentjev (2003).

NLCE are promising materials for studying and controlling surface pattern formation White and Broer (2015); Mostajeran et al. (2016); Giminez-Pinto et al. (2017); Aharoni et al. (2018); Warner and Mostajeran (2018); Kowalski et al. (2018); Lavrentovich (2018); Babakhanova et al. (2018) due to unique effects that arise from the coupling between strain and nematic order Küpfer and Finkelmann (1994); Urayama et al. (2006); Rogez et al. (2018). For example, NLCE can display elastic soft modes, in which they can undergo deformation with negligible energy cost Kundler and Finkelmann (1995); Ware et al. (2016). These materials, therefore, represent an extension of traditional anisotropic plates Rand and Rovenski (2004); Lekhnitsky (1957), where the angle of anisotropy is now free to vary while the plate remains fixed. It has been shown that, under tension, thin sheets of NLCE can also develop intricate microstructures which suppress pattern formation at larger scales Plucinsky and Bhattacharya (2017). Because the nematic phase is easily affected by electromagnetic fields, strain-order coupling also allows for nuanced actuation of these materials via light DeSimone et al. (2016); Tiwari et al. (2018); Guin et al. (2018); Auguste et al. (2018), and one can also measure and record their deformation electronically Rihani et al. (2018). They are also an excellent candidate for artificial muscle tissue Martella and Parmeggiani (2018).

Although the mechanics of NLCE has been considered by several authors Warner and Mahadevan (2004); Warner et al. (2010a); Aharoni et al. (2014); Cirak et al. (2014); Agostiniani and Desimone (2017); Keip and Nadgir (2017); Nguyen and Selinger (2017), a first principle derivation of a dimensionally reduced model from an effective theory for these materials in the limit of thin plates is still lacking. In this article we derive a Föppl-von Kármán-like plate theory for NLCE and use this theory to study the wrinkling of such thin materials under compressive loads. The specific problem to be addressed is when a thin sheet of NLCE is placed on a soft isotropic elastic foundation, or a fluid sub-phase. Our problem is thus an inversion of already existent experiments and theory on wrinkling of an isotropic plate atop a thick nematic elastomer foundation Agrawal et al. (2012); Soni et al. (2016). As depicted in Fig. 1, the wrinkling of the plate atop the foundation is induced by compressing both the foundation and the plate as a unit. In the absence of a foundation, the plate will choose an out-of-plane deformation which is a single arch; if the plate is bonded to the foundation, however, this deformation has a high energetic cost, due either to deformation of an elastic foundation or the gravitational energy stored in a fluid substrate.

In Sections II and III, we outline the dimensional reduction. We begin with a three-dimensional phenomenological energy and extend the kinematic and dynamic Kirchhoff-Love approximations with compatible approximations for the behavior of the nematic director. These approximations allow us to integrate across the thickness of the plate, thus reducing the three-dimensional energy to a two-dimensional energy comprised of two terms—a stretching-like term proportional to the thickness $t$, and a bending-like term proportional to $t^{3}$. These terms couple the geometry of the macroscopic sheet to the director orientation. We explore the ramifications of our model by considering a common wrinkling ansatz, in which a base state of stress for a compressed but planar plate is energetically coupled to the strain resulting from out-of-plane displacement. The base state of stress is derived in Section IV. The approximations for the bending and stretching energies due to the ansatz are discussed in Section V, where the scaling of the wrinkle wavelength is also considered. The scaling has the same functional form as for an isotropic plate, but includes information from the nematogen-elastomer coupling strengths. We summarize our results in Section VI.

II Three-dimensional energy

We assume that the elastomer is prepared in such a way that the material was deep in the nematic phase at the time of cross-linking, and denote by $\bar{\mathbf{n}}$ the unit vector director field in this initial configuration, also known as the “reference” configuration. This assumption is important in NLCE because there is an extra term in the energy for the elastomer depending on whether the material was prepared and cross-linked with the liquid crystal in the isotropic phase Lubensky et al. (2002); Warner and Terentjev (2003). When the NLCE is subject to external forces or torques, there exists an energetic cost to deform the elastomer and rotate the director field to a new configuration $\mathbf{n}$. To describe this energy, we review standard notions from continuum mechanics.

The reference configuration of the elastomer is a thin slab. We define an orthonormal basis $\{\mathbf{e}_{i}\}$ with $\mathbf{e}_{1}$ and $\mathbf{e}_{2}$ parallel to the plane of the slab, and $\mathbf{e}_{3}$ perpendicular to the slab. Cartesian coordinates $x_{1}$, $x_{2}$, $x_{3}$ determine the positions of points relative to this frame. The $x_{i}$ will also serve as Lagrangian coordinates for the deformed configuration, where the position of a material point $p$, originally at $\mathbf{x}$, is given by the deformation map $\mathbf{Y}(\mathbf{x})$. The deformation map defines another natural basis $\{\mathbf{g}_{i}\}$, where $\mathbf{g}_{i}=\partial_{i}\mathbf{Y}\equiv\partial\mathbf{Y}/\partial x_{i}.$ Due to the deformation, the distance between material points changes; the distance between material points is encoded in the metric tensor $\boldmath{\mathsf{g}}=g_{ij}\mathbf{e}^{i}\otimes\mathbf{e}^{j}$, where $g_{ij}=\mathbf{g}_{i}\cdot\mathbf{g}_{j}$. Denoting the inverse of the metric tensor by $g^{ij}$, we can define the dual basis $\{\mathbf{g}^{i}\}$ via $\mathbf{g}^{i}=g^{ij}\mathbf{g}_{j}$ so that $\mathbf{g}^{i}\cdot\mathbf{g}_{j}=\delta^{i}_{j}$. (Note that we use the Einstein convention of summing repeated indices.) For a deformed configuration, $\mathbf{g}^{i}\neq\mathbf{g}_{i}$, but $\mathbf{e}^{i}=\mathbf{e}_{i}$ for the Cartesian basis. The deformation tensor $\mathsf{F}$ maps the Cartesian basis vectors $\mathbf{e}_{i}$ onto the deformed basis vectors $\mathbf{g}_{i}$, i.e. $\mathsf{F}:\mathbf{e}_{i}\mapsto\mathbf{g}_{i}$. Thus, $\mathsf{F}=\mathbf{g}_{i}\otimes\mathbf{e}^{i}$, where $\otimes$ denotes the outer product. Note that $\mathsf{F}=\nabla\mathbf{Y}=\partial_{j}Y_{i}\mathbf{e}^{i}\otimes\mathbf{e}^{j}$, or $F_{ij}=\partial_{j}Y_{i}$. We can define the Lagrangian or Green strain tensor $\boldsymbol{\varepsilon}$ in terms of reference-space quantities $\mathrm{d}Y_{i}\mathrm{d}Y_{i}-\mathrm{d}x_{i}\mathrm{d}x_{j}=2\varepsilon_{ij}\mathrm{d}x_{i}\mathrm{d}x_{j}$, and the Eulerian or Almansi strain tensor $\boldsymbol{\varepsilon}^{*}$ in terms of deformed- or “target”-space quantities, $\mathrm{d}Y_{i}\mathrm{d}Y_{i}-\mathrm{d}x_{i}\mathrm{d}x_{j}=2\varepsilon^{*}_{ij}\mathrm{d}Y_{i}\mathrm{d}Y_{j}$:

[TABLE]

where $\mathsf{1}=\mathbf{e}_{i}\otimes\mathbf{e}^{i}$ is the identity on the reference space and $\mathsf{1}^{*}=\mathbf{g}_{i}\otimes\mathbf{g}^{i}$ is the identity on the target space.

The resistance of the cross-linked network to deformation is described by an energy density which is quadratic in the strain. The coupling between the orientational order and the strain of the polymer network is described by de Gennes’ phenomenological energy deGennes (1975); Mbanga et al. (2010); Sawa et al. (2011):

[TABLE]

where $\lambda$ and $\mu$ are the Lamé coefficients, and $\alpha$ and $h$ are coupling constants. The term proportional to $\alpha$ is the lowest-order term involving the nematic directors and strains. Since the reference director $\bar{\mathbf{n}}$ is a vector under the rotation group in the reference space, it must couple with the Green strain tensor in order to make the energy density a scalar Lubensky et al. (2002), which gives the term $\bar{\mathbf{n}}\cdot\boldsymbol{\varepsilon}\cdot\bar{\mathbf{n}}$. Similarly, the deformed director $\mathbf{n}$ must couple with the Almansi strain tensor, thus $\mathbf{n}\cdot\boldsymbol{\varepsilon}^{*}\cdot\mathbf{n}$. The thermomechanical history of the sample is captured by the term proportional to $h$ Urayama et al. (2009). For materials prepared deep in the nematic phase, the cross-links themselves resist rotating away from the reference state and, therefore, they “remember” the orientation of the nematic at the time of cross-linking Warner and Terentjev (2003). A low value of $h$ corresponds to a soft nematic elastomer, which has a low cross-link density, whereas a high value of $h$ corresponds to a nematic glass with high cross-link density. The function $f(\bar{\mathbf{n}},\mathbf{n})$ captures the energetic cost of rotating the director relative to the polymer matrix background. The symmetries also restrict the possible forms of the function $f(\bar{\mathbf{n}},\mathbf{n})$, since it must behave as a scalar under rotations in both the reference and deformed spaces. Since $\bar{\mathbf{n}}=\bar{n}_{i}\mathbf{e}^{i}$ and $\mathbf{n}=n_{i}\mathbf{g}^{i}$, we choose to construct $f(\bar{\mathbf{n}},\mathbf{n})$ by using the deformation map in such a way that $\mathbf{n}$ is mapped back to the reference space. Therefore, we define $f(\bar{\mathbf{n}},\mathbf{n})\equiv\mathrm{Tr}\Delta\mathsf{Q}^{2}$, where $\Delta\mathsf{Q}\equiv\left(\mathsf{F}^{-1}\mathbf{n}\right)\otimes\left(\mathsf{F}^{-1}\mathbf{n}\right)-\bar{\mathbf{n}}\otimes\bar{\mathbf{n}}$.

A derivation of the coupling between strain and orientation order starting with the neoclassical NLCE energy Warner and Terentjev (2003) reveals that $0<\alpha<3\mu/2$ Uchida and Onuki (1999). The only restriction on $h$ is that it must be positive.

III Dimensional Reduction

The dimensional reduction procedure considered here follows standard methods available for isotropic thin plates or shells Koiter (1970); Ciarlet (2005). We assume that the volume element may be decomposed into $\mathrm{d}V=\mathrm{d}\mathcal{A}^{\mbox{\tiny(0)}}\,\mathrm{d}{x_{3}}$, where $\mathrm{d}\mathcal{A}^{\mbox{\tiny(0)}}=\mathrm{d}{x_{1}}\mathrm{d}{x_{2}}$ is the area element of the middle reference-plane and $\mathrm{d}{x_{3}}$ is the integration measure through the thickness, ${x_{3}}\in[-t/2,t/2]$. Furthermore, simplifications are made if we assume two a priori conditions John (1965); Koiter (1970); Ciarlet (2005): (i) a mechanical condition, which states that the state of stress inside the body as purely parallel to the mid-surface; and (ii) a kinematic condition, known as the Kirchhoff-Love condition, which states that points located along the normal to the middle plane remain along the normal after the deformation to an arbitrary surface, while their distance with respect to this mid-surface does not change. In order to complement (i) and (ii) in light of the problem at hand, we impose a third condition, which is a kinematic constraint on the distribution of nematic directors: (iii) on every surface of constant $x_{3}$, we shall assume that both the reference and deformed director field remain tangent to the mid-surface of the film, i.e. $\bar{n}_{3}=0$ and $n_{3}=0$. Note that this last equality applies to the director in the deformed space.

The embedding of the plate in the deformed configuration is explicitly given in normal coordinates by assuming that a point $p$ (such that $\mathbf{Y}:p=(x_{1},x_{2},x_{3})\mapsto\mathbb{R}^{3}$) in the body is written in terms of a point in the mid-surface $p^{\mbox{\tiny(0)}}$ (such that $\mathbf{S}:p^{\mbox{\tiny(0)}}=(x_{1},x_{2})\mapsto\mathbb{R}^{3}$) through the relationship $p=p^{\mbox{\tiny(0)}}+x_{3}\mathbf{e}_{3}$, where $\mathbf{Y}(p^{\mbox{\tiny(0)}})=\mathbf{S}(p^{(0)})$. It is, therefore, a critical step in the dimensional reduction to express the embedding evaluated at $p$ as a Taylor series about points in the mid-surface:

[TABLE]

were we have used the definition $\mathsf{F}^{\mbox{\tiny(0)}}\equiv\left.\nabla\mathbf{Y}\right|_{p^{\mbox{\tiny(0)}}}$. We note that $\mathsf{F}^{\mbox{\tiny(0)}}$ is the deformation gradient of the mid-surface and, given the corresponding basis set on the mid-surface $\mathbf{a}_{i}=\left\{\mathbf{a}_{\alpha},\mathbf{a}_{3}\right\}$, we may write $\mathsf{F}^{\mbox{\tiny(0)}}=\mathbf{a}_{\alpha}\otimes\mathbf{e}^{\alpha}+\mathbf{a}_{3}\otimes\mathbf{e}^{3}$, where $\{\alpha,\beta,\cdots\}$ take values in the set $\{1,2\}$. An additional clarification is that the basis set is consistently defined using the mid-surface embedding, i.e. $\mathbf{a}_{\alpha}\equiv\partial_{\alpha}\mathbf{S}$ and $\mathbf{a}^{3}=\mathbf{a}_{3}=\mathbf{a}_{1}\times\mathbf{a}_{2}/|\mathbf{a}_{1}\times\mathbf{a}_{2}|$. We can relate basis vectors $\mathbf{g}_{i}\equiv\partial_{i}\mathbf{Y}$ on surfaces of constant $x_{3}$ to the corresponding basis vectors $\mathbf{a}_{i}$ on the mid surface, and vice-versa, via a translation tensor $\mathsf{T}:\mathbf{a}_{i}\mapsto\mathbf{g}_{i}$ which is given by $\mathsf{T}=\mathbf{g}_{i}\otimes\mathbf{a}^{i}$ (and $\mathsf{T}^{-1}=\mathbf{a}_{i}\otimes\mathbf{g}^{i}$) Pietraszkiewicz (1980). Through the definition of the mid-surface metric tensor, $\mathsf{a}=a_{\alpha\beta}\mathbf{a}^{\alpha}\otimes\mathbf{a}^{\beta}=\mathbf{a}_{\alpha}\otimes\mathbf{a}^{\alpha}$, and curvature tensor, $\mathsf{b}=b_{\alpha\beta}\mathbf{a}^{\alpha}\otimes\mathbf{a}^{\beta}=-\partial_{\alpha}\mathbf{a}^{3}\otimes\mathbf{a}^{\alpha}$, we may write $\mathsf{T}=\mathsf{1}-x_{3}\mathsf{b}$ ($\mathsf{T}^{-1}=\mathsf{1}+x_{3}\mathsf{b}+\cdots$). From the above definition, we may also apply the translation tensor to the deformation gradient, which results in the relation $\mathsf{F}=\mathsf{T}\mathsf{F}^{\mbox{\tiny(0)}}$. The above definitions allow us to write both strain tensors, in terms of their components, as follows:

[TABLE]

where

[TABLE]

with $c_{\alpha\beta}=\partial_{\alpha}\mathbf{a}^{3}\cdot\partial_{\beta}\mathbf{a}^{3}$ are the components of the third fundamental form. Also,

[TABLE]

where higher-order terms in $x_{3}$ have been neglected, and

[TABLE]

From the above expressions, it is clear that $\boldsymbol{\varepsilon}^{*}=\mathsf{F}^{\mbox{\tiny(0)}\mathrm{T}}\boldsymbol{\varepsilon}\mathsf{F}^{\mbox{\tiny(0)}}$ is satisfied. Despite differences in the Eulerian and Lagrangian formulations Pietraszkiewicz (1980), Eqs. (4) and (6) show that these two approaches yield identical tensorial components, in terms of the components of the fundamental forms in the deformed mid-surface do Carmo (1976), expressed with different bases. In practical terms, the mixed contribution $\mathbf{n}\cdot\boldsymbol{\varepsilon}^{*}\cdot\mathbf{n}-\bar{\mathbf{n}}\cdot\boldsymbol{\varepsilon}\cdot\bar{\mathbf{n}}$ results in contractions with the same components for the strain measure. Therefore, if we define $\varepsilon_{\alpha\beta}\equiv\left(a_{\alpha\beta}-\delta_{\alpha\beta}\right)/2-x_{3}b_{\alpha\beta}+x_{3}^{2}c_{\alpha\beta}/2$, the coupling terms become $\left(n^{\alpha}n^{\beta}-\bar{n}^{\alpha}\bar{n}^{\beta}\right)\varepsilon_{\alpha\beta}$.

Note that the curvature with a raised index is defined with the inverse of the deformed metric, $b^{\beta}_{\alpha}\equiv\left(\mathsf{a}^{-1}\right)^{\beta\delta}b_{\delta\alpha}$, which satisfies the property $\left(\mathsf{a}^{-1}\right)^{\alpha\gamma}a_{\gamma\beta}=\delta^{\alpha}_{\beta}$. However, by definition, we may write $a_{\alpha\beta}=\delta_{\alpha\beta}+2\varepsilon_{\alpha\beta}$, which through expanding in orders of the stain results in $\left(\mathsf{a}^{-1}\right)^{\alpha\beta}=\delta_{\alpha\beta}+\mathcal{O}\left(|\varepsilon|\right)$. Since the energy is already second-order in strain, we here assume the approximation that the Cartesian metric components $\delta_{\alpha\beta}$ raises and lowers the indices for the contractions in the energy, thus allowing us to use either covariant or contravariant notation, knowing that their difference only yields higher order contributions.

From the definition of the components of the stress tensor as

[TABLE]

the condition (i) reads $\sigma_{33}=\sigma_{23}=\sigma_{13}=0$ Landau and Lifshitz (1986). This step automatically implies that the third component of the three-dimensional strain tensor follows the constraint given by $\varepsilon_{33}=-\varepsilon_{\alpha\alpha}\lambda/(\lambda+2\mu)$, which allows us to rewrite the energy given in Eq. (2) entirely in terms of components on the surfaces of constant $x_{3}$.

The final assumption (iii) tells us that the director field only has components parallel to the surfaces of constant ${x_{3}}$. We can therefore uniquely express the field $\mathbf{n}$ in terms of its projection onto the mid-surface, here defined by $n^{\mbox{\tiny(0)}}_{\alpha}$, via the relationship

[TABLE]

which gives $n_{\alpha}=T_{\alpha}^{\beta}n^{\mbox{\tiny(0)}}_{\beta}$.

The dependence on $x_{3}$ coming from the nematic contribution to the energy is explicitly written though the following expression,

[TABLE]

where Eq. (9) has been used, $B_{\alpha\beta}\equiv n^{\mbox{\tiny(0)}}_{\alpha}n^{\mbox{\tiny(0)}}_{\delta}b_{\delta\beta}+b_{\alpha\delta}n^{\mbox{\tiny(0)}}_{\delta}n^{\mbox{\tiny(0)}}_{\beta}$, and $C_{\alpha\beta}\equiv 2b_{\alpha\delta}n^{\mbox{\tiny(0)}}_{\delta}n^{\mbox{\tiny(0)}}_{\gamma}b_{\gamma\beta}$ (the approximation $b^{\beta}_{\alpha}\approx b_{\alpha\beta}$ has been employed in these definitions). Substituting these expansions into Eq. (2) yields an energy density with an explicit functional dependence on ${x_{3}}$. Therefore, we may eliminate this degree of freedom by integrating along the body’s thickness. The non-zero contributions to the two-dimensional energy come only from even powers of ${x_{3}}$, because of the choice of symmetric limits of integration, $[-t/2,t/2]$. In this derivation, the terms proportional to $t^{2}\varepsilon^{\mbox{\tiny(0)}}_{\alpha\beta}c_{\alpha\beta}$ and $t^{2}\varepsilon^{\mbox{\tiny(0)}}_{\alpha\beta}C_{\alpha\beta}$ shall be neglected, as they are $\mathcal{O}\left(t^{2}\zeta^{4}/L^{6}\right)$ or smaller, where $\zeta$ is the scale of the typical out-of-plane deflection of material points in the $x_{3}$ direction and $L$ is a lateral length-scale of the plate. We write the two-dimensional energy in the usual thin-plate sense, ${E}=(1/2)\int\mathrm{d}\mathcal{A}^{\mbox{\tiny(0)}}\left[t\,\mathcal{E}_{\mathrm{s}}+(t^{3}/12)\mathcal{E}_{\mathrm{b}}\right]$, where the modified energy densities are written as follows:

[TABLE]

where $\Delta{Q}^{\mbox{\tiny(0)}}_{\alpha\beta}\equiv n^{\mbox{\tiny(0)}}_{\alpha}n^{\mbox{\tiny(0)}}_{\beta}-\bar{n}^{\mbox{\tiny(0)}}_{\alpha}\bar{n}^{\mbox{\tiny(0)}}_{\beta}$ and $\varepsilon^{\mbox{\tiny(0)}}_{\alpha\beta}\equiv\left(a_{\alpha\beta}-\delta_{\alpha\beta}\right)/2$. This energy is expressed entirely in terms of quantities defined on the mid-surface. Since the bulk modulus is typically much larger than the shear modulus in rubbery materials, we have assumed that the elastomer is nearly incompressible, i.e. $\lambda/(\lambda+2\mu)\rightarrow 1$ (or $\varepsilon_{33}=-\varepsilon_{\alpha\alpha}$). The first terms in Eqs. (11a)-(11b) are the usual stretching and bending, respectively, arising in isotropic plate theory. The last two terms in Eq. (11a) are readily seen to be inherited from the three-dimensional form of the bulk energy (2), while the last two in Eq. (11b) couple bending to the orientational order.

It is noteworthy that in the limit $h\rightarrow\infty$, a limit associated with nematic glasses Biggins et al. (2008); Cirak et al. (2014), our energy recovers the strong coupling between nematic defects and elastomer curvature Warner et al. (2010b); Modes and Warner (2011). In particular, for the energy to be bounded in this limit, we must have $B^{2}_{\alpha\beta}+C_{\alpha\beta}\Delta Q^{(0)}_{\alpha\beta}=\mathcal{O}(1/h)$, which implies that products of the director field and second fundamental form components are small at every point. Because $C_{\alpha\beta}\equiv 2b_{\alpha\delta}n^{\mbox{\tiny(0)}}_{\delta}n^{\mbox{\tiny(0)}}_{\gamma}b_{\gamma\beta}$, in a defect-free texture where the nematic director $n^{\mbox{\tiny(0)}}$ never vanishes, the only way for the energy to be bounded is to have the second fundamental form $b_{\alpha\delta}$ vanish. This implies that there is no curvature of the sheet for large $h$, which has large implications for wrinkling (since wrinkles induce curvature in the sheet), which we discuss in Section V.

Next, we use our derived plate theory to address the particular example of a thin NLCE plate bonded to an elastic or fluid substrate and placed under compression. Following previous works, we first derive the in-plane stress in the absence of buckling. For an isotropic elastomer, this step is crucial to understand how the stress and therefore the stretching energy scales in a near threshold regime. For NLCE, it is in addition necessary to determine the director orientation, which affects both the in-plane stress and the effective bending modulus, and therefore “tunes” the wavelength and critical compression for wrinkles. Having derived these quantities, we then consider a standard ansatz for the buckled shape, and solve for the amplitude, wavelength, and critical stress for buckling in both the cases of a fluid and elastic substrate.

IV In-plane stress and director orientation

First, we describe the state of stress in a flat configuration. We now adopt a global frame $(\hat{x},\hat{y},\hat{z})$ with $z$ being the thin direction, see Fig. 1. The thin elastomer is of infinite width $0<y<\infty$ and finite length $0<x<L$, with initial director $\bar{\mathbf{n}}^{(0)}=\hat{\mathbf{x}}$ (we perform the analogous finite-width calculation in the Appendix). The plate is subjected to a compressive strain $\varepsilon^{\mbox{\tiny(0)}}_{xx}=-\gamma$. We will assume the deformations are small and that the out-of-plane deflection $\zeta$ is identically zero. Therefore, we write $\varepsilon^{\mbox{\tiny(0)}}_{\alpha\beta}=(\partial_{\alpha}u_{\beta}+\partial_{\beta}u_{\alpha})/2$. We will often plot quantities against strains up to $25\%$, but the theory is expected to be numerically accurate only for small strains. To determine the base state of stress, we must solve the first Föppl-von Kármán equation $\nabla\cdot\boldsymbol{\sigma}^{\mbox{\tiny(0)}}=0$, where the stress is given by $\sigma^{\mbox{\tiny(0)}}_{\alpha\beta}=\partial\mathcal{E}_{\mathrm{s}}/\partial\varepsilon^{\mbox{\tiny(0)}}_{\alpha\beta}$ using the stretching energy density in Eq. (11a).

If $\phi$ is the angle between $\mathbf{n}^{\mbox{\tiny(0)}}$ and $\bar{\mathbf{n}}^{\mbox{\tiny(0)}}$, then these solutions are given by

[TABLE]

Note that due to the use of the rotation from the deformed space when defining $\mathbf{n}^{\mbox{\tiny(0)}}$, the angle $\phi$ may be different from what is measured in the laboratory frame. To convert to the angle measured in the lab frame, it is necessary to apply a local rotation obtained by taking the polar decomposition of the deformation tensor, which includes local metric information Mbanga (2012). However, the observables of interest in this work are macroscopic geometric quantities, such as the wavelength of wrinkles in the deformed state. Therefore, we do not convert from the measure $\phi$ to the director orientation measured in the laboratory frame. This orientation $\phi$ is determined by the balance equations $\partial\mathcal{E}_{\mathrm{s}}/\partial\phi=0$, which yields

[TABLE]

We suppose that due to the free boundary, $\sigma^{\mbox{\tiny(0)}}_{yy}=\sigma^{\mbox{\tiny(0)}}_{xy}=0$. Since the $\varepsilon^{\mbox{\tiny(0)}}_{xx}$ component of the strain is simply $-\gamma$, the imposed compression, equations (12b–12c) can be solved to give the remaining components of strain,

[TABLE]

Substituting these expressions into the balance equation for the nematic director, Eq. (13), allows us to calculate the director angle in terms of compressive strain and parameters of the NLCE. This equation admits the trivial solutions $\phi=0$, $\phi=\pi/2$, corresponding respectively to the states $\mathbf{n}^{\mbox{\tiny(0)}}=\hat{\mathbf{x}}$ and $\mathbf{n}^{\mbox{\tiny(0)}}=\hat{\mathbf{y}}$, but also admits a third intermediate solution given by

[TABLE]

By considering the regions of existence ($-1\leq\mathrm{arg}(\phi_{\mathrm{c}})\leq 1$) for the three critical points $\phi=\left(0,\phi_{\mathrm{c}},\pi/2\right)$, as well as their stability conditions,

[TABLE]

we are able to construct a phase diagram of the base nematic state. The different phases (unrotated nematic director, rotating nematic director leading to an elastic soft mode, fully-rotated nematic director) are separated by two critical strains:

[TABLE]

where for $\gamma<\gamma_{1}$ the director is unrotated, and for $\gamma>\gamma_{2}$ the director is fully-rotated. We show examples from this phase diagram in Fig. 2. Depending on the values of $(\alpha/\mu,h/\mu)$ at which the sheet was prepared, there are a few possible distinct behaviors of the director under small amounts of applied compression: no rotation, rotation only after a moderate critical strain, and an immediate jump in rotation under any amount of applied compression. In all cases, the director angle increases monotonically under the applied compression $\gamma$ and decreases monotonically in $h/\mu$. The behavior as a function of $\alpha/\mu$ is more complicated. For low values of $h/\mu$ (roughly $h/\mu<0.1$), there is a phenomenon of re-entrance as $\alpha/\mu$ varies, meaning that $\phi$ is not a monotonic function of $\alpha/\mu$, and holding $(\gamma,h/\mu)$ fixed whilst increasing $\alpha/\mu$ can take the director out of the rotating phase into the fully rotated phase, and then back into the rotating phase. The mechanism of this re-entrance is remarkably similar to that previously described for the nematic-smectic A re-entrant phase transition Pershan and Prost (1979), insofar as Eq. (16) expanded in a Taylor series gives a fourth-order term in $\phi$ with coefficient $47\alpha^{2}+8h\mu-12\alpha\gamma\mu$, which for particular fixed values of $(h,\gamma,\mu)$ (the analog of temperature in the nematic-smectic A) can change signs as a function of $\alpha$ (the analog of pressure), which leads to re-entrance as the transition changes from second-order to first-order.

When $\phi=\phi_{c}$, that is, the nematogen is neither unrotated nor fully rotated, we see the presence of a soft mode. This can be seen most easily by substituting Eq. (15) into Eq. (12a), yielding

[TABLE]

which shows that the in-plane stress is effectively constant until the nematic finishes its rotation, even if the amount of compressive strain is increased. Additionally, this stress is not necessarily negative, an unusual feature in comparison to a classical elastomer. This gives a criterion for the suppression of wrinkling of a compressed plate which is completely independent of any substrate, because there can be no buckling if the stress is positive. In fact, such a plate would not be expected to buckle into a single arch in the absence of any substrate. We call this condition in which the plate would remain flat “strong suppression” of wrinkles: plates prepared with certain values of $(\alpha,h)$ would begin the soft mode (and also be in the regime where the soft mode leads to $\sigma^{\mbox{\tiny(0)}}_{xx}>0$) under any amount of compressive strain, and would not complete the soft mode (i.e., rotate a total of $\pi/2$ from the base state) until a certain amount of strain was imposed, after which time buckling may be possible. Some examples are plotted in Fig. 3.

V Out-of-plane displacement, critical buckling strain, and scaling of wrinkles

We now assume an out-of-plane displacement which depends only on the in-plane coordinate $x$, $\zeta(x)=A\cos(qx)$ (we perform an analogous calculation for a plate of finite width $W$ and out-of-plane displacement $\zeta(x,y)$ in the Appendix, in which we demonstrate that identical results hold in comparison to this infinite, one-dimensional problem). Our particular choice of geometry allows us to simplify the bending energy (11b). The bending energy per unit wavelength (and per unit width, so effectively an energy per unit area) reads

[TABLE]

where the renormalized Young’s modulus is given by

[TABLE]

Following a standard calculation Cerda and Mahadevan (2003); Davidovitch et al. (2011); Huang et al. (2005); Audoly and Boudaoud (2008), we write the stretching energy density as

[TABLE]

where $\sigma^{\mbox{\tiny(0)}}_{xx}$ is given by Eq. (12a), and can be written in the form $\sigma^{\mbox{\tiny(0)}}_{xx}=-3\mu\gamma+(3/2)\alpha\sin^{2}\phi$. While there are other contributions to the stretching energy, they are independent of the wavelength of wrinkling and will therefore not contribute to the variational problem which gives this wavelength. For instance, following Huang et al. (2005), the total strain perpendicular to the wrinkles can be written as $\varepsilon^{\mbox{\tiny(0)}}_{xx}=-\gamma+\frac{1}{4}q^{2}A^{2}$, so that by Eq. (12a) the addition to the stress is simply $\mu q^{2}A^{2}$. The stretching energy per unit area (keeping only terms proportional to $q$) is therefore

[TABLE]

Finally, we assume that the thin NLCE film is bonded to an underlying foundation Landau and Lifshitz (1986); Cerda and Mahadevan (2003), which can be either fluid or elastic. Because the energy of the foundation depends slightly differently on the wrinkling in each case—the elastic foundation energy depends both on the amplitude and the wavelength of wrinkles, whereas the fluid foundation energy depends only on the amplitude—we will treat each case separately below.

Elastic substrate

Following previous derivations of the energetics for an elastic sheet atop an infinitely deep substrate of Young’s modulus $K$ Chen and Hutchinson (2004); Huang et al. (2005), the wrinkling induces a normal stress $\bar{\sigma}_{zz}=Kq\zeta$, so that the energy per area stored in the foundation can be written as

[TABLE]

and the total energy in the system (again, neglecting terms not proportional to $q$) is, collecting in terms of powers of $q$,

[TABLE]

When viewed as a function of the wrinkling amplitude $A$, it is clear that wrinkling can only occur when the term inside the parentheses is negative, or else the optimal amplitude will always be zero. It is convenient to write this term as $t\sigma^{\mbox{\tiny(0)}}_{xx}+\mathcal{B}(\phi)t\xi_{\mathrm{e}}$, where

[TABLE]

When the thin plate is wrinkled, the energy is then minimized at an amplitude of

[TABLE]

Therefore the plate can only wrinkle if $|\sigma^{\mbox{\tiny(0)}}_{xx}|>\mathcal{B}(\phi)\xi_{\mathrm{e}}$, and the stress must also be negative (compressive), which is not always the case, as seen in Sec. IV. Substituting Eq. (26) into Eq. (24) and minimizing the energy gives the wavelength $\lambda_{\mathrm{e}}$ of wrinkling on an elastic substrate,

[TABLE]

where $\ell_{\mathrm{e}}\equiv t[\mu/(6K)]^{1/3}$ is the natural length emerging in this problem. Eq. (27) allows us to write the critical stress for wrinkling in the case where the stress remains compressive,

[TABLE]

On their face, Eqs. (27) and (28) are identical to those found in the literature for a thin, non-nematic elastomer plate of infinite extent across an elastic foundation of infinite depth Chen and Hutchinson (2004); Huang et al. (2005); Audoly and Boudaoud (2008). However, the renormalized Young’s modulus $\mathcal{B}(\phi)$ as well as the in-plane stress $\sigma^{\mbox{\tiny(0)}}_{xx}$ now carry information from the embedded nematic phase — therefore, tuning parameters relevant to the nematic can both change the wavelength of wrinkles at a fixed amount of compression, $\gamma$, and can also change the threshold for buckling.

To demonstrate the interesting ways in which information from the nematic phase can affect the critical buckling strain of the entire sheet, we have plotted information about four representative elastomer plates in Fig. 4. These four plates were chosen because they represent the four essential behaviors of the director angle which are seen in Fig. 2 if the applied compressive strain never exceeds $\gamma=0.25$: no rotation (high $h$, low $\alpha$); a large soft mode characterized by rotation across the whole range of $\gamma$ (low $h$, high $\alpha$); rotation which begins only after a moderate applied strain (high $h$, high $\alpha$); and a plate which both begins and ends its rotation within the range of $\gamma$ used (low $\alpha$, low $h$). The curved lines in Fig. 4(a) delimit the regions where the plate will never rotate for $\gamma<0.25$ (blue-purple) and where the soft mode induces positive stress (purple-orange). The flat line (blue-red) roughly delimits plates which begin and end rotation before $\gamma=0.25$. During the rotation of the director, the critical buckling strain becomes almost completely independent of the Young’s modulus of the underlying substrate, which acts to suppress wrinkling. For instance, the red plate bonded to a substrate with $K/\mu=0.01$ has a critical buckling strain of roughly $\gamma=0.22$, whereas a blue plate atop the same substrate has a critical buckling strain of roughly $\gamma=0.03$.

Fluid substrate

Atop a fluid substrate, the increase in the energy density is due to hydrostatic pressure, which has an energetic cost given by

[TABLE]

where $\rho$ is the density of the fluid and $g$ is the gravitational field. The analogue of the quadratic term in Eq. (25) is therefore

[TABLE]

Eq. (26) and the critical stress for wrinkling are therefore unchanged (save for the substitution of $\xi_{\mathrm{f}}$ for $\xi_{\mathrm{e}}$); however, because of the extra factor of $q^{-1}$ in $\xi_{\mathrm{f}}$, the critical wavelength is slightly different:

[TABLE]

where $\ell_{\mathrm{f}}\equiv[t^{3}\mu/(12\rho g)]^{1/4}$ is the natural length associated with the cost of deforming the fluid substrate. This expression is again identical to the equivalent expression for a regular, isotropic elastomer, aside from the fact that $\mathcal{B}(\phi)$ carries information from the nematic orientation, which is dictated by the in-plane problem. Therefore, the wavelength of the wrinkles or even their existence is highly dependent upon the base state of stress, much unlike the case of an isotropic elastomer. The critical stress, again in the case where the stress remains compressive, is

[TABLE]

much of which can still be inferred from Fig. 3, though now the critical stress is dependent on the thickness $t$ of the elastomer plate.

We note that in addition to the region described in Fig. 3 where wrinkling is forbidden due to a change from compressive to extensile stress, wrinkling is also suppressed for large $h$ due to effective stiffening of the sheet. We show this in Fig. 5 for the case of an elastic substrate, although a nearly-identical result holds for the fluid substrate as well. Because this limit corresponds to a nematic glass Cirak et al. (2014); Lubensky et al. (2002), the mechanism for instability suppression in this region is the strong coupling between curvature and defects in the nematic texture Modes and Warner (2011). In Eq. (11b), we see that for $h\gg\alpha$ the bending energy is very large unless $B_{\alpha\beta}^{2}$ and $C_{\alpha\beta}\Delta{Q}^{\mbox{\tiny(0)}}_{\alpha\beta}$ are very small. Since our textures are free of defects, the presence of any curvature is highly penalized for arbitrarily-increasing $h$.

VI Summary

In this work we derived a Föppl-von Kármán type of plate theory for a thin NLCE which was cross-linked deep in the nematic phase. This was accomplished following the standard techniques of dimensional reduction for thin elastic bodies, employing the additional assumption that the nematic director was tangent to the mid-plane before and after the deformation. In the limit of small and large $h$, the coupling parameter for the “memory” or “anchoring” term corresponding to the fact that the elastomer was cross-linked in the nematic phase, our model corresponds to other objects studied in the literature: respectively, an elastomer cross-linked in the isotropic phase ($h\rightarrow 0$) Lubensky et al. (2002); Warner and Terentjev (2003) or a nematic glass ($h\rightarrow\infty$) Biggins et al. (2008); Cirak et al. (2014); Warner et al. (2010b); Modes and Warner (2011). When both this coupling and the director-strain coupling parameter vanish, we recover the equations for the isotropic plate.

The model is fairly simple and solutions to certain geometries and boundary conditions can be found analytically. As a first step in this direction, we calculated the wavelength of wrinkles for a compressed NLCE plate atop a fluid or elastic foundation. The wave number is found to be non-monotonic in the compressive strain $\gamma$ and highly dependent on the nematic-elastomer coupling parameters, as it is shown in Fig. 5. For certain sets of coupling parameters the plate does not wrinkle until a very high strain threshold is reached. We therefore believe it will be possible to design NLCE which are not subject to the elastic instabilities of classical rubbers until a prescribed strain threshold is met.

In summary, NLCE plates present a new venue to study nontrivial pattern formation in thin sheet elasticity and also present an example of robust control of mechanical stability through the coupling of geometry and microstructure. We hope that our simplified model will be useful in further study of liquid crystal elastomers.

Acknowledgements.

We are grateful to Thomas Powers for his guidance and the many discussions had concerning this project. We would also like to thank Badel Mbanga, Timothy Atherton, James Hanna, and Jonathan Selinger for helpful conversations.

Appendix A Wrinkling of a plate with finite width

When we consider a plate of finite width, the stress is the same as found in (12a)-(12c), but we allow for out-of-plane deformations which depend on both spatial coordinates, $\zeta(x,y)$. For the sake of convenience, we consider here the case of a sheet which has one small lateral dimension: $0<y<W$, where $W\ll L$. Following previous work on wrinkling in a similar geometry Cerda and Mahadevan (2003) we reason that because the out-of-plane displacement varies much less in $y$ than in $x$ the respective terms in the energy are also smaller: $\left(\partial_{y}^{2}\zeta\right)^{2}<\left(\partial_{x}\partial_{y}\zeta\right)\left(\partial_{y}^{2}\zeta\right)\ll\left(\partial_{x}^{2}\zeta\right)^{2}$; we will consequentially choose to keep only the largest terms. The total bending energy reads

[TABLE]

where the respective bending stiffnesses are given by

[TABLE]

Note that $\mathcal{B}^{xxxx}$ corresponds to $\mathcal{B}(\phi)$ in (19). Because the stress and foundation energy remain unchanged, we have now the balance equation

[TABLE]

We assume that the out-of-plane displacement takes the form $\zeta=e^{(iqx)}Y(y)$. That is, the wrinkles are cylindrical as before, but the function $Y(y)$ will serve to implement the free boundary conditions at the edges. This gives a Sturm-Liouville equation for $Y(y)$:

[TABLE]

where

[TABLE]

Upon substitution of the ansatz $Y(y)=Ae^{\kappa y}$, the characteristic equation is obtained:

[TABLE]

where $\omega^{2}=A_{1}/A_{3}$ and $2D=iA_{2}/A_{3}$. This equation has the roots $\kappa=-D\pm\sqrt{D^{2}-\omega^{2}}$. As a consequence, $Y(y)$ is either an exponential function or a sinusoid; the former must be disbarred by symmetry, so that $Y(y)=0$ whenever $\kappa\in\mathbb{R}$. Because $A_{2},\ A_{3}$ are real, $D$ satisfies $D^{2}<0$. Therefore, wrinkling solutions are only valid when $\omega^{2}\geq D^{2}$.

The scaling for the finite two-dimensional sheet is the same as in the one-dimensional cylindrical approximation. To show this, we consider the sinusoidal solution for $Y(y)$:

[TABLE]

where $\kappa^{+},\ \kappa^{-}$ indicates which sign of the $\pm$ in the root is being considered. Minimizing the bending energy allows us to simplify this expression to $Y_{n}(y)=\sin\left(\pi y/W\right)$, because additional modes increase the bending energy prohibitively. Integrating the total energy and keeping only terms algebraic in $q$ yields

[TABLE]

Because the new contributions to the bending energy which are proportional to other bending moduli do not scale with $q$ and the foundation and stretching energies remain unaltered, the scalings given in (31), (27) are still applicable for the fully finite 2D sheet with saddle-shaped wrinkles.