Foldable structures made of hydrogel bilayers

1.   Introduction

The study of shape programming and morphing of surfaces is receiving considerable attention, especially in the field of active (or smart) materials, i.e. materials that deform in response to non-mechanical stimuli. Of particular interest is the problem of exploiting material heterogeneities to induce complex shape changes—for instance, to produce curved configurations from an initially flat state ^[2,3,4,5,6]. In many natural systems, such as plants ^[7,8], shape control is usually accomplished by growth, remodelling, or swelling, in response to simple external stimuli (e.g. a uniform change in the ambient temperature or humidity). To mimic such behaviors, synthetic soft active materials appear as suitable candidates. Hydrogels are examples of active materials where spontaneous deformations are induced by swelling due to the absorption of a liquid. In the form of bilayers, they can be employed to produce curved shapes ^[9].

Self-folding is a widespread phenomenon that occurs in natural systems, such as in the opening and closing of flowers. The understanding of such mechanisms offers interesting challenges and opportunities, not only from the point of view of biology, but also from the point of view of applications, for instance, in the fabrication of actuation systems and origami-like structures ^[10]. Self-folding results from spatially heterogeneous deformations that are (in most cases) either induced by a spatial variation of the external stimulus or by modulations of the material properties imprinted during the fabrication process. We refer to ^{[11,12,13,14,15,16]} for some examples of the approaches to folding and the numerical models that might be used in the study of folding mechanisms.

In this paper, we discuss self-folding of thin sheets by using hydrogel bilayers, which act as hinges connecting flat faces. Folding is actuated by the heterogeneous swelling due to different stiffness (i.e. cross-linking density) of the polymer network across two layers. More precisely, we will focus on the case where the cross-linking density in the top and the bottom layers is a perturbation (of magnitude comparable to the total thickness of the plate, which is a small parameter) of an average constant density $\overline N$ , see formula (2.1). We show that this structure allows to endow a thin gel sheet with a controlled curvature localized at the hinges (see Figure 1), which can be realized upon swelling at low energy cost. Furthermore, such a curvature can be expressed as a function of the material parameters of the layers. Specifically, we demonstrate the feasibility of the proposed folding mechanism with two examples, corresponding to specific patterns of flat faces and hinges (see Figure 2).

Our approach to the analysis of hydrogel-based folding sheets exploits the two-dimensional (nonlinear) plate model presented in Subsection 2.1, which is obtained by applying the rigorously derived theory ^[1] for heterogeneous thin elastic plates based on $\Gamma$ -convergence arguments. The value of employing such a rigorously derived theory is twofold: first, we avoid inconsistencies often present in ad hoc, formally-deduced models by performing a mathematically rigorous limit; second, the material and geometric parameters of the reduced 2D model emerge from (i.e. are computed from) those of the 3D one, which are the only ones that are accessible from the point of view of the sample preparation in a laboratory, hence ensuring that the minimizers of the 2D model faithfully reproduce the true behavior of sufficiently thin (yet, three-dimensional) bilayers. In Section 3 we study the pointwise energy minimizers of the plate model, which describe the configurations of a plate with bilayer-like hinges and thus provide a theoretical justification for the effectiveness of the folding mechanism. In Section 4, with reference to an appropriate approximate variant of our Flory-Rehner-type model, we provide the explicit relation between the target curvature and the material parameters, see (4.25). Such relation is then used in the design of folding bilayer sheets that morph into cubes or pyramids. We end the paper with an appendix, where we provide a detailed analysis of Flory-Rehner-type energy densities and show how the dimension-reduction framework studied in ^[1] applies to the model considered here.

2.   Self-folding using bilayer gels

Let us call $\Omega_{h_0} = \omega\times (-h_0/2, h_0/2)$ the reference configuration of our foldable structure, where the reference thickness $h_0$ is much smaller than the in-plane dimensions. The planar domain $\omega$ , which represents the mid-plane of the plate, is hereafter assumed to be a union $\omega = \bigcup_{i = 1}^n\omega_i$ of polygons. We denote by $x = (x_1, x_2, x_3) = (x', x_3)$ an arbitrary point in $\Omega_{h_0}$ . We say that the interfaces which delimit the polygons form a pattern on $\omega$ .

To model our system, we consider a family $\{W_h\}$ of 3D free-energy densities, from which we can derive the corresponding 2D Kirchhoff bending model (i.e. the energies scale as $h^3$ ) in the vanishing thickness limit, by using the rigorous derivation recently provided in ^[1]. At the same time, $W_{h_0}$ represents the actual energy density corresponding to the finite (small) thickness system. A thin sheet with ideal small thickness $h$ is supposed to be a heterogeneous system where the heterogeneity, both lateral and vertical, is due to a patterned bilayer structure made of hydrogels. This structure is engineered through an $x$ -dependent density of cross-links in the polymer network, which in turn determines an $x$ -dependent density $N_h(x)$ of polymer chains. More precisely, we suppose that

for some constants $\overline N>0$ and $M_i\geq0$ , $i = 1, \ldots, n$ . Note that the average density along the thickness $\int_{-1/2}^{1/2}N_h(x', x_3){\rm d}x_3$ is supposed to be constant in $x'$ and equal to $\overline N$ .

The system is a composite of flat faces (where $M_i = 0$ ) connected by hinges (where $M_i>0$ ). Starting from an initially dry state, folding will be accomplished by putting the gel in contact with a solvent. Upon swelling, the difference in $N_h$ in the two layers of the hinges will induce bending. The prototype of this actuation mechanism is presented in Figure 1, where the larger (resp. lower) number of dots in the hinge (red patch) corresponds to a higher (resp. lower) cross-linking density in each layer.

For the constitutive model of the hydrogel, we use an energy density of Flory-Rehner type ^[17], which takes into account both the elastic and the mixing energy contributions. Namely,

where $x\in\Omega_h$ and $F$ , which represents a deformation gradient, belongs to the class ${\mathbb{R}}_1^{3 \times 3}$ of all the $3\times 3$ matrices with $\det F\geq 1$ . The last term in the above expression is the energy of the solvent within the gel, whereas the mixing energy $W_{vol}^{\chi}$ is defined in $(1, +\infty)$ by

Moreover, the parameters appearing in (2.2) are:

● $k$ - Boltzmann's constant;

● $T$ - absolute temperature of the environment;

● $v$ - volume per solvent molecule;

● $\chi\in(0, 1/2]$ - dimensionless measure of the enthalpy of mixing;

● $\mu\leq0$ - chemical potential of the solvent molecules within the gel.

Note that the range $\mu\leq0$ corresponds to having a gel in contact with an external fluid that is either a vapor ($\mu < 0$ ) or a pure liquid in equilibrium with its own vapor ($\mu = 0$ ). Finally, the range $\chi\in(0, 1/2]$ , which corresponds to a "good" solvent, implies in particular that $W_{vol}^{\chi}$ is a strictly decreasing function fulfilling

As a consequence of the heterogeneity of $W_h$ , which is due to the heterogeneous density of polymer chains $N_h$ , the energy well will in turn be $x$ -dependent. More precisely, from the computations in Appendix A (see (A.35) and (A.38)–(A.39)) we have that there exist two constants

which depend on $k$ , on the (fixed) material parameters of the gel $v$ , $\chi$ , and average density $N$ , and on the environmental conditions ($\mu$ , $T$ ), such that

where, for every $x\in\Omega_h$ , $U_h(x)$ is the scalar multiple of the identity matrix given by

The map $x\mapsto U_h(x)$ is the spontaneous stretch field that characterizes the relaxed (zero-stress) state the system would try to attain spontaneously, at each point $x$ , through a deformation with gradient $U_h(x)$ . The two constant matrices $\alpha{\rm I}_3$ and $\frac{h}{h_0}\beta\, M_i{\rm I}_3$ have been written separately in definition (2.6) to stress their different contributions as stretching and bending terms, respectively (see equations (2.8)–(2.9) defining the reduced model, where $\alpha$ enters the definition of the class on which such model is constrained and $\beta$ in that of the target curvature).

For every $0 < h\ll 1$ , the total energy of the system associated with a deformation $v_h\in {\rm{W}}^{1, 2}(\Omega_h, {\mathbb{R}}^3)$ in the absence of external loads is given by the quantity

The choice of the functional space ${\rm{W}}^{1, 2}(\Omega_h, {\mathbb{R}}^3)$ is natural since a finite energy $\mathcal E_h(v_h)$ implies that $v_h\in {\rm{W}}^{1, 2}(\Omega, {\mathbb{R}}^3)$ (see the lower bound (A.49) in Appendix A).

2.1.   Two-dimensional model of a thin gel sheet

The 3D model presented in Section 2 falls within the general framework studied in ^[1] for the derivation of the corresponding 2D Kirchhoff model, via $\Gamma$ -convergence arguments. This is proved explicitly in Appendix A. A direct consequence of this fact is that any low-energy sequence $\{v_h\}_h$ of the scaled 3D energies $\frac{1}{h^3}\mathcal E_h$ converges (in the ${\rm{W}}^{1, 2}$ topology) to a deformation $y$ that minimizes the 2D energy defined, up to constant terms, as

on the class ${\rm{W}}^{2, 2}_{\alpha, \rm iso}(\omega)$ of all $u\in {\rm{W}}^{2, 2}(\omega, {\mathbb{R}}^3)$ such that $(\nabla u)^{{\rm{T}}}\nabla u = \alpha^2\, {\rm I}_2$ . We recall that, up to an $\alpha$ -dilation, any $u\in{\rm{W}}^{2, 2}_{\alpha, \rm iso}(\omega)$ is a $C^1$ -isometric immersion of $\omega$ into ${\mathbb{R}}^3$ , which admits jumps of curvature.

In formula (2.8), ${\rm A}_u$ is the second fundamental form associated with $u$ and expressed in terms of the "undeformed" in-plane variable $x'$ . Namely, with the following choice $\nu(x')$ of the unit normal to the surface $u(\omega)$ at $u(x')$ ,

The target curvature tensor $\overline A$ is a piecewise constant map, which is defined on each subdomain $\omega_i$ of $\omega$ by

Moreover, the fourth-order tensor $\mathbb C$ is given by

where $G$ and $\Lambda(\alpha)$ are non-dimensional positive Lamé constants, which can be expressed in terms of the material parameters as

The positivity of $\lambda(\alpha)$ is not self evident from the above formula and due to the fact that $\alpha$ depends in particular on $v$ , $\overline N$ , and $\chi$ (see formulas (A.35) and (A.47)–(A.48) with θ = 1 in Appendix A).

Another direct consequence of the $\Gamma$ -convergence result mentioned at the beginning of this subsection, which we now emphasize in terms of the physical quantities and the finite thickness $h_0$ (small with respect to the in-plane characteristic size of the plate), the following asymptotic approximate identity for the low-energy values $\mathcal E_{h_0}\big(v_{h_0}\big) = \inf\mathcal E_{h_0}+o(1)$ :

Here, the additional term "${\rm ad.t.}$ " is given by the formula

Hence, the minimizers of the 2D bending model (2.8) provide reliable estimates for the "almost minimal" values of the 3D energy given by (2.2) and (2.7). Also, the positivity of the additional term (due to (2.11)) tells us that the 2D energy cannot be minimized at the value zero: this fact originates in the incompatibility of the spontaneous strain distribution and hence in the presence of residual stresses in the 3D reference configuration.

3.   Minimal energy configurations: connection between shape parameters and gel properties

In this section, we study the minimizers of the limiting bending energy $\mathcal E_0$ given by (2.8). More precisely, we focus our attention on the class of pointwise minimizers, namely, on those deformations $y\in {\rm{W}}^{2, 2}_{\alpha, \rm iso}(\omega)$ whose second fundamental form ${\rm A}_y$ minimizes pointwise the integrand in (2.8). Importantly, this analysis employs a simple structure from the fabrication viewpoint and provides the theoretical basis for successful and robust actuation using the proposed folding mechanism.

Let us first operate an affine change of variable by associating to each $y\in {\rm{W}}^{2, 2}_{\alpha, \rm iso}(\omega)$ a corresponding isometry $\tilde y\in{\rm{W}}^{2, 2}_{1, \rm iso}(\alpha\omega)$ , via

The second fundamental forms associated with $y$ and $\tilde y$ are related via the formula

In turn, we get that

According to (2.9) and (3.16), and recalling that each isometry $\tilde y$ is such that $\det {\rm A}_{\tilde y} = 0$ a.e. in $\alpha\omega$ , our minimization problem reduces to finding those isometries $\tilde y$ whose second fundamental form fulfills

for every $i = 1, \ldots, n$ , where $\bar a_i$ is defined as in (2.9). A necessary and sufficient condition for (3.17) to hold is that ${\rm A}_{\tilde y}$ is actually equal to a constant $A_i\in \mathcal N_i$ on each $i$ -th subdomain of $\alpha\omega$ (for the necessity of such condition we refer the reader to ^[18]).

The solution of the (finite dimensional) minimization problem in (3.17), as done in [1,Lemma 3.1], yields the following explicit representation for the set $\mathcal N_i$ :

In other words, the second fundamental form ${\rm A}_{\tilde y}$ of an energy minimizing isometry $\tilde y$ , when restricted to the $i$ -th subdomain of $\alpha\omega$ (that we refer to as the $i$ -th patch) with $\kappa_i\neq 0$ , corresponds to the second fundamental form of a cylindrical surface whose non-zero principal curvature equals $\kappa_i/\alpha^2$ , while the associated principal curvature direction might be (a priori) given by any unit vector in ${\mathbb{R}}^2$ . However, the isometry constraint forces a precise choice of the principal curvature direction associated with the non-zero principal curvature $\kappa_i/\alpha^2$ – it must be orthogonal to the interfaces between $i$ -th patch and all the neighboring ones. Equivalently, a cylindrical surface with $\kappa_i\neq 0$ can be glued to all the neighboring cylindrical surface patches (or even plane patches) only if its rulings are parallel to the interface between them. It is thus clear that the existence of an isometry with the above described properties heavily depends on the compatibility between the pattern on $\omega$ and the target curvature, as proved in [1,Theorem 3.9].

The theoretical results that we have just discussed provide the foundations for a successful self-folding strategy in patterned, bilayer thin sheets made of hydrogels (or other active materials). In particular, the planar domain $\omega$ is patterned in such a way that the heterogeneity in swelling of the gel due to variations in the cross-linking density induces piecewise constant target curvatures. Many interesting and feasible patterns on $\omega$ along with the induced target curvatures at the hinges satisfy the previously mentioned compatibility property, thus guaranteeing the existence of a pointwise minimizer of the 2D bending energy. To be concrete, we now consider the pyramid-type and the cube-type domains, sketched in the Figure 2(A) and (B), respectively. Both of these two patterned domains consist of two different types of patches: hinges and flat faces.

Each hinge, denoted by $\omega_{\sf hinge}$ , has a nontrivial bilayer structure, characterized by

with $M>0$ . On each flat face, denoted by $\omega_{\sf flat}$ , the density of polymer chains $N_h$ is constantly equal to $\overline N$ , i.e. $M = 0$ . Note that both in the pyramid and the cube case the interfaces between an $\omega_{\sf hinge}$ and the (at most two) neighboring $\omega_{\sf flat}$ are mutually parallel (recall that the set of the nonzero principal curvature directions of a pointwise minimizer restricted to $\omega_{\sf hinge}$ may be any unit vector of ${\mathbb{R}}^2$ , see (3.18)). Hence the existence of a pointwise minimizer of $\mathcal E_0$ is guaranteed. More in detail, such minimizer deforms the mid-plane $\omega$ in the following way: it dilates the domain $\omega$ by the factor $\alpha$ and maps each $\alpha$ -dilated hinge into a cylindrical surface with radius $r = \alpha^2/|\kappa|$ , where

the cylinder's rulings are parallel to the interfaces which delimit the hinge from the flat faces. At the same time, each $\alpha$ -dilated flat face remains flat.

Looking at the expression (3.20) for the curvature $\kappa$ that arises on the hinges, it is clear that, maintaining the other physical constants fixed, an appropriate choice of the values $\overline N$ and $M$ of the bilayer structure at the hinges (recall that $\Lambda(\alpha)$ is $\overline N$ - dependent as well, see formula (2.11) and the subsequent paragraph) is needed to induce precise self-folding as in Figure 2(C) and (D). We address the problem of designing the structure of each bilayer to produce the desired shape upon folding in the following section.

4.   Bilayer design problem

In this section we address the "inverse" problem of finding the physical properties of a hydrogel that allow us to realize a given target shape upon self-folding. First, we introduce an approximate variant of the Flory-Rehner energies $W_h$ , for which the dependence of the target curvature $\overline A$ on the fixed physical parameters of the 3D system becomes explicit. Further, we restrict our attention to gels characterized by

● $\mu = 0$

● $v$ , $\overline N$ and $\chi$ satisfying the relation $\overline N < \frac{1-2\chi}{2v}.$

The first of these two conditions means that the gel is in contact with a pure liquid in equilibrium with its own vapor. The second condition is always satisfied whenever a hydrogel with a large swelling ratio is considered, as it occurs for the typical values $\overline Nv\sim 10^{-3}$ and $\chi\sim 10^{-1}$ .

To present the approximate model we work with, let us go back to the general setting of a thin gel sheet occupying the domain $\Omega_h = \omega\times (-h/2, h/2)$ , with $\omega = \bigcup_{i = 1}^n\omega_i$ being a general patchwork of polygons. First of all, in the case $\mu = 0$ the Flory-Rehner energy density (2.2) reduces to

Then, by Taylor's expansion, $W_{vol}^{\chi}(t) = \frac{1-2\chi}{2t}+o\left(\frac{1}{t}\right)+\chi-1,$ for every $t\gg 1$ . Thus, when deformation gradients with large determinants are considered, one may approximate $W_h$ (discarding the constant term $\chi-1$ ) by

This approximate model is adequate from the point of view of applications. It has been considered for the first time in ^[17], and afterward used, for instance, in ^[19,20].

The functional form (4.21) for the energy densities $\widehat W_h$ allows us to determine the associated spontaneous stretch distribution explicitly in terms of the fixed physical parameters of the model. Indeed, from Remark 1 in Appendix A, for every $x\in \omega_i\times (-h/2, h/2)$ , $i = 1, \ldots n$ , we have

where, for every $x\in\Omega_h$ , $\widehat U_h(x)$ is the scalar multiple of the identity matrix given by

and the constants $\alpha$ and $\beta$ are explicitly given by

In this case, the constant $\Lambda(\alpha)$ appearing in (2.11) can be explicitly derived, as well. Thus, by Remark 3 (in particular, formula (A.53)) in Appendix A we obtain that the tensor $\mathbb C$ in (2.8), within this approximate 2D model, reads

Finally, as a consequence of (4.24), the target curvature $\overline A$ is given on each subdomain $\omega_i$ of $\omega$ by

Summarizing, the total energy (2.8) in this case is

for every $u\in {\rm{W}}^{2, 2}_{\alpha, \rm iso}(\omega)$ . In particular, whenever a compatible pattern is considered, from the discussion at the end of Section 3 and from formula (3.18), it follows that the isometry $\tilde y$ associated with a pointwise minimizer $y$ of $\mathcal E_0$ via (3.14), maps each $i$ -th subdomain of $\alpha\omega$ into a cylindrical surface of radius

and with rulings parallel to the interface with each neighboring patch.

We are now ready to design a bilayer structure on each $\omega_{\sf hinge}$ (i.e. to determine the values of $\overline N$ and $M$ on such subdomain, see (3.19)), in order to induce precise self-folding of a pyramid and a cube (Figure 2(C) and (D)), at the minimum energy cost.

Pyramid    The target shape that we refer to as a (precisely) folded pyramid is characterized by the following set of parameters (see Figure 3):

● $\ell$ – initial length of the pyramid basis,

● $\phi\in(0, \pi/2)$ – vertex angle,

● $H$ – height of the pyramid,

● $\alpha>1$ – in-plane swelling factor.

By 'precisely folded pyramid' we mean that the four external vertices of the unfolded pyramid meet at a single point in the folded configuration. The above four parameters completely determine a (unique) 'precisely folded pyramid' shape, since from then one can derive

● $r$ – radius of curvature of the deformed ($\alpha$ -dilated) hinges,

● $\ell_1$ – initial hinge width,

● $\ell_2$ – initial height of the pyramid side,

via the formulas

By using the first formula in (4.24) and expression (4.27), we deduce that the correct material properties ($\overline N$ and $M$ ) to be imprinted on the hinges in order to accomplish the desired self-folding are determined in terms of the swelling factor $\alpha$ and the curvature radius $r$ by

Cube    We consider as the target shape the (precisely) folded cube characterized by the following parameters (see Figure 4):

● $\ell_1$ – initial width of the hinge

● $\alpha>1$ – swelling factor

● $\phi = \pi/2$ – hinge closing angle

The radius $r$ of a cylindrical surface representing a deformed hinge, which is needed to obtain the desired folded cube from a flat, cube-type pattern, must satisfy the relation $r = \frac{2\alpha\ell_1}{\pi}.$ As in the pyramid case, formulas (4.24) and (4.27) allow us to design the hinges in order to get a (precisely) folded cube. Namely, the material parameters that determine the cross-linking density in the bilayers that constitute the hinges are

As it is clear from the above formulas, the final "size" of a folded cube or pyramid is controlled by the parameter $\overline N$ through its relation with the swelling factor $\alpha$ . Specifically, in order to produce larger folded cubes or pyramids, which correspond to a bigger value of $\alpha$ , a smaller value of $\overline N$ is needed. The curvature of the hinges $1/r$ is controlled by the parameter $M$ , in such a way that, as one might expect, the value of $r$ increases as the value of $M$ decreases. This means that, on one side, flat plates (seen as cylindrical surfaces with radius $r = +\infty$ ) correspond to the isometric deformation of patches with no through-the-thickness variation of the cross linking density, i.e. $M = 0$ . On the other side, shapes with sharp folds (i.e. with $r\to 0$ ) require large variations $M$ in the cross-linking density along the bilayer that constitutes the hinge. In turn, since $r$ and $\ell_1$ are proportional, manufacturing a folded structure with a small hinge width $\ell_1$ requires a large value of $M$ .

5.   Conclusion

We have shown how patterned, thin hydrogel bilayers can be employed to realize self-folding structures. Specifically, folding can be achieved upon swelling by inserting hydrogel bilayers with varying cross-linking density across the thickness as hinges between flat faces, where the material is homogeneous.

Using rigorous dimension reduction techniques, we have obtained a two-dimensional plate model that describes the mechanics of such a system. Then, the study of the pointwise minimizers of the 2D model has suggested a successful strategy for the realization of foldable structures with programmable shape. In particular, folded structures arise as pointwise minimizers of the 2D energy whenever the pattern of hinges and flat faces satisfy the compatibility conditions determined in ^[1]. These conditions provide the theoretical basis for successful and robust actuation.

A major result of the dimension reduction procedure is to provide us with the connection between material properties (a feature of the 3D model) and the induced curvature upon swelling (a feature of the 2D model). We have then exploited such a connection as a tool for the design of the bilayer structure able to produce target folded shapes (specifically, cubes and pyramids).

The richness of emergent shapes found in nature that one would like to mimic using synthetic systems is, of course, not limited to configurations obtainable as pointwise minimizers. This strongly motivates further, more general, studies of the geometry of the minimal energy configurations of the 2D bending functional (2.8), as well as the derivation of new dimensionally reduced models allowing, for instance, for the presence of both bending and stretching contributions in the 2D energy.

A.   Appendix

This section is devoted to the analysis of the Flory-Rehner-type energy density functions $W_h$ introduced in (2.2). Letting all the physical parameters ($v, \chi, \mu, k$ and $T$ ) be fixed and letting the mid-plane of the plate $\omega$ be a union of polygons $\omega_i$ , we recall that

for every small thickness parameter $0 < h\ll 1$ . We refer to (2.3) for the definition of $W_{vol}^{\chi}$ . For later use it is convenient to introduce the following equivalent normalized expression for the density of polymer chains:

Note that for every $F\in\mathbb{R}^{3\times 3}_{+} : = \{F\in\mathbb{R}^{3\times 3}: {\rm det}F>0\}$ , by a standard algebraic inequality we have that

where $\lambda_{1}^{2}, \lambda_{2}^{2}, \lambda_{3}^{2}$ are the eigenvalues of the matrix $F^{{\rm{T}}}F$ , and where the equality holds if and only if $\lambda_{1}^{2} = \lambda_{2}^{2} = \lambda_{3}^{2}$ . In particular, we have that

where $p:(0, +\infty)\times [1, +\infty)\to \mathbb{R}$ is defined as

The above inequality, and the related rigidity characterizing the equality case, says in particular that $F_{\rm min}$ is a minimizer of $W_h(x, \cdot)$ if and only if $F_{\rm min} = t_{\rm min}^{1/3}R$ , for some $R\in{\rm{SO}}(3)$ , where $t_{\rm min}$ is a minimizer of $p\big(f_h(x), \cdot\big)$ . This argument is detailed in the proof of Lemma 3 below.

In view of the previous discussion and since $f_h(x)$ is uniformly close to $1$ for every $h$ sufficiently small, in what follows we analyze the minimizers of $p(\theta, \cdot)$ , for $\theta$ varying in some interval $\mathcal I$ containing $1$ . To do this, it is useful to recall the following sharp logarithmic estimate (which can be found in ^[21])

This can be checked just observing that the function defined as the left-hand side minus the right-hand side of (A.34) is null at zero and has a non-positive derivative. We also recall the following 1-dimensional real analysis result, which will be useful as well.

Lemma 1. Let $\mathcal J\subset \mathbb{R}$ be an open interval. Let $f\in C^{1}(\mathcal J)$ and assume that whenever $f'(t) = 0$ , for some $t\in \mathcal J$ , the second derivative of $f$ at $t$ exists and is such that $f^{"}(t)>0$ . Then the function $f$ has at most one stationary point, which is a global minimizer.

With the following lemma we show that $p(\theta, \cdot)$ has a unique global minimum in $[1, +\infty)$ .

Lemma 2. There exists an open interval $\mathcal I\ni 1$ and a unique smooth function $\varphi:\mathcal I\to (1, +\infty)$ such that

Proof. Set $\overline{p}: = p(1, \cdot)$ . First, note that the derivative of $\overline{p}$ reads as

Let $h(t): = \log\left(\frac{t}{t-1}\right)-\frac{1}{t}-\frac{1}{2t^{2}}$ for every $t\in (1, +\infty)$ . Note that $h(t)\geq 0$ for every $t\in (1, +\infty)$ and that $t^{2}h(t)\to 0$ as $t\to +\infty$ . It is straightforward to check that $\lim_{t\to 1}\overline{p}'(t) = -\infty$ and $\lim_{t\to +\infty }\overline{p}'(t) = -\frac{\mu}{kT}\geq 0.$ Moreover, one has that

Since $\lim_{t\to +\infty}v\overline N t^{5/3}-(1/2-\chi)-t^{2}h(t)-\frac{\mu}{kT} t^{2} = +\infty$ , the existence of $t\in (1, +\infty)$ such that $\overline{p}'(t) = 0$ is guaranteed. The uniqueness of such $t\in (1, +\infty)$ follows from the above Lemma 1. Indeed, by using (A.34), one finds that $(t^{2}\overline{p}'(t))'>0$ for every $t\in (1, +\infty)$ . Since $(t^{2}\overline{p}'(t))' = 2t\overline{p}'(t)+t^{2}\overline{p}^{"}(t)$ , we get that $\overline{p}'(t) = 0$ implies $\overline{p}^{"}(t)>0$ . Therefore, by Lemma 1, there is a unique stationary point $t_{\rm min}\in (1, +\infty)$ , which is also the global minimizer of $\overline{p}$ . Recall that

By applying the Implicit Function Theorem we get that there exist an open neighborhood $\mathcal I\subseteq (0, +\infty)$ of $1$ , an open neighborhood $\mathcal J\subseteq (1, +\infty)$ of $t_1$ and a unique smooth function $\varphi: \mathcal I\to \mathcal J$ , such that

At the same time, one can check, as above, that $p_t\big(\theta, t\big) = 0$ implies $p_{tt}\big(\theta, t\big)>0$ , for every $\theta\in \mathcal I$ . Hence, again by Lemma 1, the following smooth map is well-defined

This concludes the proof of the lemma.

For the sake of brevity, we introduce the following notation. Let the interval ${\mathbb{R}}\supset\mathcal I\ni1$ and the function $\varphi$ be given by Lemma 2. For any $\theta\in \mathcal I$ , set

and define the function

According to this notation, for every $x\in \Omega_h$ and every $h$ sufficiently small we have that

Properties of the family $\{W_h\}$     Hereafter we gather some properties satisfied by the family of 3D energy densities $\{W_h\}$ and by its approximate variant $\{\widehat W_h\}$ defined in (4.21).

Energy wells

Lemma 3. Let the interval ${\mathbb{R}}\supset\mathcal I\ni1$ and the function $\varphi$ be given by Lemma 2. For every $\theta\in \mathcal I$ , the function $W_{\theta}$ is minimized precisely on the set $\alpha_{\theta}\, {\rm{SO}}(3)$ , with $\alpha_{\theta}$ defined as in (A.35).

Proof. Fix $\theta\in \mathcal I$ and denote $m_{\theta}: = \min_{[1, +\infty)}p(\theta, \cdot) = p(\theta, \alpha_{\theta}^3)$ . From inequality (A.32) and by definition of $p$ , we have that $W_{\theta}(F)\geq p(\theta, \det F)\geq m_{\theta}$ , for every $F\in {\mathbb{R}}_1^{3 \times 3}$ . Suppose that $W_{\theta}(F) = m_{\theta}$ , for some $F\in {\mathbb{R}}_1^{3 \times 3}$ . By the definition of $p$ given by (A.33), $F$ must satisfy $|F|^2 = 3(\det F)^{2/3}$ . This implies that all eigenvalues of the positive symmetric matrix $F^{{\rm{T}}}F$ are equal to some $\lambda^2\in \mathbb{R}\setminus \{0\}$ , by (A.32), and accordingly $F\in \lambda {\rm{SO}}(3)$ . On the other hand, it must be $\lambda^3 = \det F = \alpha_{\theta}^3$ (by uniqueness of the minimizer of $p(\theta, \cdot)$ ), implying that indeed $F\in \alpha{\rm{SO}}(3)$ .

Conversely, let $F = \alpha_{\theta}R$ for some $R\in {\rm{SO}}(3)$ . Then $|F|^{2} = 3\alpha_{\theta}^{2} = 3(\alpha_{\theta}^{3})^{2/3} = 3(\det F)^{2/3}$ and thus $W_{\theta}(F) = m_{\theta}$ .

Taking into account (A.37) and the above lemma, we have that

Finally, from (A.35) and from definition (A.31) of $f_h(x)$ , we have that for every $x\in \omega_i\times(-h/2, 0]$ (similarly, for $x\in \omega_i\times(0, h/2)$ ), $i = 1, \ldots n$ ,

by Taylor expansion. This gives the energy well $U_h(x){\rm{SO}}(3)$ defined in (2.4)–(2.6). Let us stress the fact that a crucial property of $U_h(x)$ which allows us to use [1,Theorem 2.6] is that such spontaneous stretch field has the structure $U_h(x) = \alpha{\rm I}_3+hB(x)$ , with $\int_{-h/2}^{h/2}B(x', x_3){\rm d}x_3 = 0$ . Note that this fact is connected with the structure of the chosen density $N_h(x)$ of polymer chains (see (2.1) and the subsequent sentence).

Remark 1. The above analysis on the energy wells of $W_h$ can be rephrased for the the variant $\widehat W_h$ introduced in Section 4, namely to

Indeed, defining the auxiliary function

for $t\in [1, +\infty)$ , we have that $\widehat W_h(x, F)\geq \hat p(\det F, f_h(x))$ , again thanks to inequality (A.32). Also, the very same argument as in the proof of Lemma 3 give that

Note that in this case there is no need of the Implicit Function Theorem to determine the minimizer of $p\big(\cdot, f_h(x)\big)$ : a straightforward computation yields the explicit formula

Using now definition (A.31) of $f_h(x)$ and Taylor expanding, from (A.41)–(A.42) we deduce (4.22)–(4.24).

Regularity    The function $W_{\theta}$ is of class $C^{\infty}$ on the set $\{F\in {\mathbb{R}}^{3 \times 3}: \, \det F>1\}$ and is continuous on ${\mathbb{R}}_1^{3 \times 3}$ . By direct computations, we have for all $F\in {\mathbb{R}}^{3 \times 3}$ with $\det F>1$ and $M, N\in {\mathbb{R}}^{3 \times 3}$ that

and

Uniform convergence    By the definition of $W_{\theta}$ , it is clear that for any sequence $\theta_n\to 1$ the sequence $\{W_{\theta_n}\}$ uniformly converges to $W = W_1$ , as $n\to +\infty$ , on any compact subset of $\{F\in {\mathbb{R}}^{3 \times 3}:\, \det F>1\}$ . Moreover, in light of (A.43) and (A.44), one can show that there exists an open neighborhood $\mathcal U$ of $\alpha{\rm{SO}}(3)$ such that

Taking into account (A.37), this convergence implies in particular the following uniform convergence:

Remark 2. A simple and interesting observation, which relates $W_h$ with the limiting (homogeneous) density $W$ , is that $W_h$ cannot be rewritten in the "prestretch" form

where $U_h(x)$ minimizes $W_h(x, F)$ at each $x\in\Omega_h$ . This is instead the case, for instance, in ^[2,3,22,23] or ^[18].

Quadratic growth    Let $\mathcal I$ be given by Lemma 2 and let $\theta\in \mathcal I$ . By plugging $F = \alpha_{\theta}{\rm I}_3$ into the expression (A.44), we get that

Note that $\lambda(\alpha_{\theta})>0$ as well. Indeed, by definition (A.35) of $\alpha_{\theta}$ and by Lemma 2, we have that $p_t(\theta, \alpha_\theta^3) = 0$ (see (A.33) for the definition of $p$ ), or, equivalently, that

As a consequence, $\lambda(\alpha_{\theta})$ can be rewritten as

Then (A.34) yields

where the last strict inequality is due to the fact that $0 < \chi\leq1/2$ and $\mu\leq0$ .

Now, note that all the differentials of $W_{\theta}$ of order greater than or equal to $3$ do not depend on $\theta$ . Hence, by Taylor expansion of $W_\theta$ around $\alpha{\rm I}_3$ and by the positivity of $G_{\theta}$ and $\lambda(\alpha_{\theta})$ we have that

for every $F$ in a neighborhood of $\alpha_\theta{\rm{SO}}(3)$ . This fact, together with Lemma 3 and the quadratic growth of $W_\theta$ at $\infty$ , shows in particular that

with $U_h(x)$ given by (2.6).

The listed properties show that the family $\{W_h\}$ of 3D energy densities – up to a multiplicative factor $\alpha$ , which is only a matter of the appropriate change of variable (as explained below)– satisfies the hypothesis of [1,Theorem 2.6]. We can then rigorously derive the (new) corresponding 2D model (2.12).

The 2D model    We recall that the 3D-to-2D model reduction of ^[1] which we aim to use is derived therein for 3D spontaneous stretch fields ${\rm I}_3+hB(x)$ , whereas in the present setting such fields are of the form $\alpha{\rm I}_3+hB(x)$ , for some $\alpha>1$ . To reduce to the setting of ^[1], it is sufficient to observe that, setting $\widetilde\Omega: = \alpha\omega\times(-1/2, 1/2)$ ,

where $\tilde v:\widetilde\Omega\to{\mathbb{R}}^3$ and $\tilde v(\tilde x', \tilde x_3) : = v(\tilde x'/\alpha, h\tilde x_3)$ . Now, setting

and, accordingly, defining

we have from definition (A.30) that $W_h\big((\tilde x'/\alpha, h\tilde x_3), \alpha F) = \widetilde W^{\alpha}_{\alpha h}(\tilde x, F)$ . In turn, we obtain from (A.50) that

We have thus reduced to compute the $\Gamma$ -limit, as $\varepsilon \to0$ , of $(1/\varepsilon ^2)\int_{\widetilde\Omega}\widetilde W_\varepsilon ^{\alpha}(\tilde x, \nabla_\varepsilon \tilde v){\rm d}\tilde x$ , where

with

Note also that $\widetilde W^{\alpha}_\varepsilon$ converges uniformly, as $\varepsilon \to0$ , to $\widetilde W^{\alpha}$ , where $\widetilde W^{\alpha}(F) = W(\alpha F)$ and $W$ is defined as in (A.36) with $\theta = 1$ . Hence, using ^[1], we can compute such $\Gamma$ -limit $\widetilde{\mathcal E}_0$ , which is given by

In this expression, the quadratic form $\widetilde Q_2^{\alpha}$ is defined, for every $F\in{\mathbb{R}}^{2\times 2}$ , as

where in the third equality the minimum problem has been solved by using the explicit expression (A.46) (with $\theta = 1$ ), and where $\Lambda(\alpha_1) : = 2G_1\lambda(\alpha_1)/(2G_1+\lambda(\alpha_1))$ , with $G_1$ and $\lambda(\alpha_1)$ given by (A.47). Moreover, in the above expression for $\widetilde{\mathcal E}_0$ , we have that

and

Finally, note that via the correspondence $\tilde u(\tilde x') = u(\tilde x'/\alpha)$ between $\tilde u\in{\rm{W}}^{2, 2}_{1, \rm iso}(\alpha\omega)$ and $u\in{\rm{W}}^{2, 2}_{\alpha, \rm iso}(\omega)$ , and by using that ${\bar A}(x') = \alpha^2 \widetilde A^{\alpha}(\alpha x')$ (see (2.9) for the definition of ${\bar A}$ ), we have that

By the very definition of the tensor $\mathbb C$ in (2.10) it is straightforward to see that the quantity $\mathcal E_0(y)$ in (2.8) with the residual term (2.8) with the residual term (2.13) added is precisely the quantity obtained in (A.51)–(A.52).

Remark 3. We here make a couple of comments about the properties of the approximate energy densities $\widehat W_h$ given in (4.21), for we will use the equivalent expression (A.40).

We start with the simple observation that the dependence on the thickness variable $x_3$ and the thickness parameter $h$ in the approximate model energies $\widehat W_h$ remains unchanged with respect to the original $W_h$ . Thus, taking into account the energy well structure showed in Remark 1, it is straightforward to verify that $\widehat W_h$ shares the same regularity, uniform convergence and quadratic growth properties as $W_h$ .

In particular, the 2D model (obtained by repeating the same procedure as in the general case provided above) will be governed by the functional $\widetilde {\mathcal E}_0$ as in (A.52), with the corresponding quadratic form $Q_2$ obtained via (A.51) from the second differential of the (approximate) homogeneous density $\widehat W$ corresponding to $\widehat W_h$ with $f_h\equiv 1$ in (A.40). Moreover, $\widehat W$ is minimized at $\alpha{\rm{SO}}(3)$ with $\alpha = (\frac{1-2\chi}{2{\rm v}\overline N})^{1/5}$ explicitly determined in (A.41)–(A.42). A direct computation yields

so that $Q_2$ in (A.51) equals $\widehat Q_2$ , defined as

Acknowledgements

This work has been funded by the European Research Council through the ERC Advanced Grant 340685-MicroMotility.

Conflict of interest

All authors declare no conflicts of interest in this paper.