## Abstract

We develop a three-dimensional model for capillary origami systems in which a rectangular plate has finite thickness, is allowed to stretch and undergoes small deflections. This latter constraint limits our description of the encapsulation process to its initial folding phase. We first simplify the resulting system of equations to two dimensions by assuming that the plate has infinite aspect ratio, which allows us to compare our approach to known two-dimensional capillary origami models for inextensible plates. Moreover, as this two-dimensional model is exactly solvable, we give an expression for its solution in terms of its parameters. We then turn to the full three-dimensional model in the limit of small drop volume and provide numerical simulations showing how the plate and the drop deform due to the effect of capillary forces.

## 1. Introduction

Elasto-capillary effects, which arise due to the interaction between elastic and capillary forces, are common in nature and everyday life, but often go largely unnoticed. This is because in our familiar macroscale environment, elastic components are usually constructed to withstand gravity and are often too rigid to be notably deformed by capillary forces. Exceptions of course include the coalescence (or clumping) of wet hair or of paintbrush bristles [1], or the adhesion of a thin, flexible sheet (e.g. made of mylar) to a soap bubble. As systems shrink in size, surface effects overtake volume effects, gravity becomes negligible and surface tension gains importance. Balancing capillary forces can then cause large elastic deformations and lead to unexpected phenomena, such as the ejection of fungal spores [2] and the buckling of pulmonary airways [3].

With a general move towards miniaturization in technology, understanding the principles of elasto-capillarity is crucial. For example, standard fabrication techniques for building micro-electromechanical systems involve wet etching, and in the drying process flexible parts become susceptible to stiction [4]. This can deform or even break components, thereby compromising the functionality of the device. On the positive side, capillarity can be used in the design of microstructures, as illustrated by the development of passive pipetting techniques [5] and by alternative approaches to self-assembly [6]. Other applications include the formation of foams by evaporating water from a vertical array of multiwalled nanotubes [7], the engineering of highly hydrophobic surfaces in the form of coated carbon nanotubes [8], and the fact that locomotion of robots and insects on interfaces is facilitated by means of flexible appendages [9,10].

One particular promising aspect of elasto-capillarity is the ability of surface tension to fold planar objects, possibly with hinges, into three-dimensional structures [6]. For elastic sheets, this process is known as capillary origami [11,12]. It could potentially lead to a new self-assembling approach to microscale manufacturing that is an alternative to standard methods, such as photolithography, which involves the deposition and etching of thin films [13] and is essentially two-dimensional. Applications have already been proposed to efficiently produce three-dimensional photovoltaic cells [14] and to form graphene nanostructures [15].

The concept of capillary origami was developed by Py *et al.* [11,12], who demonstrated that a variety of three-dimensional structures, such as tubes, pyramids and spheres, could be created by this method. Their work also included numerical simulations for predicting the shapes formed by a two-dimensional, pinned-contact line model and has led to a number of other studies. Extensions involve the addition of drop impact dynamics to speed up fabrication time [16] and electric fields for control [17]. Additionally, other approaches have focused on alternative setups to better control folding [18,19].

To the best of our knowledge, all of the theoretical descriptions of capillary origami available in the literature are based on two-dimensional set-ups or on simplified three-dimensional scenarios where reductions can be made. In the original articles, Py and co-workers wrote ad hoc governing equations in two dimensions based on Euler's elastica and assumed that the liquid–air interface of the drop is a circular arc whose contact line is pinned to the boundary of the plate. Their work was extended to three dimensions in [20], although in a less realistic situation due to assumed geometries for the plate and interface, and also further simplified in [19,21]. Other articles have explored the effect of plate geometry and fluid wetting properties on folding [22]. Recently, the authors have focused on putting two-dimensional capillary origami in a more rigorous variational framework [23–25].

The benefit of two-dimensional or simplified three-dimensional models is that they provide good qualitative information about the features of the system, capturing the spontaneous wrapping phenomenon [11] and predicting the instant encapsulation threshold [23] that is also present, yet not explicitly mentioned, in other similar systems [5]. However, the downside is that none of the reduced set-ups appropriately capture three-dimensional scenarios, resulting in inaccurate quantitative predictions. As a consequence, there is a need for the development of three-dimensional models that take into account a full description of plate elasticity.

In this article, we begin the process of developing such a model. Neglecting the effect of gravity, we assume capillary origami configurations are minimizers of an appropriate energy functional containing only the interfacial energies of the drop and the deformation energy of the sheet. Interfacial energies form the basis of the theory of capillary surfaces and come with a well-developed set of tools [26]. The handling of the deformation energy is more problematic and has been the focus of numerous recent studies (e.g. [27–31]). The difficulty lies in the fact that the energy is composed of two parts: a non-convex in-plane stretching term penalizing deviations from a flat metric and a convex out-of-plane bending term involving both the mean and Gaussian curvatures. When the thickness of the sheet limits to zero, the resulting equilibrium shape involves only bending from the initial configuration [32], and the plate can thus be locally described as a developable surface. Imposing the isometric constraint becomes difficult due to its degenerate and point-wise nature. Our present approach consists of avoiding this problem by allowing for bending using a Föppl–von Kármán ansatz, and should thus be viewed as a first step towards the description of the three-dimensional encapsulation process. Indeed, although the ansatz discussed here is not expected to be valid for the large deformations that are needed to study the requisite encapsulation, it nevertheless captures the early stages of the process, which are often inherently three-dimensional and may be crucial for determining the final structure [12].

Little work has been done on modelling three-dimensional elasto-capillary systems involving the deflection of thin plates. The first studies were carried out by Shanahan [33,34], then modified by Olives [35,36] and recently adapted in [37,38] to understand the wrinkling patterns seen in the experiments of [39]. All of them focused on radially symmetric arrangements with partially wetted plates that have a fixed boundary or a prescribed far-field behaviour. Motivated by the capillary origami set-up, our work is concerned with a very different situation: we assume that the plate is completely wetted, has a free boundary, and does not become radially symmetric, except maybe in rare cases where it is much more beneficial for the plate to stretch than bend. In the conditions of this study, boundary effects are important, and affect the overall qualitative shape of the sheet.

The rest of this article is organized as follows. The general model is presented in §2. In §3, we focus our attention on a rectangular plate in two specific limits: infinite length, which allows us to connect the present work to previously established results for two-dimensional capillary origami; and small liquid volume, which leads to a genuine, albeit simple, tri-dimensional problem that explains how the shape of the fluid interface is affected by the finite extent of the plate. Section 4 summarizes our work and discusses potential applications and extensions of the present model.

## 2. General model

To formulate our general model for capillary induced elastic deflections, we consider a thin planar sheet of thickness *h* and characteristic size , with , resting on top of a rigid, hydrophobic substrate, as illustrated in figure 1. Centred at the origin, the mid-plane of the plate traces out the domain *Ω*, which is equipped with the standard Euclidean metric. When a drop of liquid *Ξ* with volume is placed on top of the plate, the latter deforms around the drop to a static configuration with mid-plane *Σ* described by the immersion and corresponding metric *g*_{Σ} induced from the parametrization of the plate. We assume that the liquid completely wets the sheet and is thus pinned to the boundary of the plate. If we represent the liquid–air interface by the surface *Γ* with parametrization , this boundary requirement can be expressed as ∂*Σ*=∂*Γ* or equivalently, **x**|_{∂M}=** ξ**|

_{∂Ω}.

We assume that the plate is a homogeneous, isotropic material with vanishing Poisson ratio and that equilibrium configurations (**x**,** ξ**) are minimizers of the total dimensionless energy of the system
2.1
over the set of functions in that are equal on the boundary,

**x**|

_{∂M}=

**|**

*ξ*_{∂Ω}and preserve the volume of the drop. Here,

*H*

_{Σ}and

*K*

_{Σ}denote the mean and Gaussian curvatures of

*Σ*, ∥

*g*∥

^{2}=tr(

*g*

^{T}

*g*), and

*ς*is a known surface energy density of the plate. Since (2.1) is in dimensionless form, it is from now on implicitly assumed that all lengths have been scaled by a characteristic length of the plate, .

The first two terms in (2.1) come from the interfacial energies of the drop and the plate. The surface energy density *ς* is obtained by noting that the surface elements d*Σ*_{+} and d*Σ*_{−} of, respectively, the liquid–solid and air–solid interfaces of the plate can be expanded as
where is a dimensionless parameter that measures the thinness of the plate and by assumption satisfies *δ*≫1. As a consequence, the interfacial energy terms are all expressed relative to the mid-plane *Σ*. With this approximation, the density *ς* is given by
see appendix A. Here, *γ* is the liquid–air interfacial (or surface) tension of the liquid, and *γ*_{SV} and *γ*_{LS} are the adhesion coefficients of the plate between the solid–vapour and liquid–solid interfaces, respectively [26]. The perturbation *ρ* is non-zero because of the finite thickness of the plate: when the plate bends, one side stretches and the other compresses.

The last two terms of (2.1) are related to the deformation of the plate. In particular, they are an asymptotic expression of the full three-dimensional elastic energy functional for small thickness [28] and contain the dimensionless parameter λ, given by for the plate's bending stiffness *B*, capturing the relative strength of surface tension effects versus plate rigidity. These two deformation terms, respectively, measure stretching and bending of the mid-plane *Σ*, and both vanish only when ** ξ** is a flat, isometric immersion. For -immersions, it is known via

*Γ*-convergence that, as , the deformation energy reduces to minimizing the bending energy over

*g*

_{Σ}=

*I*, implying that

**is an isometry [32] and no stretching occurs.**

*ξ*The inclusion of the effect of stretching in the deformed plate of finite thinness *δ* can be achieved through a Föppl–von Kármán small slopes approximation by assuming that the parametrization ** ξ** of the deflected plate

*Σ*is of the form 2.2 for (

*u*,

*v*) in

*Ω*, where

*χ*:=(

*χ*

_{1},

*χ*

_{2}) and

*ζ*are, respectively, the in-plane and out-of-plane components. As a consequence, the strain of the elastic plate is assumed to be small. When

*χ*is of order and

*ζ*is of order , the deformation energy at leading order becomes 2.3 where

*D*

^{2}

*ζ*is the second derivative of

*ζ*, i.e. its Hessian matrix in a prescribed frame,

*ε*is the in-plane strain tensor 2.4 for the rank one 2×2 tensor ∇

*ζ*⊗∇

*ζ*:=(∇

*ζ*)(∇

*ζ*)

^{T}and the vector-valued function

*χ*=(

*χ*

_{1},

*χ*

_{2})

^{T}, with gradient matrix

*Dχ*=(∇

*χ*

_{1},∇

*χ*

_{2}) [40]. Also in this ansatz, the surface element reads , and the surface energy density

*ς*can be approximated by the effective adhesion parameter

*τ*=−(

*γ*

_{SV}+

*γ*

_{LS})/

*γ*, which simplifies the corresponding interfacial energy. Therefore, the reduced total energy becomes 2.5 up to an additive constant. We look for minimizers over the set of functions (

**x**,

*χ*,

*ζ*) in the space of functions that keep ∂

*Γ*=∂

*Σ*and the volume of the drop fixed. It is then necessary that the first variation of vanish relative to the two aforementioned constraints. By the method of Lagrange multipliers, this may be written as 2.6 for , where is the volume functional of the drop and denotes the Fréchet derivative of the functional . System (2.6) forms our governing model.

To compute the variations of and , we introduce the variations of **x** and ** ξ** defined, respectively, by the sufficiently smooth one-parameter families of immersions and , for

*μ*small, |

*μ*|<

*μ*

_{0}, satisfying 2.7 for the smooth functions

*u*

_{Γ},

*X*

_{1},

*X*

_{2}and

*Z*. Here,

**n**

_{Γ}is the Gauss map of

*Γ*, directed out of the fluid and

**T**

_{Γ}is a given tangential vector field. Ensuring that these variations stay attached on their boundaries requires 2.8 hence, the boundary variations of the drop are completely determined by those of the plate.

To handle the terms involving the plate, we first parametrize the boundary ∂*Ω* of the undeflected state by its arc-length using the function ** ω**, oriented counter-clockwise; hence its tangent becomes

**=**

*t***′ with outward directed normal**

*ω***. Denoting the angle of inclination between the**

*n**u*-axis and the normal

**by**

*n**θ*, these two vectors can be rewritten and . The first variations of the bending and stretching terms can be then be calculated in a standard manner; for details, see [41]. The results are 2.9 and 2.10 where d

*l*refers to the line element on ∂

*Ω*. For the plate's adhesion term, the perturbed strain tensor is obtained by substituting the second equation of (2.7) in the expression for the strain tensor given in (2.4) and expanding the result to order

*μ*. We obtain

*ε*

_{μ}=

*ε*+

*μ*(

*DX*+(

*DX*)

^{T}+∇

*ζ*⊗∇

*Z*+∇

*Z*⊗∇

*ζ*)+…, where

*X*=(

*X*

_{1},

*X*

_{2}) and the term in

*μ*is the lowest order contribution to the strain tensor resulting from changing the position of the plate by

*μ*(

*X*

_{1},

*X*

_{2},

*Z*), as described by the second equation of (2.7). Sequentially taking the trace, differentiating the result with respect to

*μ*and integrating over

*Ω*gives 2.11 where the divergence theorem was additionally applied to remove the derivatives from

*X*.

From work on capillary surfaces, it is well known that the first variation of the surface area functional, with moving boundary, is
2.12
where *ν*_{Γ} is the co-normal of the boundary ∂*Γ* and *H*_{Γ} is the mean curvature of *Γ* (e.g. [26,42]). Using equation (2.8), we can rewrite the product **T** _{Γ}⋅*ν*_{Γ} on the boundary of the drop ∂*Γ* as
for the standard Euclidean basis {**e**_{1},**e**_{2},**e**_{3}} of . Also since the boundaries of the drop and plate are equal, the contour integral in (2.12) can be taken along ∂*Σ*. Then evaluating the mapping ** ξ** on the parametrization

**of the boundary curve ∂**

*ω**Ω*gives a parametrization

**(**

*γ**t*):=

**(**

*ξ***(**

*ω**t*)) of ∂

*Σ*, and the first variation (2.12) can be rewritten as 2.13 For plates of the form (2.2), we approximate the volume of the drop by for

*δ*≫1, where

*Ξ*

_{0}is the region occupied by the drop when the plate is in an undeflected configuration. Its variation then becomes 2.14 after differentiating under the integral sign [43,44].

Finally, by adding equations (2.9)–(2.11), (2.13) and (2.14) after multiplying them by their appropriate coefficients, we obtain from (2.6) that *Γ* satisfies
2.15a
and the plate terms (*χ*,*ζ*) satisfy the partial differential equations
2.15b
and
2.15c
with natural boundary conditions
2.15d
2.15e
2.15f
Problem (2.15) forms our general model for elasto-capillary deformations of a thin sheet undergoing small deflections, including stretching. Since *β* is fixed, equation (2.15a) shows that, as is common in capillary problems, the liquid–air interface of the drop is necessarily a surface of constant mean curvature. Equations (2.15b,*c*) indicate that the plate satisfies a variant of the Föppl–von Kármán equations [41]. In particular, equation (2.15c) governs the in-plane stretching, and (2.15b) determines the out-of-plane bending, where adhesion effects due to the presence of the drop are captured by the third term on its left-hand side. Although it is standard to introduce an Airy stress potential to remove the in-plane strain tensor, we avoid this approach since we would like to explicitly keep the drop attached to the boundary of the plate. These PDEs are coupled through the boundary conditions given in (2.15d–*g*), as well as through the Lagrange multiplier *β* used to enforce the volume of the drop. The value of *β* captures the pressure difference between the drop and the surrounding air, which depends on the volume of liquid in the drop, and can be treated as a bifurcation parameter.

## 3. Models for a rectangular plate

To continue our study of elasto-capillarity, we next look at two specific situations in which model (2.15) can be simplified. Both arise from a rectangular initial plate *Ω*. In the first case, we assume that one of the orthogonal directions is much longer than the other, namely, *Ω*=[−1,1]×[−*L*,*L*], for *L*≫1, and suppose that the plate is approximated by an infinite strip that bends in the short direction. As a consequence, drop-plate configurations are uniform by translation in the long direction of the plate. System (2.15) can thus be reduced to governing equations for the cross-sectional curves, yielding a new two-dimensional model, which includes stretching, and which can be compared with previous two-dimensional descriptions of capillary origami. In the second situation, we do not assume any specific aspect ratio for the plate, but restrict ourselves to small volumes of liquid. This allows us to simplify the partial differential equations for the plate in (2.15b,*c*), while keeping the interface of the drop inherently three-dimensional. The result is a more physically realistic description of the system, in which the full effects of the boundary may be studied.

When the initial plate configuration corresponds to *Ω*=[−1,1]×[−*L*,*L*], the boundary ∂*Ω* can be decomposed into four components so that conditions (2.15d–*g*) may be simplified. On the right component, when *u*=1 and |*v*|≤*L*, the angle of inclination is *θ*=0, giving the normal vector ** n**=(1,0) and the tangent vector

**=(0,1). Likewise, on the left component (**

*t**u*=−1 and |

*v*|≤

*L*), the angle of inclination is

*θ*=

*π*, yielding

**=(−1,0) and**

*n***=(0,−1). Consequently, the right and left boundary conditions in equations (2.15d–**

*t**g*) reduce to 3.1a and 3.1b where ∂

*Ω*

_{±}:={(±1,

*v*):|

*v*|≤

*L*} and the curve

**on ∂**

*γ**Ω*

_{±}is given by

**(**

*γ**v*):=

**(±1,**

*ξ**v*). Similarly the top and bottom boundary conditions on

*v*=±

*L*become 3.1c and 3.1d where ∂

*Ω*

^{±}:={(

*u*,±

*L*):|

*u*|≤1} and the curve

**on ∂**

*γ**Ω*

_{±}is given by

**(**

*γ**u*):=

**(**

*ξ**u*,

*L*). Note that the natural boundary conditions given in (2.15d–

*g*) are derived in a context where parametrization of the normal and tangent vectors of ∂

*Ω*is given in terms of a well-defined angle of inclination, and that the method breaks down at the corner points of a domain. However, by using integration by parts locally, it can be shown that an additional term arises from the Hessian, which requires the mixed partial derivative

*ζ*

_{uv}to vanish at corners. We thus have 3.1e Including these boundary conditions with those given above, the model for the rectangular plate becomes (2.15a–

*c*) coupled with (3.1).

### (a) Infinite strip

The simplest situation in which to study our rectangular model is to assume that the bottom plate becomes the infinite dimensional strip (i.e. ) and that deformations are invariant by translation in the infinite direction. Hence
and
where we have assumed that cross sections of the liquid–air interface *Γ* in a plane perpendicular to *e*_{2} have length 2ℓ and arc length *s*. As a consequence, the drop has infinite volume in three dimensions, but finite cross-sectional area in two dimensions. The surfaces *Γ* and *Σ* are completely determined by the two-dimensional curves (*x*,*z*) and (*χ*_{1},*ζ*), and the model on the rectangle becomes
3.2a
for the liquid–air interface, where *κ*_{Γ} is the curvature of the (*x*,*z*) curve, and
3.2b
and
3.2c
for the plate, with the boundary conditions
3.2d
and
3.2e
In the above, a prime denotes differentiation with respect to *u* and a dot differentiation with respect to *s*. For the boundary conditions, we have used the facts that in the three-dimensional set-up ∂*Γ* is parallel to **e**_{2}, that its co-normal is given by , and that ** γ**(

*v*):=

**(±1+**

*ξ**χ*

_{1}(±1),

*v*).

It is informative to compare equations (3.2) to other models for two-dimensional capillary origami, which typically assume inextensible plates, in situations where strong contact line pinning is expected. Specifically, such models read [11,23]
3.3a
for the curvatures *κ*_{Γ} and *κ*_{Σ} of the two-dimensional mappings (*x*,*z*) and (*ξ*,*η*), respectively, representing the air–liquid interface and the plate, and where prime denotes differentiation with respect to arc-length. Equations (3.3a) may be obtained by describing two-dimensional cross sections of *Σ* as an elastica subject to uniform pressure and to surface tension forces at its endpoints [11], or equivalently [23] by minimizing the total energy of the two-dimensional drop-plate system, assuming that the curve (*ξ*,*η*), which represents the plate, has constant length, and that its endpoints are pinned to the endpoints of (*x*,*z*), which represents the air–liquid interface. These equations are to be solved with the following boundary conditions written at the endpoints of *Γ* and *Σ*:
3.3b
where **n**_{Σ} is the unit normal of the plate pointing into the liquid, *ν*_{Γ} and *ν*_{Σ} are the boundary co-normals of *Γ* and *Σ* (i.e. the outward pointing tangent vectors of (*x*,*z*) and (*ξ*,*η*)), and *α* is a Lagrange multiplier used to fix the length of the plate. In this model, the second to last boundary condition captures a transverse shear force caused by the surface tension of the liquid–air interface [45], and the last boundary condition, which enforces a contact angle, reduces to the standard condition [26] when dewetting occurs [24]. By assuming that (*ξ*,*η*) is a perturbation of the flat state of the form (*ξ*,*η*)=(*u*+*χ*_{1}(*u*),*ζ*(*u*)) for *u*∈(−1,1), where *χ*_{1} and *ζ*^{2} are of order for *δ*≫1, the differential equation in (3.3) for the curvature *κ*_{Σ} simplifies to
3.4a
as . The last three boundary conditions in (3.3b) become
3.4b
and
3.4c
where (*x*,*z*) is assumed to be parametrized in terms of its arc-length with unit tangent . In this form, it is easy to compare the two models and see that (3.2) is a modification of (3.4), valid for small plate deflections, that includes stretching through (3.2c). Indeed, from (3.2c) we can set *ζ*′^{2}+2*χ*_{1}^{′}=*λC*, where *C* is constant. Then, equations (3.4) are the same as (3.2), provided and . Given that the last equation of (3.3b) relates *α* to the angle made by the two-dimensional plate and air–liquid interfaces at their endpoints, equation (3.2e) can be viewed as a contact angle condition. Similarly, the equations of (3.2d) indicate, respectively, that the plate is torque free on its boundary and also subject to a transverse shear that accounts for the strain present in the plate due to stretching [45].

Informed by this reduction and the overall form of (3.2), we do not expect solutions to exist for all scalings of the curve (*x*,*z*) representing the liquid–air interface. Observe that since and , the left-hand side of equation (3.2b) is of order . Balancing then requires *β*—a free parameter with unknown scaling enforcing the volume condition—to be of order , and from (3.2a), we see that the curvature of the interface must then be small. In physical terms, this means that the two-dimensional cross section of the drop limits either to an infinitely large or infinitely small cap. We however do not make specific assumptions on the size of the known parameters, other than being of order one.

A convenient aspect of problem (3.2) is that it is exactly solvable. First, observe that since the interface *Γ* has constant curvature, it traces out the arc of a circle parametrized by
where we have applied symmetry conditions to remove rotations and horizontal translations (*x*(0)=0). Assuming and implies the expansions
for *δ*≫1. Next, we set (*χ*_{1},*ζ*)=(*δ*^{−2}*χ*_{10},*δ*^{−1}*ζ*_{0})+⋯ . Note that this scaling implies , which as discussed above is necessary for the equations in (3.2) to match the two-dimensional model (3.4). At leading order, (3.2b–*e*) become

The above equations are exactly solvable since the nonlinear term in (3.5a) can be removed by integrating (3.5b). Namely, setting *ζ*_{0}′^{2}+2*χ*_{10}′=*C*_{1} reduces the problem to
3.6a
and
3.6b
with *C*_{1}=λ(*τ*−1)/2. This last conditions makes (*τλ*−2*C*_{1})=λ; hence by solving (3.6a) and using the symmetry conditions *ζ*_{0}(0)=*ζ*_{0}′(0)=0 to remove rotation and translation invariances, we find
Consequently, *χ*_{10} can be completely determined as well by imposing *χ*_{10}(0)=0 and is given by
The above expression can become quite large for specific values of λ and is therefore valid only if it remains of order one. In expanding the arc-length of the interface as ℓ=ℓ_{0}+*δ*^{−2}ℓ_{1}+⋯ , we deduce that ℓ_{0}=1, ℓ_{1}=*β*_{0}^{2}/6+*χ*_{10}(1), and *z*_{0}=*ζ*_{0}(1)−*β*_{0}/2 from the boundary condition (3.6b) and the last two equations of (3.2a), respectively. The expression for *z*_{0} reads
and will remain of order one as long as and 1/λ are of order one. One can check that for *τ*=−0.5 and λ between 4 and 12, both |*χ*_{10}(1)| and *z*_{0} are less than 5. Finally, integrating (*ξ*′(*u*)^{2}+*η*′(*u*)^{2})^{1/2} over [−1,1] gives the length of the plate:
3.7
which shows that the plate compresses depending on the value of *τ*.

In summary, the above description gives a lowest order approximation of two-dimensional capillary origami configurations when the plate is allowed to stretch, and when the air–liquid interface has small curvature. For fixed values of λ, *δ* and *τ*, there is only one free parameter: *β*_{0}. The latter may be set by imposing the volume of liquid in a cross section of the system, that is by choosing the area of the region enclosed by the two curves (*ξ*,*η*) and (*x*,*z*), as long of course as the condition remains satisfied. Figure 2 shows the two-dimensional configurations predicted by equations (3.6) for three different values of λ, with *β*_{0}=1, *δ*=100 and *τ*=−0.5.

### (b) Small volume of liquid

The previous approach reduces the model to two dimensions by assuming translational invariance in the long direction of the strip. While this significantly lessens the difficulty of the problem, it is physically unrealistic. Therefore, in this section, we investigate the behaviour of the unconstrained deformations of the rectangular plate using model (2.15a–*c*) coupled with (3.1) under a small-height scaling of the liquid–air interface. This assumption is valid for small volumes and produces a problem that includes three-dimensional effects from the interface and the boundary.

Motivated by the work on the infinite strip, we again do not expect solutions of (2.15a–*c*) to exist for all scalings of **x**. Instead, by explicitly expanding the components of the plate as
3.8a
we see that balancing equation (2.15b) requires the mean curvature *β* of *Γ* to be of order . That is,
3.8b
and consequently, the liquid–air interface of the drop is nearly flat, corresponding to essentially a very large or very small cap. Our focus will be on the latter. In particular, we assume that the shape of the liquid–air interface has the non-parametric form **x**(*x*,*y*)=(*x*,*y*,*z*(*x*,*y*)) with vertical coordinate admitting the asymptotic expansion
3.8c
This reduces the constant mean curvature problem (2.15a) to
where the Laplacian Δ is taken with respect to the standard flat metric and the boundary conditions enforce ∂*Γ*=∂*Σ*. Note that as a consequence of our small-volume assumption, we have also set *M*=*Ω*. Using (3.7c), the Gauss map of *Γ* becomes , yielding the exterior pointing co-normals
on the right/left and top/bottom boundaries of the plate, respectively. Note that these boundaries are parametrized by ** ξ**(±1,

*v*) and

**(**

*ξ**u*,±

*L*) so that their corresponding tangents

**′ become either (±1,0,0) or (0,±1,0) at leading order. Therefore, substitution of these expressions into (3.1), along with plate expansion in (3.7a), gives the new right/left boundary conditions on**

*γ**u*=±1,|

*v*|≤

*L*, and the new top/bottom boundary conditions on |

*u*|≤1,

*v*=±

*L*, where (

*x*,

*y*)=(

*u*,

*v*), the notation

*ζ*

_{0,u}denotes the partial derivative of

*ζ*

_{0}with respect

*u*, and

*ε*

_{0}is the rescaled strain tensor for

*χ*

_{0}:=(

*χ*

_{10},

*χ*

_{20}). Finally by substituting the expansions (3.7a) into the PDEs for the plate (2.15a–

*c*), the model of the rectangular plate at leading order gives 3.9a for the interface of the drop (

*x*,

*y*,

*δ*

^{−1}

*z*(

*x*,

*y*)), and 3.9b 3.9c 3.9d 3.9e for the out-of-plane deformation of the plate, and 3.9f 3.9g 3.9h for the in-plane deformation of the plate, where ∂

*Ω*

_{±}and ∂

*Ω*

^{±}denote the right/left and the top/bottom boundaries, respectively, defined by ∂

*Ω*

_{±}:={(

*x*,

*y*,

*u*,

*v*):

*x*=

*u*=±1,

*y*=

*v*∈(−

*L*,

*L*)} and ∂

*Ω*

^{±}:={(

*x*,

*y*,

*u*,

*v*):

*x*=

*u*∈(−1,1),

*y*=

*v*=±

*L*}. Here for convenience, we have dropped the previously introduced order-zero subscripts. Recall that in our assumed scaling, the metric of the plate

*g*

_{Σ}expands as , and therefore its deviations, i.e. the in-plane strain tensor

*ε*, are completely determined by problem (3.9f–h). Notice that as in the infinite strip case, it depends directly on the value of the constant λ(

*τ*−1).

In the above model, the parameter *β* is free, set by fixing the volume of the drop, and the values of the parameters *δ*, *τ*, *L* and λ are known, where the latter captures the flexibility of the sheet: λ≫1 is a highly bendable regime and λ≪1 is the rigid one. Since (3.8) was derived under a small-deflection ansatz, λ must be of order for the solutions to remain physically valid.

Model (3.8) is nonlinear due to the inclusion of the stretching terms, which couple the in-plane and out-of-plane deformations of the plate. We solve these equations numerically by means of an iterative Newton–Raphson numerical scheme.

In carrying out numerical simulations (figures 3 and 4), we observe that initially for small volumes of liquid the corners of the plate are pulled up evenly by the drop, as expected in the early stages of capillary origami folding and noted in the initial experiments of [11]. By evenly, we mean that the vertical deflection, or the bending energy density, remains essentially radially symmetric about the centre of the plate. This is true for both square (figure 3) and rectangular (figure 4) plates. Most of the stretching in the deformation of the plate occurs near its boundary and is localized at the four points located in the middle of each edge. This suggests that cusps are likely to form at these points when the plate is inextensible, thereby legitimizing the method of using rigid components connected by hinges [6], or decomposing the plate into developable pieces. In the case of a rectangular plate, bending preferentially occurs along its longer direction, as can be seen in figure 4.

When the volume of the drop is sufficiently small, the configuration described above is the unique configuration. However, when the volume is continued past some threshold, denoted , two more equilibrium solutions appear. This regime corresponds to values of *β*_{0} that are of order 40 and are therefore outside of the strict range of validity of the model, since we assumed *β*_{0} to be of order one. This behaviour is however interesting because it illustrates a symmetry-breaking phenomenon that is observed in experiments. Figure 5 shows the deformation energy for a nearly square plate with (*L*=1.1), for . For small volumes there is one solution, as described above, whose corners are evenly deflected. As the volume is increased, the corners continue to deflect evenly until . At that point, the plate starts to bend more in its short direction and the deformation energy of the plate, computed via (2.3), increases more rapidly. At , the two solutions that appear have lower deformation energy: the one of lowest energy bends approximately equally in both directions, while the other bends more in the long direction of the plate.

A pitchfork bifurcation occurs at the threshold when *L*=1 and unfolds if the aspect ratio of the plate is not equal to 1, as shown in figure 5. This (imperfect) bifurcation captures a symmetry-breaking transition in the shape of the plate, from mode-4 (evenly deflected corners) to mode-2 (quasi-cylindrical), which occurs as the volume of the drop, as well as the deformation of the plate, increase. This is reminiscent of observations noted in the experiments of Py *et al.* [11]. The branch with lower deformation energy, which is the continuation past the bifurcation threshold of the unique branch that exists below the bifurcation threshold, is expected to be unstable for . We have checked that when interfacial energies are taken into account, this branch has indeed higher energy than the other two, quasi-cylindrical, branches for . Figure 5 shows the deformation energy as opposed to the total energy because the bifurcation is difficult to visualize in the latter representation, due the fact that all of the solutions have comparable total energy.

To summarize, the three-dimensional model discussed here, although limited to small deflections of the plate, is able to reproduce the initial configurations seen in capillary origami experiments, as the plate begins to wrap around the drop. Although it cannot capture folding, it nevertheless provides a physically valid reduction, which can be used to understand the three-dimensional effects that are of paramount importance in the early stages of the encapsulation process. When pushed beyond its strict limit of validity, it also captures a symmetry-breaking bifurcation that favours quasi-cylindrical bending.

## 4. Further remarks

In this article, we took a first pass at modelling three-dimensional capillary origami for systems in which the plate has finite thickness and is allowed to stretch. We assumed that the plate was rectangular and underwent small deflections, and considered two simplifying situations: infinite aspect ratio and small liquid volume.

The first simplification led to an exactly solvable model, which allowed us to conclude that the compression of the plate is directly related to the dimensionless parameters λ and *τ*. In particular, the plate shrinks depending on the relative strength of adhesive and capillary forces, as measured by *τ*. This model also illustrates how existing two-dimensional capillary origami descriptions may be modified to include stretching of the plate.

The second simplification includes three-dimensional effects in the small-volume limit where the mean curvature of the liquid–air interface scales with the out-of-plane deflection of the plate. This model removes the inherent limitations of two-dimensional models of capillary origami, which cannot capture the early stages of the folding process. It also shows that when stretching is allowed, it tends to concentrate on regions that bisect the edges of the plate. This justifies the approach of decomposing isometrically constrained plates into patches that connect with each other in a *C*^{1} fashion. This new three-dimensional model captures the seemingly inherent transition in the deformation of a rectangular plate from a shape with evenly deflected corners to one that is quasi-cylindrical. Such a transition, which is observed in experiments, is described here as being associated with the onset of a pitchfork bifurcation, which is imperfect if the plate is not square. We therefore believe that understanding the nature of this bifurcation may hold the key to predicting spontaneous encapsulation for rectangular plates.

The limitations of the present approach are that it is only physically valid for small deformations and cannot accurately capture fully folded capillary origami stages. Modelling the latter would require a complete description of the geometry of the plate, together with the boundary effects due to strong contact line pinning with the drop. This is a highly challenging problem, not only analytically but also from a numerical perspective.

## Authors' contributions

N.B. and J.L. designed the study, analysed the results and wrote the manuscript. N.B. did the calculations and the three-dimensional simulations shown in figures 3–5. J.L. reviewed the calculations and produced figure 2. Both authors gave final approval for publication.

## Competing interests

We declare we have no competing interests.

## Funding

This material is based upon work supported by the National Science Foundation under award no. 1304090.

## Appendix A. Plate adhesion energy

To derive the reduced two-dimensional interfacial energy of the plate that is given in (2.1), we assume that the three-dimensional configuration of the plate can be represented by

where is the parametrization of mid-plane with normal *n*_{Σ} and *z*∈(−*h*/2,*h*/2) is the thin direction for *h*≪1 (cf. [28]). Then the upper and lower surfaces of the plate are parametrized by
respectively, and have corresponding surface elements , where *H*_{Σ} is the mean curvature of the mid-plane. Assuming the top surface is in contact with the liquid, while the bottom surface is in contact with the vapour, the interfacial energy of the plate becomes
By rescaling *u* and *v* by the characteristic length and then dividing by the factor , the right-hand side reads
Then the dimensionless interfacial energy of the plate can be rewritten as an integral over the mid-surface:
where the integrand
for the dimensionless parameter .

## Footnotes

One contribution of 11 to a theme issue ‘Trends and challenges in the mechanics of complex materials’.

- Accepted November 13, 2015.

- © 2016 The Author(s)

Published by the Royal Society. All rights reserved.