Homogenization of architected materials incorporating shearable beams

Homogenization of architected materials incorporating shearable beams Matteo Franzoi1, Davide Bigoni∗1, and Andrea Piccolroaz1 1Instabilities Lab, University of Trento, Via Mesiano 77, 38123 Trento, Italy December 4, 2025 Abstract Two-dimensional architected materials are often realized as periodic grids of elastic beams. Conventional homogenization methods represent these structures as equivalent elastic solids but neglect shear deformation in the constituent beams. This article addresses this limitation by incorporating shear deformability through Timoshenko beam theory, enabling accurate modeling of stubby beams. Moreover, shearable beams with extreme mechanical characteristics can be obtained through the design of appropriate microstructures. Introducing shearable beams into the grid expands the design space, allowing, for instance, the control of the effective Poisson’s ratio beyond the limits achievable with slender beams. Keywords Architected materials ·Timoshenko beam ·Homogenization ·Auxetic materials 1 Introduction Homogenization of composites is a key tool in the design of architected materials with predetermined mechanical properties. In this context, elasticity is particularly appealing because architected materials often must withstand cyclic loading and vibrations without undergoing permanent deformation. Although homogenization theory is well established within the framework of linear elasticity [1–8], explicit closed-form solutions for equivalent materials are scarce, and several aspects still merit investigation. This article extends the results obtained for two-dimensional periodic grids of elastic beams introduced in [9,10] by incorporating the Timoshenko theory of shear-deformable beams. The interest in the modeling of grids composed of Timoshenko beams arises from three distinct sources. (i) Domokos [11,12], followed by Challamel and co-workers [13–15], and Paradiso et al. [16], have shown that homogenization of a discrete and periodic chain of rigid elements with lumped deformability leads to a continuous rod model that can be represented by either the unshearable Euler model or the Timoshenko model, depending on the mechanical characteristics of the junctions in the chain. In particular, the Timoshenko model is obtained with a chain such as that shown in Fig. 1, equipped with hinges (stiffness kr), sliders (stiffness kt), and axial springs (stiffness ka). In the homogenized continuous model derived from the discrete structure, the bending stiffness ‘EJ’ is given by kra, while the shear stiffness ‘κ GA’ and axial stiffness ‘EA’ are replaced by ktaand kaa, respectively, where ais the length of the unit cell. This result, detailed in Section 2, shows that designing and manufacturing elastic beams in which the three stiffnesses are mutually independent is possible. This opens the possibility of creating slender beams whose shear compliance can range from zero to infinity, potentially dominating the axial and bending responses. In such a configuration, the beam slenderness is given by Λ = lpka/krand may reach very low values, approaching zero. ∗Corresponding author: e-mail: [email protected]; phone: +39 0461 282507. 1 Fig. 1. A microstructured chain model, represented from node i−1to node i+ 2. The chain is constructed by periodically repeating a unit cell of width a; this cell comprises rigid bars connected by axial springs with stiffness ka, elastic hinges with stiffness kr, and elastic sliders (that permit only transverse relative displacement) with stiffness kt. The i-th node undergoes a displacement uiand the overall behavior of the chain can be homogenized, resulting in the Timoshenko beam model, with bending, shear, and axial stiffnesses equal to kra,kta, and kaa, respectively; these stiffnesses remain mutually independent. (ii) Beyond the approach described above, another option exists. A grid of stubby elastic rods can be realized by superimposing beams on different planes so that they interact only through their end cross-sections. This configuration enables levels of stubbiness that cannot be achieved by simply perforating a plate. In particular, the slenderness can approach the theoretical minimum value of 2 (for a hexagonal grid). This arrangement yields a homogenized response with mechanical properties beyond the capabilities of grids of slender beams. This is explained in detail in Section 2and illustrated in Fig. 3. (iii) Finally, for grids of beams lying in the same plane (i.e., obtained by perforating a plate), there is a range of beam stubbiness for which the Euler model provides unrealistic values, while the Timoshenko model maintains a correct behavior for slenderness Λ<12, as shown in Fig. 30. Fig. 2highlights the advantages of the Timoshenko beam model by plotting the effective elastic modulus E(left) and Poisson’s ratio ν(right) for an isotropic solid equivalent to a two-dimensional grid of elastic rods of length land slenderness Λ, arranged in a hexagonal geometry. The rods are either microstructured as in Fig. 1or elastic beams with Young’s modulus Es, Poisson’s ratio νs, and rectangular cross section b×h. While microstructured beams can achieve any slenderness Λ, the method of superposition of elastic beams on different planes cannot permit slenderness Λsmaller than 2 (region highlighted in gray in the figure). The various curves in the figure correspond to the Timoshenko beam model for different values of ka/kt= 2(1 + νs)/κ, whereas the case ka/kt= 0 represents the Euler-Bernoulli model. For the range of slenderness Λshown in the figure, the Euler–Bernoulli theory tends to overestimate the effective elastic modulus and may provide inaccurate values for the effective Poisson’s ratio. The Euler–Bernoulli model also predicts a negative effective Poisson’s ratio, suggesting auxetic behavior, whereas the correct values are positive. Moreover, in the Euler-Bernoulli model, the effective Poisson’s ratio νis unaffected by the constituent stiffness ratio ka/kt(or νs), whereas the Timoshenko theory shows a strong dependency on it. This example reveals a valuable design opportunity. Using microstructured or stubby beams with high shear deformability, the equivalent material can achieve high stiffness and Poisson’s ratio values significantly different from those of slender-beam assemblies. The present article shows that proper design makes grids of microstructured or stubby beams feasible (Section 2) and that employing the Timoshenko model (Section 3) combined with homogenization (Section 4) leads to mechanical characteristics unachievable with assemblies of slender beams (Section 5). These characteristics can involve all elastic parameters of the equivalent solid and highlight the role that the grid constituents play in its overall behavior as a solid. It is noticed in closure that the approach to homogenization developed here is based on the assumption that the characteristic size of the microstructure is much smaller than the specimen size. No attempt is made to develop a higher-order approach, as in [17]. 2 0 2 4 6 8 10 12 14 0.0 0.2 0.4 0.6 0.8 1.0 0 2 4 6 8 10 12 14 -0.4 -0.2 0.0 0.2 0.4 0.6 0.8 1.0 Fig. 2. Effective elastic modulus Eand Poisson’s ratio νcharacterizing an elastic solid of out-of-plane thickness h, equivalent to a two-dimensional grid of beams arranged in a hexagonal lattice (inset). Two models are considered for the constituent beams: (i) a homogenized discrete chain, characterized by the slenderness Λ = lpka/kr, where ka, ktand krare the stiffnesses of the axial spring, slider spring, and elastic hinge, respectively, and lis the total chain length, see Fig. 1; (ii) a continuous elastic beam of length l, cross-section area A=bh, slenderness Λ = 2√3l/b, and shear correction factor κ, made of an isotropic elastic material with Young’s modulus Esand Poisson’s ratio νs(the grid is obtained by superimposing beams on different planes, as explained in Section 2). Different curves represent various ratios ka/kt(EsA/(κ GsA)for the continuous model), with dashed lines indicating the limiting case kt→ ∞ (or κ GsA→ ∞), which corresponds to the Euler–Bernoulli theory. 2 Timoshenko model for grids of microstructured beams and for grids of stubby beams This section elaborates on points (i)–(iii) presented in the introduction. First, it briefly demonstrates how a discrete chain of rigid elements with lumped deformation sources can be converted into a continuous Timoshenko beam model in which axial, shear, and bending deformations are defined by independent parameters. Second, it explains, at least in principle, how to realize a junction between stubby beams working in different planes. Third, it derives the homogenization scheme for grids of Timoshenko shearable beams. 2.1 Continualization of a microstructured chain into a beam model An axially and shear-deformable rod model is obtained from the discrete chain shown in Fig. 1. The elastic energy stored in a discrete chain of nelements with lumped deformability, starting at node 0 and element 1, and ending at node nand element n, is E=1 2 n X i=1 hka(ui−ui−1)2+kt(vi−vi−1−aθi)2i+kr 2 n−1 X i=1 (θi+1 −θi)2,(1) where θiis the rotation of the i–th element, and uiand viare the horizontal and vertical components of the displacement of the i–th node. In the limit as a→0and n→ ∞, the sums become integrals and the differences become derivatives of continuous functions of the axial variable x, yielding E=a 2Zl 0hkau′2+kt(v′−θ)2+krθ′2idx. (2) Eq. (2) models the mechanics of a Timoshenko beam in which the bending stiffness EsJ, shear stiffness κ GsA, and axial stiffness EsAcorrespond to kra,kta, and kaa, respectively. Importantly, the microstructured chain enables a beam whose bending, shear, and axial stiffnesses can be adjusted independently over the entire range from zero to infinity. 3 2.2 A design strategy for grids with stubby beams Two-dimensional grids of beams can be realized by perforating a plate; however, this is not the only possible approach. In particular, Fig. 3shows a proof-of-concept realization of a node connecting four stubby beams. The key idea for the practical realization of the junction is to use thin elements operating in different planes, but preserving the out-of-plane symmetry. Fig. 3. Concept model of a nodal connection where four stubby beams are jointed. The connection is designed for two orthogonal rods as a multi-layer rigid joint (colored in purple), where the out-of-plane symmetry is preserved: (left) exploded view; (right) perspective view. The junction imposes full continuity of displacement between connected elements, which are reinforced with stiff coatings and linked by orthogonal ‘pivots’ with rectangular cross-section to prevent relative rotation, in addition to relative displacement. In this arrangement, in-plane axial and shear forces, as well as bending moments, are fully transmitted through the node, while out-ofplane forces and moments do not develop, provided that the load is symmetric in the plane. Such junctions must be rigid (made of a material stiffer than the beams) and allow the construction of grids of elastic stubby beams that are mechanically equivalent to a uniform material under plane stress conditions. Special junctions connecting elements located in different planes are required because, when a grid is cut directly from a plate, the spaces between beams shrink as the beams become stubby, invalidating the beam model; in this case, a model of a plate with holes becomes more appropriate, as shown in the insets of Fig. 30. Figure 3illustrates a simple case in which the beams converging at the node have a rectangular cross-section, for which the shear correction factor is 5/6≈0.833. Under these conditions, there exists a lower limit for the slenderness, corresponding to the situation where parallel coplanar beams touch each other, thus occupying all the space between two consecutive lines. For a square grid, the rectangular cross-section assumption yields a minimum slenderness Λof 2√3≈3.46 [see Eq. (19)]. The minimum slenderness is equal to 6 for an equilateral triangular lattice and 2 for a hexagonal grid. Smaller values can be achieved by adopting different shapes of the crosssection. For example, a double-T section allows a slenderness of 2 for the square lattice; however, in this case, the shear correction factor drops to approximately 0.296. For the sake of simplicity, this article considers only continuous elastic beams with rectangular cross-sections, although using microstructured beams would enable achieving arbitrarily small slenderness values. Therefore, results are presented starting from zero slenderness to provide a more comprehensive view of the situation. 2.3 A validation of the homogenization scheme for grids of Timoshenko beams Although it is well known that the Timoshenko beam model can effectively predict the quasistatic and dynamic responses of stubby beams [18–22], a validation of the homogenization scheme proposed in the present article, specifically through the use of the Timoshenko model, has not been previously presented. Therefore, although not a primary focus of this study, finite element computations are included in Appendix A, referring to an equilateral triangular grid of beams. By considering a finite yet sufficiently large grid of elastic rods –such that end effects can be neglected– the equivalent solid obtained through homogenization is shown to accurately capture the overall 4 behavior of the ensemble. At the same time, the limitations of the approach when a perforated plate is modeled as a beam network are also highlighted. 3 The Timoshenko beam model The Timoshenko beam theory relaxes the Euler-Bernoulli assumption that the cross sections of a deflected beam remain perpendicular to its axis. In particular, for a rectilinear beam belonging to the x–axis (parallel and orthogonal to the unit vectors e1and e2, respectively) in its reference configuration, the position vector of the displaced points is r(x) = x+u(x)e1+v(x)e2,(3) where u(x)and v(x)are the axial and transverse displacement components. The tangent to the deformed axis is given by the derivative of rwith respect to x, r′(x) = 1 + u′(x)e1+v′(x)e2,(4) so that the stretch of the axis is λ(x) = |r′(x)|. The cross-section of the beam rotates at an angle γwith respect to the tangent r′, while the normal is inclined at the angle θ, as sketched in Fig. 4. Introducing the unit vectors aand b, respectively normal and parallel to the cross section, a= cos θe1+ sin θe2,b=−sin θe1+ cos θe2,(5) the tangent r′can be understood to condensate the effects of the axial ε=λ−1and shear γ strains through the equation e1+u′= (1 + ε)(cos γa+ sin γb),(6) which upon linearization in εand γbecomes e1+u′= (1 + ε)a+γb,(7) where u(x) = u(x)e1+v(x)e2is the vector collecting the axial and transverse displacements. Under the assumption of small deflection and small rotations, Eq. (7) leads to the linearized equations governing the kinematics of a planar Timoshenko beam ε(x) = u′(x), γ(x) = v′(x)−θ(x), χ(x) = θ′(x),(8) where the last equation defines the curvature χof the beam axis associated with the rotation of the cross-section θ. The linearized equilibrium equations in the absence of distributed forces are N′= 0, V ′= 0, M′=−V, (9) where N,V, and Mare the normal and shear forces and the bending moment. In addition, linear constitutive equations are assumed N=EsA ε, V =κ GsA γ, M =EsJ χ, (10) where Esand Gsare the elastic and shear moduli, Aand Jthe cross-section area and second moment of area, while κis the shear correction factor, introduced by Timoshenko [23,24] to account for the non-uniform distribution of the shear stress across the cross-section. For a rectangular crosssection κ= 5/6; however, this value is underestimated, and more accurate relationships have been provided [25–28]. In particular, Cowper [25] obtained an analytical relation for the shear correction factor as a function of the Poisson’s ratio νsof the material, κ(νs), by integrating the equations of three-dimensional elasticity, which is suitable for any uniform beam. Assuming, for simplicity, that all the elastic coefficients are constant and combining equations (8)–(10), the equilibrium equations in terms of the kinematic variables become u′′ = 0, v′′ =θ′, θ′′ =κ GsA EsJ(θ−v′),(11) 5 Fig. 4. Deformed configuration of the planar Timoshenko beam at coordinate x. The material point x(x) = xe1 is mapped to its deformed position r(x). Vectors aand bare, respectively, orthogonal and tangent to the crosssection at r(x), as defined in Eq. (5). The figure also illustrates the slope of the beam axis, v′(x), the rotation of the cross-section, θ(x), and the shear deformation, γ(x) = v′(x)−θ(x). and the elastic strain energy of a Timoshenko beam of length lis defined as E=1 2Zl 0EsA u′2+κ GsA(v′−θ)2+EsJ θ′2dx . (12) A comparison between Eqs. (12) and (2) reveals that a microstructured rod behaves as a Timoshenko beam, but with EsA,κ GsA, and EsJreplaced by kaa,kta, and kra, respectively. It is important to note that these stiffness components can vary independently due to the microstructural assumptions. The solution to the system (11) can be sought by assuming exponential functions for axial and transverse displacements u(x) = 2 X j=1 Cjeiαjx, v(x) = 4 X k=1 Dkeiβkx,(13) where i=√−1is the imaginary unit, {C1, C2, D1, D2, D3, D4}are six arbitrary constants, and {α1, α2, β1, β2, β3, β4}are the six characteristic roots. In the present case, all characteristic roots vanish. As a result, the axial and transverse displacements reduce to a linear and a cubic polynomial, respectively: u(x) = C1+C2x, v(x) = D1+D2x+D3x2+D4x3,(14) where x∈[0, l]. In addition, the rotation θcan be obtained by imposing θ′′ =v′′′ and solving the equation (11)3for θ, leading to the general solution for the cross-section rotation θ(x) = D2+ 2D3x+ 3D4x2+ 2 EsJ κ GsA.(15) The six above constants can be determined by imposing the boundary conditions at the ends of the beam u(0) = u0, v(0) = v0, θ(0) = θ0, u(l) = ul, v(l) = vl, θ(l) = θl,(16) 6 which will be collected in the vector q= [u0, v0, θ0, ul, vl, θl]T. Consequently, the displacement field u= [u(x), v(x), θ(x)]Tcan be written as a function of the vector of nodal degrees of freedom qas u(x) = N(x)q,(17) where N(x) = e1⊗Nu(x)+e2⊗Nv(x)+e3⊗Nθ(x)is a 3×6matrix determining the exact shape functions of a Timoshenko beam, while e3=e1×e2is the unit vector denoting the out-of-plane axis. The three rows of the matrix Ndefine the shape functions for the axial and transverse displacement and the rotation of the cross-section, respectively, Nu(x) = h1−x l0 0 x l0 0i, Nv(x) = 12x1−x l 1 + 12α 01 + 12α 12x+l−2x 12l2 l(1 + 6α)−x 12l0 α+x 4l−x2 6l2 l−x−6α l +x 12l , Nθ(x) = 61−x l 1 + 12α"0−x l2 1 6+ 2α−x 2l0x l2 x12α+ 3x l−2 6(l−x)#, (18) where the following dimensionless parameters are introduced α=EsJ κ GsA l2=2(1 + νs) κΛ2,Λ = lrA J,(19) in which Λis the slenderness of the beam. The elastic energy, Eq. (12), can be expressed in a matrix form by employing Eq. (17) as E=1 2qTKbq,(20) where the stiffness matrix Kbfor a Timoshenko beam of length lcan be expressed as Kb=Zl 0 B(x)TE B(x)dx, (21) in which B(x)denotes the strain-displacement matrix, while Edenotes the matrix containing the axial, shear, and bending stiffnesses B(x) =   d dx 0 0 0d dx −1 0 0 d dx  N(x),E=  EsA0 0 0κ GsA0 0 0 EsJ .(22) The stiffness matrix can explicitly be written as Kb=α κGsA l(1 + 12α)                1 + 12α αϕ20 0 −1 + 12α αϕ20 0 ·12 6l0−12 6l · · 4l2(1 + 3α) 0 −6l2l2(1 −6α) · · · 1 + 12α αϕ20 0 · · · · 12 −6l · · · · · 4l2(1 + 3α)                ,(23) where ϕis defined as ϕ=rEsA κ GsA=r2(1 + νs) κ,(24) 7 and denotes the shear coefficient, a dimensionless parameter that quantifies the ratio between the axial and shear stiffnesses of the Timoshenko beam. For the homogenized version of the discrete chain shown in Fig. 1, the slenderness and shear coefficient can be equivalently expressed in terms of the axial, shear, and rotational spring stiffnesses as follows Λ = lrka kr , ϕ =rka kt ,(25) where lis the length of the equivalent beam. 4 Homogenization of grids of Timoshenko beams The static homogenization procedure matches the strain energies to derive an elastic continuum equivalent to a periodic grid of Timoshenko beams, following an approach similar to that applied to Euler-Bernoulli rods [29]. To this end, a unit cell composed of Timoshenko beams is introduced, as shown in Fig. 5. This cell tessellates a two-dimensional periodic lattice along the directions defined by the direct basis vectors a1and a2. Although illustrated as a square for simplicity, the unit cell can have an arbitrary shape. The unit cell shown in Fig. 5for the square lattice is chosen to minimize the number of degrees of freedom in the system. Fig. 5. Undeformed (left) and deformed (right) configurations for a square periodic lattice and its unit cell, made up of Timoshenko beams, represented thick to emphasize the effects of shear deformation. In the right panel, the gray beams represent the undeformed configuration for reference. Elastic strain-energy density of the unit cell. On introducing the vectors qand f, collecting the nodal displacements and the dual nodal forces, respectively, the equilibrium of the unit cell can be expressed as K q =f,(26) where Kis the symmetric stiffness matrix of the unit cell, obtained by assembling the stiffness matrices Kb, Eq. (23). Periodicity of the lattice requires for each pair of nodes {p, q}such that xq− xp=n1a1+n2a2,n1, n2∈ {0,1}that the components of displacement ˜ q(j)=˜u(j)˜v(j)˜ θ(j)T satisfy ˜u(q)= ˜u(p),˜v(q)= ˜v(p),˜ θ(q)=˜ θ(p).(27) Multiple-scale asymptotic expansion of displacement field for periodic composites [30,31] suggests that, in addition to a periodic field (generating a vanishing macroscopic strain), an affine deformation ε x(j)is considered [4], so that the total displacement field can be written as q(j)=˜u(j)+ (ε x(j))·e1˜v(j)+ (ε x(j))·e2˜ θ(j)T,(28) where εis a given uniform strain tensor (representing a homogeneous strain for the unit cell) and x(j)is the vector collecting the coordinates of the j-th node. With reference to the unit cell reported in Fig. 5, vectors qand fcan be partitioned as q=q(lb)q(rb)q(lt)T,f=f(lb)f(rb)f(lt)T.(29) 8 Periodicity implies that the vectors ˜ q(rb),˜ q(lt)can be expressed through the ‘left-bottom’ term ˜ q(lb)as ˜ q=  ˜ q(lb) ˜ q(rb) ˜ q(lt) =  I I I ˜ q(lb)=Z0˜ q∗.(30) Consequently, the generalized displacement vector qcan be written as the sum of the periodic component, ˜ q, and the homogeneous component, ˆ q(ε), as follows q(˜ q∗,ε) = Z0˜ q∗+ˆ q(ε),(31) where the homogeneous displacement vector of the j-th node is written as ˆ q(j)(ε) =   (ε x(j))·e1 (ε x(j))·e2 0 .(32) Finally, Eq. (31) can be substituted into the equilibrium equation (26) for the unit cell, so that a pre-multiplication by ZT 0leads to the following expression ZT 0K Z0˜ q∗+ZT 0Kˆ q(ε) = ZT 0f.(33) The right-hand side term can be represented as the sum of three nodal force vectors ZT 0f=f(lb)+f(rb)+f(lt),(34) which vanishes in the absence of external loads, as the internal forces shown in Fig. 5satisfy the equilibrium condition, f(lb)=−f(rb)−f(lt). Consequently, Eq. (33) reduces to ZT 0K Z0˜ q∗=−ZT 0Kˆ q(ε).(35) The stiffness matrix Kis rank-deficient because the unit cell is unconstrained, allowing three independent rigid-body motions. Applying the periodicity conditions, represented by Z0, suppresses the rigid-body rotation while still permitting two rigid-body translations. As a result, for any vector trepresenting a pure rigid-body translation, the following relation holds K Z0t=0,(36) which shows that dimker(ZT 0K Z0)≥2. The nullspace of ZT 0K Z0therefore contains at least two zero-energy modes associated with rigid-body translations. Any additional zero-energy mode in ker (ZT 0K Z0)is referred to as an internal mechanism, or ‘floppy mode’. Floppy modes are excluded in the present treatment. The Fredholm alternative theorem states that the system (35) admits a solution if and only if the right-hand side, −ZT 0Kˆ q(ε), lies in the orthogonal complement of the nullspace of the coefficient matrix ZT 0K Z0, that is −ZT 0Kˆ q(ε)∈ker (ZT 0K Z0)⊥.(37) This condition is always satisfied, since any vector tin ker (ZT 0K Z0)represents a rigid-body translation. As a result, t·ZT 0Kˆ q(ε) = K Z0t·ˆ q(ε)=0,(38) as implied by Eq. (36). The solution of system (35) can therefore obtained as a linear function of ε ˜ q∗=˜ q∗(ε).(39) Consequently, the elastic strain-energy density of the unit cell takes the form of a quadratic function in εexpressed as E(ε) = 1 2|C| q(˜ q∗(ε),ε)·K q(˜ q∗(ε),ε),(40) where |C| denotes the area of the unit cell. 9 Fig. 11. Hexagonal grid (left) generated by tessellating a unit cell composed of five elastic Timoshenko beams (right). Assuming that all beams have identical cross-sectional area A, moment of inertia J, and slenderness Λ, the elasticity tensor of the equivalent continuum can be expressed in Voigt contracted form as   ˜ E1111 ˜ E1122 ˜ E1112 ·˜ E2222 ˜ E2212 · · ˜ E1212 =Λ2+ 12 (ϕ2+ 3) 2√3 (Λ2+ 12 (ϕ2+ 1))    1Λ2+12 (ϕ2−1) Λ2+12 (ϕ2+3) 0 ·1 0 · · 24 Λ2+12 (ϕ2+3)   ,(59) where only two of the three components, made dimensionless through Eq. (48), are independent, in accordance with the isotropy conditions (50). The dimensionless Young’s modulus and Poisson’s ratio are obtained by inverting the elasticity tensor, and are given by ˜ E=16√3 Λ2+ 12(ϕ2+ 3), ν =Λ2+ 12(ϕ2−1) Λ2+ 12(ϕ2+ 3),(60) whereas the shear modulus is ˜ G=˜ E/(2(1+ν)), as expected for an isotropic linear elastic material. Assuming the beams are made of an isotropic material and have a rectangular cross-section, so that κ= 5/6, the Young’s modulus ˜ Eand Poisson’s ratio νcan be expressed as functions of the slenderness Λand the Poisson’s ratio of the base material νs, by using relation (24). Under these assumptions, the engineering constants are shown in Fig. 12 as functions of νs(left) and Λ(right). -1-0.75 -0.5 -0.25 0 0.25 0.5 -0.5 - 0.25 0 0.25 0.5 0.75 1 1.1 1.55 1.90 0 10 20 30 40 50 Fig. 12. Effective Young’s modulus, ˜ E, and Poisson’s ratio, νplotted as functions of the Poisson’s ratio νsof the material composing the beams (left), and the slenderness Λ(right), for κ= 5/6. The engineering constants of the effective material can be tuned by properly selecting Λand νs. It is worth noting that the engineering constants of the material equivalent to the hexagonal grid, given in Eq. (60), can be tuned by appropriately selecting Λand νs. In particular, the 16 equivalent Young’s modulus, ˜ E, increases as both the slenderness and the Poisson’s ratio of the beams decrease. Conversely, the equivalent Poisson’s ratio νdecreases as Λand νsdecrease, reaching its minimum value for νs=−1, which also corresponds to the Euler-Bernoulli beam model, and in the limit Λ→0. The fact that bending moments are essential for the equilibrium of the hexagonal grid under analysis (in contrast to the triangular grid discussed in the previous section) is reflected in the observation that smaller values of Λare required, compared to the triangular lattice, to obtain a sufficiently stiff equivalent continuum. In particular, for Λ>50, the effective Young’s modulus becomes negligible, rendering the equivalent continuum extremely compliant. Simultaneously, the effective Poisson’s ratio approaches its upper limit of 1for both the Euler-Bernoulli and Timoshenko models, causing the effective shear modulus ˜ Gto vanish. In fact, for the grid design shown in Fig. 11, bending stiffness is essential. The beams cannot be connected through hinges, as such a configuration would render the grid incapable of sustaining any stress state other than isotropic. This behavior can be demonstrated by considering the limit in which the bending stiffness of the beams vanishes and observing how the equivalent elastic properties change. Since ˜ Eand ˜ Gdepend on the slenderness, and according to Eq. (19), reducing the bending stiffness of the beams corresponds to increasing their slenderness. Consequently, in the limit Λ→ ∞ one obtains ˜ Ehex ∞= lim Λ→∞ ˜ Ehex = 0,˜ Ghex ∞=˜ Ehex ∞ 2(1 + νhex ∞)= 0, νhex ∞= lim Λ→∞ νhex = 1,(61) that is, both the effective Young’s and shear moduli vanish, allowing rigid-body displacements, while the Poisson’s ratio reaches its 2D upper limit of 1, as shown in Fig. 12 (right) as slenderness approaches 50. Conversely, it can be noted that the Young’s and shear moduli for the triangular grid shown in Fig. 8do not vanish as the bending stiffness approaches zero, but instead converge to the finite values ˜ Etri ∞=2 √3,˜ Gtri ∞=√3 4, νtri ∞=1 3,(62) so that the effective Poisson’s ratio attains the 2D plane stress value corresponding to the pinjointed triangular lattice [4,41]. The response of the hexagonal lattice under uniaxial horizontal stress is shown in Fig. 13 for stubby beams, illustrating the difference between the Euler-Bernoulli beam model in panel (a) and the Timoshenko beam model with νs= 0.3in panel (b). As expected, the Euler-Bernoulli beam model exhibits false auxeticity, causing the lattice to expand along direction e2when stretched along direction e1. In contrast, the Timoshenko beam model correctly predicts contraction along direction e2, resulting in an effective Poisson’s ratio greater than that of the base material, ν= 0.353 > νs. (a) Euler-Bernoulli beam model with Λ=1, yielding ν=−0.297 (b) Timoshenko beam model with Λ=1and νs= 0.3, yielding ν= 0.353 Fig. 13. As in Fig. 10, except that the lattice is now hexagonal. The Timoshenko model (b) corrects the false auxeticity exhibited by the Euler-Bernoulli model (a). 17 5.2.3 Triangular lattice with hexagonal rigid inclusions The mechanical response of a triangular lattice can be enhanced by embedding rigid hexagonal inclusions, which constrain the deformation of the surrounding beams and increase the overall stiffness, as illustrated in Fig. 14. Fig. 14. The unit cell composed of six elastic Timoshenko beams connected to a rigid hexagonal element (right) yields an equilateral triangular lattice, stiffened with rigid inclusions (left). The unit cell has side length land consists of six identical elastic Timoshenko beams emanating from the vertices of a rigid hexagonal inclusion, as shown in Fig. 14 on the right. The length and slenderness of the beams are given by lb=l(1 −ζ),Λ = lbrA J=l(1 −ζ)rA J,(63) where ζ∈(0,1) is a dimensionless parameter representing the ratio between the inclusion side length and the unit cell side length. The beams are rigidly connected to the hexagonal inclusion, which imposes kinematic constraints at the junctions. The effective elasticity tensor ˜ E, whose components have been made dimensionless according to Eq. (48), is computed as   ˜ E1111 ˜ E1122 ˜ E1112 ·˜ E2222 ˜ E2212 · · ˜ E1212 =3√3Λ2+ 12ϕ2+ 4 4(1 −ζ) (Λ2+ 12ϕ2)   1Λ2+12(ϕ2−1) 3(Λ2+12ϕ2+4) 0 ·1 0 · · Λ2+12(ϕ2+1) 3(Λ2+12ϕ2+4)   ,(64) which corresponds to an isotropic solid, since ˜ E1111 =˜ E2222,˜ E1112 =˜ E2212 = 0, and ˜ E1212 = (˜ E1111 −˜ E1122)/2. The dimensionless compliance tensor ˜ E−1is fully characterized by two independent engineering constants, expressed as functions of the parameters ζ,Λ, and ϕ, given by ˜ E=2Λ2+ 12(ϕ2+ 1) √3(1 −ζ) (Λ2+ 12ϕ2+ 4), ν =Λ2+ 12(ϕ2−1) 3 (Λ2+ 12ϕ2+ 4),(65) while the effective shear modulus is ˜ G=˜ E/(2(1 + ν)). Note that the effective Poisson’s ratio νis independent of the rigid inclusion size, measured by ζ. Consequently, similarly to the previously analyzed triangular lattice, the Poisson’s ratio of the equivalent material lies within the range −1< ν < 1/3. Compared to the triangular lattice response shown in Fig. 10, the presence of rigid inclusions results in a stiffer response under uniaxial horizontal uniform stress, as illustrated in Fig. 15. While the inclusions reduce the overall displacement magnitudes, the ratio between axial and transverse deformations remains unchanged for the two grids, as evidenced by a comparison of Figs. 10 and 15. 5.3 Re-entrant lattice and true auxeticity Two examples of re-entrant grids composed of Timoshenko beams are examined below, as these geometries are well-known to produce auxetic response in the equivalent material [36,38,39,42– 44]. In these lattices, flexural stiffness plays a crucial role, similarly to the hexagonal geometry 18 (a) Euler-Bernoulli beam model with Λ=1, yielding ν=−11/15 ≈ −0.733 (b) Timoshenko beam model with Λ=1and νs= 0.3, yielding ν= 0.207 Fig. 15. As in Fig. 10, but for a hexagonal lattice with rigid inclusions (ζ= 1/2). The Timoshenko model (b) corrects the false auxetic behavior shown by the Euler-Bernoulli model (a). discussed in Section 5.2.2, so that taking the limit Λ→ ∞ results in a collapse mechanism for the re-entrant lattice, corresponding to a singular constitutive tensor for the effective material. To maintain generality, the shear stiffness will be incorporated through the parameter ϕ, as defined in Eqs. (24) and (25). This parameter ranges from zero (corresponding to the EulerBernoulli beam model or to the Timoshenko model with νs=−1) to infinity (representing a beam with extremely high shear deformability, such as the microstructured beam illustrated in Fig. 1). Accordingly, the Poisson’s ratio of the beam constituent material νswill be replaced by the shear coefficient ϕ. The particular case of beams made from an isotropic material with νs= 0.3and a rectangular cross-section (κ= 5/6) corresponds to ϕ= 1.76. The first example concerns a re-entrant hexagonal unit cell, shown in Fig. 16, consisting of two vertically aligned, ‘arrow-shaped’ beams connected by two additional vertical beams along the lateral edges of the cell. The geometry of the unit cell is governed by the re-entrant angle α, which is restricted to the interval (π/3, π/2). At the upper bound α=π/2, the internal cavity assumes a rectangular shape. For values of α > π/2, the configuration ceases to be re-entrant, and the auxetic behavior is lost. At the lower bound α=π/3, the arrow tips make contact, rendering further reduction of αgeometrically infeasible due to beam interference—unless the lattice joints are designed to allow the beams to operate on separate planes. Fig. 16. A lattice (left) composed of re-entrant unit cells (right), each consisting of eight Timoshenko beams. All beams are assumed to be identical, with length l, cross-sectional area A, and slenderness Λ. Although the side lengths of the unit cell depend on land the re-entrant angle α, the unit cell remains rectangular for all values of (l, α), and the two vectors a1,a2, defining the direct basis, remain mutually orthogonal. The lattice is generated by periodically repeating the unit cell along the directions a1,a2, resulting in the geometry shown in Fig. 16. 19 According to Eq. (48), the dimensionless elasticity tensor can be expressed in Voigt form as   ˜ E1111 ˜ E1122 ˜ E1112 ·˜ E2222 ˜ E2212 · · ˜ E1212 = =k0sin α 2         (Θ −12)(1 −cos 2α) + 72 sin αtan α 2−2(Θ −12) cos α0 ·(Θ −12)(1 + cos 2α) + 24 sin αcot α 2 0 · · 8k1 k0         ,(66) where the coefficients k0and k1are defined as k0=1 Θ(2 + cos 2α) + 12(1 −cos 2α), k1=1 Θ(1 + cos α) + 4(1 −cos α),(67) both depending on the parameter Θ, that combines Λand ϕas Θ = Λ2+ 12ϕ2.(68) The analytical expressions for the dimensionless engineering constants, corresponding to the re-entrant hexagonal lattice, are ˜ E1=12 cot α 2 (Θ + 12 tan2α) cos2α,˜ E2=12 tan α 2 (Θ + 12 cot2α) sin2α+ 24,˜ G12 = 4k1sin α, ν12 =−Θ−12 Θ + 12 tan2α(1 + sec α), ν21 =−2Θ−12 (Θ + 12 cot2α) sin2α+ 24 cos αsin2α 2. (69) Tensor (66) satisfies the condition E1112 =E2212 = 0, implying that the effective material exhibits orthotropic symmetry. This symmetry class is characterized by four independent elasticity components, since the two Poisson’s ratios, ν12 and ν21, are related through the symmetry of the compliance tensor E−1 ν12 E1 =ν21 E2 .(70) 0 10 20 30 40 50 0 0.25 0.5 0.75 1 1.25 1.5 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0 10 20 30 40 50 -5 -4 -3 -2 -1 0 1 -0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 Fig. 17. Effective Young’s and shear moduli (left), and Poisson’s ratios (right) plotted as functions of the slenderness Λ, for a fixed re-entrant angle α= 75◦. Three values of shear parameter ϕare considered: Euler-Bernoulli theory (red lines), Timoshenko theory for an isotropic material with a rectangular cross-section (blue lines), and Timoshenko theory for a highly shear-deformable material (yellow lines). Note that there are two different scales for the elastic properties on the left and the right of the plots. 20 The engineering constants are plotted in Fig. 17 as functions of the slenderness Λ, for both the Euler-Bernoulli beam theory (red lines) and the Timoshenko beam theory, with a fixed re-entrant angle α= 75◦. The shear deformability of the beams reduces the stiffness of the effective material compared to the Euler-Bernoulli beam model, as shown in Fig. 17 (left), with the effect becoming more pronounced as beams become stubby. This reduction is moderate for ˜ E1(with values ˜ E1= 1.36 and ˜ E1= 1.12 for ϕ= 0 and ϕ= 1.76, respectively), but more significant for ˜ E2(decreasing from ˜ E2= 0.32 to ˜ E2= 0.15) and the shear modulus ˜ G12 (decreasing from ˜ G12 = 0.48 to ˜ G12 = 0.07), for Λ=2,α= 75◦. However, the influence of shear deformability diminishes rapidly as the beams become slender, with all curves converging to the same limiting values ( ˜ E1= 0.09,˜ E2= 0.004, and ˜ G12 = 0.001 at Λ = 50 and α= 75◦) and vanishing in the limit Λ→ ∞. In other words, for the re-entrant cell geometry, compared to the previously analyzed grids, shear deformability leads to a significant reduction in the elastic properties of the equivalent material when the beams are stubby. It is worth noting that, for the chosen re-entrant angle α= 75◦, auxetic behavior is observed for all values of ϕ, except in the range Λ<2√3, where the Euler-Bernoulli model predicts positive Poisson’s ratios. However, when Λ≲4, geometric interference occurs between the beams. In this regime, the modeling approaches described in Section 2must be adopted to avoid physical interpenetration and to preserve model consistency. The mechanical response of the re-entrant lattice under horizontal uniaxial stress is shown in Fig. 18. The beam slenderness is Λ=1, a low value chosen to emphasize the effect. The lattice on the left, constructed with Euler-Bernoulli beams, exhibits an incorrect non-auxetic response, while the lattice on the right, modeled with Timoshenko beams, correctly displays auxetic behavior. (a) Euler-Bernoulli beam model with Λ=1,α= 75◦, yielding ν12 = 0.318 (b) Timoshenko beam model with Λ=1,ϕ= 1.76, α= 75◦, yielding ν12 =−0.62 Fig. 18. As in Fig. 10, except that the lattice is now re-entrant (α= 75◦), the Timoshenko model (b) correctly exhibits auxetic behavior, in contrast to the Euler-Bernoulli model (a). The second example considers a re-entrant grid of Timoshenko beams incorporating square rigid inclusions. The unit cell, illustrated in Fig. 19 (right), consists of four rigid inclusions connected at their corners by eight Timoshenko beams. It is important to note that the beams shown in red and those shown in blue in the figure are positioned on separate planes and connected only through the rigid inclusions. The symmetry of the unit cell ensures that the equivalent continuum behaves as a cubic material. As a result, the components of the elasticity tensor satisfy the following conditions E1111 =E2222,E1112 =E2212 = 0.(71) The length of the elastic beams can be expressed as a function of αand ηby lb=l 2p(2ηcot α+ 2η−1)2+ 1, where lis the dimension of the unit cell, as illustrated in Fig. 19 (right). The parameters ηand α must satisfy the geometric conditions 0≤η≤1 2,0≤α≤π, −1 2≤ηcot α≤1 2,(72) 21 Fig. 19. The re-entrant lattice with square rigid inclusions (left) is formed by tassellating the square unit cell (right) along the directions e1,e2. The rigid inclusions (gray areas surrounded by black thick lines) are fully characterized by the parameters ηand αthat satisfy condition (72). Consequently, the slenderness of the beams is given by Λ = lbrA J=l 2rA Jp(2ηcot α+ 2η−1)2+ 1.(73) Two interesting limiting cases occur: when α=π/3and η= 1/2, the rigid inclusions intersect at the center of the unit cell; and when α=π/2and η= 0.3, the rigid inclusions flatten, effectively becoming rigid rods. These limiting configurations are shown in Fig. 20. Fig. 20. Limit geometries of the re-entrant square lattice with rigid inclusions. Elastic beams overlap the boundary of the rigid inclusions for α=π/3and η= 1/2(left); rigid inclusions degenerate into rigid rods for α=π/2and η= 0.3(right). Keeping in mind conditions (71) and using the dimensionless variables defined with Eq. (48), the expressions for the non-zero elastic constants are ˜ E1111 =√2(Θ + 12) Θp1+2η(1 + cot α)[η(1 + cot α)−1], ˜ E1212 =4k3pk2 2+ 1(k2−cot α)2 (1 + cot α)2{Λ2(k2+ 1)2(12k2 2+ Θ)(1 −sin 2α) + 4Θ(k2 2+ 1)2(1 + sin 2α)} ˜ E1122 =−√2(Θ −12)[2η(1 + cot α)−1] Θ{1+2η(1 + cot α)[η(1 + cot α)−1]}3/2, (74) which lead to the following three dimensionless engineering constants ˜ E=√2 4 (Θ + 12)2(k2 2+ 1)2−4k2 2(Θ −12)2 Θ(Θ + 12)h(k2 2+ 1)/2i5/2,˜ G=˜ E1212, ν =−2k2(Θ −12) (Θ + 12)(k2 2+ 1),(75) where Θis defined in Eq. (68), while k2=r4l2 b l2−1, k3=12(k2+ 1)2+ Θ(k2−1)2(1 + sin 2α) + 3Λ2(k2+ 1)2(1 −sin 2α),(76) 22 depend on the geometry of the unit cell and on Poisson’s ratio of the beam material, included in the parameter ϕ. The elastic constants given in Eqs. (75) are plotted in Fig. 21 as functions of ϕ(left) and Λ (right), for η= 0.3and α=π/4, which satisfy condition (72). The figure shows that, compared to the Euler-Bernoulli model, the Timoshenko beam model reduces the elastic properties of the effective material, particularly at low slenderness. It can also be observed that all the engineering constants defined in Eq. (75) decrease as ϕand Λincrease, and eventually approach constant values as {ϕ, Λ} → ∞, Fig. 21. However, in the left plot of this figure, the equivalent Poisson’s ratio remains positive for 0<ϕ<1even with the Timoshenko beam model, when the slenderness is reasonably small (blue and red dotted lines). 0 0.5 1 8 13 18 23 28 0.8 1.3 1.8 2.3 2.8 0 2 4 6 8 10 - 4 - 2 0 2 4 6 8 - 0.4 - 0.3 - 0.2 - 0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0 5 10 15 20 - 4 - 2 0 2 4 6 8 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 Fig. 21. Effective Young’s modulus, shear modulus, and Poisson’s ratio as functions of the shear coefficient ϕ(left) and slenderness Λ(right), for fixed geometric parameters η= 0.3and α=π/4. Note the two different scales: one for ˜ Eon the left, and another for ˜ Gand νon the right of the frames. Unlike the case of the triangular lattice, bending stiffness is essential for the grid shown in Fig. 19. In the limit of vanishing bending stiffness, the following expressions are obtained ˜ E∞=2(k2 2−1)2 (k2 2+ 1)5/2,˜ G∞= 0, ν∞=−2k2 k2 2+ 1,(77) showing that the shear modulus ˜ G∞vanishes, while ˜ E∞>0, and ν∞<0for any combination of ηand αsatisfying the geometric condition (72). The two re-entrant grids shown in Figs. 16 and 19 are compared in the polar plots of Figs. 22 and 23, respectively. In these figures, the effective Poisson’s ratios are plotted as functions of the re-entrant angle αfor both the Euler-Bernoulli model (left panel) (equivalent to the Timoshenko model with νs=−1), and the Timoshenko model with νs= 0.3and κ= 5/6(right panel). In particular, since the equivalent material of the grid shown in Fig. 16 is anisotropic, the Poisson’s ratios can become unbounded and arbitrarily vary over the real numbers [45,46]. This behavior is clearly illustrated in Fig. 22, where ν12 and ν21 are plotted as functions of αfor a highly stubby beam (Λ = 1). For this slenderness value, the Euler-Bernoulli beam theory (left) predicts an unexpected non-auxetic response. In contrast, when shear deformability is taken into account (right), the equivalent material exhibits auxetic behavior for any re-entrant angle α < 90◦. When the re-entrant mechanism is provided by the square rigid inclusions (resulting in a cubic symmetry, Fig. 19), the Poisson’s ratio remains bounded within the range (−1,1) and the auxetic effect becomes strongly dependent on the geometric parameters, αand η. In particular, the colored zones in Fig. 23 represent the combinations of (α, η)that satisfy the geometric compatibility conditions (72). Consequently, outside these regions, the solution is not physically meaningful. Table 1reports the corresponding ranges of the angle αfor the re-entrant lattice with square rigid inclusions at various values of η. 23 0.0. 0.0. Fig. 22. Polar plots of the effective Poisson’s ratios for the re-entrant lattice without rigid inclusions, shown as a function of the re-entrant angle α∈(0,180◦): Euler-Bernoulli model, ϕ= 0, (left) and Timoshenko model, ϕ= 1.76, (right). Extremely stubby beams are considered with slenderness Λ=1. The gray region marks the range of angles corresponding to re-entrant geometries (α < 90◦) and no interference between beams (α > 60◦). Solid curves indicate positive Poisson’s ratios, while dashed curves denote negative (auxetic) values. The Timoshenko model (right) predicts consistent auxetic behavior throughout the re-entrant regime, whereas the Euler-Bernoulli model (left) incorrectly suggests non-auxetic behavior in the same range. 0.0. 0.0. Fig. 23. As in Fig. 22, but for the re-entrant lattice with rigid inclusions. The effective Poisson’s ratio νis shown for three different values of the geometric parameter η. Colored regions denote the angular intervals that satisfy both the re-entrant condition α < 90◦and the geometric compatibility conditions (72). The ranges of angles αfor each value of ηare reported in Table 1. Table 1. Range of compatible angles αfor varying semi-diagonal length of the inclusion, η, according to Eq. (72). Note that the upper limit, αmax = 90◦is independent of ηand corresponds to the case where the square inclusions collapse into lines (Fig. 20, right). Additionally, when η= 0.5, the rigid inclusions intersect at the center of the unit cell, causing the elastic beams to overlap the boundaries of the rigid inclusions (see Fig. 20, left). η αmin [◦]αmax [◦] 0.1 11.31 90 0.2 21.80 90 0.3 30.96 90 0.4 38.66 90 0.5 45 90 24 Deformed shapes of the grid with rigid inclusions are shown in Fig. 24 for both the EulerBernoulli model (left) and the Timoshenko model (right), with parameters Λ=1,α=π/4, and η= 0.3. The Euler-Bernoulli model incorrectly predicts a positive Poisson’s ratio, in agreement with the polar plot in Fig. 23 (left, yellow line). In contrast, the Timoshenko beam grid exhibits lateral expansion along the direction e2under horizontal uniform stress, resulting in a negative Poisson’s ratio of the equivalent material, as predicted by the yellow line in the right panel of the same figure. (a) Euler-Bernoulli beam model with Λ=1,ϕ= 0 yielding ν= 0.325 (b) Timoshenko beam model with Λ=1,ϕ= 1.76 yielding ν=−0.2 Fig. 24. As in Fig. 10, except that the lattice is now re-entrant (α=π/4and η= 0.3) and includes square rigid inclusions. The Timoshenko model (b) correctly exhibits auxetic behavior, in contrast to the Euler-Bernoulli model (a). 5.4 An oblique lattice An oblique lattice is now analyzed, as illustrated in Fig. 25. The grid is constructed by tiling a parallelogram-shaped unit cell with side lengths land l/2, oriented along the a1and a2directions, respectively. The a2direction is inclined at an angle of π/3with respect to the e1axis. The horizontal and oblique beams have rectangular cross-sections with dimensions b×hand b/2×2h, respectively, such that all beams share the same slenderness, Λ = √12 l/b, and cross-sectional area, A=bh. Fig. 25. Oblique lattice (left) composed of elastic Timoshenko beams arranged on the plane to form an L-shape pattern inclined at an angle of π/3(right). This condition, together with the chosen beam lengths, is sufficient to ensure that the equivalent material is anisotropic, with an elasticity tensor given by   ˜ E1111 ˜ E1122 ˜ E1112 ·˜ E2222 ˜ E2212 · · ˜ E1212 =1 8√3Θ   33Θ + 24 3Θ −24 √3Θ + 24√3 ·9Θ + 24 3√3Θ −24√3 · · 3Θ + 72  ,(78) where dimensionless quantities are defined according to Eq. (48), and Θdepends on both Λand ϕ, as given by Eq. (68). It is worth noting that the elasticity tensor in Eq. (78) is fully populated, implying a coupling between axial and shear deformations. This coupling is better understood by 25 [29] G. Bordiga, D. Bigoni, and A. Piccolroaz. “Tensile material instabilities in elastic beam lattices lead to a bounded stability domain”. In: Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 380.2231 (Sept. 5, 2022), p. 20210388. doi:10.1098/rsta.2021.0388. [30] A. Bensoussan, J.-L. Lions, and G. Papanicolaou. Asymptotic analysis for periodic structures. North-Holland Pub. Co., 1978. [31] D. Cioranescu and S. J. S. Paulin. Homogenization of reticulated structures. Springer, 1999. [32] D. Bigoni. Nonlinear Solid Mechanics: Bifurcation Theory and Material Instability. 1st ed. Cambridge University Press, July 30, 2012. doi:10.1017/CBO9781139178938. [33] M. E. Gurtin. An introduction to continuum mechanics. Academic press, 1982. [34] S. G. Lekhnitskii, P. Fern, J. J. Brandstatter, and E. H. Dill. “Theory of Elasticity of an Anisotropic Elastic Body”. In: Physics Today 17 (Jan. 1964), pp. 84–84. doi:10.1063/1. 3051394. [35] T. C. Hales. “The Honeycomb Conjecture”. In: Discrete & Computational Geometry 25.1 (Jan. 1, 2001), pp. 1–22. doi:10.1007/s004540010071. [36] K. E. Evans. “Auxetic polymers: a new range of materials”. In: Endeavour 15.4 (Jan. 1991), pp. 170–174. doi:10.1016/0160-9327(91)90123-S. [37] K. E. Evans. “The design of doubly curved sandwich panels with honeycomb cores”. In: Composite Structures 17.2 (Jan. 1, 1991), pp. 95–111. doi:10.1016/0263-8223(91)900646. [38] L. J. Gibson and M. F. Ashby. Cellular Solids: Structure and Properties. 2nd ed. Cambridge University Press, May 1, 1997. doi:10.1017/CBO9781139878326. [39] I. Masters and K. Evans. “Models for the elastic deformation of honeycombs”. In: Composite Structures 35.4 (Aug. 1996), pp. 403–422. doi:10.1016/S0263-8223(96)00054-2. [40] I. J. Gibson and M. F. Ashby. “The mechanics of three-dimensional cellular materials”. In: Proceedings of the Royal Society of London. A. Mathematical and Physical Sciences 382.1782 (July 8, 1982), pp. 43–59. doi:10.1098/rspa.1982.0088. [41] L. Slepyan. “Feeding and Dissipative Waves in Fracture and Phase Transition. III. Triangularcell Lattice”. In: Journal of the Mechanics and Physics of Solids 49.12 (Dec. 2001), pp. 2839– 2875. doi:10.1016/S0022-5096(01)00053-9. [42] K. E. Evans, M. A. Nkansah, I. J. Hutchinson, and S. C. Rogers. “Molecular network design”. In: Nature 353.6340 (Sept. 1991), pp. 124–124. doi:10.1038/353124a0. [43] R. Lakes. “Foam Structures with a Negative Poisson’s Ratio”. In: Science 235.4792 (Feb. 27, 1987), pp. 1038–1040. doi:10.1126/science.235.4792.1038. [44] R. Lakes. “Advances in negative Poisson’s ratio materials”. In: Advanced Materials 5.4 (Apr. 1993), pp. 293–296. doi:10.1002/adma.19930050416. [45] T. C. T. Ting and T. Chen. “Poisson’s Ratio for Anisotropic Elastic Materials Can Have No Bounds”. In: The Quarterly Journal of Mechanics and Applied Mathematics 58.1 (Feb. 1, 2005), pp. 73–82. doi:10.1093/qjmamj/hbh021. [46] P. Boulanger and M. Hayes. “Poisson’s Ratio for Orthorhombic Materials”. In: Journal of Elasticity 50 (1998), pp. 87–9. doi:https://doi.org/10.1023/A:1007468812050. 32