Folds in thin strips are extremely localized curvatures that occur over a small length-scale. In the limit of such length to zero, the shape of the strip transforms into one with a sharp corner (Lechenault et al., 2014). Such mechanical instabilities in thin strips are often the consequence of an over-constraining environment and material mismatch (Tanaka et al., 1987; Bowden et al., 1998), such as films on foundations.

Examples of patterns generated in confined conditions (Wang and Zhao, 2015) include wrinkles (Huang et al., 2004, 2005; Guvendiren et al., 2010), creases (Hong et al., 2009; Lechenault et al., 2014), blisters, folds (Conti and Maggi, 2007; Pocivavsek et al., 2008; Bayart et al., 2014), crumples (Vliegenthart and Gompper, 2006; Deboeuf et al., 2013), and ridges (Lobkovsky et al., 1995; Lobkovsky, 1996; Lobkovsky and Witten, 1997; Audoly and Pomeau, 2010). The transitions from one to the other manifest as a consequence of the change in the boundary conditions (Holmes and Crosby, 2010; Jin et al., 2015).

Several researchers explored the mechanical snap-back and snap-through instabilities of thin structures constrained by controlled boundary conditions, both on displacements and rotations.

In Beharic et al. (2014), the authors examined the bi-stable snap-through behavior of a beam with both positive and negative curvatures. To achieve such configuration, they confined the beam with applied symmetrical rotations at both ends and end-to-end horizontal span smaller than the arc-length. By deriving a strain energy landscape, they were able to determine a critical angle at which the deformation becomes mono-stable. Such critical angle depends monotonically only on the (compression) ratio between the arc-length and the horizontal span and not on the material or thickness, (Gomez et al., 2017) by studying a very similar system, demonstrated that the speed of the snapping is lower than the prediction based on the speed of sound; in fact, even without dissipation, they discovered the existence of a slow-down near criticality.

Cazzolli and Dal Corso (2018) studied the snapping mechanisms for elastic strips without folds, generalizing the work of Beharic et al. (2014). They derived an energy release map of all the possible rotations (symmetrical and non-symmetrical) for which the strip exhibits a snap mechanism. Other recent works on snapping in over-confined thin strips include (Yu and Hanna, 2019), where they used the Kirchhoff equations for anisotropic rods to explain the stability under compression, shear, and symmetric clamping: they revealed a series of boundary conditions that can generate stable configurations. Sano and Wada (2018) instead explored a range of controlled asymmetrical boundary conditions (hinged-clamped) that generate snapping.

From the literature survey, there appears to be a considerable amount of research in confined thin sheets, especially in determining their stability. Very little yet exists on confined folded sheet. By fold, in this paper, we mean a sharp corner in the deformation (Lobkovsky, 1996; Lechenault et al., 2014).

Despite the abundance of solutions of the Euler's Elastica (Bigoni, 2012; Manning, 2014), there seems to be a lack of closed-form (or even numerical) solutions of the Elastica containing folds as we intend in this paper. Related works are (Dado et al., 2004) and (Phungpaingam and Chucheepsakul, 2018), where they studied the post-buckling response of a cantilevered column (fixed-free ends) with a rotational spring in the middle, subjected to a concentrated horizontal force at the free end. However, there is no specific connection to folds, nor the spring contains a residual dihedral angle.

A sharp corner translates into a discontinuity of the rotation of the strip. Using a numerical meshfree method initially developed for cracks, our recently published work (Barbieri et al., 2019) presented a discretization method to model thick plates with multiple folds with infinite rotational stiffness.

In this paper, we are particularly interested in how constraining a thin strip with a fold with finite rotational stiffness leads to non-zero curvatures. Also, we study how such curvatures can be tuned by modifying the applied rotations at the boundaries (Figure 1). The question of whether such shapes are stable or unstable remains a topic that requires further investigation.

To examine which boundary conditions create curvatures in folded strips, we carried out semi-analytical simulations. Assuming invariance in the width, we modeled the strip as a plate under cylindrical bending. Therefore, instead of the Föppl-von Kármán equations, we modeled the strip as a planar, linearly elastic, unshearable, and inextensible rod, according to the Euler's theory of the Elastica.

The equilibrium equations for the Elastica are (Landau and Lifshitz, 1959; Frisch-Fay, 1962)

with boundary conditions

where (·)′ = d/dS, with 0 ≤ S ≤ L being the curvilinear coordinate, μ is the curvature, θ is the rotation, x and y the Cartesian coordinates of the deformation, L the length of the rod, B = EI/(1 − ν²) the bending stiffness, E the Young modulus, I the second moment of area of the cross-section, ν the Poisson ratio, H₀ the applied horizontal force at S = 0, V₀ the applied vertical force at S = 0, M₀ the bending moment at S = 0, and θ₀ the rotation at S = 0, with positions x₀, y₀ at S = 0 (Figure 2A).

We pass to the dimensionless form using the following normalizations:

resulting in

and the following system of non-linear ordinary differential equations [the (·̄) is removed for ease of notation]

with boundary conditions

where (·)′=d/dS̄, with 0≤S̄≤1. Following standard techniques (Bigoni, 2012) of solution for Equation (5), we get

with θ_F0 such that

and

where am is the Jacobi amplitude function and

with

with F is the elliptic integral of the first kind and

The Cartesian coordinates of the deformation are

with

where dn is the Jacobi dn function

where

where E is the elliptic integral of the second kind, I_max = I(ϕ_max) and nϕ′↑↓ is the number of times from 0 to S when ϕ′ changes sign. Furthermore,

and

The strain energy of the rod is

Let us now assume that a fold exists at 0 < S_f < 1. The fold is a discontinuity 〚·〛 on the rotation, with the following constitutive model (Lechenault et al., 2014)

with Sf+ the curvilinear coordinate immediately to the right of S_f and Sf- the one immediately on the left (Figure 2B). The quantity κ > 0 is the rotational stiffness of the fold, and M(S_f) is the bending moment, assumed continuous in the absence of concentrated moments applied at S_f. The angle φ₀ is the rest opening angle. When φ₀ = 0, there is no fold; when φ₀ < 0 the fold is called a mountain and when φ₀ > 0 the fold is called a valley (Figure 1). The linear relationship between bending moment and dihedral angle was showed experimentally in Lechenault et al. (2014) and Jules et al. (2019) to hold for a wide range of dihedral angles. In their experiments, Lechenault and co-workers used the setup in Figure 3. The sample is clamped at one end, with a dynamometer attached to the other end. Both ends lie in the same plane. The right end can translate as the force pulls the specimen. The sample is illuminated by a laser, and the successive deformations are captured by a digital camera. In this way, it is possible to measure the dihedral angle β₀. The moment at the fold is calculated as the product of the force measured by the dynamometer and the maximum height of the profile y(S_f).

The geometrical meaning of φ₀ is shown in Figure 2B. The rest opening angle is related to the dihedral angle β₀ by the simple relation

The rest opening angle is a mechanical property related to the yield stress of the material of the strip (Lechenault et al., 2014): to create a sharp fold, one needs to apply an extreme deformation that generates localized and irreversible strains. Such strains, in turn, lead to high stress concentrations beyond the yield stress of the material. An estimate of the dihedral angle is (Lechenault et al., 2014)

where ϵ_p is the plastic strain, σ_Y the yield stress and E the Young modulus of the material composing the strip. In addition, φ₀ appears to be independent from the thickness of the strip (Lechenault et al., 2014).

To tune the curvatures, we will further assume that at least one of the following quantities is assigned:

which renders one or all of H₀, V₀, M₀ unknowns to be determined. The procedure to compute such unknowns is a shooting method (Press et al., 2007). Firstly, the solutions (7), (15), (16), and (22) are applied from 0 ≤ S ≤ S_f: let us call this solution μ⁻, θ⁻, x⁻, and y⁻. Then, the values θ(1), x(1), and y(1) are computed by applying the solution (7), (15), and (16) with the following boundary conditions:

the procedure is iterated with a non-linear solver until θ(1), x(1), and y(1) computed in this way match the conditions (27), within a certain tolerance.

In this section we report the results of the Elastica model for strips containing folds. We examine both isostatic and hyperstatic boundary conditions. In all the numerical calculations, we set the origin of the Cartesian system to the origin of the rod, therefore x₀ = y₀ = 0.

In this section, we derive a simplified analytical model to explain the numerical results obtained in section 3.

A thin sheet containing a fold assumes a roof-like undeformed state when confined isostatically. Instead, if confined hyperstatically, it will possess one or more curvatures. Such curvatures change according to the rotations applied at the ends of the sheet. Using the Euler's Elastica theory, we calculated the values of the applied rotations corresponding to a switch in the curvatures' signs.

We considered a fold as a discontinuity in the rotations, equal to the supplementary angle to its dihedral angle. The mechanical behavior of the fold depends upon two material properties: the rest opening angle and the rotational stiffness. The rest opening angle is a material property closely related to the yield stress; the rotational stiffness contributes to increasing the dihedral angle of the fold when subjected to a bending moment.

For infinite rotational stiffness, the dihedral angle depends only on the rest opening angle. For the signs of the curvatures before and after the fold, we can identify three regions: one where both curvatures are positive, one where they are both negative and one when they are opposite in sign. We presented a map of such occurrences according to the applied rotations; for small rotations and displacements, the numerical results agree with a simplified model based on the Euler-Bernoulli beam. Also, this map confirmed the presence of values of applied rotations for which the curvature is uniform throughout the sheet; in particular, there is one couple of applied rotations when such uniform curvature is zero. These rotations correspond to the isostatic confinement.

For folds with finite rotational stiffness, the two regions with same-sign curvatures shrink. Interestingly, there exists a critical value of the rotational stiffness when the fold disappears, meaning that the discontinuity on the rotation becomes zero. In particular, when the applied rotations are opposite in signs and the fold is a mountain, the critical value of the rotational stiffness is different from zero, and represents a transition from mountain to valley. In the limit of zero rotational stiffness, the sheet transforms into an undeformed valley with rest angle equal to twice the applied rotation at the end.

All datasets generated for this study are included in the manuscript and/or the supplementary files.

EB designed the research, conducted the analysis, and wrote the manuscript.