# Faceted wrinkling by contracting a curved boundary

Anshuman S. Pal, Luka Pocivavsek, Thomas A. Witten  James Franck Institute,
Dept. of Physics
University of Chicago

###### Abstract

Single-mode deformations of two-dimensional materials, such as the Miura-ori
fold, are important to the design of deployable structures because of their
robustness, but usually require careful pre-patterning of the material. Here,
we show that inward contraction of a curved boundary produces a novel single-
mode deformation without any pre-patterning. Using finite-element simulations
of the contraction of a thin circular annular sheet, we show that these sheets
wrinkle into a structure with negligible stretching energy, in which the
contracted boundary forms spontaneous facets. We construct a strictly
isometric wrinkled surface formed of triangles and cones that matches
geometric and energy features closely, suggesting that this class of partly-
faceted wrinkled deformations is isometric. Isometry favours the restriction
of such deformations to a robust low-bending energy channel that avoids
stretching. This class of buckling also offers a novel way to manipulate sheet
morphology via boundary forces. Finally, it serves as a minimal model for
illustrating the strong constraints imposed by geometry in elastic pattern
formation.

††preprint: APS/123-QED

##  I Introduction

The flourishing field of extreme mechanics [1, 2, 3, 4] seeks ways in which
large, spontaneous material deformations yield robust, reproducible structures
and motions. Many of these structures exploit the selective deformations shown
by thin sheets, which bend or fold easily but resist internal stretching
strongly [5, 6, 7, 8, 9, 10, 11, 12, 13]. The large resistance to stretching
creates strong, non-local constraints that restrict the sheet to selected
modes of deformation. A classic example is the Miura-ori folding pattern [14,
15, 16]. This periodic pattern enables to contract a spread-out flat sheet of
paper into a compactly folded state with one single mode of deformation
dictated by the pattern (see Fig. 1a). Here we report similar robust
deformation of a flat sheet without the position-dependent processing of an
origami fold. Instead, we specify an initial annular shape and a simple mode
of edge deformation: contracting the inner boundary radially. The result is a
strongly reproducible three-dimensional pattern of high-amplitude wrinkles. We
attribute the reproducibility to a qualitative feature of the simulated
deformation – it involves virtually no stretching; in the limit of zero
thickness, it is “isometric”. This unstretched morphology is made possible by
a distinct pattern of flat, triangular facets alternating with smoothly-
curving conical segments.

Figure 1: Geometry-controlled robust radial wrinkling. a) Miura-ori patterns
guide a flat sheet to contract into a folded state with a single, continuous
motion (adapted from [16]). b) Smooth radial (i.e., Lamé) wrinkling in an
annulus under both inner and outer tension, indicated by arrows (adapted from
[17]). In this paper, we study Miura-ori like robust deformation of an annulus
in a radial Lamé setting. This deformation can be generated in multiple ways.
c) Finite-element simulation under uniform radial contraction of the inner
boundary, described in Methods (inset shows zoom-in of emergent zigzag inner
boundary). d) Cone-triangle geometric construction defined in Sec III, with
triangles shown in blue and cones in red. e) Ruff collar garment (reproduced
from [18]) made by manually pleating the inner boundary of an annular piece of
cloth.

#####

This paper studies a variant of the classic Lamé radial wrinkling deformation
(see Fig. 1b), in which an annular sheet is subjected to uniform tensile
forcing at the inner and outer boundaries [19, 17, 20, 21]. Our “inner Lamé”
variant imposes an inward radial contraction at the inner boundary only.
Below, we show how the requirement of isometricity quantitatively accounts for
observed features of our simulated inner Lamé annulus, in both its shape and
its elastic energy. Specifically, we show that the simulated deformation
matches an explicitly isometric deformation of the original annulus into a
cone-triangle pattern. We describe the robustness of the structure using a
variety of boundary conditions for the contracted inner boundary. We argue
that the observed robustness of the cone-triangle structure is a natural
consequence of its isometric property. We note how the cone-triangle structure
can be induced via different paths, as illustrated in Figs. 1d, e, and f, and
discuss generalisations to a broader range of structures.

#####

Sections II and III define the inner Lamé deformation that we simulate, and
the isometric model shape used for comparison. In Section IV, we show that the
shape and energy of the simulated sheet are consistent with the model.
Finally, Section V addresses our assumptins and discusses the value of such
deformations for generating spontaneous structures.

##  II Inner Lamé deformation

Given the circular annulus of Fig. 2a, we define the inner Lamé deformation as
the equilibrium shape induced by drawing the inner boundary inward by a
distance $\Delta$, so that it is forced to live on a cylinder whose radius is
reduced by $\Delta$. We define the initial inner radius as our unit of length.
The pull distance $\Delta$ thus acts as the sole loading parameter in the
system. The system can thus be conveniently defined using only three geometric
parameters: thickness $t$, the width $w$, and displacement ‘loading’ $\Delta$.
The internal forces determining the shape arise from the bending modulus $B$,
and the in-plane stretching modulus $Y$. In thin sheets with thickness $t$
much smaller than the inner radius and the width $w$, the sheet may be taken
as virtually unstretchable, and many features of the shape are independent of
$Y$ [17].

Figure 2: a) Sketch of the inner Lamé deformation of a horizontal annulus.
(Top) Top view showing the undeformed annulus, of inner radius $1$ and width
$w$. (Bottom) Vertical cross-section of the deformed annulus through a line
where the vertical displacement of the wrinkle is maximal. The undeformed
state is shown in light shading; the dark lines are a section of the deformed
annulus, contracted horizontally by a distance $\Delta$. These lines are also
displaced upward and tilted downward. b) Sketch of the cone-triangle
construction in its initial, flat state. Here, we choose wavenumber $m=6$,
hence it has 12 triangles (shown in colour). The mid-line and height $w$ of
one triangle are indicated. The white sectors bridging the triangles are shown
in white. When the bases of the triangles are moved inward, they must tilt to
retain their shape and connectivity. The bridging sectors deform into cones. A
typical deformed example is shown in Fig. 1b.

#####

We may anticipate that this inner Lamé deformation is strictly isometric by
using the post-buckling far-from-threshold (FT) analysis of ref. [17] for the
classic Lamé problem, where tension is applied at both the inner and outer
boundaries. In the FT regime, the wrinkling slopes grow to the order of unity,
and the associated azimuthal stress becomes negligible compared to the pre-
buckling stress. Only the bending stress, owing to wrinkling and proportional
to the bending modulus, remains. The same is expected for our inner Lamé case:
consider the radial stress in a narrow sector of the annulus. Equilibrium
requires that the net radial force on the sector be zero. In the unbuckled
state, the inner tension driving the deformation is balanced by the net
outward component of the azimuthal stresses on the sides of the sector. But in
the FT regime, this stress is negligible, as argued above. In addition, there
is no radial tension at the outer boundary. We conclude that the only forces
to balance the inner tension are the weak bending stresses.

#####

The absence of any deformation stress except for bending stress would imply
that the deformed sheet deforms only by bending in the thin limit. That is, it
is isometric to the undeformed flat annulus. The limiting shape must thus be
‘developable’ (c.f. Gauss’ Theorema Egregium; [22]). This means there must be
a straight ‘generator’ line through each point of the surface that extends to
the boundary without curving in space. In the next sections, we show that the
contracted annulus is indeed consistent with a specific developable surface,
suggested by the numerically determined shape.

##  III A developable cone-triangle shape

The faceted wrinkling seen in our numerical simulations – made using the
commercial software Abaqus (SIMULIA, Dassault Systèmes) – clearly shows the
presence of straight directors at the zigzag inner boundary (see Fig. 1c,
inset). Moreover, it suggests a further hypothesis: that the wrinkled
morphology can be seen as a piecewise union of alternating flat triangles and
curved cones. We take this as our starting point. For a configuration with
wavenumber $m$, we consider a flat configuration (see Fig. 2b) where the inner
boundary is a $2m$-sided polygon. Each segment of this boundary is the base of
an isosceles triangle of height $w$. To complete the annulus, the apexes of
adjacent triangles are joined with an arc. The constructed surface then
consists of these triangles and the sectors joining them.

#####

One may readily specify an isometric deformation that allows the inner
boundary of this structure to contract. We deform the initial surface by
translating the base of each triangle inward toward the centre, specifying
that the triangles remain rigid. In order to maintain their shape and
connectivity, the triangles must tilt (i.e. rotate about their mid-lines)
increasingly as they are drawn inward; adjacent triangles tilt in opposite
directions. The sectors connecting two triangles must then also bend in order
to remain connected to their triangles. A ready isometric choice for this
deformation is to bend the flat sectors into circular conical arcs. These arcs
need not, however, join smoothly with their adjacent triangles, as further
discussed in Sec. IV.3.

##  IV Comparing simulated deformation with model

In this section, we test whether the isometric/developable cone-triangle
construction is able to predict the quantitative features of our simulated
inner Lamé deformation from Sec. II. Our commercial Abaqus software uses
standard finite-element methods to represent the deformed equilibrium shape as
the inner boundary is gradually contracted to move along the axis of the
cylinder (as detailed in Methods). In comparing this simulation to the cone-
triangle construction, we must supply two free parameters: the boundary
displacement $\Delta$, as well as the wavenumber $m$. Thus, we measure the
$\Delta$ and $m$ of the simulated deformation, and then construct the
corresponding geometric solution with the same $\Delta$ and $m$. The first
comparison we make is purely qualitative – we superpose the two solutions and
observe how well they match up. Fig. 3 shows such a comparison, where we have
rotated the geometric solution about the $z$-axis to obtain the best
superposition. As the figure shows, one may align the two figures so that they
overlay well over a significant part of the annulus, with matching zigzags at
the inner boundary and matching waves at the outer boundary. Between the inner
and outer boundaries, the surfaces match closely enough that the two surfaces
alternately hide each other.

#####

The qualitative resemblance seen in the figure encourages a more quantitative
comparison. Accurate comparison is hindered by the variability of the
simulated wrinkles. Though regular in general, the wrinkle pattern is not
strictly periodic on the annulus. To optimise the comparison, we perturb the
initial state of the simulations to promote periodicity. Specifically, we
deform the initial state with a low-amplitude, smooth, periodic wrinkle with a
wavelength similar to the original simulations. 111The initial wrinkle shape
chosen is an eigenmode obtained from linear stability analysis performed in
Abaqus.

Figure 3: Qualitative comparison between Abaqus simulation and our cone-
triangle model, for a representative sample ($w=1.67$, $t=1.33\times
10^{-3}$). (Left) Given $w$, $\Delta$ and the wavenumber $m$ from the Abaqus
solution, we construct a model cone-triangle solution (cones in red, triangles
in blue) using the procedure described in Section III. (Right) For a
qualitative comparison, we superpose the two solutions, with azimuthal
rotation as the only degree of freedom. We find close agreement between them;
with this top-down view, areas with colour indicate where the model solution
is above the Abaqus solution, while those in grey indicate the opposite. In
particular, the zoomed-in boxed portion shows close matching at both inner and
outer boundaries. Zooming out to the whole annulus, we see that there is a
‘phase mismatch’ between the two solutions, with part of the solution in-phase
(the boxed region) and part of it out-of-phase. This is due to the non-uniform
size of the wrinkles in the Abaqus solution.

###  IV.1 Inextensibility and bending strains

We first verify the anticipated isometry of the sheet as announced in the
Introduction. To this end, we consider a material circle at radius $r$ in the
annulus; its arc length $L(r)$ is necessarily $2\pi r$, and independent of the
displacement $\Delta$. However, since the position of this circle is displaced
inward by an amount $\Delta$, the horizontal projection of the circle has a
length $L_{xy}(r)=2\pi(r-\Delta)$ 222Here we treat the radial lines in the
material as nearly horizontal, as they are for large wavenumber $m$.. Thus the
relative change in $L_{xy}$ on deformation is given by

$(L_{xy}-L)/L=-\Delta/r.$ (1)

The “contraction factor” $\Delta/r\equiv\tilde{\Delta}(r)$ is useful for the
comparisons below.

#####

We test the validity of 1 by explicitly measuring $L_{xy}$ for the inner and
outer boundary circles in the simulation. In the simulation, these boundaries
consist of a sequence of straight segments, one for each finite element. We
add up the segment lengths and their $xy$ planar projections to obtain
$L_{xy}/L$ as a function of $\Delta$. Fig. 4b shows representative data
testing Eq. 1 at both inner ($r=1$) and outer ($r=1+w$) boundaries. The data
points (left side of Eq. 1) seem to line up perfectly with the linear
prediction (right side of Eq. 1). However, closer inspection shows a tiny
deviation between two two. This is because Eq. 1 requires exact isometry,
whereas the simulation cannot be an exact isometry due to its non-zero
thickness. The difference $L_{xy}/L-1+\tilde{\Delta}$ between prediction and
data in Eq. 1 is thus caused by a residual azimuthal strain, which we denote
$\epsilon_{\theta\theta}^{\rm resid}$. At the inner boundary, this strain is
caused by localised stretching at the cone tips, and can thus be neglected. At
the outer boundary, however, this strain should result from the isometric
bending stress described in Section II; measuring it thus provides a further
test of isometry.

#####

We know that Elastica bending at scale $\lambda$ induces a residual bending
stress, $\sigma_{\theta\theta}^{\rm resid}\sim B/\lambda^{2}$, where $B$ is
the bending modulus. Thus, for wrinkling of wavenumber $m$, we expect the
outer boundary to have an azimuthal bending strain given by:

$\epsilon_{\theta\theta}^{\rm
resid}\sim\frac{1}{Y}\frac{B}{\lambda^{2}}\sim\frac{m^{2}t^{2}}{(1+w)^{2}},$
(2)

where we have used $\lambda=2\pi(1+w)/m$ and $B/Y\sim t^{2}$ ($t$ is the
thickness). The inset of Fig. 4b is a log-log plot of
$\epsilon_{\theta\theta}^{\rm resid}$ versus $\frac{m^{2}t^{2}}{(1+w)^{2}}$,
at fixed $\Delta$, for multiple samples with differing values of $m$, $t$ and
$w$. We get data collapse over the range of a decade, with the predicted slope
of unity. Thus, Fig. 4b provides strong evidence for inextensibility of the
numerical solutions at both boundaries.

Figure 4: Quantitative comparisons between theory and Abaqus simulations
(using the explicitly periodic simulations described at the beginning of Sec.
IV). (a) Sketches of the wrinkle profile at the inner ($r=1$)and outer
($r=1+w$) boundaries, showing the definition of slope angles $\alpha_{\rm in}$
and $\alpha_{\rm out}$. (b) Testing inextensibility at inner and outer
boundaries (see Sec. IV.1). The main graph shows a typical example of testing
Eq. 1 (for $w=0.33$, $m=10$, $t=1.33\times 10^{-3}$); Abaqus data points
capture the left side of the equation while the black theory line is the right
side. (inset) The residual strain at $r=1+w$, $\epsilon_{\theta\theta}^{\rm
resid}\equiv L_{xy}/L-1+\tilde{\Delta}(1+w)$, plotted logarithmically against
the bending strain $m^{2}t^{2}/(1+w)^{2}$ at fixed $\Delta$, for different
values of $m$, $t$ and $w$ (see Eq. (2)). Legend: $w=0.33$ (blue circles),
$w=0.67$ (red diamonds), $w=1.0$ (yellow squares). The black line shows the
linear fit: $\epsilon_{\theta\theta}^{\rm
resid}=10.8\frac{m^{2}t^{2}}{(1+w)^{2}}+10^{-5}$. (c) Comparing wrinkle slope
angle vs. applied contraction $\tilde{\Delta}$ (see Sec. IV.2). (Top row) We
separately test Eq. (3) for $\alpha_{\rm in}$ (in red) and Eq. (4) for
$\alpha_{\rm out}$ (in black). The Abaqus data points are plotted alongside
the model predictions. Abaqus samples of differing width $w$ and wavenumber
$m$ are coloured using the dimensionless cone splay parameter,
$\eta\equiv\frac{\pi(1+w)}{mw}$. The error bars are standard deviations from
averaging over all the wrinkles. Markers: $w=0.33$ (circles), $w=0.67$
(squares), $w=1.0$ (diamonds), $w=1.67$ (triangles). (Bottom row) Example
wrinkle profiles for small and large $\eta$: (left) for $\eta=0.25$ ($w=1.0$,
$m=25$), and (right) for $\eta=1.8$ ($w=0.33$, $m=7$). Colouring denotes
height. (d) Comparing elastic energies (see Sec. IV.3): $U_{\rm elastic}^{\rm
Abaqus}(\Delta)$ (in red) with $U_{\rm elastic}^{\rm model}(\Delta)$ (in blue)
and $U_{\rm\langle\kappa^{2}\rangle}^{\rm Abaqus}(\Delta)$ (in yellow). (Left)
Annulus, indicating the region where energies were compared (in red), omitting
the high-curvature regions near $r=1$. (Right) We plot the results for a
representative sample ($w=1.0$, $m=25$, $t=6.67\times 10^{-4}$). All three
curves show Elastica-like linear scaling (solid lines are fits), but the two
theoretical curves overestimate the Abaqus energy by roughly $15\%$ and
$30\%$. Other $w$ and $m$ values also generally over-estimated the energies,
as discussed in Sec. IV.3.

###  IV.2 Evolution of wrinkle amplitude

We now compare specific features of the wrinkle shape to the predictions of
the cone-triangle construction. In Fig.4c, we consider the amplitude of the
observed wrinkles. We measure this amplitude by determining the average slope
angle $\alpha$ of the wrinkles at the inner and outer boundaries. Though the
cone-triangle construction can be implemented for a wide range of wavenumbers
$m$ for any given $\Delta$ and $w$, we considered only the narrow range of $m$
that appeared spontaneously in the unbiased Abaqus simulations. In practice,
these had roughly $7\leq m\leq 35$.

#####

We define the average slope angle $\alpha$ to be the slope angle of the line
joining the zero-crossing of the wrinkle with the adjacent maximum
(illustrated in Fig. 4a). For angle $\alpha_{\rm in}$ at the inner boundary,
the figure explicitly relates $\alpha_{\rm in}$ to $L_{xy}/L$. Using Eq. 1, we
have:

$\cos\alpha_{\rm in}=1-\Delta\equiv 1-\tilde{\Delta}(1).$ (3)

At the outer edge, $L_{xy}/L$ is the ratio of the horizontal line to the arc
length of the circular sector in Fig. 4a. For our circular cone model, it is
related to the slope angle $\alpha_{\rm out}$ by

$\sin 2\alpha_{\rm out}/(2\alpha_{\rm out})=L_{xy}/L\equiv
1-\tilde{\Delta}(1+w).$ (4)

In Fig. 4c, we compare these cone-triangle slope angles to the simulated ones
as a function of $\tilde{\Delta}(r)$ for both inner ($r=1$) and outer
($r=1+w$) boundaries. The upper red curve and data points compare $\alpha_{\rm
in}$. The agreement of the data points with the curve confirms that each
triangle segment tilts to preserve its length as it is drawn towards the axis.
The solid black curve gives the cone-triangle prediction of Eq. 4 for
$\alpha_{\rm out}$. The data points on and beneath the curve represent annuli
of different widths $w$ and wavenumber $m$, coloured by their splay:
$\eta\equiv\pi(1+w)/(mw)$. In general, our circular cone model is consistent
with the prediction for narrowly-splayed wrinkles, but over-estimates the
slopes for widely-splayed ones. The pictures below the graph show the annulus
shape for different amounts of splay $\eta$. They suggest the source of the
over-estimated slopes: for the widely splayed example with the greatest
discrepancy, the boundary curve is significantly flatter than a circular arc.
Thus the simplifying assumptions of our triangle-cone construction are not
satisfied in this widely-splayed regime.

###  IV.3 Comparing energies

Perhaps the main significance of this isometric deformation is that its energy
behaves so differently from that of a generically deformed sheet. It is
different in two respects. Firstly, the energy depends only on the bending
modulus $B$, and is independent of the stretching modulus $Y$ for given
bending modulus. Secondly, the energy in this limit is qualitatively smaller
than for a generic deformation with comparable curvature and thickness. Below
we verify that the energy of the contracted annulus has these distinctive
features.

#####

The elastic energy is a sum of strain energy and bending energy:

$U_{\rm elastic}=U_{\rm strain}+U_{\rm bend}$

Thus, for a strain-free isometric configuration, we have $U_{\rm
elastic}\approx U_{\rm bend}$. For a conical deformation extending between
$L_{\rm min}$ and $L_{\rm max}$ – as measured along the cone axis – the
bending energy is given by:

$U_{\rm bend}^{\rm cone}=\frac{B}{2}\log(\frac{L_{\rm max}}{L_{\rm min}})\int
ds\;\kappa^{2},$ (5)

where $\kappa(s)$ is the normal curvature of the outermost arc of the cone,
parameterised by the arc length $s$. For a circular cone as in our model,
$\kappa(s)=\rm{const.}=1/\rho_{c}$, where $\rho_{c}$ is the radius of the
outermost circle. Then, adding up the contribution of $2m$ cones for a
wrinkled solution of wavenumber $m$ (see SI for details), we get:

$U_{\rm elastic}^{\rm model}(\Delta)=U_{\rm bend}^{\rm model}(\Delta)\approx
B\log(\frac{L_{\rm max}}{L_{\rm min}})\frac{\pi(1+w)w}{\rho_{c}^{2}(\Delta)},$
(6)

where the approximation is for large-$m$ solutions. Here, the bending modulus
$B=\frac{Et^{3}}{12(1-\nu^{2})}$ (with $E$ the Young’s modulus, $\nu$ the
Poisson ratio, and $t$ the thickness) can be calculated directly from the
material properties of the Abaqus solution. Also, in order to avoid high-
curvature regions near $r=1$ where the thin-sheet approximation $\kappa t\ll
1$ becomes questionable, we choose to compare energies for only the outer
$3/4^{th}$ part of the annulus, so that $\log(\frac{L_{\rm max}}{L_{\rm
min}})=\log(4)$.

#####

In Eq. (6), the only dynamic variable is $\rho_{c}(\Delta)$. Our circular cone
construction (see SI) permits us to get an expression for $\rho_{c}(\Delta)$
in terms of the outer slope angle $\alpha_{\rm out}$:

$\rho_{c}(\Delta)\approx\frac{\pi(1+w)}{4m}\frac{1}{\alpha_{\rm
out}(\Delta)}.$ (7)

We recall that $\alpha_{\rm out}$ depends only on $w$ and $\Delta$, and not on
$m$. Thus $\rho_{c}\sim 1/m$, and so $U_{\rm elastic}^{\rm model}\sim m^{2}$,
which is the quadratic scaling with wavenumber (at fixed amplitude) expected
for a bending energy.

#####

Our primary subject of interest, however, is the scaling of $U_{\rm
elastic}^{\rm model}$ with $\Delta$ for fixed $m$. Since $\alpha_{\rm
out}(\Delta)\sim\sqrt{\Delta}$ approximately (see Fig. 4c), we expect the
squared curvature and the bending energy to scale linearly with $\Delta$:
$U_{\rm elastic}^{\rm model}\sim 1/\rho_{c}^{2}\sim\Delta$. We note that this
linear form is the same as for the Elastica energy of a rod under weak Euler
buckling [25], which is expected since each transverse arc of the cone
undergoes Euler buckling with this same $\Delta$ dependence.

#####

In Fig. 4d, we plot three energy curves for a typical sample. These curves
compare the simulated energy for a typical example, $U_{\rm elastic}^{\rm
Abaqus}(\Delta)$ (in red), with two theoretical estimates: $U_{\rm
elastic}^{\rm model}(\Delta)$ (in blue; see Eq. 6) calculated from our cone-
triangle model, and $U_{\rm\langle\kappa^{2}\rangle}^{\rm Abaqus}(\Delta)$ (in
yellow) calculated indirectly from the Abaqus solution through measurement of
the mean squared curvature at the outer boundary:
$\langle\kappa^{2}\rangle_{\rm out}\equiv\frac{1}{2\pi(1+w)}\int
ds\;\kappa^{2}$, and by using Eq. 5. We see that the two theoretical estimates
are larger than the simulation energy by a constant factor of order unity.
This pattern was seen in all our simulations. The extra energy seen from the
circular-cone model (blue curve) is a consequence of that very approximation –
our circular cones are constructed to fit into the narrowed sector angle of
the contracted sheet without compressing or stretching. However, the actual
bent sectors must satisfy a further condition, even in the asymptotic limit of
zero strain, viz. continuity of the first derivative. Thus, each sector must
meet the adjacent flat triangle smoothly, with no change in the normal vector.
In the regime we studied, the slope-angle at the cone-triangle boundary is too
large to match the adjacent triangle, giving excess bending energy in the
cone. The observed discrepancy in this slope is compatible with the
discrepancy in bending energy.

#####

The $U_{\rm\langle\kappa^{2}\rangle}^{\rm Abaqus}(\Delta)$ curve infers the
conical bending energy using the (non-constant) observed surface curvature at
the simulated outer boundary. The observed curvature was smaller than the
circular arc at the cone-triangle boundary and greater at the mid-line. This
estimate gave an energy much closer than the simulated energy, but still not
in agreement with it. This indicates that the surface is not of the assumed
conical shape throughout the curved surface measured. The surface could attain
a smaller energy while matching the outer boundary shape if its directors did
not converge to a point but instead converged more weakly. We have not settled
whether this remaining discrepancy is owing to a) discretisation error in the
simulation, or b) a residual effect of nonzero sheet thickness, which would
become smaller with decreasing thickness.

#####

We conclude that the elastic energy of deformation for the wrinkles
encountered in this study is strongly consistent with our claim of an
isometric conical deformation, with the accuracy increasing for narrow
wrinkles (having small parameter $\eta$). Furthermore, the explicit cone-
triangle construction gives a workable quantitative estimate of the full
elastic energy, provided we supply the observed wrinkle wavenumber $m$.

##  V Discussion

###  V.1 A developable wrinkling morphology

As shown above, our variant of the Lamé radial contraction system is
qualitatively altered by removing the outer tension from that system. The
mechanics in our system is governed instead by the constraints of isometric
deformation, dictating a developable surface in the limit of zero thickness.
While developability has been hypothesised and used as tool to study buckling
in several contexts [26, 27, 28, 29], in our favourable case, we have been
able to account for all the details of the structure with a simple periodic
cone-triangle ansatz.

#####

To the best of our knowledge, such reconstruction of a wrinkled configuration
through fundamental developable building blocks is novel in the field of thin
sheet elasticity. That an a priori complicated wrinkled surface can be broken
down into such a two-part macroscopic unit cell is indicative of the strong
geometric simplifications brought about by developability and Gauss’ Theorema
Egregium. Furthermore, unlike previously-known faceted structures – which are
composed of only flat facets [30, 31, 32] – faceted wrinkling has both flat
triangles and curved cones. This cone-triangle decomposition also provides an
instructive example of a buckled shape in a radially symmetric system that
does not have the conventional separable form – i.e. the buckling displacement
is not a product of a radial factor and azimuthal factor.

#####

The developability of the faceted wrinkling structure makes it fundamentally
geometric, so that it can be created easily on a computer. This provides a way
to programmatically construct smooth radially wrinkled patterns in materials
using intrinsically flat components, much in the spirit of origami and other
designed materials [4]. In fact, a well-known real world example of such
programmed radial wrinkling is the Elizabethan ruff collar. One way of making
a ruff is exactly using the inner Lamé boundary conditions: we take a cloth
annulus, and manually make pleats in it such that the inner boundary of the
annulus lives on a contracted radius. Figs. 1e and 1d show, resp. , a ruff
collar and its reconstruction using our cone-triangle model; the visual
similarity between the two is striking. The conical outer solution seems to
perfectly capture the ruff wrinkles, in particular their radial straightness
and their shape at the outer boundary.

###  V.2 Wavelength selection

As is clear from the ruff example, the faceting and developability discussed
in this paper is completely independent of the number of wrinkles (viz. the
wavenumber $m$). That is because these features are set independently.
Wavelength determination is controlled by the minimisation of the bending
energy, whereas we are concerned in this paper almost entirely by the
minimisation of the strain energy. This is in keeping with the standard
paradigm of understanding wrinkled morphologies in thin sheets [17], where
separation of energy scales means that the gross features of the morphology
are determined by minimising the dominant strain energy, while the finer
features (like wavelength) are controlled by the sub-dominant bending energy.
Wavenumber determination will be discussed in a future paper. For now, we note
that $m$ seems to grow without bound as the sheet becomes thinner, as in the
tensile classical Lamé case. Yet the reason must be entirely different, since
the tensile case is governed by non-bending work, which is absent in the
present system.

###  V.3 Assumptions

####  V.3.1 Circular cone ansatz

The key modelling assumption used in this paper is the circular ansatz for the
bent cones. As we have shown, this ansatz serves as a good approximation to
the true energy-minimising Elastica solution, while being simple enough to
allow for explicit calculation of several key observables, viz. the slope
$\alpha_{\rm out}$ and the energy $U_{\rm elastic}$. Sec. IV shows that this
approximation is accurate only for relatively narrow cones (small $\eta$) and
large $w$. This can be partially explained as follows. We expect the circular
ansatz to hold true only for weak Euler buckling. ‘Weak’ here means that the
applied contraction, and the resultant slope, are both small. In the SI, we
show that the contraction experienced by a conical sector increases sharply
with decreasing $w$, such that for our narrowest annuli, the contraction can
be as high as $15-20\%$. At such large values, one no longer expects the weak
Euler regime to be valid, and the resultant shape manifests its non-linearity
through the flattening apparent in Fig. 4c. We note, however, that even at
such large contractions, both the circular ansatz and the numerical
simulations give a linear energy profile, although their slopes are no longer
close.

####  V.3.2 Inner boundary approximation

In the physical problem of Sec.II, the inner boundary is constrained to lie on
a cylinder. Our cone-triangle surface clearly does not satisfy this condition
exactly; its inner boundary consists of straight lines, and does not match the
cylindrical surface exactly. Instead, the cone-triangle surface solves a
closely related problem in which the cylinder is replaced by a polygonal
cylinder with $2m$ vertical facets. For the thin sheets we studied
numerically, the number of wrinkles $m$ grows without bound for very thin
sheets, so that the mismatch between the real cylinder and the faceted
cylinder becomes arbitrarily small. In the numerical samples, we looked for
strains along the inner boundary owing to this inner boundary approximation.
We detected no change in the length of the inner zigzag line as the
contraction $\Delta$ increases.

#####

Still, the inexactness of the triangle construction at the inner surface makes
it clear that the construction is not a realistic representation of a circular
annulus contracted onto a cylinder for arbitrary wavenumber $m$. Likewise, we
have no evidence that an annulus constrained to have a given small $m$ remains
developable under contraction. Our only evidence for this developable shape
applies to the range of $m$ that the system spontaneously chooses, or larger.

####  V.3.3 Lateral constraint boundary condition

A further remark concerns the boundary conditions used in our Abaqus
simulations. As noted in Methods, we prohibit any lateral movement of the
nodes on the inner boundary, in order to stabilise the faceted wrinkled regime
against folding instabilities. However, if we do allow lateral movement, we
see the same faceted wrinkling at early times, before azimuthal movement of
the boundary quickly leads to folding instabilities. This suggests that
faceting can be seen for weak buckling, even in the absence of lateral
constraints. But our evidence is still only suggestive and needs to be
confirmed through physical experiments.

###  V.4 Role of thickness

The role of the thickness $t$ (or its absence) in our work also merits some
attention. Our developable cone-triangle construction lives in the asymptotic
zero-thickness limit, and thus thicker sheets might be expected to deviate
from it significantly. However, our analytical construction accurately models
the wrinkled morphology for even moderately thick sheets (our thickest
simulated sheet had $t=2.67\times 10^{-2}$). In fact, the example of the ruff
collar (Fig. 1e) suggests that the cone-triangle asymptotic solution is valid
practically even for macroscopically thick sheets. This suggests that
deviations from our asymptotic solution are restricted to the localised
boundary layers of stress in the otherwise isometric solution – the areas of
so-called ‘stress concentration’ [30]. These are zones where where the radius
of curvature becomes comparable to the thickness, leading to divergent bending
energy density. To prevent this, such zones are ‘smoothened out’ at the cost
of some localised strain; in our wrinkled solution, we see such point-like
boundary layers at the tips of the cones. These boundary layers are expected
to shrink to zero radius as a sheet approaches the zero-thickness limit.
However, since such boundary layers are localised and have negligible effect
on the overall solution (see SI), we did not investigate them any further in
the present study.

###  V.5 Generalisation

Can our model be generalised? The key modification applied in this paper to
the traditional Lamé problem is the freeing up of the outer boundary. This
allows the radial tension in the sheet to vanish, and also allows the outer
boundary to deform freely. In the absence of any tension or other non-elastic
sources of work [33, 34], and following our reasoning, we should see analogous
faceted wrinkling with any contraction of the inner curved boundary of a
sheet. The isosceles triangles and right-circular cones associated with the
axisymmetric Lamé geometry will have to be correspondingly generalised to
arbitrarily shaped triangles and cones, but the rest should follow as
described in this paper. Separately, one can also imagine attaching the inner
boundary to a constraining cone instead of a cylinder. This will impose an
angle instead of a displacement at the inner boundary, thereby breaking the
up-down symmetry of Lamé wrinkling. Thus, our work can act as a starting point
for analysing an entire class of boundary-controlled wrinkling morphologies,
whose prototype is an annulus contracted at its inner boundary.

###  V.6 Application

The deformation geometry treated here amounts to a strategy for avoiding the
typical fate of deformed sheets. Such sheets adopt structures that require
significant elastic energy beyond the needed bending energy. Often these
structures produce uncontrolled disorder, e.g. , as seen in a crumpled sheet.
By contrast, the radial contraction used here guides the sheet into a low-
energy channel resulting in only bending energy. This low-energy state is
moreover attainable by continuous deformation from a flat reference state.
Since these deformations lack the additional deformation costs of competing
structures, they are intrinsically stable. They provide a general means to
fold and compact a sheet continuously and reversibly with minimal external
guidance. We note that the curved inner boundary appears important in
stabilising this structure; straight boundaries with no radial contraction [9]
do not produce analogous regular wrinkling.

###  V.7 Conclusion

While the emergence of developability in wrinkled systems has been studied
intensively over the past decade, the simplicity and symmetry of the Lamé
geometry makes it a particularly clean realisation of such a phenomenon. In
fact, the inner Lamé system may be considered a direct 2d equivalent to the
prototypical 1d Euler Elastica system. Both cases consider pure bending of a
system whose stretching energy has been set to zero. The faceted, developable
morphology of the Lamé case in particular highlights the strong geometric
constraints entailed by isometry in 2d. This faceted morphology is
qualitatively different from previously known forms of 2d wrinkling.
Significantly, it creates sharply dictated, strongly distorted features
without stretching and without a separable radial/azimuthal shape. Moreover,
it has many variants to be explored, with real world applications. There
remain also fundamental mysteries like the wavelength selection to be
explained; this problem will be the focus of future work.

##  VI Acknowledgements

The authors would like to thank Benny Davidovitch, Dominic Vella and Enrique
Cerda for many insightful discussions, and Nhung Nguyen and George
Papazafeiropoulos for help with numerical analysis. This work was partially
funded by the National Science Foundation (NSF)-Materials Research Science and
Engineering Center (MRSEC) grants DMR-1420709 and DMR-2011854.

##  VII SI: Construction of circular conical sections

In the main text, we describe how we divide the flat annulus into two regions
– an inner region composed of isosceles triangles, and an outer region
composed of arc sectors. The main text describes how the inner triangles can
isometrically deform under an inner displacement (or contraction) of length
$\Delta$ by tilting about the normal to the inner boundary – this gives an
inner ‘buckled’ solution. Here, we describe how to construct the outer
solution whereby the arc sectors isometrically bend into sections of right-
circular cones. In other words, each arc bends into a circular profile.

Figure 5: Isometrically constructing a sector of a cone from a flat arc
sector. (a) The flat configuration. (b) The curved configuration; the actual
conical sector and its variables are depicted in black, the triangle
$\bigtriangleup AOB$ formed by the edges of the two adjacent tilted triangles,
along with its variables, is in blue, and the putative cone with its variables
is in grey.

#####

The requisite variables for our geometric construction are described in Fig.
5. In the flat state, each arc sector, lying between two adjacent flat
triangles, subtends an angle $\phi_{f}$ at the common vertex of these two
triangles. If the edges of the triangles have length $R_{f}$, then the length
of the largest arc is $L=R_{f}\phi_{f}$. On tilting, the edges of these same
triangles approach each other, so that they subtend a reduced angle $\phi_{t}$
(see SI Fig. 6a). This reduction in angle ($\Delta\phi=\phi_{t}-\phi_{f}\leq
0$) constitutes the effective compressive strain $\epsilon_{\rm
cone}=\Delta\phi/\phi_{f}$ for the arc sector. SI Fig. 6a shows the evolution
of $\epsilon_{\rm cone}$ with boundary contraction $\Delta$. We see clearly
that narrower sheets feel a greater $\epsilon_{\rm cone}$ for the same
$\Delta$.

#####

The excess material can be accommodated by bending the sector into a cone. The
exact shape of this cone can be determined by solving the Elastica equation
with the appropriate boundary conditions. However, for simplicity, in this
paper we choose to approximate the exact shape by a circular cone. The aim
then is to find the circular cone that fits in between the two adjacent
triangles while exactly accommodating the excess angle $\Delta\phi$. Fig. 5b
shows all the variables needed for this operation. The isometrically bent
conical section has to fit inside the triangle $\bigtriangleup AOB$ (drawn in
teal) defined by the apex $0$ and the tilted edges $OA$ and $OB$ (see also
Fig. 3 in main text). The conical section itself is drawn in bold black. Let
the circular cone have angular extent $\phi_{c}$, and let its maximum radius
be $\rho_{c}$, centred at point $C$. Then the conical shape can be determined
from the following set of geometric constraints.

#####

The first constraint is that of inextensibility, i.e. length conservation.
This gives us:

$L=\rho_{c}\phi_{c}=R_{f}\phi_{f}.$ (8)

The second constraint involves specifying the end-to-end distance $L_{1}$
between points A and B. Thus, using the blue triangle $\bigtriangleup AOB$ and
the grey triangle $\bigtriangleup ACB$ in Fig. 5b, we have:

$L_{1}=2\rho_{c}\sin(\phi_{c}/2)=2R_{f}\sin(\phi_{t}/2).$ (9)

Eqns. (8) and (9) constitute a set of simultaneous equations for the cone
variables $\rho_{c}$ and $\phi_{c}$, since $R_{f}$ and $\phi_{f}$ are fixed by
the flat geometry, while $\phi_{t}$ is set by the tilt of the neighbouring
angles. Finally, using basic trigonometry, we have for the cone’s apex angle:

$\sin\beta=\rho_{c}/R_{f}.$ (10)

In sum, we have three independent equations for three unknowns: $\phi_{c}$,
$\rho_{c}$ and $\beta$ (all positive-definite), which fully specify the cone
in space (including the centre $C$). SI Fig. 6b shows a representative
solution of such a constructed cone. SI Fig. 6c shows the evolution of this
construction with increasing $\Delta$.

Figure 6: Constructing a single cone (see Sec. III and SI for details). (a)
Two adjacent flat triangles (in yellow) subtending an angle $\phi_{f}$, on
tilting (in deep blue), subtend a smaller angle $\phi_{t}$. Thus, the region
between the triangles (shown approximately as a violet triangle) gets squeezed
by an angle $\Delta\phi=\phi_{t}-\phi_{f}$. (inset) Evolution of the conical
strain $\epsilon_{\rm cone}=\Delta\phi/\phi_{f}$ with $\Delta$ for three
different widths $w$. Model measurements are plotted alongside the predicted
lines of slope $=-1/(1+w)$. (b) To fit within the contracted boundary defined
by $OA$ and $OB$, the region must bend into a circular cone of radius
$\rho_{c}$, angular extent $\phi_{c}$, axial length $L_{c}$ and tilt angle
$\beta$. The only ingredients needed for the construction are $\triangle AOB$
and the contraction $\Delta$. (c) The evolution of the conical solution with
$\Delta$ (here, $0.05\leq\Delta\leq 0.6$). Greater opacity means larger
$\Delta$; the cones are translated with respect to each other for viewing
clarity. (inset) As $\Delta$ increases, the radius of curvature $\rho_{c}$
decreases, and for $\Delta\to 1$, attains an asymptotic value of
$\rho_{c}\approx\frac{\pi w}{2m}$ (dotted line).

###  VII.1 Correcting for scalloped arc sectors

The flat arc sectors as defined above – circular arcs with their origin at the
inner boundary – lead to the creation of a flat shape that is in fact a
scalloped annulus, whose outer circumference $2mR_{f}\phi_{f}$ is greater than
the expected $2\pi(1+w)$.

We can correct for this length discrepancy when constructing the bent cones,
by changing $R_{f}\to\chi R_{f}$ in Eq. (8) where
$\chi=\frac{\pi(1+w)}{mR_{f}\phi_{f}}\leq 1$. Eq. (9) remains unchanged since
it specifies the end-to-end-distance, which is set by the neighbouring tilted
triangles. We note that this factor $\chi$ is smallest for samples with small
$w$ and $m$ (wide splay $\eta$), and becomes $\approx 1$ for large $w$ and $m$
(narrow splay).

###  VII.2 Limits of the cone-triangle construction

To the best of our knowledge, the above conical construction works as long as
the initial (flat) angle $\phi_{f}<\pi$, i.e. as long as the edges of two
adjacent triangles define a triangle. The value of $\phi_{f}$ depends on the
width $w$ and the wavenumber $m$, and increases as $m$ decreases. This defines
a minimum wavenumber $m_{\rm min}^{\rm cone}(w)$ below which a conical
solution is invalid. For a sufficiently wide annulus, we find that $m_{\rm
min}^{\rm cone}=2$, which is the minimum possible value for any wrinkled
solution. However, for very narrow annuli, this value goes up. Thus, for
$w=0.2$, we find $m_{\rm min}^{\rm cone}=4$. Fig. 7 shows two such contrasting
geometries. $m_{\rm min}^{\rm cone}(w)$ defines a geometric limit beyond which
we expect our cone-triangle model to fail. However, the main text shows that
we already see significant deviations from our model for Abaqus solutions with
wavenumber $m$ significantly higher than $m_{\rm min}^{\rm cone}(w)$.

Figure 7: Two different flat state geometries for $w=0.2$. (Left) For $m=4$,
the flat angle $\phi_{f}<\pi$, which allows the conical construction described
in Sec. VII.2. (Right) For $m=2$, the flat angle $\phi_{f}>\pi$, which means
that our conical construction is not valid here. Figure 8: Diagram describing
the geometric problem for determining $\phi_{t}$. (Left) Top view of the flat
isosceles triangle $\bigtriangleup PQR$ lying in the x-y plane, having height
$w$ and base $h$. $OQ$ is the perpendicular bisector of $PR$. $\phi_{f}$ is
the angle between $QR$ and the dotted (vertical) blue plane. (Right) End view
of the same triangle rotated clockwise by angle $\alpha$ about the y-axis.

##  VIII SI: Calculating the conical strain $\epsilon_{\rm cone}(\Delta)$

Fig. 6a shows the evolution of the conical strain $\epsilon_{\rm cone}$ with
contraction $\Delta$. The observations (plotted with dots) show that
$\epsilon_{\rm cone}(\Delta)$ is linear, with the slope depending
significantly on the width $w$ of the sheet, but independent of the wavenumber
$m$. To understand this behaviour, we analytically derive an expression for
$\epsilon_{\rm cone}(\Delta)=(\phi_{t}-\phi_{f})/\phi_{f}$ .

#####

Since $\phi_{f}$ is fixed, to find $\epsilon_{\rm cone}$, we only need to
calculate $\phi_{t}(\Delta)$. The angle $\phi_{t}(\Delta)$ can be calculated
given any one of the inner triangles, the flat arc angle $\phi_{f}$, and the
tilt angle $\alpha$. The geometric problem is defined in Fig. 8. Here,
$\phi_{t}$ is the angle between the rotated vector
$\overrightarrow{QR^{\prime}}$ and the vertical blue plane defined by
$\overrightarrow{RZ}$. Also, $w$ is the annular width, and $h=\pi/m$. Then our
steps are as follows:

  1. 1.

Find the 3-vector
$\overrightarrow{OR^{\prime}}=\mathbb{R}_{-\alpha}\cdot\overrightarrow{OR}$,
where $\mathbb{R}_{\alpha}$ is the 2x2 rotation matrix of angle $\alpha$ about
the y-axis (viz. $\overrightarrow{OQ}$).

  2. 2.

Then
$\overrightarrow{R^{\prime}Q}=\overrightarrow{OQ}-\overrightarrow{OR^{\prime}}=(-h/2\cos\alpha,w,h/2\sin\alpha)$.

  3. 3.

We only needs its orthogonal projection on to the x-y plane: the 2-vector
$\overrightarrow{R^{\prime}Q_{\perp}}=(-h/2\cos{\alpha},w)$.

  4. 4.

We need the vector angle between this and the vertical blue plane defined by
$\overrightarrow{RZ}$. Here,
$\overrightarrow{RZ}=\mathbb{R}_{-\phi_{f}/2}\cdot\overrightarrow{RQ}$, where
$\overrightarrow{RQ}=(-h/2,w)$ .

Finally, we have
$\phi_{t}=2\mathrm{VectorAngle[}\overrightarrow{R^{\prime}Q_{\perp}},\overrightarrow{RZ}\rm{]}$.
The final expression is complex and composed of trigonometric functions. But
using $\cos\alpha=1-\Delta$ (recall that the inner radius $=1$), and expanding
to linear order in $\Delta$, we have:

$\Delta\phi\equiv\phi_{t}-\phi_{f}\approx-\frac{h/w}{(h/2w)^{2}+1}\Delta+\dots,$
(11)

where $h/w=\frac{\pi}{mw}\ll 1$ for most of our samples. Finally, we also have
$\phi_{f}\approx\frac{\pi(1+w)}{mw}$, where $1+w$ is the outer radius.
Therefore, we get for the dimensionless conical contraction:

$\epsilon_{\rm
cone}\equiv\frac{\Delta\phi}{\phi_{f}}\approx\frac{h}{w\phi_{f}}\Delta\approx-\frac{1}{1+w}\Delta.$
(12)

The prediction of Eq. 12 is plotted alongside measurements in the inset of SI
Fig. 6a. The two are fully consistent. This explains why wider cones should
feel less contraction for a given $\Delta$.

##  IX SI: Calculating bending energy for a cone

In this section, we continue with the variables introduced in section SI-I. A
right-circular cone is made up of a series of circles whose radii $\rho$
increase linearly with distance (say, $z$) along the cone axis, measured from
the cone tip. For a cone of axial length $L_{c}$ (see SI Fig. 6a) and maximum
radius $\rho_{c}$, we thus have: $\rho(z)=\frac{\rho_{c}}{L_{c}}z$. We can
also write $\frac{\rho_{c}}{L_{c}}=\tan\beta$, where $\beta$ is the cone’s
vertex angle.

#####

Now, an infinitesimally wide circular band of radius $\rho$, angular extent
$\phi_{c}$, and width $\frac{dz}{\cos\beta}$, has bending energy
$=B\frac{1}{\rho^{2}}\times(\phi_{c}\rho)\times\frac{dz}{\cos\beta}$, where
$B$ is the bending modulus. Integrating over an entire conical sector, we get:

$\displaystyle U_{\rm bend}^{\rm cone}$
$\displaystyle=B\frac{\phi_{c}}{\cos\beta}\int_{0}^{L_{c}}dz\;\frac{1}{\rho(z)}$
(13)
$\displaystyle=B\phi_{c}\frac{1}{\cos\beta\tan\beta}\int_{0}^{L_{c}}\frac{dz}{z}$
(14)

Any physical cone must have a finite bending energy, and thus must depart from
the conical shape near $z=0$. Thus, the conical shape can only be realistic
beyond a minimum axial distance $L_{\rm core}$. For this outer region, the
energy becomes:

$U_{\rm bend}^{\rm cone}=B\frac{\phi_{c}}{\sin\beta}\log(\frac{L_{c}}{L_{\rm
core}}),$ (15)

where $\beta=\tan^{-1}\frac{\rho_{c}}{L_{c}}$. For small angle $\beta$ (i.e.
large $\Delta$), Eq. (15) can be simplified to:

$U_{\rm bend}^{\rm cone}\approx
B\phi_{c}\frac{L_{c}}{\rho_{c}}\log(\frac{L_{c}}{L_{\rm core}}),$ (16)

Here, $\phi_{c}$, $\rho_{c}$ and $L_{c}$ are all functions of the contraction
$\Delta$. But $\phi_{c}$ and $\rho_{c}$ are related through the constraint:
$\phi_{c}\rho_{c}=L=\frac{\pi(1+w)}{m}$ (see Eq. 8), and $L_{c}\approx w$ (the
width of the annulus), so $U_{\rm bend}^{\rm cone}$ can be reduced to a
function of the single variable $\rho_{c}$. For the entire annulus, we need to
multiply this by the wavenumber $m$. Thus, we get:

$U_{\rm bend}^{\rm theory}\approx
B\frac{\pi(1+w)w}{\rho_{c}^{2}}\log(\frac{w}{L_{\rm core}}),$ (17)

More generally, for a cone extending between axial limits $L_{\rm min}$ and
$L_{\rm max}$, we have:

$U_{\rm bend}^{\rm theory}(\Delta)\approx
B\frac{\pi(1+w)w}{\rho_{c}^{2}(\Delta)}\log(\frac{L_{\rm max}}{L_{\rm min}}),$
(18)

This is the expression given in the main text.

#####

While Eqs. (17) and (18) seem independent of $m$, it is not so. Eq. (7) in the
main text shows that $\rho_{c}\sim(1+w)/m\implies U_{\rm bend}^{\rm
theory}\sim m^{2}$, as expected of a bending energy. We note that this
$\rho_{c}\sim 1/m$ scaling could have been predicted in another way. In the
limit of maximum contraction ($\Delta\to 1$), the point C approaches the xy-
plane, and so the diameters of the $2m$ circles must approximately equal the
reduced projected outer perimeter: $2\pi(1+w-\Delta)\to 2\pi w$. This gives
us: $4m\rho_{c}\to 2\pi w\implies\rho_{c}\to\frac{\pi w}{2m}$ (the dotted line
in the inset of SI Fig. 6c). Indeed, even if C is off the x-y plane (e.g. for
smaller $\Delta)$, the diameters of the $2m$ circles must still equal the
reduced perimeter $2\pi(1+w-\Delta)$ up to some factor. Thus, we have the
scaling relation $\rho_{c}\sim 1/m$ as expected.

#####

Finally, we note that the data presented in the paper represents the full
expression 15, without any approximation.

##  X SI: Estimating the strain energy concentrated at the cone tips

The elastic energy of an isometric system is given purely by the bending
energy: $\rm U_{elastic}=U_{bend}$. However, as mentioned in the main text,
the isometric solution is not valid at the tips of the conical sections
located on the inner boundary. The expected divergence in bending energy at
these points is regularised by localised stretching - the so-called zones of
stress (or strain) concentration (see [30]). In some systems (e.g. [28, 35,
36]), such concentrated stress energy is found to be comparable to the bending
energy $\rm U_{bend}$. While section IV.3 (alongside Fig.4d) of the main text
convincingly argues that the elastic energy of our isometric annulus is given
entirely by the bending energy, we would nevertheless like to independently
verify that the strain energy localised at cone tips is negligible compared to
$\rm U_{bend}$.

#####

For this, we estimate the strain energy concentrated at the cone tips. For a
single cone tip, this is given by $U_{\rm strain}\sim Y\epsilon_{\rm
core}^{2}R_{\rm core}^{2}$, where $Y=Et$ is the stretching modulus ($E$ is the
Young’s modulus and $t$ the thickness), $R_{\rm core}$ is the approximate size
of the strain localisation zone (the ’core region’), and $\epsilon_{\rm core}$
is the localised strain. Through measurements of the stress and strain from
the Abaqus samples, we estimate $\epsilon_{\rm core}\sim 10^{-3}$ and $R_{\rm
core}\sim 10^{-2}$ (compare to system size of roughly 1). Furthermore, we find
that these two variables are independent of system size, and only depend
(weakly) on the thickness $t$. To estimate $U_{\rm strain}$, we consider two
thicknesses. For the thinner case, $t=1\times 10^{-4}$, we have $Y\sim 100$,
$\epsilon_{\rm core}\sim 10^{-3}$ and $R_{\rm core}\sim 10^{-2}$. This gives
$U_{\rm strain}\sim 10^{-8}$, which is entirely negligible with respect to the
bending energy ($\rm U_{bend}\sim 10^{-4}-10^{-3}$). Again, for $t=1\times
10^{-3}$, we have $Y\sim 1000$, $\epsilon_{\rm core}\sim 5\times 10^{-3}$ and
$R_{\rm core}\sim 10^{-2}$, giving $U_{\rm strain}\sim 10^{-6}$, which is
again negligible before the bending energy.

##  XI SI: FEA Simulation Methods

For our simulations, we used the commercial finite-element software Abaqus
2018 (Simulia, Dassault-Systemes, Providence, RI). We employed both the
dynamic explicit and dynamic implicit methods in the quasi-static limit;
‘dynamic‘ means that the solution is found by forward integration of the
equations of motion, while ‘explicit‘ and ‘implicit‘ refer to the integration
scheme (resp. they correspond to variations of forward and backward Euler
integration). The mixture of these two methods was done partly as a
consistency check, partly for convenience, and party by necessity. While the
implicit method in the quasi-static limit was faster for most jobs, the
explicit solver was indispensable for the thinnest samples ($t=2\times
10^{-4},1\times 10^{-4}$), where the static solver ran into convergence
problems. Ensuring the quasi-static limit is also easier for dynamic implicit
than for explicit. In dynamic implicit, the quasi-static option is in-built,
but for dynamic explicit, it has to be ensured manually by applying the
loading slowly enough so that that the ratio of kinetic to strain energy stays
strictly below 10%. For the thinnest samples, e.g. this involved applying the
contraction $\Delta$ over a time period $T=120$ (system units). As a result,
the slowest simulations, for the thinnest and widest sheets, lasted $\approx
30$ hours with 4 cores; the average simulation however, lasted between 2-8
hours on 4 cores.

#####

The assembly consisted of only a single annulus, with inner radius fixed and
taken to be unity, and with varying width $w$ and thickness $t$ in order to
test the geometric Lamé system over a wide range of system parameters. For
width, we used values $w=0.20,0.33,0.67,1.0,1.67$ (a factor of roughly 10,
ranging from very narrow to moderately wide), and for thickness, we used
values $t=2.67\times 10^{-3},1.33\times 10^{-3},6.67\times 10^{-4},2.67\times
10^{-4},1.33\times 10^{-4}$ (a factor of 20, ranging from moderately thick to
very thin). While the thickness varies over a decade, the values still fall
well within the thin sheet limit. The annular part was made of 2d shell quad
(S4R) elements. This choice was made mainly to optimise speed, since we used a
very fine mesh to ensure that the results be independent of meshing. For
comparison, the coarsest mesh we used was for the $w=1.67$ annulus, with 60
elements across the radius and 1400 elements across a circle, giving a maximum
linear size for an element $\approx 0.008$. For consistency checks however, we
also ran some simulations with annuli made of 3d cubic (C3D8R) elements (which
run far more slowly).

#####

For the material properties, we used a neo-Hookean hyperelastic model, with
coefficients $\rm C_{10}=1.5375\times 10^{5}$, $\rm D_{1}=3.2520\times
10^{-7}$. These coefficients are related to the more well-known linear
elasticity moduli by the relations: $C_{10}=\rm G/2$ (where $\rm G$ is the
shear modulus), and $\rm D_{1}=2/K$ (where $\rm K$ is the material bulk
modulus). The corresponding Poisson ratio is given by
$\nu=\frac{3/(C_{10}D_{1})-2}{6/(C_{10}D_{1})+2}=0.475$. The choice of modulus
(and Poisson ratio) was made to optimise computational performance for the
dynamic explicit simulations, since the material modulus affects the wave
speed and hence the stable time increment within a dynamic formulation. Too
soft a modulus makes it difficult to attain the quasi-static limit, while a
very large modulus or perfect incompressibility ($\nu\to 0.5$) leads to
divergence in the stable time increment and very long simulation times. For
the implicit method, the choice of modulus is not so important.

#####

To generate wrinkling, we applied ‘$\theta$-fixed’ boundary conditions at the
inner boundary: $u_{r}(r=1)=-\Delta,\;u_{\theta}(r=1)=0$ , while keeping the
outer boundary at $r=1+w$ free. We usually used $\Delta\leq 0.267$. However,
applying such an in-plane loading to a planar sheet is not enough to generate
buckling automatically. To do so, we need to first perturb the flat state with
a small imperfection; for this, we first solved the linear stability prob in
Abaqus with $z$-pinned BC at $r=1$:
$u_{r}(r=1)=-\Delta,u_{\theta}(r=1)=0,u_{z}(r=1)=0$ , and then used a small
amplitude combination of calculated eigenmodes as a perturbation. This is
standard procedure for doing post-buckling analysis in Abaqus.

—————————————————

## References

  * Krieger [2012] K. Krieger, Extreme mechanics: Buckling down, Nature 488, 146 (2012).
  * Bertoldi  _et al._ [2017] K. Bertoldi, V. Vitelli, J. Christensen, and M. van Hecke, Flexible mechanical metamaterials, Nature Reviews Materials 2, 10.1038/natrevmats.2017.66 (2017).
  * Santangelo [2017] C. D. Santangelo, Extreme mechanics: Self-folding origami, Annual Review of Condensed Matter Physics 8, 165 (2017).
  * Callens and Zadpoor [2018] S. J. Callens and A. A. Zadpoor, From flat sheets to curved geometries: Origami and kirigami approaches, Materials Today 21, 241 (2018).
  * Cerda and Mahadevan [2003] E. Cerda and L. Mahadevan, Geometry and Physics of Wrinkling, Physical Review Letters 90, 074302 (2003).
  * Cerda  _et al._ [2004] E. Cerda, L. Mahadevan, and J. M. Pasini, The elements of draping, Proceedings of the National Academy of Sciences 101, 1806 (2004).
  * Mora and Boudaoud [2006] T. Mora and A. Boudaoud, Buckling of swelling gels, Eur. Phys. J. E 20, 119 (2006).
  * Klein  _et al._ [2007] Y. Klein, E. Efrati, and E. Sharon, Shaping of elastic sheets by prescription of non-euclidean metrics, Science 315, 1116 (2007).
  * Vandeparre  _et al._ [2011] H. Vandeparre, M. Piñeirua, F. Brau, B. Roman, J. Bico, C. Gay, W. Bao, C. N. Lau, P. M. Reis, and P. Damman, Wrinkling Hierarchy in Constrained Thin Sheets from Suspended Graphene to Curtains, Physical Review Letters 106, 224301 (2011).
  * Hure  _et al._ [2012] J. Hure, B. Roman, and J. Bico, Stamping and wrinkling of elastic plates, Physical Review Letters 109, 054302 (2012).
  * Vella  _et al._ [2015] D. Vella, J. Huang, N. Menon, T. P. Russell, and B. Davidovitch, Indentation of Ultrathin Elastic Films and the Emergence of Asymptotic Isometry, Physical Review Letters 114, 014301 (2015).
  * Davidovitch  _et al._ [2019] B. Davidovitch, Y. Sun, and G. M. Grason, Geometrically incompatible confinement of solids, Proceedings of the National Academy of Sciences 116, 1483 (2019).
  * Paulsen [2019] J. D. Paulsen, Wrapping liquids, solids, and gases in thin sheets, Annual Review of Condensed Matter Physics 10, 431 (2019).
  * MIURA [1985] K. MIURA, Method of packaging and deployment of large membranes in space, The Institute of Space and Astronautical Science report , 1 (1985).
  * Mahadevan and Rica [2005] L. Mahadevan and S. Rica, Self-organized origami, Science 307, 1740 (2005).
  * Schenk and Guest [2013] M. Schenk and S. D. Guest, Geometry of miura-folded metamaterials, Proceedings of the National Academy of Sciences 110, 3276 (2013).
  * Davidovitch  _et al._ [2011] B. Davidovitch, R. D. Schroll, D. Vella, M. Adda-Bedia, and E. A. Cerda, Prototypical model for tensional wrinkling in thin sheets, Proceedings of the National Academy of Sciences 108, 18227 (2011).
  * [18] Kristian schuller photography, https://www.trendhunter.com/trends/hdr-photography.
  * Huang  _et al._ [2007] J. Huang, M. Juszkiewicz, W. H. de Jeu, E. Cerda, T. Emrick, N. Menon, and T. P. Russell, Capillary wrinkling of floating thin polymer films, Science 317, 650 (2007).
  * Piñeirua  _et al._ [2013] M. Piñeirua, N. Tanaka, B. Roman, and J. Bico, Capillary buckling of a floating annulus, Soft Matter 9, 10985 (2013).
  * Paulsen  _et al._ [2017] J. D. Paulsen, V. Démery, K. B. Toga, Z. Qiu, T. P. Russell, B. Davidovitch, and N. Menon, Geometry-Driven Folding of a Floating Annular Sheet, Physical Review Letters 118, 048004 (2017).
  * do Carmo [1976] M. P. do Carmo,  _Differential Geometry of Curves & Surfaces_ (Prentice Hall, 1976).
  * Note [1] The initial wrinkle shape chosen is an eigenmode obtained from linear stability analysis performed in Abaqus.
  * Note [2] Here we treat the radial lines in the material as nearly horizontal, as they are for large wavenumber $m$.
  * Audoly and Pomeau [2010] B. Audoly and Y. Pomeau,  _Elasticity and geometry_ (OUP Oxford, 2010).
  * Cerda and Mahadevan [1998] E. Cerda and L. Mahadevan, Conical Surfaces and Crescent Singularities in Crumpled Sheets, Physical Review Letters 80, 2358 (1998).
  * Ben Amar and Pomeau [1997] M. Ben Amar and Y. Pomeau, Crumpled paper, Proceedings of the Royal Society of London. Series A: Mathematical, Physical and Engineering Sciences 453, 729 (1997).
  * Lobkovsky  _et al._ [1995] A. Lobkovsky, S. Gentges, H. Li, D. Morse, and T. A. Witten, Scaling Properties of Stretching Ridges in a Crumpled Elastic Sheet, Science 270, 1482 (1995).
  * Dias  _et al._ [2012] M. A. Dias, L. H. Dudte, L. Mahadevan, and C. D. Santangelo, Geometric Mechanics of Curved Crease Origami, Physical Review Letters 109, 114301 (2012).
  * Witten [2007] T. A. Witten, Stress focusing in elastic sheets, Reviews of Modern Physics 79, 643 (2007).
  * Seffen and Stott [2014] K. A. Seffen and S. V. Stott, Surface texturing through cylinder buckling, Journal of Applied Mechanics 81, 10.1115/1.4026331 (2014).
  * Pham Dinh  _et al._ [2016] H. Pham Dinh, V. Démery, B. Davidovitch, F. Brau, and P. Damman, From Cylindrical to Stretching Ridges and Wrinkles in Twisted Ribbons, Physical Review Letters 117, 104301 (2016).
  * Huang  _et al._ [2010] J. Huang, B. Davidovitch, C. D. Santangelo, T. P. Russell, and N. Menon, Smooth Cascade of Wrinkles at the Edge of a Floating Elastic Film, Physical Review Letters 105, 038302 (2010).
  * Paulsen  _et al._ [2016] J. D. Paulsen, E. Hohlfeld, H. King, J. Huang, Z. Qiu, T. P. Russell, N. Menon, D. Vella, and B. Davidovitch, Curvature-induced stiffness and the spatial variation of wavelength in wrinkled sheets, Proceedings of the National Academy of Sciences 113, 1144 (2016).
  * Kim  _et al._ [2012a] J. Kim, J. A. Hanna, R. C. Hayward, and C. D. Santangelo, Thermally responsive rolling of thin gel strips with discrete variations in swelling, Soft Matter 8, 2375 (2012a).
  * Kim  _et al._ [2012b] J. Kim, J. A. Hanna, M. Byun, C. D. Santangelo, and R. C. Hayward, Designing responsive buckled surfaces by halftone gel lithography, Science 335, 1201 (2012b).


# Multimodal Emotion Recognition from Raw Audio with Sinc-convolution

###### Abstract

Speech Emotion Recognition (SER) is still a complex task for computers with
average recall rates usually about 70% on the most realistic datasets. Most
SER systems use hand-crafted features extracted from audio signal such as
energy, zero crossing rate, spectral information, prosodic, mel frequency
cepstral coefficient (MFCC), and so on. More recently, using raw waveform for
training neural network is becoming an emerging trend. This approach is
advantageous as it eliminates the feature extraction pipeline. Learning from
time-domain signal has shown good results for tasks such as speech
recognition, speaker verification etc. In this paper, we utilize Sinc-
convolution layer, which is an efficient architecture for preprocessing raw
speech waveform for emotion recognition, to extract acoustic features from raw
audio signals followed by a long short-term memory (LSTM). We also incorporate
linguistic features and append a dialogical emotion decoding (DED) strategy.
Our approach achieves a weighted accuracy of 85.1% in four class emotion on
the Interactive Emotional Dyadic Motion Capture (IEMOCAP) dataset.

Index Terms—  multimodal emotion recognition, raw audio signal, Sinc-
convolution layer, feature-level fusion

##  1 INTRODUCTION

Speech emotion recognition (SER) is a major research area. Most of the
acoustic features used in the fusing system are hand-crafted features [1],
such as MFCC, pitch and voice quality features extracted from the raw audio or
features learned automatically from spectrograms [2], [3]. However, the hand-
crafted features contain more limited information than the raw audio signal,
which makes more and more researchers pay attention to extracting features
from raw audio more flexibly. For example, Guizzo et el [4] present a multi-
time-scale (MTS) convolution layer which does not increase the number of
parameters but increases the temporal flexibility compared to standard CNNs.
Xu et el [5] propose a hierarchical grained and feature model (HGFM) which
includes a frame-level representation module with before and after
information, a utterance-level representation module with context information,
and a different level acoustic feature fusing module. In this paper, we append
the sentences of raw audio as acoustic features to our model.

Due to the great success of CNN in the field of speech emotion recognition,
most systems use the convolution layer to process acoustic features. However,
the standard CNN convolution layer has many filter parameters, which makes the
convolution layer only capable of learning the low dimensional and compacted
features [6, 7]. To help the input layer discover more meaningful filters,
Ravanelli el et [8] proposes SincNet, a neural architecture of directly
processing waveform audio. The Sinc-convolution layer of SincNet is not only
faster in convergence speed than a standard CNN but also more computationally
efficient due to the exploitation of filter symmetry. For this reason, we
utilize the Sinc-convolution layer as the input layer to extract acoustic
information from raw audio signals, and we utilize LSTM followed by the Sinc-
convolution layer to extract information between sentences. Compared with
conventional convolution layers, Sinc-convolution layers can learn the
parameters of filters from raw audio signals, have fast convergence and higher
interpretability with a smaller number of parameters. Humans express emotions
through multi-modal ways, which suggests that we can process multi-modal
features by fusing model to improve the accuracy of predictions [9]. At
present, the feature-level fusion of acoustic features and linguistic features
is more commonly used. Pepino et el [2] present different fusing models for
SER by combining acoustic and linguistic features. The result shows that
multi-modal system leads to significant improvements of approximately $16\%$
on Interactive Emotional Dyadic Motion Capture (IEMOCAP) [10]. In this paper,
we study a multimodal emotion recognition model with combined acoustic and
linguistic information and utilize dialogical emotion decoder (DED) [11] to
process the prediction result of pre-trained classifier. DED decodes each
utterance into one of the four emotion categories at inference stage.
Experiments are performed on IEMOCAP dataset. Our approach achieves a weighted
accuracy of $85.1\%$ in four class emotion on IEMOCAP.

##  2 METHODOLOGY

###  2.1 Acoustic Model

To learn more useful information from raw signals, we utilize Sinc-convolution
filter layers [8] to learn custom filter banks tuned for emotion recognition
from speech audio. Neural architecture for processing raw audio samples as
shown in Figure 1.

Fig. 1: Architecture of Sinc-conv layer.

The first layer of a standard CNN performs a set of time-domain convolutions
between the input waveform and some Finite Impulse Response (FIR) filters.
Each convolution is defined as follows:

$y[n]=x[n]*h[n]=\sum\limits_{l=0}^{L-1}x[l]\cdot h[n-l]$ (1)

where $x[n]$ is a chunk of speech signal. $h[n]$ is the filter of length $L$,
and $y[n]$ is the filtered output. In conventional convolution layer, all the
$L$ elements of each filter are learned from data. Conversely, the proposed
Sinc-convolution layer performs the convolution with a pre-defined function
$g$ that depends on fewer learnable parameters $\theta$ only [8]

$y[n]=x[n]*g[n,\theta]$ (2)

Each defined filter-bank is composed of rectangular band-pass filters. In the
frequency domain, the magnitude of a generic band-pass filter can be
represented by two low-pass filters.

$G[f,f_{1},f_{2}]={\rm rect}\left(\frac{f}{2f_{2}}\right)-{\rm
rect}\left(\frac{f}{2f_{1}}\right)$ (3)

where $f_{1},f_{2}$ refers to low and high cutoff frequencies and $\rm rect$
function is the rectangular function in the magnitude frequency domain. The
time-domain representation of the function $g$ can be derived as follows

$g[n,f_{1},f_{2}]=2f_{2}{\rm sinc}(2\pi f_{2}n)-2f_{1}{\rm sinc}(2\pi f_{1}n)$
(4)

where the $\rm sinc$ function is defined as ${\rm sinc}(x)=\sin(x)/x$. We
compare Sinc-convolution layers followed by Deep Neural Network (DNN) and
LSTM. We also compare conventional convolution layers of CNN fed by
spectrogram computed from raw audio and Sinc-convolution layers of Sinc-DNN
and Sinc-LSTM fed by 250ms chunk selected from raw audio randomly. Though CNN
is based on the same architecture as Sinc-DNN, it replaces the sinc-based
convolution with a standard one. Both Sinc-convolution layers and convolution
layers of CNN use Batch Normalization. The performance is evaluated with
Sentence Error Rate (SER), which represents the average of error rate of each
sentence.

The learning curves of Sinc-DNN and CNN can be observed from Figure 2. Both
convolution layer of Sinc-DNN and standard CNN followed by a Deep Neural
Network (DNN), which contains several fully connected layers and Batch
Normalization. Figure 2 shows that both Sinc-DNN and Sinc-LSTM have faster
convergence than CNN beacuse they reduces the number of parameters in the
first convolutional layer and the function $g$ is symmetric that means we only
need to consider one side of the filter and inherit the result for the other
half.

Fig. 2: Sentence Error Rate of CNN, Sinc-DNN and Sinc-LSTM over various
training epochs

However, the CNN converges to a better performance leading to a Sentence Error
Rate (SER) of 46% against a SER of 51% achieved with the Sinc-DNN. One reason
for this result is that Sinc-DNN does not pay attention to the information
between sentences, which is helpful in the emotion recognition. We learn this
neglected information by replacing DNN with LSTM. Another reason is that some
of the extracted features may be meaningless noise due to our random chunk
extractions from raw audio. To reduce the noise interference caused by
randomly extracting features, we extract the features after pre-processing the
signals in the database, such as extracting the parts with higher power to
effectively avoid the noise interference. As observed in the Figure 2, the SER
of Sinc-LSTM dropped to 43%, which is 8% lower than Sinc-DNN and 3% lower than
CNN.

###  2.2 Linguistic Model

Recurrent Neural Network (RNN) is widely used in natural language processing
tasks. However, RNNs cannot remember longer sentences and sequences due to the
vanishing/exploding gradient problem. It can only remember the parts which it
has just seen. The emergence of LSTM has solved this problem well [12].
Although an LSTM is supposed to capture the long-range dependency better than
the RNN, it is difficult to give large weights to the important words for
estimation. To this end, we deploy a self-attention layer on top of the LSTM.
The self-attention mechanism can focus on the important parts in the sentence
[13]. The more relevant the semantic relation between words and emotion is,
the greater the weight of the connection between them, which can assist the
model to determine more accurate classification.

###  2.3 Feature Fusion

Feature-level fusion is used to fuse linguistic and acoustic features obtained
from individual networks. A 2048-D feature vector from the acoustic network
and a 4800-D feature vector from the linguistic network are concatenated. We
use attention before fusion where attention is applied on individual feature
vectors. Finally, the utterance emotion is classified with the use of a
“Softmax” activation and a dialogical emotion decoder (DED) which is a post-
processing over the final dense layer of the fusing network.

In this paper, we study two fusing models without DED. One is the acoustic
feature of Sinc-DNN fusing with the linguistic feature of LSTM (M1), and the
other is the acoustic feature of Sinc-LSTM fusing with the linguistic feature
of LSTM (M2). We also utilize DED to classify the prediction of M2, which
makes the performance better than M2.

###  2.4 Dialogical Emotion Decoder

The decoding process of Dialogical Emotion Decoder (DED) is built on top of a
SER model which is an approximate inference algorithm which decodes each
utterance into one of the four emotion categories at inference stage. This
decoder is built on three core ideas: the emotion which occurs more frequently
in dialog history is more likely to show up again; while not all utterances
have consensus labels; and the posterior distributions capturing affective
information would enable us to decode utterances in sequence. DED can model
the emotion flows in a dialog consecutively by combining emotion
classification, emotion shift and emotion assignment process together. Yeh et
el [11] achieves $70.1\%$ unweighted accuracy on four emotion class in the
IEMOCAP.

##  3 EXPERIMENTS

###  3.1 Experiment Setup

In this work, we used the IEMOCAP, a benchmark dataset containing about 12
hours of audio and video data, as well as text transcriptions. The dataset
contains five sessions, each of which involves two distinct professional
actors conversing with one another in both scripted and improvised manners. In
this work, we utilize data from both scripted and improvised conversations, as
well as mere audio data to stay consistent with the vast majority of prior
work. We also train and evaluate our model on four emotions: happy, neutral,
angry, and sad, resulting in a total of 5531 utterances (happy: 29.5%,
neutral: 30.8%, angry: 19.9%, sad: 19.5%).

###  3.2 Result and Analysis

We compare the performance of three models: 1) a standard CNN; 2) a Sinc-DNN;
and 3) a Sinc-LSTM.

Table 1: Performance of CNN, Sinc-DNN and Sinc-LSTM. Model | Training time($epoch$) | WA($\%$) | UA($\%$)
---|---|---|---
CNN | $40$ | $54.3$ | $53.2$
Sinc-DNN | $25$ | $49.5$ | $47.8$
Sinc-LSTM | $25$ | $57.1$ | $54.5$

It can be seen from Table 1 that Sinc-LSTM outperforms CNN with a $2.8\%$
improvement that is also obtained with faster convergence (25 vs 40 epochs).
Sinc-convolution layer is specifically designed to implement rectangular
bandpass filters, leading to more meaningful CNN filters, and LSTM can learn
more contextual information than DNN does. This makes Sinc-LSTM perform
better. The performance of multi-modal systems is indicated in Table 2.
CNN+LSTM is the our baseline with acoustic model CNN and linguistic model
LSTM. MDRE [14] using dual RNNs encode acoustic and linguistic information and
predict using a feed-forward neural model. MHA-2 [15] propose a multi-hop
attention model which is designed to process the contextual information from
transcripts and audio. Table 2 shows that MHA-2 has outperformed CNN+LSTM by
$5.7\%$ in the terms of WA when the model is applied to the IEMOCAP dataset
and MHA-2 has outperformed MDRE by $4.7\%$ that demonstrates the attention
framework can increase performance.

Table 2: Performance of different models. Model | Modality | WA($\%$) | UA($\%$)
---|---|---|---
CNN+LSTM | Audio+Text | $70.8$ | $69.7$
MDRE [14] | Audio+Text | $71.8$ | $-$
MHA-2 [15] | Audio+Text | $76.5$ | $77.6$
M1=Sinc-DNN+LSTM | Audio+Text | $72.1$ | $71.1$
M2=Sinc-LSTM+LSTM | Audio+Text | $75.9$ | $76.3$
IAAN [16] | Audio | $65.2$ | $67.1$
IANN+DED [11] | Audio | $69.0$ | $70.1$
M2+DED | Audio+Text | $85.3$ | $85.1$

Fig. 3: The Architecture of fusing model combining Sinc-LSTM and LSTM with DED
post-processing

As observed from Table 2, Both M1 and M2 perform better than CNN+LSTM, which
indicates that the multi-modal model using Sinc-convolution layer not only has
faster convergence than that using CNN, but also has better verification
accuracy. The unweighted accuracy of Sinc-DNN is lower than CNN but M1 has
outperformed CNN+LSTM by 1.3%. The main reason is that the information between
sentences can learned from textual features by LSTM. The M2 system using LSTM
performs better than M1. Specifically, M2 achieves 75% accuracy, with relative
3% improvement over M1, which verifies that LSTM is more suitable for
processing contextual information on IEMOCAP than DNN. Figure 3 shows the
architecture of M2+DED. Based on the M2 system, we use DED to classify the
prediction of M2, and finally achieves an accuracy rate of 85%, which improves
the accuracy rate by 8.8% compared with the MHA-2. However, as shown in Table
2, the effectiveness of IAAN+DED is not as significant as M2+DED. First, DED
relies on a well-performing classifier [11]. Due to the poorer performance of
IAAN in the IEMOCAP, the effect of DED is more limited. Furthermore, as the
IAAN+DED only extract acoustic low-level descriptors (LLD), including features
such as MFCCs, the ability of DED in properly decoding through long sequence
of dialogs is not as obvious as for the raw audio. From Figure 4, We also find
that the performance of DED is more significant when accuracy of the pre-
classifier is higher than $30\%$ because DED classifies the current emotion
based on the previous and the current prediction of the pre-classifier. If the
correct account for a large proportion of all predictions reaches, for
example, $75\%$ of M2, then even if the current is wrong, it can still be
corrected by the previous correct predictions.

Fig. 4: The performance of DED with different pre-classifiers

##  4 CONCLUSIONS

In this paper, we utilize Sinc-convolution layer to extract acoustic features
from raw audio followed by an LSTM and a multi-modal emotion recognition
system based on Sinc-LSTM and LSTM combining acoustic and linguistic data
(Sinc-LSTM+LSTM). The performance of Sinc-LSTM+LSTM is better than CNN+LSTM,
which shows that not only acoustic model but also fusing model using Sinc-
convolution layer has higher verification accuracy. We do not use the hand-
crafted features that have been used in the standard CNN but use the raw audio
append to the Sinc-convolution layer, which makes Sinc-LSTM+LSTM learn more
information than the baseline using CNN and LSTM does. Finally, the prediction
of Sinc-LSTM+LSTM is appended to dialogical emotion decoding (DED) strategy, a
post-processing that classifies the current result based on the previous and
current prediction of pre-classifier. The performance of DED has achieved a
weighted accuracy of 85.29% on the IEMOCAP database, which outperforms the
Sinc-LSTM+LSTM by 9.96%.

## References

  * [1] Christos-Nikolaos Anagnostopoulos, Theodoros lliou, and loannis Giannoukos,  “Features and classifiers for emotion recognition from speech: a survey from 2000 to 2011,”  Artificial Intelligence Review, vol. 43, pp. 155–177, 2015.
  * [2] L. Pepino, P. Riera, L. Ferrer, and A. Gravano,  “Fusion approaches for emotion recognition from speech using acoustic and text-based features,”  in ICASSP 2020 - 2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2020, pp. 6484–6488.
  * [3] Aharon Satt, Shai Rozenberg, and Ron Hoory,  “Efficient emotion recognition from speech using deep learning on spectrograms,”  in Proceedings of the IEEE. IEEE, 2017, pp. 1089–1093.
  * [4] E. Guizzo, T. Weyde, and J. B. Leveson,  “Multi-time-scale convolution for emotion recognition from speech audio signals,”  in ICASSP 2020 - 2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2020, pp. 6489–6493.
  * [5] Y. Xu, H. Xu, and J. Zou,  “Hgfm : A hierarchical grained and feature model for acoustic emotion recognition,”  in ICASSP 2020 - 2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2020, pp. 6499–6503.
  * [6] Xiaohui Zhang, Jiangyan Yi, Chenglong Wang, Chuyuan Zhang, Siding Zeng, and Jianhua Tao,  “What to remember: Self-adaptive continual learning for audio deepfake detection,”  arXiv preprint arXiv:2312.09651, 2023.
  * [7] Xiaohui Zhang, Jiangyan Yi, Jianhua Tao, Chenglong Wang, and Chu Yuan Zhang,  “Do you remember? Overcoming catastrophic forgetting for fake audio detection,”  in Proceedings of the 40th International Conference on Machine Learning, Andreas Krause, Emma Brunskill, Kyunghyun Cho, Barbara Engelhardt, Sivan Sabato, and Jonathan Scarlett, Eds. 23–29 Jul 2023, vol. 202 of Proceedings of Machine Learning Research, pp. 41819–41831, PMLR.
  * [8] M. Ravanelli and Y. Bengio,  “Speaker recognition from raw waveform with sincnet,”  in 2018 IEEE Spoken Language Technology Workshop (SLT). IEEE, 2018, pp. 1021–1028.
  * [9] Xiaohui Zhang, Jaehong Yoon, Mohit Bansal, and Huaxiu Yao,  “Multimodal representation learning by alternating unimodal adaptation,” 2023.
  * [10] C. Busso, M. Bulut, C-C. Lee, Abe Kazemzadeh, Emily Mower, Samuel Kim, Jeannette N.Chang, Sungbok Lee, and Shrikanth S.Narayanan,  “Iemocap: Interactive emotional dyadic motion capture database,”  Language resources and evaluation, vol. 42, no. 4, pp. 335, 2008.
  * [11] Sung-Lin Yeh, Yun-Shao Lin, and Chi-Chun Lee,  “A dialogical emotion decoder for speech emotion recognition in spoken dialog,”  in ICASSP 2020 - 2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). IEEE, 2020, pp. 6479–6483.
  * [12] Q. Huang, R. Chen, X. Zheng, and Z. Dong,  “Deep sentiment representation based on cnn and lstm,”  in 2017 International Conference on Green Informatics (ICGI), 2017, pp. 30–33.
  * [13] J. Xie, B. Chen, X. Gu, F. Liang, and X. Xu,  “Self-attention-based bilstm model for short text fine-grained sentiment classification,”  IEEE Access, vol. 7, pp. 180558–180570, 2019.
  * [14] S. Yoon, S. Byun, S. Dey, and K. Jung,  “Multimodal speech emotion recognition using audio and text,”  2018 IEEE Spoken Language Technology Workshop (SLT), pp. 112–118, 2018.
  * [15] S. Yoon, S. Byun, S. Dey, and K. Jung,  “Speech emotion recognition using multi-hop attention mechanism,”  ICASSP 2019-2019 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), pp. 2822–2826, 2019.
  * [16] Sung-Lin Yeh, Yun-Shao Lin, and Chi-Chun Lee,  “An interaction-aware attention network for speech emotion recognition in spoken dialogs,”  ICASSP 2019-2019 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), pp. 6685–6689, 2019.