PROBABILISTIC FINITE ELEMENT SIMULATION OF SHEET MOLDING COMPOUND COMPOSITES WITH MULTIMODAL WEIBULL DISTRIBUTION By Sakib Iqbal A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of Mechanical Engineering—Doctor of Philosophy 2023 ABSTRACT Sheet Molding Compound (SMC) is a type of ready to mold composites material. The most com- mon SMC consists of glass fiber bundles about one inch long distributed randomly in a B-stage polyester resin. SMC possesses good mechanical properties and manufacturing flexibility in form- ing complex shaped parts and is relatively low cost, making it one of the attractive choices to re- place metallic parts in automotive industry. Nevertheless, SMC composites have not been utilized in critical automobile components owing to the lack of a satisfactory predictive model, especially for crashworthiness simulations. The main challenges in analysis of SMC structures are the large scatter in mechanical properties and the large difference in strengths under different stress distri- bution or loading conditions. For example, SMC demonstrates 1.5-2 times higher strength under 3-point (3-pt) bending in comparison to uniaxial tension strength. This phenomenon is known as the size effect on strength and can be explained by Weibull’s statistical strength theory. For materi- als with large size effect such as SMC, simulations carried out with the mean mechanical properties (i.e., tension, compression, and shear data) would result in a significant underprediction of flexural responses of the structure. To improve the predictions, the statistical distributions of the mechan- ical properties need to be considered and the size effect should be examined. Although statistical analysis has long been considered in composite designs, probabilistic finite element (PFE) analysis based on statistical strength models has also been employed to consider uncertainties and design re- liability at every scale in composites, little work has been done to examine the size effect of strength in FE simulations. This work aims to incorporate the size effect in probabilistic simulations of SMC composite structures. First, we extended the unimodal Weibull strength model into multimodal one by com- bining the tensile and flexural Weibull strength models. This approach was examined with a glass fiber SMC composite. A randomization algorithm was developed to incorporate the strength dis- tribution in PFE models. The strength distribution model was discretized into a limited number of segments and the values of the average strength for each segment and their probabilities were determined. The strength values were then randomly assigned across the integration points in the PFE model according to their probabilities. This approach successfully reproduced the tensile and flexural responses with the mean peak load, post peak behavior, and energy absorption similar to experimental results within ten iterations. Next, in addition to the tensile strength, the statistical distributions of the elastic modulus and compressive strength were also considered. The tensile strength and compressive strength were modeled by bimodal Weibull strength distributions corre- sponding to the uniaxial and 3-pt bending experiments. To determine the mixture weight fraction of the bimodal models and some difficult to measure parameters in the damage mechanics based composite material model, model optimization was explored using two techniques: (1) Artificial Neural Network (ANN)-based machine learning (ML) and (2) Random Search. It was observed that although computationally inexpensive, ANN-ML was rather complicated for a general-purpose re- gression. On the other hand, RS is easy to implement. Its high computational cost is acceptable as the optimization has to be done only once for any specific material model. The PFE models optimized with RS were examined with four verification cases including tensile, compression, 3-pt and 4-pt bending, and three validation cases including open hole tension, disk bending with a fixed boundary and with a simply supported boundary conditions. The PFE predictions agreed well with the experimental results across these load cases. Copyright by SAKIB IQBAL 2023 This dissertation is dedicated to my parents, family, and friends, who have been my pillars of strength and unwavering support throughout my academic journey. v ACKNOWLEDGEMENTS I would like to express my deepest gratitude to my thesis supervisor, Dr. Xinran Xiao, for her invaluable guidance, encouragement, and unwavering support throughout my research journey. Her extensive knowledge and expertise in the field have been instrumental in shaping my understanding and developing my research skills. I am also grateful to Dr. Yang Guo, Dr. Weiyi Lu, and Siva Nadimpalli for their insightful comments, suggestions, and feedback that have significantly improved the quality of my work. I am deeply indebted to Dr. Shutian Yan, Royal Ihuaenyi for their constant support and inspira- tion. I am thankful for the effort made by Beckett Clifford and Dr. Erik Stitt in CVRC lab to train, maintain, and repair the machines. I would also like to thank the faculty and staff of Mechanical Engineering, who have provided me with access to resources, facilities, and opportunities that have enriched my academic experience. Finally, I would like to express my heartfelt appreciation to my family for their unconditional love, encouragement, and support throughout my academic journey. Their unwavering faith in me has been a constant source of strength and inspiration. This work is partially supported by General Motors LLC. I would also like to thank ME De- partment, and Graduate School for the funding support. vi TABLE OF CONTENTS LIST OF ABBREVIATIONS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix LIST OF SYMBOLS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xi CHAPTER 1 INTRODUCTION . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.1 Introduction to Sheet Molding Compound (SMC) . . . . . . . . . . . . . . . . 1 1.2 Motivations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.3 Objectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.4 Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 CHAPTER 2 LITERATURE REVIEW . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.1 SMC Manufacturing and Its Morphology . . . . . . . . . . . . . . . . . . . . 5 2.2 Micromechanics and Multiscale Models of SMC . . . . . . . . . . . . . . . . 8 2.3 Statistical Strength Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.4 Weibull Theory of Strength . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.5 Multimodal Weibull Distribution . . . . . . . . . . . . . . . . . . . . . . . . . 17 2.6 Machine Learning and Statistical Analysis . . . . . . . . . . . . . . . . . . . . 18 2.7 Random Search Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 2.8 Simulations with Statistical Consideration . . . . . . . . . . . . . . . . . . . . 23 CHAPTER 3 MECHANICAL CHARACTERIZATION OF SMC . . . . . . . . . . . 25 3.1 Uniaxial Tension Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 3.2 Uniaxial Compression Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 3.3 3-Pt Bending Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 3.4 4-Pt Bending Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 3.5 In-Plane Shear Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38 3.6 Cyclic Shear Test . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 3.7 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 CHAPTER 4 PFE MODEL WITH RANDOM TENSILE STRENGTH . . . . . . . . . 47 4.1 Mechanical Properties of SMC . . . . . . . . . . . . . . . . . . . . . . . . . . 47 4.2 Statistical Model of SMC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 4.3 Considering Randomness in PFE Models . . . . . . . . . . . . . . . . . . . . 53 4.4 Results and Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 4.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 CHAPTER 5 OPTIMIZATION OF PFE MODEL WITH MACHINE LEARNING . . . . . . . . . . . . . . . . . . . . . . . . . . 70 5.1 Mechanical Properties of SMC . . . . . . . . . . . . . . . . . . . . . . . . . . 70 5.2 Statistical Model & PFE Simulation . . . . . . . . . . . . . . . . . . . . . . . 70 5.3 Optimization of Model Parameter Through Machine Learning . . . . . . . . . 79 5.4 Results and Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 83 5.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 vii CHAPTER 6 OPTIMIZATION OF PFE MODEL WITH RANDOM SEARCH . . . . 87 6.1 Model Description . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 6.2 Model Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 6.3 Optimization strategy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89 6.4 Iterations and Computation Time . . . . . . . . . . . . . . . . . . . . . . . . . 93 6.5 Results and Discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95 6.6 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112 CHAPTER 7 SUMMARY AND CONCLUSION . . . . . . . . . . . . . . . . . . . . 113 7.1 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113 7.2 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114 7.3 Future Work Opportunity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 BIBLIOGRAPHY . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117 APPENDIX A SUPPLEMENTARY EXPERIMENTAL RESULTS . . . . . . . . . . . 129 APPENDIX B PYTHON SCRIPTS . . . . . . . . . . . . . . . . . . . . . . . . . . . . 149 APPENDIX C PUBLICATIONS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 161 viii LIST OF ABBREVIATIONS ANN Artificial neural network AUCE Area under the curve error CoV Coefficient of variance CT Computed tomography DET Deterministic finite element model DFRP Discontinuous fiber reinforced polymer DIC Digital image correlation ESDM Extended strength distribution model FE Finite element GSDM General strength distribution model HPC High-performance computer MAE Mean absolute error MC Monte Carlo ML Machine learning MPFE Mean peak force error MPFFE Mean post-failure force error MRFE Mean residual force error OF Objective function PFE Probabilistic finite element PFSE Peak force or stress error REG Regular unimodal Weibull distribution RMSE Root mean squared error RES1T Residual tensile strength RES1C Residual compressive strength ix RS Random search RVE Representative volume element SDE Strain or displacement error at peak force/stress SF Surrogate function SMC Sheet molding compound x LIST OF SYMBOLS 𝑎 Activation function in ANN 𝑏 Bias in ANN 𝐷 Distance between V-notches for shear test 𝑑 Displacement of specimen 𝐹 Force 𝑓 Frequency 𝐻 Hyperparameter ℎ Intermediate variable for a hidden layer in ANN 𝐼 Number of integration points in FE model 𝐿 Length of a specimen 𝑙 Loading span for 4-pt bending test 𝑚 Weibull modulus or shape parameter 𝑛𝑚𝑎𝑡 Number of discretization or material models from Weibull distribution. 𝑛𝑖𝑡𝑒𝑟 Number of iterations in PFE simulation. 𝑃 Probability 𝑝 Mixture weight 𝑟 Distance between support and loading roller 𝑆 Survival probability 𝑠2 Mean variance 𝑇 Thickness or height of a specimen 𝑡 Time 𝑉 The volume of the specimen under stress 𝑊 Width of a specimen 𝑤 Weights in ANN 𝑋 Strength or Peak stress 𝑋𝑇 Tensile strength xi 𝑋0 Characteristics strength 𝑋𝑢 The threshold stress below which the material will not fail. 𝑥 Input variables in ANN 𝑦 Output variable in ANN 𝛾 Shear strain 𝜖 𝑚𝑎𝑥 Strain limit in FE model 𝜖 Normal strain 𝛿 Mid-span deflection 𝜎 Normal stress 𝜇 Mean standard deviation in experimental results 𝜏 Shear stress xii CHAPTER 1 INTRODUCTION Improvement of fuel efficiency in automobiles has become increasingly significant in recent years due to international regulations on carbon footprint reduction. To meet this, weight reduction is a key factor, as it directly affects the energy required to move a vehicle. While advanced metallic al- loys have been used to reduce weight, their poor specific strength limits their effectiveness, and their high manufacturing costs make them less attractive to the consumer market. Composite materials offer a superior alternative in this regard. The use of composite materials in the aerospace industry and high-performance race cars has been long successful. This has led to their integration into au- tomobiles and mass transportation industries. This is an area where discontinuous fiber-reinforced polymer (DFRP) composites, such as Sheet Molding Compound (SMC) composites, are frequently used [1–3]. 1.1 Introduction to Sheet Molding Compound (SMC) Sheet Molding Compound composites or SMC composites are a specific branch of thermoset- ting DFRP composite materials. It has gained popularity in the automotive industry due to its excel- lent mechanical properties, superior corrosion resistance, and low manufacturing cost. Moreover, SMC composites can be molded into complex shapes and sizes [4, 5]. Its manufacturing process is relatively simple and efficient, and the speed of production can be compared to that of metal stamp- ing. These properties make SMC composites a highly attractive option for part-supplying industries in the automotive sector. The use of SMC composites is expected to increase in the coming years as manufacturers continue to seek lightweight and cost-effective alternatives to traditional materials such as continuous fiber-reinforced composites, metallic alloys, etc. SMC panels consist of the matrix material, such as unsaturated polyester or epoxy resin, re- inforced with chopped glass fiber bundles. The fiber bundles in SMC are randomly distributed throughout the matrix and typically range from 15-45 wt% of the composite. The length of the fiber bundles is typically around 20-30mm, with approximately 200 fibers in each bundle [6]. The 1 manufacturing process of SMC is characterized by mold flow, which can result in microstructural changes such as the evolution of spatial and orientation distributions of fiber bundles and the forma- tion of voids. The mold flow pattern has a significant impact on the final physical and mechanical properties of the molded parts, as well as their surface characteristics, such as pinholes, craters, and surface waviness [7–9]. Further details about SMC’s manufacturing technique and morphology can be found in section 2.1.4. The mechanical properties of the Sheet Molding Compound (SMC) exhibit significant varia- tion in different loading conditions, regardless of the manufacturing process. For instance, SMC demonstrates significantly higher strength in 3-point bending tests than in uniaxial tension tests, as reported in [7–9]. The randomness of fiber distribution and the inclusion of flaws have been cited as reasons for this variation by Thomason et al. [8]. Additionally, Knight and Hahn [10] have ex- plained that any Discontinuous Fiber-Reinforced Polymer (DFRP) composite under uniaxial tension always fails at the largest flaw, and thus tensile strength only represents the lowest possible strength in a specimen. SMC also shows tension-compression asymmetry, similar to other polymers and their composites, which means it exhibits higher strength under compressive loading [11–13]. The variation between tensile strength, compressive strength, and flexural strength is commonly known as the size effect of strength [11]. Understanding these characteristics is critical for the development of optimized SMC materials and processes for the production of automotive components. 1.2 Motivations Although having long been used, the usage of SMC or any DFRP parts has endured limited use in structural components owing to the lack of a satisfactory crashworthiness model. The main difficulty in fully exploring the capabilities of the DFRP composite lies in the apparent uncertainty regarding their safety and predictability. In order for automotive industries to use DFRP in making vital components, an algorithm needs to be developed which can effectively model the randomness of strength under different loading conditions [10, 14, 15]. Furthermore, the model needs to be easy to implement and computationally efficient. While extensive research has been focused on developing models for composite materials, most 2 of them do not work well with DFRP. For instance, classical laminated plate theory [15], which is a reliable and computationally efficient approach for continuous fiber-reinforced composites, per- forms poorly with DFRP. Other approaches, such as quasi-classical laminated composite theory or laminated random strand method, can statistically assess sample size effects through moving- window analyses but cannot predict the variation of strength under different loading conditions. Statistical analysis and statistical strength theory has long been considered in composite designs [16, 17]. Probabilistic finite element (PFE) analysis based on statistical strength models has been employed to consider uncertainties and design reliability at every scale in composites. PFE has been used to predict the variations in ply-level material properties [18], component-level responses with uncertainties in material and geometrical properties [18–21], the damage and failure propagation of structures with open holes [22, 23], and the variability of the final product from manufacturing processes [24]. Computational micromechanics techniques have been used to predict the mesoscale variations of mechanical properties for SMC composites with chopped carbon fibers [25–30]. The predicted variations have been incorporated in numerical models and demonstrated in simulations of tensile experiments [25, 26, 28, 29]. However, little work has been done to examine the size effect of strength in FE simulations of structural components. 1.3 Objectives The objective of this work is to develop a methodology for crashworthiness simulations of SMC composite structures with the consideration of randomness in the mechanical properties and size effect on strength. The proposed model will cover a broad spectrum of stress distributions from uniform to deep gradients, as it might be found in a structure. This procedure should be practical, easy to calibrate, and computationally efficient for vehicle-level simulations. The feasibility of this approach will be examined with experimental results, probabilistic simulations, verification, and validation. Thus, the first objective is to investigate the uniqueness of material strength in the conventional orthogonal direction under tension, compression, and bending and to perform statistical analysis to determine the correlation between them. 3 The second objective is to develop a General Strength Distribution Model (GSDM) for SMC composites. For this purpose, multimodal Weibull distributions will be developed from statisti- cal models of tension, compression, and bending. This model would cover a broad spectrum of stress distributions from uniform to deep gradients, and thereby, the size effect is accounted for in simulations. 3rd objective is to develop an algorithm to implement GSDM in the FE model and optimize model parameters to reduce error. This will be done by creating randomization in FE according to GSDM and fitting PFE simulations with experimental results. Finally, optimized models will be verified and validated against different mechanical charac- terization tests. This includes PFE simulations of standardized and in-house tests, i.e., tension, compression, axial bending, open-hole tension, disk bending, etc. 1.4 Outline This thesis is organized as follows. Ch.1 introduces the problem and defines the scope of this thesis. Ch. 2 provides an overview of SMC manufacturing and its morphology. This chapter will also provide a literature review regarding the mechanics of DFRP and different modeling approaches. Chapter 3 presents an overview of mechanical characterization techniques and experimental results for SMC plaques. Chapter 4 presents an Extended Strength Distribution Model (ESDM) and probabilistic simu- lation techniques for SMC. This technique was cross-examined with standardized characterization tests. Chapter 5 presents a General Strength Distribution Model (GSDM) that includes tension and compression from uniaxial and bending. Additionally, a machine learning based optimization tech- nique is developed to estimate parameters for GSDM and material damage mechanics. Chapter 6 presents a Random Search based optimization technique to further improve accuracy from PFE simulations. Additionally, results from different optimization are compared for different mechanical tests. 4 CHAPTER 2 LITERATURE REVIEW This chapter presents an overview of current understandings relevant to Discontinuous Fiber Re- inforced (DFRP) Composites, especially Sheet Molding Compound (SMC). SMC’s manufacturing processes, typical morphology are reviewed in section 2.1. Literature studies on various modeling techniques for DFRP and SMC are discussed through section 2.2-2.4. Furthermore, a review on multimodal Weibull distribution and optimization techniques are presented through section 2.5-2.7. Recent studies regarding simulation with statistical consideration are reviewed in section 2.8. 2.1 SMC Manufacturing and Its Morphology The manufacturing process of SMC can be divided into two primary stages: manufacturing of semi-finished SMC sheets or prepreg and producing parts by compression molding of SMC prepreg. 2.1.1 Manufacturing of SMC Prepreg The first step in the manufacturing process of SMC prepreg is the mixing of resin, curing agent, and filler. The resin is typically a thermosetting resin such as polyester, vinyl ester, or epoxy. The filler is a material that provides bulk to the mixture and reduces the overall cost. Common fillers include calcium carbonate, talc, and silica. Afterward, chopped glass fiber bundles are distributed heterogeneously with a resin mixture. Once the reinforcement fibers are distributed, the entire mixture is sandwiched between two release films, as shown in Figure 2.1. The release films are used to prevent the mixture from sticking to the rollers during the compression process. Finally, it undergoes multiple roller compression to remove any trapped air. Then it is coiled up and stored for later stages. At this stage, it is called SMC Prepreg which can be stored safely for a few weeks to a couple of months, depending on storage temperature. For further detailed information regarding the SMC ingredients and a detailed description of the manufacturing process of the semi- finished SMC can be found in the book chapter of “Sheet Molding Compound” by Nicolais et al. [31], and in the book by Böhlke et al. [32]. 5 Figure 2.1 Manufacturing process of semi-finished SMC sheet or SMC Prepreg (Böhlke et al. [32]). 2.1.2 Manufacturing Parts from SMC Prepreg The process of manufacturing parts from SMC prepreg involves the compression molding of the semi-finished sheets. Firstly, a molding tool is prepared based on the desired shape and size of the part to be manufactured. The tool is heated to a specific temperature to ensure that the SMC prepreg material flows and takes the shape of the mold. Later, SMC prepregs are cut into smaller manageable sizes and stacked up, which is known as SMC charge. The size of the charge and the number of prepreg layers in the charge depends on the size and thickness of the final product. In reality, multiple charges are placed at different locations and sometimes at multiple layers depending on the complexity of shape and variation of thickness at different regions. Finally, stacked-up charges undergo compression at high temperatures and pressure. The heating and pressurization process is carried out for a specific time. 2.1.3 Typical Arrangement of Charge For Molding Parts of intricate design requires stacking up of SMC charges in various arrangement. These arrangements can be simplified into three basic categories, i.e., unidirectional flow, heterogeneous flow, and merging flow. 6 • If the charge is placed in the middle, as shown in Figure 2.2a, then the resin-fiber mixture would flow in every direction. Fiber orientation can still be considered as random uniform in this case. Plaques manufactured from this manufacturing technique is known as standard plaque. Experimental results of standard plaques demonstrated macroscale planar isotropic behavior, as to be discussed in Chapter 3. • When the charge is placed at either side of the mold, as shown in Figure 2.2b, the resin- fiber mixture flows predominantly in direction-1. This results in the reorientation of fiber and macroscale in-plane anisotropy, which is experimentally shown in Chapter 3. Plaques manufactured with this technique is named flow plaque. • Merging flow is obtained when two charges are placed at opposite sides, mold flows toward the center, and merges as shown in Figure 2.2c. The resulting part is called knit plaque. This arrangement of charge typically yields very poor mechanical properties at the flow merging section or knit line. Knit plaque is considered the least desirable SMC. 2.1.4 Morphology of SMC Several research works have been conducted to determine fiber orientation and microstructures of SMC plaques using destructive and non-destructive characterization techniques. Charring and weighing techniques are some of the commonly used destructive methods to determine average fiber mass fraction, but they can’t be used to obtain fiber orientation since the fiber-bundle network gets destroyed [33, 34]. Optical or electron microscopy is a common form of non-destructive technique. It is, however, only useful to analyze the morphology of SMC at fractured surfaces as it can’t provide much in- sight into fiber distribution. Micro-computed tomography (μCT) is a more suitable non-destructive characterization method to obtain microstructure and inner fiber architecture [34–36]. Le et al. [37] performed a detailed microstructure analysis. They reported that glass fiber bun- dles in the SMC are elliptically cross-sectioned (major axis ≈0.66mm and minor axis 0.06mm). Each fiber bundle contains around 200 fibers of diameter ≈15µm which is shown in Figure 2.3a and 2.3b. Their 2D micro-imaging shows typical fiber bundle distribution and a random presence of 7 Charge Charge 2 2 Mold flow direction Mold flow direction 1 1 (a) (b) Charge 90° 2 Mold flow 2 0° direction 1 1 (c) (d) Figure 2.2 Schematic view of typical SMC molding from stack of semi-finished sheet (a) standard plaques, (b) flow plaques, and (c) knit plaques. (d) Sampling plan from plaque at orthogonal direc- tion. 0° is treated as longitudinal direction and 90° as transverse direction. voids. Microstructure analysis of SMC performed by Kehrer et al. [6] is shown in Figure 2.3c and 2.3d. Their 3D Micro-computed tomography (μCT) exhibited in-plane randomness of fiber bun- dle distribution (Fig 2.3c). Additionally, fiber orientation varies through the thickness, as shown in Figure 2.3d. Fiber bundles at the top and bottom surface were mostly curled up, while in the middle section, they were mostly straight. Kehrer et al. [6] and Le et al. [37] concluded that SMC would be in-plane anisotropic when mold flow is dominant along a particular direction. However, when a fiber-resin mixture flows unidirectionally, it results in in-plane isotropic mechanical properties. 2.2 Micromechanics and Multiscale Models of SMC As shown in the previous section, SMC exhibits complex microstructure and random orienta- tion of fibers and voids in the matrix phase. Its mechanical characteristics are influenced by the interactions between these components. Micromechanical and multiscale have been employed in this regard to predict these characteristics and to optimize their properties for various applications. 8 (a) (b) (c) (d) Figure 2.3 (a), (b) 2D micrographs of SMC showing fiber bundle distribution and random presence of voids in the matrix (Le et al. [37]) (c) 3D μCT scan of SMC composite performed by Kehrer et al. [6] exhibits in plane randomness of fiber bundles (d) Fiber Orientation at different thickness [6]. 9 Micromechanical models of SMC aim to describe the mechanical behavior of the material at the microscale level. These models typically use numerical techniques, such as Finite Element (FE) analysis, to simulate the deformation and failure mechanisms of the composite. The matrix and fiber phases, as well as the interfacial regions, are modeled by a representative volume element (RVE) [11, 12, 38]. Micromechanical models have been used to study the effect of various parameters, such as fiber content, fiber orientation, and matrix properties, on the mechanical behavior of SMC. Multiscale approaches of SMC aim to describe the mechanical behavior of the material at mul- tiple length scales, from the microscale to the macroscale. Typically, these models use homoge- nization techniques to link the microstructure of the composite to its macroscopic behavior. The properties of the microscale constituents are used to calculate effective properties at the macroscale. Some of the more commonly used homogenization methods include the Voigt, Reuss, and Mori- Tanaka homogenization schemes. The Voigt homogenization scheme assumes a constant mean strain, while the Reuss homogenization scheme assumes constant mean stress, and they represent the upper and lower bounds of the solution [39]. The Mori-Tanaka homogenization scheme as- sumes that the mean stress in the matrix of the composite is constant [40]. These homogenization schemes are suitable for composites with continuously aligned fiber structures. Extensions have been made to these approaches for misaligned fiber composites [39, 41]. While micromechanics and Multiscale Models models demonstrated promising results for uni- directional composites with easily identifiable RVE, their usage in SMC is very limited due to its intricate, random microstructures. Moreover, modeling the micromechanical structure of the entire component using this approach is not feasible for large components such as automobiles as it would require enormous elements in FE mesh. Performing a FE simulation on such a structure will require immense computational resources. Finally, these approaches are typically more computationally expensive compared to macro-scale phenomenological models. 10 2.3 Statistical Strength Theory 2.3.1 Origin of the Statistical Strength Theory Classical mechanics has been long used to describe engineering materials by internal stress, young’s modulus, density, etc., assuming that materials are homogeneous. In a nutshell, these parameters are assumed to be a definite number for a certain material and less likely to be influenced by its structure. This model has demonstrated sufficient accuracy to determine the characteristics of ductile materials which have an ordered lattice structure. However, brittle materials don’t exhibit a simple, ordered lattice structure like ductile materials; their structures are composed of numerous defects, i.e., dislocations, secondary phases, pores, fissures, and inclusions. Due to the complexity of the microstructure and the uncertainty of the property of the constituents, a set of nominally identical specimens tested under the same condition would fail at different strengths, which is very difficult to estimate with classical mechanical theories [42]. Statistical strength theory, on the other hand, considers the inhomogeneity of structures and links the material strength to its imperfections. According to this theory, materials are assumed to be composed of discrete volume elements whose characteristics are related to porosity, cracks, and other imperfections. Imperfections are randomly distributed throughout the body, which renders the microscopic stress concentration, propagation of cracks, and overall failure. Unlike the contin- uum theory, where the ultimate failure is deterministic, statistical strength theory describes it using probability due to the randomness of structural flaws. By relating the strength to the material struc- ture, this theory attempts to bridge the gap between the microscopic and continuum approaches to fracture mechanics. 2.3.2 Different Types of Statistical Theories Several statistical models such as Weibull, normal (Gaussian), log-normal, gamma, and gen- eralized exponential distributions have been investigated by researchers to describe the strength distribution in materials. Basu et al. [43, 44] summarized the most common statistical theories and compared their accuracy to represent the strength distribution of brittle materials. They analyzed different materials ranging from extremely brittle such as solid glass, to relatively tougher mate- 11 rials, like ceramics and its composites, to determine strength reliability using probability models. Based on their findings, the Weibull distribution represents the closest match to the experimental results, followed by the normal distribution. Lu et al. [45] also presented a relative comparison between Weibull and normal distributions. They stated that the two-parameter Weibull and normal distribution fits better than the three-parameter Weibull distribution for brittle materials like ceramics. They also noted that normal distribution could potentially introduce non-physical negative tensile strengths. Elimination of these irrelevant data points would disrupt the symmetry of strength distribution around the mean value [16, 46]. This limitation makes the normal distribution less versatile, and special attention is needed to at- tenuate the biased results in different cases. In contrast, the implementation of the Weibull model is much straight forward. Additionally, the Weibull model can accurately estimate flexure strength from tensile strength [10]. Other statistical models were not employed in this regard. 2.4 Weibull Theory of Strength 2.4.1 Weibull Theory Derivation Weibull strength theory [47] dictates that any brittle material can be considered as a chain, and it will always fail at its weakest link. So, the strength of a brittle material represents the strength of the weakest link. The Weibull model assumes that defects are randomly distributed throughout the body. This random flaw distribution creates the provision of merging multiple flaws into a larger flaw. This results in defects of random sizes throughout the specimen (can be seen in the micrograph of brittle materials in Figure 2.3a and 2.3b). As the volume of a material increases, the probability of having multiple weak links or having a larger flaw also increases. Based on this theory, the probability (S) of a specimen’s strength would exceed its threshold strength can be described as,  ∫  𝑚  𝑋 − 𝑋𝑢 𝑆 = 𝑒𝑥 𝑝 − 𝑑𝑉 (2.1) 𝑣 𝑋0 12 where, 𝑋0 = Characteristics strength 𝑚 = Weibull modulus or shape parameter 𝑉 = The volume of the specimen under stress 𝑋𝑢 = The threshold stress below which the material will not fail. 2.4.2 Application of Weibull Theory A major application of Weibull strength theory is in the field of brittle fracture. The Weibull theory has been applied to different materials, including ferritic steels [48], ceramics [11], rocks [45], dental materials [45], fibers [49], and composites [50]. The Weibull shape parameter is closely related to the flaw size distribution, and hence it can define the brittle fracture of a material and be used to estimate the crack propagation under load [51, 52]. Xie and Gao [42] investigated the damage evolution in rocks under load using video microscopy and scanning electron microscopy and correlated its strength and crack propagation with fracture statistics which fits properly with the Weibull model. They argued that the Weibull shape parameter could represent the structural properties and crack distribution in a rock. Another major application of the Weibull theory is in the probabilistic design and reliability of engineering structures [53–55]. Engineering structures contain defects. Under service loads, these defects may grow into micro cracks and initiate failure. The design based on classical mechanical theories with deterministic properties cannot link the material nonuniformity to the reliability of the structures. Besides large structures, Weibull theory is also effective in predicting the strength of micro-electro-mechanical system design [56]. These microscale devices demonstrated significantly higher strength compared to bulk counterparts, as predicted by the Weibull size effect. 2.4.3 Application of Weibull Theory in DFRP DFRP consists of discontinuous fibers in a polymer matrix. Unlike continuous fiber-reinforced composites, where the load is carried directly by fibers, DFRP transfers load from one fiber bundle to another through the matrix and fiber/matrix interfaces. Therefore, its mechanical property is greatly 13 influenced by the fiber length, fiber distribution, and flaw distribution in the matrix. Fiber bundles and flaws in the matrix are randomly distributed throughout the whole part as well as through- thickness [40, 41]. Consequently, the scatter in the strength of DFRP composites is much larger than that of continuous fiber composites, and statistical analysis has to be considered. The Weibull model has been successful regarding this by describing its strength with probability distribution [10]. Another issue with DFRP composites is the large difference between the tensile and flexural strengths. Its flexural strength is often 60-80% higher than the tensile strength. The statistical strength theory based on the Weibull model can explain this difference [10, 57]. In uniaxial ten- sion, the entire gage section of the specimen is under uniform stress. The tensile strength will be determined by the largest flaw in the gauge section. In the case of 3-point flexure tests, the stress is concentrated in a small effective volume under the middle roller on the surface. In this case, the flaw has to be exactly under the middle roller and on the surface. The possibility of having a flaw in this small region is infinitesimally low, and hence it would result in significantly higher strength. Mathematically, this phenomenon can be explained by the Weibull model size effect on strength. Under uniaxial tension, the specimen is under uniform stress and threshold strength 𝑋𝑢 = 0. Equation 2.1 can be written as,   𝑚 𝑋 − 𝑋𝑢 𝑆 = exp −𝑉 (2.2) 𝑋0 If the tensile test is carried out for two specimens of the same material but different volumes V1 and V2, the two specimens would have the same probability of survival. According to Equation 2.2, we have:   𝑚   𝑚 𝑋1 − 𝑋𝑢 𝑋2 − 𝑋𝑢 𝑆 = exp −𝑉1 = exp −𝑉2 (2.3) 𝑋0 𝑋0   𝑚1 𝑋1 𝑉2 =⇒ = (2.4) 𝑋2 𝑉1 The effective volume for the uniaxial tension test is the entire specimen volume (except the grip). However, according to the Weibull model, the effective volume for 3-point bending and 14 4-point bending test can be described by Equation 2.5 and 2.6 respectively as below: 𝑊 𝐿𝑇 𝑣𝑒 = (2.5) 2(𝑚 + 1) 2   𝑙 𝑟 𝑣 𝑒 = 𝑊𝑇 + (2.6) 2(𝑚 + 1) (𝑚 + 1) 2 Where, 𝑊= width of specimen 𝐿= length of specimen 𝑇= thickness of specimen 𝑙= loading span for 4-pt bending 𝑟= distance between support and loading roller Equation 2.4 is the so-called Weibull size effect equation.According to Equation 2.4, specimens with different effective volumes would have different tensile strengths. It has been used to predict the tensile and flexural strengths of composite specimens of different sizes [14, 58]. Wisnom et al. [59] and Cattell et al. [14] demonstrated that the tensile strength, 3-pt flexure, and 4-pt flexure strength fit well with the Weibull size effect when considering effective volume as shown by Figure 2.4. Bullock [57] developed a general method to determine the tensile strength of DFRP composite from inexpensive tests such as the fiber tow test or the 3-point flexure test by utilizing the Weibull size effect. They examined two materials with different Weibull shape parameters, and the predicted tensile strengths from the two tests were within 5% accuracy. Nader et al. [21] also performed similar studies for multiphase materials and reported that the Weibull model was successful in accounting for the size effect in tension specimens and variation of strength between tension and flexural tests. 2.4.4 Significance of Parameters in Weibull Model In Weibull distribution, the shape parameter is considered the most important measurement, which determines the randomness of the material properties. A smaller shape parameter specifies a larger scatter in material properties and vice versa. The shape parameter can be determined 15 Figure 2.4 Relationship between strength and the volume under stress for composites (Cattell and Kibble [14]). from experimental results and utilized to determine a material’s behavior under different loading conditions. Jayatilaka et al. Jayatilaka and Trustrum [51], Trustrum and Jayatilaka [52] considered the shape parameter of the Weibull model as a material property, and Xie et al. [42] argued that the shape parameter is related to microstructure or its flaw distribution. They successfully correlated crack propagation in materials with their shape parameter. Their investigation suggests that the shape parameter can account for material properties as well as its microstructure. The shape parameter, however, depends on how sampling and tests were performed. Xu et al. [60] investigated the influence of sample size on the shape parameter, and confidence level and summarized that the sample size should be more than 20 to minimize the error in the estimated parameter. Similar works can be found in [61–63], which confirmed that a sample size smaller than 20 would not accurately determine the shape parameter. In order to establish a Weibull model that properly describes the material property, it is very important to have sufficient experimental data points. 16 2.5 Multimodal Weibull Distribution The Weibull distribution is a common form of statistical analysis to assess reliability in engi- neering and medical studies [64–68]. However, its assumption of a unimodal distribution is not always sufficient for real-world cases such as materials with complex microstructures. Multimodal Weibull distribution provides an alternative solution to this problem. The multimodal Weibull distribution is an extension to unimodal Weibull distribution. It allows any possibility of having multiple modes in a sample population. In the simplest form, multimodal Weibull distribution can be considered as a mixture of several unimodal Weibull distributions. It can be defined as below: Õ𝑛 𝑆(𝑋) = 𝑝𝑖 .𝑆𝑖 (2.7) 𝑖=1 Õ𝑛 𝑝𝑖 = 1 (2.8) 𝑖=1 "   # Õ𝑛 𝑋𝑖 𝑚(𝑖) 𝑆(𝑋) = 𝑝𝑖 . exp − (2.9) 𝑖=1 𝑋0 where, 𝑆(𝑋) = survival probability for multimodal distribution 𝑝𝑖 = mixture weight 𝑆𝑖 = survival probability for unimodal distribution 𝑋0 = characteristics strength 𝑚 = Weibull modulus or shape parameter Mixture weight (𝑝𝑖 ) is determined from diversity in a population. Several methods have been developed to estimate this parameter. They can be broadly classified in three categories [69] i.e. • Graphical methods • Statistical methods • Combination of the above two methods 17 Graphical methods are useful when an initial rough estimation is required. But, it has limited use in most real-life cases as there is no well-developed statistical theory for determining asymptotic properties. Statistical methods are better suited as a general-purpose estimation technique [70, 71]. It works well with all kinds of models and data. Statistical methods can be done by following techniques: • Bayesian approach • Pearson’s method of moments • Maximum likelihood • General curve fitting • Hybrid methods etc. In most cases, model parameters for multimodal Weibull distribution can’t be directly mea- sured from experimental results as the diversity could be unpredictable. Instead, it is estimated from an iterative process [72]. These approaches are relatively easy to apply for bimodal Weibull distributions but become complicated for trimodal Weibull distributions. Multimodal Weibull distribution covers a broad range of phenomena, including material sci- ence, mechanical system, wind energy, financial studies, etc. Bourahli et al. analyzed the strength data of natural Diss fibers and investigated the compatibility of various probability distribution functions with their experimental results [66]. Their investigation showed that the bimodal Weibull distribution provides the least error in comparison to other distributions. Watanabe et al., A re- search group from the Institute of Power Engineering, provided an overview of statistical analysis for wind power [73]. Their review paper shows that a mixture of two Weibull distributions provides fewer relative errors in determining the annual wind power density and energy production. Similar results are found for brittle materials [74], plasma coating [75], SiC fibers [76], carbon fibers [77], titanium alloy [78], composite structures [79], etc. 2.6 Machine Learning and Statistical Analysis Machine Learning (ML) is a field of artificial intelligence. In simple terms, it can be considered a statistical model that can learn from a given data set and automatically improve its accuracy 18 Inputs 𝑥1 Hidden layer Output 𝑤1 𝑥2 𝑤2 Í 𝑎 𝑦 .. .. . . 𝑤𝑛 𝑤 𝑖 𝑥𝑖 + 𝑏 𝑥𝑛 Figure 2.5 Schematic representation of simple artificial neural network (perceptron). through a rigorous algorithm. It has broad applications ranging from simple regression analysis to autonomous driving. One of the common subsets of ML algorithms is Artificial Neural Networks (ANN). It is in- spired by biological neural networks that constitute animal brains. Its simplest form is known as a perceptron, composed of a single neuron where the linear transformation of input variable occurs [80–82]. It can be described by Equation 2.11 and illustrated as Figure 2.5. 𝑦 = 𝑓 (𝑥) = 𝑎(𝑤 𝑖 𝑥𝑖 + 𝑏) (2.10) 𝑦 = 𝑎(𝑤 1 𝑥 1 + 𝑤 2 𝑥 2 + ... + 𝑤 𝑛 𝑥 𝑛 + 𝑏) (2.11) where, 𝑎= activation function 𝑏 = bias 𝑤 = weights ANN with 3 hidden layers can be illustrated as Figure 2.6. ANN with 𝑛 hidden layers can be represented by following equations. 19 Input layer Hidden layers 𝑥1 Output layer ℎ1(1) ℎ1(2) ℎ1(3) 𝑥2 𝑦1 ℎ2(1) ℎ2(2) ℎ2(3) .. 𝑥3 .. .. .. . .. . . . 𝑦𝑖 . ℎ𝑚(1) ℎ𝑚(2) ℎ𝑚(3) 𝑥𝑛 Figure 2.6 Schematic diagram of ANN with 3 hidden layers. ℎ𝑖(1) = 𝑎 (1) (𝑤 𝑖(1) (1) (1) 𝑗 𝑥 𝑗 + 𝑏𝑖 ) ℎ𝑖(2) = 𝑎 (2) (𝑤 𝑖(2) (1) (2) 𝑗 ℎ 𝑗 + 𝑏𝑖 ) . . 𝑦 = 𝑎 (𝑛) (𝑤 𝑖(𝑛) ℎ𝑖(𝑛−1) + 𝑏 (𝑛) ) Where, ℎ𝑖 =intermediate variable for a hidden layer 𝑥 𝑗 =input variables 𝑦 =output variable ANN can be used for both classification and regression tasks. In regression, the ANN is trained to predict a continuous output variable based on a set of input features. Training in ML is a process where weights 𝑤 𝑖 𝑗 are estimated from an input data set. This is typically done using an optimization algorithm such as gradient descent. Once ANN is trained, it can be used to make predictions on a new data set. 20 (a) (b) Figure 2.7 Example of regression analysis with ANN with (a) one hidden layer (b) 3 hidden layers. ANN with one hidden layer, as shown in Equation 2.11 can fit with linear data as well as some simple nonlinear data. However, it doesn’t work well with every nonlinear data. One of the examples is shown in Figure 2.7a. Here, the blue line represents the data ANN was trained on, and the red line represents the data predicted from ANN. The model was able to capture the average trend of the input data, but the prediction wasn’t precise. When the same data was used to train an ANN with three hidden layers, it was able to predict very accurately, as shown in Figure 2.7b. In engineering applications, ANN is being adopted for manufacturing processes, optimization of FE models, structural engineering, estimation of parameters for stochastic processes, etc. Luo et al. [83] used an ANN-based optimization strategy to minimize errors in FE analysis. Their opti- mization techniques were able to achieve less than 1% error. In a study by Gholizadeh et al. [84], a machine-learning approach was used to predict the behavior of steel structures under extreme load- ing conditions. The authors trained a neural network on a dataset of simulations of steel structures under different loading scenarios and used the network to predict the behavior of new structures under similar conditions. The results showed that the neural network was able to accurately predict the behavior of new structures and outperformed traditional FEA models in terms of computational time and accuracy. Liu et al. [85] provided a comprehensive overview of ANNs in the constitutive modeling of composite materials. Their review paper enlisted some promising results regarding the improvement of accuracy and efficiency of FE modeling for composites. 21 Figure 2.8 Methodology of the blind random search algorithm (Ranazzi and Pinto [87]). 2.7 Random Search Optimization Random Search (RS) is a common optimization technique in statistical analysis to estimate the optimal values of a set of parameters [86]. The approach involves randomly selecting values for the parameters within a predefined range and then assessing their performance using a designated evaluation metric, such as accuracy or error rate. The process is then repeated with different ran- domly sampled parameter values until the best combination of parameters that yielded the highest performance metric is identified. Generic RS algorithm can be illustrated [87] as Figure 2.8. One of the key advantages of RS optimization is its simplicity and versatility. It is easy to implement and does not require any prior knowledge of the parameter space. Additionally, RS op- 22 timization isn’t susceptible to converge at a local optima as data sampling is random independent [88]. This makes it an ideal technique for exploratory research or when working with complex models that have many parameters. However, RS can be slow and computationally expensive, par- ticularly for large parameter spaces. Despite its drawbacks, random search optimization is powerful and widely used for finding the optimal set of parameters for complex models [89]. Keith et al. compared the performance of random search optimization to other popular optimiza- tion techniques, such as genetic algorithms and particle swarm optimization, for parameter tuning in support vector machines. The authors found that random search optimization outperformed the other methods in terms of accuracy. Similar studies have been done by a research group from the School of Information and Communication Technology, Australia. They remarked that RS opti- mization is a reliable and accurate optimization scheme when computation cost is not concerned. Andradóttir investigated the potentiality of RS in the optimization of simulation and concluded that it can always lead to global optimum system design [90–92]. 2.8 Simulations with Statistical Consideration Several types of research have been conducted to develop phenomenological modeling tech- niques to simulate the stochastic nature found in a wide range of fields such as material science, engineering, finance, particle physics, computational biology, etc. These simulation techniques are generally referred to as probabilistic simulation. Among them, Monte Carlo (MC) simulation is commonly used [93]. Monte Carlo simulation is a computational technique used to estimate the outcome of a stochas- tic event in a system by generating from its probabilistic model [94]. It does so by generating random samples of inputs and calculating the resulting outputs. Predictions from MC can be reliable and accurate as long as the underlying probability distribution is modeled accurately. Nonetheless, a broad range of problems has been effectively resolved with MC. In composite mechanics, MC has been used since the 1970s MC technique has been used suc- cessfully to incorporate statistical strength simulation [38, 95–100]. But, they were mostly limited to micromechanics-based models. One of the earliest applications of Monte Carlo simulation in FE 23 analysis of mesoscale composites demonstrated by Xia et al. [101–104]. The authors developed the methods to study the quasi-static and dynamic failure process of unidirectional composites in tension, and the predicted results agree well with the experiments. In recent years, Nader et al. [19–21] applied MC simulation in PFE analysis for different tests on multidirectional woven E-glass/vinyl ester polymer matrix composites. The authors treated the elastic properties and strength of the materials as random variables and generated sample sets of data parameters using the corresponding probability distribution. Their models predicted results met sufficient reliability with 500 iterations. This could be helpful for the analysis of smaller specimens but not ideal for automotive crash simulation applications. Similar studies have been performed by Emil et al. The authors compared the MC and quasi- MC method to determine the fracture of composite laminates [22]. However, their approach was slightly different from Nader et al. [19–21]. They used random sampling for every element, while Nader et al. used the same material properties for all elements. The authors demonstrated that the quasi-MC technique is significantly faster than the MC method. Overall, these models are based on unimodal distribution derived from a particular mechanical property. Their usage is mostly limited to continuous fiber-reinforced composites where the size effect on strength isn’t as significant as DFRP. In the case of DFRP composites such as SMC, the size effect is inviolable, as to be discussed in Chapter 3. Such as SMC has 60-80% higher flexural strength in comparison to its tensile strength. This phenomenon can’t be properly predicted from models that are based on unimodal distribution. 24 CHAPTER 3 MECHANICAL CHARACTERIZATION OF SMC This chapter presents the development of experimental methods for the characterization of SMC composite. It contained chopped glass fiber strands roughly 25 mm long in a thermoset resin. The nominal fiber weight fraction determined by the burn-out method [105] was about 43%. The received SMC plaques had a nominal dimension of 300x600x3.25 mm. Mechanical properties were measured with specimens cut along the longitudinal and transverse directions of the plaques as shown in Figure 2.2d. Three types of plaques of glass fiber SMC were investigated, i.e., standard, flow, and knit plaques. However, the investigation in this chapter is mainly focused on standard plaques. For statistical analysis, mechanical tests were performed under tension, 3-point (3-pt) bending, and 4-point (4-pt) bending, which is tabulated in Table 3.1. Size effect on strength [10, 14, 57– 59, 106] is investigated by comparing different categories of tests. It wasn’t tested under the same loading configuration. In other words, the specimen size was kept constant for each type of test. All tests were performed under quasi-static displacement-controlled mode and at ambient temperature (∼22°C). For all specimens and all kinds of tests, the modulus of elasticity was calculated by taking the slope of the stress-strain curves ranging from 5-30% of the ultimate stress. For example, if a specimen has ultimate stress of 160MPa, then the modulus would be determined by taking a slope between 8 MPa to 48 MPa. Table 3.1 Mechanical characterization tests performed on the SMC. Test name Standard Nominal Dimensions Uniaxial Tension ASTM D3039 [107] 20x150x3.25 Uniaxial Compression In-house 20x75x3.25 3-pt bending ASTM D790-03 [108] 12.7x52x3.25 4-pt bending ASTM D7264 [109] 13x104x3.25 In-plane shear ASTM D7078 [110] 56x56x3.25 Cyclic shear In-house 56x56x3.25 25 3.25 25 100 25 20 Metal tab Specimen Figure 3.1 Schematic diagram for uni-axial tension test. 3.1 Uniaxial Tension Test Tensile tests were performed according to ASTM D 3039/D 3039M [107]. A schematic diagram of the uniaxial tension test is shown in Figure 3.1. The test was carried out with a MTS servo- hydraulic 810 material testing system with a 22.5 KN load cell. The specimen was gripped with hydraulic grips with a pressure of 9 MPa. Specimens were loaded to the machine by MTS’s built-in specimen protection system to eliminate unexpected preloading of the specimen. A laser extensometer was primarily used to measure the displacement, in addition to machine displacement. The laser extensometer’s output (analog voltage) was calibrated from a predeter- mined length. Output from the laser extensometer was connected to MTS 810 to record data syn- chronously with force and machine displacement. The test was performed at a quasi-static rate of 1.5 mm/min or strain rate of 0.0005/sec under the displacement-controlled model. 3.1.1 Sampling For each type of plaque, 6-8 specimens were taken along the longitudinal and transverse direc- tions. Locations of these specimens on the plaques are shown in Appendix A on Figure A.1 at page 129. The specimens were marked with an identification number and photographed to keep proper documentation. The specimens were cut with a wet diamond saw, cleaned, and dried before test- ing. After specimen preparation, their dimensions were measured. The thickness and width were measured with a micrometer and slide caliper, respectively. For each specimen, four dimension measurements were taken, and their average was recorded. The tensile specimen had a nominal 26 width of 20 mm and a thickness of 3.25 mm. The gage length for the laser extensometer was about 45 mm. Exact dimensions are given in Appendix A on Table A.1 through Table A.4 at page 129. 3.1.2 Calibration of Machine Displacement Crosshead displacement of the testing machine is typically very inaccurate due to slippage of the specimen in the grips. A laser extensometer can provide reliable and precise measurements in this regard. However, the laser extensometer’s measurement becomes unreliable near the peak load of the specimen as the laser tags tend to become loose. Additionally, failure doesn’t always happen in between laser tags. To avoid these limitations, calibrated machine displacement is used. The calibration was done by correlating the machine displacement with a laser extensometer. The machine displacement vs. the laser extensometer displacement for uniaxial tension is shown in Appendix A on Figure A.2 at page 130. The slope of the graph provides the calibration factor. Since the specimen’s slippage varies from test to test, the calibration factor was determined for each specimen. From Figure A.2 it is seen that the laser displacement provided a reliable reading until machine displacement reached about 1.3mm. After that, the laser tag tends to peel off due to its poor adhesion to the specimen. Therefore, a correlation factor was determined between machine displacement and laser displacement for the 0-1 mm range. The results are given in Table 3.2. Table 3.2 Calibration factor for machine displacement during tensile test of standard plaques. Specimen ID TT1 TT2 TT3 TT4 TT11 TT12 TT13 TT14 Calibration factor (mm/mm) 0.313 0.072 0.505 0.460 0.451 0.514 0.459 0.460 3.1.3 Digital Image Correlation The local strain in a SMC composite is influenced by the fiber orientation. Vice versa, the strain field can provide useful information about the fiber orientation of the SMC. To investigate this phenomenon, the digital image correlation (DIC) measurement was performed for standard plaques during the tensile test. For DIC measurement, the surface of the specimen was sprayed with fine speckles. SMC spec- imens are naturally grey to black in color. Therefore, they were coated with a flat black paint to make a uniform non-reflecting black surface. Later, white speckles were created on top of black 27 Figure 3.2 Virtual extensometer arrangement on specimen to determine poison’s ratio of SMC using DIC. paint using flat white spray paint. Speckles’ size ranged between 0.2-0.5 mm. A Stingray F-201B monochrome CCD camera was used to capture the image during tensile tests. The camera has a resolution of 1624x1234 pixels. All images were taken with a Zeiss 50 mm macro lens. A DC incandescent light source was used during imaging to improve contrast and reduce noise. Images from the camera were synchronized with the tensile testing machine. The synchroniza- tion was done by using the same data acquisition frequency and matching rupture. For example, the tensile test data acquisition frequency was set to 10 Hz, and the image capturing system was 28 175 TT1 TT2 TT3 TT4 TT11 TT12 TT13 TT14 TL11 150 TL12 TL13 TL14 TL27 TL28 125 Stress (MPa) 100 75 50 25 0 0 0.005 0.01 0.015 0.02 0.025 Strain (mm/mm) Figure 3.3 SMC’s uniaxial tension test stress-strain results for Standard plaques. set to 10 FPS. This technique enabled us to obtain data at the same time interval. The last step towards synchronization is to match tensile testing data with camera images. It can be done if the two systems can be started exactly at the same time, but that couldn’t be done as two processes were done on a separate system. An alternative route is to match a specific point, such as the point of fail- ure. The force measurement from the testing machine and the strain measurement from DIC would provide a sharp point when the specimen fails at the ultimate stress. By matching the failure point, we can correct the phase lag between the camera and the testing machine. This synchronization process was performed in the post-processing of the tensile test experiment. The same process was used to determine the Poisson’s ratio of the material. In this setup, eight virtual extensometers were created on the specimen in the longitudinal and transverse directions to obtain an average strain measurement along those axes, as shown in Figure 3.2. Results from DIC are shown in Appendix A on Figure A.9 at page 134. 3.1.4 Uniaxial Tension Test Results Figure 3.3a presents all tensile stress-strain curves measured in SMC standard plaques. Speci- mens TL11 to TL28 were along the longitudinal direction, whereas TT1 through TT14 were along 29 the transverse direction. From these results, the following are observed: (1) tensile stress-strain curves were highly scattered, and (2) tensile properties for standard plaques in the longitudinal di- rection and transverse directions were randomly distributed. In other words, tensile properties were indistinguishable between directions. Figure 3.4 demonstrates the strain distribution on the specimen (SMC standard plaques) during the uniaxial tension test. The strain field was obtained from digital image correlation. This figure clearly shows (1) the strain field was almost uniform at the lower strain level. (2) strain became gradually non-uniform as the specimens were stressed further. (3) Near peak load strain was highly concentrated near the area of failure. The strain field from DIC provided a good prediction of where the specimen might fail. Uniaxial tension test results for Flow plaques and knit plaques are given in Appendix A Figure A.3 to Figure A.6 at page 131. 3.2 Uniaxial Compression Test The MTS 810 material testing system was used to perform uniaxial compression tests. The test protocol was developed in-house. A schematic diagram of this test protocol is shown in Figure 3.5. In this protocol, the specimen was gripped directly by the hydraulic grips at its two ends. The length of the ungripped portion was kept sufficiently short to prevent global buckling of the specimen under compression. The ASTM standard compression test usually utilizes an anti-buckling fixture for this purpose. However, the standard setup prohibits the measurement of the post-peak behavior, which is critical to crash simulations. For each type of plaque, 6-10 specimens were taken in both longitudinal and transverse direc- tions. Figure 135 shows the locations of the specimens. The specimens had a nominal cross-section dimension of 15x3.25 mm. The measured dimensions are given in Appendix A on Table A.5 to Table A.7 at page 135. The laser extensometer and the machine crosshead displacement were used for displacement measurement. The machine displacement was calibrated the same way as in the tensile experiment. Calibration curves and results are given in Appendix A on Figure A.11 at page A.11. It indicates 30 Figure 3.4 Uniaxial tension strain field for standard plaques measured by DIC. 25 15 25 3.25 10 Grip area Specimen Figure 3.5 Schematic diagram for uniaxial compression test. 31 that the laser extensometer provided reliable readings until machine displacement was around 0.6 mm. Therefore, the correlation factor was determined for a displacement range of 0-0.4 mm. Figure 3.6a plots the compressive stress-strain curves for the standard plaques. Specimens 1 through 11 were along the longitudinal direction, whereas 12 to 22 were along the transverse direc- tion. As seen, the curves were randomly distributed in the stress-strain plots. The results suggest that the standard plaque can be considered isotropic under compression. The degree of anisotropy observed in these plaques is similar to the uniaxial tension cases. Figure 3.6b demonstrates the failure mode of SMC standard plaques during uniaxial compression tests. Compressive properties of flow plaques and knit plaques are shown in Appendix A on Figure A.12 to Figure A.14 at page 136. 3.3 3-Pt Bending Test The flexural tests under 3-point bending were performed according to ASTM D790-03 [108]. A Schematic diagram and nominal dimensions for the 3-pt bending test are shown in Figure 3.7. Seven specimens were taken from standard plaques along longitudinal, and six were along trans- verse directions. Specimens were prepared using a water-jet CNC cutting machine. The exact dimensions of the specimens were given in Appendix A on Table A.8 to Table A.10 at page 138. All Tests were performed with a MTS Insight material testing system. WTF IO-108 three-point bending fixture was used to mount specimens on the MTS Insight testing machine. This fixture had three rollers of 0.25-inch diameter to provide smooth translation motion of specimen. A quasi-static displacement rate of 1.5 mm/min was used for all tests. The reaction force from specimens was measured by a 10 kN load cell. Crosshead displacement was used to measure the specimen’s mid-span deflection. The flexural stress of SMC was determined using the following equation: 3𝐹 𝐿 𝜎𝑓 = (3.1) 2𝑊𝑇 2 32 250 0 2 3 6 8 9 10 11 200 13 14 15 16 17 18 19 20 Compressive Strength (MPa) 150 100 50 0 0 0.005 0.01 0.015 0.02 0.025 Strain (mm/mm) (a) (b) Figure 3.6 Uniaxial compression test for SMC standard plaques (a) Stress-strain plots (b) Failure mode. 33 3.25 12.7 𝐹𝑜𝑟𝑐𝑒 Specimen 𝐹𝑜𝑟𝑐𝑒 𝐹𝑜𝑟𝑐𝑒 2 2 26 26 52 Figure 3.7 3-pt bending test schematic diagram (ASTM D790). where, 𝜎 𝑓 = stress in the outer fibers at midpoint [MPa] 𝐹 = applied force [N] 𝐿 = Support span [mm] 𝑊 = width of the specimen [mm] 𝑇 = Height or thickness of the specimen [mm] Figure 3.8a presents the flexural stress-strain curves obtained under 3-point bending for standard plaques. Flexural test results showed a similar scatter as uniaxial tension results. For the standard plaque, specimens B1 through B8 were along the longitudinal direction, whereas B74 to B80 were along the transverse direction. Differences in flexural strength between longitudinal and transverse directions were indistinguishable similar to uniaxial tension or compression test results. However, 3-pt bending tests of standard plaques demonstrated two key differences in compar- ison to uniaxial tests, i.e. (1) Peak stress was significantly higher compared to its tensile or com- pressive strength, and (2) Specimens failed gradually during bending, whereas in uniaxial testing specimens reached to near zero stress after peak stress. In other words, specimens had very high 34 350 B1 B2 B3 B4 B6 300 B8 B74 B75 B76 B77 Flexure stress [MPa] B78 B79 B80 250 200 150 100 50 0 0 0.02 0.04 0.06 0.08 0.1 0.12 Flexure strain [mm/mm] (a) (b) Figure 3.8 3-pt bending test results for SMC standard plaques (a) Stress-strain plots (b) Specimen deformation and failure mode at different strains. 35 strain at rupture during 3-pt bending tests. Gradual failure of SMC under 3-pt bending can be clearly seen in Figure 3.8b. These sequential images were taken during the test. They demonstrate the progression of failure of the specimens at different strains. 2nd image in the sequence (𝜖 = 0.003) represents conditions at peak stress where the failure of the lower surface was visible. The last image in the sequence represents conditions at rupture. One important key point of failure mode is: failure started at the bottom surface where the fibers were at tension. This failure then gradually progresses toward the top surface. Even at rupture, there wasn’t any sign of failure on the top surface of the specimen. A similar failure mode was seen in all 3-pt bending cases. Therefore, it can be concluded that SMC’s failure during 3-pt bending tests was dictated by tensile failure. Flow plaques and knit plaques exhibited similar results. Results for Flow plaques and knit plaques 3-pt bending tests are shown in Appendix A on Figure A.15 to A.18 at page 139. 3.4 4-Pt Bending Test 4-pt bending tests were performed according to ASTM D7264 [109]. The schematic diagram and nominal specimen dimensions for the 3-pt bending test are shown in Figure 3.9. Tests were carried out with a MTS Insight material testing system and WTF IO-108 fixture. It was the same fixture as 3-pt bending setup, but it used four rollers instead of three rollers. A total of fourteen specimens were prepared from SMC standard plaques along the longitudinal and transverse directions. All tests were carried out in displacement-controlled mode at a displace- ment rate of 1.5 mm/min. The reaction force from the specimen was measured with a 10 kN load cell. Flexural stress was determined from the following equation: 3𝐹 𝐿 𝜎= (3.2) 4𝑊𝑇 2 36 3.25 13 50 Specimen 25 25 100 Figure 3.9 Schematic diagram for 4-pt bending test. Where, 𝜎 = stress in the outer surface in the load span [MPa] 𝐹 = applied force [N] 𝐿 = Support span [mm] 𝑊 = width of the specimen [mm] 𝑇 = thickness of the specimen [mm] ASTM D7264 standard requires mid-span deflection measurement to determine strain at the outer surface of the specimen. This experiment was done with a Laser extensometer. Mid-span deflection measurement was done by putting one laser tag on the top part of the fixture and another at the mid-section of the specimen. Later output from the laser extensometer was synchronously recorded with force measurement from the MTS machine. Finally, strain at the outer surface was determined from Equation 3.3. The correlation between laser displacement and machine displace- ment is shown in Appendix A on Figure A.19 at page 142. 48𝛿𝑇 𝜖= (3.3) 11𝐿 2 37 Where, 𝜖 = strain at the outer surface of the specimen [mm/mm] 𝛿 = mid-span deflection [mm] 𝐿 = Support span [mm] 𝑇 = thickness of the specimen [mm] Figure 3.10a shows the stress-strain curves for 4-point bending tests. These test results were obtained from SMC standard plaques. Specimen B67-B73 were cut along the longitudinal direc- tion from standard plaques, whereas B50 to B56 were along the transverse directions. 4-pt tests demonstrated similar anisotropy and scatter as other characterization tests. 4-pt bending tests of standard plaques had two major differences in comparison to two 3-pt bending tests: (1) Peak stress was lower compared, yet it was higher than its tensile strength, and (2) Specimens failed more rapidly after peak stress. Rapid failure of peak stress during 4-pt bending can be clearly seen in Figure 3.10b. These sequential images were taken during the test. They demonstrate the progression of failure of the specimens at different strains. 5th image in the sequence (𝜖 = 0.003) represents conditions at peak stress where the failure of the lower surface was identified. The last image in the sequence represents conditions at rupture. Similar to 3-pt bending, failure started at the bottom surface where the fibers were at tension. Then the failure rapidly progressed through the thickness towards the top surface. This failure then gradually progresses toward the top surface. There wasn’t any noticeable failure on the top surface due to compression. 3.5 In-Plane Shear Test The in-plane shear properties of the SMC were measured using the V-Notched rail shear method according to ASTM D7078/D7078M -12 [110]. Figure 3.11 shows a schematic diagram of the V- Notched specimen and its nominal dimensions. 38 250 B67 B68 B69 B70 B71 B72 B73 B50 B51 B52 B53 B54 B55 B56 200 Stress [MPa] 150 100 50 0 0 0.01 0.02 0.03 0.04 Strain [mm/mm] (a) (b) Figure 3.10 4-pt bending test results for SMC standard plaques (a) Stress-strain plots (b) Specimen deformation and failure mode at different strains. 39 T=3.25 90° 56 D=31 12.5 56 Figure 3.11 Schematic diagram for V-notched shear test. 3.5.1 Sampling For standard plaque, a total of 15 samples have been taken. Seven test specimens were taken in the transverse direction, and eight specimens were taken along the longitudinal direction. Their locations are shown in Appendix A in Figure A.20 at 143. Specimens were cut from the SMC plaque using a CNC water jet. The actual dimensions of the V-Notched specimens are given in Appendix A on Table A.11 to Table A.13 at page 143. 3.5.2 Test setup The test setup for the in-plane shear test is shown in Figure 3.12. The specimen was clamped using a V-Notched rail fixture, as shown on the left side of the experiment setup. This specific was manufactured by Wyoming Test Fixtures, INC. The right image in Figure 3.12 shows the entire experiment setup. A DC incandescent light source with an endoscope was used to illuminate the surface of the specimen to obtain noise-free images from the camera. All images were taken with Zeiss 50 mm macro lens. The camera was set on a tripod to get the same field of view. The specimen was clamped to the fixture according to ASTM specifications, i.e., each bolt of the fixture was tightened to 50Nm torque using a digital torque wrench. Tests were performed on a MTS-810 material testing machine with a displacement-controlled mode similar to the uniaxial tension test. 40 Figure 3.12 Experimental setup for Shear test of SMC. The shear stress 𝜏 is calculated using the following formula: 𝑃 𝜏= (3.4) 𝐻𝐷 where, 𝐹 = force obtain from load cell 𝑇 = thickness of the specimen 𝐷 = distance between V-notches as shown in Figure 3.11 3.5.3 Strain Measurement According to ASTM D7078 [110], the shear strain can be measured through two single strain gages oriented at ±45º directions. The shear strain 𝛾 is computed as 𝛾 = 𝜖+45 + 𝜖−45 (3.5) In this work, the strain measurement was performed using DIC. To be consistent with ASTM D7078, virtual strain gages oriented at ±45º directions were defined. As shown in Figure 3.13a total of 12 virtual strain gages were set between the notch tips in the gage section of the specimen. 41 (a) (b) Figure 3.13 (a)Virtual strain gage setup and (b) the measured strain history curves by these gages. 6 gages were oriented 45°, the other 6 at -45°. The strain history curves measured by these strain gages are presented in Figure 3.13b. The shear strain was calculated using the averaged strain value measured from the 45° and -45° gages. 3.5.4 Synchronization of DIC and Force During the failure process, some speckles on the specimen surface fell off. This affected the strain measurement in the final stage. To estimate the strain beyond the peak load, the machine displacement was used. For this purpose, the machine displacement was correlated with the DIC strain, as shown in Figure A.21 on Appendix A at page 144. The synchronization of the two systems was similar to that in the tensile experiment. The force data were recorded at 2 Hz, and images were taken at 2 FPS for all specimens. Figure A.21 provides a linear correlation between DIC strain and machine displacement until the specimen failure point. 42 Figure 3.14 In-plane shear test results for SMC standard plaques. 3.5.5 Results for In-plane Shear Test Figure 3.14 presents the measured in-plane shear stress-strain curves for SMC standard plaques. Results showed that SMC had significantly low shear strength in comparison to other mechanical tests. However, SMC demonstrated high residual shear strength (approximately 20%). A typical shear-strain field obtained from DIC for an in-plane shear test is presented in Figure 3.15. Results are shown for five different time steps during the tests. The image taken at 47s represents the strain field at peak stress. Here, high levels of strain (≥ 2.5%) were concentrated randomly between the notches. This strain concentration path was similar to its failure mods as shown in Figure A.22 on Appendix A at page 144. 3.6 Cyclic Shear Test MAT 297 has two plasticity parameters which are determined through a cyclic shear test [2]. In the current work, the cyclic test was performed for standard plaques. This test was carried out 43 Figure 3.15 Strain filed of in-plane shear test obtained by DIC. Table 3.3 The load levels in cyclic shear test. Cycle no. 1 2 3 4 5 Load (kN) 2 4 6 8 9.5 the same way as mentioned in Section 3.5.1 and 3.5.2. Only the loading condition was varied. In the cyclic test, the specimen was loaded by five consecutive loading and unloading cycles with incremental load levels, as shown in Table 3.3. The strain measurement for the cyclic shear test was carried out with DIC using the same method as described in Section 3.5.3 and 3.5.4. Stress-strain results from the cyclic shear test are plotted in Figure 3.16 and plasticity parameters are graphed in Figure A.26 on Appendix A at page 45. 3.7 Summary Mechanical properties of SMC standard plaques are summarized in Table 3.4. SMC standard plaques demonstrated significantly high flexural strength in comparison to their tensile strength or compressive strength. Flexural strength in 3-pt bending was 1.9 times higher, and 4-pt bending was 1.5 times higher compared to uniaxial tension. Flow plaques and knit plaques exhibited similar characteristics. Their mechanical properties are summarized in Appendix A on Table A.14 through Table A.17 at page 147. 44 100 cycle 1 Cycle 2 Cycle 3 Cycle 4 Cycle 5 Shear Stress (MPa) 80 60 40 20 0 0 0.005 0.01 0.015 0.02 0.025 Shear strain Figure 3.16 Cyclic shear stress-strain results for SMC standard plaques. Table 3.4 Summary of Mechanical Properties for SMC standard plaques. Tension Compression 3-pt 4-pt Shear Strength (MPa) 135.6±15.5 174.2±25.1 262.9±39.1 189.8±17.4 105.5±7.4 Modulus (GPa) 10.14±1.1 13.4±1.0 11.7±1.6 10.3±0.44 5.1±0.5 Strain at peak load 0.016±0.003 0.018±0.002 0.023±0.002 0.024±0.003 0.034±0.004 Additionally, SMC standard plaques exhibited in-plane isotropic mechanical properties for all characterization tests. In contrast, flow plaques and knit plaques were anisotropic. Typical stress-strain curves of SMC standard plaques for different mechanical characterization tests are plotted in Figure 3.17a. SMC specimens exhibited gradual failure during flexural load, whereas specimens failed abruptly after peak load under uniaxial tension or compression. Never- theless, all tests exhibited a large scatter in strength as shown in Figure 3.17b. 45 250 3-pt Bending Uniaxial Tension 200 Uniaxial compression 4-pt Bending In-plane Shear Stress [MPa] 150 100 50 0 0.00 0.02 0.04 0.06 0.08 0.10 Strain [mm/mm] (a) 350 300 250 Peak Stress [MPa] 200 150 100 50 0 Tension Compression 3pt Bending 4pt Bending Shear (b) Figure 3.17 (a)Typical stress-strain curves, and (b) Scatter of peak stress for standard plaques in various mechanical characterization test. 46 CHAPTER 4 PFE MODEL WITH RANDOM TENSILE STRENGTH This chapter presents the development of the Probabilistic Finite Element (PFE) model for SMC standard plaques. The PFE model is based on an extended strength distribution model (ESDM) in the form of a trimodal Weibull distribution. ESDM was developed from unimodal Weibull dis- tributions of uniaxial tension, 3-pt bending, and 4-pt bending. This strength distribution doesn’t consider the randomness of compressive and shear strength. ESDM was subsequently tested in PFE simulations of SMC tensile and flexural experiments following a randomization algorithm that considers the tensile strength as a random variable and assigns it according to ESDM across finite element models. PFE Model described in this chapter followed two major steps, i.e., the development of a statis- tical model for SMC and the implementation of randomness in the FE model. It is to be noted that this model didn’t use any optimization to fit the simulation results with experimental results. The work presented in this chapter has been summarized and published [111]. 4.1 Mechanical Properties of SMC Mechanical characterization of SMC standard plaques is described in Chapter 3. Among them, only uniaxial tension, 3-pt bending, and 4-pt bending test cases were considered for the PFE model described in this chapter. A summary of these characterization tests is given in Table 4.1. The ratios of the flexural strength to the tensile strength measured in the 3-pt bending and 4-pt flexural tests were 1.9 and 1.4, respectively. Table 4.1 Summary of mechanical properties for SMC standard plaques used in this chapter. Tension 3-pt 4-pt Strength (MPa) 135.6±15.5 262.9±39.1 189.8±17.4 Modulus (GPa) 10.14±1.1 11.7±1.6 10.3±0.44 Strain at peak load 0.016±0.003 0.023±0.002 0.024±0.003 Figure 4.1a compares the typical stress-strain curves, and Figure 4.1b presents the scatter in strengths for each type of test. As shown, the 3-pt bending yielded the highest strength, followed 47 350 250 3pt Bending Uniaxial Tension 4pt bending 300 200 250 Peak Stress [MPa] Stress [MPa] 150 200 100 150 50 100 0 50 0.00 0.02 0.04 0.06 0.08 0.10 0 Strain [mm/mm] Tension 4pt Bending 3pt Bending (a) (b) Figure 4.1 (a) The typical stress-strain curves of the SMC material obtained by uniaxial tension, 3-pt, and 4-pt bending tests. (b) The scatter in strengths. by 4-pt bending and uniaxial tension. The shapes of the stress-strain curves varied with the loading modes. During the tensile experiments, the specimens failed abruptly at the peak stress with very little residual strength. On the other hand, the flexural stress-strain curves displayed a gradual failure pattern, i.e., the specimens maintained some post-peak load-carrying capabilities after the peak stress. Furthermore, the specimens exhibited a more gradual failure under 3-pt bending than that under 4-pt bending. Figure 3.8b and 3.10b present a sequence of photos taken during 3-pt and 4-pt flexural experiments. As shown, the failure occurred in the region of the specimen under tension. There was minimal to no failure due to compression. The same failure mode was observed in all flexural tests. 4.2 Statistical Model of SMC Brittle materials display a much higher strength under compression than under tensile load [11–13]. In flexural experiments, brittle materials tend to fail at the region of the specimen under tension. Therefore, the flexural strength of a brittle material is related to its tensile strength but with a much higher value. According to Weibull’s statistical strength theory [112], the strength of a material is determined by its weakest point, i.e., the flaws. As the volume of the specimen increases, the probability of having a larger size flaw also increases. Assuming that the flaws are randomly distributed throughout the specimen, the probability of survival (S) of a specimen under stress 𝑋 is 48 described by the Weibull model as:  ∫  𝑚  𝑋 𝑆 = exp − 𝑑𝑉 (4.1) 𝑣 𝑋0 where 𝑋0 is the characteristic stress, 𝑚 is the shape parameter, is the threshold stress below which the material will never fail, and 𝑉 is the volume of the specimen under stress. 𝑚 is a measure of the scatter of the data. A smaller value of m means that the material has a large variation in strength and vice versa. The Weibull theory has been used to explain the discrepancy between the flexural and tensile strengths observed in composite materials [10, 57]. In uniaxial tension, the material is under uni- form stress. With the assumption of 𝑋𝑢 = 0, Equation 4.1 becomes the two-parameter Weibull’s model expressed as:   𝑚 𝑋 𝑆 = exp −𝑉 (4.2) 𝑋0 where 𝑉 is the volume of the tensile specimen and 𝑋 is the tensile strength. If the tensile test is carried out using specimens with two different volumes of 𝑉1 and 𝑉2 of the same material and assuming that the two have the same probability of survival, we have:   𝑚   𝑚 𝑋1 𝑋2 𝑆 = exp −𝑉1 = exp −𝑉2 (4.3) 𝑋0 𝑋0   𝑚1 𝑋1 𝑉2 = (4.4) 𝑋2 𝑉1 Equation 4.4 is the so-called Weibull size effect model, which provides a relationship between the volume under tensile stress, its strength, and the scatter of the strength. This equation essentially shows that a material with a large standard deviation (smaller 𝑚) would have higher strength when the volume of the specimen under higher stresses is smaller. In flexural tests, due to the presence of a stress gradient, only a portion of the specimen is under high tensile stress. Therefore, an effective volume is used. The equations for the effective volumes of the 3-pt and 4-pt bending configurations are provided in Refs. [10, 57]. From the effective volume, one can estimate the flexural strength 49 from the tensile strength. Vice versa, the known value of the apparent tensile strength can be used to calculate the effective volume. 4.2.1 Unimodal Weibull Strength Model The Weibull statistical analysis starts by making the survival probability plot, i.e., the Weibull plot. There are a number of ways to determine the survival probability. In this work, it was deter- mined using the median rank (M.R.) as follows [10, 113]: 𝑖 − 0.3 𝑆𝑖 = 1 − M.R. = 1 − (4.5) 𝑛 + 0.4 Using Equation 4.5, the survival probability was determined for each data point. The Weibull plots for strength were subsequently constructed for tensile, 3-pt, and 4-pt flexural tests, as presented in Figure 4.2. The Weibull plots were then fitted with the two-parameter Weibull model (Equation 4.6) using the Least Square method. The obtained characteristic strength (𝑋0 ) and shape parameter (𝑚) are summarized in Table 4.2. The Weibull models obtained from a single type of experimental data are labeled as the regular Weibull models.   𝑚 𝑋𝑖 𝑆𝑖 = exp − (4.6) 𝑋0 Figure 4.2 shows that the data from different tests follow different Weibull distributions. The tensile strength curve is on the left, followed by 4-pt bending and then 3-pt bending. Although some of the data points from different tests overlap at their lower and higher ends, an individual category of tests can cover only a portion of the strength spectrum. As observed in the experimental procedure, flexural specimens failed at the region under tension. Theoretically, to cover load cases with stress distributions either uniformly as in tensile tests or with deep gradients as in flexural experiments, we need to perform tensile tests with specimens of different sizes. To measure the tensile strength in a small volume, as at the stress concentration point in 3-pt bending, the specimen size will have to be infinitesimal. Such specimens will require a small-scale test. 50 1 Tension 0.8 3-Pt bending Probability of survival (S) 4-Pt Bending 0.6 0.4 0.2 0 50 100 150 200 250 300 350 Ultimate Stress (MPa) Figure 4.2 The tensile strength (red) and flexural strength measured by 3-pt bending (blue) and 4-pt bending (black) followed three different Weibull distributions. Table 4.2 Summary of Weibull parameters for SMC plaques (Regular Weibull Model). Data set for the model Tensile 3-Pt Bending 4-Pt Bending Weibull modulus, 𝑚 9.68 7.31 12.04 Characteristics strength, 𝑋0 143 280 260 4.2.2 Multimodal Weibull Distribution The Weibull statistical strength model has been successful in describing the characteristics of brittle materials and estimating the size effect on strength. However, unimodal Weibull distribu- tion has its own limitations. For example, manufacturing process variability, fiber heat treatment temperature, and variation of strain rate create a significant shift in strength distribution which does not fit with a unimodal distribution properly [114–116]. To overcome these limitations, multimodal Weibull distributions have been used to analyze the strength of fibers [66, 117], composite struc- tures [79, 118, 119], and ceramic coating [120, 121]. They have also been used to analyze the fatigue rating of alloys [78]. One of the common ways of modeling multimodal data is by using 51 a mixture of unimodal distributions [72, 122]. Trimodal mixture Weibull distribution is a special case of multimodal distribution that has been used to model heterogenous mechanical properties [117, 121]. Mathematically, a trimodal Weibull distribution can be represented as:   𝑚    𝑚    𝑚  𝑋𝑖 1 𝑋𝑖 2 𝑋𝑖 3 𝑆𝑖 = 𝑝 1 . exp − + 𝑝 2 . exp − + 𝑝 3 . exp − (4.7) 𝑋01 𝑋02 𝑋03 Õ3 with, 𝑝𝑗 = 1 (4.8) 1 where 𝑚 𝑗 , 𝑋0 𝑗 and 𝑝 𝑗 ( 𝑗 = 1, 2, 3) are the shape parameter, characteristic strength, and mixture weight, respectively. The mixture weight 𝑝 𝑗 is a parameter to be determined with an iterative pro- cess. A number of approaches have been used to determine the parameters in multimodal Weibull distribution [72], such as Pearson’s method of moments, method of quantiles, general curve fitting, Bayesian approach, maximum likelihood, etc. These approaches are relatively easy to apply for bimodal Weibull distributions but become complicated for trimodal Weibull distributions. 4.2.3 Extended Strength Distribution Model As discussed in section 4.1, SMC exhibits higher strength in flexure tests than in uniaxial ten- sion tests. Therefore, a model based on only tensile strength distribution would always underpredict flexure strength. Similarly, if the model is based on only flexural strength distribution, it would al- ways overpredict the tensile strength. Tensile strength distribution and flexure strength distribution can’t individually represent the overall properties of SMC material. To predict material behavior in different mechanical tests, it is necessary to combine tensile test results with 3-pt flexure test results and 4-pt flexure test results in the material model. For this purpose, we hypothesize that the flexural strength is equivalent to the tensile strength measured with small tensile specimens, as depicted in Figure 4.3. It is proposed to use the flexural strength data to expand the tensile strength spectrum. This may be done by adding the flexural strength data to the tensile strength data and fitting the combined data set with a unimodal Weibull model as in [123]. In the current work, the model was built as a trimodal Weibull distribution by combining the tension, 4-pt flexure, and 3-pt flexure Weibull 52 Moment Length Figure 4.3 In flexural test of brittle materials, the failure starts from the tensile side. The volume under higher tensile stresses is relatively small, which is equivalent to tensile test with a small size specimen. models with equal mixture weights, as described by Equation 4.10. 1 𝑝 = 𝑝1 = 𝑝2 = 𝑝3 = (4.9)     𝑚1    3 𝑚 2    𝑚  1 𝑋𝑖 𝑋𝑖 𝑋𝑖 3 𝑆 𝐸 𝑆𝐷 𝑀 (𝑖) = exp − + exp − + exp − (4.10) 3 𝑋01 𝑋02 𝑋03 The model in Equation 4.10 is called the Extended Strength Distribution Model (ESDM). Figure 4.4 compares the combined data set with the model predictions using the ESDM and the unimodal Weibull distribution. It can be seen that the ESDM matches the experimental data better in com- parison to the unimodal approximation. 4.3 Considering Randomness in PFE Models For an in-plane isotropic composite material modeled using shell elements, the inputs for ma- terial models include the in-plane tensile, compression, and shear properties. In the current work, only the tensile strength of the SMC was considered as the random variable, while other material properties were kept constant. A randomization algorithm was developed to assign the random variable in the PFE models. This included two steps: (1) discretization and (2) randomization. 53 1 0.9 Experiment Probability of survival (S) 0.8 Unimodal 0.7 ESDM 0.6 0.5 0.4 0.3 0.2 0.1 0 50 100 150 200 250 300 350 Ultimate Stress (MPa) Figure 4.4 The ESDM (red) was created by combining Weibull models of tension, 4-pt bending and 3-pt bending. 4.3.1 Discretization To obtain discretized tensile strength values, the Weibull model was divided into segments and the probability for each segment 𝑃𝑖 was determined by: 𝑃𝑖 = 𝑆 𝐸 𝑆𝐷 𝑀 (𝑖) − 𝑆 𝐸 𝑆𝐷 𝑀 (𝑖+1) (4.11) This is illustrated in Figure 4.5. In this demonstration, the probability Pi distribution for the ESDM was discretized into ten segments for strength values between 104 MPa and 346 MPa. Values outside this range were not considered, as no experimental data fell outside this range. This range covers 95% of ESDM’s spectrum. The tensile strength (𝑋𝑇 (𝑖) ) for a given segment is the mean strength of the segment and is calculated using Equation 4.12. The discretized probability (𝑃𝑖 ) is shown in Figure 4.5 as a function of the tensile strength (𝑋𝑇 (𝑖) ). 1 𝑋𝑇 (𝑖) = (𝑋𝑖 + 𝑋𝑖+1 ) (4.12) 2 54 0.2 0.16 Probability (P) 0.12 0.08 0.04 0 116 140 164 188.5 213 237 261 285.5 310 334 Tensile Strength [MPa] Figure 4.5 The discretized probability distribution of the tensile strength for ESDM with 𝑛𝑚𝑎𝑡 = 10. This represents probabilistic materials models. Each segment can be considered as a unique material and hence is represented with a unique material input set in simulations. The number of material input sets, i.e., the number of segments denoted as 𝑛𝑚𝑎𝑡 . In the current work, the material input sets differ only in tensile strength. The material input sets were assigned randomly to the elements in the PFE model following the ESDM’s probability distribution (P). This was implemented by controlling the frequency of the assignment. In SMC, fiber orientations and distributions vary both in-plane and through-thickness. For FE models with shell elements, the through-thickness integration points can be utilized to assign randomness in the through-thickness direction. In this case, the frequency 𝑓 of each material input sets to be assigned is calculated as: 𝑁 𝑃𝑖 𝐼 𝑓𝑖 = Í (4.13) 𝑃 where N is the total number of shell elements, I is the number of integration points through the thickness of the shell element, and Pi is the probability of the segment determined by Equation 55 4.12. Excluding the values generated by the ESDM outside the range of 104-346 MPa, we would Í have 𝑃 < 1. This effect is adjusted in Equation 4.13. 4.3.2 Approximation of Element Size Weibull size effect Equation 4.4 provides a correlation between strength and volume. This equation can be used to require element size that can acquire flexure strength for a small test volume from tensile strength. Later, volume can be used to approximate element size for a given specimen width. Therefore, at survival probability S=0.5 we can determine 𝑋𝑇 (at S=0.5) = 138.6 and 𝑋𝐹 (at S=0.5) = 269. Using these values and putting them in Equation 4.4, we get effective volume or element volume for the 3-point flexure test as:  𝑚 𝑋𝑇 𝑉2 = 𝑊 𝐿𝑇 = 𝑉1 = 7.906𝑚𝑚 3 (4.14) 𝑋𝐹 s   𝑚 𝑉1 𝑋𝑇 ⇒𝑊 =𝐿 = ≈ 1.6𝑚𝑚 (4.15) 𝑇 𝑋𝐹 where, 𝑊 =Wdith of element for effective volume 𝑇 =Height/thickness of SMC material (in our case which is 3.25mm) 𝐿 =length of element (in FEM analysis element is assumed to be square meaning W=L) 𝑉1 =Tensile test specimen volume 𝑚 =Weibull modulus for uniaxial tension (given in Table 4.2) 4.3.3 Randomization in FE Model The randomization process is illustrated using the FE model of a 3-pt bending specimen ran- domized with 𝑛𝑚𝑎𝑡 = 10, as shown in Figure 4.6a. The specimen was modeled with 320 shell elements. Each shell element had six integration points through the thickness. Therefore, the pop- ulation (total integration points) was 1920. A Python program was written and used to randomly as- sign these integration points to 10 pins according to the frequency determined by Equation 4.13 and then associate each pin to one material input set. In this way, a probabilistic FE model with a ran- 56 (a) (b) (c) Figure 4.6 Randomization of a (a) 3-pt bending specimen (b) Uniaxial tension specimen and, (c) 4-pt bending specimen with 𝑛𝑚𝑎𝑡 = 10. Elements with different colors have different material input sets. The element size is 1.5 mm in all cases. domly distributed property following a given statistical distribution can be generated. This method, in essence, is similar to Monte Carlo (MC) or quasi-MC method [22, 72, 122, 124]. However, in- stead of considering numerous material variables and using a refined discretization, the current method considers only one critical material parameter and uses a relatively coarse discretization. As to be presented in the next section, this approach can achieve sufficient accuracy with a limited number of iterations. Similar approach was taken for modeling specimens for uniaxial tension and 4-pt bending, as shown in Figure 4.6b and Figure 4.6c. 57 300 XT=115MPa XT=220MPa 250 XT=310MPa 200 Stress [MPa] 150 100 50 0 0.00 0.01 0.02 0.03 0.04 0.05 0.06 Strain [mm/mm] Figure 4.7 Tensile responses by one element simulations with the ECDM model with three different tensile strength values. 4.3.4 Material Model The material model used in this investigation was MAT 297, an Enhanced Continuum Damage Mechanics (ECDM) model for composite materials [93] implemented for shell elements in LS- DYNA. Figure 4.7 presents the tensile stress-strain responses produced by simulations using one element FE model with the ECDM model for the SMC. Three simulation runs were performed with three different tensile strength values while all other parameters were kept the same. The ECDM uses the Hashin failure criterion, the same as that in MAT 58 of LS-DYNA [125], which is a widely used composite damage model. The main features of the ECDM are: (1) it employs two different damage evolution laws in the pre-and post-failure regions, and hence it can describe the stress-strain response in both regions more accurately than MAT 58 [126], and (2) it is a coupled damage-plasticity model and hence can capture the energy absorption related to irreversible deformation in crashworthiness simulations [127, 128]. It should be mentioned that 58 Table 4.3 Three types of simulations. Symbol Type 𝑛𝑚𝑎𝑡 𝑋𝑇 (MPa) 𝑛𝑖𝑡𝑒𝑟 DET Deterministic 1 135.6 1 REG PFE 10 114-170 10 ESDM PFE 10 116-334 10 although MAT 297 was used, the method presented here can be applied to other material models and other finite element codes. 4.3.5 PFE Models All simulations were performed with LS-DYNA explicit solver. In explicit FE analysis, the time step is controlled by the size of the smallest element, and therefore, a uniform element size is preferred. The PFE model for the 3-pt bending simulation (Figure 4.6a) had a uniform element size of 1.5 mm. This size is tolerable in the current crashworthiness FE simulations. The effect of element size will be examined in the next section. 4.4 Results and Discussions 4.4.1 Simulations With a Randomly Distributed Tensile Strength Simulations were performed for uniaxial tension, 3-pt, and 4-pt bending cases. For each case, three types of simulations were performed: one was deterministic (DET) with the mean tensile strength value of 135.6MPa, the other two were probabilistic with a tensile strength distribution determined by either the regular Weibull model (REG) or the ESDM. PFE models were random- ized with 𝑛𝑚𝑎𝑡 = 10, following the procedure discussed in the previous section. Simulations were performed in 10 independent iterations, each with new randomization. As to be discussed in the Number of Randomization Iterations, the simulations with the current approach can converge within ten iterations. Table 4.3 summarises three types of simulations. An element size of 1.5 mm was used in all FE models. Figure 4.8, 4.9, and 4.10 compares the responses from three types of simulations with the exper- imental results. For the probabilistic simulations, the response was the average of ten independent simulations. Table 4.4 compares the mean experimental and predicted values for the peak force (or stress) and the area under the load-displacement (or stress-strain) curve. 59 EXP Results 0.15 REG ESDM DET Stress (GPa) 0.1 0.05 0 0 0.005 0.01 0.015 0.02 Strain (mm/mm) Figure 4.8 Comparison of experimental results with three types of simulations. Uniaxial tensile stress-strain curves. Simulation responses were the average of ten iterations. For uniaxial tension, the results were compared in stress-strain responses. In the tensile ex- periments, the strain was measured using a laser extensometer. While in the simulation, a virtual extensometer of the same gage length was defined on each specimen to determine the strain. Fig- ure 4.8 shows that DET and REG predicted peak stress within the scatter of the experimental peak stresses but a failure strain lower than that of the experimental values. The simulation by ESDM agreed best in terms of failure strain values. The predicted mean peak stress was relatively high in Table 4.4 Comparison in Ratios of Model Prediction Average to Experimental Average Values. Ratios to Experiment Average Experiment Average DET REG ESDM Tensile strength 136.5 MPa 1.003 0.990 1.103 Area under stress-strain 1.260 MPa.mm/mm 0.804 0.875 1.175 3-pt bending peak force 423 N 0.787 0.898 1.045 Area under force-disp 2040 N.mm 0.229 0.577 0.923 4-pt bending peak force 342 N 0.906 1.019 1.026 Area under force-disp 3445 N.mm 0.703 0.589 0.943 60 0.6 Exp Results DET 0.5 REG ESDM 0.4 Force (kN) 0.3 0.2 0.1 0 0 2 4 6 8 10 12 Displacement (mm) Figure 4.9 Comparison of experimental results with three types of simulations 3-pt bending force- displacement responses. Simulation responses were the average of ten iterations. comparison to the mean experimental value but still within the scatter. The area under the stress- strain curve was calculated. The mean experimental value and the ratio of the predicted values to the experimental values are presented in Table 4.4 for the tensile strength and the area. As shown, DET and REG predicted the mean peak stress within 1% but underpredicted the area by about 12 to 20%. On the other hand, ESDM overpredicted the peak stress by 10% and the area by 17%. The difference between the three simulations is more drastic for the 3-pt bending case. In Figure 4.9, DET predicted a peak force about 21% lower than the average peak force and an abrupt failure without any post-peak response. REG yielded a better prediction of the peak force and post-peak response. ESDM provided the best prediction. The predicted peak force was within the experimen- tal scatter, and the post-peak response was similar to the experimental curves. As shown in Table 4.4, in terms of the area under the force-deflection curve, the values predicted by DET, REG, and ESDM were about 23%, 58%, and 92% of the mean experimental value, respectively. It should be 61 0.5 EXP Results REG 0.4 DET ESDM 0.3 Force (kN) 0.2 0.1 0 0 4 8 12 16 20 Displacement (mm) Figure 4.10 Comparison of experimental results with three types of simulations for 4-pt bending force-displacement responses. Simulation responses were the average of ten iterations. noted that the area under the force-deflection curve represents the energy absorbed by the structure. DET and REG significantly underpredicted this amount. In crashworthiness simulations, these two methods would result in erroneous predictions for the crash performance of the composite structures [127, 128]. In the 4-pt bending case (shown in 4.10), the force-displacement curves predicted by three types of simulations all displayed a gradual failure pattern. However, DET and REG Weibull underpre- dicted the area under the force-displacement curve by about 30 to 40%. In summary, the predictions follow the order DET1 𝑛𝑖𝑡𝑒𝑟 1. 2. 3. 4. Generate Run FE 5. Develop Random- random Model as CV of SF with surrogate search for Hyperparameter objective OF function Errors function Repeat CV⩽1 Probabilistic FE Optimum hyperparameter Machine learning Figure 5.5 Flowchart for PFE model optimization with Machine learing. The goal of the optimization is to find a set of optimum hyperparameter values that would result in the lowest error between the average of the experiments and probabilistic FE simulations for the 3-pt bending experiment. Figure 5.5 presents an overview of this process. There are five steps: (1) generating random hyperparameters; (2) performing probabilistic simulation and obtaining the objective function (OF) for the run; (3) establishing the surrogate functions (SFs); (4) grid-searching for hyperparameters; and (5) cross-validation. Steps 1-2 are the same as the regular probabilistic simulations. These two steps will repeat until a predetermined number of probabilistic simulations 79 Figure 5.6 Two sets of probabilistic simulations with similar MAE values. have been completed. Steps 3-5 will use the ML process to analyze the results. The details of Step 1 can be found in Section 3. Steps 2-5 will be explained below. 5.3.1 Objective Function The Objective Function (OF) is defined as an error function that will be minimized. A common way to determine the error between experiment and simulation is to use the mean absolute error (MAE). MAE is a good measure for comparison between models, but it can be very ambiguous when optimizing hyperparameters. Being an absolute value, MAE does not provide information about whether the model is overpredicting or underpredicting, and it becomes very difficult to con- verge MAE to zero. For example, in Figure 5.6, the two simulations had similar MAE values, but the results were very different. Instead of MAE, the simulation and experimental results may be compared in multiple cate- gories. For example, the 3-pt bending force-displacement curve can be divided into several seg- ments, and the error can be evaluated in each segment. In the current work, the OF was evaluated for the following three categories: mean peak force error (MPFE), mean post-failure force error (MPFFE), and the mean residual force error (MRFE), i.e., OF=[MPFE, MPFFE, MRFE]. Figure 5.7 illustrates the mean simulation response of the 3-pt bending test with a random set of hyper- parameters and the location and ranges where these three errors are evaluated. Figure 5.8a-5.8c presents the distributions of MPFE, MPFFE, and MRFE obtained by probabilistic simulations. 80 mean simu 500 Peak Force mean exp 400 Post Failure Force Force [N] 300 200 100 Residual Force 0 0 2 4 6 8 10 12 14 Displacement [mm] Figure 5.7 The three types of errors and their locations with various 𝐻 = [ 𝑝, 𝜖 𝑚𝑎𝑥 , 𝑅𝐸 𝑆1𝑇]. 5.3.2 Surrogate Function The large amount of information generated in probabilistic simulations is analyzed with the aid of a surrogate function (SF). In the current work, SF is established using artificial neural network (ANN), one of the machine learning algorithms. A simple ANN consists of inputs, a single hidden layer, and an output. The hidden layer works like a perceptron which receives inputs, multiplies them by weight, and then passes them into an activation function to produce an output. Further details about ANN are provided in section 2.6 on chapter 2 at page 18. In the current work, ANN was developed using Tensorflow and Keras library in Python. ANN with three hidden layers and a unit activation function was used. The three hyperparameters [ 𝑝, 𝜖 𝑚𝑎𝑥 , 𝑅𝐸 𝑆1𝑇] are the input variables. The three errors [PFE, MPPFE, MRFE] are the outputs. As ANN can only have one output, the analysis was performed with three inputs for each output, and three SFs were established. 81 (a) (b) (c) Figure 5.8 The three types of errors. (a) MPFE (b) MPFFE and (c) MRFE obtained from proba- bilistic simulations with various 𝐻 = [ 𝑝, 𝜖 𝑚𝑎𝑥 , 𝑅𝐸 𝑆1𝑇]. 82 Table 5.4 Example of Inputs and Outputs with Surrogate Function. 𝑝 𝜖 𝑚𝑎𝑥 𝑅𝐸 𝑆1𝑇 𝑀 𝑃𝐹𝐸 𝑀 𝑅𝐹𝐸 𝑀 𝑃𝐹𝐹𝐸 𝑅𝑀𝑆𝐸 0.29 0.28 0.38 4.9 0.14 3.78 3.57 0.3 0.28 0.38 -0.69 -2.97 -5.48 3.62 0.27 0.29 0.34 6.36 -3.21 -4.52 4.87 0.32 0.29 0.39 -4.11 3.82 -7.75 5.52 0.32 0.28 0.4 -6.58 2.82 -7.94 6.17 0.27 0.31 0.35 7.43 7.96 -2.51 6.45 0.29 0.26 0.4 6.35 -2.9 9.46 6.79 5.3.3 Estimation of Optimum Parameter With SFs, the errors of simulations to experimental results corresponding to any set of hyper- parameter values can be evaluated more quickly. In Step 4, a search for a set of errors is performed with a random set of hyperparameters at a large sample volume covering its entire range. Further- more, the root mean squared error (RMSE) can be evaluated by these errors as follows: r 𝑀 𝑃𝐹𝐸 2 + 𝑀 𝑃𝐹𝐹𝐸 2 + 𝑀 𝑅𝐹𝐸 2 𝑅𝑀𝑆𝐸 = (5.13) 3 The grid search will determine the optimum set of hyperparameter values corresponding to a min- imum RMSE and a 5% error range around the optimum values. Table 5.4 presents several sets of hyperparameters with relatively low RMSE computed by a random search for error with the SFs. 5.3.4 Cross Validation Based on the results of Grid Search, a second cycle optimization will be performed. Starting at Step 1, the hyperparameters will be generated randomly within the bound of the 5% range of the optimum established in Grid Search. In Step 2, probabilistic simulations will be run, and OF will be obtained. In Step 3, new SFs will be established. In Step 4, a new optimum set of hyperparameter values will be determined. This process will end when the RMSE meets the targeted value. 5.4 Results and Discussions Following the ML optimization algorithm discussed above, the optimum set of hyperparameters was determined as 𝐻 = [ 𝑝, 𝜖 𝑚𝑎𝑥 , 𝑅𝐸 𝑆1𝑇] = [0.28, 0.2, 0.46]. This set of parameters was examined in probabilistic simulations for tension, compression, 3-pt bending, and 4-pt bending cases. For each case, the simulations were performed with the optimum 𝐻 values for 20 iterations following 83 160 mean simu mean simu 140 simu 95% confidence simu 95% confidence mean exp 200 mean exp 120 exp 95% confidence exp 95% confidence 100 150 stress [MPa] stress [MPa] 80 60 100 40 50 20 0 0 0.000 0.005 0.010 0.015 0.020 0.025 0.000 0.005 0.010 0.015 0.020 0.025 0.030 strain [mm/mm] strain [mm/mm] (a) (b) Figure 5.9 (a) Uniaxial tension and (b) Uniaxial compression. The mean tress-strain curves and their scatter at 95% confidence level obtained by experiment, and by probabilistic simulations with H = [0.28, 0.2, 0.46] with 20 randomization iterations. the randomization procedure in Chapter 4 on Section 4.3.3 at page 56. Figure 5.9a presents the case of uniaxial tension. The FE model was developed with 𝑛𝑚𝑎𝑡 = 40 randomization. The stress-strain curves obtained in experiments and probabilistic FE simulations are compared. The results are shown by the mean curve and their scatters corresponding to a 95% confidence level. Simulation results agreed well with the experimental stress-strain curve. Figure 5.9b presents the results of uniaxial compression. The same as in Figure 5.9a, the stress- strain mean curves and their scatters corresponding to 95% confidence level obtained by experiment and simulations are compared. The experimental curve shows a sudden failure, whereas the curve generated by simulations exhibits a more gradual failure process. The simulations were performed with a residual compressive strength 𝑅𝐸 𝑆1𝐶 = 0.4. This value is a legacy value used in crashwor- thiness simulations of braided composites [127]. It seems to be too high for the SMC material. Figure 5.10a, and 5.10b compare the force-displacement curves obtained by experiment average and by probabilistic simulations with 𝐻 = [0.28, 0.2, 0.46] for 3-pt bending and 4-pt bending, respectively. For 3-pt bending, it is not surprising that the simulations agreed very well with the average experimental curve since the 𝐻 values were determined in this case. For 4-pt bending, the average curve generated by simulation agreed well with the experimental average up to the peak load but was higher in the post-peak region. As for the case of compression, the high 𝑅𝐸 𝑆1𝐶 value 84 (a) (b) Figure 5.10 The mean force-displacement curves and their scatter at 95% confidence level obtained by experiment average and by probabilistic simulations with H = [0.28, 0.2, 0.46] with 20 random- ization iterations. (a) 3-pt bending. (b) 4-pt bending. 85 could be the cause of more gradual post-peak responses in simulations. In future work, 𝑅𝐸 𝑆1𝐶 will be added to the set of hyperparameters in model parameter optimization. 5.5 Summary This work investigated the probabilistic simulations of SMC structures. The scatter of the ten- sile and compression strengths and moduli were considered. To account for the size effect on the strength, a General Strength Distribution Model (GSDM) was employed. To determine the mixture weight fraction of the bimodal Weibull models and some difficult-to-measure parameters in the damage mechanics-based composite material model, an optimization procedure for probabilistic FE models using an artificial neural network (ANN)-based machine learning (ML) algorithm has been developed. This procedure allows one to obtain a set of optimum parameters that can produce satisfactory simulation results for different load cases. In the current work, the optimization was performed using a 3-pt bending case for three hyperparameters, and the results were examined for uniaxial tension, uniaxial compression, and 3-pt and 4-pt bending cases. Simulations with this ap- proach yielded low errors for all 4 cases, and the results reached a steady state within 20 iterations. 86 CHAPTER 6 OPTIMIZATION OF PFE MODEL WITH RANDOM SEARCH This chapter presents a Random Search (RS) based optimization technique to estimate the hyper- parameters for Probabilistic Finite Element (PFE) model that was developed in Chapter 5. The optimum set of parameters was subsequently cross-examined in probabilistic simulations of ten- sion, compression, and 4-pt bending experiments. Later, this model was validated with open hole tension and disk bending with a fixed boundary and with simply supported boundary conditions. 6.1 Model Description Probabilistic FE model used in this chapter is identical to the model described in Chapter 5, Section 5.2 at page 70. Briefly, the model used in this chapter is based on General Strength Distri- bution Model (GSDM) as described by Equation 5.10 at page 74. This GSDM was implemented in MAT 297 with 𝑛𝑚𝑎𝑡 = 40 to perform PFE simulations which is described in Section 5.2.4.3 at page 77. 6.2 Model Optimization A few parameters in the GSDM and MAT 297 are difficult to assess. These parameters can be estimated by optimization techniques. In the model optimization process, they are called hyperpa- rameters. The goal of model optimization is to search for an optimal set of hyperparameters that leads to the least error in the objective function. Commonly used methods for model optimization include Bayesian [130], Random Search (RS) [86, 131], Multivariate Regression [129], MC [132, 133], etc. Model optimization is an iterative process and tends to be computationally expensive. In our previous attempt, machine learning- based multivariate regression was used for model optimization [129]. Although it was compu- tationally inexpensive, it was deemed to be rather complicated for a general-purpose regression. Additionally, machine learning has its own hyperparameters, which can utterly affect the optimiza- tion. MC is another commonly used technique. It is effective when the number of random variables is low. However, the results may converge to a local minima instead of global minima, as all itera- 87 Start 𝑛𝑖𝑡𝑒𝑟 Generate random Run FE Model or argmin(MAE) Hyperparameter objective function Repeat Optimum hyperparameter Figure 6.1 Graphical representation of probabilistic FE model optimization algorithm using Ran- dom Search. tions are dependent [88, 134]. In this work, RS is used. RS is the simplest approach and extremely easy to implement for general purposes. Additionally, it reliably converges to the global minima as all iterations are independent. Although RS is computationally expensive, this cost is acceptable as the optimization has to be done only once for any specific material model. 6.2.1 Optimization Flow-Chart Figure 6.1 presents the flowchart of an RS-based model optimization procedure. There are three steps: (i) generating random hyperparameters; (ii) performing probabilistic FE simulation and obtaining the objective function (OF) for the run; (iii) finding the least error and corresponding hyperparameters. These steps are explained below: i. Generate Random Hyperparameters: Hyperparameters (𝐻) are defined as arbitrary pa- rameters in the model. These parameters are typically determined through an iterative process because they cannot be determined from experimental results. In this work, we have three hyperparameters. One is the mixture weight fraction (𝑝) in GSDM (Equation 5.10). The 88 other two are MAT 297 parameters: the strain at element deletion (𝜖 𝑚𝑎𝑥 ) and the residual tensile strength (𝑅𝐸 𝑆1𝑇). Therefore, 𝐻 can be described as: 𝐻 = [ 𝑝, 𝜖 𝑚𝑎𝑥 , 𝑅𝐸 𝑆1𝑇] (6.1) These three hyperparameters are random independent variables ranging from 0-1 with two decimal precisions. A pseudo-random number generator was used for this purpose. ii. Objective Function (OF): Objective Function (OF) is an error function that will be mini- mized through optimization. It performs probabilistic FE simulation to determine the error between simulation and experiment. In model optimization, OF is often chosen as the dif- ference between the responses for a particular load case generated by FE simulations and by experiment. In this work, OF is defined as a Mean Absolute Error (𝑀 𝐴𝐸) determined by comparing the mean force-displacement curve from probabilistic simulations and exper- iments. For each set of H, a predetermined number of probabilistic simulation runs (𝑛𝑖𝑡𝑒𝑟 ) are performed to obtain MAE. It then goes back to step (𝑖) for the next set of 𝐻. This pro- cess repeats until the entire spectrum for H is covered, and the MAE as a function of H is established 𝑀 𝐴𝐸 = 𝑂𝐹 (𝐻) (6.2) iii. Estimate Optimum Hyperparameter: This step determines the lowest value of the OF and the corresponding H, the optimum hyperparameter 𝐻𝑜 𝑝𝑡 . 𝐻𝑜 𝑝𝑡 = 𝑎𝑟𝑔𝑚𝑖𝑛(𝑀 𝐴𝐸) (6.3) The optimization procedure is driven by Python and shell script. 6.3 Optimization strategy The model optimization may be performed with a selected load case or multiple load cases. In this work, three strategies were examined. These strategies use the same optimization procedure as depicted in Figure 6.1 but with different load cases or load case combinations. They are labeled as Opt1-3. The model developed in Chapter 4 is served as the baseline, labeled as Base. Table 6.1 provides a summary of these model optimization strategies. 89 Table 6.1 Summary of model optimization plans and optimized hyperparameters (𝐻𝑜 𝑝𝑡 ). ID Load Case(s) in OF 𝐻𝑜 𝑝𝑡 𝑀 𝐴𝐸 𝑝 𝜖 𝑚𝑎𝑥 𝑅𝐸 𝑆1𝑇 Base Baseline, no optimization 0.33 0.4 0.1 - Opt1 3-pt bending 0.28 0.2 0.46 9.71 Opt2 4-pt bending 0.31 0.37 0.19 10.81 Opt3 3-pt 4-pt bending 0.29 0.23 0.37 13.90 Figure 6.2 Trimodal Weibull Strength distribution used in Base optimization strategy. 6.3.1 Base It is identical to the model described in Chapter 4 but with 𝑛𝑚𝑎𝑡 = 40. This model considered the randomness in tensile strength by a trimodal Weibull distribution of tensile, 3-pt, and 4-pt strength as shown in Figure 6.2. The contribution of the strengths from these three cases was assumed to be equal, i.e., 𝑝=0.33. The model parameters 𝜖 𝑚𝑎𝑥 , 𝑅𝐸 𝑆1𝑇 were assumed based on previous experience without optimization. Furthermore, the compressive strength was kept constant at 300 MPa so that the compressive failure was suppressed. Nonetheless, this model performed fairly well for standardized tests, including uniaxial tension and 3-pt and 4-pt bending, as shown in Chapter 4. 90 (a) (b) Figure 6.3 Opt1 optimization strategy (a) 3-pt bending force-displacement responses with various H (b) The mean absolute error (MAE) versus random H. 91 6.3.2 Opt1 Model optimization was performed upon a 3-pt bending load case. The optimization follows the procedure described in section 6.2.1. It takes a random set of 𝐻 and runs 3-pt bending proba- bilistic simulations for 𝑛𝑖𝑡𝑒𝑟 = 20. The mean force-displacement curve from these 20 iterations is determined, and the 𝑀 𝐴𝐸 for this 𝐻 set is computed. The process repeats for the next random set of 𝐻. This process is repeated 300 times to obtain sufficient resolution for the spectrum of H. Figure 6.3a shows the mean force-displacement curves of 3-pt bending simulations with dif- ferent 𝐻. Each dotted line represents a mean curve generated from 20 probabilistic simulations with a specific set of H, and the solid black line represents the mean experimental response. As shown, some mean curves significantly underpredicted/overpredicted the entire force-displacement response. A few simulations predicted the peak load properly but not at the post-failure region. Figure 6.3b demonstrates the MAE as a function of H. It should be noted that 𝑝, 𝜖 𝑚𝑎𝑥 , 𝑅𝐸 𝑆1𝑇 are varied simultaneously, i.e., it would be a 4D plot but is represented in a 2D plot. The relation between H and MAE is not obvious from this 2D plot, and the optimum parameters cannot be estimated graphically. Therefore, 𝐻𝑜 𝑝𝑡 was determined from the array of hyperparameters and their corresponding MAE using Equation 6.2. 𝐻𝑜 𝑝𝑡 from Opt1 optimization are shown in Table 6.1. 𝐻𝑜 𝑝𝑡 from Opt1 optimization plan can be used to plot the GSDM described by Equation 5.10 and 5.11 at page 74. As mentioned in section 5.2.2.3, GSDM is a 3D vector of tensile strength, compressive strength, and probability. To make visualization easier, the strength vector is repre- sented by a 2D plot with two 𝑌 axes, as shown in Figure 6.4. Each discretized strength in this plot has two 𝑌 values, i.e., 𝑌 1 = compressive strength, which is denoted by dots, and 𝑌 2 = probability, which is denoted by bars. 6.3.3 Opt2 The optimization was performed upon 4-pt bending, i.e., 𝐻𝑜 𝑝𝑡 was determined in such a way that would lead to the lowest mean absolute error in 4-pt bending simulation. Hence, the error is denoted as 𝑀 𝐴𝐸 4𝑝𝑡 . 𝐻𝑜 𝑝𝑡 from Opt2 are shown in Table 6.1. Strength distribution from Opt2 is 92 Figure 6.4 Bimodal Weibull distributions with 𝑝 = 0.28 obtained from model optimization Opt1. shown in Figure 6.5a. 6.3.4 Opt3 This was a combination of Opt1 and Opt2. 𝐻𝑜 𝑝𝑡 was determined in such a way that would lead to the least mean absolute error in 3-pt bending as well as 4-pt bending simulation. 3-pt and 4-pt bending simulation was run 𝑛𝑖𝑡𝑒𝑟 = 20 times for the same random set of Hyperparameter to determine error as defined by 6.4. This process was repeated 300 times as described in section 6.2.1. 1 𝑀 𝐴𝐸 = (𝑀 𝐴𝐸 3𝑝𝑡 + 𝑀 𝐴𝐸 4𝑝𝑡 ) (6.4) 2 6.4 Iterations and Computation Time Probabilistic FE simulations are independent iterations. As a result, iterations can be performed simultaneously. The computation time almost linearly reduces with the number of parallelization. In other words, the number of available CPU cores and threads would be the only limiting factor regarding the computation time [135–137]. Model optimization can be computationally expensive as the number of iterations is usually 93 (a) (b) Figure 6.5 Bimodal Weibull distributions with p obtained by model optimization: (a) Opt2, 𝑝 = 0.31. (b) Opt3, 𝑝 = 0.29. 94 Table 6.2 Summary of computation time requirement for different analysis. Total number Number of Time for each Total time Analysis type of iterations threads iteration [min] [min] Opt1 and Opt2 6000 100 20-40 1200-2400 Optimization Opt3 12000 100 20-40 2400-4800 Ten, Comp 20 20 10-15 10-15 Probabilistic 3-pt, 4-pt 20 20 20-40 20-40 Simulation Hole, Disk1, Disk2 20 20 40-120 40-120 high. In the current study, an optimization usually takes 6000-12000 iterations. On the other hand, the subsequent probabilistic FE simulations with an optimized model do not take significant compu- tation time as it requires rather limited iterations. In the current study, all probabilistic simulations were performed with 20 iterations. Table 6.2 is a summary of the computation time in various analyses in this work. All computa- tions were performed using a High-Performance Computing (HPC) system with 112 Virtual CPUs or threads (2 x Intel® Xeon® Platinum 8180, 56 cores) and 768GB of RAM. The load cases listed in the last row will be discussed in section 6.5. 6.5 Results and Discussions The optimal sets of 𝐻 obtained from the three optimized models were examined with proba- bilistic simulations for seven different load cases. These load cases are summarized in Table 6.3. The first four are standard mechanical tests, including uniaxial tension, compression, 3-pt bending, and 4-pt bending. They serve the purpose of model verification. The next three tests are specially designed validation tests representing several common structural loading scenarios. They serve the purpose of model validation. For comparison, probabilistic simulations were also performed with the baseline model. For each test, probabilistic FE models with corresponding 𝐻𝑜 𝑝𝑡 values were randomized with 𝑛𝑚𝑎𝑡 = 40 following the procedure discussed in chapter 4 on Section 4.3.3 at page 56. This process was repeated 20 times, i.e., 𝑛𝑖𝑡𝑒𝑟 = 20. Each iteration is run with a new randomization. Finally, the responses from 20 iterations were averaged into the arithmetic mean. An element size of 1.5 mm was used in all FE models. 95 Table 6.3 Seven test cases examined in probabilistic simulations with optimized parameters. ID Description Standard Ten Uniaxial tension ASTM D3039 [107] Comp Uniaxial compression In-house 3-pt 3-point bending ASTM D790-03 [108] 4-pt 4-point bending ASTM D7264 [109] Hole Open hole tension ASTM D5766 [138] Disk1 Disk bending with fixed boundary In-house Disk2 Disk bending with simply supported boundary In-house 6.5.1 Model Verification Figure 6.6 through 6.8 compares the mean simulation responses with the experimental responses for the four verification tests. The mean experimental response is shown in a solid black line. The area in gray represents the 95% confidence level for the experimental results. To quantitatively measure how well a model performed for different tests, the error for different quantities such as the peak load/stress, the area under the curve, and the strain at peak load/stress have been calculated by the following equation: 𝑀𝑒𝑎𝑛𝑠𝑖𝑚𝑢𝑙𝑎𝑡𝑖𝑜𝑛𝑣𝑎𝑙𝑢𝑒 − 𝑀𝑒𝑎𝑛𝑒𝑥 𝑝𝑒𝑟𝑖𝑚𝑒𝑛𝑡𝑎𝑙𝑣𝑎𝑙𝑢𝑒 𝐸𝑟𝑟𝑜𝑟 = 𝑥100 (6.5) 𝑀𝑒𝑎𝑛𝑒𝑥 𝑝𝑒𝑟𝑖𝑚𝑒𝑛𝑡𝑎𝑙𝑣𝑎𝑙𝑢𝑒 The following errors were computed for the Base and Opt1-3 models: the Peak Force or Stress Error (PFSE), Strain/Displacement Error (SDE) at the peak force/stress, and the Area Under the Curve Error (AUCE). These values are summarized in Table 6.4. Table 6.4 Error (%) for Peak Force/Stress Error (PFSE), Strain/Displacement Error (SDE) at the peak force/stress, and the Area Under the Curve Error (AUCE). ID Base Opt1 Opt2 Opt3 PFSE SDE AUCE PFSE SDE AUCE PFSE SDE AUCE PFSE SDE AUCE Ten 16.7 1.2 5.4 7.6 0.74 9 4.2 7.2 8.3 1.4 11.5 5.9 Comp 31.2 71.3 276.8 7.8 6.4 31.6 15.4 7.8 132.6 6.9 6.4 50.4 3-pt 5.8 7.5 10.1 2.5 0.9 1.4 3.7 0.9 0.1 1.3 3.9 3.7 4-pt 2.6 13.3 5.8 6.9 2.2 20.9 1.6 4.3 1.9 7.1 4.3 6.5 Figure 6.6a presents the case of uniaxial tension. In the tensile experiments, the strain was measured using a laser extensometer. In simulations, a virtual extensometer of the same gage length was defined on each specimen to determine the strain. As shown, the prediction with Opt3 was most 96 160 mean expt 140 expt 95% 120 mean base mean opt1 100 mean opt2 Stress [MPa] 80 mean opt3 60 40 20 0 0.000 0.005 0.010 0.015 0.020 0.025 Strain [mm/mm] (a) mean expt expt 95% 200 mean base mean opt1 150 mean opt2 Stress [MPa] mean opt3 100 50 0 0.00 0.01 0.02 0.03 0.04 0.05 Strain [mm/mm] (b) Figure 6.6 The experimental mean response curves, their 95% confidence band, and the mean re- sponse curves generated by probabilistic simulations with Base, and Opt1-3 for (a) uniaxial tension, (b) uniaxial compression. 97 600 mean expt expt 95% 500 mean base 400 mean opt1 mean opt2 Force [N] 300 mean opt3 200 100 0 0 2 4 6 8 10 12 14 Displacement [mm] Figure 6.7 The experimental mean response curves, their 95% confidence band, and the mean re- sponse curves generated by probabilistic simulations with Base and Opt1-3 for 3-pt bending. accurate in terms of the peak stress and the strain at peak stress. Opt1 and Opt2 overpredicted the peak stress by 5% and 3%, respectively. The baseline model overpredicted the peak stress by 10%. In general, all simulation results agreed fairly well with the experimental stress-strain curve. This shows that the uniaxial tension is not sensitive to the variation of the three hyperparameters. Figure 6.6b presents the results for uniaxial compression. As seen, simulations varied dras- tically with the H values obtained with different optimization strategies. The experimental curve shows a sudden failure, whereas the curve generated by simulations exhibits a gradual failure pro- cess. As mentioned above, to suppress the compressive failure, the baseline model was assigned a compressive strength higher than the value measured in the experiment. This resulted in a com- puted peak stress 31% higher than the experiment value. Additionally, the baseline model was not optimized for element erosion criteria (𝜖 𝑚𝑎𝑥 ), which resulted in a 276% larger area under the curve. Opt2 overpredicted this area by 133%. This may be due to the relatively high 𝜖 𝑚𝑎𝑥 value of 0.37. Opt1 and Opt3 had a 𝜖 𝑚𝑎𝑥 value at 0.2 and 0.23, respectively. These resulted in a smaller error for the area under the curve. Overall, Opt1 and Opt3 predicted the compression failure better than the 98 mean expt 500 expt 95% mean base 400 mean opt1 mean opt2 Force [N] 300 mean opt3 200 100 0 0.0 2.5 5.0 7.5 10.0 12.5 15.0 17.5 Displacement [mm] Figure 6.8 The experimental mean response curves, their 95% confidence band, and the mean re- sponse curves generated by probabilistic simulations with Base and Opt1-3 for 4-pt bending. two other simulations. Figure 6.7 compares the force-displacement curves for 3-pt bending between the simulations and experiments. It is not surprising that all simulations agreed well with the experiment since both the trimodal and bimodal Weibull strength models considered the 3-pt flexural strength. Further- more, Opt1 and Opt3 were optimized with a 3-pt bending test and hence resulted in the lowest error for the area under the curve, i.e., 1.4% and 3.7%, respectively. The baseline model and Opt2 were not optimized for 3-pt bending and hence a slightly larger error. Figure 6.8 compare the force-displacement curves for 4-pt bending. As shown, all simulations except Opt1 matched the experiment mean force displacement. Opt1 overpredicted the area under the curve by 21% as it was only optimized for the 3-pt bending case. Opt2 and Opt3 had the least error (1.9% and 6.5% respectively) as they were optimized for or with a 4-pt bending case. The baseline model, though not optimized for any cases, resulted in only a 5.8% error for the area under the curve. Figure 6.9 compares the errors of Base and Opt1-3 models for the four verification cases in the 99 Base Opt1 Opt2 Opt3 300 Ten 3-pt Comp 4-pt 250 200 Error [%] 150 100 50 0 PFSE PFSE PFSE PFSE SDE SDE SDE SDE AUCE AUCE AUCE AUCE Figure 6.9 Comparison of errors of four different models for four verification load cases. form of cumulated error. For example, the first column is the Peak Force or Stress Error (PFSE) of the baseline model from all four cases, where the error for each case is represented by a specific color block. For Base, the errors were large for compression but relatively small for other load cases. This is due to the fact that the compression failure strength was set high in the baseline model. The errors were reduced significantly with optimized models. Furthermore, a model performed better for the load case it was optimized for but was less reliable for other cases. For example, Opt1 performed best for 3-pt bending but not so well for 4-pt bending. It was vice versa for Opt2. Opt3 exhibited moderate error in all cases. 6.5.2 Model Validation 6.5.2.1 Open Hole Tension (Hole) The open hole tensile strength is one of the limiting factors in composite structural design. The open hole strength is often used to examine the capability of failure prediction for composite 100 3.25 25 150 25 Laser Tag 𝜙12.5 ± 0.25 50 Metal tab Specimen Figure 6.10 Schematic diagram for open-hole tension test. structures. [139–142]. This load case is selected as one of the validation cases in the current work. Figure 6.10 shows the schematic of the open-hole specimen. The setup of the open-hole tension experiment follows the guideline of ASTM D5766 [138]. It is shown in Figure A.27 at page 148. The specimens had a nominal hole diameter of 12.5±0.25 mm and a width of 50 mm. Aluminum end tabs were placed on the two ends of the specimen to prevent damage under the grips. The experiment was performed with a servo-hydraulic testing frame under displacement-controlled mode. A laser extensometer with a gauge length of 25 mm was used to measure the displacement of the specimen. The laser tag positions are shown in Figure 6.10. Eight specimens were tested, and their mean force-displacement curve is shown by the black line in Figure 6.11. The grey area represents the 95% confidence band in experiment results. The results exhibited large variations from the mean response. This could be due to the scatter in the hole diameter, the imperfections near the hole, in addition to the randomness in properties. Probabilistic FE simulations were performed for the open-hole tension case. Simulations were run with FE models with a hole diameter of 12.5 mm, and the variation of the hole diameter was neglected. A virtual extensometer of the same gage length as in experiments was defined on each specimen to determine the displacement. Simulations were performed for 𝑛𝑖𝑡𝑒𝑟 = 20 with the Base and Opt1-3 models. Figure 6.11a compares the mean force-displacement curves with the experimental results. As shown, the baseline and Opt2 predicted peak stress and post-peak response within the 95% band but at the lower side. The prediction by Opt1 was at the higher end. Opt3 101 20.0 mean expt mean opt1 17.5 expt 95% mean opt2 15.0 mean det mean opt3 mean base 12.5 Force [kN] 10.0 7.5 5.0 2.5 0.0 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 Displacement [mm] (a) ε1 [%] DIC 1.50 1.500 1.350 DIC 1.20 1.200 OPT3 1.050 OPT3 0.90 0.900 0.750 0.60 0.600 DIC OPT3 0.450 0.30 0.300 0.150 0.00 0.000 (b) Figure 6.11 (a) The comparison of simulation results of different optimization plans with experi- ment results (b) The comparison of the mean response and its 95% confidence band of Opt3 with experiment results. The inserts show the major strain field from the experiment (DIC) and simula- tions. 102 [MPa] 150 120 Y Stress 90 60 30 0 1.0 0.8 Y Damage 0.6 0.4 0.2 0.0 Y Displacement=0.3 Y Displacement=0.6 Y Displacement=0.8 Y Displacement=1.2 Figure 6.12 Typical Y stress, Y damage parameter distributions at the front surface under tensile load for open hole specimen at different displacements, obtained with Opt3. yielded the best prediction. Figure 6.11b examines the simulations by Opt3. The red shaded area represents the 95% con- fidence band of the simulations. The scatter in simulations is significantly lower than the experi- mental scatter. This may be attributed to the constant hole diameter used in the FE model. In the open-hole tension experiment, Digital Image Correlation (DIC) was used to measure the strain field evolution. The inserts in Figure 6.11b compare the major strain fields obtained from DIC and Opt3 at displacement d=0.3 mm, 0.8 mm (the peak load), and 1.2 mm (near rupture). The major strain distributions in simulations resembled DIC measurement closely at these instants. Figure 6.12 presents the typical y-stress and y-damage parameter distributions at the front sur- face of the open hole specimen during simulations. The damage parameter ranges from 0-1, repre- senting no damage to total damage. As shown in Figure 6.12, at d=0.3 mm, there was no damage. At d=0.6 mm, a few elements in the damage field changed their color from dark blue to light blue, corresponding to damage =0.2. At d=0.8 mm (the peak load), some elements near the edge of the hole had reached damage =1, and the stress field had been redistributed. At d=1.2 mm, the ele- ments across the center line had severely deformed, and the stress reduced to nearly zero for the entire specimen, indicating the specimen had ruptured. 103 Force 𝜙75 𝜙50 𝜙50 3.25 Specimen 𝜙50 Specimen (a) Front View (b) Top View Figure 6.13 Schematic diagram for disk bending test with fixed boundary (Disk1). 6.5.2.2 Disk Bending with Fixed Boundary (Disk1) The lateral bending of a disk creates a biaxial stress field in the specimen. In this work, disk bending was tested under two boundary conditions: fixed and simply supported. For fixed boundary, the specimen was sandwiched between a pair of steel frames with a circular opening of 75 mm diameter. The specimen and the steel frames were secured together with eight bolts. The force was applied to the specimen by a hemispherical nose of 50 mm diameter made of aluminum. Figure 6.13 presents a schematic diagram of the experimental setup. In the side view, the grey area represents the specimen, the hatched lines represent the fixed boundary, and the white area represents the hemispherical nose for loading. The experiment was performed with an electromechanical MTS machine. The load and displacement data were measured from the load cell and machine displacement. Figure 6.14 presents the photos of the experiment. The side view shows the experimental setup with the specimen sandwiched between the fixture plates. The bottom-front view shows a failed specimen before being removed from the fixture. Four specimens were tested with this configura- tion. 104 Nose Load Cell Specimen Nose Fixture Front-Top View Specimen Fixture Fixture Specimen (Sandwiched in between fixtures) Nose Front View Front-Bottom View Figure 6.14 Experiment setup for disk bending with a fixed boundary (Disk1). Probabilistic FE simulations were performed with 𝑛𝑖𝑡𝑒𝑟 = 20 using Base and Opt1-3 models. The mean force-displacement curves from these simulations are compared with the experimental curve in Figure 6.15a. The difference between the baseline and optimized models was drastic in this load case. The baseline underpredicted the peak load. In the post-peak region, the predicted re- sponse was below the experimental 95% confidence band. The predictions of the optimized models were acceptable for the peak load and pos-peak region till the displacement reached 12 mm. The differences between Opt1-3 were relatively small. Opt1 and Opt3 matched the experimental peak load slightly better than Opt2. Figure 6.15b compares simulations with Opt3 and experiment results. The red shaded area represents the 95% confidence band of the simulations. With 𝑛𝑖𝑡𝑒𝑟 = 20, the scatter in simulations in the post-peak region was much greater than the experiment, particularly after 12 mm displacement. This error may be attributed to the limitation of the shell element. In the experiment, the area of 105 7 mean expt 6 expt 95% mean det 5 mean base mean opt1 Force [kN] 4 mean opt2 3 mean opt3 2 1 0 0 5 10 15 20 Displacement [mm] (a) ε1 [%] 7 mean expt mean opt3 5.00 expt 95% opt3 95% 6 4.10 5 Force [kN] 4 3.20 3 2.30 2 1 1.40 0 0.50 0 5 10 15 20 Displacement [mm] (b) Figure 6.15 Disk bending with a fixed boundary. The mean force-displacement curves were ob- tained by experiment and by probabilistic simulations. (a) Comparison of simulations with different models and experimental results. (b) Comparison of the mean and its 95% confidence level by Opt3 with experiment results. Inserts compare the failure modes. 106 Figure 6.16 Disk bending with fixed boundary. Typical major stress, damage distributions at the bottom surface of the disk at different displacements in a simulation with Opt3. d=6 and d=17 correspond to the peak load and near final rupture, respectively. the specimen under the loading nose broke apart, and yet the specimen still carried some load. This scenario is very difficult to simulate with shell elements. The inserts in Figure 6.15b compare the failure patterns in the simulation and in the experiment. The insert for simulation shows the major strain distribution. To reveal the failure, the elements with a strain exceeding 5% were shown in white. The failure patterns in simulations had some similarities to that observed in experiments. Figure 6.16 presents the typical major stress and damage distributions at the lower surface ob- tained with Opt3. At d=3 mm, the damage started to appear, which corresponds to a small bend in the slope of the load-displacement curve in Figure 13b. In the experiment, however, this bend happened earlier, roughly at d = 2 mm. At d=6 mm, the peak load, an array of elements failed at the disk center, forming a cross-liked pattern. The major stress field revealed that the stress in these elements had reduced to nearly zero. At d=11 mm, the damage spread to a larger area, and finally, at d=17 mm, most of the elements failed. 107 𝜙100 Force 𝜙80 ° 120 𝜙50 𝜙50 𝜙50 Specimen 3.25 𝜙5 𝜙5 Specimen 𝐹𝑜𝑟𝑐𝑒 𝐹𝑜𝑟𝑐𝑒 𝐹𝑜𝑟𝑐𝑒 3 3 3 (a) Front View (b) Top View Figure 6.17 Schematic diagram for disk bending with a simply supported boundary (Disk2). 6.5.2.3 Disk Bending with Simply Supported Boundary (Disk2) For disk bending with a simply supported boundary, a disk-shaped specimen of 100 mm diam- eter was supported at three points by supporting wedges with a 5 mm diameter round top. Three wedges were placed 120º angle apart and 40 mm away from the center of the specimen. The force was applied to the specimen using the same hemispherical nose with the same testing frame as in the fixed boundary experiment. Figure 6.17 shows the schematic diagram of this setup. The ex- periment setup is shown in Figure A.28 at page 148. Eight specimens were tested, and their mean force-displacement curve, along with the 95% confidence band, are presented in Figure 6.18. Probabilistic simulations were performed as in the previous cases. The mean force-displacement curves from simulations are presented in Figure 6.18a. The difference between simulations from Base and Opt1-3 models was less drastic as compared with the previous test. However, Base re- sulted in the largest error among the four models, i.e., 6.7% for the peak force and 16% for the area under the curve. The force-displacement curve exhibited a bend at d=4 mm, an indication of damage initiation, which was not shown on the mean experimental curve. On the other hand, all optimized models performed well in this test. Overall, Opt3 yielded the most reliable results with moderate errors for the peak force and area under the curve. 108 (a) (b) Figure 6.18 Disk bending with a simply supported boundary. The mean force-displacement curve and 95% confidence band by experiment and by probabilistic simulations. (a) Compares the experi- ment results with simulations of different optimization models. (b) Compares the mean and its 95% confidence band of Opt3 with experimental results. Inserts shows the failure modes in experiment and simulation. 109 Figure 6.19 Disk bending with a simply supported boundary. The typical major stress, damage distributions at the bottom surface at different displacements predicted by Opt3. Figure 6.18b compares the results of Opt3 with the experiment. The red shaded area represents the 95% confidence band of Opt3. Opt3 yielded only a 4% error for the peak force and a 9% error for the area under the curve. As in the fixed boundary case, Opt3 exhibited a larger uncertainty at a larger displacement range. Finally, the inserts in Figure 6.18b show that the predicted failure mode was similar to that seen in the experiment. Figure 6.19 presents the typical major stress and damage distributions at the lower surface ob- tained in simulations with Opt3 at different displacements. At d=5 mm, there were a few damaged elements, and a few elements had major stress above 240MPa. The areas over the bottom supports had very low stress, whereas the areas between the supports carried stress ranging from 60-120MPa. At the peak load (d=11 mm), the damage started to appear at the disk center and along three radially extending lines between the supports. However, most of the elements still carried stress at the 120-200MPa range. With increasing the displacement to d=15 mm, the damage became focused, and the stress was reduced in more elements. At d=18 mm displacement, most elements between the supports reached the maximum damage. 110 Base Opt1 Opt2 Opt3 80 Hole Disk2 70 Disk1 60 50 Error [%] 40 30 20 10 0 PFSE PFSE PFSE PFSE SDE SDE SDE SDE AUCE AUCE AUCE AUCE Figure 6.20 Comparison of errors of Base and Opt1-3 models for three validation load cases. 6.5.3 Model Validation Summary Table 6.5 summarizes the errors for Base and Opt1-3 models for all three validation cases. Figure 6.20 compares these errors in the form of a bar chart. Table 6.5 Error (%) for the Peak Force/Stress Error (PFSE), Strain/Displacement Error (SDE) at the peak force/stress, and the Area Under the Curve Error (AUCE). ID Base Opt1 Opt2 Opt3 PFSE SDE AUCE PFSE SDE AUCE PFSE SDE AUCE PFSE SDE AUCE Hole 15.1 15.9 16.5 8.6 17.8 22.4 10.2 6.9 11.9 1.2 7.9 10.5 Disk1 21.4 0.3 48.6 9.8 0.7 5.9 0.3 8.6 18.1 3.6 1.6 11.8 Disk2 6.7 11.1 15.9 8.1 15 3.1 0.6 8.2 17.1 4.6 13 9.1 Figure 6.20 clearly shows that optimized models outperformed the baseline model in every validation case and in every category of errors. Among the optimized models, Opt3 had the lowest cumulative errors. In summary, the optimization significantly improved the model prediction. The model optimized for multiple cases performed better than models optimized with one load case. 111 6.6 Summary This work investigated the probabilistic simulations of SMC structures. The scatter of the tensile and compression strengths and elastic moduli have been considered. The size effect on the strength was modeled with a General Strength Distribution Model (GSDM). It is a mixture of the uniaxial strength and flexural strength Weibull models from tension and compression. Subsequently, the PFE model was optimized with a random search procedure to estimate hyperparameters from GSDM and MAT 297. The optimization was performed through correlating simulations with experimental results, i.e., (Opt1) optimized for the 3-pt bending case, (Opt2) for the 4-pt bending case, and (Opt3) for both 3-pt and 4-pt cases. The performances of the optimized models were verified and validated in probabilistic FE sim- ulation for seven mechanical characterization tests. The results showed that the probabilistic FE simulations with the randomization technique and the GSDM were successful in predicting a va- riety of structural responses of the SMC material. Furthermore, the optimization reduced errors in every category and provided much better predictions in every validation. Among the optimized models, Opt3 resulted in the lowest cumulative errors. In conclusion, the optimization significantly improved the model prediction. The model opti- mized for multiple cases performed better than the model optimized with one load case. 112 CHAPTER 7 SUMMARY AND CONCLUSION 7.1 Summary In this research, a Probabilistic Finite Element (PFE) simulation algorithm has been devel- oped to analyze Sheet Molding Compound (SMC). The stochastic nature of the material has been characterized using multimodal Weibull strength distribution models. The proposed algorithm is designed to be compatible with any finite element software and can be utilized in structural analysis for automobile crash simulation. The key findings of this study are as follows: 7.1.1 Mechanical Characteristics of SMC SMC composites manufactured from three molding processes, i.e., standard plaques, flow plaques, and knit plaques, have been characterized with standardized mechanical tests. Their key character- istics can be listed below: • Standard plaques demonstrated in-plane isotropic mechanical properties. In contrast, flow plaques and knit plaques were orthotropic. • All plaques exhibited substantial uncertainty in strength. Standard deviation for strength raged from 5-10% in all characterization tests. • SMC specimens showed notable size effect on strength, i.e., its strength largely depended on stress distribution. For example, SMC’s strength was 90% higher in 3-point bending and 50% higher in 4-point bending compared to the strength obtained in uniaxial tension. • SMC failed abruptly in the uniaxial tension and uniaxial compression test. In comparison, it failed progressively during the 3-point, 4-point bending test. 7.1.2 PFE Model With Random Tensile Strength The stochastic nature of SMC and the size effect in strength has been modeled with an Extended Strength Distribution Model (ESDM). The ESDM is built upon a trimodal tensile strength distribu- tion that covers a broader distribution than unimodal Weibull models. In addition, a randomization technique was developed to incorporate the statistical strength distribution in FE models. The pro- 113 posed method was tested in PFE simulations of tensile and flexural experiments. This strategy ob- tained desired results with as few as ten iterations, indicating outstanding computational efficiency over existing PFE simulation techniques. 7.1.3 Optimization of PFE Model With Machine Learning The General Strength Distribution Model (GSDM) has been introduced as a means of assimi- lating the uncertainty of tensile and compression strengths in SMC materials. GSDM is based on two bimodal Weibull distributions that are correlated by Weibull survival probability. One distri- bution was developed for tension, and the other for compression, with each derived from unimodal Weibull models of corresponding uniaxial and flexural strength. In addition, an optimization procedure has been developed for probabilistic FE models using an artificial neural network (ANN)-based algorithm. The goal of model optimization was to esti- mate parameters for GSDM and damage mechanics, rather than relying on manual guesswork. The validation results indicated that PFE simulation with GSDM and optimization could accurately assimilate experimental results. 7.1.4 Optimization of PFE Model With Random Search An optimization scheme based on random search (RS) has been investigated to minimize error and improve efficiency in PFE simulations. Optimization was performed using randomization from GSDM for three different cases of PFE simulations. The performance of the optimized models was verified and validated through simulation of seven mechanical characterization tests. The results showed that model optimization led to a reduction in error for every test compared to the non- optimized model. Additionally, it was found that optimization was more reliable when multiple cases were optimized simultaneously. 7.2 Conclusion Sheet Molding Compound (SMC) is a Discontinuous Fiber Reinforced Polymeric (DFRP) com- posite material with remarkable potential in the automotive industry. SMC is known for its out- standing strength-to-weight ratio and remarkable cost-effectiveness. Additionally, SMC has ex- 114 cellent resistance to corrosion, impact, and thermal expansion, making it ideal for use in harsh environments. Using SMC in automobiles can reduce weight, improve fuel efficiency, and reduce emissions. Nonetheless, SMC’s mechanical characteristics have always been challenging to predict in sim- ulation. One of the crucial challenges is size effect on strength. Weibull’s statistical strength theory well explains this phenomenon. However, it is not administered in the FE model for structural analysis. Therefore, this research has been dedicated to incorporating well-known size effects in FE analysis. In conclusion, the following remarks can be made: • Uncertainty of strength and size effect in SMC can be described by multimodal strength distribution. In this regard, two models have been developed, i.e., ESDM and GSDM. • ESDM is the simplest form among the models developed in this work. Additionally, it is very efficient computationally as fewer iterations are involved and requires no optimization. This model works well when a crash simulation does not involve compression failure. However, its prediction could be unreliable as the accuracy is inconsistent for various test cases. Overall, ESDM can be helpful when a computation cost is prioritized over accuracy. • GSDM model outperforms ESDM in terms of accuracy when appropriately optimized. Fur- thermore, PFE simulations with GSDM are more reliable and consistent despite increased computation costs. In the case of optimization, random search is the easiest to implement and faultless. Generally, GSDM would be preferred over ESDM for automobile crash simu- lation. 7.3 Future Work Opportunity Despite success in various verification and validation tests, the PFE simulation algorithm is not entirely settled. Several opportunities exist for model improvement, extension and adaptation to serve as a general-purpose crashworthiness model for discontinuous continuous fiber-reinforced polymeric (DFRP) composites. 115 7.3.1 Model Improvement • Optimize some other hyperparameters that are not addressed in this investigation, i.e. residual compressive strength, randomness in fiber angle etc. • Incorporate trimodal Weibull distribution in general strength distribution model. • Furthermore, improve the algorithm to balance between accuracy and computation cost. 7.3.2 Model Extension The model described in the research is based on in-plane isotropic SMC materials. It can be extended to anisotropic SMC materials. The following works are required: • Develop a GSDM for each orthogonal direction of anisotropic SMC materials. • Develop techniques to control the degree of anisotropy in the strength distribution models so that they can serve all SMC materials. • Estimate hyperparameters and perform validation with mechanical characterization tests. • Develop techniques to incorporate isotropic and anisotropic models in a given structure as real-life SMC materials contain both. • Implement, test, optimize and validate this algorithm in other finite element software pack- ages. 7.3.3 Adaptation to Other DFRP Composites Currently, the model is established for SMC. However, it may be applicable for various types of DFRP composites, which are unpredictable through deterministic simulations, have similar anisotropic characteristics; demonstrate similar stochastic nature; have extensive size effect on strength. 116 BIBLIOGRAPHY [1] P Beardmore and CF Johnson. The potential for composites in structural automotive appli- cations. Composites science and Technology, 26(4):251–281, 1986. [2] Ahmed Elmarakbi. Advanced composite materials for automotive applications: Structural integrity and crashworthiness. John Wiley & Sons, 2013. [3] Zhao Liu, Jiahai Lu, and Ping Zhu. Lightweight design of automotive composite bumper system using modified particle swarm optimizer. Composite Structures, 140:630–643, 2016. [4] George C Jacob, John F Fellers, Srdan Simunovic, and J Michael Starbuck. Energy absorp- tion in polymer composites for automotive crashworthiness. Journal of composite materials, 36(7):813–850, 2002. [5] George C Jacob, J Michael Starbuck, John F Fellers, Srdan Simunovic, and Raymond G Boeman. Crashworthiness of various random chopped carbon fiber reinforced epoxy com- posite materials and their strain rate dependence. Journal of applied polymer science, 101 (3):1477–1486, 2006. [6] Loredana Kehrer, Jeffrey T Wood, and Thomas Böhlke. Mean-field homogenization of ther- moelastic material properties of a long fiber-reinforced thermoset and experimental investi- gation. Journal of Composite Materials, 54(25):3777–3799, 2020. [7] Tseng S.-C Osswald T. Flow and rheology in polymer composites manufacturing. Compos- ites materials series, 10(10):477–484, 1994. [8] JL Thomason, MA Vlug, G Schipper, and HGLT Krikor. Influence of fibre length and concentration on the properties of glass fibre-reinforced polypropylene: Part 3. strength and strain at failure. Composites Part A: Applied Science and Manufacturing, 27(11):1075–1084, 1996. [9] V Feuillade, Anne Bergeret, J-C Quantin, and A Crespy. Relationships between the glass fibre sizing composition and the surface quality of sheet moulding compounds (smc) body panels. Composites Science and Technology, 66(1):115–127, 2006. [10] Marvin Knight and HT Hahn. Strength and elastic modulus of a randomly-distributed short fiber composite. Journal of Composite Materials, 9(1):77–90, 1975. [11] AM Hartl, M Jerabek, and RW Lang. Anisotropy and compression/tension asymmetry of pp containing soft and hard particles and short glass fibers. Express Polymer Letters, 9(7): 658–670, 2015. [12] AM Hartl, M Jerabek, P Freudenthaler, and RW Lang. Orientation-dependent compres- sion/tension asymmetry of short glass fiber reinforced polypropylene: Deformation, damage 117 and failure. Composites Part A: Applied Science and Manufacturing, 79:14–22, 2015. [13] Justin D Littell, Charles R Ruggeri, Robert K Goldberg, Gary D Roberts, William A Arnold, and Wieslaw K Binienda. Measurement of epoxy resin tension, compression, and shear stress–strain curves over a wide range of strain rates using small test specimens. Journal of Aerospace Engineering, 21(3):162–173, 2008. [14] MK Cattell and Kevin A Kibble. Determination of the relationship between strength and test method for glass fibre epoxy composite coupons using weibull analysis. Materials & Design, 22(4):245–250, 2001. [15] A Ionita and YJ Weitsman. On the mechanical response of randomly reinforced chopped- fibers composites: Data and model. Composites science and technology, 66(14):2566–2579, 2006. [16] Christos C Chamis. Probabilistic assessment of composite structures, volume 106024. NASA, 1993. [17] MW Long and JD Narciso. Probabilistic design methodology for composite aircraft struc- tures. Technical report, NORTHROP GRUMMAN CORP DALLAS TX COMMERCIAL AIRCRAFT DIV, 1999. [18] Srinivas Sriramula and Marios K Chryssanthopoulos. Quantification of uncertainty mod- elling in stochastic analysis of frp composites. Composites Part A: Applied Science and Manufacturing, 40(11):1673–1684, 2009. [19] Jacques W Nader, Habib J Dagher, Roberto Lopez-Anido, Fadi El Chiti, Ghassan N Fayad, and Lawrence Thomson. Probabilistic finite element analysis of modified astm d3039 ten- sion test for marine grade polymer matrix composites. Journal of reinforced plastics and composites, 27(6):583–597, 2008. [20] Jacques W Nader, Habib J Dagher, Fadi El Chiti, and Roberto Lopez-Anido. Probabilistic finite element analysis of astm d6641 compression test for marine grade polymer matrix composites. Journal of reinforced plastics and composites, 28(8):897–911, 2009. [21] Jacques W Nader, Habib J Dagher, and Roberto Lopez-Anido. Size effects on the bending strength of fiber-reinforced polymer matrix composites. Journal of reinforced plastics and composites, 30(4):309–317, 2011. [22] EMIL J Pitz and KISHORE V Pochiraju. Quasi monte carlo simulations for stochastic fail- ure analysis in composites. In Proceedings of the American Society for Composites-Thirty- Fourth Technical Conference, 2019. [23] Hoil Choi, Seungwoo Jung, Chao Zhang, and Gun Jin Yun. A three-dimensional stochastic progressive damage simulation model for polymer matrix-based laminate composites. Me- 118 chanics of Advanced Materials and Structures, 29(5):633–650, 2022. [24] TS Mesogitis, Alexandros A Skordos, and AC Long. Uncertainty in the manufacturing of fibrous thermosetting composites: A review. Composites Part A: Applied Science and Man- ufacturing, 57:67–75, 2014. [25] Paolo Feraboli, Tyler Cleveland, Patrick Stickler, and John Halpin. Stochastic laminate anal- ogy for simulating the variability in modulus of discontinuous composite materials. Com- posites Part A: Applied Science and Manufacturing, 41(4):557–570, 2010. [26] Marina Selezneva, Steven Roy, Sean Meldrum, Larry Lessard, and Ali Yousefpour. Mod- elling of mechanical properties of randomly oriented strand thermoplastic composites. Jour- nal of Composite Materials, 51(6):831–845, 2017. [27] Zhangxing Chen, Tianyu Huang, Yimin Shao, Yang Li, Hongyi Xu, Katherine Avery, Danielle Zeng, Wei Chen, and Xuming Su. Multiscale finite element modeling of sheet mold- ing compound (smc) composite structure based on stochastic mesostructure reconstruction. Composite Structures, 188:25–38, 2018. [28] Zhangxing Chen, Haibin Tang, Yimin Shao, Qingping Sun, Guowei Zhou, Yang Li, Hongyi Xu, Danielle Zeng, and Xuming Su. Failure of chopped carbon fiber sheet molding com- pound (smc) composites under uniaxial tensile loading: Computational prediction and ex- perimental analysis. Composites Part A: Applied Science and Manufacturing, 118:117–130, 2019. [29] Sergii G Kravchenko, Drew E Sommer, Benjamin R Denos, Anthony J Favaloro, Clark M Tow, William B Avery, and R Byron Pipes. Tensile properties of a stochastic prepreg platelet molded composite. Composites Part A: Applied Science and Manufacturing, 124:105507, 2019. [30] Yizhuo Li and Soraia Pimenta. Development and assessment of modelling strategies to predict failure in tow-based discontinuous composites. Composite Structures, 209:1005– 1021, 2019. [31] Luigi Nicolais, Assunta Borzacchiello, and Stuart M Lee. Wiley encyclopedia of composites. Wiley Online Library, 2012. [32] Thomas Böhlke, Frank Henning, Andrew N Hrymak, Luise Kärger, Kay Weidenmann, and Jeffrey T Wood. Continuous–discontinuous fiber-reinforced polymers: an integrated engi- neering approach. Carl Hanser Verlag GmbH Co KG, 2019. [33] Sari K Christensen, Bryan Hutchinson, Esther M Sun, Tim A Osswald, and BA Davis. Fiber- matrix separation in ribbed smc and bmc parts. In Technical Papers of the Annual Technical Conference-Society of Plastics Engineers Incorporated, volume 1, pages 782–787. Society of Plastics Engineers Inc, 1997. 119 [34] P Dumont, Laurent Orgéas, Denis Favier, Patrick Pizette, and Cécile Venet. Compression moulding of smc: In situ experiments, modelling and simulation. Composites Part A: Ap- plied Science and Manufacturing, 38(2):353–368, 2007. [35] Colin Servais, André Luciani, and Jan-Anders E Månson. Fiber–fiber interaction in concen- trated suspensions: dispersed fiber bundles. Journal of rheology, 43(4):1005–1018, 1999. [36] Eric Maire and Philip John Withers. Quantitative x-ray tomography. International materials reviews, 59(1):1–43, 2014. [37] T-H Le, Pierre JJ Dumont, L Orgéas, D Favier, Luc Salvo, and Elodie Boller. X-ray phase contrast microtomography for the analysis of the fibrous microstructure of smc composites. Composites Part A: Applied Science and Manufacturing, 39(1):91–103, 2008. [38] Sivasambu Mahesh, S Leigh Phoenix, and Irene J Beyerlein. Strength distributions and size effects for 2d and 3d composites with weibull fibers in an elastic matrix. International Journal of Fracture, 115:41–85, 2002. [39] Mark R Garnich and Ghodrat Karami. Finite element micromechanics for stiffness and strength of wavy fiber composites. Journal of Composite Materials, 38(4):273–292, 2004. [40] Shao-Yun Fu and Bernd Lauke. Effects of fiber length and fiber orientation distributions on the tensile strength of short-fiber-reinforced polymers. Composites Science and Technology, 56(10):1179–1190, 1996. [41] Seyyedvahid Mortazavian and Ali Fatemi. Effects of fiber orientation and anisotropy on ten- sile strength and elastic modulus of short fiber reinforced polymer composites. Composites part B: engineering, 72:116–129, 2015. [42] Heping Xie and Feng Gao. The mechanics of cracks and a statistical strength theory for rocks. International Journal of Rock Mechanics and Mining Sciences, 37(3):477–488, 2000. [43] Bikramjit Basu, Devesh Tiwari, Debasis Kundu, and Rajesh Prasad. Is weibull distribution the most appropriate statistical strength distribution for brittle materials? Ceramics Interna- tional, 35(1):237–246, 2009. [44] D Dierickx, Bikramjit Basu, Jef Vleugels, and Omer Van der Biest. Statistical extreme value modeling of particle size distributions: experimental grain size distribution type estimation and parameterization of sintered zirconia. Materials characterization, 45(1):61–70, 2000. [45] Chunsheng Lu, Robert Danzer, and Franz Dieter Fischer. Fracture statistics of brittle mate- rials: Weibull or normal distribution. Physical Review E, 65(6):067102, 2002. [46] Eckhard Limpert and Werner A Stahel. Problems with using the normal distribution–and ways to improve quality and efficiency of data analysis. PloS one, 6(7):e21403, 2011. 120 [47] Waloddi Weibull. A statistical theory of strength of materials. IVB-Handl., 1939. [48] Xiaosheng Gao and Robert H Dodds Jr. Constraint effects on the ductile-to-brittle transition temperature of ferritic steels: a weibull stress model. International Journal of Fracture, 102 (1):43–69, 2000. [49] J Andersons, Roberts Joffe, M Hojo, and S Ochiai. Glass fibre strength distribution de- termined by common experimental methods. Composites Science and Technology, 62(1): 131–145, 2002. [50] Enrique Barbero, José Fernández-Sáez, and C Navarro. Statistical analysis of the mechanical properties of composite materials. Composites Part B: Engineering, 31(5):375–381, 2000. [51] A de S Jayatilaka and K Trustrum. Statistical approach to brittle fracture. Journal of Mate- rials Science, 12:1426–1430, 1977. [52] K Trustrum and A De S Jayatilaka. Applicability of weibull analysis for brittle materials. Journal of Materials Science, 18:2765–2770, 1983. [53] J Margetson and NR Cooper. Brittle material design using three parameter weibull distri- butions. In Probabilistic Methods in the Mechanics of Solids and Structures: Symposium Stockholm, Sweden June 19–21, 1984 To the Memory of Waloddi Weibull, pages 253–262. Springer, 1985. [54] Gabriel J DeSalvo. Theory and structural design applications of weibull statistics. Technical report, Westinghouse Electric Corp., Pittsburgh, Pa. Astronuclear Lab., 1970. [55] Kailash Chander Kapur and Leonard R Lamberson. Reliability in engineering design. New York, 1977. [56] Osama M Jadaan, Noel N Nemeth, Joerg Bagdahn, and WN Sharpe. Probabilistic weibull behavior and mechanical properties of mems brittle materials. Journal of materials science, 38:4087–4113, 2003. [57] RE Bullock. Strength ratios of composite materials in flexure and in tension. Journal of Composite Materials, 8(2):200–206, 1974. [58] Karen E Jackson, Sotiris Kellas, and John Morton. Scale effects in the response and failure of fiber reinforced composite laminates loaded in tension and in flexure. Journal of composite materials, 26(18):2674–2705, 1992. [59] MR Wisnom and JW Atkinson. Reduction in tensile and flexural strength of unidirectional glass fibre-epoxy with increasing specimen size. Composite Structures, 38(1-4):405–411, 1997. 121 [60] Yongdong Xu, Laifei Cheng, Litong Zhang, Dantao Yan, and Chang You. Optimization of sample number for weibull function of brittle materials strength. Ceramics International, 27 (2):239–241, 2001. [61] JD Sullivan and PH Lauzon. Experimental probability estimators for weibull plots. Journal of Materials Science Letters, 5(12):1245–1247, 1986. [62] A Khalili and K Kromp. Statistical properties of weibull estimators. Journal of materials science, 26:6741–6752, 1991. [63] PC Gope. Determination of sample size for estimation of fatigue life by using weibull or log-normal distribution. International journal of fatigue, 21(8):745–752, 1999. [64] Koichi Goda and Hideharu Fukunaga. The evaluation of the strength distribution of sili- con carbide and alumina fibres by a multi-modal weibull distribution. Journal of Materials Science, 21:4475–4480, 1986. [65] Fuchen Mu, Changhua Tan, and Mingzhen Xu. Proportional difference estimate method of determining the characteristic parameters of monomodal and multimodal weibull distribu- tions of time-dependent dielectric breakdown. Solid-State Electronics, 44(8):1419–1424, 2000. [66] MEH Bourahli. Uni-and bimodal weibull distribution for analyzing the tensile strength of diss fibers. Journal of Natural Fibers, 15(6):843–852, 2018. [67] Jie Zhang, Souma Chowdhury, Achille Messac, and Luciano Castillo. Multivariate and multimodal wind distribution model based on kernel density estimation. In Energy Sustain- ability, volume 54686, pages 2125–2135, 2011. [68] Philippe Zinck, Jean-Franįois Gerard, and HD Wagner. On the significance and description of the size effect in multimodal fracture behavior. experimental assessment on e-glass fibers. Engineering fracture mechanics, 69(9):1049–1055, 2002. [69] CD Lai, DNP Murthy, and M Xie. Weibull distributions. Wiley Interdisciplinary Reviews: Computational Statistics, 3(3):282–287, 2011. [70] Jerald F Lawless. Statistical models and methods for lifetime data. John Wiley & Sons, 2011. [71] William Q Meeker, Luis A Escobar, and Francis G Pascual. Statistical methods for reliability data. John Wiley & Sons, 2022. [72] Aydin Karakoca, Ulku Erisoglu, and Murat Erisoglu. A comparison of the parameter es- timation methods for bimodal mixture weibull distribution with complete data. Journal of Applied Statistics, 42(7):1472–1489, 2015. 122 [73] Piotr Wais. A review of weibull functions in wind sector. Renewable and Sustainable Energy Reviews, 70:1099–1107, 2017. [74] Kanji Ono. A simple estimation method of weibull modulus and verification with strength data. Applied Sciences, 9(8):1575, 2019. [75] Peng Wang, Guozheng Ma, Fenghua Su, Weiling Guo, Shuying Chen, Haichao Zhao, Ming Liu, and Haidou Wang. Microstructure and micromechanical properties of in situ synthe- sized tio2- x coatings by plasma spraying of bimodal feedstocks. Journal of Thermal Spray Technology, 31(8):2300–2313, 2022. [76] TS Creasy. Modeling analysis of tensile tests of bundled filaments with a bimodal weibull survival function. Journal of composite materials, 36(2):183–194, 2002. [77] LH Peebles. Carbon fibres: Structure and mechanical properties. International materials reviews, 39(2):75–92, 1994. [78] WANG Tianshuai, HE Xiaofan, WANG Jinyu, and LI Yuhai. Detail fatigue rating method based on bimodal weibull distribution for ded ti-6.5 al-2zr-1mo-1v titanium alloy. Chinese Journal of Aeronautics, 35(4):281–291, 2022. [79] David Cohen. Application of reliability and fiber probabilistic strength distribution concepts to composite vessel burst strength design. Journal of composite materials, 26(13):1984– 2014, 1992. [80] Stephen I Gallant et al. Perceptron-based learning algorithms. IEEE Transactions on neural networks, 1(2):179–191, 1990. [81] Matt W Gardner and SR Dorling. Artificial neural networks (the multilayer perceptron)—a review of applications in the atmospheric sciences. Atmospheric environment, 32(14-15): 2627–2636, 1998. [82] Frank Rosenblatt. The perceptron: a probabilistic model for information storage and orga- nization in the brain. Psychological review, 65(6):386, 1958. [83] Aoyi Luo, Hang Zhang, and Kevin T Turner. Machine learning-based optimization of the design of composite pillars for dry adhesives. Extreme Mechanics Letters, 54:101695, 2022. [84] Saeed Gholizadeh, Akbar Pirmoz, and Reza Attarnejad. Assessment of load carrying capac- ity of castellated steel beams by neural networks. Journal of Constructional Steel Research, 67(5):770–779, 2011. [85] Xin Liu, Su Tian, Fei Tao, and Wenbin Yu. A review of artificial neural networks in the constitutive modeling of composite materials. Composites Part B: Engineering, 224:109152, 2021. 123 [86] Dean C Karnopp. Random search techniques for optimization problems. Automatica, 1(2-3): 111–121, 1963. [87] Paulo Henrique Ranazzi and Marcio Augusto Sampaio Pinto. Assisted history matching us- ing combined optimization methods. In Proceedings of the Iberian Latin-American Congress on Computational Methods in Engineering, Florianópolis, Brazil, pages 5–8, 2017. [88] WL1551847 Price. Global optimization by controlled random search. Journal of optimiza- tion theory and applications, 40:333–348, 1983. [89] Sigrún Andradóttir. A review of random search methods. Handbook of simulation optimiza- tion, pages 277–292, 2014. [90] Keith Seymour, Haihang You, and Jack Dongarra. A comparison of search heuristics for empirical code optimization. In 2008 IEEE International Conference on Cluster Computing, pages 421–429. IEEE, 2008. [91] Seyedali Mirjalili and Andrew Lewis. Obstacles and difficulties for robust benchmark prob- lems: A novel penalty-based robust optimisation method. Information Sciences, 328:485– 509, 2016. [92] Sigrún Andradóttir. An overview of simulation optimization via random search. Handbooks in operations research and management science, 13:617–631, 2006. [93] David Landau and Kurt Binder. A guide to Monte Carlo simulations in statistical physics. Cambridge university press, 2021. [94] Robert L Harrison. Introduction to monte carlo simulation. In AIP conference proceedings, volume 1204, pages 17–21. American Institute of Physics, 2010. [95] Hiroshi Fukuda and Kozo Kawata. On the strength distribution of unidirectional fibre com- posites. Fibre Science and Technology, 10(1):53–63, 1977. [96] I Kimpara and T Ocaki. Dynamic simulations on stochastic tensile failure process of unidi- rectional fiber reinforced composites. Proceedings of the ISCMS (1st), Beijing, pages 682– 92, 1986. [97] Koichi Goda. The role of interfacial debonding in increasing the strength and reliability of unidirectional fibrous composites. Composites science and technology, 59(12):1871–1879, 1999. [98] M Ibnabdeljalil and SL Phoenix. Scalings in the statistical failure of brittle matrix composites with discontinuous fibers—i. analysis and monte carlo simulations. Acta metallurgica et materialia, 43(8):2975–2983, 1995. 124 [99] Sivasambu Mahesh, Irene J Beyerlein, and S Leigh Phoenix. Size and heterogeneity effects on the strength of fibrous composites. Physica D: Nonlinear Phenomena, 133(1-4):371–389, 1999. [100] Chad M Landis, Irene J Beyerlein, and Robert M McMeeking. Micromechanical simulation of the failure of fiber reinforced composites. Journal of the Mechanics and Physics of Solids, 48(3):621–648, 2000. [101] Jianming Yuan, Yuanming Xia, and Baochang Yang. A note on the monte carlo simulation of the tensile deformation and failure process of unidirectional composites. Composites science and technology, 52(2):197–204, 1994. [102] Zhen Wang, Jianming Yuan, and Yuanming Xia. A dynamic monte-carlo simulation for unidirectional composites under tensile impact. Composites science and technology, 58(3- 4):487–495, 1998. [103] Yuanxin Zhou, Wen Huang, and Yuanming Xia. A microscopic dynamic monte carlo sim- ulation for unidirectional fiber reinforced metal matrix composites. Composites science and technology, 62(15):1935–1946, 2002. [104] Tianle Cheng, Rui Qiao, and Yuanming Xia. A monte carlo simulation of damage and failure process with crack saturation for unidirectional fiber reinforced ceramic composites. Com- posites science and technology, 64(13-14):2251–2260, 2004. [105] ASTM D2584. Standard test method for ignition loss of cured reinforced resins. In American society for testing and materials, volume 100, 2002. [106] James M Whitney and Marvin Knight. The relationship between tensile strength and flexure strength in fiber reinforced composites. Technical report, AIR FORCE WRIGHT AERO- NAUTICAL LABS WRIGHT-PATTERSON AFB OH, 1980. [107] ASTM Committee D-30 on Composite Materials. Standard test method for tensile properties of polymer matrix composite materials. ASTM international, 2008. [108] Inernational ASTM. Standard test methods for flexural properties of unreinforced and rein- forced plastics and electrical insulating materials. ASTM D790-07, 2007. [109] ASTM International. Standard test method for flexural properties of polymer matrix com- posite materials. ASTM International, 2015. [110] D ASTM. Standard test method for shear properties of composite materials by the v-notched beam method. Annu. Book ASTM Stand, 15:227–1993, 2000. [111] Sakib Iqbal, Beichen Li, Kestutis Sonta, Royal Chibuzor Ihuaenyi, and Xinran Xiao. Prob- abilistic finite element analysis of sheet molding compound composites with an extended 125 strength distribution model. Finite Elements in Analysis and Design, 214:103865, 2023. [112] Waloddi Weibull. A statistical distribution function of wide applicability. Journal of applied mechanics, 1951. [113] Wayne B Nelson. Applied life data analysis. John Wiley & Sons, 2005. [114] JL Thomason, Chih-Chuan Kao, J Ure, and Liu Yang. The strength of glass fibre reinforce- ment after exposure to elevated composite processing temperatures. Journal of Materials Science, 49:153–162, 2014. [115] Dieter Loidl, Oskar Paris, H Rennhofer, Martin Müller, and Herwig Peterlik. Skin-core structure and bimodal weibull distribution of the strength of carbon fibers. Carbon, 45(14): 2801–2805, 2007. [116] Zhen Wang and Yuanming Xia. Experimental evaluation of the strength distribution of fibers under high strain rates by bimodal weibull distribution. Composites science and technology, 57(12):1599–1607, 1998. [117] Jun Liu and Paul Bowen. Strength of scs-6 fibers in a ti (beta) 21s/scs-6 composite before and after fatigue. Metallurgical and Materials Transactions, 34(5):1203, 2003. [118] Linli Zhu, Xiang Guo, Haihui Ruan, and Jian Lu. Prediction of mechanical properties in bi- modal nanotwinned metals with a composite structure. Composites Science and Technology, 123:222–231, 2016. [119] C Robin, Susanne S Scherrer, HWA Wiskott, WG De Rijk, and UC Belser. Weibull parame- ters of composite resin bond strengths to porcelain and noble alloy using the rocatec system. Dental Materials, 18(5):389–395, 2002. [120] Jia-jie Kang, Jian-long Ma, Guo-lu Li, Hai-dou Wang, Bin-shi Xu, and Cheng-biao Wang. Bimodal distribution characteristic of microstructure and mechanical properties of nanos- tructured composite ceramic coatings prepared by supersonic plasma spraying. Materials & Design, 64:755–759, 2014. [121] BA Kahl, CC Berndt, and ASM Ang. Boride-based ultra-high temperature ceramic coatings deposited via controlled atmosphere plasma spray. Surface and Coatings Technology, 416: 127128, 2021. [122] Donald A Mumford. Robust parameter estimation for the mixed weibull (seven parameter) including the method of minimum likelihood and the method of minimum distance. Tech- nical report, AIR FORCE INST OF TECH WRIGHT-PATTERSON AFB OH, 1997. [123] Sakib Iqbal, Xinran Xiao, Beichen Li, and Kestutis Sonta. Considering the randomness of mechanical properties in simulations of discontinuous fiber reinforced composites. In 126 Proceedings of the American Society for Composites Thirty-fifth Technical Conference, 2020. [124] José David Arregui-Mena, Lee Margetts, and Paul M Mummery. Practical application of the stochastic finite element method. Archives of computational methods in engineering, 23: 171–190, 2016. [125] John O Hallquist et al. Ls-dyna theory manual. Livermore software Technology corporation, 3:25–31, 2006. [126] Danghe Shi and Xinran Xiao. An enhanced continuum damage mechanics model for crash simulation of composites. Composite Structures, 185:774–785, 2018. [127] Xinran Xiao, Carla McGregor, Reza Vaziri, and Anoush Poursartip. Progress in braided composite tube crush simulation. International Journal of Impact Engineering, 36(5):711– 719, 2009. [128] Xinran Xiao, Mark E Botkin, and Nancy L Johnson. Axial crush simulation of braided carbon tubes using mat58 in ls-dyna. Thin-Walled Structures, 47(6-7):740–749, 2009. [129] Sakib Iqbal and Xinran Xiao. Optimization of fe crashworthiness model for sheet mold- ing compound (smc) with extended strength distribution model and machine learning. In Proceedings of The American Society For Composites-Thirty-Seventh Technical Conference, 2022. [130] Jasper Snoek, Hugo Larochelle, and Ryan P Adams. Practical bayesian optimization of machine learning algorithms. Advances in neural information processing systems, 25, 2012. [131] George A Bekey and Sami F Masri. Random search techniques for optimization of nonlinear systems with many parameters. Mathematics and Computers in Simulation, 25(3):210–213, 1983. [132] Nicholas Kantas, Arnaud Doucet, Sumeetpal Sindhu Singh, and Jan Marian Maciejowski. An overview of sequential monte carlo methods for parameter estimation in general state- space models. IFAC Proceedings Volumes, 42(10):774–785, 2009. [133] Dennis L Jackson. Sample size and number of parameter estimates in maximum likelihood confirmatory factor analysis: A monte carlo investigation. Structural Equation Modeling, 8 (2):205–223, 2001. [134] Elchanan Mossel and Eric Vigoda. Limitations of markov chain monte carlo algorithms for bayesian inference of phylogeny. 2006. [135] Cleve Ashcraft, Jef Dawson, Roger Grimes, Erman Guleryuz, Seid Koric, Robert Lucas, James Ong, Francois-Henry Rouet, Todd Simons, and Ting-Ting Zhu. Running jet engine models on thousands of processors with ls-dyna implicit. In 12th European LS-DYNA Con- 127 ference, 2019. [136] Francois-Henry Rouet, Cleve Ashcraft, Jef Dawson, Roger Grimes, Erman Guleryuz, Seid Koric, Robert F Lucas, James S Ong, Todd A Simons, and Ting-Ting Zhu. Scalability challenges of an industrial implicit finite element code. In 2020 IEEE International Parallel and Distributed Processing Symposium (IPDPS), pages 505–514. IEEE, 2020. [137] Ting-Ting Zhu and Jason Wang. Ls-dyna scalability analysis on cray supercomputers. In 13th International LS-DYNA Users Conference. LSTC, 2020. [138] Inernational ASTM. Standard test method for open-hole tensile strength of polymer matrix composite laminates. ASTM D5766, 2018. [139] Stephen R Hallett, Ben G Green, Wen-Guang Jiang, Kin Hei Cheung, and Michael R Wis- nom. The open hole tensile test: a challenge for virtual testing of composites. International Journal of Fracture, 158:169–181, 2009. [140] Emmanuelle Abisset, Federica Daghia, and Pierre Ladevèze. On the validation of a damage mesomodel for laminated composites by means of open-hole tensile tests on quasi-isotropic laminates. Composites Part A: Applied Science and Manufacturing, 42(10):1515–1524, 2011. [141] Alfredo Balacó de Morais. Open-hole tensile strength of quasi-isotropic laminates. Com- posites Science and Technology, 60(10):1997–2004, 2000. [142] BY Chen, TE Tay, PM Baiz, and ST Pinho. Numerical analysis of size effects on open-hole tensile composite laminates. Composites Part A: Applied Science and Manufacturing, 47: 52–62, 2013. [143] Royal Chibuzor Ihuaenyi, Shutian Yan, Jie Deng, Chulheung Bae, Iqbal Sakib, and Xinran Xiao. Orthotropic thermo-viscoelastic model for polymeric battery separators with elec- trolyte effect. Journal of The Electrochemical Society, 168(9):090536, 2021. [144] Sakib Iqbal, Beichen Li, Kestutis Sonta, Royal Chibuzor Ihuaenyi, and Xinran Xiao. Simu- lations of sheet molding compound composite structures with a randomly distributed tensile strength. 2021. 128 APPENDIX A SUPPLEMENTARY EXPERIMENTAL RESULTS Uniaxial Tension Figure A.1 Sampling plan for uniaxial tension test from standard SMC plaques. Red marked labels are specimens for tension tests. Table A.1 Dimensions of standard plaque specimens for tension tests (mm). ID TT11 TT12 TT13 TT14 TL11 TL12 TL13 TL14 Width 20.02 20.05 20.05 19.92 19.92 20.02 20.01 20.07 Thickness 3.061 3.112 3.333 3.28 2.901 2.989 3.242 3.278 Gage length 46.47 44.20 45.85 44.17 47.392 46.351 49.633 45.678 Table A.2 Dimensions of flow plaque specimens in transverse direction for tension tests (mm). ID TT15 TT16 TT17 TT18 TT19 TT20 TL13 TL14 Width 3.348 3.559 3.115 3.308 2.665 2.826 20.01 20.07 Thickness 15.71 15.55 15.14 15.45 15.06 15.15 3.242 3.278 Gage length 42.03 40.96 42.77 40.17 39.75 44.65 49.633 45.678 129 Table A.3 Dimensions of flow plaque specimens in longitudinal direction for tension tests. ID TL15 TL16 TL17 TL18 TL19 TL20 TL21 TL22 TL23 TL24 TL25 TL26 b 3.111 3.127 3.317 3.324 3.419 3.428 2.798 2.853 2.943 2.973 3.023 3.035 d 15.39 15.71 15.69 13.99 14.88 15.83 17.44 13.81 15.1 14.27 15.27 13.16 L 62.41 62.40 63.91 64.31 62.76 64.24 64.12 64.45 59.91 65.38 68.31 62.45 b= width of specimen; d= thickness/depth of specimen; L= gage length of specimen. Table A.4 Dimensions of knit plaque specimens for tension test. ID TL01 TL02 TL04 TL05 Width 19.71 19.61 22.45 20.51 Thickness 3.142 3.254 3.438 3.414 Gage length 41.97 44.09 47.38 40.99 Displacement from Laser extensometer (mm) 2 TT1 TT2 TT3 TT4 TT12 TT13 TT14 1.5 1 0.5 0 0 0.5 1 1.5 2 2.5 Machine Displacement (mm) Figure A.2 Calibration curve for machine displacement during tension test of SMC standard plaques. 130 200 180 160 140 Stress (MPa) 120 TL15 100 TL16 TL18 80 TL19 TL20 60 TL21 TL24 40 TL25 20 TL26 0 0 0.005 0.01 0.015 0.02 Strain (mm/mm) Figure A.3 SMC tensile test results of flow plaques longitudinal direction. 100 90 TT15 80 TT16 70 Stress (MPa) TT17 60 50 TT18 40 TT19 30 TT20 20 10 0 0 0.005 0.01 0.015 0.02 Strain (mm/mm) Figure A.4 SMC tensile test results of flow plaques transverse direction. 131 20 TL01 18 TL02 16 14 TL04 Stress (MPa) 12 TL05 10 8 6 4 2 0 0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 Strain (mm/mm) Figure A.5 SMC tensile test results of knit plaques. Standard Plaque Transverse Direction Flow Plaque Longitudinal Direction 150 Flow Plaque Transverse Direction Knit Plaque 120 Stress (MPa) 90 60 30 0 0 0.005 0.01 0.015 0.02 Strain (mm/mm) Figure A.6 Tensile test results of SMC plaques. 132 Figure A.7 Tensile failure mode for (left) Standard Plaques (right) Knit plaques. 133 Figure A.8 Tensile strain field measured by DIC. (a) Standard Plaques (b) Knit plaques. -0.15 Poison’s ratio = 0.2256 -0.12 Transverse strain -0.09 -0.06 -0.03 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 Longitudinal strain Figure A.9 Poisson’s ratio measurement for standard plaques from DIC. 134 Uniaxial Compression Figure A.10 Compression test specimens from standard plaques. Table A.5 Dimensions of compression specimens from standard plaques and knit plaques (mm). Standard Plaque Knit Plaque ID C1 C2 C3 C4 C5 C6 C32 C33 C34 Width 25.69 25.03 23.56 25.02 24.5 23.91 15.77 15.71 14.88 Thickness 3.404 3.256 3.248 3.144 3.291 3.195 3.267 3.119 3.206 Gage Length 11.97 19.83 19.34 19.30 22.73 21.78 20.72 21.04 24.83 Table A.6 Dimensions of specimens from flow plaques in transverse direction for compression test. ID C12 C13 C14 C20 C21 C22 C29 C30 C31 Thickness 3.441 3.522 3.299 2.991 3.186 3.238 2.482 2.637 2.716 Width 15.7 15.93 15.15 13.73 13.42 13.18 17.18 17.28 17.37 Gage Length 21.78 27.84 26.20 15.70 14.33 12.59 12.93 15.02 16.47 Table A.7 Dimensions of specimens from flow plaques in transverse direction for compression test. ID C15 C16 C17 C18 C19 C23 C24 C25 C26 C27 C28 d 3.236 3.407 3.338 3.516 3.403 3.034 2.848 3.165 2.988 3.181 2.923 b 15.43 14.24 14.36 15.04 15.04 16.15 15.74 14.36 14.01 14.91 14.88 L 23.919 19.914 14.895 15.903 14.168 16.063 14.088 13.63 14.227 14.498 17.133 135 -1 Laser extensometer output voltage -0.9 C12 C13 C14 C15 C16 C17 C18 C19 -0.8 C20 C21 C22 -0.7 -0.6 -0.5 -0.4 -0.3 -0.2 -0.1 0 0 -0.2 -0.4 -0.6 -0.8 -1 -1.2 Machine displacement Figure A.11 Calibration curve of machine displacement for compression tests of SMC standard plaques. -210 Compressive Stress (MPa) -160 Solid lines represents specimen along longitudina direction & dotted line are transverse direction -110 -60 -10 0 -0.005 -0.01 -0.015 -0.02 -0.025 -0.03 -0.035 -0.04 Strain (mm/mm) Figure A.12 Compression test results of flow plaques. 136 -90 C32 C33 C34 -80 Compressive Stress (MPa) -70 -60 -50 -40 -30 -20 -10 0 0 -0.005 -0.01 -0.015 -0.02 -0.025 -0.03 -0.035 -0.04 Strain (mm/mm) Figure A.13 Compression test results of knit plaques. -200 Standard Plaque -180 Flow Plaque Longitudinal Direction Compressive Stress (MPa) -160 Flow Plaque Transverse Direction -140 Knit Plaque -120 -100 -80 -60 -40 -20 0 0 -0.005 -0.01 -0.015 -0.02 -0.025 Strain (mm/mm) Figure A.14 Compression test results of glass fibers SMC plaques. 137 3-Pt Bending Table A.8 Dimensions of specimens from standard and knit plaques for 3-pt bending tests. Standard Plaque Knit Plaque ID Thickness Width Span ID Thickness Width Span B1 3.394 12.68 52 B01 3.474 11.67 52 B2 3.241 11.58 52 B02 3.483 10.71 52 B3 3.181 13.38 52 B03 3.493 12.73 52 B4 3.251 11.29 52 B04 3.549 11.2 52 B6 3.383 13.05 52 B05 3.524 13.72 52 B7 3.143 11.89 52 B06 3.548 14.56 52 B8 3.124 12.55 52 - - - - Table A.9 Dimensions of flow specimen for 3-pt bending tests (longitudinal direction). ID B18 B19 B20 B21 B22 B23 B27 B28 B29 B30 B31 B32 d 3.266 3.203 3.386 3.367 3.354 3.394 2.986 2.799 3.126 2.851 3.092 3.028 b 15.54 15.18 14.66 15.02 15.39 15.47 13.6 13.82 14.61 14.4 15.38 15.46 L 52 52 52 52 52 52 52 52 52 52 52 52 b= width of specimen; d= thickness/depth of specimen; L= span of test fixture. Table A.10 Dimensions of flow specimen for 3-pt bending tests (Transverse direction). ID B15 B16 B17 B24 B25 B26 B33 B34 B35 d 3.31 3.419 3.54 3.026 3.25 3.271 2.554 2.702 2.808 b 13.93 13.92 13.87 15.75 15.69 15.63 13.98 13.98 14.05 L 52 52 52 52 52 52 52 52 52 138 450 B15 B16 B17 B18 B19 B20 400 B21 B22 B23 B24 B25 B26 Flexure strength (MPa) 350 B27 B28 B29 300 B30 B31 B32 B33 B34 B35 250 200 150 100 50 0 0 0.02 0.04 0.06 0.08 0.1 Flexure strain (mm/mm) Figure A.15 3-pt bending test results of Flow Plaques. 100 B01 B02 B03 90 80 B04 B05 B06 Flexure strength (MPa) 70 60 50 40 30 20 10 0 0 0.01 0.02 0.03 0.04 0.05 Flexure strain (mm/mm) Figure A.16 3-pt bending test results of Knit Plaques. 139 350 Standard Plaque 300 Knit Plaque Flexure Strength (MPa) 250 Flow Plaque longitudinal direct 200 Flow Plaque Transverse Directio 150 100 50 0 0 0.02 0.04 0.06 0.08 0.1 Flexure Strain (mm/mm) Figure A.17 3-pt bending test results of glass fibers SMC plaques. 140 Figure A.18 Failure mode during 3-pt bending tests for (left) Standard plaques (right) Knit plaques. 141 4-pt Bending 20 Laser displacement [mm] 15 B67 B68 B69 10 B70 B71 B72 B73 5 0 0 5 10 15 Machine displacement [mm] Figure A.19 4-pt bending test calibration curve. 142 In-Plane Shear Figure A.20 Sampling from standard plaques and specimen dimension for shear test. Table A.11 The dimensions of V-Notched specimens from standard plaques. ID S1 S2 S3 S4 S14 S15 Thickness, h 2.722 2.853 2.976 2.972 2.941 2.947 Width, d1 30.77 32.42 32.25 31.84 32.28 31.22 Table A.12 The dimensions of V-Notched specimens from standard plaques. ID S16 S17 S18 Thickness, h 3.176 3.33 3.369 Width, d1 29.08 29.94 34.41 Table A.13 The dimensions of V-Notched specimens from standard plaques. ID S19 S20 S21 S22 S23 S24 Thickness, h 3.65 3.587 3.437 3.437 3.148 3.128 Width, d1 31.56 28.63 30.54 30.62 30.63 30.69 143 0.05 S2 S3 S4 S14 S15 0.045 0.04 DIC Strain (mm/mm) 0.035 0.03 0.025 0.02 0.015 0.01 0.005 0 0 0.5 1 1.5 2 2.5 Machine displacement in mm Figure A.21 Correlation between DIC strain and machine displacement during shear test. Figure A.22 Failure mode of SMC standard plaques during in-plane shear test. 144 140 S19 S20 S21 S22 S23 S24 120 Shear Stress (MPa) 100 80 60 40 20 0 0 0.01 0.02 0.03 0.04 0.05 Shear strain Figure A.23 Shear stress-strain graphs for flow plaques. 35 S16 30 S17 Shear stress (MPa) 25 S18 20 15 10 5 0 0 0.005 0.01 0.015 0.02 Shear strain Figure A.24 Shear stress-strain graphs for knit Plaques. 145 120 Standard Plaque 100 Knit Plaque Shear stress (MPa) 80 Flow Plaque Longitudina Direction 60 40 20 0 0 0.01 0.02 0.03 0.04 0.05 0.06 Shear strain Figure A.25 Shear test results of glass fibers SMC. 0.12 0.1 Applied Stress [GPa] 0.08 0.06 0.4861 y = 2.1737x 0.04 0.02 0 0 0.0005 0.001 0.0015 0.002 0.0025 Plastic deformation Figure A.26 Plasticity parameters estimation from cyclic shear test. 146 Summary of Mechanical Properties Table A.14 Mechanical properties of SMC Flow Plaques in longitudinal direction. Tensile Compressive Flexure Shear Strength (MPa) 159.2 ± 28.0 184.9 ± 15.3 297.7 ± 70.5 106.8 ± 8.5 Modulus (GPa) 18.7 ± 2.3 14.6 ± 2.9 13.3 ± 1.8 5.1 ± 0.5 Strain at peak load 0.014 ± 0.002 0.018 ± 0.005 0.028 ± 0.003 0.033 ± 0.004 Table A.15 Mechanical properties of SMC flow plaques in transverse direction. Tension Compressive Flexure Strength (MPa) 87.4 ± 8.0 137.1 ± 29.6 176.7 ± 33.4 Modulus (GPa) 11.6 ± 1.7 11.1 ± 3.0 9.5 ± 1.1 Strain at peak load 0.011 ± .001 0.0174 ± 0.004 0.027 ± 0.006 Table A.16 Mechanical properties of SMC knit plaques. Tension Compressive Flexure Shear Strength (MPa) 16.1 ± 4.5 69.3 ± 8.5 68.7 ± 23.2 30.3 ± 0.9 Modulus (GPa) 7.3 ± 0.9 7.4 ± 1.0 7.9 ± 0.5 3.8 ± 1.0 Strain at peak load 0.0035 ± 0.002 0.012 ± 0.001 0.011 ± 0.004 0.010 ± 0.003 Table A.17 Comparison of mechanical properties of three types of SMC plaques. Ratio to Standard Properties Standard Flow Trans Flow Knit Flow Trans Flow Knit Tensile Strength (MPa) 151.1 159.1 87.4 16.1 1.05 0.58 0.11 Compressive Strength (MPa) 182.0 184.9 137.1 69.3 1.02 0.75 0.38 Shear Strength (MPa) 105.5 106.8 - 30.3 1.01 - 0.29 Flexural Strength (MPa) 253.9 297.7 176.9 68.7 1.17 0.70 0.27 Tensile Modulus (GPa) 14.7 18.7 11.6 7.3 1.27 0.79 0.50 Compressive Modulus (GPa) 13.4 14.6 11.1 7.4 1.09 0.83 0.55 Shear Modulus (GPa) 4.1 5.1 - 3.8 1.24 - 0.91 Flexural Modulus (GPa) 11.3 13.3 9.5 7.9 1.18 0.84 0.70 147 Model Validation Figure A.27 Experiment setup for open hole tension test. Nose Specimen Fixture Figure A.28 Experiment setup for disk bending with simply supported boundary. 148 APPENDIX B PYTHON SCRIPTS Python Script for Bimodal Weibull Strength Distribution This following python script generates 𝑛𝑚𝑎𝑡 = 40 material cart based on bimodal Weibull dis- tribution as described in in Chapter 5, Section 5.2 at page 70. In this script 𝐻 = [0.28, 0.29, 0.2] was used. ””” @author : S a k i b I q b a l comment : s t r e n g t h s f o l l o w B i m o d a l W e i b u l l d i s t r i b u t i o n . P=p∗P ( X_T )+(1 − p )∗P ( X_3P ) ””” import numpy a s np import m a t p l o t l i b . p y p l o t a s p l t import p a n d a s a s pd # ###################################################### # ##### h y p e r p a r a m e t e r s f o r m a t e r i a l m o d e l s ############ # ###################################################### n = 40 # number o f s e g m e n t o f d i s t r i b u t i o n ( N_mat ) p , EPSMAX, r e s 1 t = 0 . 2 8 , 0 . 2 9 , 0 . 2 # ###################### # E x p e r i m e n t r e s u l t s ## # ###################### # T e n s i l e s t r e n g t h ( same f o r d i r e c t i o n 1 o r 2 ) 149 X1T = [ 1 0 6 , 1 1 3 , 1 2 0 , 1 2 7 , 1 2 9 , 1 3 2 , 1 3 5 , 1 3 9 , 1 4 0 , 1 4 1 , 150 , 151 , 153 , 160] #3 p t F l e x u r e s t r e n g t h ( same f o r d i r e c t i o n 1 o r 2 ) X3F = [ 1 8 5 , 2 2 4 , 2 3 5 , 2 5 0 , 2 5 0 , 2 5 1 , 2 5 3 , 2 5 6 , 2 5 7 , 2 8 1 , 289 , 298 , 304 , 346] # C o m p r e s s i v e s t r e n g t h ( same f o r d i r e c t i o n 1 o r 2 ) X1C= [ 1 7 3 , 1 7 0 , 1 5 3 , 1 4 7 , 2 0 1 , 1 8 9 , 1 5 6 , 1 5 7 , 1 5 0 , 1 5 6 , 176 , 154 , 191 , 164 , 176 , 199 , 193 , 188 , 231 , 224 , 232 , 212 , 191 , 196] ## C o m p r e s s i v e s t r e n g t h d u r i n g 3 p t f l e x u r e can ’ t be d e t e r m i n e d from e x p e r i m e n t i t w i l l be d e t e r m i n e d from w e i b u l l model o f s i z e e f f e c t # T e n s i o n E l a s t i c m o d u l u s ( same f o r d i r e c t i o n 1 o r 2 ) E1T = [ 8 . 4 8 , 8 . 9 0 , 9 . 1 3 , 9 . 6 1 , 9 . 6 2 , 9 . 6 4 , 9 . 9 7 , 1 0 . 1 7 , 1 0 . 2 7 , 1 0 . 4 8 , 10.75 ,10.80 ,11.05 ,13.14] # C o m p r e s s i o n ( same f o r d i r e c t i o n 1 o r 2 ) E1C = [ 7 . 0 4 , 1 3 . 9 4 3 , 1 0 . 6 1 7 , 9.837 , 8.546 , 12.135 , 12.176 , 10.647 , 12.535 , 9.374 , 13.445 , 11.031 , 11.195 , 12.735 , 10.703 , 8.868 , 13.662 , 13.899 , 11.239 , 12.01 , 12.98 , 13.31 , 13.98] # ################################################################ # umat47v − New ECDM I n p u t Order − S a k i b I q b a l 08−24−2021 # # ################################################################ 150 # ### a i s an a r r a y f o r t h e m a t e r i a l c a r t ##### # ### a a r r a y d o e s n ’ t i n c l u d e t e x t s . t e x t s w i l l be added when writing output a=np . z e r o s ( ( 1 0 , 8 ) , d t y p e = f l o a t ) ## m a t e r i a l p a r a m e t e r a r r a y # f i b e r a n g l e s w i l l be d e t e r m i n e d f r o m ∗PART_COMPOSITE # F i r s t two l i n e s o f m a t e r i a l c a r t h a s t o be i n t e g e r . They w i l l be assigned l a s t section a [ 2 , 0 ] = 3 #AOPT a [ 3 , 0 ] = 1 #V1 a [ 3 , 1 ] = 0 #V2 a [ 3 , 2 ] = 0 #V3 a [ 4 , 4 ] = 4 . 2 4 #G12 a [ 4 , 5 ] = 2 . 0 3 #G23 a [ 4 , 6 ] = 2 . 0 3 #G13 a [ 4 , 7 ] = 1 #SOFT a [ 5 , 4 ] = 0 . 1 0 7 #XS a [ 5 , 5 ] = a [ 5 , 6 ] = a [ 5 , 7 ] = 0 . 3 3 # pr a [ 6 , 0 ] = a [ 6 , 1 ] = a [ 6 , 2 ] = a [ 6 , 3 ] = a [ 6 , 4 ] = 0 . 0 3 #E1TPK , E1TCPK , E2TPK , E2CPK , ESPK a [ 6 , 5 ] = a [ 6 , 6 ] = 4 0 # pm1t , pm1c a [ 7 , 4 ] = 0 . 0 2 2 #RESS a [ 7 , 5 ] = a [ 7 , 6 ] = a [ 7 , 7 ] = 4 0 # pm2t , pm2c , pms #RES1T and RES1T w i l l be d e t e r m i n e d l a t e r , a s a p e r c e n t a g e o f 151 ultimate strength a [ 8 , 1 ] = 2 #EROD a [ 8 , 2 ] = a [ 8 , 3 ] = a [ 8 , 4 ] = a [ 8 , 5 ] =EPSMAX #EPSMAX , e 1 d e l e , e 2 d e l e , e 4 d e l e a [ 8 , 6 ] = 6 0 #K a [ 8 , 7 ] = 1 5 #G # ############################################################ # W e i b u l l model f o r T e n s i l e , c o m p r e s s i v e s t r e n g t h & m o d u l u s # # ############################################################ # ## W e i b u l l model p a r a m e t e r d e v e l o p m e n t def weibull ( x ) : x= s o r t e d ( x ) s n =np . a r a n g e ( l e n ( x ) ) + 1 S=1 −( sn − 0 . 3 ) / ( l e n ( s n ) + 0 . 4 ) v a r =np . p o l y f i t ( np . l o g ( x ) , np . l o g ( −1∗ np . l o g ( S ) ) , 1 ) m= v a r [ 0 ] x0=np . e ∗∗( −1∗ v a r [ 1 ] / v a r [ 0 ] ) r e t u r n m, x0 # D e t e r m i n e W e i b u l l model p a r a m e t e r s m, x0 f r o m m a t e r i a l properties # # t h e s e p a r a m e t e r w i l l be w r i t t e n t o m and x0 l i s t # # P a r a m e t e r a r e w r i t t e n i n o r d e r i . e . X1T , X3F , X1C , X3FC m= [ ] x0 = [ ] m. a p p e n d ( w e i b u l l ( X1T ) [ 0 ] ) ; x0 . a p p e n d ( w e i b u l l ( X1T ) [ 1 ] ) 152 m. a p p e n d ( w e i b u l l ( X3F ) [ 0 ] ) ; x0 . a p p e n d ( w e i b u l l ( X3F ) [ 1 ] ) m. a p p e n d ( w e i b u l l ( X1C ) [ 0 ] ) ; x0 . a p p e n d ( w e i b u l l ( X1C ) [ 1 ] ) # Determine compressive s t r e n g t h f o r 3 pt f l e x u r e t e s t # e q u a t i o n 0 . 5 ∗w∗ l ∗ t / ( m+1)^2 t h r e e p t v e = 0 . 5 ∗ 1 2 . 5 ∗ 4 8 ∗ 3 . 2 4 / (m[ 2 ] + 1 ) ∗ ∗ 2 # e q u a t i o n w∗ l ∗ t compressionn_ve =20∗3.305∗100 X1C=np . a r r a y ( X1C ) X3FC=X1C∗ ( c o m p r e s s i o n n _ v e / t h r e e p t v e ) ∗ ∗ ( 1 /m[ 2 ] ) m. a p p e n d ( w e i b u l l ( X3FC ) [ 0 ] ) ; x0 . a p p e n d ( w e i b u l l ( X3FC ) [ 1 ] ) # t e s n i o n m o d u l u s o f e l a s t i c i t y w e i b u l l model m. a p p e n d ( w e i b u l l ( E1T ) [ 0 ] ) x0 . a p p e n d ( w e i b u l l ( E1T ) [ 1 ] ) # c o m p r e s s i o n m o d u l u s o f e l a s t i c i t y w e i b u l l model m. a p p e n d ( w e i b u l l ( E1C ) [ 0 ] ) x0 . a p p e n d ( w e i b u l l ( E1C ) [ 1 ] ) # D e v e l o p b i m o d a l d i s t r i b u t i o n f r o m W e i b u l l model d i s t r i b u t i o n parameters # # f u n c t i o n t o d e t e r m i n e s t r e n g t h s from bimodal w e i b u l l d i s t r i b u t i o n d e f b i m o d a l (m, x0 ) : xlow=x0 [ 0 ] ∗ ( − 1 ∗ np . l o g ( 0 . 9 6 ) ) ∗ ∗ ( 1 /m[ 0 ] ) xm=x0 [ 1 ] ∗ ( − 1 ∗ np . l o g ( 0 . 0 1 ) ) ∗ ∗ ( 1 /m[ 1 ] ) X=np . l i n s p a c e ( xlow , xm+1 , n + 1 ) S=p∗ np . e ∗ ∗ ( − 1 ∗ (X/ x0 [ 0 ] ) ∗ ∗m[ 0 ] ) + 153 (1 − p ) ∗ np . e ∗ ∗ ( − 1 ∗ (X/ x0 [ 1 ] ) ∗ ∗m[ 1 ] ) x = (X[ : − 1 ] +X [ 1 : ] ) / 2 x=np . round ( x ) return x , S stren_mod =[] # s t r e n g t h / modulus d i s t r i b u t i o n f o r i i n range ( 2 ) : s t r e n _ m o d . a p p e n d ( b i m o d a l (m[ 2 ∗ i : 2 ∗ i + 2 ] , x0 [ 2 ∗ i : 2 ∗ i + 2 ] ) [ 0 ] ) S= b i m o d a l (m[ 0 : 2 ] , x0 [ 0 : 2 ] ) [ 1 ] # E l a s t i c modulus f o l l o w s unimodal W e i b u l l d i s t r i b u t i o n # Determine e l a s t i c modulus d i s t r i b u t i o n from W e i b u l l paramters # m o d u l u s random d i s t r i b u t i o n s t r e n _ m o d . a p p e n d ( np . random . w e i b u l l (m[ 4 ] , n ) ∗ x0 [ 4 ] ) s t r e n _ m o d . a p p e n d ( np . random . w e i b u l l (m[ 5 ] , n ) ∗ x0 [ 5 ] ) s t r e n _ m o d =np . a r r a y ( s t r e n _ m o d ) ## d e t e r m i n e p r o b a b i l i t y f r o m S ( w e i b u l l s u r v i v a l p r o b a b i l i t y ) P = −1∗ np . d i f f ( S ) # ############################################################# # Determine r e s i d u a l s t r e n g t h f o r t e n s i l e , compression Failure # # ############################################################# ## r e s i d u a l c o m p r e s s i v e s t r e n g t h would be 40% o f u l t i m a t e s t r e n g t h 154 RES1T=np . round ( s t r e n _ m o d [ 0 ] ∗ r e s 1 t / 1 0 0 0 , 3 ) RES1C=np . round ( s t r e n _ m o d [ 1 ] ∗ 0 . 0 0 0 5 , 3 ) # ########################################################### # D e t e r m i n e damage p a r a m e t e r f r o m s t r e n g t h & e l a s t i c m o d u l u s # # ########################################################### # i t i s h y p o t h e s i z e d t h a t n o l i n e a r i t y s t a r t e d a t 30% o f u l t i m a t e strength # # ## Ec ( e l a s t i c m o d u l u s a t u l t i m a t e s t r e s s ) i s 40% o f e l a s t i c modulus ( t h i s i s d e t e r m i n e d from e x p e r i m e n t ) # # Damage p a r a m e t e r s f o r t e n s i o n and c o m p r e s s i o n a r e d e t e r m i n e d h e r e . Damage p a r a m e t e r s f o r s h e a r # # and o t h e r p a r a m e t e r s w i l l be i n L a d e v e z e . t x t file # # damage p a r a m e t e r d e t e r m i n a t i o n d e f damage (X, E ) : y10 = 0 . 3 ∗X/ ( 1 0 0 0 ∗ ( 2 ∗ E ) ∗ ∗ 0 . 5 ) Ec = 0 . 4 ∗ E d=1−Ec / E y1c =X/ ( 1 0 0 0 ∗ ( 1 − d ) ∗ ( 2 ∗ Ec ) ∗ ∗ 0 . 5 ) r e t u r n np . r o u n d _ ( y10 , 3 ) , np . r o u n d _ ( y1c , 3 ) y10 , y1c =damage ( s t r e n _ m o d [ 0 ] , s t r e n _ m o d [ 2 ] ) y10c , y 1 c c =damage ( s t r e n _ m o d [ 1 ] , s t r e n _ m o d [ 3 ] ) # ######################################################## # Write output # # ######################################################## b=pd . DataFrame ( [ 155 [ ’ $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ’ ] , [ ’ $ ###### umat47v −NewECDMInputOrder − S a k i b _ I q b a l ###### $ ’ ] , [ ’ $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ’ ] , [ ’ ∗MAT_USER_DEFINED_MATERIAL_MODELS ’ ] ] ) c=pd . DataFrame ( [ ’ ∗END ’ ] ) f o r i i n range ( n ) : a [ 4 , 0 ] = a [ 4 , 2 ] = np . r o u n d _ ( s t r e n _ m o d [ 2 , i ] , 1 ) a [ 4 , 1 ] = a [ 4 , 3 ] = np . r o u n d _ ( s t r e n _ m o d [ 3 , i ] , 1 ) a [ 5 , 0 ] = a [ 5 , 2 ] = np . r o u n d _ ( s t r e n _ m o d [ 0 , i ] / 1 0 0 0 , 3 ) a [ 5 , 1 ] = a [ 5 , 3 ] = np . r o u n d _ ( s t r e n _ m o d [ 1 , i ] / 1 0 0 0 , 3 ) a [ 7 , 0 ] = a [ 7 , 2 ] = RES1T [ i ] #RES1T , RES2T a [ 7 , 1 ] = a [ 7 , 3 ] = RES1C [ i ] a [ 9 , 0 ] = a [ 9 , 4 ] = y10 [ i ] a [ 9 , 1 ] = a [ 9 , 5 ] = y1c [ i ] a [ 9 , 2 ] = a [ 9 , 6 ] = y10c [ i ] a [ 9 , 3 ] = a [ 9 , 7 ] = y1cc [ i ] # conver a matrix to dataframe mat=pd . DataFrame ( a , d t y p e = ’ o b j e c t ’ ) # w r i t e down f i r s t two rows o f m a t e r i a l c a r t mat . i l o c [ 0 , 0 ] = 5 0 2 0 1 + i # m a t e r i a l ID mat . i l o c [ 0 , 1 ] = 1 . 7 2 e −6 # r o mat . i l o c [ 0 , 2 ] = 4 7 #MT mat . i l o c [ 0 , 3 ] = 4 8 #LMC mat . i l o c [ 0 , 4 ] = 7 0 #NHV mat . i l o c [ 0 , 5 ] = 1 #IORTHO mat . i l o c [ 0 , 6 ] = 3 9 #IBULK 156 mat . i l o c [ 0 , 7 ] = 4 0 #IG mat . i l o c [ 1 , 0 ] = 1 #IVECT mat . i l o c [ 1 , 1 ] = 0 # IFAIL mat [ 2 ] [ 1 ] = mat [ 3 ] [ 1 ] = mat [ 4 ] [ 1 ] = mat [ 5 ] [ 1 ] = mat [ 6 ] [ 1 ] = mat [ 7 ] [ 1 ] = 0 # output to csv f i l e b . t o _ c s v ( ’ mat /ECDM%s . k ’%( i + 1 ) , i n d e x = F a l s e , h e a d e r = F a l s e ) # append mat a r r a y t o CSV f i l e mat . t o _ c s v ( ’ mat /ECDM%s . k ’%( i + 1 ) , mode= ’ a ’ , i n d e x = F a l s e , header=False ) # append CSV f i l e c . t o _ c s v ( ’ mat /ECDM%s . k ’%( i + 1 ) , mode= ’ a ’ , i n d e x = F a l s e , header=False ) # ## w r i t e h y p e r p a r a m e t e r s t o a f i l e c=np . z e r o s ( ( 2 , n ) ) c [0]=P ; c [1 ,0]= n np . s a v e t x t ( ’ mat / p a r a m e t e r ’ , c , d e l i m i t e r = ’ , ’ , f m t = ’ %1.4 f ’ ) Python Script for Randomization The following python script creates material.k file that contains material cart distribution as- signed to various *PART_COMPOSITE. Randomization was done with 𝑛𝑚𝑎𝑡 = 40, *PART_COMPOSITE = 50. Each *PART_COMPOSITE was configured with 8 integration point. In addition to material cart number, material angle 𝛽 was randomly varied between -90° to 90°. import random import p a n d a s a s pd from math import c e i l 157 a =50 # number o f p a r t s b=8 # number o f i n t e g r a t i o n p o i n t d=round ( 3 . 2 5 / b , 3 ) # t h i c k n e s s o f e a c h i n t e g r a t i o n p o i n t r = c e i l ( a / 8 ) + 1 # s i z e o f ∗SET_PART_LIST a r r a y # ####################################################### # ### P , n ( N_mat ) a r e w r i t t e n by m o d i f i e d w e i b u l l . py #### # ######## Read p a r a m e t e r . t x t f o r P , n ( N_mat ) ########## # ####################################################### x = pd . r e a d _ c s v ( ’ mat / p a r a m e t e r ’ , h e a d e r =None ) . v a l u e s n = int (x [1 ,0]) P = ( ( x [ 0 ] ) ∗ a ∗b ∗ 1 . 0 4 ) . round ( ) . a s t y p e ( i n t ) # c r e a t e random l i s t o f m a t e r i a l m o d e l s p =[] f o r i i n range ( n ) : p . append ([50201+ i ]∗ P [ i ] ) p=sum ( p , [ ] ) random . s h u f f l e ( p ) filename= ’ material . k ’ ## C r e a t e m a t e r i a l . k f i l e w i t h t h i s i n t r o d u c t i o n g=pd . DataFrame ( [ [ ’ $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ’ ] , [ ’ $$ #### RandomAssingmentofMATmodelin ∗PART_COMPOSITE##### $$ ’ ] , [ ’ $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ ’ ] ] ) 158 g . to_csv ( filename , index=False , header= False ) ## f i s t h e d a t a f r a m e f o r d e f i n a t i o n o f ∗PART_COMPOSITE f =pd . DataFrame ( c o l u m n s = range ( 8 ) , i n d e x = range ( b / / 2 + 1 ) ) f o r i i n range ( a ) : f . i l o c [0 ,0]=100+ i f . i l o c [0 ,1]=16 f . i l o c [0 ,2]=1 e=pd . DataFrame ( [ [ ’ ∗PART_COMPOSITE ’ ] , [ ’ P a r t ␣%s ’ %(100+ i ) ] ] ) f o r j i n range ( b / / 2 ) : f . i l o c [ j + 1 , 0 ] = p [ i ∗b+ j ∗ 2 ] f . i l o c [ j + 1 , 4 ] = p [ i ∗b+ j ∗2+1] f . i l o c [ j +1 ,1]= f . i l o c [ j +1 ,5]= d f . i l o c [ j + 1 , 2 ] = random . r a n d r a n g e ( − 9 0 , 9 1 , 1 5 ) f . i l o c [ j + 1 , 6 ] = random . r a n d r a n g e ( − 9 0 , 9 1 , 1 5 ) e . t o _ c s v ( f i l e n a m e , mode= ’ a ’ , i n d e x = F a l s e , h e a d e r = F a l s e ) f . t o _ c s v ( f i l e n a m e , mode= ’ a ’ , i n d e x = F a l s e , h e a d e r = F a l s e ) # ############################################### # ############# C r e a t e s e t _ p a r t _ l i s t ############ # ############################################### f =pd . DataFrame ( c o l u m n s = range ( 8 ) , i n d e x = range ( r ) ) e=pd . DataFrame ( [ ’ ∗SET_PART_LIST ’ ] ) f . i l o c [0 ,0]=1 f . i l o c [ 0 , 5 ] = ”MECH” 159 f o r i i n range ( 1 , r ) : f o r j i n range ( 8 ) : f . i l o c [ i , j ] = 1 0 0 + 8 ∗ ( i −1)+ j i f 100+8∗( i −1)+ j ==100+ a − 1 : break e . t o _ c s v ( f i l e n a m e , mode= ’ a ’ , i n d e x = F a l s e , h e a d e r = F a l s e ) f . t o _ c s v ( f i l e n a m e , mode= ’ a ’ , i n d e x = F a l s e , h e a d e r = F a l s e ) # ############################################### # ######## add ∗INCLUDE m a t e r i a l c a r t f i l e ###### # ########### i n t o m a t e r i a l . k f i l e ############## a = pd . DataFrame ( [ ’ ∗INCLUDE ’ ] ) f o r i i n range ( n ) : b = ’ mat /ECDM%s . k ’%( i + 1 ) a = a . append ( [ b ] , i g n o r e _ i n d e x =True ) a . t o _ c s v ( f i l e n a m e , mode= ’ a ’ , i n d e x = F a l s e , h e a d e r = F a l s e ) 160 APPENDIX C PUBLICATIONS 1. Sakib Iqbal, Beichen Li, Kestutis Sonta, Royal Chibuzor Ihuaenyi, and Xinran Xiao. Prob- abilistic finite element analysis of sheet molding compound composites with an extended strength distribution model. Finite Elements in Analysis and Design, 214:103865, 2023 2. Sakib Iqbal and Xinran Xiao. Optimization of fe crashworthiness model for sheet mold- ing compound (smc) with extended strength distribution model and machine learning. In Proceedings of The American Society For Composites-Thirty-Seventh Technical Conference, 2022 3. Royal Chibuzor Ihuaenyi, Shutian Yan, Jie Deng, Chulheung Bae, Iqbal Sakib, and Xinran Xiao. Orthotropic thermo-viscoelastic model for polymeric battery separators with elec- trolyte effect. Journal of The Electrochemical Society, 168(9):090536, 2021 4. Sakib Iqbal, Beichen Li, Kestutis Sonta, Royal Chibuzor Ihuaenyi, and Xinran Xiao. Simu- lations of sheet molding compound composite structures with a randomly distributed tensile strength. 2021 5. Sakib Iqbal, Xinran Xiao, Beichen Li, and Kestutis Sonta. Considering the randomness of mechanical properties in simulations of discontinuous fiber reinforced composites. In Proceedings of the American Society for Composites Thirty-fifth Technical Conference, 2020 161