justified

# Supplementary Information
Mechanics and wrinkling patterns of pressurized bent tubes

César L. Pastrana  Physics of Complex Biosystems, Technical University of
Munich, 85748 Garching, Germany  Luyi Qiu  John A. Paulson School of
Engineering and Applied Sciences, Harvard University, Cambridge, Massachusetts
02138, USA  John W. Hutchinson  John A. Paulson School of Engineering and
Applied Sciences, Harvard University, Cambridge, Massachusetts 02138, USA
Ariel Amir  John A. Paulson School of Engineering and Applied Sciences,
Harvard University, Cambridge, Massachusetts 02138, USA  Department of Physics
of Complex Systems, Weizmann Institute of Science, 7610001 Rehovot, Israel
Ulrich Gerland  Physics of Complex Biosystems, Technical University of Munich,
85748 Garching, Germany

###### Contents

  1. I Details of the computational model
    1. A Initial configuration
    2. B Model description and energy minimization
  2. II Model
  3. III Estimation of mechanical and geometrical properties
    1. A Curvature, $\kappa$
    2. B Torque, $M$
    3. C Stress resultants, $N_{z}$
    4. D Extent of wrinkles, $\theta_{w}$
    5. E Corrections for the non-linearity of the mesh
  4. IV Mechanics and Geometry of unpressurized tubes
  5. V Finite length dependency
    1. A Non-pressurized tubes
    2. B Pressurized tubes
  6. VI Bending of pressurized cylinder: a tension field theory model
    1. A Direct approach
    2. B Energy approach
  7. VII Bending of pressurized cylinder: non-linear stress-strain relation
  8. VIII Wrinkles geometry from axial in-extensibility conditions
  9. IX Supplementary Figures

###### List of Figures

  1. S1 Estimation of torque in the computational model
  2. S2 Longitudinally aligned meshes
  3. S3 Non-linearity of the mesh
  4. S4 Determination of the extent of the wrinkles
  5. S5 Mechanical response to bending in unpressurized tubes
  6. S6 Wrinkling patterns during bending
  7. S7 Influence of pressure on the wrinkle patterns
  8. S8 Comparison of analytical models for $\kappa$ and $M$
  9. S9 Impact of tube’s length on mechanic parameters (non-pressurized)
  10. S10 Impact of tube’s length on geometric parameters (non-pressurized)
  11. S11 Impact of tube’s length on mechanic parameters (pressurized)
  12. S12 Impact of tube’s length on geometric parameters (pressurized)
  13. S13 Geometry in pressurized tubes during bending
  14. S14 Dependency of the wavelength $\lambda$ of the sinusoidal wrinkles
  15. S15 Expansion of wrinkles $\theta_{w}$ (tension field theory)
  16. S16 $(A/\lambda)^{2}$ in the sinusoidal regime during bending

## Notation

Symbol | Parameter
---|---
$p$ | Pressure
$E$ | 3D Young’s modulus
$\nu$ | Poisson’s ratio
$t$ | Shell thickness
$R_{0}$ | Undeformed radius of the tube
$R$ | Radius of the tube (pressurized)
$L_{0}$ | Undeformed length of the tube
$L$ | Length of the tube
$M$ | Torque
$N$ | Stress resultants
$\varepsilon$ | Strains
$\mathcal{K}$ | Curvature of the medial axis
$\kappa$ | Curvature of the innermost (concave) side of tube
$A$ | Amplitude of the wrinkles
$\lambda$ | Wavelength of the wrinkles
$\theta_{0}$ | Region span by the wrinkles (polar coordinate)

## Supplementary Movies

Supplementary Movie 1 |  Bending of a pressurized tube at $\hat{p}=0.026$. The surface is colored according to Gaussian curvature, $\kappa_{G}$ (Meyer  _et al._ , 2003). The onset of wrinkling occurs at $t=$ 3 s. The first substate of the multi-wavelength wrinkling mode occurs at $t=$ 10 s and the transition to the second substate with a larger wavelength occurs at $t=$ 20 s. Note that the curvature of the tube is not a linear function of time to highlight the different states.
---|---
Supplementary Movie 2 |  Bending of a pressurized tube at $\hat{p}=0.088$. The onset of wrinkling occurs at $t=$ 5 s. The first substate of the multi-wavelength wrinkling mode occurs at $t=$ 13 s and the transition to the second substate with a larger wavelength occurs at $t=$ 19 s.
Supplementary Movie 3 |  Bending of a pressurized tube at $\hat{p}=0.164$. The onset of wrinkling occurs at $t=$ 6 s. The first substate of the multi-wavelength wrinkling mode occurs at $t=$ 11 s and the transition to the second substate with a larger wavelength occurs at $t=$ 15 s.

##  I Details of the computational model

###  A Initial configuration

The main body of the tube is constructed by positioning particles on the
surface of a cylinder describing a perfectly-ordered equilateral triangular
mesh. The tube includes two lids to close the volume. The lids are constructed
by including veritices with progressively reduced radius and organized to have
an hexatic arrangement describing (almost) equilateral triangles. The
geometrical details of the initial configuration are shown in Table S1. The
mesh data is obtained from the vertices’ coordinates by calling to the
Advancing Front Surface Reconstruction routine on the Computational Geometry
Algorithms Library (CGAL).

We used a circumferentially aligned mesh because of its easiness and accuracy
for the determination of both the mechanical and geometrical parameters.
Similar results for the mechanical properties and the features of the wrinkles
are obtained by considering a mesh with a long axis alignment (Figure S2).

###  B Model description and energy minimization

##  II Model

Computational model. The cylindrical shell is modeled as triangulated mesh
with stretching and bending potentials. Torque is induced on the tube by
considering an energetic penalty for the miss-alignment of the normal vectors
of the lids with respect to a target vector.

We study a tube of unpressurized radius $R_{0}$ and length $L_{0}$, described
by a triangulated mesh with two lids that close the volume of the tube. The
total energy of the system is given by
$\mathcal{E}=\mathcal{E}_{\textrm{int}}+\mathcal{E}_{\textrm{ext}}$, where
$\mathcal{E}_{\textrm{int}}$ is the internal (mechanical) energy of the tube
and $\mathcal{E}_{\textrm{ext}}$ is an external potential employed to bent the
tube.

The internal energy of the pressurized shell is given by:

$\mathcal{E}_{\textrm{int}}=\mathcal{E}_{\textrm{b}}+\mathcal{E}_{\textrm{b}}-pV,$
(1)

where $\mathcal{E}_{\textrm{s}}$ is a stretching energy,
$\mathcal{E}_{\textrm{b}}$ is a bending energy penalizing out-of-plane
deformations of the shell, $p$ is the pressure difference and $V$ the volume
enclosed by the shell (Pastrana  _et al._ , 2023). We calculate the volume $V$
enclosed by a non-convex polyhedron using the methods described in (Lien and
Kajiya, 1984).

The stretching energy of the mesh h is described by the harmonic potential:

$\mathcal{E}_{\textrm{s}}=\frac{k_{s}}{2}\sum_{\langle
i,j\rangle}(r_{ij}-r_{ij}^{0})^{2},$ (2)

The sum runs over all pairs of connected vertices $i,j$ and
$r_{ij}=\lVert{\mathbf{r}_{i}-\mathbf{r}_{j}}\rVert$, where $\mathbf{r}$ are
the coordinates of the vertices. The surface bending potential is given by

$\mathcal{E}_{\textrm{b}}=k_{b}\sum_{\langle\alpha,\beta\rangle}(1-\mathbf{\hat{n}}_{\alpha}\cdot\mathbf{\hat{n}}_{\beta}),$
(3)

where the sum is over pairs of triangles $\alpha,\beta$ sharing an edge, and
$\mathbf{\hat{n}}$ are their respective normal vectors.

This model maps to an (equivalent) continuous elastic material with 2D Young
modulus $Y=\frac{2}{\sqrt{3}}k_{s}$ and Poisson ratio $\nu=1/3$ (Seung and
Nelson, 1988). The 2D Young modulus $Y$ is related to the 3D Young modulus $E$
by $Y=Et$, with $t$ the thickness of the shell. The discrete bending stiffness
$k_{b}$ is assigned considering its relation with the continuum flexural
rigidity $B$, given by $B=\frac{\sqrt{3}}{2}k_{b}$ (Gompper and Kroll, 2004;
Seung and Nelson, 1988). The value for $B$ is related to $E$ and $t$ via
(Landau and Lifshitz, 2016):

$B=\frac{Et^{3}}{12(1-\nu^{2})}.$ (4)

The mechanical parameters employed are shown in Table S1. We noticed that the
response of the mesh is only linear for low pressures (low strains). Then, we
correct $E$ and $\nu$ to account for the non-linear response of the mesh. See
Section IIIE for the details of the procedure.

The external potential is
$\mathcal{E}_{\mathrm{ext}}=\mathcal{E}_{\tau}+\mathcal{E}_{\mathrm{CM}}$. The
component $\mathcal{E}_{\tau}$ is included to induce the bending of the tube
and is given by:

$\mathcal{E}_{\tau}=k_{\tau}\sum_{i\in\begin{subarray}{c}\text{left}\\\
\text{lid}\end{subarray}}(1-\mathbf{\hat{n}}_{i}\cdot\mathbf{\hat{n}}^{l}_{0})+k_{\tau}\sum_{i\in\begin{subarray}{c}\text{right}\\\
\text{lid}\end{subarray}}(1-\mathbf{\hat{n}}_{i}\cdot\mathbf{\hat{n}}_{0}^{r})$
(5)

where the summation runs over all the triangles of the right (or left) lid and
$\mathbf{\hat{n}}_{i}$ is the normal vector of the triangle $i$. The term
$\mathbf{\hat{n}}^{l,r}_{0}$ is the preferred orientation vector for the right
or left lid, $r$ and $l$ respectively.

The term $\mathcal{E}_{\mathrm{CM}}$ is used to enforce the alignment of the
two lids in the plane normal to the long axis of the tube

$\mathcal{E}_{\mathrm{CM}}=\frac{k_{\mathrm{CM}}}{2}\lVert{\mathbf{R}_{\mathrm{CM},x}^{l}-\mathbf{R}_{\mathrm{CM},x}^{r}}\rVert^{2}+\frac{k_{\mathrm{CM}}}{2}\lVert{\mathbf{R}_{\mathrm{CM},z}^{l}-\mathbf{R}_{\mathrm{CM},z}^{r}}\rVert^{2},$
(6)

where $\mathbf{R}_{\mathrm{CM},x}^{l,r}$ is the center of mass of the
left/right lid and $x,z$ indicates the Cartesian coordinate. We found that
this potential has no influence in the mechanical and geometrical response.
However, in its absence a single kink is formed anywhere along the distance of
the tube.

The approach described resembles the action of _virtual hands_ to induce
bending. The force exerted on the lids is adjusted with the constant
$k_{\tau}$, given by $k_{\tau}:=K_{\tau}k_{s}l_{0}^{2}$, where $l_{0}$ is the
average bond length of the mesh. $K_{\tau}$ is a free parameter determined
empirically to achieve the target curvatures of the tube. Similarly,
$k_{\mathrm{CM}}:=k_{s}K_{\mathrm{CM}}$. Notably, this procedure to bend the
tube does not impose any constraint along the axial direction. The routine to
bend the tube is as follows. The simulations start by first pressurizing the
tube to the target pressure $p$ by progressively increasing the pressure in
small steps from $p=0$. Next, to induce the bending of the tube, the preferred
orientation vectors of the lids $\mathbf{\hat{n}}^{l,r}_{0}$ are tilted in a
step-wise fashion. For instance,
$\mathbf{\hat{n}}^{r}_{0}=0\mathbf{\hat{i}}-\sin(s\delta\theta)\mathbf{\hat{j}}+\cos(s\delta\theta)\mathbf{\hat{k}}$,
where $s$ is the step number and $\delta\theta$ is a defined change in the
angular position. Each step, pressurization and bending, is followed by energy
minimization.

We use a non-linear conjugate gradient algorithm to find the configuration of
the system with minimum total energy $\mathcal{E}$ using the the secant-method
for line search (Nocedal and Wright, 2006). The tolerance criterion to escape
the minimization loop relies on the average per-vertex force and the precise
value is adjusted empirically for the different pressures. We consider
explicit analytical expressions for all the gradients. The source code is
fully written in C++. We use OpenMP to parallelize the calculations of the
per-vertex gradients as well as vector operations involved in the minimization
routine. This significantly improves the execution run-time for the large
number of vertices employed, a requirement to appropriately characterize the
wrinkles. Importantly, our mesh model is _idealized_ because it does not
include long range interactions. Since the surface can intersect with itself
near the kinks, further bending of the kinked state would not be physical, and
we exclude them from any analysis.

The compiled code, the sources and instructions of use are freely available on
the public repository in (git, 2023).

Parameter | Value
---|---
Tube configuration |
Tube radius (unpressurized), $R_{0}$ | 0.250 $\mu$m
Tube length (unpressurized), $L_{0}$ | 2.0-12.0 $\mu$m
Number of vertices, $N$ | 65,000 - 113,117
Mechanical parameters |
Young modulus (3D), $E$ | 50 MPa
Shell thickness, $t$ | 2.0 nm
Pressure, $p$ | 0.025-0.6 atm
Table S1: Simulation parameters in real units. The mechanical parameters were
selected to be within the ranges observed for Gram negative bacteria (Yao  _et
al._ , 1999; Deng  _et al._ , 2011).

##  III Estimation of mechanical and geometrical properties

###  A Curvature, $\kappa$

We find the curvature at the bottom side of the tube $\kappa$ from the medial
axis curvature $\mathcal{K}$. First, we find a set of particles defining the
medial axis of the tube $\mathbf{x}$ by taking the center of mass of
consecutive sections along the tube. Then, the $\mathbf{x}$ coordinates are
smoothed considering a moving average between consecutive points.

A discrete local tangent vector is given by
$\mathbf{t}_{i}=\mathbf{x}_{i+1}-\mathbf{x}_{i}$, where the nodes are
separated by $s_{i}=\lVert\mathbf{x}_{i+1}-\mathbf{x}_{i}\rVert$. From the
tangent vectors a local (unsigned) curvature of the medial axis is obtained
as:

$\mathcal{K}_{i}=\left\lVert\frac{\mathbf{\hat{t}}_{i+1}-\mathbf{\hat{t}}_{i}}{s_{i+1}-s_{i}}\right\rVert.$
(7)

Considering the local minimum radius of the tube $R_{y}$, the local curvature
at the inner most side of the tube, is then

$\kappa(s)=\frac{\mathcal{K}}{1-\mathcal{K}R_{y}},$ (8)

where $\mathcal{K}=\mathcal{K}(s)$ and $R_{y}=R_{y}(s)$ since in absence of
pressure, the curvature is not homogeneous along the contour length of the
tube (see Section VA and Figure S10). Then, we consider the average value of
$\kappa$ on the 10% central region of the tube, where the curvature is
maximal.

In the pressurized configuration, where the curvature is mostly constant
everywhere along the contour of the tube (Section VB; Figure S10 and Figure
S13), we used an alternative method. From the the geometry of a curved tube
with an angle $\gamma$ between the lids and the $z$ axis
($\cos\gamma=\langle\mathbf{\hat{n}}_{i}\cdot\mathbf{\hat{k}}\rangle_{N}$,
where the average is over the $N$ particles of the lids), it can be shown
that:

$\kappa=\frac{2\gamma}{L-2R\gamma},$ (9)

where $R$ and $L$ are the pressurized radius and lenght, respectively. This
method to calculate an average $\kappa$ is reliable after the occurrence of
the wrinkling transition, where the highly curved regions near the wrinkles
results in a large bias in the estimation of the mean curvature using the
method based on Equation 8 . Therefore, we opted for Equation 9 for the
analysis of the complete response of the tube, including the wrinkled state
(Figures 4 and 5 in the main text). This method, however, slightly
overestimates the curvature at the highest pressures, as a result of
deformation of the lids. Thus, we resort to the more precise values obtained
from Equation 8 for the determination of $\kappa_{w}$ (Figure 6a in the main
text).

###  B Torque, $M$

The external energy $\mathcal{E}_{\textrm{ext}}$ employed to bent the tube
(Equation 5) can be rewritten as:

$\mathcal{E}_{\textrm{ext}}=k_{\tau}\sum_{i\in\begin{subarray}{c}\text{left}\\\
\text{lid}\end{subarray}}(1-\cos\phi_{l})+k_{\tau}\sum_{i\in\begin{subarray}{c}\text{right}\\\
\text{lid}\end{subarray}}(1-\cos\phi_{r}),$ (10)

where $\phi$ are the angles formed between the normal vector of the lids and
the target normal vector. The torque is, per definition, $M=\frac{\partial
E}{\partial\phi}$. Then, considering both lids:

$\begin{split}M&=\frac{1}{2}\left(\frac{\partial\mathcal{E}_{\textrm{ext}}}{\partial\phi_{l}}+\frac{\partial\mathcal{E}_{\textrm{ext}}}{\partial\phi_{r}}\right)\\\
&=\frac{k_{\tau}}{2}\left(\sin\phi_{l}+\sin\phi_{r}\right).\end{split}$ (11)

Alternatively, the torque can be determined from the mechanical energy of the
tube:

$M=\frac{\partial(\Delta\tilde{E}_{int})}{\partial\kappa},$ (12)

where $\Delta\tilde{E}_{int}=\tilde{E}_{int}(\kappa)-\tilde{E}_{int}(0)$ is
the difference of the internal energy _per length_ between the bent and the
unbent pressurized configurations. We calculate the internal energy of the
shell $\mathcal{E}_{\textrm{int}}=E_{s}+E_{b}-pV$, excluding from the terminal
regions of the tube with the objective of minimizing artifacts associated to
the lids. The torque is then calculated by obtaining the numerical derivative
of the energy with respect to the curvature (Equation 12), with $\kappa$
determined as described in the previous section.

Both methods provide equivalent results (Figure S1). Nevertheless, calculating
the torque from Equation 12 is computationally expensive and prone to noise,
specially after the wrinkling transition, due to the numerical derivative over
$\kappa$. Consequently, all the results shown in this work used Equation 11
for the determination of $M$.

###  C Stress resultants, $N_{z}$

The stresses shown in Figure 7 of the main text are the average of the
stresses determined on the bottom-most part of the tube. The method described
in this section is nonetheless valid to find per-vertex stress resultants of
triangulated surfaces at any location.

First, we determine the normal vector $\mathbf{\hat{n}}$ of the given vertex
of interest. The normal vector is obtained from the triangles in the adjacent
neighbour to the main vertex:

$\mathbf{n}_{i}=\sum_{t=1}^{Q}(\mathbf{r}_{1,t}-\mathbf{r}_{i})\times(\mathbf{r}_{2,t}-\mathbf{r}_{i}),$
(13)

where $Q$ is the coordination number, $\mathbf{r}_{i}$ are the coordinates of
the main particle and $\mathbf{r}_{1,t}$ and $\mathbf{r}_{2,t}$ are the
coordinates of the other two particles defining the triangle $t$. More precise
estimates of the local normal vector can be obtained by weighting the normal
vector of each triangle by the corresponding local Voronoi area

$S_{t}=\frac{1}{8}(\lVert{\mathbf{r}_{1}-\mathbf{r}_{i}}\rVert^{2}\cot\angle\mathbf{r}_{2}+\lVert{\mathbf{r}_{2}-\mathbf{r}_{i}}\rVert^{2}\cot\angle\mathbf{r}_{1}),$
(14)

where $\angle$ is the angle at the indicated vertex. For additional details
see Meyer  _et al._ (2003).

Having defined the local normal vector, we project the coordinates of the $Q$
vertices of the first-ring neighbourhood of the particle $\mathbf{r}_{i}$ on
the local plane
$\mathbf{r}_{j}^{\text{proj}}=\mathbf{r}_{j,i}-(\mathbf{r}_{j,i}\cdot\mathbf{\hat{n}}_{i})\mathbf{\hat{n}}_{i}$.
From now on, and with the objective of keeping the notation brief, the
reference to the coordinates $\mathbf{r}$ states for the coordinates after the
plane projection.

We express the coordinates of all the particles in terms of a local and
arbitrary basis. The first basis vector is given by an edge connecting one of
the vertices with the main vertex:
$\mathbf{e}_{1}=\mathbf{r}_{i}-\mathbf{r}_{0}$. A second basis vector can be
obtained as $\mathbf{e}_{2}=\mathbf{e}_{1}\times\mathbf{n}$. Combining these
two basis vectors with the normal vector allows to obtain the local coordinate
basis, defined as
$\mathbf{B}=[\mathbf{\hat{e}}_{1},\mathbf{\hat{e}}_{2},\mathbf{\hat{n}}]$.
Since the coordinate system here defined is orthogonal, we can express the
coordinates of a given point on the surface in terms of this basis as
$\mathbf{r}^{\prime}=\mathbf{B}^{\text{T}}\mathbf{r}$.

Finally, the full stress tensor of particle $i$ is calculated as (Tsai, 2008):

$N_{i}=\frac{1}{S}\sum_{\left<i,j\right>}\mathbf{F}_{i,j}\otimes\mathbf{r}_{i,j}^{\prime},$
(15)

where $\otimes$ indicates the dyadic product and $\mathbf{F}$ is the resulting
force on the edge connecting vertices $i,j$. The total area of the patch is
given by the sum over the Voronoi areas of the $n_{t}$ triangles around
particle $i$: $S=\sum_{t=1}^{n_{t}}S_{t}$. The range of indexes of the tensor
are $\alpha=(1,2,3)$ and $\beta=(1,2,3)$. Plane projection and the change of
basis results in $N^{\alpha\beta}=0$ for every component of tensor having
$\alpha=3$ or $\beta=3$.

The stress tensor obtained is defined on the arbitrary basis given by the
connectivity. Thus, we need to change to a target basis of choice
$\mathbf{B}=[\mathbf{m}_{1},\mathbf{m}_{2},\mathbf{m}_{3}]$ following the
tensor transformation law,
$\mathbf{N}^{\prime}=\mathbf{Q}\mathbf{N}\mathbf{Q}^{\textrm{T}}$, where the
matrix $\mathbf{Q}$ is given by

$\mathbf{Q}=\begin{bmatrix}\mathbf{m_{1}}\cdot\mathbf{\hat{e}_{1}}&\mathbf{m_{1}}\cdot\mathbf{\hat{e}_{2}}&\mathbf{m_{1}}\cdot\mathbf{\hat{n}}\\\
\mathbf{m_{2}}\cdot\mathbf{\hat{e}_{1}}&\mathbf{m_{2}}\cdot\mathbf{\hat{e}_{2}}&\mathbf{m_{2}}\cdot\mathbf{\hat{n}}\\\
\mathbf{m_{3}}\cdot\mathbf{\hat{e}_{1}}&\mathbf{m_{3}}\cdot\mathbf{\hat{e}_{2}}&\mathbf{m_{3}}\cdot\mathbf{\hat{n}}\\\
\end{bmatrix}.$ (16)

In our particular case, the target basis are those giving the orientation
along the axial and circumferential directions of the tube and
$\mathbf{m}_{3}\equiv\mathbf{\hat{n}}$. The basis vector in the axial
direction is determined considering the tangent vectors of the vertices
defining the bottom-most part of the tube:
$\mathbf{m}_{1}=\mathbf{r}_{i+1}-\mathbf{r}_{i}$, excluding the component
along the $x$ axis, which arises as a consequence of the triangular geometry
of the mesh. The basis vector in the circumferential direction can be obtained
as $\mathbf{m}_{2}=\mathbf{m}_{1}\times\mathbf{n}_{i}$.

In Figure 4 we show the stress tensor in the axial direction as the average
over vertices on the bottom-most part of the tube and excluding regions near
the terminal sides. Note that the tensor obtained is of mixed components
$N_{z}^{z}$. In other words, we use lower indexes $N_{zz}$ for convenience of
notation, without implying a covariant form.

###  D Extent of wrinkles, $\theta_{w}$

We use the curvature along the circumferential direction
$\kappa_{\theta\theta}$ as a proxy for the extent of the wrinkles. The
curvature $\kappa_{\theta\theta}$ is calculated from Gaussian and mean
curvatures, as described in Meyer  _et al._ (2003).

The bending of the tube results in a deformation of the cross-section to an
elliptic shape. Though the eccentricity is very small (Figure S13), as a
consequence of the elliptic cross-section of tube, the curvature along the
circumferential direction has two maxima: a peak at $\sim\pi/2$ and other
caused by the wrinkled region prior to decay to non-wrinkled regions (Figure
S4a and S4b). However, particularly at the largest pressures where the
eccentricity of the tube is larger, the position of the first peak associated
to wrinkled area is hard to resolve. To make the identification of the extent
of the wrinkles easier, we subtract the response immediately before the onset
of wrinkling. In mathematical terms,
$\kappa_{\theta\theta}^{off}(\theta,\kappa)=\kappa_{\theta\theta}(\theta,\kappa)-\kappa_{\theta\theta}(\theta,\kappa_{w})$.
The angular extent of the wrinkles is given by:
$\theta_{w}(\kappa)=\operatorname*{arg\,max}_{\theta}\\{\kappa_{\theta\theta}^{off}(\theta,\kappa)\\}$.
Examples for $\kappa_{\theta\theta}^{off}$ at different post-wrinkling
curvatures are shown in Figure S4a-c.

###  E Corrections for the non-linearity of the mesh

The mapping between the mesh model and linear elasticity only applies for low
strains. At the pressures employed in this work the strains are moderately
small ($\approx$ 0.17 at the highest pressure applied), yet significant to
manifest a deviation from linearity (Figure S3). We correct the Young’s
modulus $E$ and the Poisson’s ratio $\nu$ to account for the aforementioned
non-linearity of the mesh model.

The starting point are the expected strains in a pressurized hollow tube:

$\varepsilon_{\theta\theta}=\frac{pR}{Et}\left(1-\frac{\nu}{2}\right)$ (17)
$\varepsilon_{zz}=\frac{pR}{Et}\left(\frac{1}{2}-\nu\right),$ (18)

where $\varepsilon_{\theta\theta}$ and $\varepsilon_{zz}$ are the strains in
the circumferential and axial direction, respectively. These equations arise
after considering that the stresses in a pressurized cylinder are
$2N_{zz}=N_{\theta\theta}=pR$. Taking the derivatives with respect to pressure
of Equations 17 and 18, we can solve for an effective Young’s modulus
$E_{\textrm{eff}}$ and Poisson’s ratio $\nu_{\textrm{eff}}$,

$E_{\textrm{eff}}=\frac{3R_{0}\partial_{p}\varepsilon_{\theta\theta}}{2t(2\partial_{p}\varepsilon_{\theta\theta}-\partial_{p}\varepsilon_{zz})},$
(19)
$\nu_{\textrm{eff}}=\frac{2\partial_{p}\varepsilon_{zz}-\partial_{p}\varepsilon_{\theta\theta}}{\partial_{p}\varepsilon_{zz}-2\partial_{p}\varepsilon_{\theta\theta}}.$
(20)

The average strains on the tube are given by
$\varepsilon_{\theta\theta}=(R-R_{0})/R_{0}$ and
$\varepsilon_{zz}=(L-L_{0})/L_{0}$. To determine the $E_{\textrm{eff}}$ and
$\nu_{\textrm{eff}}$ from the simulation, we find average radius $R$ and
length $L$ in a region away from the terminal sides to prevent artifacts
associated to the lids (Figure S3). The values
$\partial_{p}\varepsilon_{\theta\theta}$ and $\partial_{p}\varepsilon_{zz}$,
are obtained from the derivatives of a 5th order polynomial fit to the strains
at different pressures in absence of torque. The resulting
$\partial_{p}\varepsilon_{\theta\theta}$ and $\partial_{p}\varepsilon_{zz}$
for each pressure is inserted in Equation 19 and 20 to determine the effective
Young’s modulus and Poisson’s ratio.

##  IV Mechanics and Geometry of unpressurized tubes

The buckling of a tube in the absence of pressure has been widely studied
during the last century (Calladine, 1983; Brazier and Southwell, 1927).
Therefore, we use the unpressurized condition as a proof-of-principle of the
ability of the computational model to produce well-established results. First,
we introduce the non-dimensional curvature and torque in the so-called _shell_
normalization, $\bar{\kappa}:=\frac{\kappa R}{\alpha}$ and torque
$\bar{M}:=\frac{M}{\alpha EtR^{2}}$, where
$\alpha:=\frac{t}{R}\sqrt{\frac{1}{12(1-\nu^{2})}}$. For thin shells
$\alpha\ll 1$ and in our case, and unless otherwise stated, $\alpha=$ 0.0024
in the unpressurized configuration. In light of the indicated normalization,
the relation between torque and curvature is, up to second order approximation
(Qiu  _et al._ , 2022; Karamanos, 2002):

$\bar{M}=\pi\left(\bar{\kappa}-\frac{\bar{\kappa}^{3}}{8}\right).$ (21)

We find that, for sufficiently long tubes ($L\gg R$), the torque-curvature
response is very well described by the numerical model (Figure S5a). An
extensive analysis of length dependencies is shown in Sec. VA. We find
quantitative agreement between the model and the analytical theory for the
full $\bar{\kappa}$ range. Moreover, our simulations show that the curvature
leading to the buckling instability is in accord with the theoretical
estimation, which predicts $\bar{\kappa}_{w}=1.320$ for the wrinkling
instability (Qiu  _et al._ , 2022). We observe that the simulation shows a
slightly earlier instability, occurring at a $\sim$3% lower $\bar{\kappa}$
than the predicted $\bar{\kappa}_{w}$. Adjusting the simulation parameters to
produce smaller changes in $\kappa$ per iteration (e.g. reducing $k_{\tau}$)
results in very similar $\kappa_{w}$ values and in every case occurring
slightly earlier than the theoretical prediction.

Furthermore, we analyze the geometry of the tube during bending. It is
predicted that the cross-section of the tube acquires an elliptical shape
(Calladine, 1983; Qiu  _et al._ , 2022). Simulations confirm the formation of
an elliptical cross-section during bending and the eccentricity at the
instability is in very good agreement with the theoretical predictions for the
semi-major and semi-minor radii (Figure S5b, (Qiu  _et al._ , 2022)).

##  V Finite length dependency

The analytical approximations to understand the buckling phenomena rely on the
assumption of an infinitely long tube with isotropic curvature and radius. In
the real world, as well as in the simulations, an infinite length is not
feasible. In this section we study the finite-size implications on the
mechanical and geometrical properties of tubes subjected to bending by using
tubes of different lengths.

###  A Non-pressurized tubes

First, we investigate the mechanical response for different lengths in an
unpressurized situation. From the torque-curvature curves we find that shorter
tubes show a linear torque-curvature dependency (Figure S9a), resembling the
mechanics of a solid tube (Landau and Lifshitz, 2016). In addition, the
curvature at buckling (and thus the torque) is higher in short tubes than in
the longer ones. As the tube length increases, the $\kappa_{w}$ and $M_{w}$
required for buckling progressively decrease, ultimately converging to the
values expected from analytical theory (Figure S9b and Figure S9c).

Next, we observe that shorter tubes do not undergo elliptical deformations of
the cross-section during the bending process (Figure S10a). Similarly as we
observe with the mechanical properties, the semimajor and semiminor radii
(measured at a central region and accounting for a length of a 10% of the
total length) steadily converge to the expected radii predicted by analytical
theory as the length of the tube increases (Figure S10b). Notably, we find
that the radius of the tube is not constant along the contour length. We
observe that, for long tubes, the resulting eccentricity is maximum at the
central region of the tube, while the terminal sides maintain a circular
cross-section (Figure S10c). The curvature of the tube displays a similar
profile: short tubes exhibit a contour-independent curvature and longer tubes
are maximally curved on the middle region and minimal at the lids (Figure
S10d). These results are consistent with the finite-length effects described
when $\Omega:=\sqrt{L^{2}h/R_{0}^{3}}\ll 1$, where
$h:=\frac{t}{\sqrt{1-\nu^{2}}}$ (Calladine, 1983).

In light of these findings, we decided to use tubes with $L/R_{0}=48$
(corresponding to $\Omega=4.42$) for the analysis of the mechanical and
geometrical response to bending in the unpressurized configuration (Figure
S5).

###  B Pressurized tubes

In contrast to the unpressurized case, the mechanical response of the
pressurized tube is length-independent for the range of pressures studied. We
do not find significant differences or length dependencies in the wrinkling
and torque curvatures, even at the lowest pressures (Figure S11). We find a
negligible effect of the length on the geometrical properties of the tube:
even at the lowest pressure, the dependency of the radii and the curvature on
the position along the tube is suppressed by the internal pressure (Figure
S12). Therefore, we resort to tubes with $L/R_{0}=24$ for the analysis of the
wrinkles shown in the main text. This length allowed us to appropriately
measure the amplitude and wavelength of the wrinkles, keeping the number of
vertices within a computationally tractable limit, $N=113,117$ vertices. For
the determination of the $A^{2}/\lambda^{2}$ relation (Fig. 4d of the main
text) we increased the refinement to $N=230,534$ vertices. This was a
requirement to precisely quantify $A$ near the onset of wrinkling, where the
amplitude is minimal.

Remarkably, we observe that for the pressures employed in this work, the
cross-section remains mostly circular during bending (Figure S13a). Pressure
prevents the deformation of the cross-section of the tube and then higher
pressures are less eccentric for a given $\kappa$. However, since the
curvature necessary to buckle the tube grows with $p$, the eccentricity at the
wrinkling curvature also grows with pressure (Figure S13b). We find that even
at very small pressures (radial strains post pressurization below 1%), the
radii of tube and the curvature are independent on the position along the
contour length (Figure S13b and Figure S13c).

##  VI Bending of pressurized cylinder: a tension field theory model

Here we consider the pure bending of a long-capped pressurized cylindrical
membrane accounting for wrinkling where axial compression occurs.

Geometry of bent tubes: notation. The long axis is defined along the length of
the tube.

The long-capped cylindrical membrane is first pressurized to $p$ which is then
held constant and the cylinder is subject to pure bending about the $x$-axis.
In the modeling analysis conducted here, the strains are assumed to be small
and the membrane is modeled by plane stress linear elasticity. In the region
of the cylindrical membrane tube where the axial stress would be negative for
the wrinkling to occur, the axial stress is taken to be zero as an
approximation for the convenience of calculation. This is the source of the
non-linearity and the buckling-like behavior. Prior to bending, the resultant
stress in the cylindrical membrane tube subject to pressure $p$ are

$N_{z}=\frac{pR}{2}\;,N_{\theta}=pR\;,$ (22)

where $R$ is the radius of the tube. During pure bending at constant $p$, the
axial strain satisfies $\varepsilon_{z}=\varepsilon_{0}+\mathcal{K}z$, where
we will regard the axial curvature of the cylinder center-line as imposed and
its axial strain $\varepsilon_{0}$ as an unknown. In the wrinkled region,
$\varepsilon_{z}$ is the axial strain averaged over the wrinkles and not the
axial strain in the membrane itself, however, $\varepsilon_{z}$ in the
wrinkled region does not play a role in the analysis.

In the unwrinkled region of the tube the linear relations relating the
resultant membrane forces per length to the strains are

$N_{z}=\frac{Et}{1-\nu^{2}}(\varepsilon_{z}+\nu\varepsilon_{\theta})\;,N_{\theta}=\frac{Et}{1-\nu^{2}}(\varepsilon_{\theta}+\nu\varepsilon_{z})\;.$
(23)

Here we limit our attention to isotropic linear elasticity. In this analysis
we assume (as an approximation) the pressure is large enough such that each
cross-section remains circular and, thus, circumferential equilibrium requires
$N_{\theta}=pR$ around the entire circumference whether wrinkling occurs or
not.

The energy per unit length of the tube is the sum of the strain energy and the
potential energy of the pressure:

$\Phi=\int_{0}^{\theta_{0}}N_{z}\varepsilon_{z}Rd\theta+\int_{0}^{\pi}N_{\theta}\varepsilon_{\theta}Rd\theta-p\Delta
V\;,$ (24)

with, for $|\theta|<\theta_{0}$,

$\displaystyle N_{z}=Et(\varepsilon_{0}+\mathcal{K}R\cos\theta)+\nu
pR\;,N_{\theta}=pR\;,$ (25)
$\displaystyle\varepsilon_{z}=\varepsilon_{0}+\mathcal{K}R\cos{\theta}\;,\varepsilon_{\theta}=(1-\nu^{2})pR/Et-\nu(\varepsilon_{0}+\mathcal{K}R\cos{\theta})\;.$
(26)

And, for $|\theta|>\theta_{0}$,

$N_{z}=0\;,N_{\theta}=pR\;,\mbox{and}\;\varepsilon_{\theta}=\frac{N_{\theta}}{Et}\;.$
(27)

And the expression for $\Delta V$ is

$\Delta
V=\int_{0}^{\pi}2(1+\varepsilon_{z})R^{2}\sin{\theta}^{2}d\theta+\int_{0}^{\pi}2R^{2}\varepsilon_{\theta}d\theta\;.$
(28)

Given that the condition for boundary $\theta_{0}$ between wrinkled and
unwrinkled regions is $N_{z}=0$, we have the equation for $\theta_{0}$:

$\varepsilon_{0}+\mathcal{K}R\cos\theta_{0}=-\nu\frac{pR}{Et}\;.$ (29)

We then show two equivalent ways to determine the strain $\varepsilon_{0}$ and
then the state of the system: direct approach using the force-balance equation
and energy minimization with regard to $\varepsilon_{0}$.

###  A Direct approach

The overall equilibrium in the axial direction gives:

$\int_{0}^{\theta_{0}}2N_{z}Rd\theta=\pi R^{2}p\;,$ (30)

which results in

$2\mathcal{K}R(\sin{\theta_{0}}-\theta_{0}\cos{\theta_{0}})=\pi\frac{pR}{Et}\;.$
(31)

after the wrinkling has happened ($\theta_{0}<\pi$). The expansion of the
wrinkles with the increase in curvature is shown in Figure S15.

The moment about the x-axis (the cross-section center) gives,

$M=\int_{0}^{\theta_{0}}2N_{z}R^{2}\cos{\theta}d\theta\
=EtR^{2}\left[\theta_{0}-\frac{\sin{(2\theta_{0})}}{2}\right]\mathcal{K}R\;.$
(32)

Again this is the result after the onset of wrinkling. Thus, the moment-
curvature relation for fixed $p$ is determined by varying $\theta_{0}$ and
plotting the above two relations.

Prior to the onset of wrinkling, we won’t find a point where $N_{z}=0$ so
Equation 29 won’t be applicable for this case. Then the overall equilibrium in
the axial direction (Equation 31) should become:

$\int_{0}^{\pi}2N_{z}Rd\theta=\pi R^{2}p\;,$ (33)

which results in

$Et\varepsilon_{0}=\frac{pR}{2}-\nu pR.$ (34)

And thus we have the expression for $N_{z}$ and $M$ before the onset of
wrinkling:

$N_{z}=\frac{pR}{2}+Et\mathcal{K}R\cos{\theta}$ (35)
$M=\pi\mathcal{K}EtR^{3}.$ (36)

Since we assume $N_{z}=0$ for $\theta_{0}$, we can use Equation 35 to
determine the curvature at wrinkling

$\hat{\mathcal{K}}_{w}=\frac{\hat{p}}{2},$ (37)

where $\hat{\mathcal{K}}:=\mathcal{K}R$ and $\hat{p}:=pR/(Et)$.

###  B Energy approach

With $p$ and $\mathcal{K}$ regarded as prescribed, $\varepsilon_{0}$ must
minimize $\phi$ :

$0=\frac{\partial\Phi}{\partial\varepsilon_{0}}=\frac{d}{d\varepsilon_{0}}\left(\int_{0}^{\theta_{0}}N_{z}\varepsilon_{z}Rd\theta\right)+\int_{0}^{\pi}\frac{d(N_{\theta}\varepsilon_{\theta})}{d\varepsilon_{0}}Rd\theta-p\int_{0}^{\pi}\left(2|x|R\sin\theta+2R^{2}\frac{d\varepsilon_{\theta}}{d\varepsilon_{0}}\right)d\theta\;.$
(38)

And we can calculate each term separately,

$\displaystyle\frac{d}{d\varepsilon_{0}}\left(\int_{0}^{\theta_{0}}N_{z}\varepsilon_{z}Rd\theta\right)$
$\displaystyle=\int_{0}^{\theta_{0}}(2N_{z}-\nu N_{\theta})Rd\theta+0\;,$ (39)
$\displaystyle\int_{0}^{\pi}\frac{d(N_{\theta}\varepsilon_{\theta})}{d\varepsilon_{0}}Rd\theta$
$\displaystyle=\int_{0}^{\theta_{0}}-\nu
N_{\theta}Rd\theta+N_{\theta}\varepsilon_{\theta}2R\frac{d\varepsilon_{\theta}}{d\varepsilon_{0}}\Bigr{|}_{\theta=\theta_{0}}\;,$
(40) $\displaystyle-p\int_{0}^{\pi}2|x|R\sin\theta d\theta$
$\displaystyle=-p\pi R^{2}\;,$ (41)
$\displaystyle-p\int_{0}^{\pi}2R^{2}\frac{d\varepsilon_{\theta}}{d\varepsilon_{0}}d\theta$
$\displaystyle=-p\int_{0}^{\theta_{0}}2R^{2}(-\nu)d\theta-2R^{2}p\varepsilon_{\theta}\frac{d\theta_{0}}{d\varepsilon_{0}}\Bigr{|}_{\theta=\theta_{0}}\;.$
(42)

Adding all terms together, we have

$\int_{0}^{\theta_{0}}2N_{z}Rd\theta=p\pi R^{2}\;,$ (44)

which is the same as the equation derived from force-balance equation.

The moment $M$ is given by,

$M=\frac{\partial\Phi}{\partial\mathcal{K}}\;.$ (45)

Using the equation specifying the wrinkling boundary and the expression for
$\varepsilon_{0}$ above, we have:

$\frac{M}{EtR^{2}}=\left[\theta_{0}-\frac{\sin{(2\theta_{0})}}{2}\right]\mathcal{K}R\;.$
(46)

The resulting expression for $M$ agrees exactly with the relation derived by
the direct method. The pre-wrinkling results are the same too and these can be
derived by setting $\theta_{0}=\pi$ in the system energy $\Phi$ .

##  VII Bending of pressurized cylinder: non-linear stress-strain relation

Here, we summarize the main findings of the model presented in (Qiu  _et al._
, 2022). Compared to the TFT model described in the main text, this model took
cross-section deformations into account and used the exact expression for the
critical compressive stress in Equation 7 in the main text.

In particular, the theoretical expression for the shape of the cross-section,
namely the semi-major and semi-minor radii, $R_{x}$ and $R_{y}$, is given by:

$R_{x}=\frac{R}{2}\int_{-\pi/2}^{\pi/2}[\cos\mu+a_{2}\sin(2\mu)]d\mu$ (47)
$R_{y}=\frac{R}{2}\int_{0}^{\pi}[\sin\mu+a_{2}\sin(2\mu)]d\mu$ (48)

The factor $a_{2}$ is defined in Equation 7 in Qiu  _et al._ (2022):

$a_{2}=-\frac{\hat{\kappa}^{2}}{8(\alpha^{2}+\hat{p}/3)}\;,$ (49)

where $\alpha:=\frac{t}{R}\sqrt{\frac{1}{12(1-\nu^{2})}}$. Taylor expansion of
the trigonometric functions in Eqs. 47 and 48 and considering up to first
order allows to write:

$R_{x}\approx R\left(1-\frac{2}{3}a_{2}\right),$ (50) $R_{y}\approx
R\left(1+\frac{2}{3}a_{2}\right),$ (51)

The eccentricity is then given by:

$\frac{R_{x}}{R_{y}}\approx\frac{3-2a_{2}}{3+2a_{2}}\;,$ (52)

Equations 50 and 51 are derived from Equations S26 and S27 in Qiu  _et al._
(2022), which denotes the shape of the cross-section from the coefficients
$a_{n}$ (for $n=1,2,...,\infty$). For the zero pressure case, it was proven
that $a_{2}$ is the leading term. For the pressurized case, we proved that
only even terms are allowed and we assume $a_{2}$ is the leading term while we
approximate all the other $a_{n}$ terms to be zero. All together leads to Eqs.
50-51. Interestingly, Eq. 52 predicts a non-monotonic dependency of the cross-
section shape with pressure, reaching almost no-deformation in the range of
$\hat{p}$ from 10-3 to 10-1. Our numerical results are in good agreement with
the predictions given by Eq. 52 (Figure S13b).

Accounting for the deformation of the cross-section, we obtain the following
expression for the torque-curvature relation:

$\hat{M}=\pi\left[\hat{\kappa}+\frac{\hat{\kappa}^{3}}{\alpha^{2}+\hat{p}/3}\right].$
(53)

In Figure S8a we show the simulation results together with the predictions of
53 for different pressures.

We found that the dimensionless critical wrinkling curvature
$\hat{\kappa}_{w}$ follows the equation:

$8\left(\alpha^{2}+\frac{\hat{p}}{3}\right)\left(\frac{\hat{p}}{2}+2\alpha-\hat{\kappa}_{w}\right)+\left(\frac{2\hat{\kappa}_{w}}{3}-4\alpha\right)\hat{\kappa}_{w}^{2}=0\,.$
(54)

For pressurized thin shells such that $\alpha\ll\hat{p}$, Equation 54 can be
written, to first order as, $\hat{\kappa}_{w}\approx\hat{p}/2$. We compare the
results from computational model with the predicted $\kappa_{w}$ given by 54
and by the tension field theory model described in the previous section
(Equations 54 and 27) in Figure S8b and Figure S8c. We find good agreement
between the simulation results and the theoretical predictions for both
$\hat{\kappa}_{w}$ and $M_{w}$. We remark the quantitative agreement at the
lowest pressures between the model and the theory.

##  VIII Wrinkles geometry from axial in-extensibility conditions

Here we derive the expression for the relation between the amplitude $A$ and
wavelength $\lambda$ of the wrinkles in the sinusoidal regime by considering a
geometric constrain due to axial in-extensibility. We will follow an approach
similar to that described by Mahadevan et al. (Mahadevan  _et al._ , 2004). We
assume that the radial displacement of the wrinkles is sinusoidal in the
longitudinal and in the circumferential directions. Thus, in cylindrical polar
coordinates $(r,\theta,z)$ (as defined in Fig. 1 in the main text), we can
write down the radial displacement of the wrinkles as
$u(r,\theta,z)=A\sin(2\pi z/\lambda)\cos((\pi-\theta)/\theta_{w}\pi/2)$. The
term $\theta_{w}=\pi-\theta_{0}$ is the extent of the wrinkles.

The axial strain of wrinkles forming on a bent tube with imposed medial axis
curvature $\mathcal{K}$, compared to that of $\mathcal{K}_{w}$, has two
contributions. One is the strain caused by the sinusoidal wrinkles and it
bears the form
$\frac{\sqrt{dz^{2}+du^{2}}}{dz}-1\approx\frac{1}{2}\left(\frac{\partial
u}{\partial z}\right)^{2}$. The other comes from the imposed curvature
$\mathcal{K}$ and it bears the form
$R(\mathcal{K}-\mathcal{K}_{w})\cos(\pi-\theta)$ under the assumption of
circular cross-section with radius $R$. Then, at certain position
$(r,\theta,z)$, the total axial strain $\varepsilon$ is approximately:

$\varepsilon=\frac{1}{2}\left(\frac{\partial u}{\partial
z}\right)^{2}-R(\mathcal{K}-\mathcal{K}_{w})\cos(\pi-\theta)\;.$ (55)

We express the axial in-extensibility constraint in the following form:

$0=2\int_{\pi-\theta_{w}}^{\pi}\int_{0}^{\lambda}\left[\frac{1}{2}\left(\frac{\partial
u}{\partial
z}\right)^{2}-(\hat{\mathcal{K}}-\hat{\mathcal{K}}_{w})\cos(\pi-\theta)\right]dzd\theta\;,$
(56)

which leads to

$\frac{A^{2}}{\lambda^{2}}=\frac{2\sin(\theta_{w})}{\pi^{2}\theta_{w}}(\hat{\mathcal{K}}-\hat{\mathcal{K}}_{w})\;.$
(57)

Since $\theta_{w}$ is independent of pressure (Fig. 4c of the main text) the
$A/\lambda$ relation is preserved for all pressures. The functional form of
equation 57 is different to that by Mahadevan et al. (Mahadevan  _et al._ ,
2004). However, we recover the same expression for $\theta_{w}=\pi/2$, which
leads to
$\frac{A^{2}}{\lambda^{2}}=\frac{4}{\pi^{3}}(\hat{\mathcal{K}}-\hat{\mathcal{K}}_{w})$.

To better quantify $A$ near the onset of wrinkling, where the amplitude is
negligible compared with the radius of the tube, we increased the resolution
the mesh to $N=230,534$ vertices to produce the data in Fig. 4d of the main
text. Still, we observe a progressive deviation between the numerical results
and the theory as the pressure increases (Fig. S16). We attribute this
deviation to the implicit loss of refinement with the expansion of the area
tube $p$, reducing the accuracy on the determination of both $A$ and $\kappa$
beyond the wrinkling transition.

##  IX Supplementary Figures

Figure S1: Estimation of torque in the computational model. Comparison of
methodologies. a) $M-\kappa$ curves in the unpressurized configuration. b)
$M-\kappa$ curves in the pressurized configuration for $\hat{p}=0.126$. Large
noise is observed after wrinkling when the torque is calculated from the
numerical derivative of the energy-curvature
($M=\partial_{\kappa}\Delta\tilde{E}_{int}$).

Figure S2: Longitudinally aligned meshes: summary of mechanical and
geometrical properties a) Comparison of the different mesh geometries. b)
Torque-curvature ($M-\kappa$) curves in the circumferentially aligned mesh vs.
the longitudinally aligned mesh. c) Wavelength of the wrinkles $\lambda$ in
circumferentially and longitudinally aligned meshes. The rescaled pressures in
panels b) and c) are, from left to right, $\hat{p}=$ 0.026, 0.088 and 0.164.
The correction in the effective Young’s modulus $E_{eff}$ and Poisson’s ratio
$\nu_{eff}$ for the rescaling were obtained considering pressurization data
from the circumferentially aligned configuration. The subtle differences in
both $M-\kappa$ and $\lambda$ as pressure increases can consequently be
attributed to not having considered specific $E_{eff}$ and $\nu_{eff}$ for the
longitudinally aligned mesh.

Figure S3: Non-linearity of the mesh. a) Radius as a function of pressure. b)
Length as a function of pressure. The red line is the expected response from
linear elasticity.

Figure S4: Determination of the extent of the wrinkles, $\theta_{w}$. a)
Curvature along the circumferential direction $\kappa_{\theta}$ immediately
before wrinkling ($\kappa/\kappa_{w}=0.99$). b) $\kappa_{\theta}$ at the onset
of wrinkling ($\kappa/\kappa_{w}=1.00$). The extent of the wrinkles is
obtained from the peak in
$\kappa_{\theta}^{off}=\kappa_{\theta\theta}(\theta,\kappa)-\kappa_{\theta\theta}(\theta,\kappa_{w})$.
Examples for $\kappa/\kappa_{w}=$ 1.02, 1.13 and 1.4 are shown in panels c),
d) and e) at $\hat{p}=0.026$. The position of the peak, indicated with a
vertical dashed line, serves as a signature for the extent of the wrinkles,
$\theta_{w}$.

Figure S5: Mechanical response to bending in unpressurized tubes. a)
Simulation results at progressively increasing induced curvatures. b) Torque-
curvature response. The theoretical predictions considering the full
expression and second order approximation (Equation 21) are shown in dotted
and continuous black lines, respectively. The expected buckling curvature is
indicated by the vertical line ($\bar{\kappa}_{w}=1.320$). c) Evolution of the
semi-major and semi-minor axes of the cross-section of tube during bending.
The dots are the simulation results. Continuous lines are the theoretical
result determined from the numerical integration of Equations 47 and 48. The
horizontal dashed lines indicate the expected radii at the wrinkling curvature
($R_{x}=$ 1.1318; $R_{y}$ = 0.8430).

Figure S6: Wrinkling patterns during bending. a) Resulting configurations
after the wrinkling transition for $\hat{p}\equiv pR/(Et)=0.026$ at
progressively increasing curvatures $\kappa/\kappa_{w}$ from top to bottom,
where $\kappa_{w}$ is the curvature at the onset of wrinkling. The surface is
colored according to Gaussian curvature, $\kappa_{G}$ (Meyer  _et al._ ,
2003). Radius at the bottom of the tube $R_{y}$ along the arc length of the
tube $s$ (normalized to the range -1 to 1). The pictures at
$\kappa/\kappa_{w}=1.69$ and $\kappa/\kappa_{w}=2.09$ are representative
examples of the system in the multi-wavelength state (MWS) in the low and high
wavelength substates, respectively. b) Configurations after the wrinkling
transition for $\hat{p}=0.164$ and $R_{y}$ as a function of $s$. The pictures
at $\kappa/\kappa_{w}=1.13$ and $\kappa/\kappa_{w}=1.20$ are representative
examples of the system in the MWS at the low and high wavelength substates,
respectively.

Figure S7: Influence of pressure on the wrinkle patterns. a) Radii from the
center of mass of the tube to the bottom side of the tube (_left_). Power
spectral density of the $R(s)$, _right_. Data for $\hat{p}=$ 0.026. The
curvatures are, from top to bottom, $\kappa/\kappa_{w}$ = 1.01, 1.13, 1.69,
2.09. b) Radii, _left_ , and power spectral density of the $R(s)$ function,
_right_. Data for $\hat{p}=0.164$. The curvatures are, from top to bottom,
$\kappa/\kappa_{w}$ = 1.00, 1.09, 1.13, 1.20. The power spectral density is
defined as $S_{RR}(f_{s})=|\tilde{R}(f_{s})|^{2}$, with
$\tilde{R}=\sum_{s=0}^{L}R(s)\exp(isf_{s})$ the discrete Fourier transform of
the radius.

Figure S8: Comparison analytical models for $\hat{\kappa}$ and $M$. a) Torque-
curvature response. Model considering non-linear stress-strain relations. b)
Torque-curvature response. Detail at low pressures and comparison between
tension field theory (dashed lines) and the model considering non-linear
stress strains relations (continuous lines). c) Curvature $\kappa_{w}$ as a
function of pressure. Continuous lines are the predictions given by Equations
37 and 54. d) Torque $\hat{M}_{w}$ at the wrinkling transition as a function
of pressure. Continuous lines are the predictions given by Equations 36 and
53. The rescaled medial axis curvature is related with the inner most
curvature by $\hat{\mathcal{K}}=\hat{\kappa}/(1+\hat{\kappa})$.

Figure S9: Impact of tube’s length on mechanic parameters (non-pressurized).
a) Torque-curvature curves for different lengths $L$ of the tube. b) Curvature
at wrinkling $\bar{\kappa}_{w}$ depending on lengths. c) Torque at the onset
of wrinkling $\bar{M}_{w}$ for different lengths. The continuous lines are the
theoretical expectations (see unpressurized section in the main text).

Figure S10: Impact of tube’s length on geometric parameters (non-pressurized).
a) Maximum and minimum cross-section radii ($R_{x}$ and $R_{y}$) of the cross
section as a function of the curvature for different lengths $L$ of the tube
(ratio considering the undeformed shape). b) Ratio between $R_{x}$ and $R_{y}$
at the onset of wrinkling for different lengths. Radii obtained in the 10%
central region. c) Dependency of the maximum and minimum radii $R_{x}$ and
$R_{y}$ (dashed lines) as a function of $s$, the arc length distance to the
center of the tube, (scaled to the interval -1 to 1), for different lengths c)
Curvature of the tube $\kappa$ as a function of $s$. Data in panels c) and d)
taken from to immediate step before the wrinkling transition.

Figure S11: Impact of tube’s length on mechanic parameters (pressurized). a)
Curvature at wrinkling $\bar{\kappa}_{w}$ depending on length. b) Torque at
the onset of wrinkling $\bar{M}_{w}$ for different lengths. The dashed lines
are the theoretical expectation (Equations 54 and 53). The insets are a detail
of the lowest pressures regime.

Figure S12: Impact of tube’s length on geometric parameters (pressurized). a)
Dependency of the semimajor ($R_{x}$) and semiminor ($R_{y}$) radii
immediately before the onset of wrinkling $\bar{\kappa}_{w}$ for different
lengths $L$ depending on the position along the contour of the tube, $s$. b)
Local curvature $\hat{\kappa}$ as a function of the contour length.

Figure S13: Geometry of pressurized tubes during bending. a) and b) Evolution
of the semi-major and semi-minor axes of the cross-section of tube during
bending. The dots are simulation results. Continuous lines are the theoretical
result determined from the numerical integration of Eqs. 47 and 48. The
horizontal dashed lines indicate the expected radii at the wrinkling
curvature. Results for $\hat{p}=0.013$ (a) and $\hat{p}=0.164$ (b). c) Cross-
section dimensions for the pressures shown in the main document. The deviation
is $<$ 2% for all the pressures (compare with Figure S10a). Low pressures
suppress the dependency on the radius and the curvature on the arc length. d)
Eccentricity at the onset of wrinkling as a function of pressure for the
pressures studied on the main document. The grey line is corresponds to the
theoretical predictions obtained from Equation 52. The inset shows the
theoretical prediction over a wide range of pressures, highlighting the non-
monotonic response of the deformations with pressure. The pressures studied in
the main documents are shown with coloured dots. Very low pressures are shown
with black dots. e) Deformations of the cross-section as a function of the
normalized arc length of the tube, $s$ at the onset of wrinkling at low
pressures (corresponding to the black dots in panel d). The continuous lines
is the semi-major radius $R_{x}$ and the dotted line is the semi-minor radius
$R_{y}$. f) Curvature as a function of $s$ at the onset of wrinkling at low
pressures. The pressure in the _shell_ normalization is defined as
$\bar{p}:=pR^{2}$. The relation with the _balloon_ normalization is
$\hat{p}=\alpha^{2}\bar{p}$. Here $\hat{p}$ are on the order $10^{-5}$ \-
$10^{-4}$. Data for $L/R_{0}=$ 48.

Figure S14: Dependency of the wavelength $\lambda$ of the sinusoidal wrinkles.
Values determined at the onset of wrinkling. The dashed line is the predicted
$\lambda_{\mathrm{sin}}$ of an unpressurized bent tube at the onset of
wrinkling (Seide and Weingarten, 1961):
$\lambda=2\pi\left(\frac{B}{Y}R_{\theta}^{2}\right)^{1/4}$. This simple
formula predicts an implicit increase of $\lambda$ with pressure as a result
of the expansion of $R_{\theta}$ . In contrast, we observe that pressurization
leads to a decrease in $\lambda$.

Figure S15: Expansion of wrinkles $\theta_{w}$ from tension field theory as a
function of curvature $\kappa_{w}$ for different pressures (Equation 31). The
area expanded by the wrinkles is defined as $\theta_{w}:=\pi-\theta_{0}$.

Figure S16: $(A/\lambda)^{2}$ in the sinusoidal regime as a function of
curvature. The continuous lines are the predictions from Eq. 57.

## References

  * Meyer  _et al._ (2003) M. Meyer, M. Desbrun, P. Schröder, and A. H. Barr, in  _Visualization and Mathematics III_ , edited by H.-C. Hege and K. Polthier (Springer Berlin Heidelberg, Berlin, Heidelberg, 2003) pp. 35–57.
  * Pastrana  _et al._ (2023) C. Pastrana, L. Qiu, S. Armon, U. Gerland, and A. Amir, Soft Matter 19, 2224 (2023).
  * Lien and Kajiya (1984) S.-l. Lien and J. T. Kajiya, IEEE Computer Graphics and Applications 4, 35 (1984).
  * Seung and Nelson (1988) H. S. Seung and D. R. Nelson, Phys. Rev. A 38, 1005 (1988).
  * Gompper and Kroll (2004) G. Gompper and D. M. Kroll, “Triangulated-surface models of fluctuating membranes,” in  _Statistical Mechanics of Membranes and Surfaces_ (World Scientific, 2004) Chap. 12, pp. 359–426.
  * Landau and Lifshitz (2016) L. D. Landau and E. M. Lifshitz,  _Theory of Elasticity_ (Elsevier, 2016).
  * Nocedal and Wright (2006) J. Nocedal and S. J. Wright,  _Numerical Optimization_ (Springer, 2006).
  * git (2023) “brazier_buckling,” https://github.com/gerland-group/brazier_buckling (2023).
  * Yao  _et al._ (1999) X. Yao, M. Jericho, D. Pink, and T. Beveridge, Journal of Bacteriology 181, 6865 (1999).
  * Deng  _et al._ (2011) Y. Deng, M. Sun, and J. W. Shaevitz, Phys. Rev. Lett. 107, 158101 (2011).
  * Tsai (2008) D. H. Tsai, The Journal of Chemical Physics 70, 1375 (2008).
  * Calladine (1983) C. R. Calladine, “The brazier effect in the buckling of bent tubes,” in  _Theory of Shell Structures_ (Cambridge University Press, 1983) p. 595–625.
  * Brazier and Southwell (1927) L. G. Brazier and R. V. Southwell, Proceedings of the Royal Society of London. Series A 116, 104 (1927).
  * Qiu  _et al._ (2022) L. Qiu, J. W. Hutchinson, and A. Amir, Phys. Rev. Lett. 128, 058101 (2022).
  * Karamanos (2002) S. A. Karamanos, International Journal of Solids and Structures 39, 2059 (2002).
  * Mahadevan  _et al._ (2004) L. Mahadevan, J. Bico, and G. McKinley, Europhys. Lett 65, 323 (2004).
  * Seide and Weingarten (1961) P. Seide and V. I. Weingarten, Journal of Applied Mechanics 28, 112 (1961).