A Refined Finite Element Formulation for the Microstructure-Dependent Analysis of Two-Dimensional (2D) Lattice Materials.

A finite element approximation is proposed for the dynamic analysis of two-dimensional (2D) lattice materials. The unit cell is modeled by means of a defined number of shear deformable micro-beams. The main innovative feature concerns the presence of a microstructure-dependent scale length, which allows the consideration of the so called size-effect that can be highly relevant, due to the characteristics of the lattice at the local scale. Some numerical results show the influence of the microstructure parameter on the dynamic behavior of two-dimensional lattice materials.

1. Introduction

At the local scale, many periodic materials can be seen as made of interconnected beams. Such materials, also called lattice materials [1,2], have received a wide interest in the scientific and industrial community over the last years, mainly due to their high strength to weight and stiffness to weight ratios, which, of course, are an important advantage from a static point of view. Lattice materials, in fact, can be designed to be stretching-dominated, thus providing higher structural performances. Many scientific papers dealing with the characterization of their effective elastic properties are available in the literature [3,4,5,6,7]. In comparison, only a few papers focus on their dynamic behavior, which is of an identical practical interest in the engineering field. Innovative applications are, in fact, related to many aspects concerning the elastic wave propagation mechanism within them. A small number of recently published works [8,9,10,11,12] has pioneered the study of the wave propagation phenomena in two-dimensional (2D) periodic lattice materials, with the main scope of detecting the existence of phononic band gaps, i.e. frequency ranges, within which the propagation of elastic/acoustic waves is prevented regardless of the incident wave direction. These studies look at a new generation of devices for energy absorption, noise and vibration control. Moreover, very interesting results have also been established within the field of civil engineering, where a new seismic mitigation strategy has been proposed that exhibits a close analogy with the idea of a lattice material [13].

In the context of lattice materials, the periodic topology plays a relevant role, considering the strict relation with the bending mechanism, which can be significantly dominant at the micro scale, where both the displacements and rotations can be similarly important [14]. This, in addition to experimental evidence reported by many authors [15,16,17,18,19], suggests an overcoming of the classical continuum theory, which assumes only displacements as kinematic quantities.

According to the classical continuum theory, stress-transfer locally occurs only through a force vector per unit area. As a consequence, stresses and strains can be represented by symmetric second-order real tensors. Nevertheless, for the accuracy of any prediction given via the classical theory, the condition that changes in stresses and strains (with wavelength λ) can be considered as uniform at the local scale is pivotal: D > λ ≫ l, with D denoting the global (structural) size, while l is a characteristic length representative of the microstructure. On the contrary, when dealing with a number of situations [20,21,22,23,24] (thin films, adhesive interfaces, notches, crack tips, localized deformations, boundary layers), it can be observed that the scale lengths of the problem agree with a new different condition (D > λ > l), rather than the previous one (D > λ ≫ l). More precisely, the lower the ratio λ/l > 1.0, the greater the importance of considering non-uniform stresses and strains at the local scale.

Unfortunately, due to the lack of an internal characteristic length, classical models are unable to capture the microstructure-dependent size-effect and, therefore, need to be extended by using higher order, non-local continuum theories.

Many continuum theories have been introduced over time in order to account for the size-effect [25,26,27,28,29,30,31,32]. Recently, non-local theories for the Bernoulli–Euler [33] and Timoshenko [34] beams have been proposed according to a modified couple stress theory [35] by using the non-local constitutive equation introduced in [30].

In line with these last developments, the present work proposes the analysis of the microstructure-dependent size-effect on the dynamic behavior of 2D lattice materials.

Due to the very important conclusions outlined in [33,34], this influence is expected to be relevant also on the local resonance phenomena, which determine the existence of frequency band gaps within the low frequency region.

To this scope, a refined 1D finite element approximation, accounting for a micro-scale length parameter, as well as the shear deformability, has been proposed.

2. Frequency Gaps

Due to the complexity of the equations governing the wave propagation within a periodic material, there is no general rule that can be used to predict whether phononic band gaps exist before a full analysis has been performed.

The existence of phononic band gaps is commonly related to the destructive interference occurring in the multiple scattering and reflection of elastic waves, which propagate throughout a periodic material. It is well known that the position and width of the phononic band gaps depend on the lattice topology, as well as the elastic moduli of the material and the mass density. Nevertheless, the influence of the micro-scale constitutive parameters has yet to be analyzed in detail.

Recent experimental and numerical studies [36,37,38], however, show the existence of phononic band gaps in frequency ranges lower than those evaluable from the above cited destructive interference theory, revealing the link between these low frequency band gaps and the local resonance phenomena of the inner microstructures. Moreover, in a recent work [39], a modified 2D lattice material with square topology, obtained by connecting auxiliary cantilever beams to the primary micro-structure, was investigated. The frequency band diagrams, numerically identified, show as a result the relevant role played by the physical mechanism of local resonance.

Based on the previous considerations, it is the opinion of the authors that the study of the dynamic response of 2D lattice materials should require a more refined continuum approach, which accounts for the stiffness modification related to the size-effect. To this aim, a micro-scale length parameter has been proposed according to the hypotheses of the modified couple stress theory developed in [35].

3. Free Waves Propagation in 2D Periodic Lattice Materials

The study of the wave propagation within a periodic material generally involves the Bloch’s theorem, which allows one to reduce the analysis of an infinite lattice to that of a single unit cell. An example of a simple two-dimensional square lattice material and the corresponding unit cell is depicted in Figure 1, where e1 and e2 are the unit vectors along the directions of spatial periodicity. The square topology implies in this case that e1 and e2 are normal to each other.

Figure 1

Example of two-dimensional (2D) square lattice material.

In the example above (Figure 1), the microstructure is composed of interconnected beams arranged along the directions e1 and e2 only. In general, however, more complex microstructures can be investigated. For instance, in [39], a modified 2D lattice material has been considered, which exhibits a more elaborated unit cell.

According to the Bloch’s theorem, if a plane elastic wave with angular frequency ω propagates in the lattice material (Figure 1), the displacement of an arbitrary point P is given by:

In Equation (1), the symbol r indicates the position of the point P, k denotes the Bloch wave vector and u (r) is the amplitude, which is characterized by the same spatial periodicity as the point lattice. As a consequence, it is possible to express the position of any point P as a function of the position of the corresponding point O, which is located in a reference unit cell. Once the reference cell has been selected (see Figure 1), it results: r = r0 + n1e1 + n2e2, where (n1, n2) is an integer pair and r0 is the position vector of O. It is then possible to update Equation (1) as follows:

In brief, changes in wave amplitude from a generic cell to the adjacent one do not depend on the cell location within the point lattice. By virtue of this, it is possible to restrict the study of an infinite lattice system to the analysis of a unit cell.

It is easy to verify that the reciprocal unit vectors and assume the following form: where e3 is the unit vector normal to the plane of the problem (Figure 1).

If the Bloch wave vector is expressed in the reciprocal space as follows:

Then, from Equation (2), the following equation is obtained, which provides a periodic boundary condition for the dynamic analysis of the unit cell:

4. Finite Element Approximation

The analysis of the free wave propagation in 2D lattice materials can be performed by modeling a unit cell only. Rigid connections are usually considered between the micro-beams of which the unit cell is composed. Examples of the finite element (FE) mesh over the unit cell are shown in Figure 2. In this figure, both a simple 2D square lattice materials (I) and a more elaborated one (II) are sketched, the last one being composed of auxiliary cantilever beams interconnected to the primary microstructure. It is worth noting that this second example is substantially the same as reported in [39].

Figure 2

Examples of the FE mesh over the unit cell of a 2D square lattice material (a being the lattice constant).

Regardless of how many micro-beams are interconnected to form the unit cell, the mesh here proposed is composed of finite elements characterized by six degrees of freedom (d.o.f.), including three d.o.f per each node p (p = node i or node j): the displacements , along the directions , and the rotation about the axis , where (, , ) are the local co-ordinates depicted in Figure 3, with aligned from node i to node j and = . It is important to highlight that, without losing generality, the first four nodes of the mesh are supposed to be located at the connections with the adjacent cells, while the fifth node is located at the center of the unit cell under consideration (Figure 2).

Figure 3

Finite element.

Let be the numeric vector collecting all degrees of freedom related to a generic finite element (Figure 3):

Accounting for a general orientation of the finite element, it is useful to express the numeric vector , which refers to the local co-ordinates system (, , ) as a function of the nodal displacements referred, instead, to the global co-ordinates system (, , = ): with = and = . It results: being the angle from to

By virtue of the next position: where η indicates the following dimensionless quantity accounting for the influence of the shear deformations the displacements field , expressed in the local co-ordinates, can be interpolated over the generic finite element as follows:

In Equation (10), the symbols E, G, I and indicate the Young modulus of the bulk material, the shear modulus, the inertia of the cross-section and the area of the cross-section incorporating the shear correction (). With reference to the well-known plane strain hypothesis, both I and are evaluated per unit length along . Furthermore, denotes the length of the finite element.

Moreover, in Equation (11), the symbols and denote the classical linear Lagrange polynomials, while the symbols , , and denote four enhanced shape functions, which allow the incorporation of the effect of shear deformations. It results: with indicating the normalized axial coordinate (, Figure 3) and being four appropriate interpolating cubic functions (p = 1,2), (q = 0,1):

The finite element proposed is able to account for shear strains and is locking-free. These features, as it is well-known, are essential for short transient and wave propagation analysis. A Bernoulli-Euler element, in fact, exhibiting an infinite phase velocity, because the parabolic order of the equation of motion, is useless for simulating wave propagation.

According to the non-local beam model presented in [35] (see Appendix), the generalized stresses within the generic finite element can be expressed as a function of the nodal displacements: where N, M, Y and V denote the generalized stresses of the beam model, C denotes the matrix given in Equation (A30) and B is a numeric matrix with the following entries: with reference to the global coordinates, the element stiffness matrix assumes the following form:

On the other hand, the kinetic energy can be expressed as follows: where the consistent mass matrix is computed over the generic finite element by means of:

By standard procedures, the equations of motion can be assembled in the following matrix form: where and are the overall mass and stiffness matrices of the unit cell and and denote the nodal displacements and external forces vectors in the global reference system (, , ), respectively, being the number of nodes over the unit cell. If a plane elastic wave with angular frequency is considered, Equation (27) assumes the new form:

The last form of the equations of motion (Equation (28)) must be coupled with the periodic boundary conditions given by Equation (6): where , , and are the nodal displacements relative to node 1 to 4, indicated in Figure 2 (connections with the four adjacent cells).

By means of Equations (29) and (30), the number of degrees of freedom is equal to . The global displacements vector , in fact, can be expressed as a function of a reduced displacements vector , according to the following relationship: where the matrix works as a transfer operator: with:

Finally, the equations of motion of the unit cell can be expressed in the following form:

The condition (free wave motion) implies that the study of the wave transmission within the periodic square lattice material reduces to the eigenvalue problem of Equation (35). For fixed values of and , to be expressed in the reciprocal space (, ) according to Equation (5), the frequencies of the wave propagation coincide with the eigenvalues of the problem formulated in Equation (35). Varying and along the boundary of the irreducible part of the first Brillouin zone, the band structure of the lattice material can be detected. In Figure 4, the filled region (triangle OAB) indicates the irreducible part of the first Brillouin zone concerning the square topology under consideration.

Figure 4

Irreducible part of the first Brillouin zone.

5. Numerical Results

In order to elucidate the investigation approach, the finite element model presented in the previous section has been applied for studying the dynamic behavior of a 2D square lattice material. More in detail, a parametric analysis has been carried out by varying the role of the microstructure length l. The example under consideration concerns the scheme (I) shown in Figure 2, the geometry and the mechanical parameters being fixed, as indicated in Table 1.

Table 1

Geometry and mechanical parameters.

lattice constant:	a=1.0 mm
cross-section area per unit length along i3	A=1.000×10^−1 mm^2/mm
cross-section shear area per unit length along i3	A_s=8.333×10^−2 mm^2/mm
flexural inertia per unit length along i3	I=8.333×10^−5 mm^4/mm
Young modulus	E=2.100×10^5 N/mm^2
shear modulus	G=8.077×10^4 N/mm^2
Poisson ratio	ν=0.3
mass density	ρ=7.850×10^-6 Kg/mm^3

The mesh used for the numerical simulations has been obtained by dividing both the horizontal beam and the vertical one shown in Figure 2 (I) by 50 finite elements, each one. The overall number of finite element results equal to 100, while the overall number of nodes (N) is equal to 101. The final number of degrees of freedom (Ndof) is thereby equal to 297. Many tests have been carried out in order to check the convergence rate and to assess the accuracy of the numerical solutions.

The following values of the dimensionless microstructure parameter have been considered in the analysis: {0.00, 0.01, 0.05, 0.10}. It is important to remark that the first value () implies that the size-effect is not present.

The eigenvalue problem given by Equation (35) has been solved via a call to a dedicated routine belonging to the IMSL Math Library. Thus, the frequency values have been identified and represented in a dimensionless form by means of: where denotes the first bending resonance frequency of a pinned-pinned beam with the same properties given in Table 1, its length being equal to the lattice constant a:

In Figure 5, a typical diagram of versus the modal number n is presented, with . The figure refers to fixed values of k1 and k2 (i.e., a well-defined point on the boundary of the irreducible part of the first Brillouin zone).

Figure 5

Dimensionless frequencies versus n (fixed point on the boundary OABO).

As it is easy to argue, the size-effect is found to be strongly relevant within the high frequency region when dealing with λ > 0.01.

For what concerns the dynamic behavior within the low frequency region, Table 2 indicates also in this case a considerable influence of the parameter .

Table 2

Ten lowest frequency values (dimensionless).

	λ = 0.00			λ = 0.01			λ = 0.05			λ = 0.10
	O	A	B	O	A	B	O	A	B	O	A	B
1	0.00	0.15	1.79	0.00	0.15	1.80	0.00	0.19	1.85	0.00	0.29	1.89
2	0.00	1.83	1.81	0.00	1.83	1.81	0.00	1.88	1.88	0.00	1.91	1.94
3	4.21	4.04	4.35	4.27	4.09	4.40	5.40	5.18	5.41	7.17	6.91	7.06
4	4.21	4.50	4.35	4.27	4.56	4.41	5.40	5.62	5.43	7.17	7.34	7.11
5	6.07	6.17	6.27	6.16	6.26	6.37	7.88	7.99	8.12	10.74	10.87	11.02
6	8.60	8.70	8.80	8.70	8.81	8.91	10.53	10.65	10.75	12.94	13.06	13.17
7	15.62	15.32	15.33	15.77	15.46	15.48	18.10	17.76	17.77	20.66	20.30	20.27
8	15.62	15.64	15.34	15.77	15.78	15.48	18.10	18.10	17.77	20.66	20.63	20.27
9	20.49	20.60	20.72	20.72	20.83	20.95	24.62	24.71	24.82	29.29	29.31	29.31
10	22.97	23.08	23.18	23.18	23.28	23.38	26.45	26.54	26.61	30.14	30.12	30.06

With respect to the corresponding predictions obtained by assuming λ = 0.00, relevant variations emerge, being the maximum percentage increase (+93%) found at point A when λ = 0.10.

6. Conclusions

As highlighted in a number of papers, the dynamic behavior of a square lattice material is significantly affected by the underlying microstructure. Due to the spatial periodicity, the hypothesis of infinite lattice points allows us to reduce the analysis to the unit cell and to investigate the dynamic behavior by means of the Bloch theorem.

It has been found by many authors that the natural frequencies predicted via a non-local beam model are always higher than those evaluated by classical beam models, due to the increased bending stiffness related to the so-called size-effect.

A more refined beam model is therefore required.

For the scope of a reliable evaluation of the existence, position and width of frequency band gaps, in the opinion of the authors, the study of a 2D lattice materials should account for possible size-effects.

The paper provides a simple 1D finite element, developed in accordance with a simplified couple stress theory, which is able to simulate the size-effect via a unique microstructural parameter. Shear deformations are also accounted for.

The numerical results confirm the influence of the size-effect on the dynamic behavior of 2D periodic materials.