Structured Porous Metamaterial

Abstract

A structured porous metamaterial includes a three-dimensional matrix of at least one repeating base unit. The matrix is formed from an array of at least eight base units, each base unit including a platonic solid including at least one shaped void, wherein each base unit has void geometry tailored to provide a porosity of between 0.3 and 0.97, and to provide the metamaterial with a response that includes a Poisson's ratio of 0 to 0.5 when under tension and compression, or negative linear compression (NLC), negative area compression (NAC), zero linear compression (ZLC), or zero area compression (ZAC) behaviour when under pressure.

Claims

1. A structured porous metamaterial comprising a three-dimensional matrix of at least one repeating base unit, the matrix formed from an array of at least eight base units, each base unit comprising a platonic solid including at least one shaped void, wherein the geometry of the at least one shaped void of each base unit is tailored to: provide a porosity of between 0.3 and 0.97; and provide the metamaterial with a response comprising at least one of: a Poisson's ratio of 0 to 0.5 when under tension and compression; or negative linear compression (NLC), negative area compression (NAC), zero linear compression (ZLC), or zero area compression (ZAC) behaviour when under pressure.

2. A metamaterial according to claim 1, wherein the base unit comprises at least one of a tetrahedron, cube, cuboid, parallelepiped, octahedral, dodecahedron, or icosahedron.

3. (canceled)

4. A metamaterial according to claim 1, wherein the base unit includes a geometric center, and the geometry of the void is centered about that geometric center.

5. A metamaterial according to claim 1, wherein the base unit includes a width, height and length, and the at least one dimension of the base geometric shape of the void is larger than at least one of the width, height or length of the base unit.

6. A metamaterial according to claim 1, wherein the void comprises at least one of: a truncated form of a base geometric shape; or an interconnected combination of at least two geometric shapes.

7. (canceled)

8. A metamaterial according to claim 1, wherein the void includes an opening in at least one, preferably two sides of the base unit.

9-10. (canceled)

11. A metamaterial according to claim 1, wherein the base geometric shape of the voids comprises at least one of spherical, ovoid, ellipsoid, cubic, cuboid, parallelepiped, hyperboloid, conical.

12. A metamaterial according to claim 1, wherein the void geometry of the base unit is tailored to provide a porosity of one of: between 0.69 and 0.97 for a spherical shaped void; between 0.30 and 0.90 for regular non-spherical shaped voids; or between 0.3 and 0.98 for optimised shaped voids.

13. A metamaterial according to claim 1, wherein shaped void comprises an optimised shaped void formed using optimization algorithms, preferably bi-directional evolutionary structural optimization.

14. A metamaterial according to claim 1, wherein the base unit comprises a cube and the base geometric shape of the void comprises a sphere.

15. A metamaterial according to claim 1, wherein the base geometric shape of the void comprises shape having a greater central length than central height, the shape having a central length axis, the matrix of base units being arranged such that the central length axis of the void of each base unit is perpendicular to the central length axis of the void of each adjoining base unit.

16. A metamaterial according to claim 15, wherein the void shape comprises an ovoid or an ellipsoid.

17. A metamaterial according to claim 1, wherein the base unit is cubic and the shaped void is ellipsoid and wherein the porosity is between 0.3 and 0.87.

18. A metamaterial according to claim 1, wherein the base unit includes at least two shaped voids.

19. A metamaterial according to claim 1, comprising a three-dimensional matrix of at least two different repeating base units, comprising a first base unit comprising platonic solid including a first shaped void and a second base unit comprising platonic solid including a second shaped void.

20. (canceled)

21. A metamaterial according to claim 1, wherein the voids are composed of a compressible material, preferably a compressible material having a high compressibility.

22. A metamaterial according to claim 1, wherein the voids include at least one fluid, preferably at least one liquid.

23-27. (canceled)

28. A method of determining the configuration of a structured porous metamaterial comprising a three-dimensional matrix of at least one repeating base unit, comprising: determining a base unit topology using a structural optimization algorithm, each base unit comprising a platonic solid including at least one shaped void, the geometry of the at least one shaped void of each base unit being tailored to provide a metamaterial with a porosity of between 0.3 and 0.97 and a response comprising at least one of: a Poisson's ratio of 0 to 0.5 when under tension and compression; or negative linear compression (NLC), negative area compression (NAC), zero linear compression (ZLC), or zero area compression (ZAC) behaviour when under pressure; and simplifying the configuration of the at least one shaped void of each base unit to form a structural base unit; and forming a three-dimensional matrix from an array of at least eight structural base units.

29. A method according to claim 28, wherein the configuration of the shaped voids within each base unit is derived from a bi-directional evolutionary structural optimization (BESO) model.

30. A method according to claim 28, wherein the step of simplifying the configuration of the at least one shaped void of each base unit comprises reconfiguring the topology of the shaped void or voids to have a more regular geometric shape.

31-33. (canceled)

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0073] The present invention will now be described with reference to the figures of the accompanying drawings, which illustrate particular preferred embodiments of the present invention, wherein:

[0074] FIG. 1A provides the geometric configurations for a comparative three dimensional structure porous material without negative Poisson's ratio showing (A) the base cell unit; (B) a block of the comparative material comprising an 888 matrix of the base unit; and (C) representative volume unit of the comparative material.

[0075] FIG. 1B provides the geometric configurations for a three dimensional structure porous metamaterial according to the first embodiment of the present invention showing (A) the base cell unit; (B) a block of the inventive metamaterial comprising an 888 matrix of the base unit; and (C) representative volume unit of the inventive metamaterial.

[0076] FIG. 2 provides photographs of samples of the metamaterial shown in FIG. 1B (A) with supporting material fabricated using 3D printing and (B) without supporting material fabricated using 3D printing.

[0077] FIG. 3A shows the deformation patterns and thus buckling model for materials with (A) comparative a face-centred cubic cell (volume fraction: 51.0%) and (B) a cubic building cell according to the present invention (volume fraction: 12.6%).

[0078] FIG. 3B provides a view of force-displacement response of inventive metamaterial in two different directions D1 and D2.

[0079] FIG. 3C provides a comparison of nominal stress-strain curve of comparative structure porous material with face-centred cubic cell along three different loading directions.

[0080] FIG. 4 provides a comparison of deformation pattern of inventive metamaterial (volume fraction: 12.6%, load direction: D2 (FIG. 3), strain rate: 10.sup.3 s.sup.1) between (A) experiment and (B) finite element model.

[0081] FIG. 5 provides a comparison of nominal stress-strain curve of inventive metamaterial between experiment and finite element model for spherical voids and slightly ovoid shaped voids shaped (spherical with imperfection).

[0082] FIG. 6 provides a comparison of deformation pattern of an embodiments of the inventive metamaterial including slightly ovoid shaped voids (volume fraction: 12.6%, strain rate: 10.sup.3 s.sup.1) between (A) experiment and (B) finite element model.

[0083] FIG. 7A provides the geometric configurations for a three dimensional structure porous metamaterial with tetrahedron in cube building cell, showing (A) the base cell unit; (B) a block of the inventive metamaterial comprising an 888 matrix of the base unit; and (C) an isometric view of the representative volume unit of the inventive metamaterial.

[0084] FIG. 7B provides the geometric configurations for a three dimensional structure porous metamaterial with ellipsoid in cube building cell, showing (A) the base cell unit; (B) a block of the inventive metamaterial comprising an 888 matrix of the base unit; and (C) an isometric view of the representative volume unit of the inventive metamaterial.

[0085] FIG. 8A provides the deformation pattern for the metamaterial shown in FIG. 7A under load, showing (A) deformation pattern for bulk material (888) in xz plane; (B) deformation pattern for bulk material (888) in yz plane; (C) deformation pattern for the representative volume unit (222) in xz plane; and (D) an isometric view of the deformation pattern of the representative volume unit (222).

[0086] FIG. 8B provides the deformation pattern for the metamaterial shown in FIG. 7B under load, showing (A) deformation pattern for bulk material (888) in xz plane; (B) deformation pattern for bulk material (888) in yz plane; (C) deformation pattern for the representative volume unit (222) in xz plane; and (D) an isometric view of the deformation pattern for the representative volume unit (222).

[0087] FIG. 9 provides the geometric configurations for a three dimensional structure porous metamaterial with NLC according to the second embodiment of the present invention showing (A) the optimised building cell from BESO; (B) the simplified building cell unit; and (C) a block of the comparative material comprising an 888 matrix of the building cell unit.

[0088] FIG. 10 provides a comparison of deformation pattern of inventive NC metamaterial with NLC shown in FIG. 9 between (A) experiment and (B) finite element model; and (C) Comparison of strain-pressure history between FE results and experimental data for NLC material under pressure.

[0089] FIG. 11 provides the geometric configurations for a three dimensional structure porous metamaterial with NAC according to the second embodiment of the present invention showing (A) the optimised half building cell from BESO; (B) the optimised building cell from BESO; (C) the simplified building cell unit; and (D) a block of the material comprising an 888 matrix of the building cell unit.

[0090] FIG. 12 provides the geometric configurations for a three dimensional structure porous metamaterial with ZLC according to the second embodiment of the present invention showing (A) the optimised half building cell from BESO; (B) the optimised building cell from BESO; (C) the simplified building cell unit; and (D) a block of the material comprising an 888 matrix of the building cell unit.

[0091] FIG. 13 provides the geometric configurations for a three dimensional structure porous metamaterial with ZAC according to the second embodiment of the present invention showing (A) the optimised half building cell from BESO; (B) the optimised building cell from BESO; (C) the simplified building cell unit; and (D) a block of the material comprising an 888 matrix of the building cell unit.

DETAILED DESCRIPTION

[0092] The present invention generally relates to a series of 3D structured porous metamaterial with specific deformation pattern under applied loading, and more particularly a 3D structured porous metamaterial having at least one of: [0093] a negative Poisson's ratio under uniaxial tensile or compression; and/or [0094] zero or negative compressibility under uniform pressure, such as negative linear compressibility (NLC), negative area compressibility (NAC), zero linear compressibility (ZLC) and/or zero area compressibility (ZAC).

[0095] The initial design of the microstructure of an auxetic metamaterial form of the present invention originates from using a three-dimensional repeating matrix formed from a base unit comprising a platonic solid such as a cube having a shaped void space such as a sphere or ellipsoid. The platonic solid provides a repeatable and stackable base structure, and the shaped void imparts the required characteristic to the void space and the surrounding base unit framework structure (around the void). The void geometry of each base unit is tailored to provide a porosity of between 0.3 and 0.97; and provide the metamaterial with a response under tension and compression having a Poisson's ratio of 0 to 0.5. The specific porosity depends on the type of shaped void used. Therefore the porosity is typically between 0.69 and 0.97 for a spherical shaped void; between 0.30 and 0.90 for regular non-spherical shaped voids; or between 0.3 and 0.97 for optimised shaped voids. Furthermore, as will be explained in more detailed below with reference to specific example material configurations, this structure imparts a tailored deformation character to the material, with the negative Poisson's ratios achieved through the a specific deformation characteristic of the voids (alternating opening and closing pattern of adjacent voids) in the material combined with the spatial rotation and translation of a rigid part of base unit material accompanied by the bending and stretching of the thinner or more flexible part of the base unit material.

[0096] The initial design of the microstructure of the zero or negative compressibility (NC) metamaterial form of the present invention originates from using a three-dimensional repeating matrix formed from a base unit comprising a platonic solid, such as a cube, having one or more shaped void spaces. The shape of the voids within that base unit and thus the topology of those building unit is derived from a bi-directional evolutionary structural optimization (BESO) model formed to provide the desired NC properties using the desired base unit (again for example a cube). That BESO result is then altered to simplify the topology of the void or voids to have a more regular shape. This simplified shape is typically more suitable for 3D printing construction. The platonic solid provides a repeatable and stackable base structure, and the shaped void or voids in the base unit cell (an optimised shaped void) imparts the required characteristic to the void space and the surrounding base unit framework structure (around the void). The void geometry (the optimised shape of the void or voids) of each base unit is tailored to provide a porosity of between 0.3 and 0.95; and provide the NC metamaterial with a response under uniform pressure having one of the following behaviour: NLC, NAC, ZLC and ZAC.

[0097] The material of the base unit can be polymeric including, but not limited to, unfilled or filled vulcanized rubber, natural or synthetic rubber, cross-linked elastomer, thermoplastic vulcanizate, thermoplastic elastomer, block copolymer, segmented copolymer, cross-linked polymer, thermoplastic polymer, filled or unfilled polymer, or epoxy. In other embodiments, the material of the base unit but may also be non-polymeric including, but not limited to, metallic and ceramic and composite materials. Exemplary metals include aluminium, magnesium, titanium, iron and alloys thereof.

[0098] Fabrication of 3D structures according to the present invention can be achieved through 3D printing, dissolving or melting patterned voids from a base material and sintering techniques well known in the art.

Bi-Directional Evolutionary Structural Optimization (BESO)

[0099] The optimization method used for the initial design of the microstructure of the zero or negative compressibility (NC) metamaterial form of the present invention is based on the bi-directional evolutionary structural optimization (BESO). The basic idea of BESO is that by gradually removing inefficient material from a ground structure and redistributing the material to the most critical locations, the structure evolves towards an optimum.

[0100] For a 3D continuum material the ground structure is a unit cubic cell and the material properties (e.g. elasticity matrix) is determined using the homogenization theory. For NC forms of the present invention, the BESO method was applied to the design of materials of four types, namely, NLC, NAC, zero linear compressibility (ZLC) and zero area compressibility (ZAC).

Determining Linear, Area and Volume Compressibilities of a Material by Homogenization

[0101] A cellular material consisting of a base material and voids is often modelled as a microstructure of a periodic base cell (PBC) using finite element (FE) analysis. According to the homogenization theory (Hassania, B., Hintona, E., 1998. A review of homogenization and topology optimization Ihomogenization theory for media with periodic structure. Computers & Structures 69 (6), 707-717), the effective elastic constants can be expressed as

[00001]$\begin{array}{cc}{E}_{\mathrm{ij}}^H=\underset{e=1}{\overset{\mathrm{NE}}{.Math.}}.Math.\left(\frac{1}{Y_e}.Math.{}_{Y_e}.Math.\left({\mathrm{.Math.}}_i^{0^T}-e_i^T\right).Math.E\left({\mathrm{.Math.}}_j^0-{\mathrm{.Math.}}_j\right).Math.Y_e\right).Math..Math.\left(i,j=1.Math..Math.\mathrm{to}.Math..Math.6.Math..Math.\mathrm{for}.Math..Math.3.Math.D\right)& \left(1\right)\end{array}$

where E is the elastic matrix of the base material, NE is the number of elements, .sub.i.sup.0 is the i-th unit strain field and .sub.i is the corresponding induced strain field.

[0102] For 3D materials, it involves applying six cases of periodic boundary conditions and unit strain fields. Then the 66 make up the elasticity matrix E.sup.H. The homogenized compliance matrix C.sup.H is the inverse of E.sup.H, i.e.

C.sup.H=[C.sub.ij]=E.sup.H-1(2)

[0103] As the materials studied here is orthotropic, there is no axial-shear coupling and thus the 33 sub-matrix of the axial components can be extracted as below

[00002]$\begin{array}{cc}{C}_A^H=\left[C_{\mathrm{ij}}^A\right]=\left[\begin{array}{ccc}{C}_{11}& C_{12}& C_{13}\\ C_{21}& C_{21}& C_{21}\\ C_{31}& C_{32}& C_{33}\end{array}\right]=E_A^{H^{-1}}& \left(3\right)\end{array}$

[0104] Based on the above compliance matrix, the linear compressibility in axis i (i=1, 2, 3) can be expressed as

.sub.Li=C.sub.i1+C.sub.i2+C.sub.i3(4)

[0105] which has the dimension of inverse of stress. The area compressibility in the ij plane is defined as

.sub.Aij=.sub.Li+.sub.Lj, ij(5)

and the volume compressibility as

.sub.v=.sub.L1+.sub.L2+.sub.L3(6)

[0106] It is noted that Eq. (6) is the summation of the nine constants of the compliance matrix in Eq. (3), which is numerically equivalent to twice the strain energy of the microstructure under the unit hydrostatic stress. Since the strain energy is greater than or equal to zero, it is clear that for orthotropic materials the volume compressibility can either be positive or zero.

[0107] 1. Negative Linear Compressibility

[0108] A typical optimization problem is usually defined in terms of the objective function(s) and constraints(s). Here an obvious choice of the objective function is the linear compressibility in a particular direction. For example, we may aim to minimize the compressibility in axis 3, .sub.L3=C.sub.31+C.sub.32+C.sub.33. We choose the solid material as the initial design for the optimization process. For such an initial design, C.sub.31 and C.sub.32 are both negative and therefore .sub.L3 can be re-written as .sub.L3=(|C.sub.31|+|C.sub.32|)+C.sub.33. It is noted that .sub.L3 is initially positive and one way to drive it to become negative is to increase the weighing of the two negative terms, i.e. .sub.L3=(p|C.sub.31|p|C.sub.32|)+C.sub.33 with p>1. Here p can be regarded as a stress factor or a penalty parameter: instead of the unit stress ={1,1,1}, a modified stress ={p,p,1} is applied during the optimization process. The lower bound of p is 1, which must be reached on convergence. The upper bound of p is specified by assuming the linear compressibility equal to zero, i.e.

.sub.L3=p.sup.upperC.sub.31+p.sup.upperC.sub.32+C.sub.33=0(7a)

[0109] In order to maintain the orthotropy of the material, Eq. (7a) is re-written as

.sub.L3=(p.sup.upperC.sub.31+p.sup.upperC.sub.32+p.sup.upperC.sub.13+p.sup.upperC.sub.23+2C.sub.33)=0(7b)

and p.sup.upper is found to be

[00003]$\begin{array}{cc}{p}^{\mathrm{upper}}=-\frac{2.Math.C_{33}}{C_{31}+C_{32}+C_{13}+C_{23}}& \left(7.Math.c\right)\end{array}$

[0110] With p[1, p.sup.upper] specified, the value of p is to be determined. Because of the same p value being applied to axes 1 and 2, the resulting material is to be symmetrical to the 45 degree line in plane 1-2.

[0111] Next we discuss what constraints should be included in the optimization process apart from the volume constraint. As the NLC design is likely to be very flexible, it is necessary to prevent the design from becoming singular. In other words, we need to maintain reasonable stiffness. The stiffness in axis 3 is maintained by including C.sub.33 in the objective function. The stiffness in axes 1 and 2 can be considered by specifying a constraint on C.sub.11 and C.sub.22, for example, by requiring them to be less than 1/E*, where E* is a prescribed stiffness target.

[0112] From the above discussions, the design of NLC materials can be treated as the following optimization:

[00004]$\begin{array}{cc}\mathrm{Minimize}.Math..Math.{}_{L.Math..Math.3}=\frac{1}{2}.Math.\left({\mathrm{pC}}_{31}+{\mathrm{pC}}_{32}+{\mathrm{pC}}_{13}+{\mathrm{pC}}_{23}+2.Math.C_{33}\right)& \left(8.Math.a\right)\\ \mathrm{Subject}.Math..Math.\mathrm{to}.Math..Math.C_{11}\frac{1}{E^{*}}& \left(8.Math.b\right)\\ C_{22}\frac{1}{E^{*}}& \left(8.Math.c\right)\\ C_{11}=C_{22}& \left(8.Math.d\right)\\ p=1.Math..Math.\mathrm{and}& \left(8.Math.e\right)\\ \underset{e=1}{\overset{\mathrm{NE}}{.Math.}}.Math..Math.V_e.Math.x_e=V,x_e=x_{\min }.Math..Math.\mathrm{or}.Math..Math.1& \left(8.Math.f\right)\end{array}$

where V is the prescribed volume, V.sub.e is the volume of element e, and x.sub.e is the design variable, with x.sub.e=x.sub.min for void and x.sub.e=1 for solid.

[0113] The Lagrangian function combining the objective function and constraints is

[00005]$\begin{array}{cc}{f}_L=\frac{1}{2}.Math.\left({\mathrm{pC}}_{31}+{\mathrm{pC}}_{32}+{\mathrm{pC}}_{13}+{\mathrm{pC}}_{23}+2.Math.C_{33}\right)+\left(C_{11}-\frac{1}{E^{*}}\right)+\left(C_{22}-\frac{1}{E^{*}}\right)& \left(9\right)\end{array}$

[0114] Since C.sub.11=C.sub.22, the same Lagrangian multiplier is applied to constraints (8b) and (8c).

Sensitivity Analysis of Elasticity and Compliance Constants

[0115] The sensitivity of the Lagrangian function with respect to the design variable is

[00006]$\begin{array}{cc}\frac{f_L}{x_e}=\frac{1}{2}.Math.\left(p.Math.\frac{C_{31}}{x_e}+p.Math.\frac{C_{32}}{x_e}+p.Math.\frac{C_{13}}{x_e}+p.Math.\frac{C_{23}}{x_e}+2.Math.C_{33}\right)+\left(\frac{C_{11}}{x_e}+\frac{C_{22}}{x_e}\right)& \left(10\right)\end{array}$

which calls for the sensitivity analysis of the compliance constants. To achieve this, the sensitivity of elasticity constants can be obtained by using the adjoint method (Bendsee, M. P., Sigmund, O., 2003. Topology optimization: theory, methods and applications 2nd ed. Springer, Berlin). From Eq. (1), the sensitivity of E.sub.ij.sup.H can be expressed as

[00007]$\begin{array}{cc}\frac{E_{\mathrm{ij}}^H}{x_e}=\frac{1}{Y_e}.Math.{}_{Y_e}.Math.\left({\mathrm{.Math.}}_i^{0^T}-{\mathrm{.Math.}}_i^T\right).Math.\frac{E}{x_e}.Math.\left({\mathrm{.Math.}}_j^0-{\mathrm{.Math.}}_j\right).Math.Y_e& \left(11\right)\end{array}$

[0116] The term

[00008]$\frac{E}{x_e}$

depends on the function used for interpolating the Young's modulus E. Here the interpolation scheme is based on

[00009]$\begin{array}{cc}E\left(x_e\right)=E_{b.Math..Math.1}+\frac{x_e\left(E_{b.Math..Math.2}-E_{b.Math..Math.1}\right)}{1+q\left(1-x_e\right)}& \left(12\right)\end{array}$

where E.sub.b1 and E.sub.b2 are the Young's moduli of the base materials and q acts as a penalty factor. Typical values of q are equal to or greater than 3. For the examples considered in this paper, it is found that q=6 gives the best results. The present study is focused on designing cellular materials and therefore one of the base materials is void, i.e. either E.sub.b1 or E.sub.b2 is approaching zero.

[0117] Making use of Eq. (3), the sensitivity of the mean compliance matrix C.sup.H is calculated by using the chain rule, i.e.

[00010]$\begin{array}{cc}\frac{C_{\mathrm{ij}}}{x_e}=\underset{k,l=1}{\overset{3}{.Math.}}.Math.\frac{C_{\mathrm{ij}}}{E_{\mathrm{kl}}^H}.Math.\frac{E_{\mathrm{kl}}^H}{x_e}.Math.\left(i,j=1.Math..Math.\mathrm{to}.Math..Math.3\right)& \left(13\right)\end{array}$

which can be calculated analytically by following a series of matrix operations.

Sensitivity Number

[0118] The above sensitivity analysis forms the basis of the sensitivity number which is used as the search criterion in the BESO solution process. From Eq. (10), the sensitivity number is defined as

[00011]$\begin{array}{cc}{}_e=-\frac{f_L}{x_e}& \left(14\right)\end{array}$

[0119] The sensitivity number .sub.e is then filtered through a spherical range of radius r.sub.min to obtain a weighted average, i.e.

{tilde over ()}.sub.e=(.sub.e)(15)

[0120] Taking the centre of a brick element e as reference, the neighbouring elements within the radius r.sub.min are included for the calculation of the average sensitivity of element e. The contributions from neighbouring elements depend on the sensitivity of each element and its distance to element e. Details of the filtering methodology are presented in Huang, X., Xie, Y. M., 2010. Evolutionary Topology Optimization of Continuum Structures: Methods and Applications. John Wiley & Sons, Chichester, England, the contents of which should be understood to be incorporated into this specification by this reference.

[0121] The sensitivity of the compliance matrix is filtered in the same way, i.e.

[00012]$\begin{array}{cc}\frac{{\tilde{C}}_{\mathrm{ij}}}{x_e}=f\left(\frac{C_{\mathrm{ij}}}{x_e}\right)& \left(16\right)\end{array}$

[0122] Assuming there are totally m elements modified in one iteration, the increment of C.sup.H is then

[00013]$\begin{array}{cc}.Math..Math.C_{\mathrm{ij}}\underset{e=1}{\overset{m}{.Math.}}.Math.\frac{{\tilde{C}}_{\mathrm{ij}}}{x_e}.Math..Math..Math.x_e& \left(17\right)\end{array}$

[0123] The predicted mean compliance after the modification is

C.sub.ijC.sub.ij+C.sub.ij(18)

BESO Procedure

[0124] Like most numerical methods based on sensitivity analysis, BESO performs the search for the optimal solution iteratively until certain criteria are satisfied. Details of the solution procedure are as follows:

A. Parameters

[0125] There are three parameters which control the step length of iteration, namely the evolutionary ratio ER, the maximum ratio R.sub.max and the maximum ratio of added elements AR.sub.max. Assume that there are totally NE elements in the design domain and the volume constraint (the target volume) is V. The volumes of the current and the next iterations are V.sup.k and V.sup.k+1, respectively. V.sup.k+1 is predicted as V.sup.k+1=V.sup.k(1ER) and the threshold for element modification is set as

NE.sub.thre=NEV.sup.k+1=NEV.sup.k(1ER)(19)

[0126] The modification according to the threshold is conducted as follows. First, sort the sensitivity numbers of the NE elements in a descend order. Then void elements above the threshold NE.sub.thre are switched to solid, and solid elements below the threshold are switched to void. As a result, the total numbers of elements removed and added are NR and NA, respectively.

[0127] The net number of modified element is NRNA, which is positive if the volume is approaching from the initial high value to the target. A parameter AR.sub.max is introduced to ensure that the number of added elements in one iteration is not too large, i.e. when the ratio NA/NE exceeds AR.sub.max, NA is reduced to NA.sub.max=AR.sub.maxNE.

[0128] Also, it is required that the total of NR and NA (or NA.sub.max if applicable) is not too high, that is,

[00014]$\begin{array}{cc}\frac{\mathrm{NR}+\mathrm{NA}}{\mathrm{NE}}{R}_{\max }& \left(20\right)\end{array}$

[0129] If the ratio is exceeded, the numbers of removed and added elements are reduced according to the following equations:

[00015]$\begin{array}{cc}{\mathrm{NR}}^{}=R_{\max }\mathrm{NE}\frac{\mathrm{NR}}{\mathrm{NR}+\mathrm{NA}}& \left(21\right)\\ {\mathrm{NA}}^{}=R_{\max }\mathrm{NE}\frac{\mathrm{NA}}{\mathrm{NR}+\mathrm{NA}}& \left(22\right)\end{array}$

B. The Overall Procedure

[0130] The outer loop of the BESO procedure is as follows: [0131] 1. Discretize the periodic base cell with finite elements and define the initial design. [0132] 2. Apply the periodic boundary conditions and corresponding unit strain fields. [0133] 3. For each boundary and unit strain case, conduct finite element analysis to obtain the induced strain field . [0134] 4. Calculate the elasticity matrix E.sup.H and the compliances matrix C.sup.H. [0135] 5. Determine the stress factor p and the Lagrangian multipliers (inner loops), as detailed in Section CStress factor and Lagrangian multipliers. [0136] 6. Calculate the sensitivity number {tilde over ()}.sub.e using Eqs. (1115). [0137] 7. Update the topology of the base cell according to {tilde over ()}.sub.e, using the threshold and parameters as detailed in Section AParameters. [0138] 8. Repeat Steps 2 to 7 until the objective function is stabilized between iterations.

C. Stress Factor and Lagrangian Multipliers

[0139] In Eq. (9) the Lagrangian function .sub.L has two unknowns, namely the stress factor p and the Lagrangian multiplier associated with the stiffness constraints. If the constraint is too stringent, i.e. the value of

[00016]$\frac{1}{E^{*}}$

is too small, the objective compressibility may not be reduced enough to be below zero. Therefore

[00017]$\frac{1}{E^{*}}$

should be reasonably high to allow the structure to be sufficiently flexible. In the early stage of iterations (starting from a solid structure as the initial design), the structure is quite stiff with the constraint

[00018]$C_{11}<\frac{1}{E^{*}}$

being satisfied and thus the Lagrangian multiplier is =0. Only the stress factor p needs to be solved at this stage. As the iterations continue, p will converge to unit and the stiffness will gradually reduce till C.sub.11 becomes greater than

[00019]$\frac{1}{E^{*}}.$

At this point the Lagrangian multiplier becomes activated and needs to be solved. Once p and are solved, they are averaged between the current and the last iterations, respectively.

D. Determination of the Stress Factor

[0140] The stress factor p is solved by a general bi-section method: [0141] 1. Calculate the upper bound of p using Eq. (7c). Then the search range of p is [1, p.sup.upper]. [0142] 2. Assign =0 and assign initial value of p=1. [0143] 3. Calculate the sensitivity number {tilde over ()}.sub.e using Eqs. (1115). [0144] 4. Obtain an assumed topology which has the volume equal to the constraint V. This is similar to Step 7 in Section 3.4.2. Now the threshold is NE.sub.thre=NEV. Sort the sensitivity numbers of the NE elements in a descend order. Then void elements above the threshold NE.sub.thre are switched to solid, and solid elements below the threshold are switched to void. [0145] 5. For the assumed new topology, estimate the compliance matrix C.sub.ij.sup.V (i, j=1, 3) using Eqs. (1618). Calculate the stress factor

[00020]$\begin{array}{cc}{p}^V=-\frac{2.Math.C_{33}^V}{C_{31}^V+C_{32}^V+C_{13}^V+C_{23}^V}& \left(23\right)\end{array}$ where the superscript V refers to the volume constraint V. [0146] 6. Check the convergence of p using the following criteria: [0147] a. p.sup.V=1. [0148] b. If a) is not satisfied, the change in p.sup.V between the current and the last iterations is small. [0149] 7. If the above convergence criteria are not satisfied, update p according to the bi-section rule, i.e. if p.sup.V>1, then p.sup.m+1=(p.sup.m+p.sup.upper), and the lower bound is reset as p.sup.lower=p.sup.m. [0150] 8. Proceed to iteration m+1. Repeat Steps 3-7 till convergence is reached. The final stress factor is taken as the average of p and p.sup.V. [0151] 9. At the convergence of p, also determine whether the following stiffness constraint is activated,

[00021]$\begin{array}{cc}{f}_{\mathrm{con}}=C_{11}^V-\frac{1}{E^{*}}& \left(24\right)\end{array}$ if .sub.con0, the Lagrangian multiplier is not activated and thus =0. if .sub.con>0, the Lagrangian multiplier is activated. Then proceed to next the step to determine , as described below.

E. Determination of the Lagrangian Multiplier

[0152] 1. From Eq. (24), calculate the upper bound of as follows

[00022]$\begin{array}{cc}{}^{\mathrm{upper}}={\left(\frac{C_{11}^V}{1/E^{*}}\right)}^2-1& \left(25\right)\end{array}$ It is noted that by applying a power of 2 in Eq. (25), the bound is further relaxed and thus the solution is searched in a wider range. [0153] 2. The stress factor is as per the already converged value of unit p=1. Assign initial value of =0. [0154] 3. Calculate the sensitivity number {tilde over ()}.sub.e using Eqs. (1115). [0155] 4. Obtain an assumed topology which has the volume equal to the constraint V, same as Step 4 in Section D. [0156] 5. For the assumed new topology, estimate the compliance matrix C.sub.jj.sup.V (j=1, 3) using Eqs. (1618), from which the value of C.sub.11.sup.V is obtained. [0157] 6. Check the convergence of C.sub.11.sup.V using the following criteria: [0158] 7.

[00023]$f_{\mathrm{con}}=C_{11}^V-\frac{1}{E^{*}}$

satisfies .sub.con=0. [0159] 8. If a) is not satisfied, the change in C.sub.11.sup.V between the current and the last iterations is small. [0160] 9. If the above convergence criteria are not satisfied, update A according to the bi-section rule, i.e. [0161] a) If .sub.con>0, then .sup.m+1=(.sup.m+.sup.upper), and the lower bound is reset as .sup.lower=.sup.m. [0162] b) If .sub.con<0, then .sup.m+1=(.sup.m+.sup.lower), and the upper bound is reset as .sup.upper=.sup.m. [0163] 10. Proceed to iteration m+1. Repeat Steps 3-7 till convergence is reached.

[0164] 2. Negative Area Compressibility

[0165] When addressing the NLC problem, we have introduced a stress factor p to drive the linear compressibility .sub.L3 to zero and then towards minimum. A similar strategy is used for NAC. Here the design objective is to minimize .sub.A23=.sub.L2+.sub.L3 and it is assumed that .sub.L2=.sub.L3. In order to make the material shrink more in axes 2 and 3, a larger stress is applied in axis 1 during the early stages of the optimization process. Therefore, a stress vector incorporating the stress factor p is defined as ={p,1,1}, where p1. The compressibility of the material under this stress is rewritten as

.sub.L3=(pC.sub.31+C.sub.32+pC.sub.13+C.sub.23+2C.sub.33)(27a)

.sub.L2=(pC.sub.21+C.sub.23+pC.sub.12+C.sub.32+2C.sub.22)(27b)

with C.sub.21=C.sub.31 and C.sub.33=C.sub.22(27c)

[0166] For the same reason as given in section 3.1, the upper bound of p can be obtained by setting .sub.L3=0, i.e.

[00024]$\begin{array}{cc}{p}^{\mathrm{upper}}=-\frac{C_{32}+C_{23}+2.Math.C_{33}}{C_{31}+C_{13}}& \left(28\right)\end{array}$

[0167] With p[1, p.sup.upper] specified, the value of p is to be determined using the bi-section method described in Section 3.4.3.1. After a number of iterations, p will converge to 1.

[0168] The optimization problem for designing NAC materials is stated as:

[00025]$\begin{array}{cc}\mathrm{Minimize}.Math..Math.{}_{A.Math..Math.23}={}_{L.Math..Math.2}+{}_{L.Math..Math.3}& \left(29.Math.a\right)\\ \mathrm{Subject}.Math..Math.\mathrm{to}.Math..Math.C_{11}\frac{1}{E^{*}}& \left(29.Math.b\right)\\ p=1& \left(29.Math.c\right)\\ \underset{e=1}{\overset{\mathrm{NE}}{.Math.}}.Math.V_e.Math.x_e=V,x_e=x_{\min }.Math..Math.\mathrm{or}.Math..Math.1& \left(29.Math.d\right)\end{array}$

[0169] The Lagrangian function is

[00026]$\begin{array}{cc}{f}_L=\frac{1}{2}.Math.\left({\mathrm{pC}}_{31}+C_{32}+{\mathrm{pC}}_{13}+C_{23}+2.Math.C_{33}\right)+\frac{1}{2}.Math.\left({\mathrm{pC}}_{21}+C_{23}+{\mathrm{pC}}_{12}+C_{23}+2.Math.C_{22}\right)+\left(C_{11}-\frac{1}{E^{*}}\right)& \left(30\right)\end{array}$

[0170] Similar to the Lagrangian function for NLC optimization given in Eq. (9), the above equation has two unknowns, namely the stress factor p and the Lagrangian multiplier . The same methodology as detailed for NLC is used to solve these two unknowns. Then the same overall BESO procedure described above is followed to find the optimal NAC design.

[0171] 3. Zero Linear Compressibility

[0172] For ZLC calculations, it was assumed that the material is under unit hydrostatic pressure and one way to measure its overall stiffness is the strain energy, i.e.

[00027]$=\frac{1}{2}.Math.\underset{\underset{j=1,3}{i=1,3}}{.Math.}.Math..Math.C_{\mathrm{ij}}.$

To design the stiffest material with zero linear compressibility (in axis 3), we state the optimization problem as

[00028]$\begin{array}{cc}\mathrm{Minimize}.Math..Math.=\frac{1}{2}.Math.\underset{\underset{j=1,3}{i=1,3}}{.Math.}.Math.C_{\mathrm{ij}}& \left(31.Math.a\right)\\ \mathrm{Subject}.Math..Math.\mathrm{to}.Math..Math.{}_{L.Math..Math.3}=0& \left(31.Math.b\right)\\ {}_{L.Math..Math.1}={}_{L.Math..Math.2}& \left(31.Math.c\right)\\ \underset{e=1}{\overset{\mathrm{NE}}{.Math.}}.Math..Math.V_e.Math.x_e=V,x_e=x_{\min }.Math..Math.\mathrm{or}.Math..Math.1& \left(31.Math.d\right)\end{array}$

[0173] The Lagrangian function is

[00029]$\begin{array}{cc}{f}_L=\frac{1}{2}.Math.\left(\underset{\underset{j=1,3}{i=1,3}}{.Math.}.Math.C_{\mathrm{ij}}\right)+\frac{1}{2}.Math.\left(C_{31}+C_{32}+C_{13}+C_{23}+2.Math.C_{33}\right)+{}_{12}\left({}_{L.Math..Math.1}-{}_{L.Math..Math.2}\right)& \left(32\right)\end{array}$

[0174] Due to the cubic symmetry of the initial design, .sub.L1=.sub.L2 is satisfied from the beginning and thus the last term in the above equation vanishes. The first multiplier is solved by using the bi-section method as detailed below.

[0175] The overall BESO procedure (outer-loop) is similar to that described previously. At each iteration an inner loop is conducted to solve the Lagrangian multiplier . Its value is then averaged between the current and last iterations. The procedure to determine A is as follows. [0176] 1. Assume that varies in the range of [0,1] and assign initial values of .sup.lower and .sup.upper to be equal to 0 and 1, respectively. [0177] 2. Assign the initial value of to be equal to 0. [0178] 3. Calculate the sensitivity number {tilde over ()}.sub.e using Eqs. (1115). [0179] 4. Obtain an assumed topology which has the volume equal to the constraint V, in the same way as in Step 4 in previous stress factor procedure. [0180] 5. For the assumed new topology, estimate the compliance matrix C.sub.ij.sup.V (i, j=1, 3) using Eqs. (1618). Calculate the compressibility .sub.L3.sup.V. [0181] 6. Check the convergence of .sub.L3.sup.V using the following criteria: [0182] a) .sub.L3.sup.V satisfies .sub.L3.sup.V=0. [0183] b) If a) is not satisfied, the change in .sub.L3.sup.V between the current and the last iterations is small. [0184] 7. If the above convergence criteria are not satisfied, update according to the bi-section rule, i.e. [0185] a) If .sub.L3.sup.V>0, then .sup.m+1=(.sup.m+.sup.upper), and the lower bound is reset as .sup.lower=.sup.m. [0186] b) If .sub.L3.sup.V<0, then .sup.m+1=(.sup.m+.sup.lower), and the upper bound is reset as .sup.upper=.sup.m. [0187] 8. Proceed to iteration m+1. Repeat Steps 3-7 till convergence is reached. [0188] 9. On convergence, assuming that a stress vector ={p,p,1} is applied to the material so that .sub.L3.sup.V=0, calculate the stress factor as

[00030]$\begin{array}{cc}{p}^V=-\frac{2.Math.C_{33}^V}{C_{31}^V+C_{32}^V+C_{13}^V+C_{23}^V}& \left(33\right)\end{array}$

p.sup.V is then used to modify the Lagrangian function as follows

[00031]$\begin{array}{cc}{f}_L=\frac{1}{2}.Math.\left(\underset{\underset{j=1,3}{i=1,3}}{.Math.}.Math.C_{\mathrm{ij}}\right)+\frac{1}{2}.Math.\left(p^V.Math.C_{31}+p^V.Math.C_{32}+p^V.Math.C_{13}+p^V.Math.C_{23}+2.Math.C_{33}\right)& \left(34\right)\end{array}$

[0189] This function is used to calculate the sensitivity of the subsequent iteration in the outer-loop.

[0190] 4. Zero Area Compressibility

[0191] For ZAC calculations, the following problem statement was followed:

[00032]$\begin{array}{cc}\mathrm{Minimize}.Math..Math.=\frac{1}{2}.Math.\underset{\underset{j=1,3}{i=1,3}}{.Math.}.Math.C_{\mathrm{ij}}& \left(39.Math.a\right)\\ \mathrm{Subject}.Math..Math.\mathrm{to}.Math..Math.{}_{L.Math..Math.3}=0& \left(39.Math.b\right)\\ {}_{L.Math..Math.2}=0& \left(39.Math.c\right)\\ \underset{e=1}{\overset{\mathrm{NE}}{.Math.}}.Math..Math.V_e.Math.x_e=V,x_e=x_{\min }.Math..Math.\mathrm{or}.Math..Math.1& \left(39.Math.d\right)\end{array}$

[0192] The Lagrangian function is

[00033]$\begin{array}{cc}{f}_L=\frac{1}{2}.Math.\left(\underset{\underset{j=1,3}{i=1,3}}{.Math.}.Math.C_{\mathrm{ij}}\right)+\frac{1}{2}.Math.\left(C_{31}+C_{32}+C_{13}+C_{23}+2.Math.C_{33}\right)+\frac{1}{2}.Math.\left(C_{21}+C_{23}+C_{12}+C_{32}+2.Math.C_{22}\right)& \left(40\right)\end{array}$

[0193] The procedure to solve the Lagrangian multiplier is similar to that in discussed above for NLC calculations. For Step 9 in calculating the stress factor, the stress vector assumed here is ={p,1,1}, where p1. By setting .sub.L3=0 and .sub.L2=0, the stress factor is

[00034]$\begin{array}{cc}{p}^V=-\frac{C_{32}^V+C_{23}^V+2.Math.C_{33}^V}{C_{31}^V+C_{13}^V}=-\frac{C_{23}^V+C_{32}^V+2.Math.C_{22}^V}{C_{21}^V.Math.{+}_{12}^V}& \left(37\right)\end{array}$

which is used to modify the Lagrangian function as follows

[00035]$\begin{array}{cc}{f}_L=\frac{1}{2}.Math.\left(\underset{\underset{j=1,3}{i=1,3}}{.Math.}.Math.C_{\mathrm{ij}}\right)+\frac{1}{2}.Math.\left(p^V.Math.C_{31}+C_{32}+p^V.Math.C_{13}+C_{23}+2.Math.C_{33}\right)+\frac{1}{2}.Math.\left(p^V.Math.C_{21}+C_{23}+p^V.Math.C_{12}+C_{32}+2.Math.C_{22}\right)& \left(38\right)\end{array}$

[0194] This function is used to calculate the sensitivity of the subsequent iteration in the outer-loop.

EXAMPLES

Example 1

Cubic Base Cell with Spherical Shape Void

[0195] The geometry of the base cell for this example 3D auxetic metamaterial is formed by creating a hollow spherical cavity inside a cube, as shown in FIG. 1A(A) and FIG. 1B(A). Each of the building cells was repeated to form a 3D cellular material as respectively shown in FIG. 1A(B) and FIG. 1B(B). The experimental bulk metamaterial was constructed by repeating nine building cells along three normal directions and cut half of the both end-cells in each direction. Each of the specimens of the bulk 3D material were manufactured using 3D printing (ObjetConnex350) with a silicone-based rubber material (TangoPlus) and a supporting material.

[0196] According to the deformation pattern after buckling, the Representative Volume Element (RVE) contains four building cells as shown in FIG. 1A(C) and FIG. 1B(C). According to the ratio (R) of the diameter of the sphere to the length of the cube, two resultant geometry were established: [0197] (1) a face-centred cubic cell with 0<R<1 (FIG. 1A(A)), which is used as a comparative design for the present invention; and [0198] (2) an inventive cubic cell 1<R<2 (FIG. 1B(A)) which comprises a metamaterial in accordance with an embodiment of the present invention. The porosity of this unit cell was found to be in the range of 0.69 to 0.97.

[0199] The material properties of the printed TangoPlus material were measured through standard compression test with six printed cylinders, up to the true strain of =0.70. Each of the 3D materials and their responses to strain and compression were also modelled as a linear elastic model using finite element analysis. A comparison of the deformation patterns between the experimental (A) and model (B) is provided in FIG. 4. Comparative Force-Displacement curves of the experimental (A) and model (B) are shown in FIGS. 5A and 5B.

[0200] The results indicate that the constitutive behaviour of each of the comparative face-centred cubic cell and inventive cubic cell could be accurately represented by a linear elastic model. It should be noted that the printed TangoPlus material exhibited slightly anisotropic behaviour with the Young's modulus along printing direction, 0.9250.02 MPa, being slightly lower than its lateral direction, 1.050.03 MPa. The Poisson's ratio of the face-centred cubic cell was found to be +0.47.

[0201] The performance of the inventive 3D cubic metamaterial was tested using standard compression tests similar to those commonly used for other cellular materials. To obtain a reliable homogenized material properties, the dimensions of the test specimens were selected as heightwidthdepth=100.0100.0100.0 mm. This resulted in a material built from a matrix having eight building cells in each normal direction as shown in FIG. 1B(B).

[0202] Two samples of the inventive cubic cell material are shown in FIG. 2. The left sample (A) still includes supporting material for the 3D printing. The right sample (B) has the supporting material removed. In spite of extreme care being taken during the removal of the supporting material, a few of the thinnest links in the bulk material were broken. An epoxy adhesive was used to repair that damage.

[0203] Comparative compression tests between the (1) comparative face-centred cubic cell and (2) the inventive cubic cell. The compression tests were conducted at a fixed strain rate of 10.sup.3 s.sup.1 using a Shimazu machine. Two cameras were used to capture the deformation in two lateral directions so as to determine the evolution of the Poisson's ratio of the metamaterial. The end strain were fixed at a nominal strain up to 0.3 for specimen formed from the comparative face-centred cubic building cells and 0.5 for specimens with inventive cubic building cells to avoid potential damage of the specimens. It was found that within these strain ranges, the deformation was purely elastic and totally reversible.

[0204] The bulk material composed of the comparative face-centred cubic building cells only exhibited global buckling at a very large strain of 0.25 as shown in FIG. 3A(a). Furthermore, as shown in FIG. 3B the stress-strain curve is linear before the buckling occurs. No obvious auxetic behaviour was observed. The Applicant notes that the localised buckling modes with alternating ellipsoids similar to 2D NPR materials reported in for example Overvelde et al (Compaction Through Buckling in 2D Periodic, Soft and Porous Structures: Effect of Pore Shape. Advanced Materials. 2012; 24:2337-2342) did not occur for this type of material.

[0205] The bulk material composed of the inventive cubic building cells, showed localised buckling modes with alternating ellipsoids. This material therefore deformed with clearly observable auxetic behaviour as shown in FIG. 3A(b). Furthermore, the force-deflection response of inventive metamaterial (shown in FIG. 3C) in two different directions D1 and D2 also shows observable auxetic response.

[0206] The different buckling behaviour of the materials formed from the face-centred cubic building cell and inventive cubic cell indicates that there is a critical porosity or volume fraction for the desired buckling mode. In this respect, auxetic behaviour is not possible when the porosity of the material is below 0.60, for example the face-centred cubic building cell material. The Applicants have unexpectantly found that a porosity of at least 0.6, preferably between 0.6 and 0.9 is necessary for the 3D material to display auxetic behaviour.

Example 2

Mechanism Analysis (Buckling Mode)

[0207] Numerical simulations were carried out using the commercial finite element (FE) software package ABAQUS (Simulia, Providence, R.I.) to determine the mechanisms of the auxetic behaviour observed in the inventive metamaterial discussed in Example 1.

[0208] The ABAQUS/standard solver was employed for buckling analysis and ABAQUS/explicit solver was employed for postbuckling analyses. Quadratic solid elements with secondary accuracy (element type C3D10R with a mesh sweeping seed size of 0.4 mm) were used. The analyses were performed under uniaxial compression. The buckling mode with 3D alternating ellipsoidal pattern from buckling analysis was used as the shape change or imperfection factor for non-linear (large deformation) post-buckling analysis. The finite element models were validated using experimental results.

[0209] FIG. 4 shows the comparison of deformation process of the metamaterial from numerical simulation and experimental result from one direction. Both the experimental results (A) and modelled behaviour (B) exhibit the auxetic behaviour in a similar manner. A noticeable difference is the long axis of ellipsoid of the representative volume unit (marked with dots) in the centre of the specimen. The similarity remained in the other lateral direction. According to the linear buckling analysis, these two different deformation patterns have nearly identical eigenvalues. Based on this analysis, the inventors consider that the actual deformation pattern after buckling is determined by the imperfection of the initial geometry.

[0210] It was found that the buckling mode was influenced by the boundary conditions of the FE model. Two boundary conditions were examined. One constrains all freedoms of the nodes on top and bottom surface except for the freedom on loading direction on the top surface and the other constrains only the freedom of the nodes bottom surface along loading direction. For the former boundary condition, the first buckling mode from the numerical simulation exhibited local buckling with alternating ellipsoids. For the latter boundary condition, the first buckling mode exhibited a planar pattern which was similar to the deformation patterns observed previously by Willshaw and Mullin (Soft Matter. 2012, 8, 1747). The 3D buckling pattern occurred as the fifth buckling mode.

[0211] It can also be observed from FIG. 4 that the deformation of the specimen is uniform at the early stage of the compression test before buckling occurs. The material behaves like conventional material with a positive Poisson's ratio. Only after the buckling occurred, the auxetic behaviour become evident, which indicated the value of negative Poisson's ratio is changing during the deformation process. This could be a disadvantage for applications with a required negative Poisson's ratio.

Example 3

Cubic Base Cell with Ovoid Shaped Void

[0212] To overcome the buckling disadvantages of Example 1 and 2, the geometry of the base cell for this example 3D auxetic metamaterial is formed by creating a hollow ovoid cavity inside a cube, as shown in FIG. 6. The designed ovoid comprised an 8% imperfection in the shape of the spherical void used in the material discussed in Examples 1 and 2. In addition, the matrix of base units in the material was arranged such that the central length axis of the ovoid void of each base unit was perpendicular to the central length axis of the ovoid void of each adjoining base unit. This, in effect, introduced the pattern of the buckling mode seen in Examples 1 and 2 into the void pattern of this embodiment of the metamaterial. The porosity of this unit cell was found to be 87.4% for Example 1 and 87.2% for Example 2.

[0213] A direct comparison of nominal stress-strain curves between experimental and numerical results is shown in FIG. 5. Both curves exhibit a similar trend and the corresponding stress are at similar levels. This demonstrates general agreement between experimental results and the finite element model. It is noted that the lower stress level in the experimental results could be attributed to the broken links during the removal process of the supporting material. The stress-strain curves of the inventive metamaterial are similar to the other cellular materials undergoing plastic deformation, with the difference that deformation of the inventive metamaterials is purely elastic and fully reversible. This appears to be an attribute of the properties of the base material used.

[0214] The overall deformation patterns for the proposed metamaterial with 8% imperfection (of the spherical shape of the voids) are shown in FIG. 6, and it is clear that the auxetic behaviour starts from the very beginning.

[0215] The Applicant observes that if the magnitude of the imperfection in the spherical shape of the void is increased (and thus the shape of the ovoid void altered or flattened), the Poisson's ratio of material could be altered, and thus effectively tailored to a desired value. This would produce a series of inventive cubic 3D metamaterials with prescribed initial negative Poisson's ratio value. This approach provides a fundamentally new way for generating a serial of 3D materials with a desired initial value of negative Poisson's ratio.

[0216] It should be noted that the volume fraction for the base cell and representative volume element of the inventive metamaterial varies with different imperfection magnitude. A combination of this approach with the initial geometry design can therefore be considered to design metamaterials with a desired volume fraction.

[0217] The above shows that the configuration of the base unit, void geometry and pattern of the matrix formed from the base units can be tailored using a buckling mode obtained through finite element analysis. The introduction of the buckling pattern into the matrix of the material and varying the magnitude of the imperfection in the spherical shape of the void enables so that it provides a mean to tailor the initial value of Poisson's ratio in a range from 0 to 0.5.

Example 4

Cubic Base Cell with Tetrahedron or Ellipsoid Shaped Void

[0218] A metamaterial of the present invention can also be formed using a cubic base cell with other void shapes, such as tetrahedron, or ellipsoid.

[0219] FIG. 7A provides the geometric configurations for a three dimensional structure porous metamaterial with tetrahedron in cube building cell. FIG. 7B provides the geometric configurations for a three dimensional structure porous metamaterial with ellipsoid in cube building cell. The geometry of the base cell for this example 3D auxetic metamaterial is formed by creating a hollow tetrahedron or ellipsoid cavity inside a cube, as shown in FIG. 7A(A) and FIG. 7B(A). Each of the building cells was repeated to form a 3D cellular material as respectively shown in FIG. 7A(B) and FIG. 7B(B). FIG. 7A(C) and FIG. 7B(C) illustrate a representative volume unit of the inventive metamaterial.

[0220] The porosity of this type of unit cell was found to be 0.63 in FIG. 7A (in the range of 0.5 to 0.91) and 0.69 in FIG. 7B (in the range of 0.6 to 0.97).

[0221] Tests have shown that this material has similar deformation behaviour as with previous cubic base cell with spherical voids. FIG. 8A provides the deformation pattern of the metamaterial shown in FIG. 7A under load. FIG. 8B provides the deformation pattern of the metamaterial shown in FIG. 7B under load. The deformation pattern shown in FIGS. 8A and 8B illustrates similar behaviour with previous cubic base cell with spherical voids examples.

Example 5

Metamaterial with Negative Linear Compression (NLC) Under Uniformed Pressure

[0222] A negative compression (NC) metamaterial of the present invention can be formed using a frame work similar to the topology resulting from bi-directional evolutionary structural optimization (BESO).

[0223] FIG. 9 provides the geometric configurations for the resulting three dimensional structure porous NC metamaterial with NLC. FIG. 9(A) provides the topology obtained from BESO. The geometry of the building cell for this example 3D NC metamaterial is formed by simplifying the irregular members in FIG. 9(A) to a truss with varied cross-section maximized at the middle span. The simplified building cell is shown in FIG. 9(B). Each of the building cells was repeated to form a 3D cellular material as respectively shown in FIG. 9(C). The porosity of this unit cell was found to be 0.902.

[0224] As shown in FIG. 9(A), the basic form of the NLC metamaterial is derived from BESO calculations, as discussed previously in relation to NLC optimisation.

[0225] In these calculations for this and subsequent examples, the finite element analysis is conducted by using ABAQUS version 10.1. Due to symmetry in three directions for orthotropic materials, only one eighth of the unit cell needs to be modelled. The one eighth model is divided into a mesh of 303030 brick elements (element type: C3D8). The resulting topology is smoothened based on curve and surface fitting. The target volume V is 30%. The unit for the compressibility is Pa.sup.1.

[0226] For the BESO calculations, the base materials are E.sub.b1=10.sup.15 (void) and E.sub.b2=1 (solid). To represent an impressible base material (such as silicon rubber) the Poisson's ratio v.sub.b is assumed to be 0.49. It is noted that for a given target volume V (which is the same as the volume fraction as the total volume of the unit cell is 1), the maximum achievable stiffness along a single axis is E.sub.max=VE.sub.b2. The stiffness target is then specified as E*=aE.sub.max=aVE.sub.b2 where a is the prescribed stiffness ratio. The stiffness ratio a is equal to 0.10. Therefore E*=aVE.sub.b2=0.100.31=0.030.

[0227] The result is shown in FIG. 9, with the unit cell having a topology similar to a truss-like system. The linear compressibility .sub.L3 is 17.53. Furthermore:

[00036]$\begin{array}{cc}E.Math..Math.\mathrm{matrix}& C.Math..Math.\mathrm{matrix}\\ \left[\begin{array}{ccc}0.068& 0.008& 0.063\\ 0.008& 0.068& 0.063\\ 0.063& 0.063& 0.141\end{array}\right]& \left[\begin{array}{ccc}33.302& 16.538& -22.171\\ 16.538& 33.302& -22.171\\ -22.171& -22.171& 26.807\end{array}\right]\end{array}$${}_{L.Math..Math.1}={}_{L.Math..Math.2}=27.67,{}_{L.Math..Math.3}=-17.53,E_1=E_2=0.030,E_3=0.037,.Math.v_{12}=-0.497,v_{23}=v_{13}=0.666$

[0228] The procedure has designed an optimised shaped void comprising a regular but complex shape, providing a cutout aperture in the truss structure, and an open end.

[0229] To verify the above material properties, a numerical simulation of a stress test was conducted on a model constructed from 888 unit cells of the above topology. The model is resized to 100 mm100 mm100 mm and was meshed with 7424 quadratic tetrahedral elements (ABAQUS element type C3D10I). A hydrostatic pressure P=1.4410.sup.3 is applied through rigid plates attached to the six faces. Displacements at the rigid plates are extracted and then strains are calculated, which result in .sub.1=.sub.2=41.3910.sup.3 and .sub.3=24.6910.sup.3, respectively. Normalizing these strains by the pressure P gives the following values of linear compressibility: .sub.L1=.sub.L2=28.74 and .sub.L3=17.15, which are very close to calculated values, with differences being less than 4%. The discrepancies are attributed to the different finite element models used for the unit cell and the array of 888 cells.

[0230] Tests have shown that this metamaterial expand in one direction while shrinking in the other two directions under pressure. In these tests, the NLC design shown in FIG. 9 was obtained by assuming a nearly incompressible base material with a Poisson's ratio of 0.49. A bulk material model is constructed from 888 unit cells of this topology. The model is resized to 100 mm100 mm100 mm. Using a 3D printer (Object Connex350), a prototype model was fabricated with a silicone-based rubber material (TangoPlus) and a support material in lower density. After the support material has been carefully removed, the NLC material design as shown in FIG. 10(A) was obtained. In the following discussions, X, Y and Z directions correspond to axes 1, 2 and 3, respectively. The material properties of the TangoPlus material are measured through standard compression test on three printed cylindrical samples, with the true strain up to 0.70. The results indicate that the constitutive behaviour of the base material can be accurately represented by a linear elastic model. It is found that the Young's modulus is 1.05 MPa and the Poisson's ratio is 0.48. These values are used in the FE simulations described below.

Uniaxial Compression Test

[0231] Uniaxial compression tests were conducted in the X, Y or Z direction separately. From these experiments, the effective (average) compliance matrix for the bulk material can be obtained. Also, linear elastic finite element analyses of the bulk material model with 888 cells are performed by applying unidirectional pressures through two rigid plates. From the FE results, the effective compliance matrix of the material can be calculated as well. The compliance matrix for the unit cell and the effective compliance matrices for the model with 888 cells from both experimental and FE results (all normalized with respect to Young's moduli) are given in Table 1.

TABLE-US-00001 TABLE 1 C Matrix Unit cell (FIG. 2) [00037]$\left[\begin{array}{ccc}33.302& 16.538& -22.171\\ 16.538& 33.302& -22.171\\ -22.171& -22.171& 26.807\end{array}\right]$ 8 8 8 cells, FE results [00038]$\left[\begin{array}{ccc}36.843& 15.431& -23.282\\ 15.431& 36.835& -23.281\\ -23.217& -23.175& 29.632\end{array}\right]$ 8 8 8 cells, experimental results [00039]$\left[\begin{array}{ccc}39.68& 13.17& -21.03\\ 13.17& 39.68& -21.03\\ -22.79& -22.79& 31.43\end{array}\right]$

[0232] It is seen that the FE results agree reasonably well with the experimental data.

[0233] It is noted that the C matrix for the unit cell and that of the 888 cells (FE results) are also similar. The discrepancies are mainly attributed to different boundary conditions. For the unit cell, periodic boundary conditions are applied; while for the bulk material model, all nodes on top and bottom surfaces are only allowed to move in the loading direction.

Triaxial Compression Test

[0234] In order to examine the behaviour of the NLC design under uniform pressure, a tri-axial pressure test was performed using a standard tri-axial test machine commonly used for soil testing. Firstly the prototype was put inside a sealed plastic bag, to which a plastic tube of 2 mm in diameter is connected. During the application of the uniform pressure on the outside surface of the plastic bag, the air inside the bulk material is pressed out through the plastic tube. The uniform pressure is gradually increased from 0 to 5 kPa. The final deformed shape of the material at 5 kPa is given in FIG. 10(A). It is seen that the original cube has become narrower and taller under uniform pressurea clear sign of the NLC effect.

[0235] A finite element simulation of the tri-axial test has been conducted. To capture the large deformation observed in the experiment, a nonlinear finite element analysis considering large deformation is carried out. The base material is assumed linear elastic, with Young's modulus E.sub.0=1.05 MPa and .sub.0=0.48. The plastic bag is modelled using membrane elements, with thickness t=0.2 mm, Young's modulus E.sub.m=6 MPa and .sub.m=0.48. The deformed shape of the model at 5 kPa from the FE simulation is given in FIG. 10(B), which is very similar to the experimental result shown in FIG. 10(A). Furthermore, the average strains in X and Z directions from the experiment and the FE simulation are given in FIG. 10(C). The experimental data agree reasonably well with the FE results.

Example 6

Metamaterial with Negative Area Compression (NAC) Under Uniformed Pressure

[0236] FIG. 11 provides the geometric configurations for a three dimensional structure porous negative compression metamaterial with NAC. FIGS. 11(Ahalf cell) and 11(Bfull unit cell) provides the topology obtained from bi-directional evolutionary structural optimization (BESO) as discussed previously in relation to NAC optimisation calculations.

[0237] The porosity of this unit cell was found to be 0.696.

[0238] In the BESO calculations, the base materials were assumed to have E.sub.b1=10.sup.15 and a common Poisson's ratio and a stiffness ratio a=0.05. The calculated parameters were:

[00040]$\begin{array}{cc}E.Math..Math.\mathrm{matrix}& C.Math..Math.\mathrm{matrix}\\ \left[\begin{array}{ccc}0.038& 0.060& 0.060\\ 0.060& 0.227& 0.082\\ 0.082& 0.082& 0.227\end{array}\right]& \left[\begin{array}{ccc}67.248& -13.121& -13.121\\ -13.121& 7.610& 0.744\\ -13.121& 0.744& 7.610\end{array}\right]\end{array}$${}_{A.Math..Math.23}=-9.534,{}_{L.Math..Math.2}={}_{L.Math..Math.3}=-4.767,E_1=0.015,.Math.E_2=E_3=0.131,v_{12}=v_{13}=0.195,v_{23}=-0.098$

[0239] The resulting topology is shown in FIGS. 11(A) and 11(B). The topology is symmetrical with respect to the 45 degree plane perpendicular to plane 2-3. The stiffness in axis 1 is E.sub.1=0.015 which is the same as that of the NLC design in Example 5 because both designs have the same stiffness constraint. The compressibility of the NLC design .sub.L3 is equal to 44.21 as compared to here .sub.L3=.sub.L2=4.767. It is noted that the absolute value of the compressibility of the NAC material is significantly smaller than that of the NLC material discussed in the previous example, even though they have the same stiffness in one direction. The procedure has designed an optimised shaped void comprising multiple voids within the cubic base unit forming a complex shape. In this respect, the optimised shape void includes two internal voids and three external voids forming the topology of the base building cell.

[0240] The geometry of the building cell for this example 3D NC metamaterial is formed by simplifying the irregular members in FIG. 11(B) to a truss with varied cross-section maximized at the middle span. The simplified building cell is shown in FIG. 11(C). Each of the building cells was repeated to form a 3D cellular material as respectively shown in FIG. 11(D).

Example 7

Metamaterial with Zero Linear Compressibility (ZLC) Under Uniformed Pressure

[0241] FIG. 12 provides the geometric configurations for a three dimensional structure porous negative compressibility metamaterial with ZLC. FIG. 12(Ahalf cell) and FIG. 12(Bfull cell) provides the topology obtained from bi-directional evolutionary structural optimization (BESO).

[0242] In the BESO calculations, the procedure provided in Example 5 was used to determine a ZLC design. The material was designed with constraint on the linear compressibility, i.e. .sub.L3=0, is shown in FIG. 7 and a strain energy of 6.33. The linear compressibility .sub.L3 equals 0.002, which is very close to zero. The procedure has designed an optimised shaped void comprising multiple voids within the cubic base unit forming a complex shape. In this respect, the optimised shape void includes two internal voids and at least three external voids (sides) forming the topology of the base building cell.

[0243] The calculated parameters were:

[00041]$\begin{array}{cc}E.Math..Math.\mathrm{matrix}& .Math.C.Math..Math.\mathrm{matrix}\\ \left[\begin{array}{ccc}0.124& 0.034& 0.079\\ 0.034& 0.126& 0.079\\ 0.079& 0.079& 0.173\end{array}\right]& \left[\begin{array}{ccc}11.438& 0.228& -5.334\\ 0.228& 11.438& -5.334\\ -5.334& -5.334& 10.665\end{array}\right]\end{array}$${}_{L.Math..Math.3}=-0.002,E_1=E_2=0.087,E_3=0.094,v_{12}=-0.020,.Math.v_{23}=v_{13}=0.466$

[0244] The geometry of the building cell for this example 3D NC metamaterial is formed by simplifying the irregular members in FIG. 12(B) to a truss with varied cross-section maximized at the middle span. The simplified building cell is shown in FIG. 12(C). Each of the building cells was repeated to form a 3D cellular material as respectively shown in FIG. 12(D).

[0245] The porosity of this unit cell was found to be 0.854.

Example 8

Metamaterial with Zero Area Compressibility (ZAC) Under Uniformed Pressure

[0246] FIG. 13 provides the geometric configurations for a 3D structure porous negative compressibility metamaterial with ZAC. FIGS. 13(Ahalf cell) and 13(Bfull unit cell) provides the topology obtained from BESO.

[0247] In the BESO calculations, the material is designed to the ZAC criterion following a procedure corresponding to the ZLC example (Example 7). The result is shown in FIGS. 13(A) and (B). The calculated parameters were:

[00042]$\begin{array}{cc}E.Math..Math.\mathrm{matrix}& C.Math..Math.\mathrm{matrix}\\ \left[\begin{array}{ccc}0.071& 0.071& 0.071\\ 0.071& 0.173& 0.090\\ 0.071& 0.090& 0.173\end{array}\right]& \left[\begin{array}{ccc}30.681& -8.335& -8.335\\ -8.335& 14.935& -1.905\\ -8.335& -1.905& 10.238\end{array}\right]\end{array}$${}_{A.Math..Math.23}=-0.002,{}_{L.Math..Math.2}={}_{L.Math..Math.3}=-0.001,.Math.E_1=0.033,E_2=E_3=0.098,v_{12}=v_{13}=0.272,v_{23}=0.186$

[0248] The strain energy is 7.00, which is higher than that of ZLC (6.33). This is because of the additional constraint on .sub.L2 compared to the ZLC design. The area compressibility .sub.A23 is equal to 0.002, which is negligibly small (in terms of its absolute value) compared to that of the NAC design shown in FIG. 11 (Example 6) (25.40). The procedure has designed an optimised shaped void comprising multiple voids within the cubic base unit forming a complex shape. In this respect, the optimised shape void includes an internal void and at least two external voids (sides) forming the topology of the base building cell.

[0249] The porosity of this unit cell was found to be 0.893.

[0250] The geometry of the building cell for this example 3D NC metamaterial is formed by simplifying the irregular members in FIG. 13(B) to a truss with varied cross-section maximized at the middle span. The simplified building cell is shown in FIG. 13(C). Each of the building cells was repeated to form a 3D cellular material as respectively shown in FIG. 12(D).

[0251] There are several special features about the performance of the metamaterial of the present invention: [0252] The deformation of embodiments of the inventive metamaterials are purely elastic and fully reversible like other elastomers, but the stress-strain curves of our NPR metamaterials exhibits plateau feature as other cellular material undergo plastic deformation. The negative Poisson's ratio of the inventive metamaterials is retained over a wide range of applied strain and the range can be altered by the initial volume fraction and the magnitude of imperfection. [0253] The proposed design approach can be applied to any length-scale. It can be extended for tuning other properties of a material from the smallest scale. [0254] The inventive metamaterial can also be combined with stimuli responsive material to switch between different deformation patterns.

[0255] The material of the present invention can be used to fabricate sensors, actuators, prosthetics, surgical implants, anchors, (as for sutures, tendons, ligaments, or muscle), fasteners, seals, corks, filters, sieves, shock absorbers, impact-mitigating materials, hybrids, or structures, impact absorption or cushioning materials, hybrids, or structures, wave propagation control materials, hybrids, or structures, blast-resistant materials, hybrids, or structures, micro-electro-mechanical systems (MEMS) components, and/or stents.

[0256] Applications of this invention directed at the biomedical field include uses relating to prosthetic materials, surgical implants, and anchors for sutures and tendons, endoscopy, and stents.

[0257] Applications of this invention directed the mechanical/electrical field include uses in piezoelectric sensors and actuators, armours, cushioning, and impact and blast resistant materials, as deployable material and defence materials for infrastructures, the filter and sieve field, the fastener field, the sealing and cork fields, and the field of micro-electro-mechanical systems (MEMS).

[0258] In one exemplary embodiment, the inventive metamaterial can be formed as a compressible biocompatible polymer for use in intervertebral disc replacement. In some forms, the configuration and patterning the voids can be configured to allow the flow of fluid. The fluid can be used as a dampening mechanism within the material.

[0259] An immediate application of NLC/NAC metamaterials is the optical component in interferometric pressure sensors due to the higher sensitivity achieved by a combination of large volume compressibility with negative linear compressibility.

[0260] One significant application of the NC metamaterials is to be used as inserted foam for the OA treatment surgery using a NPWT system. The NC metamaterial will maintain their height but contract laterally under negative pressure and thereby enable the OA wound to close directly without using invasive mechanical devices.

[0261] With further understanding the mechanisms of negative compressibility, NLC/NAC materials also have potential to be used as efficient biological structures, nanofluidic actuators or as compensators for undesirable moisture-induced swelling of concrete/clay-based engineering materials (Cairns et al., 2013).

[0262] In exemplary embodiments, the inventive metamaterials can be used in a new type of smart amour for defence engineering or in blast control from explosive devices and projectiles. In one embodiment, the inventive material is formed from a Titanium or titanium alloy base unit matrix. The material can be used to compresses to the point of impact thereby providing lightweight armour plating.

[0263] In yet another exemplary application, the material can be used as lightweight cellular materials with enhanced energy absorption for motor vehicles.

[0264] Those skilled in the art will appreciate that the invention described herein is susceptible to variations and modifications other than those specifically described. It is understood that the invention includes all such variations and modifications which fall within the spirit and scope of the present invention.

[0265] Where the terms comprise, comprises, comprised or comprising are used in this specification (including the claims) they are to be interpreted as specifying the presence of the stated features, integers, steps or components, but not precluding the presence of one or more other feature, integer, step, component or group thereof.