Acta Polytechnica CTU Proceedings https://doi.org/10.14311/APP.2022.34.0127 Acta Polytechnica CTU Proceedings 34:127–132, 2022 © 2022 The Author(s). Licensed under a CC-BY 4.0 licence Published by the Czech Technical University in Prague INFLUENCE OF MICROSTRUCTURE ON RATE-DEPENDENT RESPONSE OF UNIDIRECTIONAL FIBROUS COMPOSITES Soňa Valentová, Michal Šejnoha∗, Jan Vorel, Zdeněk Prošek Czech Technical University in Prague, Faculty of Civil Engineering, Thákurova 7, 166 29 Prague 6, Czech Republic ∗ corresponding author: sejnom@fsv.cvut.cz Abstract. This paper outlines prediction of the macroscopic response of unidirectional fibrous composites made either from basalt or carbon fibers impregnated by a polymeric epoxy matrix. The viscoelastic response of the matrix phase was represented by the Maxwell chain model. A series of creep tests performed at several stress levels served as a stepping stone for the model calibration. The macroscopic behavior of both composites was first examined via computational homogenization. Attention was accorded to computational cells with variable size extracted from large representative images. We observed that selecting the computational model as a sufficiently large test window should be approached with caution. Because initial designs often need a large parametric study to test various material and geometrical patterns, this study was then supported by computationally much more effective Mori-Tanaka averaging scheme, clearly showing its potential even if loading the composite beyond its elastic limit. Keywords: Viscoelasticity, fibrous composites, Mori-Tanaka method, periodic unit cell, generalized Leonov model, Maxwell chain model, creep. 1. Introduction Fibrous composites with a polymeric matrix are used in many engineering applications. This type of ma- terial is often favored owing to, e.g., high strength, stiffness, thermal or chemical resistance, while the weight and volume are relatively low. Presently, also an ecological point of view plays a significant role and both basalt and carbon fibers, used in the examined composites, appear convenient. Thus the resulting system may well substitute more traditional materials, such as steel or concrete, in many industrial fields. On the other hand, prediction of the macroscopic behavior of such a material, where the matrix is vis- coelastic, calls for combining a comprehensive, rela- tively time demanding, laboratory program and com- putational analysis to calibrate the selected nonlinear viscoelastic model. Having the calibrated model, the determination of macroscopic response then builds upon the applica- tion of either computational homogenization [1, 2] or application of one of the macromechanical model [3, 4] such as the Mori-Tanaka method [5]. We begin with the first-order computational ho- mogenization, but unlike in [6, 7] we abandon the approach based on the statistically equivalent peri- odic unit cell (SEPUC) and formulate the representa- tive volume element (RVE) as a certain test window taken from randomly from a large image of a real mi- crostructure [8]. Despite random distribution of fibers, for unidirectional fibrous composite a periodic unit cell (PUC) with the periodic hexagonal arrangement (PHA), shown in the Fig. 1, is often assumed suffi- cient. The comparative results shows some drawbacks of both models when loading the present composites beyond elasticity. (a) (b) Figure 1. a) PUC with PHA, b) example of mesh. As already pointed out the finite element element simulations are generally computationally expensive. It is therefore useful to have more efficient, yet reliable, method at hand especially if the averages of local fields are sufficient. In this regard we attempt to reproduce the finite element simulations with the modified Mori- Tanaka method, here presented in the framework of Dvorak’s transformation field analysis [3, 9]. 2. Theoretical background 2.1. Generalized Leonov model As the epoxy matrix influences the macroscopic re- sponse time and rate-dependently, a reliable model describing its nonlinear viscoelastic behavior is indis- pensable. Herein the generalized Leonov model is used. A brief description is provided next while the details can be found in [4, 10, 11]. The volumetric response is assumed elastic σm = Kεv, (1) 127 https://doi.org/10.14311/APP.2022.34.0127 https://creativecommons.org/licenses/by/4.0/ https://www.cvut.cz/en S. Valentová, M. Šejnoha, J. Vorel, Z. Prošek Acta Polytechnica CTU Proceedings whereas the deviatoric response is described by the generalized Maxwell chain model as dsij dt = M∑ µ=1 2Gµ(deij dt − dep,µij dt ), (2) sij = M∑ µ=1 sµij . (3) The Eyring flow model written as dep dt = 1 2A sinh(τ/τ0) , (4) represents the plastic shear rate of deformation and determines the creep strain. The model parameters A and τ0 are obtained from laboratory experiments, see [10–12]. The description of the required creep tests to cali- brate the Maxwell model for the adopted epoxy ma- trix is available in [13]. Because of space limitation we present in Table 1 only the resulting parameters (Jµ, τmu) of individual chains necessary for the con- struction of creep compliance function of the matrix phase. The Laplace transform is used next to get the necessary stiffnesses Gµ entering Eq. (2). 2.2. Mori-Tanaka method The classical micromechanical Mori-Tanaka method belongs to frequently used two-point averaging schemes. Unlike the finite element approach exploit- ing periodic unit cells, this method brings a significant computational simplification, as just the knowledge of the shape and orientation of reinforcements, and the volume fraction and material properties of both phases, fibers and matrix, are needed. We start the description from the local stress and strain increments provided by ∆σf = Lf∆εf , ∆σm = L̂m(∆εm −∆µm), (5) where L̂m represents the matrix stiffness dependent on the viscoelastic modulus. According to Dvorak’s transformation field analysis, see [9], the local strain increments are written as ∆εf = Âf∆E + D̂fm∆µm, (6) ∆εm = Âm∆E + D̂mm∆µm. (7) The strain localization factors Âr (r = f,m stand for the fiber and matrix phase, respectively), attain the form Âm = [ cmI + cf T̂f ]−1 , Âf = T̂f Âm. (8) With reference to [5], the transformation influence matrices D̂rm are given by D̂rm = ( I− Âr )( L̂m − Lf )−1 L̂m. (9) It was pointed out in [14] that if the nonlinear matrix influences the macroscopic behavior, the Mori-Tanaka method gives a noticeably stiffer response in compar- ison to numerical predictions provided by the finite element method. Thus in our recent paper, see [13], we introduced a certain modification to the original Mori-Tanaka formulation which showed potential in arriving at results being reasonably close to those provided by FEM. 2.3. Homogenization using FEM The first-order homogenization assumes periodic fields. The periodic hexagonal array, seen in the Fig. 1, can be selected as a suitable candidate. In the present study, however, other RVEs, generally non-periodic, will also be examined while enforcing the periodicity on the fluctuation part of the displacement field u∗(x). The local displacement and strain increments can be then split into homogeneous and fluctuation parts as ∆u(x) = ∆E · x + ∆u∗(x), (10) ∆ε(x) = ∆E + ∆ε∗(x). (11) If the macroscopic strain E is prescribed, with arrive at a particular form of the Hill lemma〈 δεT∆σ 〉 = 0, (12) where 〈·〉 stands for the volume averaging. The solu- tion of the resulting system of algebraic equations [4] is searched in terms of fluctuations u∗. The local stress increment can be expressed as ∆σ(x) = L(x)(∆ε(x)−∆µ(x)). (13) Their volume averages then allow us to construct the desired macroscopic stress-strain diagrams. Further details can be found for example in [14]. 3. Results Two types of the unidirectional fibrous composites re- inforced by basalt and carbon fibers were investigated. For both systems, the same epoxy resin was used. Figure 2 presents binary images of large cross-sections taken from corresponding yarns clearly suggesting a random nature of fiber distributions. The fiber and matrix volume fractions were estimated via standard image analysis. The corresponding values are listed in Table 2 together with the material properties of fibers taken from [15]. The two images in Fig. 2 also served to construct the computational models of variable complexity while at- tempting to match the corresponding volume fractions as close as possible. Twenty such models, henceforth denoted as unit cells (UC), were generated by taking a random cut from the original images allowing us to cover a sufficiently large range of cross-sectional areas. The unit cells are numbered in the ascending order from the lowest to the highest area. Two particular examples identifying the cells with the lowest and 128 vol. 34/2022 Viscoelastic response of UD fibrous composites µ τµ [s] Jµ [MPa−1] θµ [MPa·s] Eµ [MPa] 1 0.001 2.606512×10−4 9.927397×10−3 2.787166×101 2 0.01 1.905071×10−6 9.966502×10−2 1.278184×101 3 0.1 8.808431×10−7 9.815126×10−1 7.056602×101 4 1 4.934025×10−6 9.543319×10+0 1.711529×102 5 10 1.276165×10−5 9.344254×10+1 2.334448×101 6 100 1.969419×10−5 9.580883×10+2 1.418353×102 7 1000 1.290521×10−5 8.275395×10+3 5.659977×102 8 10000 6.291266×10−5 9.647045×10+4 1.586346×102 9 100000 7.887707×10−6 2.005373×10+5 1.944645×103 10 1000000 1.577867×10−3 4.168654×10+5 5.096147×102 Table 1. Parameters of Maxwell chain model. EA ET GA GT νA cf [GPa] [GPa] [GPa] [GPa] [-] [-] Carbon fibers 294 13 12 5 0.24 0.568 Basalt fibers 69.68 64.82 28.10 26.14 0.4 0.563 Table 2. Material properties of carbon and basalt fibers and its volume fraction in yarn. (a) (b) Figure 2. Cross-sections of examined fibrous com- posites. a) Carbon fiber composite, b) basalt fiber composite. (a) (b) Figure 3. a) Basalt cell with the lowest volume, b) Basalt cell with the highest volume. (a) (b) Figure 4. a) Carbon cell with the lowest volume, b) Carbon cell with the highest volume. highest area plotted for each system in Figs. 3 and 4 for illustration. Both materials were subjected to transverse shear loading in the strain control regime and the macro- scopic response was recorded. At the initial stage the composite was subjected to the prescribed shear strain rate of 0.001 s−1 until the maximum value of 0.1 was reached. Then the shear strain was held constant for the next 100 s to monitor relaxation. The resulting macroscopic response is plotted in Figs. 5 - 9. Figures 5 and 6 show the shear stress- strain diagram and the time variation of the shear stress, respectively, for the basalt composite (Bi UC). Similar graphs are displayed for the carbon composite (Ci UC) in Figs. 7 and 8. The bold lines correspond to the cells with the maximum and minimum area (VMIN, VMAX), stiffness (SMIN, SMAX), and volume fraction (ciMIN, ciMAX), respectively. 129 S. Valentová, M. Šejnoha, J. Vorel, Z. Prošek Acta Polytechnica CTU Proceedings The results corresponding to the PHA model and those derived from the Mori-Tanaka method are also provided for the sake of comparison. It should be pointed out that the MT results correspond to the standard formulation with no particular adjustments as presented in [13] to match the results of the PHA model. Figure 5. Shear stress - strain diagram for basalt composites. Figure 6. Time evolution of shear stress for basalt composites. Figure 9 then compares the overall behavior of the two systems. The graphs corresponding to the basalt and carbon UC represent averages obtained from all unit cells. For both systems we notice a relatively wide spread of results with no particular trend from the UC area point of view. The expected convergence with increas- ing size of UC as shown, e.g., in [8] for the nonlinear response of random masonry walls has not been ob- served for the viscoelastic response of the examined material systems. But it is interesting to point out that essentially all curves fall within the limits set by SMIN and SMAX cells. The response of the PHA model appears considerably more compliant when compared Figure 7. Shear stress - strain diagram for carbon composites. Figure 8. Time evolution of shear stress for carbon composites. Figure 9. Time evolution of shear stress for carbon and basalt composites. to all unit cells. Surprisingly, this is also true even for the original format of the MT method. This may suggest the need for a rigorously defined statistically equivalent periodic unit cell as proposed, e.g., in [6]. 130 vol. 34/2022 Viscoelastic response of UD fibrous composites But this goes beyond the present scope. It is also clear that at some applications the PHA model rec- ommended for random composites with a sufficiently high volume fraction of fibers may not always provide reliable results. (a) (b) Figure 10. a) Basalt cell with minimal stiffness, b) basalt cell with maximal stiffness. (a) (b) Figure 11. Basalt composite with minimal stiffness, time = 50 s: a) Shear stress [MPa], b) creep strain. (a) (b) Figure 12. Basalt composite with minimal stiffness, time = 100 s: a) Shear stress [MPa], b) creep strain. (a) (b) Figure 13. Basalt composite with maximal stiffness, time = 50 s: a) Shear stress [MPa], b) creep strain. The evolution of stresses and creep strains is further demonstrated on surface plots in Figs. 11 - 14 and (a) (b) Figure 14. Basalt composite with maximal stiffness, time = 100 s: a) Shear stress [MPa], b) creep strain. (a) (b) Figure 15. Carbon cell: a) minimal stiffness, b) maximal stiffness. (a) (b) Figure 16. Carbon composite with minimal stiffness, time = 50 s: a) Shear stress [MPa], b) creep strain. (a) (b) Figure 17. Carbon composite with minimal stiffness, time = 100 s: a) Shear stress [MPa], b) creep strain. Figs. 16 - 19 for the basalt and carbon fibers, respec- tively. Given the results in Figs. 5 - 8 we limited our attention to cells with the lowest and highest stiff- ness illustrated in Fig. 10 for the basalt fiber based composite and in Fig. 15 for the carbon fiber based system. 131 S. Valentová, M. Šejnoha, J. Vorel, Z. Prošek Acta Polytechnica CTU Proceedings (a) (b) Figure 18. Carbon composite with maximal stiffness, time = 50 s: a) Shear stress [MPa], b) creep strain. (a) (b) Figure 19. Carbon composite with maximal stiffness, time = 100 s: a) Shear stress [MPa], b) creep strain. We see that for the SMIN cell the creep strain is more localized when compared to the SMAX cell where the corresponding distribution is more diffuse. This par- tially explains deviations of the macroscopic response pertinent to individual unit cells. 4. Conclusion Viscoelastic response of the epoxy matrix based basalt and carbon fiber composites loaded in in-plane shear was examined. The option of generating the com- putational models applicable with the first-order ho- mogenization as test windows, taken randomly from large binary images, in place of statistically equivalent periodic unit cell was examined. This approach, how- ever, seems inappropriate given the large scatter of the predicted stress-strain curves. The present results also suggested that application of the PHA model [1], generally accepted for random fibrous composites with high volume fraction of fibers, should be approched with caution. The same holds for the Mori-Tanaka method when applied to complex nonlienar systems. Acknowledgements The financial support provided by the Czech Technical University in Prague within SGS project with the applica- tion registered under the No. SGS21/037/OHK1/1T/11 and by the GAČR grant No. 21-28525S is gratefully ac- knowledged. References [1] J. L. Teplý, G. J. Dvorak. Bound on overall instantaneous properties of elastic-plastic composites. Journal of the Mechanics and Physics of Solids 36(1):29–58, 1988. [2] J. Fish, Q. Yu, K. Shek. Computational damage mechanics for composite materials based on mathematical homogenization. International Journal for Numerical Methods in Engineering 45(11):1657–1679, 1999. [3] G. Dvorak. Micromechanics of composite materials. Springer Dordrecht Heidelberg New York London, 2013. [4] M. Šejnoha, J. Zeman. Micromechanics in Practice. WIT Press, Southampton, Boston, 2013. [5] Y. Benveniste. A new approach to the application of Mori-Tanaka theory in composite materials. Mechanics of Materials 6:147–157, 1987. [6] J. Zeman, M. Šejnoha. Numerical evaluation of effective properties of graphite fiber tow impregnated by polymer matrix. Journal of the Mechanics and Physics of Solids 49(1):69–90, 2001. [7] J. Zeman, M. Šejnoha. From random microstructures to representative volume elements. Modelling and Simulation in Materials Science and Engineering 15(4):S325–S335, 2007. [8] S. Tiberti, G. Milani. 3d homogenized limit analysis of non-periodic multi-leaf masonry walls. Computers and Structures 234:106253, 2020. https://doi.org/10.1016/j.compstruc.2020.106253. [9] G. J. Dvorak, Y. Benveniste. On transformation strains and uniform fields in multiphase elastic media. Proceedings of the Royal Society of London Series A - Mathematical, Physical and Engineering Sciences 437(1900):291–310, 1992. [10] T. A. Tervoort. Constitutive modeling of polymer glasses: Finite, nonlinear visocelastic behaviour of polycarbonate. Ph.D. thesis, Eindhoven University of Technology, Eindhoven, 1996. [11] R. Valenta, M. Šejnoha, J. Zeman. Macroscopic constitutive law for mastic asphalt mixtures from multiscale modeling. International Journal for Multiscale Computational Engineering 8(1):131–149, 2010. [12] R. Valenta, M. Šejnoha. Hierarchical modeling of mastic asphalt in layered road structures based on the Mori-Tanaka method. Acta Polytechnica 52(6):48–58, 2012. [13] S. Valentová, M. Šejnoha, J. Vorel, et al. Application of the mori-tanak method to describe the rate- dependent behavior of unidirectional fibrous composites. Acta Polytechnica CTU Proceedings 30:114–120, 2021. [14] S. Valentová, M. Šejnoha, J. Vorel. Comparing Mori-Tanaka method and first-order homogenization scheme in the viscoelastic modeling of unidirectional fibrous composites. Acta Polytechnica CTU Proceedings 26:133–138, 2020. https://doi.org/10.14311/APP.2020.26.0133. [15] J. Vorel, E. Grippon, M. Šejnoha. Effective thermoelastic properties of polysiloxane matrix based plain weave textile composites. International Journal for Multiscale Computational Engineering 13(3):181–200, 2015. 132 https://doi.org/10.1016/j.compstruc.2020.106253 https://doi.org/10.14311/APP.2020.26.0133 Acta Polytechnica CTU Proceedings 34:127–132, 2022 1 Introduction 2 Theoretical background 2.1 Generalized Leonov model 2.2 Mori-Tanaka method 2.3 Homogenization using FEM 3 Results 4 Conclusion Acknowledgements References