Mechanical characterization of composite materials with rectangular microstructure and voids

The purpose of this work is to study the mechanical behavior of microstructured materials, in particular porous media. We consider a detailed description of the material through a discrete model, considered as the benchmark of the problem. Two continuous models, one micropolar and one classic, obtained through a homogenization procedure of the material, are studied both in static and dynamic conditions. Furthermore, the internal characteristics of the material, such as the internal scale of the microstructure and the percentage of the voids, are made to vary in order to investigate the mechanical response and to have an exhaustive comparison among the models.


Introduction
Micromechanical modeling of materials has sparked a lot of attention in recent years, as it is well known that many materials contain heterogeneous microstructures that play a dominating role in defining macromechanical behavior. Multiphase fiber and particle composites, ceramic and metal composites, poly-crystals (tungsten-carbide, zinc oxide, aluminia, zirconia), rocks, concrete, masonry structures and granular materials are examples of materials where this happens. These materials show microstructures that occur at a variety of length scales, from meters to nanometers. The micromechanical behaviors across distinct material phases have a profound influence on the response of such heterogeneous solids. A detailed discrete model, considering these materials as Lagrangian systems, is characterized by the expensive computational cost of the numerical model [1,2] and for this, an alternative approach is preferable such as to adopt multiscale procedures for deriving equivalent homogenized continua. The capability of classical continuum mechanics theories to predict such phenomena is restricted, prompting the creation of a slew of novel micromechanical solids theories. There have been several instances where significant differences between theory and experimentation have been noticed. This refers, first and foremost, to stress situations with significant stress gradients. Stress concentrations in the vicinity of holes, notches, and fractures might be instances of such states. The disparity between the classical theory of elasticity and the tests is also noticeable in dynamical situations, such as elastic vibrations with high frequency and short wave lengths, such as ultrasonic waves [3][4][5][6][7][8]. More recently, different kind of continuum models have been proposed, in which the material response is dependent on certain microscale length parameters linked to the existence of inner degrees of freedom and non-local continuum behavior. 'Explicit' non-local behavior refers to the fact that the stress at a given location is influenced not only by the strain at that location, but also by the stresses at nearby locations. This necessitates non-local descriptions that incorporate internal length parameters, 'implicit' non-local description, with internal length and dispersion properties in wave propagation [3,6,9,10], can be also obtained by adding extra degrees of freedom [11][12][13][14]. To explain the micromechanical behavior of materials, non classical and non-local theories based on higher-order gradients have been developed as well as micromorphic continua [9,15]: the micropolar theory [16][17][18] is an example of 'implicit' non-local descriptions, in which the rotation of the material point (microrotation), different from the classical local rigid rotation (macrorotation), is introduced as an additional kinematic parameter. This theory has been successfully used to a variety of composites with both periodic and random microstructure [19][20][21] as well as the couple stress theory [22,23], in which the microrotation and the macrorotation coincide. Non-local 'explicit' theories have been adopted for the modeling of heterogeneous materials [24,25], such as composites with periodic microstructure [26] or random composites [27][28][29]. Given the extensive application of composite materials in many engineering structures, strong is the interest of researcher about the thermomechanical behavior of anisotropic materials [30,31]. In particular, techniques based on homogenization have been widely exploited to study the mechanical behavior of materials: for example, the effects of micro-fractures and contact simulations on the macroscopic response of deformed composites, fiber-reinforced materials [32] and the failure and the damage of periodic masonry assemblies [33,34] and of ceramic and advanced materials [35][36][37]. There are also many examples in which materials with voids [38][39][40] have been studied referring to Cosserat theory, such as foams [41][42][43][44], granular and geo-materials [45][46][47][48][49][50][51], cellular solids [52][53][54] and different kind of structures such as rods [55] and shells [56][57][58] or even multifield models with affine microstructure [5,7]. In order to extend the numerical investigations about the microstructured composite materials made of rigid particles and elastic thin interfaces [59][60][61], the goal is to test the efficiency of the homogenization technique adopted [62,63], to describe porous materials, not yet analyzed with this approach by the authors themselves [64][65][66]. The study is conducted both in static and dynamic conditions, varying the internal characteristics and then proposing a comparison between a discrete model and two continuous models, one micropolar and one classic. This is how the article is structured: we start from a brief description of the micropolar theory, especially referring to the 2D case (Sect. 2), then the homogenization technique used for the homogenization of the materials is mentioned (Sect. 3), the simulations will be illustrated both as regards the numerical implementation and both as regards the results (Sect. 4), and, finally, the conclusions will be remarked (Sect. 5).

Micropolar theory
Extra degrees of freedom are considered in micropolar theory, in particular, the microrotation of the material particles is considered together with the displacements. Therefore, the displacement field of a material particle in a micropolar continuum is characterized by displacements and rotations. In a 3D framework, there are three displacement and three rotation components, which reduce to two displacement and one rotation components in the 2D case. In the present study, we refer to 2D problems in a linearized kinematical framework, where, the displacement components being gathered in the vector u = u 1 u 2 ω , u 1 and u 2 being the displacement components, and ω the micro-rotation. Note that the variable ω is different from the macro-rotation θ , defined as the skew-symmetric part of the gradient of displacement. The strain vector is: ε = ε 11 ε 22 ε 12 ε 21 κ 1 κ 2 , where ε i j (i, j = 1, 2) are the normal and shear strain components, and the microcurvatures are indicated by κ 1 and κ 2 . Differently from the classical continuum, the strain components are not reciprocal ε 12 = ε 21 . The vector σ = σ 11 σ 22 σ 12 σ 21 μ 1 μ 2 collects the stress components, where σ i j (i, j = 1, 2) represents the normal and shear stress components, and μ 1 , μ 2 are the microcouples. Due the presence of microcouples, differently from the classical continuum, the shear stress components are not reciprocal: σ 12 = σ 21 . The couple stress components μ 1 , μ 2 allow to satisfy the moment equilibrium of the micropolar body.
The operator D is defined below: as consequence, is possible to write the kinematic compatibility relation between the vectors u and ε: Using the Hamilton's principle, the equilibrium of the body can be expressed as: where the term K indicates the kinetic energy and is the total potential energy, which is defined by the sum of the strain energy U and the potential of external loads V : The variation of the kinetic energy is: h is the thickness of the body which can be assumed unitary and m is mass matrix defined as: where ρ is the material density and J c represents the rotary inertia of the material point. The variation of the strain energy is written in the form: and using Eq. (1): Below, the variation of the potential of external loads is reported: where the body and surface forces are represented by the vectors b and t, respectively. At this point, the micropolar anisotropic constitutive equation can be reported: where As hyperelastic media are considered, the stiffness matrix of the discrete system is symmetric, this yields that the continuum constitutive matrix C is symmetric (C ∈ Sym). Thus, and they are a fourth-order, third-order and second-order tensors, respectively [67]. On the basis of these assumptions, the Hamilton's principle can be written as: 3 Continuum model

Homogenization
In this work, materials made of a rectangular microstructure and with the presence of voids ( Fig. 1), representing to a porous medium, are taken into account: the blocks of the microstructure are assumed to be constant while the voids percentage v varies. In this case, two different size of the internal microstructure and three different percentage of the voids have been considered. In order to test the homogenization procedure for these kind of microstructured materials, we introduce a scale parameter s corresponding to the block lengths l 1 , l 2 and consider two cases: s = 1, with l 1 = 0. Three different voids size are considered, as reported in Fig. 2 where the representing volume element (RVE) is depicted. To apply the homogenization method to periodic assemblies, the RVE must be identified, as the volume element, made of the minimal number of elements and joints necessary to fully define the behavior of a discontinuous and heterogeneous material. The blocks are assumed rigid respect to the stiffness of the contact interfaces, and the elasticity of the microstructure is assumed to be concentrated in the joints between the blocks. For the interfaces, a linear elastic constitutive law is considered.
As homogenization criterion, an energy equivalence criterion is adopted, based on the a generalized version of the so-called Cauchy-Born rule [68][69][70][71][72][73][74] that is a kinematic correspondence map between discrete and continuous fields proposed in [62,63]. Non-classical continuous models of discrete systems of various kinds have been derived using revised discrete-continuum correspondence maps using this approach, always accounting for the internal size effect of the microstructure [75][76][77]. The feasibility of giving a macroscopic description for heterogeneous media has been widely studied within the framework of the homogenization or coarse graining theories [11]. The power equivalence approach provides a straightforward way to determine all the constitutive parameters of the continuous medium, the components of the matrix C in Eq. 11, by establishing the elastic constants and the geometry of the discrete system. Finally, the constitutive matrix for the Cauchy model is derived from the previous Eq. 11 using the procedure in [78] as: It should be remarked that Cauchy model does not keep memory of the internal length of the microstructure.

Finite element model
The continuum models have been solved with the finite element method that is based on the approximation of nodal displacements: The current implementation was carried out utilizing an in-house MATLAB code [79]. The kinematic displacement vector is arranged as follows: there are 12 degrees of freedom overall (3 per node). The vector N, which gathers the Lagrangian shape functions, makes up the matrix of shape functions: By replacing the above expression of the displacement vector, in the Hamilton principle, the kinetic energy becomes: The mass matrix reads: The internal work takes the form: where B = D N , thus the element stiffness matrix is: which must be integrated using a 2 × 2 Gauss integration for normal components and micro-couples, while reduced integration is used for shear components. Finally, the potential of external forces is: where F is the global vector of volume and surface forces.

Numerical simulations
It is emphasized that the constitutive matrix of the continuum model, C (Eq. 11), for all the cases here considered, with rectangular microstructure, is centrosymmetric and moreover orthotropic (Eqs. [22][23][24][25][26][27]. Only the submatrix D takes into account the internal scale of the microstructure of the material. In particular, due to the central symmetry, the submatrix B = 0 and there is no coupling between normal and shear stresses/strains with curvatures/micro-couples, moreover, due to the absence of dilatancy effect in the discrete original model, no Poisson effect is present. A rectangular panel, of width L x = 12 μm and height L y = 11.6 μm, in the 2D plane stress case, both in static and in dynamic conditions, is studied. As a comparison, a lattice model as benchmark has been carried out in ABAQUS , and two continuum models, the Cosserat model and the Cauchy model, are compared with the latter. For the static case, the panel is subjected to a vertical load p applied at the top distributed in a foot print equal to a = L x /4 (Fig. 3). In the dynamic case, free vibrations are considered.

Scale s = 1
The constitutive matrix for the microstructure with the scale parameter s = 1 and voids size v = 0.25 is: In Fig. 4, the vertical displacement component is shown and due to the symmetry of the problem only half of the panel is reported (only the vertical component of displacements is reported, as the most significant of the case study). It can be noted that the contour plot of the micropolar continuum is in line with the discrete system. There are not evident differences among the Cauchy and Cosserat models, even if, this last is able to better distribute the level curves in comparison with the classical model. These differences have been shown to strongly increase when the load tip decreases [67].
The constitutive matrix for the microstructure with the scale parameter s = 1 and voids size v = 0.5 is: An increase in the size of the voids causes a decrease in the elastic contact among the rigid blocks of the microstructure, and this translates into a reduction in the elastic constants of the submatrix A, hence, in a loss of the stiffness. Differently, for the terms of the submatrix D, which depends, not only on the interface stiffness, but also on the internal size, there is an increment of the value of the term D 11 , while the term D 22 decreases respect to the previous case. The lower stiffness led to an increase in the maximum displacements, and the classical model gives reliable results only for the static case (Fig. 4), although different from the discrete model which results are in agreement with the micropolar ones.
Finally, the constitutive matrix for the microstructure with the scale parameter s = 1 and the voids size v = 0.75 is: In the latter case, there is a further reduction in the elastic constants of the matrix A, while as regards the terms of the sub-matrix D, the term D 11 decreases with respect to the case in which the dimensions of the voids are v = 0.5, but it is greater than the case v = 0.25, the term D 22 is smaller in both cases. This implies an increase in the maximum displacement, and both continuous models provide a good characterization of the discrete model behavior (Fig. 4).
Regarding the dynamic analysis of the microstructure with a voids size v = 0.75 (Fig. 5), it can be observed as for the first three modes, both the classical and the micropolar continuum are able to detect the displacement plots. However, for this last model the maximum error, in the natural frequencies evaluation (Table 1), is around 3%, whereas for the Cauchy model the error is very high compared to the discrete system and only the axial vibration mode exhibits a good result (the errors reported are the relative errors of the continuous models computed with respect to the discrete system). Furthermore, only the micropolar model is capable of detecting the correct dynamic behavior (mode shape) of the fourth mode.
From the dynamic analysis, the growth of the voids size highlights that the Cauchy model has worse results than the micropolar one, when compared to the discrete system and also with respect to the previous simulations with smaller voids percentage. Such large differences are noticed in terms of frequency (Table 1) and mode shapes (Fig. 6). The fourth mode of the Cauchy model, even if different from the previous case, does not match with the discrete one, and only the Cosserat model satisfactory represents the behavior of the benchmark.
For the voids size v = 0.75 case, the micropolar model is still efficient as regards the representation of the displacement fields of the discrete model (Fig. 7). The estimation of frequencies for both continuous models deteriorated as the size of the voids has become larger (Tab. 1), even if the largest error threshold of the Cosserat model is about the 3.5%, and as a result, it continues to be reliable.  The scale parameter s is taken into account by the submatrix D, and, specifically, it can be seen that D s=1 ii ≈ 4D s=0.5 ii (i = 1, 2). There is an improvement in the correspondence between the results of the discrete model and the micropolar model when the internal scale of the microstructure is reduced, as shown in the contour plots of the vertical displacements for the static case (Fig. 8).
The constitutive matrix for the microstructure with the scale parameter s = 0.5 and voids size v = 0.5 is: Also in this case, the displacement field of the discrete model is well represented by the two continuous models, and the reduction in the scale favors the correspondence between the models; moreover, there is an increase in the maximum vertical displacement compared to the previous case (Fig. 8).
With the decrease in the scale, as seen in the previous two cases, there is still a better correspondence between the continuous models and the discrete system (Fig. 8).  As regards the dynamic analysis, it is clear that only the micropolar model gives satisfactory results both in terms of plots (Fig. 9) and frequencies ( Table 2). In particular, it has been noticed that the errors of the present results are smaller compared to the previous case with the scale parameter s = 1, in fact the the maximum error is around the 0.70%.
Also for case where the voids size is v = 0.5, the maximum error is reduced when the microstructure scale is reduced, which improves the natural frequencies evaluation for the micropolar model (Table 2); however, also in this case, the classical model is not able to predict the dynamic behavior of the benchmark and this is related to the orthotropic character of the material as well, except for the axial mode (Fig. 10).
Finally, the dynamic study reveals that decreasing the scale parameter improves the outcomes and in general, it can be stated that with the increase in the voids size the evaluation of the frequencies presents some worsening; however, the Cosserat model still offers a good estimate of the dynamic response, the maximum error being around the 2% (Fig. 11 and Table 2).

Conclusions
In this paper, we have dealt with the numerical study of an orthotropic material, with the presence of voids, therefore assimilable to a porous medium. The comparisons among the discrete systems, assumed to be a micromechanical model of the material and supposed as the benchmark of the problem, and two continuous models were carried out through static and dynamic analyses. When the contact area between the blocks becomes significantly smaller than their geometric size, there is a worsening of the capability of the continuous model to represent the lattice system behavior. The differences between the Cosserat and Cauchy solutions arise as the scale ratio, understood as the ratio between the characteristic length of the macrostructure and the internal length of the material microstructure, increases because the Cauchy model cannot account for scale effects [80,81]. The work highlights the possibility of describing porous media through this homogenization procedure, not yet previously detailed for these kinds of materials, and how the introduction of a new internal characteristic of the material, the size of the voids, involves a deterioration of the results compared to nonporous materials for the classical model [65]. As widely shown in earlier works [59,62,81], for orthotropic media the inadequacy of the a classical continuum arises even when the actual value of the length parameter of the internal microstructure can be considered small compared to the macroscopic characteristics. Only in the isotropic case, the classical model can provide results comparable with the results of the micropolar model when the internal length decreases. The numerical results show that the Cauchy model is not sufficient to adequately characterize problems in which material internal lengths play a significant role and the case of porous media analyzed in this work also contribute to highlight, both in the static and the dynamic case, this aspect. In this work, it has been particularly highlighted, from the dynamic analyses, not already performed in earlier papers, that the micropolar model, accounting for material size and anisotropy [67], makes an important contribution and allows us to accurately describe the mechanical behavior of the discrete system adopted to finely represent the material behavior.