© M. Schaeffer and M. Ruzzene, published by EDP Sciences, 2015

1 Introduction

In previous studies [1, 2], the authors have shown that periodic magneto-elastic structures can be reconfigured and tuned to change the way that plane waves travel through them. Though it is evident that such reconfiguration and tuning also affects the mechanical properties of the structures, these changes and the quantification thereof have not yet been studied in detail.

Assembly of conventional materials in novel ways creates new structures that can perform better than their constituents. Such structures are often composed of “metamaterials” and are called “metastructures”. Periodic and cellular structures are common products of such an endeavor as they can be extremely stiff for their weight, or exhibit novel wave propagation properties [3–8]. Oftentimes, it is convenient to investigate periodic metamaterial assemblies as equivalent continuua through homogenization. The process of homogenization has been investigated for discrete periodic structures by, for example Born [9], Gonella and Ruzzene [10], Maradudin et al. [11], Suiker et al. [12]. Similar methods are applied to find the equivalent properties of the discrete models discussed in this publication.

Related work has been presented [13] in which finite magneto-elastic structures show that their Poisson’s ratio can be tuned with the application of a non-uniform external magnetic field and through changes in constituent dimensions. The magnetic dipole moments in the 2D structures discussed therein are oriented in the plane of the system and the changes in Poisson’s ratio are a result of the re-orientation of the system’s rigid entities under the constant application of the external magnetic field. Also related is the homogenization methodology in reference [10], which is well suited for calculating equivalent mechanical properties of periodic lattices and is applied to hexagonal and re-entrant beam lattices. It is also generalized for arbitrarily complex elastic structures, but is not presented for use in systems exhibiting non-nearest neighbor interactions, such as magnetic interactions, and must therefore be extended for application to magneto-elastic lattices. In contrast to the work by Grima et al. [13], the structures discussed in this publication contain magnetic dipole moments perpendicular to the lattice plane, and mechanical properties including stiffnesses are presented. The changes Poisson’s ratio in this publication are not brought about solely by geometric orientation as in reference [13], but also through changes in interaction stiffnesses, which are caused by the included magnets. Furthermore, this publication focus specifically on the multistability and bistability of the studied lattices and how reconfiguration effects properties.

The objective of this paper is to extend the homogenization method presented in reference [10] to structures with non-nearest neighbor interactions and quantify the tunability of the equivalent mechanical properties of reconfigurable periodic magneto-elastic lattices such as those studied in reference [1]. This includes stiffness and density for 1D lattices and Poisson’s ratios and stiffnesses for 2D lattices.

The paper is organized as follows. Following this introduction the modeling approach used is discussed in Section 2, which contains the methodology for the study of periodic lattices, including the identification of equilibrium configurations (Sect. 2.4). Then the homogenization procedure applied is described in Section 3. Finally the homogenization method is applied to magneto-elastic structures to obtain solutions for 1D and 2D structures in Sections 4 and 5 respectively. Section 4 includes a comparison to numerical simulations for validation of the homogenization method. Conclusions are then summarized in Section 6.

2 Theoretical background

2.1 Lattice model

Magneto-elastic lattices are modeled as systems of permanent magnetic particles with translational and rotational degrees of freedom (DOFs). Axial and torsional springs connect the particles to form an elastic structure. The particles have a finite radius r_p, and feature both translational and rotational inertias. All particles are identical in terms of their inertial and magnetic properties. Similarly, all springs connecting the particles are massless and feature the same axial and torsional stiffness. The interactions between particles are governed by: (1) magnetic interactions, and (2) elastic interactions through axial and torsional springs. The energy functionals associated with the interactions above and the governing equations are described in the following subsections. Details can be found in reference [1].

2.2 Energy functionals

The behavior of the generic i-th particle of a 2D lattice is described by the position of its center of mass r_i = x_ii + y_ij, and by a rotation angle θ_i. Accordingly, the vector of the generalized coordinates associated with the i-th particle is given by:(1)

2.2.1 Mechanical interactions

The strain energy associated with the elastic interactions is defined in terms of axial and torsional contributions. A phenomenological, simplified model is formulated whereby the axial spring is characterized by a bilinear force-displacement relation, which is chosen as a simplified model for a contact condition between two particles. The contact stiffness is described as proportional to the stiffness k_a through a parameter α > 1. A contact distance r_c = 2r_p defines the contact conditions. When the contact conditions are satisfied the axial stiffness of the bilinear spring is αk_a. The strain energy associated with the mechanical components connecting particles i and j is therefore given by:(2)where k_τ denotes the torsional spring constant, r _ij =  r _j −  r _i, r₀ defines the inter-particle distance when the axial springs are undeformed, and H is the Heaviside unit step function.

2.2.2 Magnetic interactions

The magnetic interaction force is described by modeling the permanent magnets as magnetically rigid dipole moments [14]. Furthermore, all magnetic dipole moments are considered perpendicular to the direction of motion. The force between two dipoles is given by Yung et al. [15]. In the 1D and 2D lattices considered, the motion of the masses is constrained to the direction i and to the lattice plane i , j respectively. This simplifies the expression of the magnetic interaction force as follows:(3)where μ₀ is the permeability of free space, m _i and m _j are the dipole moments of particles i and j, n  =  r _ij/r_ij is a unit vector and ψ ^(m) = m_i m_j3μ₀/(4π) contains all the magnetic constants. Evaluating the work of the magnetic force (Eq. (3)) from a reference state at r_ij → ∞ to the current state r_ij defines a magnetic potential energy function .

The magnetic force described by equation (3) in theory acts for r → ∞. This implies that the lattice is, in general, characterized by non-nearest neighbor interactions that extend to infinity. However, the magnitude of the magnetic force rapidly decays with inter-particle distance, which for practical purposes allows the introduction of a radius of influence r_∞, beyond which magnetic interactions are neglected. As many as 3,000 particles may be included in a typical analysis of a magnetized structure.

2.3 Lattice description and unit cell energy

The considered lattices are characterized by periodic distributions of masses along the i direction in the case of the 1D lattices, and over the plane i, j for the 2D lattices. Their positions are described through a set of lattice vectors di for the 1D lattices and d₁, d₂ for the 2D lattices.

A unit cell is here defined as the smallest structure which can be repeated to describe the mass, stiffness, and magnetic characteristics of the assembly. The total energy of the unit cell is given by the sum of the potential energy U and kinetic energy T and the work of the external forces W:(4)where e = 1, …, E identifies the particles belonging to a unit cell. The magnetic interactions can be considered as external forces or, since they are magnetically rigid, they can be included in U_e. As the lattice is periodic, this energy functional is equal for all unit cells and its expression does not need to keep track of a specific unit cell identifier, which simplifies notation. The generalized coordinates for the unit cell are therefore:(5)

The set of governing equations describing the behavior of the cell and its interactions with connected units can be derived from the energy functional through application of Lagrange’s equations. The resulting set of governing equations is:(6)where are unit cell level matrices, and are unit cell level vectors of generalized coordinates and magnetic forces.

2.4 Identification of equilibrium configurations

The energy functional in equation (4), restricted to the time constant terms, is employed for the identification of the periodic equilibrium configurations and for the analysis of their stability which is discussed in this section. The stable periodic equilibrium configurations are identified through an optimization search algorithm that seeks for local energy minima [16]. The cell interacting with a reference unit are denoted with the index n = −N, …, N in 1D lattices, and the index pair n = −N, …, N; m = −M, …, M in 2D lattices.

In order to reduce the number of possible combinations only configurations consisting of four or less particles are considered in 1D lattices, and only configurations consisting of two particles are considered in 2D lattices. This is enforced by imposing an E-particle periodicity (where E = 1, 2, or 4) as a constraint to the optimization problem defining the equilibrium conditions. This is practically implemented by imposing periodic conditions prior to the minimization of the static energy functional, i.e. by letting:(7)with n = −N, …, N; m = −M, …, M, and where for simplicity q = q_0,0 denotes the DOFs of the reference unit. This implies that the unknowns of the nonlinear algebraic problem seeking for a minimum of the energy functional are: the unit cell DOFs q including the coordinates x, y of the particles within a unit cell and their rotation θ, and the lattice vectors d₁, d₂ in the case of 2D lattices. The resulting nonlinear system of equations is expressed as:(8)where x = [q, d₁, d₂]. The problem in equation (8) is solved through an unconstrained Quasi-Newton optimization scheme as described by Rao [17].

Results of this process are shown in Figures 1 and 2 for 1D and 2D lattices respectively. All 1D configurations pictured in Figure 1 include periodic systems with 1–4 particle unit cells. All 2D configurations to be discussed are pictured in Figure 2. Additional examples of equilibrium configurations are found in reference [1].

Figure 1. Periodic equilibrium examples for 1D structures with one (a), two (b–d), and four (e–g) particle unit cells. Red and blue circles denote magnets with dipole moments perpendicular to the lattice plane pointing towards and away from the reader respectively. Black lines denote elastic connections between magnets. Reprinted from M. Schaeffer and M. Ruzzene, Wave propagation in multistable magneto-elastic lattices, Int. J. Solids Struct. 56–57 (2015) 78–95, copyright (2015), with permission from Elsevier.

Figure 2. Periodic equilibrium examples for 2D structures with two particle unit cells. Red and blue circles denote magnets with dipole moments perpendicular to the lattice plane pointing towards and away from the reader respectively. Black lines denote elastic connections between magnets. Particles in the unit cell have been highlighted by placing yellow dots in their centers. Reprinted from M. Schaeffer and M. Ruzzene, Wave propagation in multistable magneto-elastic lattices, Int. J. Solids Struct. 56–57 (2015) 78–95, copyright (2015), with permission from Elsevier.

3 Homogenization procedure

3.1 Linearization about equilibrium

The considered homogenization approach assumes a linear elastic continuum. This requires linearization of the governing equations under the assumption that the lattice undergoes small displacements from a stable equilibrium.

A first order Taylor series expansion is employed to linearize the forces. Thus, the gradient of the forces with respect to the unit cell DOFs defines a stiffness matrix, and the equations of motion for a unit cell can be expressed as:(9)where q_m,n now denotes the displacement from equilibrium of the DOFs in cell m, n, and the absence of external forcing is assumed. In equation (9), K_m,n describes both the linearized elastic and magnetic interactions between the considered unit cell, here identified by the subscripts m = 0, n = 0, and the (2N + 1)(2M + 1) − 1 connected neighboring cells.

3.2 Reduction of the degrees of freedom

The homogenization method, developed by Gonella and Ruzzene [10], derives an equation of motion from the lattice structure that is compared to the continuum equations of motion in order to extract the mechanical properties. Before the method is applied the lattice must be reduced to an equivalent system with only one particle per unit cell, which ensures that the spacing between the particles makes a square grid in the lattice coordinate space. This is required to achieve cancellation between the 1st order partial derivatives that would otherwise manifest in the equation of motion derived by the homogenization method and break the analogy to the continuum equation of motion. Guyan reduction [18] accomplishes this task. Through the algebraic manipulation and substitution in Guyan reduction the n_m DOFs of the master particles m are used to define the n_s DOFs of the slave particles s. Any particle in a unit cell can be chosen as the master particle as long as the same particle is chosen in each unit cell. Partitioning the equation of motion into master and slave DOFs gives:(10)where n_dof = n_m + n_s is the number of DOFs describing the particles within the distance r_∞ from a reference unit cell.

Though the equilibrium geometry of the structure is periodic, the deformation is arbitrary and cannot be described by the reference unit cell alone. Therefore, equation (10) is the equation of motion for an entire structure, not just one cell as equation (9), and [q_m; q_s] contains the q_m,n corresponding to all m, n pairs. Having equation (10), a transformation matrix is then calculated, which, when applied, solves for the master DOFs in terms of the slave DOFs(11)where I is the identity matrix. Applying the transformation matrix produces the reduced mass and stiffness matrices , and the displacement vector .(12)

In final preparation for homogenization methodology the mass and stiffness matrices are partitioned into blocks, each of which describe the effect of cell m, n’s DOFs on the reference cell’s DOFs. For the lattices studied, equation (13) separates out three equations from the n_m reduced equations of motion, which correspond to the three DOFs of the master node in the reference unit cell.(13)

It is worth noting here that for a static problem Guyan reduction introduces no assumptions, and thus no error regarding static elastic properties.

3.3 Long-wavelength approximation

The homogenization method follows the procedure described by Gonella and Ruzzene [10] and is summarized here. First the lattice coordinate space η₁, η₂ is defined using the lattice vectors d₁, d₂:(14) (15)

Then it is assumed that the displacement of the particles is continuous in space and can be represented by a second order Taylor series in the lattice coordinate space η₁, η₂. Thus, the master DOFs of the cell md₁ + nd₂ away from the reference unit cell may be approximated as:(16)

Application of equation (16) to equation (13) produces:(17)where, although the absence of external forcing will be assumed, a force vector f continuous in x, y is included for generality. Equation (17) assumes that for m, n ≠ 0, meaning there is no inertial coupling between cells. Of note is the fact that the matrices in equation (17) are dependent not on nearest neighbors only but on the cells up to M or N indexes away, where M, N ≠ 1 in general, allowing the inclusion of non-nearest neighbor magnetic interactions. The matrices in equation (17) are defined:(18)

Equation (18) differs from that in reference [10] as the typo omitting the 1/2 from A′ and C′ has been fixed. In order to convert to the physical space x, y from the lattice coordinate space η₁, η₂ the derivatives are converted using:(19) (20)

Note that these equations differ from those in reference [10] slightly as two small typos have been fixed. The change of basis gives the continuum equation of motion in the physical space x, y:(21)where A, …, F, F_m are defined in terms of A′, …, F′, F′_m and E as follows:(22)

As a comparison will be made with the equation of motion for a continuum the equation must be expressed per unit volume V, which for a 2D structure is:(23)where b is the out-of-plane thickness of the structure. The new equation of motion is:(24)where b is a vector of body forces.

From here the procedure of [10] is continued, omitting terms that decrease with the length scale of the micro structure ε ³. Then the rotational DOF θ is algebraically removed (assuming no moments are applied) and the two resulting equations are:(25)where q_x, q_y and b_x, b_y are the components of q and b respectively, and a … f are constants. Equation (25) is compared with the equations of motion for a thin orthotropic material (plane stress) which are:(26)

Comparing like terms allows algebraic solution of the density ρ the Young’s moduli E₁ and E₂, the Poisson’s ratios ν₁₂ and ν₂₁, and the shear modulus G₁₂ in terms of a … f.

4 Homogenization of 1D magneto-elastic lattices

The procedure is first applied to the 1D lattices shown in Figures 1a and 1b. The analytical expressions of the equivalent properties of both structures derived below are used to check the convergence of the approximate solutions described in the previous section.

For both lattices the same approach is used. A variable s, used to denote the pattern type, allowing specification between the two lattices mathematically. It is −1 for alternating polarities (Figure 1b) and 1 if all polarities are the same (Figure 1a). First the equation of motion for the generic particle i is derived. It is assumed that the lattice is in equilibrium, and evenly spaced with distance h between all nearest-neighbor magnets. The equation of motion is:(27)where u_i is the displacement of particle i in the x direction, and is the equivalent linearized stiffness of the magnetic interaction with a particle distance hn away from particle i.

Assuming a continuous distribution of displacement, a second order Taylor series can be expanded about particle i for a step n, which approximates u_i+n. The i is dropped from the notation as its value is not important so u_i, u_i+n → u, u_n, yielding:(28)

For 1D lattices the unit cell volume V is:(29)where b and w are the thicknesses of the rectangular prism unit cell in the 2 directions perpendicular to the length of the lattice. Substituting equation (28) into equation (27) and dividing by the unit cell volume V simplifies to:(30)

By comparing equation (30) to the 1D wave equation(31)where ρ is the material density, it can be seen that the equivalent Young’s modulus E is:(32)

Linearizing equation (3) it can be determined that:(33)where r_n is the distance from particle i to particle n and s is defined above.

4.1 Young’s modulus for two particle cells

For convenience, the modulus of elasticity is normalized to the unmagnetized case (ψ ^(m) = 0), E_o. If all magnets are repelling (s = 1) (Figure 1a) the normalized modulus of elasticity is:(34)where h_o = h/r₀ and ζ is the Reimann zeta function defined as [19]:(35)

Since there exist methods to accurately calculate the Reimann zeta function and it is computationally reasonable to achieve a value for it to machine precision the radius of influence r_∞ is not used in this case and no modeling error is introduced by omitting magnetic interactions in this 1D structure. Furthermore, the normalized magnetization is:(36)

The normalized value can be interpreted as the ratio of the magnetic force at r = r₀ to the force required to fully compress the axial spring in the absence of contact, i.e. k_ar₀. Hence, this normalized quantity provides a reference of the strength of the magnetic interaction with respect to the spring restoring force. In general, will simply denote magnitude for a structure, since attraction and repulsion of magnets depends on the specific magnet pair in a lattice. Therefore when the terms “high” and “low” are used to refer to it is the magnitude that is being discussed [1].

If the polarity of the magnets alternates (s = −1) (Figure 1b) the normalized modulus of elasticity is:(37)

Now only h_o must still be defined. The distance between each particle h is a stable critical point of the energy functional equation (4) when reduced to its static terms. This leads to the following polynomial problem for nearest neighbors repelling and attracting respectively:(38) (39)

The variation of the equivalent Young’s modulus with is presented in Figure 3. The figure describes the tunability of the stiffness with changes in magnetization for simple structures defined by a 1 and 2 magnet unit cell. It is observed that the presence of repulsive magnetic interactions produces a hardening of the structure (solid black curve), and the presence of attractive nearest neighbor magnetic interactions has a softening effect (dashed red curve). There is a critical magnetization for the case with attracting nearest neighbors near where the structure approaches 0 stiffness as the structure becomes unstable. This is the magnetization magnitude at which the magnetic forces overpower the elastic forces from the axial springs. Beyond this point only reconfigured versions of the structure are stable.

Figure 3. Exact analytical solution of normalized modulus of elasticity for: the structure pictured in Figure 1a (black –) and the structure pictured in Figure 1b (red – –).

4.2 Young’s modulus for four particle cells

The larger a unit cell becomes, the more possible configurations there are. The 1D results in this paper only consider unit cells of up to four magnets for simplicity. Furthermore, not all possible magnetic polarity patterns are considered. The results from four magnet unit cells are now discussed, with two magnet unit cells considered as a special case of four magnet unit cells. The group of configurations considered here are those pictured in Figures 1b–1g. The configurations in this group can reconfigure into any other configurations in the group, provided that both can exist for a given magnetization. Any structures that have contacting particles require a minimum magnetization to hold the particles in contact against the compressed spring’s force. In addition, when a structure is magnetized highly enough certain configurations will not be possible as magnetic forces will overcome elastic forces. Therefore, there are values of that do not permit all the configurations. For the following results a contact stiffness of α = 40 and a contact distance of r_c = 0.5r₀ are used. The properties of the different configurations are shown in Figures 4 and 5. For this single system can take any of the six studied configurations, which have different stiffnesses and densities. The modulus of elasticity E and density ρ are normalized by the modulus of elasticity and density calculated when , which are E₀ and ρ₀. Conveniently, the normalized density is the reciprocal of the normalized equilibrium unit cell width d, i.e. d/d₀ = (ρ/ρ₀)⁻¹, where d₀ is the equilibrium distance calculated when . For the simplest cases, which are shown in Figure 3, d/d₀ = 2h/2r₀) = h₀.

Figure 4. Normalized modulus of elasticity of structures with an even number of magnets per unit cell; configuration in Figure 1b (black circles), configuration in Figure 1c (red squares), configuration in Figure 1e (red – –), configuration in Figure 1f (black –), configuration in Figure 1g (blue – ·).

Figure 5. Normalized density of structures with an even number of magnets per unit cell; configuration in Figure 1b (black circles), configuration in Figure 1c (red squares), configuration in Figure 1d (black rings), configuration in Figure 1e (red – –), configuration in Figure 1f (black –), configuration in Figure 1g (blue – ·).

Figure 6. Comparison between normalized modulus of elasticity calculated by homogenization and numerical simulation for the configuration in Figure 1a (blue – · and black squares respectively), the configuration in Figure 1b (red – – and black circles respectively), and the configuration in Figure 1c (black – and black diamonds respectively).

The Young’s modulus and its variation with is presented in Figure 4. Through reconfiguration alone there are as many values available as configurations. The structures corresponding to Figure 4 all soften in response to an increase of their magnetization. This is due to nearest neighbors that are in attraction, because the linearized model of attracting magnets exhibits negative axial stiffness. If both reconfiguration and changes in are considered there is no stiffness between the values of E/E₀ ≈ 0.2 and 2 that cannot be achieved by at least one configuration, and there are some stiffnesses which can be achieved by as many as four configurations. This shows that stiffness adaptation is not limited to a few discrete values, allowing easier integration into possible design applications. The structure pictured in Figure 1d is shown separately in Figure 7 because its stiffness is much greater than the other structures’ as the contact stiffness used is α = 40. The structure is not a full 40 times stiffer though because the attractive magnetic interactions soften the structure, even under the minimum magnetization required for it to exist. Having only contact interactions, the structure in Figure 1d is the stiffest possible configuration for a structure with alternating magnetic poles and the parameters considered here, and thus gives an upper bound. The same can be said for the density, as contacting particles provides the closest packed structure.

Figure 7. Comparison between normalized modulus of elasticity calculated by homogenization and numerical simulation for the configuration in Figure 1d (red – – and black circles respectively).

Figure 8. Normalized moduli and Poisson’s ratio for the hexagonal lattice with repelling magnets (a, c, and e) and with attracting nearest neighbors (b, d, and f). k τ /k τo  = 102 (red – –), k τ /k τo  = 100 (black triangles), k τ /k τo  = 10−1 (blue – ·), and k τ /k τo  = 10−3 (black –).

Figure 9. Poisson’s ratio for the hexagonal lattice with repelling magnets (a) and with attracting nearest neighbors (b) for values of k τ that allow changes in to change the sign of ν. k τ /k τo  = 10−0.1 (blue circles), k τ /k τo  = 10−0.2 (red squares), k τ /k τo  = 10−0.3 (black –), k τ /k τo  = 10−0.4 (red – –), k τ /k τo  = 10−0.5 (blue – ·), and k τ/k τo  = 10−0.6 (black diamonds).

The densities of the 2 and 4 magnet configurations are presented in Figure 5. The density is primarily effected by the configuration, but there is some change with . The density can be approximately doubled through reconfiguration. The configurations pictured in Figures 1c and 1f exhibit similar densities because both have the same ratio of contacting to non-contacting magnets.

The configurations pictured in Figures 1c and 1f are worthy of note because they would be considered equivalent in terms of E and ρ if only nearest neighbor interactions were considered. In fact, for low enough they have almost equivalent properties. But with higher values of their dissimilar properties show that non-nearest-neighbor interactions have a significant effect. Even so, the density does not change significantly with increasing , which means that a structure can have little change in its global size but a significant change in its stiffness by changing only its internal configuration. For example, if reconfiguring from the configuration in Figure 1c to the one in Figure 1f would produce only a 2.7% increase in density, but a 33% decrease in stiffness.

4.3 Numerical simulations

For the purpose of validating the homogenization method, as applied to the magneto-elastic lattices, numerical simulations are completed. The equations of motion are integrated using an explicit, variable-step Runge-Kutta method to obtain the time history of the positions of the particles in a finite lattice. This provides a framework for performing numerical experiments. A finite sample of each configuration to be validated, 40 particles long, is compressed quasi-statically over a strain of less than 0.001, in which the structure can be considered linear. The compression is performed through an applied displacement of one end of the system while the other end is held fixed. There is linear damping included as described in reference [1] to damp out vibrations, so slow enough forcing will ensure that dynamic forces are negligible. To obtain a strain rate slow enough to consider quasi-static a convergence study is performed. The force versus displacement information is used to calculate a stiffness simply by comparing the change in reaction force to the structure’s axial strain. By normalizing to the case when the simulation results can be compared with the homogenization results directly. Their comparison is shown in Figures 6 and 7, with curves representing the approximate analytical results and the shapes representing the numerical results.

5 Results: 2D lattices

For simplicity the results presented for 2D magneto-elastic structures are restricted to systems described by a two magnet unit cell. A magnetization pattern where all magnets repel (Figure 2a) is used to investigate and compare to the pattern where nearest neighbors attract (Figure 2b). The magnetization pattern with attracting nearest neighbors is used to showcase the drastic changes in properties that can be achieved by reconfiguration. The reconfiguration considered is the transition to the re-entrant lattice (Figure 2c) as this is possibly the easiest way to reconfigure the magneto-elastic hexagonal lattice studied based on previous research by the authors Schaeffer and Ruzzene [1, 20]. In contrast with the hexagonal patterns, the re-entrant lattice is anisotropic. Figure 10a defines the 1 and 2 directions for the re-entrant lattice, which are used to identify the orthotropic properties. Furthermore, the re-entrant lattice cannot exist for all the same values of k_τ/k_τo and as the hexagonal lattice. Namely, there is an upper limit on k_τ/k_τo near 0.1. To explore the tunability of the equivalent mechanical properties of the aforesaid lattices the homogenization procedure as described in Section 3 is applied to structures with varying torsional stiffnesses k_τ and magnetization magnitudes . As in the previous section a contact stiffness of α = 40 and a contact distance of r_c = 0.5r₀ are used. The properties are normalized with respect to a reference case defining the reference properties E_o and G_o, which is the unmagnetized () hexagonal lattice with k_τo = 1 N/rad and k_ao = 1 N/m.

Figure 10. Re-entrant unit cell with coordinate system 1, 2 pictured (a), and θ defined as the angle between coordinate axes 1, 2 and x, y (b).

The first results presented are for the two considered hexagonal lattice structures (Figures 2a and 2b). The homogenization produces properties that describe an isotropic material, i.e. E₁ = E₂ = E, ν₁₂ = ν₂₁ = ν, and G₁₂ = E/(2(1 + ν)) = G, which is to be expected for regular hexagonal geometry [4]. Variation in k_τ is seen to cause substantial changes in material stiffness so the value of k_τ is incremented exponentially. Indeed, as k_τ → 0 the quantities E, G → 0 and the lattice approaches a mechanism. In this state the structure’s global deformation is dominated by the rotational deformation of the torsional springs. In such a case magnetization has a stiffening effect on the structure, which is evident in Figures 8a–8d for the cases where k_τ/k_τo = 10⁻³ (solid black line), though it is more prominent in the lattice with only repelling magnetic interactions. The stiffening of the structure can be explained predominantly through the repulsion of the next-nearest neighbor interactions, which are repulsive in both of the lattices discussed here. As the lattices are dominated by rotational deformation, the nearest neighbor magnetic interactions do not significantly effect the global properties. The exception comes when the the effective axial stiffness is nearly 0 close to the point of instability for attracting nearest neighbors, e.g. near . Though next-nearest neighbor interactions can explain the dominant trend for both magnetization patterns, comparison of Figures 8a and 8b reveals that more distant magnetic interactions have an non-negligible effect on the structure. Thus, efforts to include the non-nearest neighbor interactions appear to be valuable. Significant tunability is seen for low k_τ when the magnetic interactions dominate the structure’s properties. For example, in the case where k_τ/k_τo = 10⁻³ in Figure 8a E increases by a factor of 48 over the range of shown.

Increasing k_τ causes the global deformation of the lattice to be dominated by the deformation of the axial springs in the lattice. In the limit as k_τ → ∞ the structure will become like the one discussed in reference [21]. The axial deformation dominated structures are stiffened by repulsive nearest neighbor interactions and softened by attractive nearest neighbor interactions, due to the equivalent linearized axial stiffness of magnetic interactions being positive for repulsion and negative for attraction. The greatest torsional stiffness shown is k_τ/k_τo = 10² (dashed red line). Increasing k_τ/k_τo above this does not have a significant effect on the mechanical properties. In addition, decreasing the torsional stiffness below k_τ/k_τo = 10⁻³ does not make a significant change. Therefore, the properties of the hexagonal lattices discussed are bounded by the solid black and red dashed curves in Figure 8.

An interesting property of both magnetization patterns is that the Poisson’s ratio can be either positive or negative depending on the torsional stiffness k_τ and the magnetization (Figures 8e and 8f). In the limit as k_τ → 0, when the compliance of the system is dominated by the torsional springs, it is observed that ν → 1, as expected for the analogous structure of the regular hexagonal honeycomb with slender beams [4]. In the limit as k_τ → ∞ it is observed that ν → −1/3, which is consistent with equations given by Grima et al. [21] for regular hexagonal lattices with the angles between bars constrained. In this second case the structure’s compliance is dominated by the axial springs. Thus, it is clear that the ratio between rotational and axial stiffness affects the sign of ν. Since magnetization can have a similar effect to changing the rotational stiffness, as discussed above, for proper values of k_τ it is possible to switch between positive and negative ν simply by changing . The ranges of k_τ that allow for such tunability of ν for both magnetization patterns discussed are shown in Figure 9.

As stated previously the magneto-elastic hexagonal lattice with attracting nearest neighbors can be reconfigured into the re-entrant lattice (Figure 2c) where strong contact forces hold the structure in equilibrium against magnetic attraction. The homogenized properties for varying values of k_τ/k_τo and are presented in Figure 11 and they are compared to the hexagonal lattices that reconfigure into them. In contrast to analyses of re-entrant lattices such as those by Gibson and Ashby [4] or Gonella and Ruzzene [10], where the geometry of a cell is explicitly defined, the re-entrant lattices discussed in this article have their geometry defined by the equilibrium configuration associated with a chosen value of . The equilibrium length of one axial spring relative to another may differ by more than 10% depending on the choice of k_τ/k_τo and , but the equilibrium lengths are approximately equal for low k_τ/k_τo and . Therefore, the tuning of properties with changing discussed hereafter are caused in part by the change in the re-entrant lattice equilibrium geometry.

Figure 11. Normalized equivalent orthotropic properties for the re-entrant lattice with α = 40 and k τ /k τo  = 10−1 (blue · –), and k τ /k τo  = 10−3 (black –).

Reconfiguration from the hexagonal to the re-entrant lattice has a drastic effect on the stiffness of the lattice. The drastic change is seen comparing in Figures 11, 8b, 8d, and 8f. Two choices of k_τ/k_τo are used. The value k_τ/k_τo = 10⁻¹ is near the maximum value of k_τ/k_τo for which the re-entrant lattice can exist, and the value k_τ/k_τo = 10⁻³ is representative of the asymptotic solution as k_τ/k_τo → 0. Values for E₁ and E₂ (Figures 11a and 11b) in the re-entrant lattices are generally much higher than E in hexagonal lattices of equal k_τ/k_τo. The internal structure of the re-entrant lattice makes strain along the one and two directions produce axially dominant local deformation, regardless of the torsional stiffness. Furthermore, the stiff contact interactions make E₁ stiffer than values of E associated with axial deformation dominant hexagonal lattices with the same and magnetization pattern. In the two direction, the increased density from the reconfiguration makes values of E₂ stiffer than associated values of E for axial deformation dominant hexagonal lattices. For the k_τ/k_τo = 10⁻³ case the Young’s moduli of the hexagonal configuration are dominated by magnetic interactions, so since the minimum magnetization for the existence of the re-entrant configuration is low in this case the Young’s moduli can be increased by more than two orders of magnitude simply through structural reconfiguration to the re-entrant lattice. Even for higher magnetizations there is over a 10-fold increase in stiffness due to reconfiguration. For k_τ/k_τo = 10⁻¹ the structure is still significantly stiffened by reconfiguration, but not as drastically as in the previous case. The k_τ/k_τo = 10⁻¹ also has a more limited range of in which it can exist in the re-entrant configuration.

As in the hexagonal lattice the shear modulus G₁₂ in the re-entrant lattice is dependent on k_τ. When k_τ/k_τo = 10⁻³ reconfiguration can account for over order of magnitude increase in G₁₂, but when k_τ/k_τo = 10⁻¹ the shear moduli for the hex and re-entrant lattices are more comparable. The Poisson’s ratios also change significantly with reconfiguration, with the re-entrant lattices exhibiting lower Poisson’s ratios than corresponding hexagonal lattices in general. Although the re-entrant lattice structure is generally known for being auxetic [4] this magneto-elastic example has positive ν₁₂ and ν₂₁ because of the stiff contact interactions that prevent the structure from folding in on itself. However, sufficiently soft contact interactions do allow negative ν₁₂ and ν₂₁.

In the re-entrant lattice the shear moduli are small compared to the axial stiffnesses, and this produces significant variation of stiffness with respect to direction. Figure 12 shows the variation of the equivalent properties in the x, y coordinate system, which is redefined here to be oriented by θ with respect to the nominal coordinate system 1, 2 as defined in Figure 10b. The properties of a re-entrant lattice (solid black line) with k_τ/k_τo < 10^−1.5 is compared to its hexagonal configuration (dashed blue line). Note that the radii of the polar plots in Figures 12a–12c are logarithmic. It is evident that the re-entrant lattice is much softer for θ = 45 degrees (Figures 12a and 12b solid black line) than θ = 0 degrees, becoming over half an order of magnitude softer. As the Young’s moduli soften with changing θ the shear modulus stiffens by almost one order of magnitude (Figure 12c). Simultaneously, there are substantial changes in the Poisson’s ratios, though they remain positive. So, even when E₁ and E₂ are nearly equal the re-entrant lattice is far from isotropic in contrast to the hexagonal lattice.

Figure 12. Variation of mechanical properties with coordinate axis orientation θ of the hexagonal (blue – –) and re-entrant lattice (black –) with kτ/kτo = 10−1.5, , α = 40. (a) log10(Ex/Eo), (b) log10(Ey/Eo), (c) log10(Gxy/Go), (d) νxy, (e) νyx.

6 Conclusion

The equivalent continuum properties of magneto-elastic lattices are calculated through the extension of an existing homogenization method. Using 1D structures it is shown that many stiffnesses can be achieved by reconfiguring the structure between different stable equilibria and by changing the lattice magnetization. Using 2D structures the isotropic mechanical properties of hexagonal lattices are shown to vary with lattice physical properties and lattice magnetization. This includes the ability to change sign of the Poisson’s ratio of the structure simply by changing the lattice magnetization magnitude and the ability to stiffen or soften the lattice. The effect of reconfiguration from the hexagonal to the re-entrant lattice on mechanical properties is presented, and the variation of the equivalent orthotropic properties of the re-entrant lattice are shown to be tunable with changing magnetization magnitude.

Acknowledgments

This research was funded by a grant from the Air Force Office of Scientific Research Monitored by Dr. David Stargel, whose input and guidance are acknowledged and appreciated. The work of William Storrs, who performed the numerical simulations for validation, is also acknowledged and appreciated.

References

All Figures

Figure 1. Periodic equilibrium examples for 1D structures with one (a), two (b–d), and four (e–g) particle unit cells. Red and blue circles denote magnets with dipole moments perpendicular to the lattice plane pointing towards and away from the reader respectively. Black lines denote elastic connections between magnets. Reprinted from M. Schaeffer and M. Ruzzene, Wave propagation in multistable magneto-elastic lattices, Int. J. Solids Struct. 56–57 (2015) 78–95, copyright (2015), with permission from Elsevier.

In the text

Figure 2. Periodic equilibrium examples for 2D structures with two particle unit cells. Red and blue circles denote magnets with dipole moments perpendicular to the lattice plane pointing towards and away from the reader respectively. Black lines denote elastic connections between magnets. Particles in the unit cell have been highlighted by placing yellow dots in their centers. Reprinted from M. Schaeffer and M. Ruzzene, Wave propagation in multistable magneto-elastic lattices, Int. J. Solids Struct. 56–57 (2015) 78–95, copyright (2015), with permission from Elsevier.

In the text

	Figure 3. Exact analytical solution of normalized modulus of elasticity for: the structure pictured in Figure 1a (black –) and the structure pictured in Figure 1b (red – –).
In the text

	Figure 4. Normalized modulus of elasticity of structures with an even number of magnets per unit cell; configuration in Figure 1b (black circles), configuration in Figure 1c (red squares), configuration in Figure 1e (red – –), configuration in Figure 1f (black –), configuration in Figure 1g (blue – ·).
In the text

	Figure 5. Normalized density of structures with an even number of magnets per unit cell; configuration in Figure 1b (black circles), configuration in Figure 1c (red squares), configuration in Figure 1d (black rings), configuration in Figure 1e (red – –), configuration in Figure 1f (black –), configuration in Figure 1g (blue – ·).
In the text

	Figure 6. Comparison between normalized modulus of elasticity calculated by homogenization and numerical simulation for the configuration in Figure 1a (blue – · and black squares respectively), the configuration in Figure 1b (red – – and black circles respectively), and the configuration in Figure 1c (black – and black diamonds respectively).
In the text

	Figure 7. Comparison between normalized modulus of elasticity calculated by homogenization and numerical simulation for the configuration in Figure 1d (red – – and black circles respectively).
In the text

	Figure 8. Normalized moduli and Poisson’s ratio for the hexagonal lattice with repelling magnets (a, c, and e) and with attracting nearest neighbors (b, d, and f). k τ /k τo  = 102 (red – –), k τ /k τo  = 100 (black triangles), k τ /k τo  = 10−1 (blue – ·), and k τ /k τo  = 10−3 (black –).
In the text

	Figure 9. Poisson’s ratio for the hexagonal lattice with repelling magnets (a) and with attracting nearest neighbors (b) for values of k τ that allow changes in to change the sign of ν. k τ /k τo  = 10−0.1 (blue circles), k τ /k τo  = 10−0.2 (red squares), k τ /k τo  = 10−0.3 (black –), k τ /k τo  = 10−0.4 (red – –), k τ /k τo  = 10−0.5 (blue – ·), and k τ/k τo  = 10−0.6 (black diamonds).
In the text

	Figure 10. Re-entrant unit cell with coordinate system 1, 2 pictured (a), and θ defined as the angle between coordinate axes 1, 2 and x, y (b).
In the text

	Figure 11. Normalized equivalent orthotropic properties for the re-entrant lattice with α = 40 and k τ /k τo  = 10−1 (blue · –), and k τ /k τo  = 10−3 (black –).
In the text

	Figure 12. Variation of mechanical properties with coordinate axis orientation θ of the hexagonal (blue – –) and re-entrant lattice (black –) with kτ/kτo = 10−1.5, , α = 40. (a) log10(Ex/Eo), (b) log10(Ey/Eo), (c) log10(Gxy/Go), (d) νxy, (e) νyx.
In the text