On the buckling of an elastic holey column

We report the results of a numerical and theoretical study of buckling in elastic columns containing a line of holes. Buckling is a common failure mode of elastic columns under compression, found over scales ranging from metres in buildings and aircraft to tens of nanometers in DNA. This failure usually occurs through lateral buckling, described for slender columns by Euler’s theory. When the column is perforated with a regular line of holes, a new buckling mode arises, in which adjacent holes collapse in orthogonal directions. In this paper, we firstly elucidate how this alternate hole buckling mode coexists and interacts with classical Euler buckling modes, using finite-element numerical calculations with bifurcation tracking. We show how the preferred buckling mode is selected by the geometry, and discuss the roles of localized (hole-scale) and global (column-scale) buckling. Secondly, we develop a novel predictive model for the buckling of columns perforated with large holes. This model is derived without arbitrary fitting parameters, and quantitatively predicts the critical strain for buckling. We extend the model to sheets perforated with a regular array of circular holes and use it to provide quantitative predictions of their buckling.


Introduction
Buckling instabilities of elastic structures subjected to deforming forces are found on all scales [1] ranging from large-scale applications such as aircraft to the engineering of DNA [2]. An everyday example can be realized using a plastic coffee stirrer, which deflects laterally if compressed with sufficient force between the forefinger and thumb. As in this example, failure due to compression is typically through a global buckling instability with a wavelength comparable to the length of the structure. Non-uniformities in internal structure introduce localized (short-wave) buckling instabilities that can compete and interact with the global (long-wave) buckling. For example, when under compression the long-wave mode transitions to wrinkling in multilayered composites made of thin interfacial layers embedded in a softer matrix [3]. This short-wavelength buckling behaviour can also be exploited in the design of novel mechanical metamaterials that contain periodic arrays of holes [4]. If an elastic sheet is perforated with a two-dimensional square array of circular holes, the sheet can exhibit pattern switching upon compression that internalizes the buckling: the circular holes deform into ellipses with adjacent holes elongated in orthogonal directions [5,6]. The resulting material properties of the sheet, including negative Poisson's ratio [7], have been applied to design of photonic [8,9] and phononic [10] cellular devices and have even been used in soft robotics [11]. Similar pattern switching also occurs in a column containing a line of equally spaced holes first studied by Pihler-Puzović et al. [12], in which the traditional lateral buckling of a column under compression can be preceded by an instability of the micro-structure between the holes (figure 1).
Here, we focus on in-plane (two-dimensional) buckling of a holey column (figure 1a) as an exemplar system to elucidate the mechanisms underlying the pattern switching in perforated elastic solids. As the column is compressed, it initially reduces in length, but preserves its vertical and horizontal reflection symmetries until it reaches a critical strain at which it buckles. We find from numerical calculations that the first unstable mode of in-plane buckling (the mode occurring at the smallest strain) is always one of two modes, which we term the alternating mode and the Euler mode (figure 1b), although other modes of buckling, including second-and higher-order Euler modes, bifurcate from the unbuckled compression branch at higher strains. In the alternating mode, the column remains straight under loading and the holes contained within the structure arrange in series of ellipses for which each major axes is orthogonal to its neighbour. The Euler mode is a classic sideways buckling, but, as we shall show, it can act in two qualitatively distinct ways. The Euler mode can act as a global mode, where the critical strain and other properties of the mode are dependent on the dimensions of the whole column. It can, alternatively, act as a localized buckling mode in which, although deformation occurs throughout the length of the column, the critical strain is independent of column length and depends only on the geometry of the holes.
Pihler-Puzović et al. [12] considered a fixed hole size and spacing relative to the column width, identified the buckling modes that occur in two-and three-hole columns, and quantified the transition between the alternating and Euler modes as the length of the column increased. They found good agreement between their finite-element numerical calculations and experiments on holey columns made of the hyperelastic material extra hard Sid AD Special (Feguramed GmbH).
In this paper, we use the validated numerical model from Pihler-Puzović et al. [12] to investigate, in detail, the effects of hole size, number and spacing on the pattern formation. We show that the presence of holes in a column introduces localized bucking. This localized buckling occurs not only in the alternating mode but also in the Euler mode, which can be qualitatively altered by the presence of holes. To analyse this localized buckling, we study two-hole columns with periodic boundary conditions, which are unstable only to the localized buckling modes under compression. In addition, we develop a new analytical model for the localized buckling of holey columns when the holes are large and separated by relatively thin ligaments of elastic material. In these columns, the deformation occurs almost entirely through buckling of the ligaments, with thicker regions of the column remaining nearly rigid. This allows the buckling of the column as a whole to be understood through the deformation occurring in the ligaments. The Each buckled state occurs through a pitchfork bifurcation, which breaks a reflection symmetry of the uncompressed column. (c) The alternating mode with an odd number of holes, which does not break a reflection symmetry. In (b),(c), the shading indicates the distribution of strain energy in the deformed state (as defined in [12]), from black (low strain) to yellow (maximum strain for that column). (Online version in colour.) asymptotic predictions from this model have no arbitrary or fitted parameters and compare very favourably with numerical solutions for a surprisingly wide range of hole sizes. Consequently, predictions of the asymptotic behaviour of the system provide insights into the observed pattern formation. In particular, we find that the transition between localized Euler and alternating mode buckling is controlled by the bending stiffness of the two types of ligament in the column: (i) those (initially horizontal) separating the holes and (ii) those (initially vertical) between the hole and the outer edge of the column. We extend the algebraic model to predict the critical stress and strain for deformation of a two-dimensional cellular sheet, which can be viewed as an array of connected holey columns. Good quantitative agreement between these algebraic predictions and numerical finite-element calculations of a cellular sheet provide a simple predictive model of buckling in sheets with square arrays of holes, indicating that the underlying behaviour is also controlled by the ligaments between holes.
The paper is organized as follows: in §2, we use numerical solutions of a plane-strain, finiteelement model to study the global and localized buckling of many different holey columns. We briefly revisit the formulation of the problem and the numerical methods from Pihler-Puzović et al. [12] in §2a. In §2b, we demonstrate how the Euler mode of buckling is qualitatively altered by the presence of holes, and discuss how the column length and hole geometry result in several different regimes of behaviour that are observed within this mode. In §2c, we show how these geometrical parameters influence buckling in the alternating mode, and in §2d, we describe the ranges of column geometries for which each different mode is the first to become unstable. In §2e, we discuss the secondary bifurcations in the finite-length holey column. Our asymptotic model for a localized buckling column is introduced in §3, and is applied to the alternating mode in §3a and to the Euler mode in §3b. Finally, in §4, we extend the asymptotic model to the twodimensional cellular sheet. The conclusions are presented in §5. geometry is, therefore, characterized by three non-dimensional parameters: N, D/W and h/W. As in [12], we investigate the behaviour of the column under compression by solving the equations for an incompressible neo-Hookean 1 hyperelastic material within this geometry, using the C++ finite-element library oomph-lib [13]. In our numerical code, these equations are formulated using the principle of virtual displacements, with fields discretized within a single finite element using quadratic interpolation for positions and linear interpolation of the solid pressure, which is continuous over element boundaries. We treat the problem as two-dimensional with no out-ofplane buckling and, assuming that the compression is quasi-static, solve for equilibrium states.

Parametric study of
We also assume that loading by gravity is insignificant. As the deformations studied here are not large enough to cause contact or self-intersection of the column, we apply zero-stress boundary conditions to the column side walls and the boundaries of the holes. We consider two cases for the boundary conditions at the top and bottom of the column. For columns of finite length, 'clamped' boundary conditions are applied, i.e. the top and bottom of the column are constrained to allow no horizontal deformation, and a uniform vertical deformation is prescribed to emulate compression. We also model columns without these boundaries, by applying periodic boundary conditions that match the deformation and stress at the top and bottom edges of a column with two holes. In this case, the distance between the upper and lower boundaries of the periodic domain is varied to apply compression. These periodic columns are not influenced by the presence of upper and lower boundaries, and so they provide a straightforward way of studying localized buckling modes.
We quantify the compression by recording the engineering strain, ε = y/(Nh), where Nh is the total length of the uncompressed column and y is the amount by which the entire column is compressed vertically. For small strains, and before buckling occurs, the engineering strain varies linearly with engineering stress, defined as the ratio of the total applied compressive force F to column width W, E eff ε = F/W. The proportionality constant E eff is the effective Young's modulus. The critical buckling strain, ε cr , is defined as the value of engineering strain at which buckling occurs; the corresponding stress is termed the critical stress. The buckling is a local bifurcation, and therefore occurs at a change in the sign of the real part of the most unstable eigenvalue. We use numerical continuation to solve for both the equilibrium configuration and the stability eigensystem [14], at successively higher compressive strains. A bifurcation tracking procedure is used to determine the compression at which this bifurcation occurs. The domain is discretized for finite element analysis using a mesh that preserves the symmetries of the undeformed configuration, namely the horizontal and vertical reflection symmetries of the column, and the permutation symmetry among the unit cells ( figure 1a). This choice of mesh facilitates accurate identification of the symmetry-breaking bifurcations. Uniform refinement of the mesh has been used to verify the convergence of the numerical results presented here.
Recall that an ideal undeformed holey column has both horizontal and vertical symmetries. The Euler buckling of a holey column is concurrent with the breaking of the horizontal symmetry (figure 1b, centre), and occurs through a pitchfork bifurcation. If the column has an even number of holes, buckling in the alternating mode breaks the vertical symmetry, through a pitchfork bifurcation (figure 1b, right). If the column has an odd number of holes (figure 1c) the alternating mode buckling does not break the vertical symmetry, and instead occurs via a transcritical bifurcation with two non-conjugate branches that meet in a limit point [12]. We restrict our study to columns with an even number of holes, in which both the Euler mode and the alternating mode are symmetry breaking, each occurring through a supercritical pitchfork bifurcation that breaks a different symmetry. Naturally, the loads at which each of these two bifurcations occurs vary with column length N and the geometrical parameters h/W and D/W, and the bifurcations may swap their order of occurrence. We now study the behaviour of each mode in turn.

(b) Euler buckling of perforated columns
The lateral buckling of long slender columns without holes is described by Euler's theory, which predicts that the critical strain for buckling scales as the inverse square of the column length [15]. The Euler buckling of sufficiently long columns with a line of holes is qualitatively similar to that of solid columns, with the critical strain for buckling scaling as the inverse square of the column length (or number of holes), cr ∼ 1/(Nh) 2 (figure 2a-c). We call this a global buckling regime, because the critical buckling strain is dependent on the column length, and buckling occurs over the length of the whole column, and not at the lengthscale of the unit cells (figure 2d, column i).
Counterintuitively, the critical buckling strain for long columns can be increased beyond the Euler prediction by the addition of holes to the column. This is shown in figure 2c, where the critical buckling strain for columns with holes has been normalized with respect to the critical strain predicted by Euler's theory for a solid column of the same dimensions, and plotted against non-dimensional hole size D/W. Particularly, for large N, the numerical predictions collapse over a range of hole sizes, supporting the scaling cr ∼ 1/(Nh) 2 , and show a trend of increasing critical strain with increasing hole size. The cause of the increase of critical buckling strain with hole diameter is that enlarging a hole not only decreases the critical buckling stress, but also decreases the effective Young's modulus of the column. This decrease in Young's modulus can increase the critical buckling strain of a long column perforated with large holes by a factor of two, relative to the critical buckling strain of a solid column of the same length and width (figure 2c). For small holes (D/W 1), the critical buckling strain approaches the prediction from Euler's theory for a solid column, and, during the initial compression phase prior to bifurcation, the stress is concentrated near the holes.
It is evident from figure 2c that the collapse of critical buckling strains on a 1/(Nh) 2 scaling no longer holds when the hole size D/W is large. Geometrical constraints mean that the diameter of the holes D must be smaller than both the column width W and the distance between hole centres h. As the hole diameter approaches this limit, two distinct behaviours are observed in the Euler mode, depending on whether the distance between hole centres is smaller or larger than the column width.
In the first case, when the distance between hole centres is smaller than the column width (h < W), large holes result in very thin ligaments of elastic material between adjacent holes, in which the strain energy is concentrated. These thin ligaments allow the column to shear with minimal resistance (figure 2d, column ii). Consequently, as the hole diameter approaches the distance between hole centres, the critical buckling strain decreases rapidly (figure 2c). This reduction in critical buckling strain occurs for a wide range of hole sizes in short columns (N 8), but is significant in long columns only as the diameter of holes approaches their spacing (D approaches h).
In the second case, when the distance between hole centres is larger than the column width (h > W), large holes instead result in very thin ligaments of elastic material separating the edges of the column from the holes, with relatively thicker regions of material between adjacent holes. In contrast with the behaviour when h < W, the critical buckling strain of such columns is nearly independent of the number of holes (figure 2b), and decreases with increasing hole size. The buckling of the whole column in this regime results from the individual buckling of the thin ligaments at the edges of the column (figure 2d, column iii). We refer to this as the 'sliding' regime of the Euler mode, as the thicker regions of the column separating the holes do not deform or rotate significantly, but instead slide left or right. An increase in the hole diameter D reduces the thickness of the ligaments at the edges of the column and causes a corresponding decrease in the critical buckling strain (figure 2b). The buckling of the ligaments at the edges of the column in this sliding regime is in contrast with the global regime of buckling ( cr ∼ 1/(Nh) 2 ) observed for long columns, in which the ligaments at the edges of the column are compressed or stretched by the bending of the column, but do not individually buckle.   [12]. (d) Post-buckling states of columns with eight holes after buckling in the Euler mode, with shading indicating local strain energy, increasing from black to yellow, where yellow indicates the maximum local strain energy for each column. Roman numerals identify the bifurcations that led to these states in subfigures (a-c). Column iv is a periodic column with period of two holes, with the same unit cell geometry and strain as column iii; the periodic extension to eight holes is shown with reduced opacity. (Online version in colour.) As the critical buckling strain for the column in the sliding regime is independent of the number of holes in the column, it depends only on the localized geometry of these thin ligaments. In the sliding regime, the deformation of each ligament is isolated from the deformation in adjoining ligaments by the thicker region of material between holes. This means the ligaments at the edges of the column buckle to the left or to the right almost independently of one another.  A consequence of this is that the critical buckling strain in this sliding regime is in quantitative agreement with the buckling strain in a two-hole periodic column with the same unit cell geometry (dashed lines in figure 2b), in which the ligaments buckle alternately left and right (figure 2d, column iv). A related consequence of this isolation of the buckling in adjacent holes is that the critical buckling strains for second-and higher-order Euler modes in the sliding regime are only slightly larger than that of the first Euler mode. This 'clustering' of many successive Euler modes around a single critical strain in the sliding regime explains why the two-hole periodic column is a good predictor of the critical strain of the first Euler mode (figure 2d, column iii), despite the fact that when extended periodically to N holes it is closest in appearance to the Euler mode of order (N − 1) (figure 2d, column iv).
The addition of a line of holes to an elastic column, therefore, significantly alters its Euler buckling behaviour. The main influence of the holes is the new sliding regime of buckling, which is a localized buckling mode that depends on the properties of the holes and not the column length. Importantly, this sliding regime is not a separate mode of instability, but is a new regime of behaviour of the existing Euler buckling mode, induced by the presence of holes. There is consequently a continuous transition between the localized behaviour of the sliding regime and the classical Euler buckling regimes, as illustrated in figure 2b.

(c) Alternating mode buckling
The critical buckling strain for the alternating mode is qualitatively different from the Euler mode, in that it tends to be a constant for long columns (figure 3a), rather than scaling as 1/N 2 . The critical buckling strain for long columns is in good agreement with that of a twohole periodic column buckling in the alternating mode, indicated by dashed lines in figure 3a. As with the sliding regime of the Euler mode, this is a localized buckling mode, characterized by independence of the critical strain on N, and agreement between predictions from finite-length and periodic columns. The post-bucking states, illustrated in figure 3b, indicate that buckling of both the ligaments separating holes and the ligaments on the column edge occurs in this localized mode. The critical buckling strain for the alternating mode decreases with increasing hole size D/W, the opposite trend from that of long columns in the Euler buckling mode.
In the alternating mode, the critical strain in short columns is larger than that in long columns. This increase in strain is caused by boundary effects from the clamped boundary conditions; these

(d) Exchange of stability between Euler and alternating modes
The number of holes N and the geometrical parameters h/W and D/W determine which of the two modes, Euler or alternating, occurs at the onset of the instability from the unbuckled compression branch.
The critical strain at which buckling occurs in the Euler mode and the critical strain at which buckling occurs in the alternating mode are compared in figure 4 as functions of the number of holes in the column N. At large N, the cr ∼ 1/N 2 behaviour of the Euler mode and cr ∼ const. behaviour of the alternating mode means that sufficiently long columns will always buckle first in the Euler mode. For short columns, the frustration of the alternating mode by the clamped boundary conditions (leading to increased cr ) has the same result; sufficiently short columns also buckle in the Euler mode. For the shortest column with an even number of holes, N = 2, we have found no column geometry where the alternating mode bifurcation (figure 1b right) occurs at a lower strain than the Euler mode (figure 1b centre). However, at intermediate column lengths, the first buckling instability of the column can be in the alternating mode.
The first mode to become unstable is illustrated for such a column of intermediate length (N = 6 holes) in figure 5, as a function of both hole size D/W and hole spacing h/W. Each point on the figure corresponds to a numerical computation of solutions along the unbuckled compression branch, in which the critical strains at which the Euler and alternating bifurcations occur are recorded. Crosses show the region in parameter space where the Euler mode bifurcation occurs at a lower strain than the alternating mode bifurcation, and dots the region where the alternating mode occurs at a lower strain. The solid black line separating these two regions indicates the column geometries for which the two bifurcations occur at the same strain, and is obtained by plotting the zero contour of the difference between the critical strains for the Euler and alternating bifurcations. The alternating mode occurs at a lower strain than the Euler mode in regions of this geometric parameter space where the ligaments are relatively thin both at the side walls of the column and between the holes. In other regions, the Euler mode occurs at a lower strain.
Comparison of figure 2a and 2b suggests that an increase in h/W promotes the localized sliding regime of the Euler mode. This is supported by the results in figure 5, in which the colours of the crosses indicate the ratio between critical strains of the second to first Euler mode bifurcations; as noted previously, these occur at very similar strains in the localized sliding regime. For a given hole size D, figure 5 indicates that the alternating mode is preferred for a limited range of hole spacings h/W. For holes spaced more widely than this range, the ligaments separating holes are thick, and a localized Euler mode in the sliding regime is preferred, where these ligaments separating holes do not deform significantly. For holes spaced more closely than this range, the ligaments at the edge of the column are thick compared to those separating holes, and a global Euler buckling mode is preferred, in which the thicker ligaments adjoining the edge of the column do not buckle, except over the length of the whole column.
Our simulations suggest that as the column length increases, the region of the alternating mode in the phase space of the unit cell geometry (dots in figure 5) shrinks and ultimately disappears at around N = 16 holes. In a sufficiently long column, the first bifurcation to occur during compression, therefore, always corresponds to the Euler mode. However, the alternating mode may still play a role in the subsequent deformation, due to the presence of secondary bifurcations.

(e) Secondary bifurcations
Our study, so far, has explored states arising at the primary bifurcation that occurs at the smallest strain on the unbuckled compression branch. We now study the secondary bifurcations that occur on the first branch to bifurcate from the unbuckled compression branch, either the alternating mode or the first Euler mode. In both cases, the primary bifurcation at the onset of the instability is supercritical. The continued compression of the column, therefore, initially results in stable solutions with a single broken symmetry that have bifurcated from the unbuckled compression branch. On these already bifurcated solution branches, we find that secondary bifurcations can exist, of the type first described by Bauer et al. [16]. These secondary bifurcations result in solutions with a second broken symmetry; in the example shown in figure 6a, a column with eight holes exhibits both the lateral buckling and the alternating mode, with both horizontal and vertical symmetries broken.
We examine the behaviour of secondary bifurcations in the neighbourhood of a point in the geometric parameter space at which the Euler and alternating modes bifurcate from the unbuckled compression branch at the same critical strain. Choosing N = 8 and D/W = 0.7525, this point occurs at h/W = r 0 ≈ 0.89426, and is marked by the vertical dashed line in figure 6b that separates the two regions of parameter space, in which either the Euler mode (h/W > r 0 ) or the alternating mode (h/W < r 0 ) occurs at the lowest strain on the trivial branch. By exploring the secondary bifurcations, we investigate the nature of the stability exchange between the two modes.
In numerical calculations, we impose increasing applied strain until we reach the first bifurcation to occur on the unbuckled compression branch, at which point we perturb the numerical solution with the eigenvector associated with the bifurcation. This bifurcation is supercritical, so by increasing the applied strain further we find solutions on the bifurcated branch until the secondary bifurcation is located. The critical strains of both the primary and secondary bifurcations are shown in figure 6b. For h/W < r 0 , the alternating mode bifurcation (A) occurs at the lowest strain, and the solutions with broken top-bottom symmetry that originate from this bifurcation are subject to the secondary bifurcation (AE) where the solutions also buckle laterally, breaking the left-right symmetry (figure 6c). For h/W > r 0 , order of the bifurcations is reversed: here the Euler mode bifurcation (E) occurs at the lowest strain, and the laterally buckled solutions that originate from this bifurcation subsequently become unstable to an alternating mode (breaking the top-bottom symmetry) at the secondary bifurcation (EA). At h/W = r 0 , the Euler and alternating mode primary bifurcations occur at the same strain, and all the bifurcation points, both primary (A and E) and secondary (AE and EA), then coalesce to meet at the same critical strain, forming a multiple primary bifurcation point on the unbuckled compression branch. This feature has been identified and explained by previous studies on secondary bifurcations [16,17]. The coalescence of critical strains of the primary and secondary bifurcations gives insight into the behaviour of columns when they are compressed beyond the initial instability; for columns with geometry where the alternating and Euler modes occur at nearly the same critical strain (i.e. those close to the black line in figure 5), only a small additional compression is required for both top-bottom and left-right symmetries of the column to be broken.

Theoretical modelling of localized buckling modes
Motivated by the existence of localized buckling modes and their importance in the overall behaviour of the column, we develop a theoretical model for the stress at the onset of bifurcation in a short section of column with two holes and with periodic boundary conditions imposed on the top and bottom edges (figure 7a). These periodic columns with a period of two holes exhibit instability to both the alternating mode (figure 3b, column ii) and to the localized sliding regime of the Euler mode (figure 2d, column iv). Furthermore, as we have shown, the critical buckling strains of the periodic columns are in quantitative agreement with the appropriate localized buckling regime of non-periodic columns of finite length: localized Euler buckling in short columns with hole spacing larger than the width (figure 2b), and alternating mode buckling in long columns (figure 3a). The theory presented in this section is, thus, applicable both to periodic columns and to the localized buckling regimes of non-periodic columns, where the buckling properties are independent of the column length.
When the holes in the column are similar in diameter to the column width or hole spacing, the deformation is localized, occurring through bending of the resulting thin regions of elastic material, or ligaments, highlighted in figure 7b. The thicker uncoloured parts of the structure remain comparatively rigid. By exploiting the large aspect ratio of the ligaments, and the comparative rigidity of other parts of the structure, we can obtain algebraic expressions for the critical buckling force for each mode. There are two types of ligament present in the system; those separating two holes, which we denote s, and those adjoining the side wall of the column, which we denote w (figure 7b). The minimum width of these ligaments is a s = h − D and a w = 1 2 (W − D), (3.1) respectively, and we develop the model in the asymptotic limit where a s D and a w D.
The width of the ligament separating two adjacent holes w s (figure 7a) is given as a function of the coordinate along its length x by When a ligament is of a similar thickness to its thinnest point (w s = O(a s )), we find from (3.2) that the ligament has lengthscale x = O( √ Ra s ) O(a s ), implying that the ligament is much longer than it is wide. We can, therefore, model the ligaments as thin Euler-Bernoulli beams of nonuniform width, with the width of the two types of ligament given by (3.2) and (3.3).

(a) Alternating mode
Motivated by the results of finite-element computations, which indicate that when a s , a w R, the strain energy is concentrated almost entirely in the thin ligaments, we assume that a column buckling in the alternating mode can be modelled as a system of rigid sections, formed from thick regions of elastic material, that are connected by the thin, flexible ligaments. Each ligament acts as a hinge with a torsion spring (figure 8a).
We start by calculating the stiffness of these torsion springs from the ligament geometry, using Euler-Bernoulli beam theory. The equation for small deformations of an Euler-Bernoulli beam [15] is where E is Young's modulus of the material, x is the distance along the beam, φ is the perpendicular displacement of the beam, where M is the bending moment. As the length of each ligament is much smaller than the column width or hole spacing, the moment M is constant (to leading order) along the length of the ligament. We consider a ligament with a general parabolic width profile where the constants a and x 0 parametrize the beam geometry in either type of ligament (see (3.2), (3.3)). Using (3.7) in (3.5) reduces (3.6) to (3.8) Integrating again and seeking a symmetric solution with dφ/dx = 0 at x = 0, we obtain where X = x/x 0 . We note that namely that the beam becomes straight at large |x|. For small deflections of the beam, we may, therefore, define a bending angle for the entire beam, as From (3.10), the moment exerted by the beam can be written as where torsion coefficient κ is Setting a = a s and x 0 = √ Ra s for the s ligaments (3.2), and a = a w and x 0 = √ 2Ra w for the w ligaments (3.3), we obtain the expressions for torsion coefficient of the spring formed by these ligaments as respectively. With these coefficients evaluated, we can now use (3.12) to relate the angle by which each ligament is bent, θ , to the moment it exerts M.
A single periodic unit of our model of the holey column consists of four w ligaments and two s ligaments (figure 8a). When this periodic unit is compressed by a force F and each of these ligaments is bent by an angle θ, the height of the periodic unit is reduced from 2h to 2h cos(θ/2). The potential energy of this unit is, therefore, and the equilibria lie at the stationary points of this energy, where The change in stability of the equilibrium branch at θ = 0, a bifurcation, occurs at a critical force, To validate this prediction, we compare it to the critical force at the bifurcation obtained from finite-element numerical calculations of a two-hole column with periodic boundary conditions applied to the top and bottom of the column. The compression of the system is controlled by altering the vertical offset between periodic units in the deformed configuration, and the force F measured by integrating the traction over the upper boundary. For these computations, an unstructured triangle mesh was used to discretize the domain, with elements concentrated in the thin ligaments. Our theoretical prediction (3.17) closely matches the numerical calculations, particularly as the hole diameter approaches the width of the column (figure 9a).

(b) Localized Euler buckling
Euler buckling of the two-hole periodic system occurs in the localized sliding regime. In this mode, only the w ligaments adjoining the edges of the column undergo deformation, whereas the s ligaments separating holes are undeformed ( figure 8b). Here, we derive the critical force required for this Euler mode to be realized. The buckling analysis for this problem is very similar to that for classical Euler buckling (e.g. [18]). Each w ligament supports half the compression force F, and so from (3.6) we have where the moment is the product of the force F/2 and the perpendicular distance φ. Evaluating I using the geometrical parameters corresponding to the w ligaments (3.3), as in (3.8), gives (3.19) where X = x/ √ 2Ra w . In the sliding regime of the Euler mode, the thicker regions of the beam between ligaments do not rotate, and so the boundary conditions are 'clamped', but able to slide in a direction perpendicular to the load, dφ dX → 0 as X → ±∞.

Asymptotic model for the onset of bifurcation in a two-dimensional cellular sheet
The asymptotic theory for bifurcation presented in the previous section extends naturally to twodimensional cellular solids perforated by a square lattice of circular voids, where the geometry is described by a single non-dimensional parameter, a/h, representing the minimum distance between holes divided by the spacing between hole centres ( figure 11). When compressed uniaxially, such two-dimensional cellular materials exhibit an alternating-mode bifurcation at a critical strain [5,7], similar to the bifurcation we have demonstrated in columns. Changing the pattern of holes in the solid, for example to a triangular lattice, leads to a range of related instabilities [19], which can be exploited in the fabrication of metamaterials with desired mechanical properties [4]. As before, we consider only two-dimensional buckling modes, neglecting out-of-plane buckling, which may occur in sufficiently thin plates (e.g. [20]). The onset of the buckling instability in cellular sheets has been studied numerically, under general in-plane loading in the case where a/h = 1/2 [21], and under uniaxial compression for more general a/h [7]. Here, we derive algebraically how the critical force F cr and strain cr vary with the hole size a/h, under uniaxial compression of the solid, when the holes are relatively large compared with the elastic regions separating them (a/h 1). Our model for alternating-mode buckling in elastic cellular sheets perforated with large holes is similar to that in holey columns (figure 8a), in that we assume that the thin ligaments bend, while thicker regions of the structure rotate without significant deformation. This approach is supported by the concentration of the strain energy within the ligaments in numerical calculations. Furthermore, the approach predicts that Poisson's ratio of the cellular sheet approaches −1 for large holes, which is consistent with numerical calculations [7].
Defining the hole radius as R = (h − a)/2, we find from (3.13) that the bending strength of a single ligament is A single periodic cell of two-by-two holes contains eight such ligaments, so the potential energy of a single periodic cell of this system with an applied compression force F is As before, the equilibria lie at the stationary points of the energy, and a bifurcation on the solution branch θ = 0 thus occurs at at a critical force F cr , where In addition to this critical force, the symmetry of the two-dimensional lattice also allows us to calculate easily the critical strain at bifurcation. We do this by evaluating the effective Young's modulus of the cellular solid under small loads, before the bifurcation occurs. The solid can be divided into two regions: a region of length and width O(h) (shaded grey in figure 11a), and a region corresponding to thin ligaments that are parallel to the direction of compression, of length O((aR) 1/2 ) and width O(a) (shaded red in figure 11). The linear deformations resulting from applying a force to these regions scale 2   parts. These scalings suggest that, for small deformation when a/h 1, the thick regions may be considered rigid and the deformation of the whole periodic unit predominantly results from compression of the thin ligaments that are parallel to the applied strain. Each periodic repeating unit of the solid consists of four such ligaments, each of which supports a force F/2. When compressed under a force F/2, the length change of a single column of the two-dimensional lattice is The total compression of the repeating periodic unit is given by 2 l, and so the relative strain is This relationship between cr and F, which holds immediately prior to a bifurcation, allows us to convert the critical force for the bifurcation (4.4) into a critical strain, cr =

9
a h 2 . (4.8) The asymptotic calculations for the critical force (4.4) and strain (4.8) are in good agreement with numerical computations of the compression of a neo-Hookean material for thin-walled two-dimensional lattices ( figure 12). In these computations, we simulate a single periodic unit of four holes, with periodic boundary conditions applied on both the upper/lower and left/right boundaries. As with the periodic columns, the amount of compression is imposed by specifying the periodic offset in the vertical direction. At the left and right periodic boundaries, we leave the horizontal periodic offset undetermined, and instead specify that the integral of the stress over the sides of the periodic unit is zero. The periodic cell is, therefore, free to expand or first bifurcation, these ligaments, therefore, bend outwards and compress, which complicates the evaluation of the effective Young's modulus.  11b). For a/h = 1/10 (as illustrated in figure 11), the asymptotic prediction of the critical stress (4.4) is within 5% of the value obtained from numerical computations, with closer agreement at smaller values of a/h. Our asymptotic prediction from (4.8) of cr = 0.03 for a/h = 0.13 is also in good agreement with the change in stress-strain behaviour observed in the periodic numerical calculations of Mullin et al. [5], and is similar to their experimental measurement of cr = 0.04 for a sheet with 10 × 10 holes.

Conclusion
We have used a previously validated plane-strain, finite-element model [12] to investigate the inplane buckling behaviour of a holey column under uniaxial compression, focusing on how this behaviour varies with the column geometry; specifically the number, size and spacing of the holes. Insight into the detailed mechanics of buckling is provided by a new asymptotic theory, which quantitatively predicts the critical stress at which a column with large holes buckles, without arbitrary or fitted parameters. We find that for all columns, the first buckling mode to occur as strain is increased, i.e. the mode with smallest critical strain, is always either (i) a Euler mode, in which the column buckles in a direction perpendicular to the applied compression; or (ii) an alternating mode, in which adjacent holes are deformed into ellipses with orthogonal major axes, but the column remains straight overall.
The alternating mode is generally the first mode to become unstable under compression in columns of moderate length (≈4-16 holes) if the holes are relatively large, diameter greater than ≈75% of both the column width and hole spacing. For other column geometries, the Euler mode is the first mode to become unstable either via global buckling on the scale of the column (favoured in long columns or columns with small holes) or localized buckling on the scale of a hole (favoured when the ligaments adjoining the column edge are thinner than those separating holes). We have also studied the exchange of stability between the two modes by exploring the secondary bifurcations in the system.
Both Euler and alternating modes exhibit regimes of localized buckling, in which the buckling of the whole column occurs through buckling of thin ligaments around each hole, and at a critical strain that is nearly independent of the column length. This localized buckling of ligaments is captured by our asymptotic model, which provides an accurate quantitative description of the buckling when the holes are large. The model illustrates how the different types of buckling arise through deformation of the ligaments. The alternating mode involves bending of both those ligaments separating holes and those adjoining the column edges; this explains why the first buckling mode observed is the alternating mode only when both types of ligament are thin. By contrast, the Euler mode in its localized 'sliding' form involves buckling only of the ligaments adjoining the column edges, whereas in its global form the Euler mode involves the compression or extension, but not buckling, of these ligaments. The buckling behaviour of the column, and its dependence on the column geometry, can be understood and predicted from ligament deformations, as we demonstrate quantitatively for localized buckling in the alternating and Euler modes, equation (3.23). One reason for the success of the theory is that the buckling is strain-dominated and occurs at such small strains that material nonlinearities, and hence details of the constitutive law do not play a prominent role.
As shown in §4, these asymptotic ideas can be extended successfully to other elastic structures. Further extensions, for example, to more general grid structures [22], could allow this type of modelling to be used as a predictive tool for buckling in a wider range of perforated structures and metamaterials [4]. Such structures may be unstable to out-of-plane buckling modes [20], in addition to the in-plane modes we describe.
We find that, counterintuitively, the critical strain of Euler buckling in a column can be significantly increased by the presence of holes in a column. One avenue for future investigation is a more detailed study and experimental validation of this potentially useful phenomenon. Although the compressed column bifurcates first from the unbuckled compression branch only through the Euler or the alternating modes, we have shown that other bifurcations can occur at only slightly greater strains, notably higher-order Euler modes of buckling in the sliding regime. This invites the possibility of custom tuning of buckling states, through the addition of small geometrical changes or other imperfections in the system, or through dynamical effects. Most studies of buckling in cellular structures, to date, have considered only quasi-static compression of samples, but the dynamical compression is still largely unexplored, and could be of more practical relevance in large structures. A step in this direction was taken by Box et al. [23], who studied buckling in elastic and plastic cellular materials under dynamic compression, showing a considerable difference in the response of the latter depending on the rate of compression. This type of dynamic compression could, therefore, be explored with the aim of accessing modes other than the first Euler mode or the alternating mode in the holey column.