- Email: [email protected]
A substructure-based homogenization approach for systems with periodic microstructures of comparable sizes
Accepted Manuscript A Substructure-Based Homogenization Approach for Systems with Periodic Microstructures of Comparable Sizes Yang Chen, Leiting Dong, Bing Wang, Yuli Chen, Zhiping Qiu, Zaoyang Guo PII: DOI: Reference:
S0263-8223(16)31940-7 http://dx.doi.org/10.1016/j.compstruct.2017.01.050 COST 8184
To appear in:
Composite Structures
Received Date: Revised Date: Accepted Date:
24 October 2016 31 December 2016 17 January 2017
Please cite this article as: Chen, Y., Dong, L., Wang, B., Chen, Y., Qiu, Z., Guo, Z., A Substructure-Based Homogenization Approach for Systems with Periodic Microstructures of Comparable Sizes, Composite Structures (2017), doi: http://dx.doi.org/10.1016/j.compstruct.2017.01.050
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
18 January 2017
Prepared for Composite Structures
A Substructure-Based Homogenization Approach for Systems with Periodic Microstructures of Comparable Sizes Yang Chen1, Leiting Dong2, Bing Wang3, Yuli Chen2, Zhiping Qiu2, Zaoyang Guo2,1* 1
Department of Engineering Mechanics, Chongqing University, Chongqing 400044, CHINA
2
School of Aeronautic Science and Engineering, Beihang University, Beijing 100191, CHINA 3
Center for composite materials, Harbin Institute of Technology, Harbin 150080, CHINA
Abstract The classical homogenization method has been widely adopted to capture the effective behaviors of heterogeneous materials. However, when the characteristic length of the microstructure of the heterogeneous material is comparable to the size of the structure, the classical homogenization method is mathematically no longer valid. In this paper, a new substructure-based homogenization approach is proposed to predict the mechanical responses of systems with periodic microstructures of comparable sizes. A substructure element is developed to reconstruct the system with periodic microstructure of comparable size. It is verified that this substructure-based
*
Corresponding author: [email protected]. Address: Institute of Solid Mechanics, Beihang University, 37
Xueyuan Road, Haidian District, Beijing 100191, CHINA.
homogenization approach can accurately predict the mechanical responses of the system. Comparing with the full finite element analysis, the computational scale is dramatically decreased. After that, a simplified substructure element is developed by using less surface nodes in the “full” substructure element. The numerical results show that, with further significantly reduced computational cost, the third-order simplified substructure element can provide a good prediction of the responses of the system with periodic microstructure of comparable size. Keywords:
substructure
approach,
periodic
microstructure,
finite
element,
homogenization
1. Introduction In recent years, the idea of homogenization has been widely adopted to estimate the effective behaviors of heterogeneous materials (e.g., composites, porous materials). One theoretical approach is to use theories of micromechanics to study the representative volume element (RVE) models and develop homogenized constitutive models to predict the effective responses of composites [1-6]. Another numerical homogenization approach is to use finite element method (FEM) to simulate the mechanical responses of the RVE models [7-10]. The RVE models are also widely used in analysis of system with periodic microstructures [11-13]. In the classical homogenization procedure, the obtained constitutive models are usually used to describe the stress-strain relation at the Gaussian integration points in finite element (FE) analysis of systems of heterogeneous materials. This implies that the classical
2 of 35 pages
homogenization approach is at the (macroscopic) material point scale. Therefore, mathematically the RVE size should be sufficiently smaller than the structure size but sufficiently larger than the characteristic length of the microstructure [14]. For example, Cricri and Luciano [15] suggested that the RVE models of cellular materials have at least 10 cells in each direction to get accurate effective parameters. Hence, the classical homogenization approaches based on RVE models cannot be adopted when the characteristic length of the microstructure of the heterogeneous material is comparable to the size of the structure (e.g., Lstruture Lmaterial < 10 ). However, due to the
novel manufacturing technologies such as 3D printing, the idea of integrated design of materials and structures leads to many structures with periodic microstructures of comparable sizes. In this paper, a substructure-based homogenization approach is developed to simulate the mechanical responses of systems with periodic microstructures of comparable sizes. In this substructure-based homogenization approach, the microstructures of the materials are considered as substructures of the structure because they are of comparable sizes. Therefore, in the substructure-based homogenization approach, the homogenization is at the structure scale rather than the material point scale. In the literature, Karpov et al. [16] adopted the substructures to analyze the static properties of finite repetitive structures, and some illustrative examples (i.e., truss bridge, clamped grid and honeycomb structures) are discussed with the substructure method. Li and Law [17] used substructures to reconstruct the system and study the response of a frame under the excitation force with wavelet domain method. Mencik [18] investigated the harmonic force response of
3 of 35 pages
one-dimensional periodic structures by analyzing the substructure with wave finite element method. Boldrin et al. [19] explored the dynamic behavior of gradient composite hexagonal honeycombs with a substructure-based component mode synthesis (CMS) method. Yu et al. [20] also adopted the CMS to study the mode shapes of the global structure with element-by-element model, which is used to update the large-scale structure. In these papers, the substructure concept is considered at the structure level, while in our approach, the substructure concept is also applied to the material level, i.e., the microstructure of the heterogeneous materials. In this paper, the substructure element is developed to describe the microstructure of the heterogeneous materials. The global stiffness matrix of the RVE model is transformed to the elemental stiffness matrix of the substructure element by eliminating the interior degrees of freedoms (DOFs). Using the substructure elements, the simulation results obtained are identical to those from full FEM simulations, but the computational cost is much lower. Comparing to classical homogenization approaches, the substructure-based approach is particularly useful when the characteristic length of the heterogeneous material is comparable to the size of the structure. Moreover, to reduce the computational scale further, the simplified substructure element is proposed and it is shown that the third-order simplified substructure element can well predict the effective behavior of the systems with periodic microstructures of comparable sizes. The structure of this paper is as follows: In Section 2, the substructure element is constructed, and the applications of the substructure-based homogenization
4 of 35 pages
approach are illustrated. Then the simplified substructure element is developed to further decrease the computational scale and it is applied to the analysis of the porous beams in Section 3. After that, some conclusion remarks are given in section 4. 2. Substructure-based homogenization approach
In this section, the substructure element is developed to describe the microstructure of the heterogeneous material and its elemental stiffness matrix is derived
using
the
concept
of
substructure.
Then
the
substructure-based
homogenization approach is applied to some examples to verify its accuracy and efficiency. 2.1. Stiffness of substructure element
Considering an RVE model of the microstructure of a heterogeneous material (an example shown in Fig. 1a), it is meshed and its mechanical response can be simulated using FEM by the following equation: Kd = P ,
(1)
where K is the system stiffness matrix of the RVE model, d denotes the vector of nodal displacements, and P represents the vector of nodal forces. The nodes are further classified as the surface nodes and the interior nodes. Then Eq. (1) can be alternatively written as K ss K is
K si d s Ps , = K ii di Pi
(2)
where subscript s and i denote the two kinds of nodes respectively. Hence, the response of the RVE can be approximated by the surface nodes only as K *ss d s = Ps* , 5 of 35 pages
(3)
where the elemental stiffness matrix of the substructure element is defined as K *ss = K ss − K si K −ii1K is ,
(4)
and the vector of equivalent surface nodal forces is computed as Ps* = Ps − K si K ii−1Pi .
(5)
Therefore, in comparison to the original RVE model, the dimension of the elemental stiffness matrix of the substructure element is reduced significantly due to the elimination of the DOFs of the interior nodes, which leads to much smaller computational scale. For instance, the RVE model of a unit square with a circular hole shown in Fig. 1(a) has 567 nodes and 987 triangular elements, which can be utilized in FEM to simulate its mechanical responses. The corresponding substructure element contains only 100 (surface) nodes (Fig. 1b), in which each edge is divided to 25 segments by 24 nodes on the edge. The number of elements also is reduced from 987 triangular elements to 1 single substructure element. It will be shown next that the substructure element can achieve exactly the same solution as the full FEM simulation of the detailed RVE model. 2.2. Applications of substructure element
In this subsection, the proposed substructure element is adopted to analyze some examples to verify the accuracy and efficiency of the substructure-based homogenization approach. 2.2.1. Beams with periodic holes
First we consider a beam, periodically stacked by 10 unit squares with circular holes, fixed at both ends and applied by a force P = 1 N m in the middle (Fig. 2).
6 of 35 pages
We assume that the thickness of the beam is much smaller than its height and length, and then plane stress is assumed in the analysis. We note that this kind of structures has been widely used in the fields of architecture and aircraft due to its light weight comparing to solid beams [21]. The material’s parameters are set as the Young’s modulus E = 210 GPa and the Poisson’s ratio v = 0.3 . The RVE model of the periodic microstructure is shown in Fig. 1. The corresponding substructure element (Fig. 1b) is developed and its elemental stiffness matrix is derived from the global stiffness matrix of the RVE model using Eq. (4). The triangular elements used in the RVE model (Fig. 1a) are normal first order triangular elements and the linear interpolation function is used for the triangular elements [22]. The beam is then modelled by 10 substructure elements and the finite element (FE) analysis is implemented using MATLAB code. The deflections of the top and bottom surfaces of the beam computed by the substructure FEM code are plotted in Fig. 3 and compared with the full FEM simulation results. The full FEM analysis uses the mesh of the RVE model for each unit. It can be found that the substructure FEM results are identical to the full FEM simulation. Here the full FEM analysis uses 5436 nodes and 9870 triangular elements, while the substructure FEM analysis has only 766 nodes and 10 elements. Therefore the substructure FEM can reduce the computational scale significantly without losing any accuracy of the results. The classical homogenization method is also attempted for this problem. To obtain the effective properties of the porous material, the RVE model in Fig. 1(a) is adopted in the commercial FEM software ABAQUS/standard 6.12 within the
7 of 35 pages
framework of linear elasticity [23]. The periodic boundary conditions, performed by “Equation” type of constraints in ABAQUS 6.12, are applied in the analysis [23]. When a uniaxial tensile load is applied to the RVE model, the effective modulus and the effective Poisson’s ratio of the porous material are computed from the FEM simulation result as Eeff E = 0.515 and ν eff ν = 0.787 respectively. Then the obtained effective properties are applied to the beam model and the computed deflections of the top and bottom surfaces are plotted in Fig. 3. In this figure, the deflections at the midpoints of both surfaces predicted by the classical homogenization method are about 24% larger than the results from the full FEM simulation as well as the substructure FEM simulation. Hence, it suggests that the classical homogenization approach is inappropriate when the characteristic length of the material microstructure is comparable to the size of the structure, while the substructure-based homogenization approach can exactly reproduce the results from full FE analysis. For the deflection of the beam fixed at both ends and applied by a concentrated force in the middle, the solution from the classical Euler beam theory is written as (assuming constant bending stiffness EI through the beam length) 1 Px3 PLx 2 − (0 ≤ x ≤ L / 2) 16 EI 12 y= 3 2 2 3 1 − Px + 3 PLx − PL x + PL EI 12 16 8 48
,
(6)
( L / 2 ≤ x ≤ L)
where L is the length of the beam, I = bh3 12 (b and h represent the thickness and height of the beam respectively) is the inertia moment of the rectangular cross section.
8 of 35 pages
Eq. (6) shows that the deflection is a third-order polynomial function of the position x (i.e., the distance from the end). The deflection of the top surface from the substructure FEM is fitted by Eq. (6) in Fig. 4(a) using the least square method and the fitted effective modulus E fit E = 0.5768 . As observed from the figure, the fitted deflection curve agrees with the one from substructure FEM simulation very well, showing an average difference of 5.57%. Here the average difference is defined as
ε ave =
1 N
2
fi − 1 , ∑ i =1 yi N
(7)
where N is the number of the data points, yi and fi are the original data and fitted data respectively. Comparing the effective modulus obtained from the classical homogenization approach Eeff = 0.515 E and the one fitted by the classical beam theory E fit = 0.5768 E , it can be easily found that the effective modulus from the classical homogenization approach Eeff is 10.71% less than E fit . We note that the classical beam theory is based on the assumption that the length of the beam is sufficiently larger than its height (i.e., L h > 10 ). To further explore the behavior of
this kind of porous beams with the classical beam theory, the beams with 20 and 30 porous units are investigated, i.e., L h =20 and
L h =30 respectively. The
deflections of the top surfaces computed from the substructure FEM are fitted by the Euler beam theory in Fig. 4 (b) and (c). The relative average differences are 3.01% and 1.84% for the beams with 20 pores and 30 pores respectively. The fitting results suggest that the deflection curves obtained from substructure FEM simulations can be well represented by the function derived from the Euler beam theory. It is found that the effective modulus E fit fitted from the two beams’ results are different from the 9 of 35 pages
one obtained from the beam with 10 holes. For the beam with 20 holes, the fitted effective modulus E (20) = 0.7709 E is 33.65% higher than E (10) fit . When the length of fit the beam increases to 30 units (i.e., with 30 holes), the fitted effective modulus (10) E (30) fit = 0.8180 E is even 41.82% higher than E fit . Both fitted effective moduli are
significantly higher than the effective modulus Eeff obtained from the classical homogenization approach, which confirms the conclusion that the classical homogenization approach is invalid when the characteristic length of the material is comparable to the size of the structure. The results also suggest that the effective modulus fitted by the Euler beam theory depends on the size of the structure. This implies that, even though the deflection of this kind of beams can be described by the same deflection function derived from the Euler beam theory, it is still difficult to use the classical beam theory to analyze this kind of beams because the effective modulus cannot be predicted correctly. To explore the application of substructure-based homogenization approach to the structure with unsymmetrical boundary, the porous beam is fixed at only one end (i.e., porous cantilever beam) and is applied to uniformly distributed load q = 1 N m 2 (Fig. 5). The deflections of the top and bottom surfaces of the beam computed by the substructure FEM code, the full FEM simulation and the classical homogenization method are plotted in Fig. 6. As can be observed from this figure, the deflections computed from the substructure FEM code are still superposed with the full FEM results, which indicates again the excellent predictions of the substructure approach. However, the classical homogenization still overestimates the deflections up to 60.88%.
10 of 35 pages
A fitted modulus E Cantilever = 0.8137 E is computed by fitting the computed deflections fit with the theoretical results of the classical Euler beam theory. The average difference is about only 3.39%, which indicates an excellent fitting. We note that the fitted modulus from the cantilever beam is 41.07% bigger than the one from the beam fixed at both ends and 58% bigger than the effective shear modulus computed from the homogenization method. This suggests that the effective modulus fitted by the Euler beam theory depends on not only the size of the structure but also the loading condition of the structure.
2.2.2. Honeycomb beam Honeycomb structures have also been widely used in many engineering fields because of their advantages such as light weight, high specific stiffness and strength [24]. To verify the accuracy and efficiency of the proposed substructure-based homogenization approach, a honeycomb beam with periodic microstructure is constructed and applied by a force P = 1 N m in the middle (Fig. 7). Both ends of the beam are fixed. In the analysis of this structure the same material parameters are assumed, i.e., E = 210 GPa and ν = 0.3 . Substructure element, developed from the RVE of honeycomb beam, is also used to reconstruct the honeycomb beam and predict the response of the beam. Fig. 8 presents the deflections of the top surface and bottom surface of the honeycomb beam, computed by the full FEM and substructure FEM. It is found that the results from the substructure-based homogenization
11 of 35 pages
approach and the full FEM simulations are identical. This suggests again that the substructure-based homogenization approach can be applied to predict the responses of honeycomb beams with periodic microstructure of comparable sizes. Similarly, the effective mechanical properties can be obtained by applying the classical homogenization method to the RVE model using ABAQUS 6.12. The computed effective properties are Eeffx E = 5.1638 × 10−3 and ν effxy ν = 3.613 . We note that only the mechanical properties along the length direction (i.e., Eeffx and
ν effxy ) are desired in our analysis according to the classical Euler beam theory, even though the RVE model of the honeycomb beam is anisotropic. The computed deflections using the classical homogenization method are also plotted in Fig. 8. It is shown that the deflections computed by the classical homogenization method are about 58% smaller than those from the full FEM simulations. Therefore we confirm that the classical homogenization approach cannot be applied to systems with microstructures of comparable sizes. The deflection curve computed from the substructure FEM simulation can be fitted by Eq. (6) from the classical Euler beam theory. The fitting result is illustrated in Fig. 9, which suggests a close fitting with an average difference of 6.72%. With the classical Euler beam theory, the fitted effective modulus is E xfit E = 2.8401× 10−3 , which is 0.45 times smaller than that computed from the classical homogenization approach. As pointed out before, the main difficulty to use Euler beam theory here is to predict the effective modulus of the heterogeneous material when the characteristic length of its microstructure is comparable to the size of the structure.
12 of 35 pages
3. Simplified substructure element
By eliminating the DOFs of the interior nodes, the substructure element approach can reduce the computational scale significantly. However, because the surface nodes of one substructure element are connected via interior nodes, the elemental stiffness matrix of a substructure element in Eq. (4) is usually a full matrix. For large scale system (with a large number of substructure elements), the system matrix will still be a sparse matrix even though every elemental stiffness matrix is full. Nevertheless, further reducing the number of DOFs of each substructure element will improve the computational efficiency significantly. This section aims to develop a simplified substructure element to further decrease the computational scale, but keeping the prediction accurate enough. In order to reduce the computational scale, the number of nodes on the surface (or boundary) of the substructure element is reduced. To implement that, some nodes are chosen at the surface of the “full” substructure element for the simplified substructure element. This can usually be done by keeping some of the surface nodes in the “full” substructure element and eliminating others. We may also choose a set of new locations on the surface as the surface nodes of the simplified substructure element as long as they can be used to interpolate the displacement field of the surface effectively. The surface nodal displacements of the “full” substructure element are interpolated from those of the simplified substructure element. To obtain an effective interpolation and a well-posed elemental stiffness matrix of the simplified substructure element, the selected surface nodes of the simplified substructure element should usually be evenly distributed on the surface
13 of 35 pages
area. The elemental stiffness matrix of the simplified substructure element can then be derived from that of the full substructure element as follows. First denote the nodal displacement vectors of the simplified and “full” substructure elements as u and u′ respectively. The vector u′ is interpolated by u as u′ = Nu ,
(8)
where N is the interpolation function matrix. The elastic strain energy of the substructure element can be represented as 1 1 W = uT Ku = u′T K ′u′ , 2 2
(9)
where K and K ' denote the elemental stiffness matrices of the simplified and the “full” substructure elements. Based on (8) and (9), the elemental stiffness matrix of the simplified substructure element is obtained as K = NT K′N .
(10)
Similarly, we can derive the equivalent nodal force vector P of the simplified substructure element from the nodal force vector P ′ of the “full” substructure element using the following principle of virtual work PT du = P′T du′ .
(11)
Here du and du′ are corresponding virtual nodal displacement vectors. The equivalent nodal force vector P of the simplified substructure element is then obtained as P = N T P′ .
(12)
After the elemental stiffness matrix of the simplified substructure element is created, the analysis can be carried out with a highly efficient method. Let’s consider
14 of 35 pages
the example of the beam with 10 pores again (Fig. 2), in accordance to (10), the “full” substructure element (originally with 100 surface nodes) is simplified to first-order, second-order and third-order elements, which have 4, 8, and 12 surface nodes respectively. Similarly, FE analysis implemented by MATLAB codes is carried out to study the behavior of the porous beam using the simplified substructure elements. Fig. 10 demonstrates the deflections of the bottom surface of the beam computed using the 3 types of the simplified substructure elements as well as the “full” substructure element. It can be found that, comparing with the result from the analysis using the “full” substructure element, the deflections predicted by the first-order and second-order simplified substructure elements underestimate significantly the beam’s deflection (the average differences are 20.51% and 17.02% respectively). The prediction by the third-order simplified substructure element, however, reveals a good approximation of the response of the beam, which only underestimates the deflection for only about 4.00%. The detailed deformed shapes of the elements are also explored. In Fig. 11, the shapes of two deformed third-order simplified substructure elements at the middle of the porous beam (the concentrated force is applied at the top joint node) are plotted against those of two “full” substructure elements (the deformation is magnified by 5 × 109 times). It is clear that the deformations computed from the two types of elements are very close, which implies the third-order simplified element can well capture the deformation of the porous beam. The first-order, second-order and third-order simplified substructure element are also employed to study the example of the porous cantilever beam, and the defections
15 of 35 pages
of bottom surface computed from the 3 types of elements as well as the “full” substructure element are plotted in Fig. 12. The deflections predicted by the first-order and second-order simplified substructure elements also underestimate the beam’s deflection, showing the average differences of 16.88% and 6.88% respectively. Again the prediction by the third-order simplified substructure element demonstrates a very good approximation of the response of the beam, which only underestimates the deflection for about 2.22%. Comparing to the “full” substructure element with 100 surface nodes, the third-order simplified substructure element only has 12 nodes, which implies nearly 90% reduction of the number of DOFs. The accuracy of the prediction by the third-order simplified substructure element is acceptable with considerably reduced computational cost. Therefore, it is shown that the simplified substructure element, which uses much less surface nodes than the “full” substructure element, has the ability to accurately approximate the responses of structures with periodic microstructures of comparable sizes using significantly reduced computational resources. 4. Conclusion Remarks
Although the classical homogenization approach has been widely used to predict the effective behaviors of heterogeneous materials in the literature, it is demonstrated in this paper that, when the characteristic length of the material is comparable to the size of the structure, it cannot correctly capture the mechanical responses of the structure. This is reasonable because the classical homogenization method is mathematically
16 of 35 pages
invalid in this situation. To predict the mechanical behavior of systems with periodic microstructures of comparable sizes, a new substructure-based homogenization approach is proposed in this paper. The substructure element is developed based on the RVE model representing the microstructure of the material. The elemental stiffness matrix of the substructure element is derived from the global stiffness matrix of the RVE model by eliminating the interior DOFs. Several numerical simulation results show that this substructure-based homogenization approach can accurately predict the mechanical responses of the system. Comparing to the full finite element analysis, the computational scale of the substructure-based approach is dramatically decreased. The classical beam theory is used to fit the deflections of the beams and it is found that the fitted effective modulus depends on the size of the structure (that is, it is not a material property). It shows the difficulty of applying the classical beam theory directly to the structure with microstructure of comparable size. After that, a simplified substructure element is developed by using less surface nodes in the “full” substructure element. The numerical results show that, with further significantly reduced computational cost, the third-order simplified substructure element can provide a good prediction of the responses of the system with periodic microstructure of comparable size. Acknowledgement
The financial support from NSFC (No. 11272362, 11572020, 11332013, CHINA), Chongqing Science and Technology Commission (cstc2013jcyjjq50003), and Beihang University startup grant is greatly appreciated. References
17 of 35 pages
[1] Hill R. CONSTITUTIVE MACRO-VARIABLES FOR HETEROGENEOUS SOLIDS AT FINITE STRAIN. Proceedings of the Royal Society of London Series a-Mathematical and Physical Sciences. 1972;326:131-47. [2] Hashin Z. Large Isotropic Elastic-Deformation of Composites and Porous-Media. International Journal of Solids and Structures. 1985;21:711-20. [3] Avazmohammadi R, Castaneda PP. Tangent Second-Order Estimates for the Large-Strain, Macroscopic Response of Particle-Reinforced Elastomers. Journal of elasticity. 2013;112:139-83. [4] Guo ZY, Peng XQ, Moran B. A composites-based hyperelastic constitutive model for soft tissue with application to the human annulus fibrosus. Journal of the mechanics and physics of solids. 2006;54:1952-71. [5] Guo ZY, Caner F, Peng XQ, Moran B. On constitutive modelling of porous neo-Hookean composites. Journal of the mechanics and physics of solids. 2008;56:2338-57. [6] Guo Z, Caner FC. MECHANICAL BEHAVIOUR OF TRANSVERSELY ISOTROPIC POROUS NEO-HOOKEAN SOLIDS. International Journal of Applied Mechanics. 2010;2:11-39. [7] Segurado J, Llorca J. A numerical approximation to the elastic properties of sphere-reinforced composites. Journal of the mechanics and physics of solids. 2002;50:2107-21. [8] Guo Z, Shi X, Chen Y, Chen H, Peng X, Harrison P. Mechanical modeling of incompressible particle-reinforced neo-Hookean composites based on numerical homogenization. Mechanics of Materials. 2014;70:1-17. [9] Guo ZY, Chen Y, Peng XQ, Shi XH, Li HT, Chen YL. Shear stiffness of neo-Hookean materials with spherical voids. Composite Structures. 2016;150:21-7. [10] Guo ZY, Chen Y, Wan Q, Li HT, Shi XH, Tang S, et al. A Hyperelastic Constitutive Model for Chain-Structured Particle Reinforced Neo-Hookean Composites. Materials & Design. 2016;95:580-90. [11] Hohe J, Becker W. Effective mechanical behavior of hyperelastic honeycombs and two-dimensional model foams at finite strain. International Journal of Mechanical Sciences. 2003;45:891-913. [12] Iltchev A, Marcadon V, Kruch S, Forest S. Computational homogenisation of periodic cellular materials: Application to structural modelling. International Journal of Mechanical Sciences. 2015;93:240-55. [13] Mihai LA, Goriely A. Finite deformation effects in cellular structures with hyperelastic cell walls. International Journal of Solids and Structures. 2015;53:107-28. [14] Drugan WJ, Willis JR. A micromechanics-based nonlocal constitutive equation and estimates of representative volume element size for elastic composites. Journal of the mechanics and physics of solids. 1996;44:497-524. [15] Cricri G, Luciano R. Homogenised properties of composite materials in large deformations. Composite Structures. 2013;103:9-17. [16] Karpov EG, Stephen NG, Dorofeev DL. On static analysis of finite repetitive structures by discrete Fourier transform. International Journal of Solids and Structures. 2002;39:4291-310. 18 of 35 pages
[17] Li J, Law SS. Substructural Response Reconstruction in Wavelet Domain. Journal of Applied Mechanics-Transactions of the Asme. 2011;78:10. [18] Mencik JM. New advances in the forced response computation of periodic structures using the wave finite element (WFE) method. Computational Mechanics. 2014;54:789-801. [19] Boldrin L, Hummel S, Scarpa F, Di Maio D, Lira C, Ruzzene M, et al. Dynamic behaviour of auxetic gradient composite hexagonal honeycombs. Composite Structures. 2016;149:114-24. [20] Yu JX, Yong X, Lin W, Zhou XQ. Element-by-element model updating of large-scale structures based on component mode synthesis method. Journal of Sound and Vibration. 2016;362:72-84. [21] Ashby MF, Evans T, Fleck NA, Hutchinson JW, Wadley HNG, Gibson LJ. Metal Foams: A Design Guide: Elsevier Science, 2000. [22] Hughes TJR. The finite element method : linear static and dynamic finite element analysis. Englewood Cliffs: Prentice-Hall ; London : Prentice-Hall International, 1987. [23] ABAQUS. Analysis User's Manual, Version 6.10: SIMULIA Inc., 2010. [24] Gibson LJ, Ashby MF. Cellular Solids: Structure and Properties: Cambridge University Press, 1999.
19 of 35 pages
(a)
(b)
Fig. 1. The schematic diagrams of an RVE model (a) and the corresponding substructure element (b).
20 of 35 pages
Fig. 2. The structure of the beam consisted of 10 unit squares with circular holes, fixed at both ends and applied by a force in the middle
21 of 35 pages
(a) Top surface
(b) Bottom surface 22 of 35 pages
Fig.3. The deflections of the porous beam computed from the substructure FEM, the full FEM and the classical homogenization method.
23 of 35 pages
(a) The beam with 10 pores
(b) The beam with 20 pores
24 of 35 pages
(c) The beam with 30 pores
Fig. 4. The deflections of the porous beam computed from the substructure FEM and fitted by the classical Euler beam theory.
25 of 35 pages
Fig. 5. The structure of the porous cantilever beam, which is subjected to a uniformly distributed load on the top surface.
26 of 35 pages
(a) Top surface
(b) Bottom surface
Fig. 6. The deflections of the porous cantilever beam computed from substructure FEM, the full FEM and the homogenization method.
27 of 35 pages
Fig. 7. The structure of the honeycomb beam, fixed at both ends and applied by a force in the middle.
29 of 35 pages
(a) Top surface
(b) Bottom surface
Fig. 8. The deflections of the beam computed from the substructure FEM, the full FEM
30 of 35 pages
and the classical homogenization method.
31 of 35 pages
Fig. 9. The deflections of the honeycomb beam computed from the substructure FEM and fitted by the classical Euler beam theory.
32 of 35 pages
Fig. 10. The deflections of the bottom surface of the porous beam computed from the “full” substructure element and the simplified elements of different orders.
33 of 35 pages
Fig. 11. The deformed shapes (magnified by 5 ×109 times) of two substructure elements at the middle of the porous beam subjected to a concentrated force. The results from the “full” substructure element simulation is plotted in solid curves, while dotted curves represent the results from the third-order simplified substructure element simulation.
34 of 35 pages
Fig. 12. The deflections of the porous cantilever beam computed from the “full” substructure element and the simplified elements of different orders.
35 of 35 pages
S0263-8223(16)31940-7 http://dx.doi.org/10.1016/j.compstruct.2017.01.050 COST 8184
To appear in:
Composite Structures
Received Date: Revised Date: Accepted Date:
24 October 2016 31 December 2016 17 January 2017
Please cite this article as: Chen, Y., Dong, L., Wang, B., Chen, Y., Qiu, Z., Guo, Z., A Substructure-Based Homogenization Approach for Systems with Periodic Microstructures of Comparable Sizes, Composite Structures (2017), doi: http://dx.doi.org/10.1016/j.compstruct.2017.01.050
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
18 January 2017
Prepared for Composite Structures
A Substructure-Based Homogenization Approach for Systems with Periodic Microstructures of Comparable Sizes Yang Chen1, Leiting Dong2, Bing Wang3, Yuli Chen2, Zhiping Qiu2, Zaoyang Guo2,1* 1
Department of Engineering Mechanics, Chongqing University, Chongqing 400044, CHINA
2
School of Aeronautic Science and Engineering, Beihang University, Beijing 100191, CHINA 3
Center for composite materials, Harbin Institute of Technology, Harbin 150080, CHINA
Abstract The classical homogenization method has been widely adopted to capture the effective behaviors of heterogeneous materials. However, when the characteristic length of the microstructure of the heterogeneous material is comparable to the size of the structure, the classical homogenization method is mathematically no longer valid. In this paper, a new substructure-based homogenization approach is proposed to predict the mechanical responses of systems with periodic microstructures of comparable sizes. A substructure element is developed to reconstruct the system with periodic microstructure of comparable size. It is verified that this substructure-based
*
Corresponding author: [email protected]. Address: Institute of Solid Mechanics, Beihang University, 37
Xueyuan Road, Haidian District, Beijing 100191, CHINA.
homogenization approach can accurately predict the mechanical responses of the system. Comparing with the full finite element analysis, the computational scale is dramatically decreased. After that, a simplified substructure element is developed by using less surface nodes in the “full” substructure element. The numerical results show that, with further significantly reduced computational cost, the third-order simplified substructure element can provide a good prediction of the responses of the system with periodic microstructure of comparable size. Keywords:
substructure
approach,
periodic
microstructure,
finite
element,
homogenization
1. Introduction In recent years, the idea of homogenization has been widely adopted to estimate the effective behaviors of heterogeneous materials (e.g., composites, porous materials). One theoretical approach is to use theories of micromechanics to study the representative volume element (RVE) models and develop homogenized constitutive models to predict the effective responses of composites [1-6]. Another numerical homogenization approach is to use finite element method (FEM) to simulate the mechanical responses of the RVE models [7-10]. The RVE models are also widely used in analysis of system with periodic microstructures [11-13]. In the classical homogenization procedure, the obtained constitutive models are usually used to describe the stress-strain relation at the Gaussian integration points in finite element (FE) analysis of systems of heterogeneous materials. This implies that the classical
2 of 35 pages
homogenization approach is at the (macroscopic) material point scale. Therefore, mathematically the RVE size should be sufficiently smaller than the structure size but sufficiently larger than the characteristic length of the microstructure [14]. For example, Cricri and Luciano [15] suggested that the RVE models of cellular materials have at least 10 cells in each direction to get accurate effective parameters. Hence, the classical homogenization approaches based on RVE models cannot be adopted when the characteristic length of the microstructure of the heterogeneous material is comparable to the size of the structure (e.g., Lstruture Lmaterial < 10 ). However, due to the
novel manufacturing technologies such as 3D printing, the idea of integrated design of materials and structures leads to many structures with periodic microstructures of comparable sizes. In this paper, a substructure-based homogenization approach is developed to simulate the mechanical responses of systems with periodic microstructures of comparable sizes. In this substructure-based homogenization approach, the microstructures of the materials are considered as substructures of the structure because they are of comparable sizes. Therefore, in the substructure-based homogenization approach, the homogenization is at the structure scale rather than the material point scale. In the literature, Karpov et al. [16] adopted the substructures to analyze the static properties of finite repetitive structures, and some illustrative examples (i.e., truss bridge, clamped grid and honeycomb structures) are discussed with the substructure method. Li and Law [17] used substructures to reconstruct the system and study the response of a frame under the excitation force with wavelet domain method. Mencik [18] investigated the harmonic force response of
3 of 35 pages
one-dimensional periodic structures by analyzing the substructure with wave finite element method. Boldrin et al. [19] explored the dynamic behavior of gradient composite hexagonal honeycombs with a substructure-based component mode synthesis (CMS) method. Yu et al. [20] also adopted the CMS to study the mode shapes of the global structure with element-by-element model, which is used to update the large-scale structure. In these papers, the substructure concept is considered at the structure level, while in our approach, the substructure concept is also applied to the material level, i.e., the microstructure of the heterogeneous materials. In this paper, the substructure element is developed to describe the microstructure of the heterogeneous materials. The global stiffness matrix of the RVE model is transformed to the elemental stiffness matrix of the substructure element by eliminating the interior degrees of freedoms (DOFs). Using the substructure elements, the simulation results obtained are identical to those from full FEM simulations, but the computational cost is much lower. Comparing to classical homogenization approaches, the substructure-based approach is particularly useful when the characteristic length of the heterogeneous material is comparable to the size of the structure. Moreover, to reduce the computational scale further, the simplified substructure element is proposed and it is shown that the third-order simplified substructure element can well predict the effective behavior of the systems with periodic microstructures of comparable sizes. The structure of this paper is as follows: In Section 2, the substructure element is constructed, and the applications of the substructure-based homogenization
4 of 35 pages
approach are illustrated. Then the simplified substructure element is developed to further decrease the computational scale and it is applied to the analysis of the porous beams in Section 3. After that, some conclusion remarks are given in section 4. 2. Substructure-based homogenization approach
In this section, the substructure element is developed to describe the microstructure of the heterogeneous material and its elemental stiffness matrix is derived
using
the
concept
of
substructure.
Then
the
substructure-based
homogenization approach is applied to some examples to verify its accuracy and efficiency. 2.1. Stiffness of substructure element
Considering an RVE model of the microstructure of a heterogeneous material (an example shown in Fig. 1a), it is meshed and its mechanical response can be simulated using FEM by the following equation: Kd = P ,
(1)
where K is the system stiffness matrix of the RVE model, d denotes the vector of nodal displacements, and P represents the vector of nodal forces. The nodes are further classified as the surface nodes and the interior nodes. Then Eq. (1) can be alternatively written as K ss K is
K si d s Ps , = K ii di Pi
(2)
where subscript s and i denote the two kinds of nodes respectively. Hence, the response of the RVE can be approximated by the surface nodes only as K *ss d s = Ps* , 5 of 35 pages
(3)
where the elemental stiffness matrix of the substructure element is defined as K *ss = K ss − K si K −ii1K is ,
(4)
and the vector of equivalent surface nodal forces is computed as Ps* = Ps − K si K ii−1Pi .
(5)
Therefore, in comparison to the original RVE model, the dimension of the elemental stiffness matrix of the substructure element is reduced significantly due to the elimination of the DOFs of the interior nodes, which leads to much smaller computational scale. For instance, the RVE model of a unit square with a circular hole shown in Fig. 1(a) has 567 nodes and 987 triangular elements, which can be utilized in FEM to simulate its mechanical responses. The corresponding substructure element contains only 100 (surface) nodes (Fig. 1b), in which each edge is divided to 25 segments by 24 nodes on the edge. The number of elements also is reduced from 987 triangular elements to 1 single substructure element. It will be shown next that the substructure element can achieve exactly the same solution as the full FEM simulation of the detailed RVE model. 2.2. Applications of substructure element
In this subsection, the proposed substructure element is adopted to analyze some examples to verify the accuracy and efficiency of the substructure-based homogenization approach. 2.2.1. Beams with periodic holes
First we consider a beam, periodically stacked by 10 unit squares with circular holes, fixed at both ends and applied by a force P = 1 N m in the middle (Fig. 2).
6 of 35 pages
We assume that the thickness of the beam is much smaller than its height and length, and then plane stress is assumed in the analysis. We note that this kind of structures has been widely used in the fields of architecture and aircraft due to its light weight comparing to solid beams [21]. The material’s parameters are set as the Young’s modulus E = 210 GPa and the Poisson’s ratio v = 0.3 . The RVE model of the periodic microstructure is shown in Fig. 1. The corresponding substructure element (Fig. 1b) is developed and its elemental stiffness matrix is derived from the global stiffness matrix of the RVE model using Eq. (4). The triangular elements used in the RVE model (Fig. 1a) are normal first order triangular elements and the linear interpolation function is used for the triangular elements [22]. The beam is then modelled by 10 substructure elements and the finite element (FE) analysis is implemented using MATLAB code. The deflections of the top and bottom surfaces of the beam computed by the substructure FEM code are plotted in Fig. 3 and compared with the full FEM simulation results. The full FEM analysis uses the mesh of the RVE model for each unit. It can be found that the substructure FEM results are identical to the full FEM simulation. Here the full FEM analysis uses 5436 nodes and 9870 triangular elements, while the substructure FEM analysis has only 766 nodes and 10 elements. Therefore the substructure FEM can reduce the computational scale significantly without losing any accuracy of the results. The classical homogenization method is also attempted for this problem. To obtain the effective properties of the porous material, the RVE model in Fig. 1(a) is adopted in the commercial FEM software ABAQUS/standard 6.12 within the
7 of 35 pages
framework of linear elasticity [23]. The periodic boundary conditions, performed by “Equation” type of constraints in ABAQUS 6.12, are applied in the analysis [23]. When a uniaxial tensile load is applied to the RVE model, the effective modulus and the effective Poisson’s ratio of the porous material are computed from the FEM simulation result as Eeff E = 0.515 and ν eff ν = 0.787 respectively. Then the obtained effective properties are applied to the beam model and the computed deflections of the top and bottom surfaces are plotted in Fig. 3. In this figure, the deflections at the midpoints of both surfaces predicted by the classical homogenization method are about 24% larger than the results from the full FEM simulation as well as the substructure FEM simulation. Hence, it suggests that the classical homogenization approach is inappropriate when the characteristic length of the material microstructure is comparable to the size of the structure, while the substructure-based homogenization approach can exactly reproduce the results from full FE analysis. For the deflection of the beam fixed at both ends and applied by a concentrated force in the middle, the solution from the classical Euler beam theory is written as (assuming constant bending stiffness EI through the beam length) 1 Px3 PLx 2 − (0 ≤ x ≤ L / 2) 16 EI 12 y= 3 2 2 3 1 − Px + 3 PLx − PL x + PL EI 12 16 8 48
,
(6)
( L / 2 ≤ x ≤ L)
where L is the length of the beam, I = bh3 12 (b and h represent the thickness and height of the beam respectively) is the inertia moment of the rectangular cross section.
8 of 35 pages
Eq. (6) shows that the deflection is a third-order polynomial function of the position x (i.e., the distance from the end). The deflection of the top surface from the substructure FEM is fitted by Eq. (6) in Fig. 4(a) using the least square method and the fitted effective modulus E fit E = 0.5768 . As observed from the figure, the fitted deflection curve agrees with the one from substructure FEM simulation very well, showing an average difference of 5.57%. Here the average difference is defined as
ε ave =
1 N
2
fi − 1 , ∑ i =1 yi N
(7)
where N is the number of the data points, yi and fi are the original data and fitted data respectively. Comparing the effective modulus obtained from the classical homogenization approach Eeff = 0.515 E and the one fitted by the classical beam theory E fit = 0.5768 E , it can be easily found that the effective modulus from the classical homogenization approach Eeff is 10.71% less than E fit . We note that the classical beam theory is based on the assumption that the length of the beam is sufficiently larger than its height (i.e., L h > 10 ). To further explore the behavior of
this kind of porous beams with the classical beam theory, the beams with 20 and 30 porous units are investigated, i.e., L h =20 and
L h =30 respectively. The
deflections of the top surfaces computed from the substructure FEM are fitted by the Euler beam theory in Fig. 4 (b) and (c). The relative average differences are 3.01% and 1.84% for the beams with 20 pores and 30 pores respectively. The fitting results suggest that the deflection curves obtained from substructure FEM simulations can be well represented by the function derived from the Euler beam theory. It is found that the effective modulus E fit fitted from the two beams’ results are different from the 9 of 35 pages
one obtained from the beam with 10 holes. For the beam with 20 holes, the fitted effective modulus E (20) = 0.7709 E is 33.65% higher than E (10) fit . When the length of fit the beam increases to 30 units (i.e., with 30 holes), the fitted effective modulus (10) E (30) fit = 0.8180 E is even 41.82% higher than E fit . Both fitted effective moduli are
significantly higher than the effective modulus Eeff obtained from the classical homogenization approach, which confirms the conclusion that the classical homogenization approach is invalid when the characteristic length of the material is comparable to the size of the structure. The results also suggest that the effective modulus fitted by the Euler beam theory depends on the size of the structure. This implies that, even though the deflection of this kind of beams can be described by the same deflection function derived from the Euler beam theory, it is still difficult to use the classical beam theory to analyze this kind of beams because the effective modulus cannot be predicted correctly. To explore the application of substructure-based homogenization approach to the structure with unsymmetrical boundary, the porous beam is fixed at only one end (i.e., porous cantilever beam) and is applied to uniformly distributed load q = 1 N m 2 (Fig. 5). The deflections of the top and bottom surfaces of the beam computed by the substructure FEM code, the full FEM simulation and the classical homogenization method are plotted in Fig. 6. As can be observed from this figure, the deflections computed from the substructure FEM code are still superposed with the full FEM results, which indicates again the excellent predictions of the substructure approach. However, the classical homogenization still overestimates the deflections up to 60.88%.
10 of 35 pages
A fitted modulus E Cantilever = 0.8137 E is computed by fitting the computed deflections fit with the theoretical results of the classical Euler beam theory. The average difference is about only 3.39%, which indicates an excellent fitting. We note that the fitted modulus from the cantilever beam is 41.07% bigger than the one from the beam fixed at both ends and 58% bigger than the effective shear modulus computed from the homogenization method. This suggests that the effective modulus fitted by the Euler beam theory depends on not only the size of the structure but also the loading condition of the structure.
2.2.2. Honeycomb beam Honeycomb structures have also been widely used in many engineering fields because of their advantages such as light weight, high specific stiffness and strength [24]. To verify the accuracy and efficiency of the proposed substructure-based homogenization approach, a honeycomb beam with periodic microstructure is constructed and applied by a force P = 1 N m in the middle (Fig. 7). Both ends of the beam are fixed. In the analysis of this structure the same material parameters are assumed, i.e., E = 210 GPa and ν = 0.3 . Substructure element, developed from the RVE of honeycomb beam, is also used to reconstruct the honeycomb beam and predict the response of the beam. Fig. 8 presents the deflections of the top surface and bottom surface of the honeycomb beam, computed by the full FEM and substructure FEM. It is found that the results from the substructure-based homogenization
11 of 35 pages
approach and the full FEM simulations are identical. This suggests again that the substructure-based homogenization approach can be applied to predict the responses of honeycomb beams with periodic microstructure of comparable sizes. Similarly, the effective mechanical properties can be obtained by applying the classical homogenization method to the RVE model using ABAQUS 6.12. The computed effective properties are Eeffx E = 5.1638 × 10−3 and ν effxy ν = 3.613 . We note that only the mechanical properties along the length direction (i.e., Eeffx and
ν effxy ) are desired in our analysis according to the classical Euler beam theory, even though the RVE model of the honeycomb beam is anisotropic. The computed deflections using the classical homogenization method are also plotted in Fig. 8. It is shown that the deflections computed by the classical homogenization method are about 58% smaller than those from the full FEM simulations. Therefore we confirm that the classical homogenization approach cannot be applied to systems with microstructures of comparable sizes. The deflection curve computed from the substructure FEM simulation can be fitted by Eq. (6) from the classical Euler beam theory. The fitting result is illustrated in Fig. 9, which suggests a close fitting with an average difference of 6.72%. With the classical Euler beam theory, the fitted effective modulus is E xfit E = 2.8401× 10−3 , which is 0.45 times smaller than that computed from the classical homogenization approach. As pointed out before, the main difficulty to use Euler beam theory here is to predict the effective modulus of the heterogeneous material when the characteristic length of its microstructure is comparable to the size of the structure.
12 of 35 pages
3. Simplified substructure element
By eliminating the DOFs of the interior nodes, the substructure element approach can reduce the computational scale significantly. However, because the surface nodes of one substructure element are connected via interior nodes, the elemental stiffness matrix of a substructure element in Eq. (4) is usually a full matrix. For large scale system (with a large number of substructure elements), the system matrix will still be a sparse matrix even though every elemental stiffness matrix is full. Nevertheless, further reducing the number of DOFs of each substructure element will improve the computational efficiency significantly. This section aims to develop a simplified substructure element to further decrease the computational scale, but keeping the prediction accurate enough. In order to reduce the computational scale, the number of nodes on the surface (or boundary) of the substructure element is reduced. To implement that, some nodes are chosen at the surface of the “full” substructure element for the simplified substructure element. This can usually be done by keeping some of the surface nodes in the “full” substructure element and eliminating others. We may also choose a set of new locations on the surface as the surface nodes of the simplified substructure element as long as they can be used to interpolate the displacement field of the surface effectively. The surface nodal displacements of the “full” substructure element are interpolated from those of the simplified substructure element. To obtain an effective interpolation and a well-posed elemental stiffness matrix of the simplified substructure element, the selected surface nodes of the simplified substructure element should usually be evenly distributed on the surface
13 of 35 pages
area. The elemental stiffness matrix of the simplified substructure element can then be derived from that of the full substructure element as follows. First denote the nodal displacement vectors of the simplified and “full” substructure elements as u and u′ respectively. The vector u′ is interpolated by u as u′ = Nu ,
(8)
where N is the interpolation function matrix. The elastic strain energy of the substructure element can be represented as 1 1 W = uT Ku = u′T K ′u′ , 2 2
(9)
where K and K ' denote the elemental stiffness matrices of the simplified and the “full” substructure elements. Based on (8) and (9), the elemental stiffness matrix of the simplified substructure element is obtained as K = NT K′N .
(10)
Similarly, we can derive the equivalent nodal force vector P of the simplified substructure element from the nodal force vector P ′ of the “full” substructure element using the following principle of virtual work PT du = P′T du′ .
(11)
Here du and du′ are corresponding virtual nodal displacement vectors. The equivalent nodal force vector P of the simplified substructure element is then obtained as P = N T P′ .
(12)
After the elemental stiffness matrix of the simplified substructure element is created, the analysis can be carried out with a highly efficient method. Let’s consider
14 of 35 pages
the example of the beam with 10 pores again (Fig. 2), in accordance to (10), the “full” substructure element (originally with 100 surface nodes) is simplified to first-order, second-order and third-order elements, which have 4, 8, and 12 surface nodes respectively. Similarly, FE analysis implemented by MATLAB codes is carried out to study the behavior of the porous beam using the simplified substructure elements. Fig. 10 demonstrates the deflections of the bottom surface of the beam computed using the 3 types of the simplified substructure elements as well as the “full” substructure element. It can be found that, comparing with the result from the analysis using the “full” substructure element, the deflections predicted by the first-order and second-order simplified substructure elements underestimate significantly the beam’s deflection (the average differences are 20.51% and 17.02% respectively). The prediction by the third-order simplified substructure element, however, reveals a good approximation of the response of the beam, which only underestimates the deflection for only about 4.00%. The detailed deformed shapes of the elements are also explored. In Fig. 11, the shapes of two deformed third-order simplified substructure elements at the middle of the porous beam (the concentrated force is applied at the top joint node) are plotted against those of two “full” substructure elements (the deformation is magnified by 5 × 109 times). It is clear that the deformations computed from the two types of elements are very close, which implies the third-order simplified element can well capture the deformation of the porous beam. The first-order, second-order and third-order simplified substructure element are also employed to study the example of the porous cantilever beam, and the defections
15 of 35 pages
of bottom surface computed from the 3 types of elements as well as the “full” substructure element are plotted in Fig. 12. The deflections predicted by the first-order and second-order simplified substructure elements also underestimate the beam’s deflection, showing the average differences of 16.88% and 6.88% respectively. Again the prediction by the third-order simplified substructure element demonstrates a very good approximation of the response of the beam, which only underestimates the deflection for about 2.22%. Comparing to the “full” substructure element with 100 surface nodes, the third-order simplified substructure element only has 12 nodes, which implies nearly 90% reduction of the number of DOFs. The accuracy of the prediction by the third-order simplified substructure element is acceptable with considerably reduced computational cost. Therefore, it is shown that the simplified substructure element, which uses much less surface nodes than the “full” substructure element, has the ability to accurately approximate the responses of structures with periodic microstructures of comparable sizes using significantly reduced computational resources. 4. Conclusion Remarks
Although the classical homogenization approach has been widely used to predict the effective behaviors of heterogeneous materials in the literature, it is demonstrated in this paper that, when the characteristic length of the material is comparable to the size of the structure, it cannot correctly capture the mechanical responses of the structure. This is reasonable because the classical homogenization method is mathematically
16 of 35 pages
invalid in this situation. To predict the mechanical behavior of systems with periodic microstructures of comparable sizes, a new substructure-based homogenization approach is proposed in this paper. The substructure element is developed based on the RVE model representing the microstructure of the material. The elemental stiffness matrix of the substructure element is derived from the global stiffness matrix of the RVE model by eliminating the interior DOFs. Several numerical simulation results show that this substructure-based homogenization approach can accurately predict the mechanical responses of the system. Comparing to the full finite element analysis, the computational scale of the substructure-based approach is dramatically decreased. The classical beam theory is used to fit the deflections of the beams and it is found that the fitted effective modulus depends on the size of the structure (that is, it is not a material property). It shows the difficulty of applying the classical beam theory directly to the structure with microstructure of comparable size. After that, a simplified substructure element is developed by using less surface nodes in the “full” substructure element. The numerical results show that, with further significantly reduced computational cost, the third-order simplified substructure element can provide a good prediction of the responses of the system with periodic microstructure of comparable size. Acknowledgement
The financial support from NSFC (No. 11272362, 11572020, 11332013, CHINA), Chongqing Science and Technology Commission (cstc2013jcyjjq50003), and Beihang University startup grant is greatly appreciated. References
17 of 35 pages
[1] Hill R. CONSTITUTIVE MACRO-VARIABLES FOR HETEROGENEOUS SOLIDS AT FINITE STRAIN. Proceedings of the Royal Society of London Series a-Mathematical and Physical Sciences. 1972;326:131-47. [2] Hashin Z. Large Isotropic Elastic-Deformation of Composites and Porous-Media. International Journal of Solids and Structures. 1985;21:711-20. [3] Avazmohammadi R, Castaneda PP. Tangent Second-Order Estimates for the Large-Strain, Macroscopic Response of Particle-Reinforced Elastomers. Journal of elasticity. 2013;112:139-83. [4] Guo ZY, Peng XQ, Moran B. A composites-based hyperelastic constitutive model for soft tissue with application to the human annulus fibrosus. Journal of the mechanics and physics of solids. 2006;54:1952-71. [5] Guo ZY, Caner F, Peng XQ, Moran B. On constitutive modelling of porous neo-Hookean composites. Journal of the mechanics and physics of solids. 2008;56:2338-57. [6] Guo Z, Caner FC. MECHANICAL BEHAVIOUR OF TRANSVERSELY ISOTROPIC POROUS NEO-HOOKEAN SOLIDS. International Journal of Applied Mechanics. 2010;2:11-39. [7] Segurado J, Llorca J. A numerical approximation to the elastic properties of sphere-reinforced composites. Journal of the mechanics and physics of solids. 2002;50:2107-21. [8] Guo Z, Shi X, Chen Y, Chen H, Peng X, Harrison P. Mechanical modeling of incompressible particle-reinforced neo-Hookean composites based on numerical homogenization. Mechanics of Materials. 2014;70:1-17. [9] Guo ZY, Chen Y, Peng XQ, Shi XH, Li HT, Chen YL. Shear stiffness of neo-Hookean materials with spherical voids. Composite Structures. 2016;150:21-7. [10] Guo ZY, Chen Y, Wan Q, Li HT, Shi XH, Tang S, et al. A Hyperelastic Constitutive Model for Chain-Structured Particle Reinforced Neo-Hookean Composites. Materials & Design. 2016;95:580-90. [11] Hohe J, Becker W. Effective mechanical behavior of hyperelastic honeycombs and two-dimensional model foams at finite strain. International Journal of Mechanical Sciences. 2003;45:891-913. [12] Iltchev A, Marcadon V, Kruch S, Forest S. Computational homogenisation of periodic cellular materials: Application to structural modelling. International Journal of Mechanical Sciences. 2015;93:240-55. [13] Mihai LA, Goriely A. Finite deformation effects in cellular structures with hyperelastic cell walls. International Journal of Solids and Structures. 2015;53:107-28. [14] Drugan WJ, Willis JR. A micromechanics-based nonlocal constitutive equation and estimates of representative volume element size for elastic composites. Journal of the mechanics and physics of solids. 1996;44:497-524. [15] Cricri G, Luciano R. Homogenised properties of composite materials in large deformations. Composite Structures. 2013;103:9-17. [16] Karpov EG, Stephen NG, Dorofeev DL. On static analysis of finite repetitive structures by discrete Fourier transform. International Journal of Solids and Structures. 2002;39:4291-310. 18 of 35 pages
[17] Li J, Law SS. Substructural Response Reconstruction in Wavelet Domain. Journal of Applied Mechanics-Transactions of the Asme. 2011;78:10. [18] Mencik JM. New advances in the forced response computation of periodic structures using the wave finite element (WFE) method. Computational Mechanics. 2014;54:789-801. [19] Boldrin L, Hummel S, Scarpa F, Di Maio D, Lira C, Ruzzene M, et al. Dynamic behaviour of auxetic gradient composite hexagonal honeycombs. Composite Structures. 2016;149:114-24. [20] Yu JX, Yong X, Lin W, Zhou XQ. Element-by-element model updating of large-scale structures based on component mode synthesis method. Journal of Sound and Vibration. 2016;362:72-84. [21] Ashby MF, Evans T, Fleck NA, Hutchinson JW, Wadley HNG, Gibson LJ. Metal Foams: A Design Guide: Elsevier Science, 2000. [22] Hughes TJR. The finite element method : linear static and dynamic finite element analysis. Englewood Cliffs: Prentice-Hall ; London : Prentice-Hall International, 1987. [23] ABAQUS. Analysis User's Manual, Version 6.10: SIMULIA Inc., 2010. [24] Gibson LJ, Ashby MF. Cellular Solids: Structure and Properties: Cambridge University Press, 1999.
19 of 35 pages
(a)
(b)
Fig. 1. The schematic diagrams of an RVE model (a) and the corresponding substructure element (b).
20 of 35 pages
Fig. 2. The structure of the beam consisted of 10 unit squares with circular holes, fixed at both ends and applied by a force in the middle
21 of 35 pages
(a) Top surface
(b) Bottom surface 22 of 35 pages
Fig.3. The deflections of the porous beam computed from the substructure FEM, the full FEM and the classical homogenization method.
23 of 35 pages
(a) The beam with 10 pores
(b) The beam with 20 pores
24 of 35 pages
(c) The beam with 30 pores
Fig. 4. The deflections of the porous beam computed from the substructure FEM and fitted by the classical Euler beam theory.
25 of 35 pages
Fig. 5. The structure of the porous cantilever beam, which is subjected to a uniformly distributed load on the top surface.
26 of 35 pages
(a) Top surface
(b) Bottom surface
Fig. 6. The deflections of the porous cantilever beam computed from substructure FEM, the full FEM and the homogenization method.
27 of 35 pages
Fig. 7. The structure of the honeycomb beam, fixed at both ends and applied by a force in the middle.
29 of 35 pages
(a) Top surface
(b) Bottom surface
Fig. 8. The deflections of the beam computed from the substructure FEM, the full FEM
30 of 35 pages
and the classical homogenization method.
31 of 35 pages
Fig. 9. The deflections of the honeycomb beam computed from the substructure FEM and fitted by the classical Euler beam theory.
32 of 35 pages
Fig. 10. The deflections of the bottom surface of the porous beam computed from the “full” substructure element and the simplified elements of different orders.
33 of 35 pages
Fig. 11. The deformed shapes (magnified by 5 ×109 times) of two substructure elements at the middle of the porous beam subjected to a concentrated force. The results from the “full” substructure element simulation is plotted in solid curves, while dotted curves represent the results from the third-order simplified substructure element simulation.
34 of 35 pages
Fig. 12. The deflections of the porous cantilever beam computed from the “full” substructure element and the simplified elements of different orders.
35 of 35 pages