Homogenization of fully nonlinear rod lattice structures: on the size of the RVE and micro structural instabilities

This work investigates the possibility of applying two-scale computational homogenization to rod lattice structures emerging, for instance, from additive manufacturing. The influence of the number of unit cells within the representative volume element (RVE), thus, the RVE’s size on the homogenized mechanical response is studied for occurring microscopic structural instabilities. Therein, the macro-scale, described in terms of three-dimensional continuum mechanics, is coupled to the micro-scale described by geometrically exact rods, enabling arbitrary large deformations and rotations. A special feature of the presented framework is that the rods building the lattice structures are not restricted to deform purely elastically but may deform inelastically. The mechanical response of lattice structures is investigated by applying the developed homogenization method to an exemplary lattice. Under special loads the structure reaches an instable state and may buckle. The appearance of instabilities depends on the geometric properties of the lattice’s underlying rods and the RVE’s size.


Introduction
In recent years additive manufacturing (AM) has experienced a tremendous increase in popularity [1]. Originally developed for rapid prototyping, AM has evolved into a tool which nowadays is also used to manufacture whole components. Besides the fact that the quality and the performance of the produced samples increased rapidly, a further aspect leading to this development is that AM opens the opportunity to fabricate meta-materials with cellular micro-structures. Such meta-materials are characterised by tailored mechanical properties such as anisotropy, spatially varying stiffness and density as well as sophisticated mechanical functionalities [2,3], e.g. auxetic properties [4]. In the case considered in this contribution the cellular meta-material is build from slender rods resulting in a rod lattice structure. By varying the arrangement, thickness and cross-sectional geometry of individual rods engineers are able to control the resultant mechanical properties of lattice structures. B Ludwig Herrnböck ludwig.herrnboeck@fau.de 1 Institute of Applied Mechanics, Friedrich-Alexander University Erlangen-Nuremberg, Egerlandstraße 5, 91058 Erlangen, Germany Common manufacturing processes used to produce lattices include electron beam melting, selective laser melting as well as selective laser sintering. Rod lattices find their applications mostly in the area of lightweight structures in the aerospace industries as well as in the medical field [1,5]. In theses applications the lattice is assumed to be stiff and a linear stress-strain relation is typically employed. However, applications in which the stress-strain relation may no longer be assumed linear are emerging more and more often. Possible nonlinearities may appear on the geometrical level as well as on the material level -for instance due to hyperelastic or inelastic material behaviour in the micro-structure. Examples of materials combining both geometric and material nonlinearities include cushioning shoe soles [6], synthetic meshes, used for example in hernia repair surgery [7] and impact absorbing materials used in defence applications such as body and object protection [8][9][10]. The latter is designed to absorb maximum energy at nearly constant stress. Geometric nonlinearities bring new challenges to the mechanical modelling of rod lattices due to finite deformations and the possibility of buckling rods in the lattice structure. Then, the use of linear rod theory is no longer possible. A rod model allowing for large deformations and rotations must be introduced.
A common procedure used to estimate the effective macro-scale mechanical response of a meta-material is twoscale computational homogenization. This approach estimates the effective mechanical response by homogenizing the mechanical behaviour of a micro-structure (micro-scale) based on a representative volume element (RVE) where the RVE consists of periodically repeating unit cells of the metamaterial. A widespread convention is to define the RVE as the smallest repeating unit cells. In this contribution different sizes of an RVE will be analysed, i.e. RVEs with different numbers of unit cells.
Two-scale homogenization of rod lattices has been addressed in numerous work [2,[11][12][13][14]. The two-scale homogenization may be used to derive a material model able to replicate the mechanical response of the meta-material for the purely elastic case. This may be done by e.g. fitting an appropriate strain energy function to the homogenized properties of an RVE consisting of the smallest repeating unit cell or by setting up an artificial neuronal network. Alternatively two-scale homogenization may also be used to obtain the constitutive behaviour of the meta-material by on-the-fly homogenization. The aim of this paper is to investigate if the number of unit cells within the RVE, i.e., the size of the RVE affects the homogenized properties. Especially, this is done for load cases leading to buckling within the RVE, since it is known from [15][16][17] that local instabilities at the level of the micro-structure may be influenced by the size of the RVE. Further, it is explored if the mechanical response converges by increasing the RVE's sizes. To this end, a computational homogenization framework that combines well known threedimensional continuum mechanics on the macro-scale with an appropriate beam formulation on the micro-scale is developed. The mechanical response of a representative lattice structure shall be analysed under the consideration of large deformations and inelastic behaviour. This involves both the pre-and post-buckling state for both the purely elastic and inelastic case. For this purpose an RVE modelled by geometrically exact elastoplastic rods is considered. In contrast to the classical beam formulation, such as Euler-Bernoulli or Timoshenko beams, the theory of geometrically exact rods (also known as the Cosserat theory of rods [18,19]) is not restricted to the geometrically linear case and is suitable for large displacements and rotations [18][19][20][21]. Hence, this theory is widely used to model long nano wires and tubes, cables, hair and other biological fibres [22][23][24][25]. An extension to the inelastic case is presented in Smriti et al. [26] and further elaborated in [27]. The bottleneck in modelling elastoplastic geometrically exact rods is the lack of knowledge of an appropriate yield-function in terms of the rod's stress resultants (internal contact forces and internal contact moments). Yield-functions have been derived for different rod crosssections in [28][29][30][31] by applying the upper and lower bound techniques of limit analysis in which the rod's cross-section is assumed to be fully plastified at the yield-limit. An approach where yield is not defined in terms of plastification but in terms of dissipated work has recently been formulated in Herrnböck et al. [32]. The method presented there is generic in terms of the rod's cross-sectional shape and underlying constitutive model.
Considering the aim of this work, the following structure is prescribed. In Sect. 2 a brief review of the theory of geometrically exact rods with rigid cross-sections is presented for both the purely elastic and inelastic case. The presentation of the inelastic case involves the derivation of the evolution equations for the plastic quantities. In Sect. 3 the constraints necessary to model networks of multiple rods are introduced, where it is assumed that the displacement and rotation of two distinct rods are equal at their intersection. Finally, in Sect. 4 a first-order homogenization framework combining classical three-dimensional continuum mechanics on the macro-scale with the geometrically exact rod theory on the micro-scale is proposed. Thus, this homogenization framework allows for large deformations and rotations of rods constituting the underlying rod lattice while allowing inelastic behaviour. The concepts introduced in Sects. 2-4 are finally used to obtain numerical results investigating the macroscopic mechanical response of an exemplary periodic rod lattice. These results are presented in Sect. 5. Investigations are performed for different uniaxial strain cases. On the one hand, the focus of the investigations is on the influence of the rods cross-sectional size and, therefore, the change of the mechanical properties of a single rod. On the other hand, an investigation into the stability of the RVE as a function of its size is also performed. The examples presented in this work serve to demonstrate the introduced framework. Imperfections emerging from the manufacturing process are not considered here.
Before setting the stage by recalling the theory of geometrically exact rods the nomenclature used throughout this contribution is introduced. As in Herrnböck et al. [32] strain measures and their energetic conjugates are denoted differently in the context of three-dimensional continuum mechanics and geometrically exact rods. The different terminology is summarized in Table 1.
All simulations carried out to generate the results shown in this contribution are based on the open source finite-element library deal.II [33].

Geometrically exact rods with rigid cross-sections
In this section the theory of geometrically exact rods with rigid cross-sections is presented. This presentation begins with the kinematics followed by the stress resultants and the balance equations for linear and angular momentum. The  [27,32] and the references therein. Figure 1 shows a typical rod deforming from its undeformed straight material configuration B 0 into its deformed spatial configuration B t . In this contribution a rod is assumed to have a constant cross-section along its arc length s. Here X denotes the position of a material point in the rod's undeformed material configuration. Straight rods are considered with the orientation of the cross-sections in the material configuration, independent of the rod's arc length s, described by the set of directors D i . The orthogonal directors D 1 and D 2 span the plane of the undeformed planar cross-section Ω 0 and are usually aligned with the cross-section's principal axis. D 3 = D 1 × D 2 is perpendicular to the rod's undeformed cross-section. The directors D i form an orthonormal basis and are related to the basis E i via

Kinematics
where R 0 denotes a three-dimensional rotation tensor, a member of the special orthogonal group SO (3). That is R 0 describes the orientation of the undeformed cross-sections with respect to the reference frame E i . One may decompose X as follows: where s(s) = s 0 +sD 3 denotes the position along the rod's arc length s by a starting point s 0 and the direction D 3 . Additionally, X C S = X C S α D α describes the material points lying in the cross-section. Here, Latin indices run from 1 to 3, whereas Greek indices run from 1 to 2. The position of the deformed centerline is denoted by r(s). The orthogonal directors d 1 (s), d 2 (s) span the plane of the deformed cross-section Ω t which is assumed to be planar. An orthonormal basis is constructed with d 3 (s) = d 1 (s) × d 2 (s). In this theory, it is not required that d 3 (s) is tangential to the deformed centerline r(s) [20].
Together the directors d i (s) describe the orientation of the deformed cross-section. The directors in the spatial configuration are related to the directors in the material configuration via The quantities r(s) and R(s) are the kinematic unknowns of this theory. Note that the rotation described by the rotation matrix R(s) can also be parametrised by a vector ϑ(s), where ϑ(s) gives the angle of the rotation, and ϑ(s) ϑ(s) describes the direction of rotation [27]. Assuming a rigid cross-section the position x(s) of any material point in the spatial configuration is described by For the sake of readability, the dependencies on s and X C S are omitted from now on. In this contribution, the cross-section of the deformed rod remains rigid. A framework to model the cross-sectional deformation depending on strain measures only is outlined in [34,35].
Strain-prescriptors are introduced as the translational strain v and the rotational strain k (also called curvature): Illustration of a geometrically exact rod undergoing a deformation from its straight material configuration B 0 into its deformed spatial configuration B t . Any point in the material configuration is denoted by X. The position of the centerline of the deformed rod is defined by r.
Material points in the spatial configuration are denoted by x where • defines the derivative of an arbitrary quantity with respect to s as The rotational pull backs v 0 and k 0 of v and k respectively are defined as: Here K and K 0 are skew-symmetric tensors whose axial vectors are k and k 0 respectively. The operator ax(A) transforms an arbitrary skew-symmetric tensor A into its axial vector a such that: We now present the definition of the deformation gradient F of the deformation map (4) The origin of the undeformed rod is set to the origin of a global coordinate system defined by the reference frame E i , i.e., s 0 = 0. Inserting equation (4) into equation (9) yields Recalling the dependencies r = r(s), R = R(s) and X = X(s; X C S ) from equation (2) and (4) respectively the following relations hold: Furthermore, Finally, equation (10) is reformulated as By extracting R in equation (13) the deformation gradient can be written in terms of strain prescriptors as follows:

Stress resultants and balance equations
The internal (contact) force and internal (contact) moment, collectively named stress resultants, are denoted by n and m respectively with: Their rotational pull backs are For the elastic case and assuming the existence of a scalar valued strain energy density function ψ rod (v 0 , k 0 ) the stress resultants follow as The second derivative of ψ rod (v 0 , k 0 ) with respect to the strain prescriptors yields the 6 × 6 tangent stiffness matrix C 0 which relates the change in the stress resultants to the change in the strain prescriptors as Remark 1 A framework to obtain the tangent stiffness matrix of strained elastic rods was presented in Arora et al. [35].
In their work, the assumption is made that the rod's crosssections need not remain plane but can warp due to imposed strain prescriptors (v 0 and k 0 ). The resulting stiffness is a function of the strain prescriptors. The framework allows for the use of arbitrary hyperelastic constitutive models on the cross-sectional level and is extended to the inelastic case in Herrnböck et al. [32]. In the contribution at hand the elastic tangent stiffness matrix is assumed to be constant. This yields a linear relation between stress resultants and strain prescriptors Nevertheless, the constant tangent stiffness matrix must be evaluated by, for example, using the framework presented in [35].
In the elastoplastic case the above definition of the stress resultants is no longer valid (see equation (17)). The key components of modelling elastoplastic rods are briefly presented here. The strain prescriptors are assumed to decompose into an elastic and a plastic part as follows: A free energy density function Ψ rod is introduced as where ψ rod denotes the strain energy density and H describes the hardening depending on the internal variables ξ . The dissipation power D is defined as the difference between the applied stress power and the material time derivative of the free energy density and is required to be greater than or equal to 0: where• := d• dt X . Performing the time derivative by applying chain rule and stating that the resulting equation has to hold for any admissible process one obtains the constitutive equations (similar to equation (17)): Furthermore, we define the hardening-stress β as Inserting the equations (24) and (25) into the dissipation inequality yields the reduced form of the dissipation inequality Equation (26) is maximised under the constraint Φ(n 0 , m 0 , β) ≤ 0. Here Φ describes a yield-function in terms of stress resultants. A framework to obtain such a function, is comprehensively presented in Herrnböck et al. [32].
Finally, the evolution equations for the plastic quantities are obtained: The Lagrange multiplierγ ensures that the resulting stress resultants are admissible and satisfy the Karush-Kuhn-Tucker (KKT) conditionṡ The changes in stress resultants are related to the changes in strain prescriptors via the 6×6 elastoplastic tangent stiffness Due to the additive nature of the strain prescriptors, well known methods from small-strain plasticity may be applied to derive an update algorithm and the algorithmic elastoplastic tangent (for three-dimensional small strain plasticity the reader is referred to [36], for the special case of inelastic rods [27] is recommended).
To conclude this section the balance equations of linear and angular momentum are presented and take the following localized forms n +n = 0, m +m + r × n = 0.
Heren andm represent the distributed external force and moment respectively. In this contribution they are assumed zero, thus: The kinematic unknowns can be determined by solving (31) using the finite element method. The framework used in this work is presented in [21]. In this contribution the rods are discretized using two-noded finite elements with three translational degrees of freedom (DOF) and three rotational DOFs assigned to each node. The kinematic unknowns (r, ϑ) within one element are approximated by linear shape functions. The discretized balance equations are evaluated using a one-point Gauss integration.

Network of rods
In this contribution not only one single rod but networks of multiple rods resulting in periodic lattices are considered. It is assumed that intersecting rods are fixed at their intersection, i.e., their kinematic unknowns r and R must coincide. Let the superscript A and B denote two distinct rods. Their material centerline is defined by Assuming that both rods intersect at s A = l A and s B = l B the spatial position and rotation are constrained via With (33) it is ensured that the connection is fixed througout the deformation and that the angle between two rods remains constant.
Perfectly rigid connections between the rods are assumed throughout this contribution. However, one could also think of pinned or semi-rigid joints. Pinned joints are easily realised by simply omitting the constraints on the rotation matrices R in equation (33). To realise semi-rigid joints the following approach can be pursued. Instead of constraining the rotations of two intersecting rods the rods are connected by a rotational spring with stiffness K . Rigid joints are obtained for K → ∞, while pinned joints result from K → 0. For more details, the reader is referred to [37]. A rheological model that links the rotations and translations of two intersecting rods is also conceivable.

Homogenization
As described in the introduction of this paper the authors aim to investigate the mechanical response of periodic rod lattices. Assuming that the length-scale of a smallest repeating unit cell is sufficiently small compared to the length-scale of the overall lattice the principle of scale-separation applies. Thus, the microand the macro-scale are introduced denoted by superscript m and M, respectively. A representative volume element (RVE) describing the smallest repeating unit of the overall lattice may be considered as a material point of the macroscale. To transfer the mechanical properties of the RVE to the macro-scale methods from multiscale homogenization are useful tools. Such methods are well known and widely discussed for all kind of problems described by nonlinear three-dimensional continuum mechanics. Key references are [15,16]. This section is devoted to a multiscale homogenization framework with the macro-scale described by classical nonlinear three-dimensional continuum mechanics and the micro-scale described by geometrically exact elastoplastic rods. First, the macro-to-micro transition is presented. A relation to transfer the kinematics between the scales is shown and appropriate boundary conditions are derived to successfully fulfil the Hill-Mandel condition. The starting point is the homogenization method from classical three-dimensional continuum theory. Here we strongly orient on Geers et al. [16]. It is then shown that the resulting boundary conditions can be transferred to the kinematic unknowns of the geometrically exact rod theory. Later, the micro-to-macro transition is presented. The rod's stress resultants are transferred to three-dimensional stress measures. Similar approaches can be found in Glaesener et al. [11,12], Weeger [13] and Gärtner et al. [14].

Macro-to-micro
In first-order homogenization a material point in the deformed configuration of the RVE is defined as: The first term in equation (35) denotes a homogeneous deformation defined by the macroscopic deformation gradient F M and the second term describes a field of fluctuations, which is a local contribution to the deformation superimposed to the homogenized deformation. The deformation gradient on the micro-scale is written as A widely used relation to transfer the kinematics between the scales is given by where V m denotes the volume of the RVE's undeformed configuration B m 0 . The RVE has to be solved under appropriate boundary conditions, which have to be derived. This is done by fulfilling the Hill-Mandel condition. We begin by stating that the increment or variation of the work performed on the macro-scale equals the volume average of the variation of work on the micro-scale.
where P denotes the Piola stress, a non-symmetric two point tensor. The Hill-Mandel condition requires This means that the fluctuations in the second term of equation (38) have to vanish. Hereinafter, the index • m is omitted for the sake of readability Applying Einstein's index notation equation (39) can be rewritten as The localized force balance in the RVE yields P i J,J = 0 i . Further, applying the Gauss divergence theorem equation (40) is reformulated as where N J denotes the outwards pointing unit vector normal to ∂B 0 . For the case of a porous RVE the boundary ∂B 0 of the body is split into the shell ∂B S 0 and the faces ∂B F 0 (shown in Fig. 2 for a cross-shaped RVE).
By splitting the integral (41) one obtains: Assuming traction free boundaries on the shell of the RVE the integral over the shell vanishes. This yields Thus far the presented concepts are not distinguished from homogenization methods applied in three-dimensional continuum mechanics [16]. Now, however, the fluctuationsx i are described in terms of the rod's kinematic unknowns. Therefore, recall that the deformed configuration of a geometrically exact rod is defined as where X C S J denotes the coordinates in the cross-section. Finally, the fluctuationsx i are defined as and their variations are Inserting equation (46) into equation (43) and splitting the boundary ∂B F 0 into periodic pairs by • + and • − one obtains Using the well known relation T i = P i J N J , where T denotes the traction vector, the above equation is reformulated as To satisfy the above equation the following conditions have to be fulfilled The first two conditions are ensured by applying periodic boundary conditions to the field of unknowns. The third and fourth ones are already satisfied since jumps on the boundary are not allowed and the normals are pointing counter-wise on the RVE's periodic boundary.

Micro-to-macro
Having shown that periodic boundary conditions in term of the rod's kinematic unknowns fulfil the Hill-Mandel condition it is now demonstrated how to obtain the averaged macro stress P M . Again, one makes use of Einstein's index notation.
In index notation, the Hill-Mandel condition yields The rod's material points X are given by equation (2). Substituting equation (2) into equation (50) one obtains The integrals are again split into their periodic pairs as Rewriting equation (52) in terms of tractions yields The integral over the periodic boundaries can be rewritten as where is a place-holder for either + or −. The upper boundary n C S of the sum indicates the number of cross-sections laying on the RVE's boundary. Finally, this yields The above statement implies TdA shown in [32,34,35], and since (56) shows that by choosing the boundary of the RVE to coincide with the rod's cross-section (i.e N = D 3 ) one can represent the averaged stress in terms of the rod's stress resultants. Further, it is noted that only the resulting forces n 0 contribute to the stresses in the macro-scale. The moments do not contribute. In the above derivations the scale-separation is assumed to be given.
In the next section use will often be made of the Piola-Kirchhoff stress, the pull-back of the Piola stress:

Numerical results
The mechanical properties discussed in this section result from multiscale homogenization of an exemplary rod lattice, which is not directly motivated by a real-life application. However, there are countless conceivable applications requiring geometrical and material nonlinearities. Consider, for example the presented lattice as an impact absorbing metamaterial in which the consideration of geometrical and material nonlinearities is crucial. The rod lattice considered here is illustrated in Fig. 3. It consists of smallest repeating unit cells with edge length 1 constructed by 8 rods coinciding with the space-diagonals of the unit cell (returning rods with ). The rods intersect at the barycentre of the unit cell. The structure of the unit cell is called body-centred To neglect the fact that multiple rods intersect, and therefore the stiffness at the intersection may vary from the stiffness of a single rod, we assume that the rods are thin, i.e.
For circular cross-sections C 0 takes following diagonal form where A = μA denotes the shearing stiffness and B = E A the tensile stiffness. The parameters C = E I and D = μJ define the bending and torsional stiffness, respectively. Further, A defines the cross-sectional area, I the planar moment of inertia and J the polar moment of inertia. The Young's modulus is denoted by E. It is special to circular crosssections that they do not show any coupled stiffnesses on the off-diagonals in the strain free state. In addition to the tangent stiffness matrix, one needs an appropriately defined yield-function in terms of stress resultants to properly model inelastic rods. Obtaining such a yield-function is not trivial. Herrnböck et al. [32] is fully devoted to the development of a framework to determine yield-functions depending on the cross-sectional shape and its underlying three-dimensional constitutive relation. Here let the inelastic behaviour of the continuous rod be described by the von-Mises yield criterion with the Kirchhoff stress τ and q being the energetic conjugate to an internal hardening variable. The yield stress σ y is set to σ y = 0.45 GPa. The introduced material parameters are taken from Simo et al. [38] and describe standard steel. Together with the circular cross-section, the following yield-function in terms of stress resultants results: where the values with subscript y represent the components of the stress resultants at yield. Note that for the sake of readability the subscript 0 is omitted. The hardening-stresses are denoted by β i (see Eq. 25). They result from a quadratic hardening potential where H is a 6 × 6 diagonal matrix with Note that nonzero hardening is favourable to ensure a stable simulation. It is shown in [32] that the stress resultants at yield can be scaled by the size of the cross-section without noticeable loss in accuracy. In the following the dimensionless exponents in equation (62) take the values α = 2.04, γ = 1.76, δ = 2.09, η = 1.73. The yield limits are set to n 1 y = n 2 y = 0.70θ 2 kN, n 3 y = 1.47θ 2 kN, m 1 y = m 2 y = 0.62θ 3 kNmm, m 3 y = 0.56θ 3 kNmm. The dimensionless scaling ratio is given by θ := d [mm] .
In the following examples, the mechanical response of the RVE undergoing small macroscopic strains (ε < 5%) will be investigated. The examples will demonstrate that, despite assuming moderate strains on the macro-scale, the deformation in the RVE can be large and geometric nonlinearities as well as structural instabilities may appear. The presented results follow from an RVE subjected to periodic boundary conditions.
The aim of this section is to find out if the mechanical response shows a dependency on the number of unit cells constituting the RVE and to find the smallest RVE whose mechanical response is no more sensitive to the RVE's size. Therefore, the size of the RVE is varied. A parameter n RVE is introduced that describes how many unit cells the RVE consists of. Choosing n RVE = 2 results in a single repetition of the unit cell in every space direction. Thus, instead of one single unit cell, the RVE consists of n 3 RVE unit cells (assuming three-dimensional space).

Remark 2 Prescribing a macroscopic deformation to an RVE consisting of thin rods may lead to structural instabilities on the microscopic level, i.e. buckling and bifurcation points may occur. To capture bifurcation points and to follow the energetically favourable post-bifurcation path methods as presented in [39-41] are used. The bifurcation point is searched by an eigenvalue-analysis of the RVE's stiffness matrix. At the instability point the solution is perturbed by the eigenvector corresponding to the smallest eigenvalue. This method ensures that the exact buckling load is found and that the deformation chooses the energetically favourable secondary path. However, this approach is expensive, since an eigenvalue problem must be solved. When considering RVEs of large size (n RVE > 2) a different approach is used.
To capture the favourable secondary path the undeformed configuration of the lattice is slightly perturbed by a linear combination of buckling modes [2,14]. The perturbed lattice does no longer face any structural instabilities. The search of the exact instability point and the correct secondary path is no longer necessary and no longer possible. This approach may be motivated by imperfections in the lattice. If the perturbation is small enough the lattice approaches the behaviour of a perfect lattice without any imperfection.

Stiffness of RVE at zero strain
Before analysing the mechanical response of the RVE subjected to uniaxial strains, its stiffness at zero strain is briefly discussed with the help of the Young's modulus surface [42]. To numerically evaluate the stiffness matrix perturbation techniques such as those presented in Miehe et al. [43] are applied. Figure 4 displays the normalised Young's modulus surface of the introduced rod lattice with varying cross-section's diameter d. The Young's modulus surface visualises the direction dependency of the Young's modulus.
The Young's modulus surfaces for BCC lattices are well known from literature [44,45]. Their interpretation is recaptured as follows. Independent of d the behaviour of the

Uniaxial tension
In this section the mechanical response of the RVE subjected to uniaxial tension is considered. Let the macroscopic deformation gradient be with the load parameter λ. In the following the tensile component of the Piola-Kirchhoff stress S 11 is plotted as a function of λ. The purely elastic case as well as the inelastic case is considered. Furthermore, the curves are shown for different diameters d of the rod's cross-section resulting in different stiffnesses. Since the course of the curves is qualitatively identical for varying d the stress-strain relation is only represented for few values of d in Fig. 5. All plots in figure 5 depict a linear relation of S 11 and λ in the elastic region. For the inelastic cases (dotted lines) the start of yield is sudden and a considerable loss in stiff- ness is clearly visible. The results differ quantitatively but not qualitatively with respect to the rod's cross-sectional dimension. The present load case does not result in any structural instabilities. Further, the size of the RVE (parametrized by n RVE = 1, 4, 6) does not affect the results. It is evident, but still mentioned for completeness, that the homogenized tensile stress S 11 increases with increasing cross-sectional diameter d.
We now investigate the onset of yield in terms of the rod's diameter d and the RVE's size n RVE (for sake of visibility we choose n RVE = 1, 4, 6). Therefore, the load parameter λ at yield is visualised as a function of d in Fig. 6.
Here, it occurs that onset of yield (subsequently marked by a cross) is not a function of the rod's diameter d for the case of uniaxial tension. The structure begins yielding at very small values of λ (less than 0.01). In this state the structure has not yet largely deformed. Considering the tensile stress S 11 at the onset of yield one observes a nonlinear increase (quadratic in d). Since the case of uniaxial tension does not show any instability in the structure the shown plots do not reveal any dependency on n RVE , which can be regarded as a validation of the homogenization framework. An RVE consisting of only one unit cell is sufficient to show the macroscopic material behaviour for the case of uniaxial tension.

RVE under compression
The case of uniaxial compression yields more noteworthy results. This is due to the appearance of buckling. It is possible to subdivide the deformation state into a pre-and with λ ∈]0, 1]. Figure 7 displays S 11 as a function of λ for the purely elastic case. Further, the number of unit cells within the RVE is varied. In all curves, one can discern a kink in S 11 .
At this point the RVE shows a structural instability. Further loading leads to buckling of the structure and tremendous loss in stiffness. The onset of buckling varies with the crosssectional diameter d. With increasing d both the compressive stress, and the loading parameter λ at buckling increases. In contrast to tensile loading, the size of the RVE (parametrized by n RVE ) clearly impacts the simulation outcome. However, with increasing n RVE the discrepancy between the curves reduces and reaches a converged state for n RVE = 6. In the pre-buckled region n RVE does not show any influence on the results. A more complex course of S 11 with respect to λ stands out for the inelastic case in Fig. 8. Here the lattice not only undergoes buckling but also starts yielding. Both the onset of buckling as well as the onset of yield are dependent on the cross-sectional size. For thin rods (≈ d < 0.05) the structure first buckles and starts yielding afterwards. Thicker rods show the opposite tendency. There the material first yields and the structure buckles afterwards. However this observation is also strongly dependent on n RVE . Throughout all cases the stiffness strongly changes when passing the onset of yield and the onset of buckling. With increasing n RVE a converged state is reached. For n RVE = 6 a sensitivity towards the RVE's Then, instantaneously S 11 remains nearly constant. The structure is buckling and can not support more load. The strain at buckling differs for different cross-sectional diameters d. In general, for increasing n RVE , λ and the corresponding stress S 11 at buckling decreases and finally converges to a constant value size is no longer discernible. Thus, the smallest RVE able to return the homogenized properties for the case of uniaxial compression consists of 6 × 6 × 6 = 216 unit cells. This applies both for the elastic and inelastic case. Recall that for n RVE ≤ 2 the undeformed lattice has no imperfections, whereas for n RVE > 2 the undeformed lattice is slightly perturbed enabling buckling (see remark 2). In perfect lattices the plastification propagates simultaneously throughout the lattice. In perturbed lattices, the plastification may begin locally in the lattice. One have to be aware that this may influence the mechanical response of the lattice, for instance by showing a smoother transition between the elastic and plastic region. Figure 9 shows the load parameter λ at onset of buckling and yield and the corresponding compressive stress S 11 as a function of the cross-sectional diameter d for different n RVE . For sake of visibility, here we choose n RVE = 1, 2, 4, 6. The values of λ and S 11 at buckling and yield reveals the following. The onset of buckling (marked by a circle) clearly shows a dependency on d and n RVE . The load parameter λ at buckling increases with d but decreases with n RVE . Besides this, S 11 at onset of buckling and yield decrease with increasing n RVE . To further discuss the results, the following distinct cases are introduced. In the first case, the compressive load parameter at onset of yield is smaller than the load parameter at onset of buckling, that is, the structure yields in the prebuckled state (e.g. d = 0.07mm, n RVE = 1). The second case is characterized by yielding in the post-buckled state (e.g. d = 0.035mm, n RVE = 1). A third case is introduced, where buckling and yield takes place simultaneously (e.g. d = 0.07mm, n RVE = 6). It strongly depends on d and n RVE which case applies. The appearing case also impacts the Even with moderate prescribed macroscopic deformation the deformation and rotation of the constituting rods are large. Additionally, it is visible that the deformation of the unit cells varies for n RVE = 1 and n RVE = 2 (the planar view clearly displays this phenomenon). The rods in the RVE buckle in a wave-like manner. The wavelength of the buckled rods in the RVE with n RVE = 1 is half the wavelength of the rods in the E2 E2 E1 E1 E3 E3 Fig. 10 Buckled configuration of an RVE with n RVE = 1 and n RVE = 2, respectively, subjected to uniaxial compression. Despite moderate macroscopic deformation, the displacements and rotations on the constituting rods is remarkable large. Furthermore, the unit cells do not deform equally when n RVE = 1 or n RVE = 2 RVE with n RVE = 2, which results into a higher curvature of the rods for n RVE = 1 than for n RVE = 2. The observation is made that the prescribed load case results in a micro-buckling wavelength which is dependent on the RVE's size.
In Fig. 11 it is shown that despite a significant difference in the mechanical response there are hardly any variations to be seen in the deformed configurations of the RVE when considering elastic or inelastic material: the only visible difference is the slightly smoother bending of the rods for the elastic case, compared to the inelastic case. To make the plastified regions of the rods visible the norm of the plastic curvature k p 0 is coloured. It is obvious that no plastified regions appear when assuming purely elastic material. When considering inelastic material the plastified regions concentrate on the most bent regions of the rod. The location of these regions depends on the choice of n RVE . The highest plastic curvature lies in the middle of each rod for n RVE = 1, while for n RVE = 2 most plastification takes place at the joints.

RVE under shear
Finally, the RVE undergoes uniaxial shear. The macroscopic deformation gradient is given by In analogy to the case of uniaxial compression, shear may lead to buckling in the micro structure. The shearing stress component S 12 is presented as a function of λ for the purely elastic case in Fig. 12. In analogy to the case of uniaxial com- pression one detects a kink in S 12 . This point denotes the start of buckling. The change in stiffness is less pronounced than for the case of uniaxial compression. This results from the fact that under shear some rods in the RVE are loaded in the tensile direction. These rods do not buckle. Whereas for the case of uniaxial compression all rods are loaded compressively, which leads to buckling in all of them. The onset of buckling is both a function of the rod's cross-sectional dimension d and of the RVE's size (parametrized by n RVE ). However, here the curves already converge for n RVE = 2. This is in contrast to the case of uniaxial compression, where convergence is only reached for n RVE = 6.
Let us now consider the case of uniaxial elastoplastic shear in Fig. 13. Again, S 12 is displayed as a function of the load parameter λ for different d and n RVE . The curves reveal two kinks that are more or less pronounced. One kink denotes the buckling at the structural instability. A second kink denotes the onset of yield within the underlying rods. In analogy to the case discussed in Sect. 5.3 yield may take place in the pre-or post-buckled region.
The plastified structure shows a considerably lower stiffness than the elastic structure. Compared to the purely elastic case (Fig. 12) convergence is only attained for n RVE = 6. Thus, the smallest RVE able to show the homogenized properties resulting uniaxial shear consists of 6 × 6 × 6 = 216 unit cells. Figure 14 depicts the load parameter λ and its corresponding shear stress S 12 at onset of buckling (marked by a circle) and yield (marked by a cross) as a function of d. Again, for sake of visibility we restrict ourselves to four different values for n RVE . In analogy to the case of uniaxial compression, the load parameter λ at onset of buckling shows a strong dependency on both d and n RVE . The shearing load at buckling increases with d and decreases with n RVE . With increasing n RVE the stress at both onset of yield and buckling decreases. Three distinct cases are introduced: pre-buckled yield, postbuckled yield and simultaneous yield and buckling. As for the case of uniaxial compression the first case is defined in that the load parameter λ at buckling is higher than the shear at yield (e.g. d = 0.05, n RVE = 1). In the second case this relation is opposite (e.g. d = 0.035, n RVE = 1). In the third case both values coincide (eg. d = 0.07, n RVE = 6). The choice of d and n RVE impacts the appearing case.

The inluence of n RVE
The results shown in previous sections revealed that the homogenized properties of the RVE are dependent on the number of unit cells within the RVE, i.e., the size of the RVE (parametrized by n RVE ). This dependency only appears if the RVE faces structural instabilities. The homogenized response in the pre-buckled state is independent of n RVE . The same occurs for load cases that do not lead to buckling, such as uniaxial tension. There the results do not depend on n RVE . These findings are consistent with those in [15,17,41,46]. One can conclude from this that the choice of an appropriate number of unit cells within the RVE is of critical importance for the determination of macroscopic mechanical properties. Recalling Fig. 10  ture the stored bending energy in high wavelength is less than in small wavelengths. Thus, buckling in high wavelengths takes place at lower stress and strain level. The reason for different wavelengths in different RVEs is that the prescribed periodic boundary conditions force the buckling wavelength to fit into the size of the RVE. By increasing the RVE's size, convergence in the homogenized properties may be reached. This brings up the following train of thoughts. The presented homogenization framework requires that the scale-separation (34) must be fulfilled. This means that the scale of the macro-problem must be larger than the scale of the micro-problem. It has been shown that here the homogenized properties of the micro-problem converge for n RVE ≥ 6, which is now defined as the smallest possible size the RVE must take to show no more sensitivity towards n RVE . This RVE can be used to derive a material model for the meta-material in a sequential homogenization framework. Alternatively it may be used in a concurrent homogenization framework such as in the FE 2 method. Then one must make sure that the size of the macroscopic problem is significantly larger than the size of the RVE, which means that the macroproblem must consist of way more that 6 × 6 × 6 unit cells. Especially in applications emerging from additive manufacturing this may not be the case. In this case, the use of an RVE with n RVE = 6 would lead to the violation of scaleseparation. Resolving the whole macro-structure with fully nonlinear beam elements may be a convenient alternative.

Summary
This contribution introduces a multiscale homogenization framework combining well known three-dimensional continuum mechanics on the macro-scale with geometrically exact elastoplastic rods on the micro-scale. This framework allows the modelling of lattice micro-structures undergoing large deformations. Additionally, the theory of geometrically exact rods has been presented briefly for both the elastic and the inelastic case. A multiscale homogenization framework was deduced from its well established three-dimensional setting. Considering the example of a specific rod lattice homogenized mechanical response for the cases of uniaxial tension, compression and shear were analysed with respect to varying cross-sectional dimension d and RVE's size. The results clearly demonstrated the appearance of particular properties resulting from the slender character of the rod structure: the appearance of structural instabilities (buckling) and yield. The influence of the rod's dimension and the RVE's size on its homogenized properties were discussed in a comprehensive way.
The results and findings of this work will be used to drive further research. In this work rods without imperfections were considered. However, this is by no means representative of practical applications. For example in the field of additive manufacturing imperfections are present in the lattice. Thus, the implementation of such shall be included in future work. These may be included for example by perturbing the geometry of the lattice (as done in Sect. 5.5 to obtain a buckled configuration), tapering the underlying rods, or by considering heterogeneous properties in the rod's cross-section. A method which would also consider the imperfections appearing during the manufacturing process is to resolve the crystalline structure and therein spatially varying material properties of the cross-section. Imperfections in the crystalline structure would then manifest as changes of the rod's stiffness and inelastic behaviour. This would, of course, require the modelling of crystal plasticity within the cross-sectional level.
The presented homogenization framework can be used to derive an effective constitutive model for different lattice structures by fitting the parameters of an appropriate material model to the data of uniaxial and biaxial strain states or by setting up an artificial neuronal network (as demonstrated in [14]). Additionally, the consideration of inelastic behaviour allows one to find an effective yield-surface. Strategies to do so are presented in [32]. If one is not interested in a sequential but concurrent homogenization the FE 2 method is a promising but expensive tool. Regardless of the chosen method (sequentially or concurrent), first the right size of the RVE must be investigated to make sure that the homogenized mechanical response does not show dependencies on the RVE's size. Further the assumption of scale-separation must be preserved. That means that the macro-scale must be a lot larger than the micro-scale. If this is not the case, the lattice shall better be fully resolved.
Acknowledgements The authors acknowledge support from the German Science Foundation (DFG) within the Collaborative Research Center 814 "Additive Manufacturing". Further, Daniel van Huyssteen is thanked for proofreading the manuscript.
Funding Open Access funding enabled and organized by Projekt DEAL.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.