Characteristics of mechanical metamaterials based on buckling elements

Accepted Manuscript

Characteristics of mechanical metamaterials based on buckling elements Claudio Findeisen, Jorg ¨ Hohe, Muamer Kadic, Peter Gumbsch PII: DOI: Reference:

S0022-5096(16)30788-8 10.1016/j.jmps.2017.02.011 MPS 3068

To appear in:

Journal of the Mechanics and Physics of Solids

Received date: Revised date: Accepted date:

30 November 2016 17 February 2017 18 February 2017

Please cite this article as: Claudio Findeisen, Jorg ¨ Hohe, Muamer Kadic, Peter Gumbsch, Characteristics of mechanical metamaterials based on buckling elements, Journal of the Mechanics and Physics of Solids (2017), doi: 10.1016/j.jmps.2017.02.011

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

Highlights • Metamaterials composed of unstable elements are studied on different length scales. (85)

• Energy dissipation relies on transformation of external work to vibrations. (77)

• Dissipation shows neither viscoelastic nor plastic charac-

CR IP T

teristics. (69)

• Effective properties can be programmed within the same graded microstructure. (77)

• Programming is triggered by loading history only with-

AC

CE

PT

ED

M

AN US

out additional manipulation. (83)

1

ACCEPTED MANUSCRIPT

Characteristics of mechanical metamaterials based on buckling elements Claudio Findeisena,b,1 , J¨org Hoheb , Muamer Kadicc,d , and Peter Gumbscha,b a Institute

CR IP T

for Applied Materials, Karlsruher Institute of Technology (KIT), 76128 Karlsruhe, Germany b Fraunhofer Institute for Mechanics of Materials IWM, 79108 Freiburg, Germany c Institute of Applied Physics, Karlsruhe Institute of Technology (KIT), 76128 Karlsruhe, Germany d Institut FEMTO-ST, CNRS, Universit´ e de Bourgogne France-Comt´e, 25044 Besancon Cedex, France

Abstract

Metamaterials are composed of structural elements and derive their properties mainly from the inner structure of the elements, rather than the properties of their constituent material. By designing an unstable structural element as the building block of a metamaterial, many interesting effective material properties can be obtained. The deformation and dissipation mechanisms of such a material built from unstable structural elements is studied in detail. To do so a combination of analytical, semi-analytical, and numerical models

AN US

are applied to a single buckling element, a periodic cell, and finite size combinations of buckling elements including gradients in the properties of the building blocks. This not only provides insight into the micromechanics and the resulting effective behavior of such metamaterials, but also makes them accessible on the different relevant length scales. A metamaterial built from these building blocks shows programmable or switchable properties and can display energy dissipation with fully reversible deformation, distinguishing it from plastic materials, and timescale independent behavior, distinguishing it from viscoelastic materials.

M

Keywords: buckling, microstructures, inhomogeneous material, stability and bifurcation, metamaterial

ED

1. Introduction

The concepts of structural design of materials have recently been widened significantly by the introduction of mechanical

PT

metamaterials. Metamaterials derive their effective material properties mainly from the structure of their elementary building

CE

blocks rather than the properties of the materials they are build from. The book of Milton [1] gives a nice introduction to the general idea of metamaterials and [2–6] give an overview of re-

AC

cent developments in this field. Structural optimization of metamaterials can for example be performed with respect to the effective elastic properties [7, 8], heat transfer capabilities [9], or energy absorption and energy dissipation [10–12]. The design space for metamaterials is currently mostly limited by available fabrication methods. However, during the last years, several new fabrication methods, such as direct laser writing, have enabled the production of complex, three dimensional, nearly ar1 Corresponding

bitrary inner structures with reasonable effort and in practical relevant dimensions. So far, most metamaterials are designed with respect to linear properties and the nonlinear regime is only partially explored. Metamaterials may, however, achieve unusual nonlinear properties from geometrically nonlinear deformation behavior or structural instabilities. Many interesting effective properties can be achieved by such unstable structures. For example the classical stiffness bounds of composites (Voigt [13] and Reuss [14]), which are derived, assuming a positive definite strain energy, can be overcome by unstable inclusions [15, 16]. The Poison ratio ν, which is limited for an isotropic elastic solid (−1 < ν < 0.5) is in fact unbounded for an unstable structure [17]. Furthermore, unstable structures can be used to generate programmable material properties. Florijn et. al. [18], for example, presented a way to program the effective stress-strain behavior induced by the pattern transformation of a periodic ar-

author

Preprint submitted to Journal of the Mechanics and Physics of Solids

February 21, 2017

ACCEPTED MANUSCRIPT

rangement of elliptical and cylindrical holes. A more recent idea is to design metamaterials for damping and energy absorption applications. Besides optimizing for the specific dissipated energy, the unstable inner structure can be designed such that the deformation is reversible and that the maz

terial can thus absorb multiple impacts without losing its func-

z

y

x

x

tionality [18–23]. This is possible by designing a metamaterial

Figure 1: Unit cell of the investigated inner structure. Hexagonal frame and

in which energy dissipation is driven by elastic buckling.

CR IP T

connection elements (grey) between the 2 × 3 parallel buckling elements (blue).

The conceptual idea of energy dissipation through instabili-

Denoted values di and hi are the diameter and height of the individual bars.

ties goes back to the work of Lakes [24] who first presented the

(with permission from [19])

capabilities of materials with negative stiffness inclusions and investigated the damping capabilities of such materials. Many

cussion and summary.

experimental works in the recent few years have used the buck-

2. Investigated unit cell

AN US

ling principles as the main design parameter [18, 20, 22, 23].

More recently Frenzel et. al. [19] designed, produced and

The unit cell of the investigated structure is shown in Fig-

characterized a lightweight unstable three dimensional micro-

ure 1. This cell is built from solid bars and has been fabri-

lattice, showing both programmable properties and repeatable

cated and experimentally characterized in [19]. The structure

and timescale independent energy dissipation. To contribute to

is three-dimensional, but buckling only occurs under compres-

a better insight into unstable metamaterials, this paper presents

sion normal to the hexagonal base structure. On the top and on

M

a detailed analytical and numerical study of the deformation and dissipation mechanism of the unstable unit cell presented

ED

in [19]. The behavior of the inner structure is investigated by

the bottom of this hexagonal basis 3 × 2 buckling elements of a

sinusoidal shape are attached. This sinusoidal shape is chosen to keep the local strain within the buckling elements at low lev-

analyzing a single buckling element, a periodic cell and the in-

els before during and after buckling. The buckling elements are

teraction between multiple unstable elements. This multistep

connected to the next row of buckling elements by a set of star-

PT

approach is intended to provide a better understanding of the

shaped connection elements. In the horizontal direction these

mechanical response of such unstable structures. Although, the

are also connected to the connection elements of the neighbor-

presented study is performed on a specific unit cell, most of the

CE

ing cells. Both the hexagonal frame and the connection ele-

results are valid for unstable unit cells in general.

ments are designed to prevent other, unintended deformation

The outline of the paper is as follows, in Section 2 the inves-

modes and instabilities such as lateral expansion, shear defor-

AC

tigated unit cell is described in detail. In Section 3 the behavior

mation, and rotation around the vertical axis of the center con-

of the unit cell is investigated by applying a simple bar model

nection of the buckling elements.

and analyzing a periodic cell. In Section 4 a semi-analytical model is used in combination with a full finite element model

3. Behavior of the inner structure

to analyze the macroscopic behavior of a finite combination of unstable cells. While the main focus is on an homogeneous

The behavior of the inner structure can be investigated using

distribution of buckling cells we also introduce a gradient into

a combination of a very simple bar model and a periodic cell

the overall structure by systematically changing the structural

model considering the real geometry of the unit cell.

characteristics of the individual cells. This is followed by a dis3

ACCEPTED MANUSCRIPT

F1 , u

ical equivalence, which is enforced by the Hill-Mandel condition [25]

k1

h0 k3

1 V

2v k2

l0

Z

Ω0

P · dFdΩ = P∗ · dF∗ ,

(4)

where F and P denote the real structural strain and stress respectively. Equation (4) is required to hold for equivalent de-

Figure 2: Simple bar model of one buckling element (as investigated in Sec-

formations [26], defined as Z 1 dFdΩ = dF∗ . V Ω0

tion 3). Denoted vales: u, v axial and vertical displacement, F1 external applied

CR IP T

force, h0 , l0 overall height and length respectively, and ki stiffness of the individual elements.

(5)

So simply speaking the incremental strain energy of the discrete

3.1. Simple bar model

inner structure needs to be the same as the energy described by

The simple bar model is composed of bars and springs as

the effective quantities F∗ and P∗ , when the deformation incre-

shown in Figure 2. The stiffness k1 is taken as the reference

ment is the same in an volume average sense.

To fulfill condition (5), an explicit relation between the ef-

contraction of the structure. The spring with stiffness k3 is intro-

fective deformation gradient F∗ and the deformation of the pe-

AN US

value, k2 imitates the resistance against lateral expansion and duced to capture the bending strain energy of the real buckling

riodic cell is required. For this, finite tangents ∆X of the refer-

elements, and of course h0 is related to the buckling element

ence configuration are mapped to finite tangents of the current

height h1 . Indeed, this is a very rough approximation of the real

configuration ∆x by

structure. However, this model can be solved analytically and

∆x = F∗ ∆X + ∆w,

M

thus lends itself to the systematic investigation of all possible

characteristics of such structures. The potential energy of the

1



application of the Green-theorem, leads to Z Z ∇∆wdΩ = ∆w ⊗ ndΩ = 0.

2k1 ∆s2 + (2v)2 k2 + k3 u2 , (1) q q with ∆s = s − s0 = (l0 − v)2 + (h0 − u)2 − l02 + h20 (2) 2

PT

E=

where ∆w (called microfluctuation) contains all higher order deformation modes. To fulfill condition (5) by this mapping,

ED

described system is given by

(6)

(7)

∂Ω0

Ω0

A common, physically motivated way to fulfill this, is to require

following equilibrium relation      ∂E  F   ∂u   1   ∂E  =   . 0

periodicity of the microfluctuation on the boundary ∂Ω (see for

AC

CE

and thus the behavior of the system is obtained by solving the

example [27]). For this, (3) ∆w|∂Ωi+0 = ∆w|∂Ωi−0

∂v

(8)

is required for corresponding points on opposite boundary sur-

3.2. Periodic cell model

i− i i faces ∂Ωi+ 0 and ∂Ω0 with normals e and −e respectively (see

To gain deeper insight into the deformation modes of the

Figure 3). By evaluation of Equation (6) for opposite surfaces

real inner structure, periodic basis cells (depicted in Figure 3)

and taking into account Equation (8), the following boundary

are analyzed. To measure the effective behavior of such a peri-

condition is obtained

odic arrangement of cells with one cell Ω0 , the effective energy conjugated first Piola-Kirchhoff stress P∗ and deformation gra-

∗ FiJ =

∗

dient F are introduced. Their values are defined by a mechan4

∆xi |∂Ω0J+ − ∆xi |∂Ω0J− lJ

,

∀ (i, J) ∈ (1, 2, 3)2 ,

(9)

ACCEPTED MANUSCRIPT

This periodic cell model is solved using finite elements,

l1 ∂Ω2+ 0 ∂Ω1− 0

∂Ω1+ 0

y x

l2 −e1

e3

∂Ω3+ 0

based on the Timoshenko hypothesis. An element size of b/15

∂Ω2+ 0

d1

h1

whereby the structure itself is modelled using beam elements

z

is used where b is the hexagonal cell size (see Figure 1). The y

−e2

e1

e2

material behavior is modeled using the Saint-Venant Kirchhoff

∂Ω2− 0 ∂Ω2− 0

∂Ω3− 0

model with material parameters E = 2GPa and ν = 0.42, which

−e3

are typical values for the linear elastic regime of the polymer

CR IP T

used in [19]. Since the structure exhibits instability phenomena, a con-

Figure 3: Periodic lattice, cut through x − y plane (left), cut through y − z plane

(right), and cut-out thereof used for the periodic cell model. Buckling elements

venient stabilization of the static boundary value problem is

are depicted as dotted lines. ∂Ωi+/i− is the surface of the periodic cell in the 0

needed in order to ensure that the potential energy is minimized.

i ∈ {1, 2, 3} direction with surface normals +ei and −ei respectively, and l J are

This is done by solving the following eigenvalue problem (See

the overall dimensions of the periodic cell.

for example [28]) or, equivalently =

∆ui |∂Ω0J+ − ∆ui |∂Ω0J− lJ

(Kl + λi Ku ) Φi = 0

,

2

∀ (i, J) ∈ (1, 2, 3) ,

AN US

HiJ∗

(10)

(13)

for the reference configuration Ω0 and applying a small initial

are the nine components of the displacement gradient

perturbation in the form of the first ten eigenvectors Φi . By

and l J are the distances between opposite surfaces (see Figure

monitoring the negative eigenvalues of the stiffness matrix it is

3).

ensured that this initial perturbation is sufficient for a stabiliza-

where

HiJ∗

ED

the boundary value problem, which reads Z P · δFdΩ = δWext = VP∗ · δF∗ ,

M

The effective stress P∗ can be obtained by the weak form of

Ω0

tion.

The periodicity assumed in condition (10), implies an equally deforming neighbourhood of Ω0 . Since an instability leads to

(11)

a localization of the deformation, this assumption is not justified beyond the first instability. This problem is related to the

∗ since the only external forces are conjugated to the applied FiJ .

assumption of scale separation in classical homogenization the-

fill condition (4). Thus, it is straightforward to obtain the ef-

ory. This assumption can be relaxed in a certain range by for

fective stresses P∗ which are simply the reaction forces of the

example including higher gradients or independent microdefor-

∗ additional degrees of freedom FiJ multiplied by the volume of

mations [29–31]. However in case of the microstructure stud-

the periodic cell.

ied here, this is questionable since scale separation is violated

CE

PT

The multiplication with the volume V is needed in order to ful-

due to the localization within the unitcell. Thus the presented

bitrary point is fixed to zero and the deformation gradient is

method should be seen only as a tool to characterize the unit cell

required to be symmetric. The latter condition suppresses rigid

and not as a method to derive some effective model, which, by

body rotations of Ω0 , as can be seen by virtue of the polar de-

the way, would also encounter difficulties due to non-unique so-

composition of the deformation gradient

lutions and subsequent mesh dependence using finite elements

AC

To suppress rigid body motions, the displacement of an ar-

 1/2 F∗ = R∗ U∗ = R∗ F∗T F∗ = R∗ F∗ ,

[32, 33]. However, the applied approach can still be used to de(12)

termine the first bifurcation point upon loading, and, assuming all cells are buckled, the bifurcation upon unloading.

and thus R∗ = 1. Furthermore the deformation gradient equals the right stretch tensor F∗ = U∗ . 5

ACCEPTED MANUSCRIPT

·10−2 k2 /k1 = 0.1 h0 /l0 = 0.5

4

stable unstable k3 /k1 = 0.01 h0 /l0 = 0.5

k2 /k1 = 2

4

[−]

1

0

F1 k1 l0

0.5

0

k2 /k1 = 0

h0 /l0 = 0.25

0

k3 /k1 = 0

2

−2

0.2 0.4 0.6 0.8 norm. displacement 2hu0 [−]

1

0

0.2 0.4 0.6 0.8 norm. displacement 2hu0 [−]

(a)

(b)

CR IP T

0

h0 /l0 = 0.6 1

norm. force

norm. force

F1 k1 l0

k3 /k1 = 0.04

2

[−]

1.5

[−] F1 k1 l0

norm. force

·10−2

stable unstable

k2 /k1 = 0.1 k3 /k1 = 0.01

3

−1

·10−2

2

stable unstable

1

0

0.2 0.4 0.6 0.8 norm. displacement 2hu0 [−]

1

(c)

Figure 4: Normalized force as a function of normalized displacement for the simple bar model shown in Figure 2. Dotted part shows the unstable (inadmissible) part of the equilibrium curve (for simplicity only depicted for the reference curve). (a) Variation of the stiffness k3 to control the potential energy in the buckled

AN US

configuration. (b) Variation of the buckling element height h0 . (c) Variation of the stiffness k2 to control the lateral expansion. Reference (red curve) is for k3 /k1 = 0.01, h0 /l0 = 0.5, and k2 /k1 = 0.1.

branches. These two branches correspond to the buckling of

3.3. Results

the two sets of buckling elements within a single unit cell. The

rized for a variation of the parameters k2 , k3 and h0 . The vari-

critical force of the second buckling event is slightly higher

ation is always performed with respect to the same reference

than the first one because of a slight lateral contraction which

(red curves). In general the force-displacement curve is non-

occurs after the first element has buckled. In the following,   only the average critical stresses Pmin := 1/2 P1min + P2min and   Pmax := 1/2 P1max + P2max are reported. The influence of the

M

In Figure 4 the behavior of the simple bar model is summa-

ED

monotonic showing a single minimum and a single maximum. Between these extreme values there is an unstable part (dotted

geometry of the periodic cell is investigated by a variation of

to as unstable equilibrium branch in the following. The two re-

the parameters d1 , d2 , d3 , and h1 in Figure 5 - 7.

PT

line) with a negative definite stiffness. This part will be referred

First, the influence of the stiffness k3 of the simple bar model

curve) will be called the primary and secondary stable equilib-

and the corresponding buckling element diameter d1 of the pe-

CE

maining stable parts with a positive definite stiffness (full red

riodic cell is discussed. For k3 = 0 the structure is symmetric

path are inadmissible, meaning, the structure will not stay at

with respect to the horizontal axis, therefore the force-displace-

such a point and will buckle to a stable equilibrium configu-

ment curve is symmetric with respect to the zero force level

ration, which fulfills the same boundary condition. However,

(See Figure 4a). This zero force level occurs when both bars are

depending on the boundary condition, such an unstable buck-

straight (u/2h0 = 0.5). For k3 > 0 an additional strain energy is

ling can be suppressed, for example by prescribed displacement

stored during the deformation, therefore the symmetry is bro-

boundary conditions. In this structure buckling therefore only

ken and the zero force level is shifted to higher displacements.

occurs under force boundary conditions.

Furthermore force minimum and maximum are raised by in-

AC

rium branches. In general the unstable parts of the equilibrium

Figure 5a shows a generic result of the periodic cell model

creasing the stiffness parameter k3 and the intersection with the

which is characterized by two stress minima and maxima (Pimin

zero force level vanishes. Above k3 /k1 = 0.03 the behavior

and Pimax , i = 1, 2), and the corresponding unstable equilibrium

switches from non-monotonic unstable to monotonic stable. 6

ACCEPTED MANUSCRIPT

2

1

0

1

0.5e-5

0

0.04

0.06

0.08

0.5 0

0.1 0.2 0.3 0.4 displacement gradient −H33 [−]

= d1 snap through = d1 snap back = 2d1 snap through = 2d1 snap back

0.6 0.4 0.2 0

Programmable 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18

0.5

di di di di

0.8

P33 E

1.5

detail

1 Pmin

2 Pmin

0

1e-5

·10−4

1

critical stress

P33 E

critical stress

norm. stress (1.PK)

3

1.2

= d1 snap through = d1 snap back = 2d1 snap through = 2d1 snap back [−]

2

di di di di

buckling diameter

(a)

d1 b

Programmable −0.2 0.12 0.14 0.16 0.18 0.2 0.22 0.24 0.26 0.28

CR IP T

4

2.5

P33 E

[−]

5

·10−4 Eq. path periodic cell (stable) Eq. path periodic cell (unstable) 2 Pmax 1 Pmax

[−]

·10−5

buckling height

[−]

(b)

h1 b

[−]

(c)

Figure 5: Results from the periodic cell model. (a) Stress-strain curve showing a non-monotonic behavior with the critical stresses P1max , P2max , P1min , and P2min . (b) Average critical stresses as a function of the buckling diameter d1 for d2 /b = d3 /b = d1 /b and h1 /b = 0.2 (solid lines), and for d2 /b = d3 /b = 2d1 /b and h1 /b = 0.2

·10−5

[−]

the periodic cell. Here, bending strain energy is stored during ter d1 , both critical stress minima and maxima increase (see

M

solid lines in Figure 5b). By further increasing the diameter

above d1 /b = 0.15 the behavior switches to monotonic stable. However d1 /b can not be chosen close to zero, since the overall

norm. crit. stress

P33 E

the deformation, thus increasing the buckling element diame-

5 4 3 2 1

0.5

1

1.5

2

ED

frame diameter

stiffness would tend to zero. Consequently there is no negative

·10−3

Pmax - av. max stress Pmin - av. min stress

max lateral expansion H22 [−]

The simple bar model helps to understand the behavior of

2.5 d2 d1

3

8 6

10

8

4 6

2 0 0.5

[−]

(a)

critical stress in the reported parameter combination because

max Lat. expansion H22 Rotation ϕmax z

1

1.5

frame diameter

2 d2 d1

Rotation ϕmax [◦ ] z

AN US

(dashed lines). (c) The same as in (b), but with a variation of the buckling element height h1 and fixed buckling element diameter d1 /b = 0.1.

2.5 [−]

(b)

Figure 6: Results from periodic the cell model for a variation of the frame

PT

the influence of the bending strain energy is simply too high.

element diameter d2 for fixed d1 /b = d3 /b = 0.1 and h1 /b = 0.2. (a) Resulting

A variation of the buckling element height h0 (see Figure

average critical stresses as a function of the frame element diameter d2 . (b)

4b) shows a completely different behavior. With increasing

Maximum lateral expansion H22 and maximum rotation ϕz of the connection point of the set of three parallel buckling elements as a function of the frame

tonic unstable. Thereby the force maximum increases and the

element diameter d2 .

CE

height the behavior changes from monotonic stable to non-mono-

AC

minimum decreases. In general, (depending on k3 ) this change

tions of the periodic cell this leads to unwanted deformation

is non-symmetric in the sense that the maximum is increasing

modes, such as shear instabilities.

faster than the minimum is decreasing. This again matches

Similarly increasing the lateral stiffness (k2 reported in Fig-

well with the behavior of the periodic cell (solid lines in Fig-

ure 4c) results in a suppression of the lateral expansion, and the

ure 5c). The critical maximum stress increases faster than the

minimum and maximum critical forces increase and decrease

minimum decreases. Furthermore the minimum approaches the

respectively. This change, however, is symmetric with respect

zero-force level in the periodic cell but does not become nega-

to the limiting case of zero lateral stiffness (k2 = 0). Interest-

tive. One might assume that further increasing the height would

ingly, after reaching a certain lateral stiffness, further increasing

also lead to negative forces; however without further modifica-

has little effect on the extreme values. In the case of the peri7

ACCEPTED MANUSCRIPT

3 2 1 0.5

1

1.5

2

connection diameter

2.5 d3 d1

3

6

ture stays in its buckled configuration, a negative critical force Lat. expansion Rotation ϕmax z

max H22

6 4 4 2 2 0 0.5

with this unit cell with very stiff connection elements d3 = 2d1 and very stiff frame elements d2 = 2d1 as shown in Figure 5b and c (dashed lines). For this combination of parameters there exists a parameter range (d1 /b ∈ [0.04, 0.08] for fixed buckling

0 1

1.5

2

connection diameter

[−]

is necessary. Such a programmable structure can be realized

8

(a)

2.5 d3 d1

3

element height h1 and h1 /b > 0.23 for fixed buckling element

[−]

diameter d1 respectively) where the critical stress minimum is

(b)

negative.

Figure 7: The same as in Figure 6, but with a variation of the connection ele-

CR IP T

4

·10−3

Rotation ϕmax [◦ ] z

8 max lateral expansion H22 [−]

norm. crit. stress

P33 E

[−]

·10−5 Pmax - av. max stress Pmin av. min stress 5

The previous investigation has demonstrated that all fea-

ment diameter d3 for fixes d1 /b = d2 /b = 0.1 and h1 /b = 0.2. (a) Resulting

tures of the simple bar model can be imitated by the investi-

average critical stresses as a function of the connection element diameter d3 .

gated unit cell. Or vice-versa: the single buckling element is

nection point of the set of three parallel buckling elements as a function of the

an appropriate model to describe the most important features of

connection element diameter d3 .

the designed unit cell. Therefore, the results from the simplified

AN US

(b) Maximum lateral expansion (H22 ) and maximum rotation ϕz of the con-

model will be used in the following to derive a semi-analytical

odic cell the same qualitative behavior is observed (Figure 6a).

model for the analysis of the behavior of a finite serial combi-

The minimum and maximum are changed by approximately the

nation of unstable elements.

same amount and for d2 /d1 > 2 further change is negligible.

To better understand the deformation mode within the periodic

4. Effective macroscopic behavior

M

max cell, the maximum lateral expansion H22 and the maximum

rotation of the connecting point of the individual buckling el-

In this section the behavior of a finite number of unstable

are shown in Figure 6b. The lateral expansion is ements ϕmax z

elements is investigated. For this purpose, the results of the simple bar model investigated in the previous section are used

the rotation is increased. Thus the buckling mode has actually

to derive a semi-analytical multi cell model for a serial com-

been changed from an expansion dominated to a rotation dom-

bination of unstable elements. Furthermore an explicit finite

inated deformation mode.

element model is used to confirm these results on a structure

PT

ED

reduced by an increased frame diameter, but at the same time

set-up by the unit cells investigated in the previous periodic cell

The influence of the connection element diameter d3 as a

model. These models are applied to combinations of geometri-

investigated in the following. In Figure 7a it is seen, that in-

cally equal cells and to structures with a gradient.

CE

further possibility to decrease the critical stress minimum is creasing the connection stiffness results in a similar effect as in-

AC

4.1. semi-analytical model for multiple unstable elements

creasing the frame stiffness, namely the critical stress minimum

With the semi-analytical model, the behavior of a fictitious

decreases and the maximum increases. However, the maximum

serial combination of n unstable elements is investigated. For

lateral expansion and maximum rotation now both decrease (see

this, the overall potential of such combination of elements is

Figure 7b). The rotation quickly converges to zero while the

given by

lateral expansion converges much slower. So far, there has been no parameter combination for the pe-

Π=

n X i=1

riodic cell, where the critical force minimum is negative. In

Ei (ui − ui+1 ) − λFext u1 ,

(14)

where ui and ui+1 are the upper and lower displacement of ele-

order to also realize programmable properties, where the struc-

ment i, Ei is the inner potential energy of element i given by 8

ACCEPTED MANUSCRIPT

Equation (1), whereby the lateral displacement v = v (u) is

4.2. Finite element multi cell model

given by the solution of Eq. (3). Furthermore Fext is an external

In addition to the semi-analytical model, the explicit finite

force applied on the upper end of the first element. Equilibrium

element method is applied to analyze finite size samples made

configurations for such a system are given as a solution of   ∂E1   for i = 1   ∂s1 − λF ext    ∂Π   ∂Ei ∂Ei−1 0= = − (15) for i = 2, . . . , n ∂si ∂si−1  ∂ui        − ∂En for i = n + 1,

of a spatial arrangement of unit cells with an internal gradi-

is noth-

To damp the high frequency oscillations initiated by the buck-

ing else but the solution of (3) with u substituted by si . Since

ling events, an artificial damping is introduced by a mass pro-

this is a system with multiple unstable elements, the standard

portional damping matrix Di j = αMi j , where Mi j is the mass

Newton procedure fails to solve this system of equations when

matrix and the coefficient α is chosen such that vibrations initi-

the resulting behavior is non-monotonic in both force and dis-

ated by the buckling decay quite fast, but without significantly

placement. Therefore a path-following procedure [28, 34–36] is

influencing the low frequency deformation modes.

ent. Similar to the periodic cell model, the structure is approximated with beam elements based on the Timoshenko hypothesis. The material is again modeled as linear elastic using the Saint-Venant Kirchhoff model with mass density of 1.18g/cm3 .

∂Ei ∂si

AN US

where the abbreviation si = ui − ui+1 is used and

CR IP T

∂sn

applied, whereby the arc-length ∆s of the resulting equilibrium

4.3. Results

branch is prescribed, instead of prescribing the usual Neumann

Figure 8a-e shows the effective behavior for selected num-

or Dirichlet boundary condition. This is secured by using the

bers of equal elements calculated using the semi-analytical model. The equilibrium path of a single element is given as a reference

(16)

M

so-called Crisfield constraint [35] q f (u, λ) = (u − um )T (u − um ) + (λ − λm ) − ∆s = 0,

(grey dashed lines). The stable (unstable) equilibrium configurations are shown as a solid (dotted) blue lines. Similar to

ternal force factor, and the superscript index m denotes values

the single element case, the unstable configurations are inad-

from the last converged solution. The resulting extended system    ∂Π(u,λ)   ∂u  (17)   = 0 f (u, λ)

missible. However, in contrast to the single element case, a

ED

where u and λ are the current displacement and the current ex-

stabilization by a displacement boundary condition is only par-

PT

tially possible. The buckling under force boundary conditions is shown for the two element case in Figure 8a but does not change

is then solved using the Newton procedure together with a block

with the number of elements. In contrast to this, the buck-

CE

elimination (See [36] for details). By doing this, the second

ling behavior under a controlled overall displacement changes

derivative of the potential Π is needed, which can be calculated

significantly. The corresponding snap through and snap back

analytical by (15). Similar to the simple bar model, the stability

AC

events for loading and unloading under controlled displacement

of an equilibrium configuration is evaluated by the definiteness

are depicted in Figure 8b-e (green and red lines respectively).

of the second derivative. A equilibrium configuration is stable,

The dissipated energy Ediss of one loading-unloading cycle is

if all eigenvalues of the stiffness matrix K are positive, whereas

given by the enclosed area. The normalized energy dissipation

it is unstable, if at least one eigenvalue is negative. The lat-

Ediss /Emax for different numbers of unstable elements is sum-

ter case implies, that there exists a deformation increment ∆u

marized in Figure 8f. In this Emax is the maximum possible

such that the incremental change in the strain energy ∆uT K∆u

energy dissipation for this combination of elements, which can

is negative. However, ∆u must be in the space of admissible

be achieved by force boundary conditions.

deformations, otherwise such a deformation is not possible, and the instability is stabilized by the boundary condition. 9

ACCEPTED MANUSCRIPT

0.5

0.5 0.25

0.2 0.4 0.6 norm. displacement

0.8 [−]

0

1

0

Pu i 2h0

1

single element n = 40 elements, unstable eq n = 40 elements, stable eq. n = 40 elements, buckling disp. b.c. n = 40 elements, unbuckling disp. b.c.

0.75

0.5

0.8 [−]

Pu i 2h0

CE

0.2 0.4 0.6 norm. displacement

(d)

1

0

0

0

0.2 0.4 0.6 norm. displacement

0.8 [−]

1

Pu i 2h0

[−]

(c)

theoretical max.

1

0.8

0.6

0.4

0.2

0.25

PT

0.25

0

0

1

M

1.25

norm. force

0.5

0

0.8 [−]

ED

norm. force

0.75

·10−2

1.5

F1 k1 l0

1

single element n = 12 elements, stable eq. n = 12 elements, unstable eq. n = 12 elements, buckling disp. b.c. n = 12 elements, unbuckling disp. b.c.

F1 k1 l0

[−]

1.25

0.5

(b)

[−]

·10−2

1.5

0.75

Pu i 2h0

[−]

(a)

0.2 0.4 0.6 norm. displacement

Ediss Emax

0

1

0.25

dissipated energy

0

0.75

AN US

0.25

1.25

norm. force

0.75

single element n = 5 elements, stable eq. n = 5 elements, unstable eq. n = 5 elements, buckling disp. b.c. n = 5 elements, unbuckling disp b.c.

F1 k1 l0

1

·10−2

1.5

CR IP T

norm. force

norm. force

1.25

single element n = 3 elements, unstable eq. n = 3 elements, stable eq. n = 3 elements, buckling disp. b.c. n = 3 elements, unbuckling disp. b.c.

F1 k1 l0

1

F1 k1 l0

[−]

1.25

single element n = 2 elements, unstable eq. n = 2 elements, stable eq. n = 2 elements, buckling force b.c. n = 2 elements, unbuckling force b.c.

·10−2

1.5

[−]

·10−2

1.5

0.2 0.4 0.6 norm. displacement

0.8 [−]

Pu i 2h0

(e)

1

0

0

10 20 30 number of elements n [−]

40

(f)

Figure 8: Behavior of a serial combination of n equal unstable elements with parameter k2 /k1 = 0.1, k3 /k1 = 0.01, and h0 /l0 = 0.5 analyzed with the semi-

AC

analytical model. (a-e) Equilibrium path (blue), buckling events for loading (green) and unloading (red), and as a reference the behavior of a single buckling element (dashed grey) for (a) n = 2, (b) n = 3, (c) n = 5, (d), n = 12, and (e) n = 40 elements. (f) Effective energy dissipation ratio within one displacement controlled loading-unloading cycle normalized with the maximum, achievable by a force controlled cycle.

10

ACCEPTED MANUSCRIPT

ity point, one element can suddenly snap while the surrounding

shows a smooth behavior. However, the unstable parts of the

elements are unloading. Following the effective stiffness argu-

equilibrium path are much steeper than in the single element

ment from above, it is obvious that the specific energy dissipa-

case. By reaching a first limit point, exactly one element is un-

tion Ediss must increase as the number of elements increases.

loading on the unstable equilibrium branch, while the remain-

The energy dissipation then approaches the theoretical maxi-

ing structure is unloading on the primary stable branch. Minor

mum (Emax ) that would be obtained by force boundary condi-

(numerical) differences decide which element takes the unsta-

tions (see Figure 8f). Already with 11 elements, 80% of the

ble branch. It is important to note that only one element unloads

maximum possible energy dissipation is achieved.

CR IP T

For n = 2 elements (Figure 8a) the force-displacement curve

on the unstable branch at any time, since this leads to mini-

Figure 9 shows results for an inhomogeneous combination

mization of the (internal) strain energy. Loading and unloading

of eleven elements with a gradient in the parameter k3 /k1 ∈

{0.0, 0.005, . . . , 0.04} . The other parameters (k2 /k1 = 4; h0 /l0 =

with a controlled displacement follows the same path for two

0.5) are chosen such that most elements have a negative force

rations can be reached and stabilized by the applied boundary

minimum and thus stay in the buckled configuration upon un-

condition.

loading. The unstable parts of the equilibrium path of this ele-

AN US

elements. Thus, in this case all unstable equilibrium configu-

For a combination of at least 3 elements in series, the ef-

ment arrangement are plotted in dotted blue and the stable parts

fective behavior under controlled displacement shows a first

are highlighted in green. These 12 stable branches correspond

buckling event at a normalized displacement of approximately

to the following configurations: the first k elements are buckled,

0.25 and 0.8 during loading, and 0.75 and 0.2 during unloading

while 11 − k elements are not buckled. Thus the total number of stable states is n + 1, where n is the number of elements. This

(See Figure 8b). Loading such a structure with a displacement

behavior appears at first sight similar to the above examples.

at these critical points. The unstable configurations between

However, the individual minima and maxima are different, and

these buckling and unbuckling events could only be reached by

more importantly, each stable equilibrium path now crosses the

prescribing more than only the displacement un+1 .

zero force level. Thus, if one unloads such a structure on any

ED

M

boundary condition thus leads to an abrupt strain localization

one of these branches, the structure stays at this zero force con-

calization are getting larger with increasing number of buckled

figuration. Upon subsequent loading, the structure follows this

elements upon loading, and increasing number of unbuckled

specific stable branch. The effective compliance of this branch

elements upon unloading. This is expected since the elements

is given by the sum of the compliances of the individual ele-

that are already buckled are then unloading on the second stable

ments, which in turn depends on whether they are buckled or

branch which is much steeper in the neighborhood of the crit-

not. By loading to a specific preload one can actually pro-

ical force than the first stable branch. Thus, if more elements

gram the effective behavior of such a structure. The effective

are unloading on the second stable branch, the force jumps are

stiffnesses of these stable branches are given in Figure 9b as a

more pronounced.

function of the corresponding stable range ∆¯u on which these

AC

CE

PT

For n > 5 it is clearly seen, that the force jumps at the lo-

effective properties can be realized.

The buckling and unbuckling occur at different total strains upon loading and unloading which leads to a hysteretic behav-

Combining three elements with a gradient in the param-

ior and energy dissipation in a cyclic process. The abrupt lo-

eter h0 ∈ {0.25, 0.435, 0.6} and investigating all equilibrium paths displays a somewhat more complicated picture (see Fig-

calization can be rationalized by the fact that the single element in a larger structure is not under displacement control but

ure 10a). Although only three elements are combined, there

rather under an effective force control. Reaching a first instabil-

are already eight stable equilibrium branches. This is because 11

ACCEPTED MANUSCRIPT

single elements n = 3 elements, unstable eq. n = 3 elements, stable eq. linear fit

0.075

F1 k1 l0

norm. force F¯ =

0.025

0

0.025

CR IP T

[−]

0.05

F1 k1 l0

norm. force F¯ =

0.05

[−]

single elements n = 11 elements, unstable eq. n = 11 elements, stable eq. linear fit

0.075

0

−0.025

−0.025

0

0.2 0.4 0.6 norm. displacement u ¯= (a)

0.8 Pu i 2h0

1 [−]

0.24

n = 11 cells

0.22

0.32

0.315

0.18 0.2 0.22 valid range ∆¯ u [−]

∆F¯1 ∆¯ u

1.2

0

(a)

n = 3 cells n = 5 cells n = 9 cells

0.2

0.18 0.16 0.14 0.12 0.1 0.05

0.24

0.1

0.15 0.2 0.25 valid range ∆¯ u [−]

0.3

(b)

CE

0.16

effective stiffness

M ED

0.325

PT

effective stiffness

∆F¯1 ∆¯ u

[−]

[−]

0.33

0.2 0.4 0.6 0.8 1 norm. displacement u ¯ = P u2hi [−]

AN US

0

−0.05

(b)

Figure 10: The same as in Figure 9 but with a variation of h0 /l0 . (a) Forcedisplacement diagram: Stable equilibrium branches (green) unstable equi-

[0.0, 0.005, . . . , 0.04], k2 /k1 = 4 and h0 /l0 = 0.5. (a) Force-displacement dia-

librium configurations (dotted blue), reference curves of the single elements

gram: stable Equilibrium branches (green), unstable equilibrium configurations

(dashed grey). h0 /l0 is varied in three equal steps from 0.25 to 0.6, other pa-

AC

Figure 9: Behavior of a serial combination of 11 different elements with k3 /k1 ∈

(dotted blue), and reference curves of the single element (dashed grey). (b) The

rameters are k2 /k1 = 0.1 and k3 /k1 = 0.1. (b) The 8 (n = 3), 32 (n = 5), and

eleven programmable stiffness values associated with the stable equilibrium

512 (n = 9) programmable stiffness values for different number of elements n

branches and the corresponding displacement range.

of the stable equilibrium branches and the corresponding displacement ranges. Values for n = 3 (blue) corresponds to results shown in Figure 9a.

12

ACCEPTED MANUSCRIPT

·10−4

the force-displacement curves of the individual elements cross each other. As a consequence, the first element that snaps back again upon unloading is not the element that has snapped down [−] norm. force

down first. Thus, with a combination of loading and unloading

F3 EA0

1

as the last element, but instead the element that has snapped all possible combinations of buckled and unbuckled elements P  can be reached and there are nk=1 nk = 2n possible stable configurations.

degressive gradient (loading) degressive gradient (unloading) progressive gradient (loading) progessive gradient (unloading)

0.5

CR IP T

0

To illustrate how this influences the possible effective stiffness, the linear fit for all stable configurations obtained for a

0

combination of three, five, and nine elements are summarized in Figure 10b, where the ratio h0 /l0 is varied in three, five, and

0.1 0.2 0.3 norm. displacement

u3 l0

0.4 [−]

0.5

Figure 11: Results from an explicit dynamic calculation of a real finite structure with 8×10×11 cells. Investigation in different gradients in the buckling element

tion in h0 /l0 has a large effect on the single elements properties,

height h1 . The corresponding values for of each layer are given in Table 1.

the stiffness can be changed in a much larger range than in the

Fixed values are d1 /b = 0.1 and d2 /b = d3 /b = 0.2.

AN US

nine equal increments from 0.25 to 0.6. Since the chosen varia-

example in Fig. 9. However, within this programmable range,

present multi cell model. Even for this gradient structure a good

much smaller relative increments are possible and in the case of many elements, the possible states are almost continuous.

agreement is obtained. This further justifies the applied periodic cell model, even if it is not rigorous.

To confirm the findings of the semi-analytical model, a struc-

M

ture of a finite combination of cells is analyzed using the finite

element model. The behavior of 8 × 10 × 11 cells loaded in the

ED

z-direction is simulated. Two different gradients, resulting from

Finally, comparing the two types of gradients in the elemen-

tal structures it becomes obvious that one can actually shape the effective macroscopic stress-strain behavior of such a structure in nearly arbitrary ways. The gradients investigated here for

a variation of the parameter h1 from layer to layer are investi-

example lead to a progressively or a degressively increasing ef-

gated. The remaining parameters are set to very stiff values to

fective stress-strain behavior.

PT

prevent a destructive failure of the weaker cells (all parameters are summarized in Table 1 and Figure 11).

CE

Figure 11 shows the resulting normalized force-displacement

5. Discussion Normally the behavior of a material with a regular micro-

compacting of the individual layers, no other unwanted defor-

structure can be fully characterized by the unit cell. In case

mation modes, like for example, shear bands or macroscopic

of a irregular microstructure, a statistically representative char-

buckling occur. A comparison of the structure corresponding

acterization is possible via so called statistically representative

to the loading and unloading paths (Figure 11) shows that un-

volume elements. Changing the number of unit cells has no

buckling occurs in the inverse order than buckling has occurred.

further effect on the effective properties. As shown in this pa-

AC

curve of one complete loading and unloading cycle. Beside the

2

Thus one can indeed reach n stable states by combining several

per, this is fundamentally different for unstable microstructures.

loading and unloading cycles.

Due to multiple equilibrium configurations for the same exter-

In Table 1 the critical stresses

per Pmax

predicted by the peri-

nal force, localization becomes possible and thus the effective

odic cell model (see dashed lines in Figure 5c) are compared

behavior becomes crucially depended on the number of cells

with the critical stresses

Pexpl max

at which buckling occurs in the

and the loading history. For the structure presented in this pa13

ACCEPTED MANUSCRIPT

Layer Gradient 1 Gradient 2 per Pmax /E

Pexpl max /E

1

2

3

4

5

6

7

8

9

10

11

0.16

0.164

0.168

0.172

0.18

0.188

0.196

0.21

0.226

0.24

0.28

0.16

0.2

0.214

0.23

0.244

0.252

0.26

0.268

0.274

0.278

0.28

5

3.4

5.3

6.1

7.1

7.9

8.4

8.9

9.3

9.6

9.8

9.9

× 105

3.5

5.0

5.7

6.5

7.6

8.4

8.9

9.2

9.3

9.4

9.5

× 10

per

Table 1: Values h1 /b for the gradient type of structure investigated with the full finite element model and the critical stresses predicted by the periodic cell Pmax

CR IP T

compared to those predicted by the explicit multi cell finite element model

expl Pmax .

formation from strain energy to the energy of local vibration is

fore been approached by successively increasing the number

initiated and driven by the inner structure, rather than the con-

of cells. It has been found that the behavior under a force con-

stituent material. This has a few interesting implications. First,

trolled boundary condition does not change in general. All cells

if the external loading speed or frequency is significant lower

are loaded by the same force and forces beyond the buckling

than the frequency of the local vibrations, the transformation of

point can only be realized if all cells are buckling simultane-

strain energy to vibrational energy should be independent of the

ously. If one controls the effective deformation instead, there

speed of external loading. Thus, energy dissipation will effec-

are multiple (stable and unstable) ways to realize a certain ef-

tively appear to be instantaneous and to occur under quasi-static

fective deformation. Thus reaching a first instability point ex-

as well as dynamic loading independent of loading velocity.

actly one element is unloading on the unstable branch, while all

This stands in sharp contrast to a viscoelastic material where

other elements are unloading on the stable branch. For increas-

the amount of energy dissipation is proportional to the speed of

M

AN US

per the transition from a structure to a (meta)material has there-

deformation and tends to zero in the case of quasistatic loading.

sis become smaller, the overall behavior becomes smoother and

The effective behavior rather mimics idealized plastic behavior

the effective material behavior asymptotically approaches the

in this respect. However, the presented unit cell can be designed

one found for a single cell under a force boundary conditions

in such a way, that the deformation of the inner structure is fully

(see Figure 8). In the literature dealing with unstable metama-

elastic and the effective deformation is reversible. Thus the ef-

PT

ED

ing numbers of cells, the discontinuities in the effective hystere-

fective behavior is also fundamentally different from plastically

missing, and thus it has not been clear how these structures be-

deforming materials where energy dissipation is directly con-

have in the case of many unstable cells.

nected to irreversible deformation. Both aspects, reversibility

CE

terials this important step from a structure to a material is often

and the independence from loading speed have been demon-

arises how energy can be dissipated in materials that only de-

strated in actual realizations of such structures [19]. The en-

AC

Irrespective of the exact mode of deformation, the question form and buckle elastically. One expects the buckling of one

ergy dissipation characteristic of such a material might be used

single element in a structure of many elements to induce severe

as an alternative concept to protect buildings. In contrast to

local vibrations in all modes belonging to acoustic and optical

other approaches [38, 39], where Rayleigh waves are reflected

branches depending on the details of the geometry of the inner

by a metamaterial, energy could be dissipated over a broad fre-

structure of the elements [37]. These elastic vibrations, in a

quency range and thus protect from different sources. In con-

real material, will be damped by the viscoelasticity of the con-

trast to other engineering applications a high dissipated energy

stituent material and thus the energy will finally be dissipated

ratio per mass is not a leading criterion.

to heat.

Some of these properties remind of on the crushing behav-

Importantly, the first step of this process, namely the trans-

ior of polymeric and metallic foams or of honeycombs used for 14

ACCEPTED MANUSCRIPT

In general, all of the discussed properties can also be found

the behavior is indeed very similar. In both cases the deforma-

in a metamaterial with an internal gradient. However there are

tion path shows a relatively high stiffness in the beginning, then

some more interesting properties that can be achieved by a gra-

a first elastic and or plastic buckling leads to a limitation of the

dient in the inner structure. First of all, such a gradient in the

effective forces and finally a densification of the microstruc-

internal structure can be directly used to shape the effective

ture leads to a stiffening behavior for high strains above 30%

stress-strain curve of loading and unloading in nearly arbitrary

(compare Fig. 11 and [40]). Even a reversible deformation,

ways (see Fig. 11). Since, minimum and maximum critical

and thus a repeatable behavior has been observed for polymeric

stresses of the single element can be influenced by appropriate

foams [46, 47]. However, the ”failure” of these structures lead-

parameters, loading and unloading can even be influenced in-

ing to the plateau in the stress-strain curve consists of a global

dependently from each other. Exploiting this, one can design

shear mode, and a local rotational mode on the cell level. Es-

metamaterials which decelerate (or accelerate) an object in a

pecially the rotational mode leads to highly localized deforma-

arbitrarily prescribed way. This could for example be used for

tions in the postbuckling regime and thus inevitable to a damage

different kinds of crash absorption applications. A possible sce-

of the microstructure if the constitutive material cannot with-

nario would be that one first wants to keep deceleration at a low

stand very high strains [41]. This failure mode is certainly trig-

level avoiding any damage to a human but if the impact is too

gered by the irregularity and randomness of the microstructure.

severe one wants to raise the deceleration accepting damage but

However, high localization and damage is also observed dur-

still avoiding fatal injury.

AN US

CR IP T

crash absorption [40–47]. From a macroscopic point of view

It is not only possible to tailor the effective properties of

ing crushing of regular and periodic microstructures like honeycombs and the different variations thereof (see [42–45]).

such a metamaterial a priori but it is also possible to change them by the loading history. This is achieved by combining

tices has been studied in detail [10, 20, 48], although these

elements with a negative force minimum and thereby giving

structures seems to almost recover from high strains of up to

up full local effective reversibility of the structure. Unload-

50%, there is a significant decrease in the hysteresis in sub-

ing such a metamaterial on one of the stable branches leaves

sequent cycles resulting from damage of the microstructure.

it in the buckled configuration. The metamaterial then shows

The metamaterial investigated here is fundamentally different

the effective stiffness of the branch it is now located on which

in that the effective deformation properties are driven by a con-

can be very different from the stiffness of the initial equilibrium

trolled buckling, while all other (unwanted) deformation modes

branch.

CE

PT

ED

M

Recently the recovery behavior of hollow Nickel microlat-

One can also take advantage of a combination of elements

fully reversible deformations for a broad range of constituent

where the stress-strain behavior curves are crossing each other.

materials and not only rubber like materials, like the ones used

Unloading such a structure, the first element that is unbuck-

for foams [46] or the microstructures presented in [23]. How-

ling, is the one that has buckled first, instead of the one that has

ever, this also directly leads to a fundamental limitation on the

buckled last. Thus, in the sense of programmability, the pro-

presented metamaterial: since only a few buckling elements are

grammable range can be increased from n + 1 configurations

active with respect to the energy dissipation, the amount of dis-

to n2 possible configurations (compare Figure 9b with Figure

sipated energy (per mass) is generally lower compared to col-

10b). It is worth noting here that, although these metamateri-

lapsing metallic foams where energy is dissipated due to the

als are not locally reversible, they may still be fully reversibly

plastic deformation of all elements (See also measured energy

brought back to the initial configuration if an appropriate nega-

dissipation ratios in [19]).

tive force is applied.

AC

are suppressed by the stabilizing basis structure. This enables

15

ACCEPTED MANUSCRIPT

Several approaches to generate similar properties are re-

thermore far away from a realization as a bulk material. Of

ported in the literature [18, 49–57]. Compared to the meta-

course this programmability does not only hold for the effec-

material presented in this paper, there are some crucial disad-

tive stiffness but can be applied to other properties as well. Even

vantages common to all these approaches: 1. The activation of

the complete stress-strain hysteresis can be changed within the

the switching between different properties needs some kind of

same structure. This can be viewed as an artificial realization

external activation [18, 50, 51], 2. most of them are compli-

of transformational or reversible plasticity [55], although, the

cated to produce, especially if one wants to go beyond only

term ”reversible plasticity” is somewhat misleading. 6. Summary

vestigated and a realization as a material, especially as a three dimensional bulk material is questionable and only rarely ad-

CR IP T

one cell [57], and 3. often, only a single mechanism is in-

The effective properties of mechanical metamaterials com-

dressed [49, 53, 54]. Some effort has, for example, been made

posed of unstable building elements is investigated here. A re-

to generate programmable effective properties through so called

cently proposed unstable microstructure [19] is studied in de-

origami based metamaterials. This is realized either by using

tail. For this particular microstructure a rich behavior has been

bi-stabilities within such structures [49], or by using folding

AN US

found including self-recovering unit cells and the possibility to

mechanisms [54]. However, these structures are rather compli-

tailor critical buckling forces upon loading and unloading.

cated to produce. Including possibilities to externally (for ex-

The consequences of this broad range of possible charac-

ample electrically, thermally or even by light [50–52]) activate

teristics of one single unstable element leads to even more di-

the switching process is even more ambitious. A much sim-

verse possible characteristics on the next higher scale which

pler and more promising approach are cellular structures with

has exemplary been investigated by analyzing combinations of

M

periodic arrangement of ellipsoidal and circular holes which undergo a pattern transformation upon loading (see [18] and

ED

references therein). This transformation process can be easily

unstable unit cells. Effective properties such as stress-strain curves, or local properties such as the differential stiffness can be tailored in almost arbitrary ways by combinations of appro-

controlled by an external lateral confinement, leading to a rich

priate unstable cells. However, they can also be changed (pro-

variety of the resulting stress-strain behavior [18]. However, to

grammed) within the same material by making use of the nega-

PT

realize such a confinement, there is still the need of some extra

tive force minimum of the individual cells or by exploiting the

external manipulation. Whereas in the metamaterial presented

multiple stable equilibrium states for the same external force or

here the reaction to an external force can be switched just by

CE

effective deformation. The complete hysteresis, including all

the same external load.

local properties, can then be changed within the same material.

There is yet another way to think about programmability of

The energy dissipation mechanism of such metamaterials

AC

metamaterials without the need of a negative force minimum:

are found to be completely different from usual materials. They

due to the order of the inverse unbuckling behavior there can

rely on converting the externally applied work into vibrational

be multiple stable configurations of buckled and unbuckled el-

modes of the metamaterial structure, effectively showing nei-

ements for the same effective strain range (compare Figure 9a

ther viscoelastic nor plastic energy dissipation.

with 10a). Thus, the metamaterial can show a different behavior, even if the effective deformation is the same. A similar

7. Acknowledgement

behavior has very recently been realized by Hame et. al. [56] with a parallel arrangement of bistable elements. However this

This project is financially supported by the Hector Fellow

mechanism also needs some extra external activation and is fur-

Academy through the project ”mechanical metamaterials”. This 16

ACCEPTED MANUSCRIPT

support is gratefully acknowledged. The authors would like to

[16] Y.C. Wang, M. Ludwigson, and R.S. Lakes, Deformation of extreme viscoelastic metals and composites, Mater. Sci. Eng. A 370 (2004), 41–49.

thank Martin Wegener, Tobias Frenzel, and Matthew Berwind

[17] W.J. Drugan, Composite Materials Having a Negative Stiffness Phase Can

for many suggestions and useful discussions.

be stable, Phys. Rev. Lett. 98 (2007), 055502. [18] B. Florijn, C. Coulais, and M. van Hecke, Programmable Mechanical

References

Metamaterials, Phys. Rev. Lett. 113 (2014), 175503 . [19] T. Frenzel, C. Findeisen, M. Kadic, P. Gumbsch, and M. Wegener, Tai-

[1] G. W. Milton, The Theory of Composites, Cambridge University Press,

lored buckling micro-lattices as reusable light-weight shock absorbers,

Cambridge, UK (2002).

Adv. Mater. 28 (2016), 5865–5870.

[2] M. Kadic, T. B¨uckmann, R. Schittny, P. Gumbsch, and M. Wegener, Pen-

mechanisms in hollow metallic microlattices, J. Mater. Res. 29 (2014),

and Mass Density, Phys. Rev. Appl. 2 (2014), 054007.

1755–1770.

[3] T. B¨uckmann, R. Schittny, M. Thiel, M. Kadic, G. W. Milton, and M.

[21] D. M. Correa, C. C. Seepersad, and M. R. Haberman, Mechanical de-

Wegener, On three-dimensional dilational elastic metamaterials, New J.

sign of negative stiffness honeycomb materials, Integrated Materials and

Phys. 16 (2014), 033032.

Manufacturing Innovation 4 (2015),1–11.

[4] M. Kadic, T. B¨uckmann, R. Schittny, and M. Wegener, Metamaterials

[22] A. Rafsanjani, A. Akbarazadeh and D. Pasini, Snapping Mechanical

beyond electromagnetism, Rep. Prog. Phys. 76 (2013), 126501.

Metamaterials under Tension, Adv. Mater. 27 (2015), 5931–5935.

[5] J. Christensen, M. Kadic, M. Wegener, and O. Kraft, Vibrant times for

AN US

mechanical metamaterials, MRS Commun. 5 (2015), 453–462.

[23] S. Shan, S. H. Kang, J. R. Raney, P. Wang, L. Fang, F. Candido, J. A. Lewis, and K. Bertoldi, Multistable Architected Materials for Trapping

[6] T. B¨uckmann, N. Stenger, M. Kadic, J. Kaschke, A. Fr¨olich, T. Ken-

Elstic Strain Energy, Adv. Mater. 27 (2015), 4296–4301.

nerknecht, C. Eberl, M. Thiel, and M. Wegener, Tailored 3D Mechani-

[24] R.S. Lakes, T. Lee, A. Bersie, and Y.C. Wang, Extreme damping in com-

cal Metamaterials Made by Dip-in Direct-Laser-Writing Optical Lithography, Adv. Mater. 24 (2012), 2710–2714.

posite materials with negative stiffness inclusions, Nature 410 (2001),

565–567.

[7] L.J. Gibson and M.F. Ashby, Cellular Solids: Structure and Properties,

M

2nd edn. Cambridge University Press, Cambridge (1997).

[25] R. Hill, On constitutive macro-variables for heterogeneous solids at finite

[8] M.F. Ashby, A.G. Evans, N.A. Fleck, L.J. Gibson, J.W. Hutchinson, and Oxford (2000).

ED

H.N.G. Wadley, Metal Foams: A Design Guide, Butterworth Heinemann,

of structured elastoviscoplastic solids at finite strains, Mech. Mater. 33

honeycombs and two-dimensional model foams at finite strain, Int. J.

PT

Mech. Sci. 45 (2003), 891–913.

[10] T.A Schaedler, A.J. Jacobsen, A. Torrents, A.E. Sorensen, J.Lian, J.R.

[28] P. Wriggers, Nonlinear Finite Element Methods, Springer-Verlag Berlin

Greer, L.Valdevit, and W.B. Carter, Ultralight Metallic Microlattices, Sci-

Heidelberg (2008)

CE

[29] T. Antonakakis, R.V. Craster, and S. Guenneau, Homogenization for elas-

[11] L. R. Meza, S. Das, and J. R. Greer, Strong, lightweight, and recoverable

tic photonic crystals and dynamic anisotropy, J. Mech. Phys. Solids 71

three-dimensional ceramic nanolattices, Science 345 (2014), 1322–1325.

(2014), 84–96.

[12] T. A. Schaedler, C. J. Ro, A. E. Sorensen, Z. Eckel, S. S. Yang, W. B.

[30] S. Forest, K. Sab, Cosserat overall modeling of heterogeneous materials,

AC

Carter and A. J. Jacobsen, Designing Metallic Microlattices for Energy

Mech. Res. Commun. 25 (1998), 449–454.

Absorber Applications, Adv. Eng. Mater. 16 (2014), 276–283. [13] W. Voigt,

[26] O. van der Sluis O., P.J.G. Schreurs, and H.E.H. Meijer , Homogenisation

[27] J. Hohe and W. Becker, Effective mechanical behavior of hyperelastic

Wadley, The topological design of multifunctional cellular metals, Prog. Mater Sci. 46 (2001), 309–327.

strain, Proc. R. Soc. Lond. A. 326 (1972), 131–147.

(2001), 499-522.

[9] A.G. Evans, J.W. Hutchinson, N.A. Fleck, M.F. Ashby, and H.N.G.

ence 334 (2011), 962–965.

CR IP T

[20] L. Salari-Sharif, T. A. Schaedler, and L. Valdevit, Energy dissipation

tamode Metameterials with Independently Tailored Bulk Modulus and

[31] V. Kouznetsova, M.G.D. Geers, and W.A.M. Brekelmans, Multi-scale

Ueber die Beziehung zwischen den beiden Elas-

constitutive modelling of heterogeneous materials with a gradient-

ticit¨atskonstanten isotroper K¨orper, Ann. Phys. Chem. 274 (1889),

enhanced computational homogenization scheme. Int. J. Numer. Methods

573–587.

Eng. 54 (2002), 1235–1260.

[14] A. Reuss, Berechnung der Fliegrenze von Mischkristallen auf Grund der

[32] R. de Borst, L.J. Sluys, H.B. M¨uhlhaus, and J. Pamin, Fundamental issues

Plastizit¨atsgrenze f¨ur Einkristalle, Z. Ang. Math. Mech. 9 (1929), 49–58.

in finite element analysis of localization of deformation, Eng. Computa-

[15] R.S. Lakes, T. Lee, A. Bersie, and Y.C. Wang, Dramatically stiffer elastic

tions 10 (1993), 99–121

composites materials due to a negative stiffness phase?, J. Mech. Phys.

[33] S. Forest, Nonlinear regularization operators as derived from the micro-

Solids 50 (2002), 979–1009.

morphic approach to gradient elasticity, viscoplasticity and damage, Proc. R. Soc. A 472 (2016), 20150755.

17

ACCEPTED MANUSCRIPT

[34] E. Riks, The application of newtons method to the problem of elastic

9 (2013), 7688–7649.

stability, J. Appl. Mech. 39 (1972), 1060–1066.

[53] M. Schenk and S.D. Guest, Goemetry of Miura-folded metamaterials,

[35] M.A. Crisfield, A fast incremental/iterative solution procedure that han-

Proc. Natl. Acad. Sci. 110 (2013), 3276-3281.

dles snap-through, Comput. Struct. 13 (1980), 55-62.

[54] C. Lv, D. Krishnaraju, G. Konjevod, H. Yu, and H. Jiang, Origami based

[36] K.H. Schweizerhof and P. Wriggers, Consistent linearization for path fol-

Mecanical Metamaterials, Sci. Rep. 4 (2014), 5979.

lowing methods in nonlinear FE analysis, Comput. Methods Appl. Mech.

[55] G. Puglisi and L. Truskinovsky, A mechansim of tranformational plastic-

Eng. 59 (1986), 261–279.

ity, Continuum Mech. Thermodyn. 14 (2002), 437-457.

[37] A.M. Balk, A.J. Cherkaev, and L.I. Slepyan, Dynamics of chains with

[56] R.L. Harne, Z.Wu, and K.W. Wang, Designing and Harnessing the

non-monotone stress-strain relations. I. Model and numerical experi-

Metastable States of a Modular Metastructure for Programmable Mechan-

ments. J. Mech. Phys. Solids 49 (2001), 131-148.

CR IP T

ical Properites Adaption, J. Mech. Des. 138 (2016), 021402.

[38] S. Brˆul´e, E.H. Javelaud, S. Enoch, and S. Guenneau, Experiments on Seis-

[57] N. Bassik, G.M. Stern, and D.H. Gracias, Microassembly based in hands

mic Metamaterials: Molding Surface Waves, Phys. Rev. Lett. 112 (2014),

free origami with bidirectional curvature, Appl. Phys. Lett. 95 (2009),

133901.

091901.

[39] P.G. Martinsson and A.B. Movchan, Vibrations of lattice structures and phononic band gaps, Q. J. Mech. Appl. Math. 56 (2003), 45–64. [40] M.F. Ashby, The properties of foams and lattices, Phil. Trans. R. Soc A 364 (2006), 15-30.

AN US

[41] L.J. Gibson, Mechanical Behaviour of Metallic Foams, Annu. Rev. Mater. Sci. 30 (2000), 191-227.

[42] M.F. Ashby and R.F.M. Medalist, The Mechanical Properties of Cellular Solids, Metall. Trans. A 14 (1983), 1755–1769.

[43] S.D. Papka and S. Kyriakides, In-plane Compressive Response and Crushing of Honeycomb, J. Mech. Phys. Solids 42 (1994), 1499-1532. eycombs, Int. J. Solids Structures Vol. 35 (1998), 239-267.

M

[44] S.D. Papka and S. Kyriakides, In-plane Crushing of Polycarbonate Hon-

[45] S.D. Papka and S. Kyriakides, Experiments and Full-Scale Numerical

ED

Simulations of In-plane crushing of a honeycomb, Acta mater. 46 (1998), 2765-2776.

[46] L. Gong, S. Kyriakides, and W.-Y. Jang, Compressive response of opentures 42 (2005), 1355-1379.

PT

cell foams. Part I: Morphology and elastic properties, Int. J. Solids Struc[47] L. Gong and S. Kyriakides, Compressive response of open-cell foams. Part II: Initiation and evolution of crushing, Int. J. Solids Structures 42

CE

(2005), 1381-1399.

[48] A. Torrents, T.A. Schaedler, A.J. Jacobson, W.B. Carter, and L. Valdevit, Characterization of nickel-based microlattice materials with structural hi-

AC

erarchy from the nanometer to the millimeter scale, Acta Materialia 60 (2012), 3511-3523.

[49] J. L. Silverberg, A. A. Evans, L. McLeod, R.C. Hayward, T. Hull, C.D. Santangelo, and I. Cohen, Using origami design principles to fold reprogrammable mechanical metamaterials, Science 345 (2014), 647–650. [50] Y. Liu, J.K. Boyles, J. Genzer, and M.D. Dickey, Self-folding of polymer sheets using local light absorption, Soft Matter 8 (2012), 1764–1769. [51] J. Ryu, M. D’Amato, X. Cui, K.N. Long, J. Qi, and M.L. Dunn, Photoorigami – Bending and folding polymers with light, Appl. Phys. Lett. 100 (2012), 161908. [52] S. M. Felton, M. T. Tolley, B. Shin, C.D. Onal, E.D. Demaine, D. Rus, and R.J. Wood, Self-folding with shape memory composites, Soft Matter

18