E F F E C T O F P R O C E S S I N G O N T H E B E H A V I O U R O F L A M I N A T E D C O M P O S I T E S T R U C T U R E S : A N U M E R I C A L A N D P R O B A B I L I S T I C A P P R O A C H by A B D U L R A H I M A H A M E D A R A F A T H B.Sc. (Civil Engineering), University of Peradeniya, Sri Lanka, 1998 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF M A S T E R OF APPLIED SCIENCE in THE F A C U L T Y OF G R A D U A T E STUDIES (Department of Civil Engineering) We accept this thesis as conforming to the required standard THE UNIVERSITY OF BRITISH C O L U M B I A February 2002 © A.R.A. Arafath, 2002 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of CL.i\/'l ^ w g g - r / ^ ^ The University of British Columbia Vancouver, Canada Date ftb j f U o-oo DE-6 (2/88) A B S T R A C T The manufacturing of composite structures made of fibre-reinforced thermoset matrix materials generally results in residual stress build-up during the curing process. As a consequence, the precise geometry and material properties of the final structure after tool (mould) removal are often difficult to control. To account for such process-induced variations, engineers either have employed empirical knockdown factors, which are determined experimentally for a range of distinct structures, materials and manufacturing processes or have incorporated the measured variations as input into sophisticated finite element analyses. The first approach may lead to an unconservative design while the second approach is not suitable for normal design practice. In this study, a combined deterministic-stochastic simulation approach is presented to demonstrate the manner in which manufacturing-induced variabilities in composite structures can be controlled to achieve targeted reliability levels in structural performance. A finite element code, COMPRO, which deterministically simulates the various physical phenomena during manufacturing of composite structures, is integrated with non-linear structural analysis and reliability analysis to compute the statistics of the parameters that must be controlled at the manufacturing level in order to result in an optimum or reliable structural performance during service. The methodology is demonstrated through a case study that examines the buckling behaviour of a composite plate in the presence of manufacturing-induced imperfections. The intrinsic and extrinsic process parameters that affect the dimensional stability of more complex composite structures such as L- and C-shaped angular laminate sections, are numerically investigated using COMPRO. The predictive capability of COMPRO is examined by comparing the simulation results with available experimental measurements. Finally, to account for all sources of residual stress build-up required for accurate assessment of failure, the micro-level stresses arising from the mismatch between the coefficient of thermal expansions of the fibre and the matrix and the matrix cure shrinkage are computed through micromechanical finite element analysis of a representative volume of the fibre and the matrix. The ply- and laminate-level residual stresses computed by COMPRO during the curing process are then transferred to the micro-level using suitably derived stress magnification factors. - i i i -T A B L E OF CONTENTS Abstract ii Table of Contents iv List of Tables vii List of Figures viii Acknowledgements xiii Chapter 1: Introduction 1 1.1 Background 1 1.2 Software Tools Used in This Study 4 1.2.1 Composite Process Modelling Software - COMPRO 4 1.2.2 Inverse Reliability Analysis Software - I R E L A N 9 1.3 Research Objective and Scope of Thesis r 15 1.4 Figures 17 Chapter 2: Probabilistic Design of a Composite Plate Against Buckling 20 2.1 Introduction 20 2.2 Problem Definition 24 2.3 Process Modelling 25 2.4 Structural Modelling 25 2.5 Prediction of the 3D State of Stress 26 2.5.1 Plane Strain 27 2.5.2 Three Dimensional Prediction 27 2.6 Parametric Study 33 2.7 Reliability-Based Determination of Process Control Parameters 34 2.8 Other Sources of Variabilities 37 2.9 Summary and Conclusions 39 2.10 Tables 40 2.11 Figures 44 Chapter 3: Spring-in of Angled Composite Laminate 55 3.1 Introduction 55 -iv-3.2 Re view of Experi merits 57 3.3 Process Modelling using COMPRO 58 3.4 Sources of residual stresses 59 3.4.1 Anisotropy 59 3.4.2 Tool-Part Interaction 64 3.5 Effect of Intrinsic and Extrinsic Parameters 68 3.5.1 Effect of Part Thickness 69 3.5.2 Effect of Tool Material 69 3.5.3 Effect of Lay-up 70 3.5.4 Effect of Part Shape 72 3.6 Summary and Conclusions 74 3.7 Tables 76 3.8 Figures , 81 Chapter 4: Process-Induced Residual Stresses 100 4.1 Introduction 100 4.2 Micromechanical Finite Element Modelling 101 4.3 Boundary Conditions 103 4.4 Material Modelling in A B A Q U S 103 4.5 Stress Analysis 106 4.5.1 Thermal Stresses 107 4.5.2 Cure Shrinkage Stresses 108 4.6 Concentric Cylindrical Model 111 4.7 Transferring macro level stresses to micro level 112 4.8 Prediction of Engineering Properties Using Micromechanics 113 4.9 Summary and Conclusions 115 4.10 Tables 116 4.11 Figures 117 Chapter 5: Conclusions and Future Work 127 5.1 Future Work 128 References 130 Appendix A: Material Properties of Hercules AS4/8552 Unidirectional Prepreg 133 -v-A . 1 Fibre Elastic Properties 133 A.2 Resin Cure Kinetics Model 133 A.3 Resin Cure Shrinkage Model 134 A.4 Resin Modulus Development Model 135 A.5 Resin Poisson's Ratio Development Model 135 A. 6 Thermal Expansion Coefficient 136 Appendix B: Ply Micromechanical Equations 137 B. l Calculation of Elastic Constants 137 B . l . l In-Plane Moduli 137 B.1.2 Shear Moduli 138 B.l .3 Poisson's ratios 138 B.2 Calculation of Thermal and Cure Shrinkage Strains 139 Appendix C: Residual Thermal Stresses in Two Concentric Transversely Isotropic Cylinders 140 Appendix D: Tool-Part Interaction 142 D . l Bi-Metallic Strip 142 D.2 Tool-Part Interaction During Heating Up 142 D.3 Tool-Part Interaction During Cooling Down 143 D.4 Tool Removal 144 D.5 Figures 145 - v i -LIST OF T A B L E S Table 2.1: Laminate material properties 40 Table 2.2: Database 40 Table 2.3: Statistics of random variables 43 Table 2.4: Inverse reliability and optimum design parameters, E X A M P L E 1 43 Table 2.5: Inverse reliability and optimum design parameters, E X A M P L E 2 43 Table 3.1: Evaluated laminate equivalent properties (the values in [ ] are those obtained using CLT) 76 Table 3.2: Spring-in of the L-shaped composite parts when constrained by the tool 76 Table 3.3: Shear layer properties and their identification labels 76 Table 3.4: Selected shear layer properties to represent the experimental tool surface conditions 77 Table 3.5: Spring-in variation with tool surface condition. COMPRO predictions are compared with experimental results 77 Table 3.6: Spring-in variation with part thickness. COMPRO predictions are compared with experimental results 78 Table 3.7: Spring- in variation with tool material for an 8-layer unidirectional part. COMPRO predictions are compared with experimental results 78 Table 3.8: Spring-in variation with tool material for a 16-layer unidirectional part. COMPRO predictions are compared with experimental results 79 Table 3.9: Variation of maximum warpage of a flat part (600mm long, 16 plies thick) with lay-up 79 Table 3.10: Spring-in variation with part lay-up. COMPRO predictions are compared with experimental results 79 Table 3.11: Spring-in variation with part shape for aluminum tool. COMPRO predictions are compared with experimental results 80 Table 3.12: Spring-in variation with part shape for steel tool. COMPRO predictions are compared with experimental results 80 Table 4.1: Material properties of isotropic fibre and matrix 116 Table 4.2: Material properties of composite for fibre volume fraction of 0.573. Both fibre and matrix are isotropic 116 Table 4.3: Material properties for an AS4/8552 lamina. Material properties are calculated at both initial (uncured) and final (cured) resin modulus 116 -vii-LIST OF FIGURES Figure 1.1: Schematic of COMPRO structure and program flow (Johnston, 1997) 17 Figure 1.2: Four-noded and eight-noded isoparametric quadrilateral element employed in COMPRO 17 Figure 1.3: Composite element (a) isometric view (b) 2D view (orientation of element local axis is also shown) (Johnston, 1997) 18 Figure 1.4: Curved region, showing element local coordinates, local axes orientation and the region reference boundary (Johnston, 1997) 18 Figure 1.5: Geometric representation of reliability calculation for F O R M (Li, 1999) 19 Figure 2.1: Basic structural elements (shell, plate and stiffener) of an aircraft 44 Figure 2.2: The air temperature and pressure cycle as applied in the autoclave during the processing of the composite plate (Twigg, 2001) 44 Figure 2.3: Non-dimensional axial load versus transverse deflection curve of the plate with various initial imperfections. Pcr is the critical Euler buckling load of the plate 45 Figure 2.4: Finite element mesh of the plate and the tool used in COMPRO to virtually manufacture the part. Shear layer is used to model the tool - part interaction 45 Figure 2.5: The warped profile of the part as it is removed from the tool. The prediction of COMPRO is compared with experimental results 46 Figure 2.6: Schematic of the composite plate used for the case study. The applied loading and boundary conditions are also shown 46 Figure 2.7: The superposition method to obtain the 3D state of stress from plane strain results for applied mechanical load (a) load P is applied in all three directions on a cube (b) 3D element: free to expand in all three directions (c) 2D plane strain element: load is applied in 1 and 2 directions, free to expand only in 1 and 2 direction, displacement in the 3rd direction is constrained to be zero (d) load is applied only in the 3rd direction, free expand in all three directions. Combination of (c) and (d) gives (b) 47 Figure 2.8: The difference between 3D element and plane strain element (a) temperature increment is applied to a cube (b) 3D element: free to expand in all three directions, no stresses developed (c) 2D plane strain element: free to expand only in 1 and 2 directions, displacement in the 3rd direction is constrained to be zero , 48 Figure 2.9: Superposing the plane strain stresses to obtain the 3D results (a) stresses from the plane strain analysis are applied as initial stresses to the cube (b) the -viii-cube is allowed to equilibrate under these initial stresses with free boundary conditions (c) the cube problem is modeled with a number of plane strain elements (d) the plane strain block in (c) is combined to form the full cube and allowed to equilibrate under the stresses 49 Figure 2.10: Plate simply supported at the corners under thermal load (finite element mesh is shown for one quarter of the plate) 50 Figure 2.11: Developing the full 3D results using the results of a row of shell elements (a) a row of shell elements (b) full 3D plate 50 Figure 2.12: Comparison of stresses from the actual 3D analysis and the modified superposition method (a) stress <TX (b) stress ay (c) stress 51 Figure 2.13: Comparison of vertical displacement from the actual 3D analysis and the modified superposition method 52 Figure 2.14: Mapping the stresses from plane strain solid element to shell elements, (a) The plate problem is modeled with plane strain solid elements and (b) the stresses from the plane strain analysis are mapped to the shell elements 52 Figure 2.15: Comparison of stresses from the actual 3D analysis and the modified superposition method (a) stress <JX (b) stress ay (c) stress 53 Figure 2.16: Comparison of vertical displacement from the actual 3D analysis and the modified superposition method 54 Figure 2.17: Sensitivity study of the processing parameters that mostly affect the warpage of the part. Sensitivity is defined as the ratio of the percentage increase in warpage to the percentage increase in the parameter 54 Figure 3.1: Tool geometry used in the experiment (Albert, 2001) 81 Figure 3.2: Part geometries used in the experiment (a) C-shaped (b) L-shaped (Albert, 2001) 81 Figure 3.3: The geometry of the L-shaped composite part and the tool used in modelling 82 Figure 3.4: Finite element representation of the angle laminate and the boundary conditions applied 82 Figure 3.5: Target cure cycles used in autoclave processing of the L-shaped part. Both one and two hold cure cycles (similar to the ones used by Albert, 2001) are used 83 Figure 3.6: Measurement of spring-in of the composite L-shape part 83 Figure 3.7: Predicted flange warpage profile of a 16-layer L-shaped part 84 Figure 3.8: Comparison of spring-in prediction between both plane strain and three-dimensional analysis 84 Figure 3.9: Numerical model to predict the equivalent properties of a composite laminate using a 3D elastic analysis, (a) Schematic diagram of a 16-layer composite element (b) & (c) identification of faces for boundary conditions 85 Figure 3.10: Temperature profile through part thickness 86 - i x -Figure 3.11: Part geometry for the analytical model. Sliding friction occurs at both the tool-part interface and at the interply region (Twigg, 2001) 86 Figure 3.12: Illustration of interfacial shear stress distribution as the debond front travels from the tool end towards the center (Twigg, 2001) 87 Figure 3.13: The variation of interfacial shear stress distribution with part thickness (600 mm long part with unidirectional layers) according to COMPRO prediction using a shear layer of modulus 6.9 kPa 87 Figure 3.14: The variation of interfacial shear stress distribution with part length (8 unidirectional layer part) according to COMPRO prediction using a shear layer of modulus 6.9 kPa 88 Figure 3.15: The variation of interfacial shear stress distribution with part length. By modifying the shear layer property with length ratio, the peak shear stresses are made equal (average shear stresses are also equal) 88 Figure 3.16: The variation of the warpage of a flat part with thickness. The prediction of COMPRO is compared with analytical solution proposed by Twigg (2001) 89 Figure 3.17: The variation of the warpage of a flat part with part length. The prediction of COMPRO is compared with analytical solution proposed by Twigg (2001) 89 Figure 3.18: Contribution by the CTE anisotropy and cure shrinkage to the total spring-in when the part is constrained by the tool 90 Figure 3.19: Spring-in of the parts with different shear layers. Only the Young's modulus is shown in here although the shear modulus is also changed by the same order 90 Figure 3.20: Comparison between the predicted and measured warpage components of an 8-layer part to select the appropriate shear layer property to represent the tool surface condition used in the experiment 91 Figure 3.21: Spring-in variation with tool surface condition. COMPRO predictions are compared with experimental results for both 8-layer and 16-layer unidirectional parts 91 Figure 3.22: Schematic showing the corner and flange components of the L-shaped part 92 Figure 3.23: Variation of spring-in with part thickness predicted by COMPRO. The results are shown for two different tool surfaces 92 Figure 3.24: Spring-in variation with part thickness. COMPRO results are compared with experimental results for two different tool surface conditions 93 Figure 3.25: Spring-in variation with tool material for an 8-layer unidirectional part. COMPRO results are compared with experimental results for two different tool surface conditions 93 Figure 3.26: Spring-in variation with tool material for a 16-layer unidirectional part. COMPRO results are compared with experimental results for two different tool surface conditions 94 Figure 3.27: The predicted variation of spring-in with part lay-up for two different thickness of the part. The results for three different shear layer properties are also shown. No-Shear indicates that there is minimal tool part interaction 94 Figure 3.28: Spring-in variation with part lay - up. COMPRO results are compared with experimental results for two different thickness of the part 95 Figure 3.29: The geometry of half of the C-shaped composite part and the tool 95 Figure 3.30: Variation of spring-in with part shapes predicted by COMPRO. The results are shown for three different tool surface conditions. No-shear indicates that there is minimal tool part interaction 96 Figure 3.31: The experimentally measured spring-in of C-shaped parts. In all cases, the warpage of the flange is negligible (Albert, 2001) 96 Figure 3.32: Schematic showing stretching of the C-shaped part due to corner-locking as the tool expands 97 Figure 3.33: Illustration of the three zones of the shear layer used in the modified method to overcome the corner-locking 97 Figure 3.34: Spring-in results by the modified method to overcome the corner-locking problem. C is the original results and C m is the modified results 98 Figure 3.35: Effect of part shape on spring-in. COMPRO predictions by both the original method (C) and the modified method (C m ) are compared with experimental results for aluminum tool 98 Figure 3.36: Effect of part shape on spring-in. COMPRO predictions by both the original method (C) and the modified method (C m ) are compared with experimental results for steel tool 99 Figure 4.1: The element selected for micro stress analysis in an L-shaped part 117 Figure 4.2: The cross section of a unidirectional ply with commonly used periodic fibre arrangements: (a) square array (b) hexagonal array 117 Figure 4.3: Cross section of graphite-reinforced composite material (Hyer, 1997) 118 Figure 4.4: The representative volume element (RVE) concept for (a) square array (b) hexagonal array 118 Figure 4.5: Finite element mesh used for the micro stress analysis: (a) square array (b) square array (fine mesh) (c) hexagonal array 119 Figure 4.6: The applied boundary condition for the R V E . Due to the symmetry only a quarter of the R V E is considered. "Sym" stands for symmetric boundary and "Coup" stands for the couple boundary 119 Figure 4.7: The temperature-time history of the element shown in Figure 4. 1 as obtained from COMPRO 120 Figure 4.8: The development of degree of cure of resin as it cures during the autoclave processing of the L-shaped composite part: the A B A Q U S results for three - X I -different time steps (10 min, 5 min and 2 min) are compared with COMPRO results. The circled portion of the graph is magnified and shown 120 Figure 4.9: The development of stresses due to CTE mismatch between the fibre and the matrix and cure shrinkage of resin as the resin cures during processing of an L-shaped composite part 121 Figure 4.10: The axial (7), normal (r), and the circumferential (0) stress components at the fibre matrix interface. The maximum stress occurs at the interface 121 Figure 4.11: The distribution of normal (ar) and shear (xre) stress components along the fibre-matrix interface due to the CTE mismatch and resin cure shrinkage 122 Figure 4.12: The distribution of circumferential (oe) stress components along the fibre-matrix interface due to the CTE mismatch and resin cure shrinkage 122 Figure 4.13: The distribution of the axial stress components, c i , along the fibre-matrix interface due to the CTE mismatch and resin cure shrinkage 123 Figure 4.14: Comparison of circumferential stresses between concentric cylinders model (elasticity approach) and finite element solution for a temperature change of -160 °C 123 Figure 4.15: Comparison of normal and shear stresses between concentric cylinders model (elasticity approach) and finite element solution for a temperature change of-160 °C 124 Figure 4.16: Comparison of axial stresses between concentric cylinders model (elasticity approach) and finite element solution for a temperature change of -160 °C 124 Figure 4.17: The variation of the stress multiplication factor for axial stress at the interface during the cure cycle. The multiplication factor is for the fibre at 9 = 0 (Figure 4. 10) for an applied axial mechanical stress 125 Figure 4.18: Comparison of transverse modulus, E2 by three different approaches for different matrix modulus. Both fibre and matrix are isotropic 125 Figure 4.19: Unit cells considered in: (a) elasticity approach (b) hexagonal array and (c) square array 126 Figure 4.20: A l l unit cells superposed on top of each other to reveal the difference between different approaches 126 Figure D . l : Bi-metallic strip under temperature change, AT 145 Figure D.2: First ply slip during the initial stage of curing process (heating up). The resin is not cured and it has very low modulus 145 Figure D.3: The normal stress distribution through part thickness during heating up 146 Figure D.4: Tool-part interaction at the final stage of the curing process (cooling down). The resin is fully cured and it has high modulus. The entire part is moved with the tool 146 Figure D.5: The normal stress distribution through part thickness during cooling down 147 Figure D.6: Resultant normal stress distribution through part thickness 147 -xn-A C K N O W L E D G E M E N T S I would firstly like to express my extreme gratitude to my supervisor Dr. Reza Vaziri for his guidance, constant encouragement and help throughout this work. Many thanks to Dr. Anoush Poursartip for sharing his expertise and offering his valuable time to discuss various problems on composite manufacturing with me. The advice and encouragement offered by Dr. Goran Fernlund was also invaluable to me during the course of my work. I would also like to thank Dr. Ricardo Foschi and Dr. Hong L i for their guidance and help in my work on reliability analysis. The assistance offered by Mr. Robert Courdji in helping me to understand and run COMPRO is highly appreciated. Many thanks to my colleagues in the Composites Group each of whom contributed in their own way to both my academic and social experiences at UBC. Specially, I would like to thank Mr. Amir Osooly for his invaluable help regarding the new version of COMPRO. -xii i-C H A P T E R 1: INTRODUCTION 1.1 BACKGROUND Fibre-reinforced plastic composite materials are replacing traditional ones in many of today's structural engineering applications. When using composites, unlike their metallic counterparts, the complete (large-scale) structure is manufactured from the raw materials in one step. Our focus is high-performance structural components made of advanced thermoset matrix composites typically employed in the aerospace industry. An autoclave process is commonly used to manufacture such structures. Broadly speaking, this process involves stacking of pre-impregnated sheets of unidirectional fibres (commonly called prepreg) at various orientations over a tool of desired shape and then subjecting the whole assembly to a controlled cycle of temperature and pressure inside an autoclave. The process results in the compaction and curing of the composite part. However, because of the residual stresses that build up during the process, the precise shape and dimensions of the final part after tool removal are often difficult to control. In thin-walled composite structures spring-back and warpage are commonly observed process-induced imperfections. In addition to the detrimental effect on the dimensional accuracy, process-induced residual stresses are also critical to the performance of the manufactured parts. Sometimes, the residual stresses are so severe that transverse matrix cracks and delaminations are introduced within the component after processing, and the parts fail earlier than expected. Therefore, the manufacturing process and its consequences on the material and geometrical properties of the final structure must be taken into account when considering the response to, and reliability under, the applied loads. Otherwise, the designer of composite structures is forced to Chapter 1: Introduction use empirical safety factors, which may result in a conservatism that, in many cases, nullifies the benefits that are offered by composite materials. There has been work done to develop more sophisticated analytical and numerical models to predict the residual stress development in cured composite parts. These models range from one-dimensional analysis (Loos et al., 1983), two-dimensional analysis (Bogetti et al., 1992), to three-dimensional analysis (Chen et al., 1988). In the same front, a comprehensive, multi-physics, 2D finite element code, COMPRO, has been developed to analyze industrial autoclave processing of composite structures of intermediate size and complexity. The model caters for a number of important processing parameters and the development of residual stress and deformation. This model advances previous work by other researchers and, in addition, accounts for the effects of tool/part interactions, which have been neglected by other investigators. The model assumptions, theoretical background, solution strategies and case studies demonstrating its predictive capabilities have been documented in a recent article by Johnston et al (2001). A brief review of the software will be given in Section 1.2. Historically, engineers have employed empirical knockdown factors to account for the large discrepancies between experimentally measured loads and analytical values due to the presence of geometric imperfections. Knockdown factors are often determined experimentally for a range of distinct structures, materials, and manufacturing processes. But these knockdown factors should be constantly updated to include new experimental results. A more sophisticated approach to determine the knockdown factors is the use of finite element analysis with experimentally measured imperfections as input (Starnes, 1999, L i et al, 1995). However, these approaches depend on availability of experimental measurements for a given processed Chapter 1: Introduction composite part and therefore cannot be used in a general manner to study the effect of imperfections on structural behaviour. Up to the present time, little work has been done on modelling of the composite processing and incorporating the process effects into the structural analysis and design. Furthermore, these manufacture-induced phenomena cannot be accounted for in the current so-called "building block approach" to design of composite structures, where the behaviour and properties of the final, full-scale structure are estimated from smaller scale coupon tests and sub-structure analyses. The objective of this work is to present a computational framework for analysing the complete sequence of events from the manufacturing stage to the in-service structural performance of a composite structure. There are many input process parameters that influence, to various degrees, the outcome of the manufacturing, influencing, in turn, the structural performance. Many of these parameters have inherent stochastic variability that result in uncertainty in the final part geometry and material properties. By coupling stochastic and deterministic computational models of the physical events one can study the effect of those variabilities as they propagate through the various stages of the analyses, and eventually influence the response of the composite structure to in-service loading. This would ultimately allow process variables to be controlled or manipulated at the outset (manufacturing level) to result in an optimum or reliable structural performance during service. In this study, we present how process-induced uncertainties can be included in the analysis of composite structures and how the processing parameters can be controlled to achieve a certain reliability level in structural performance. The deterministic process simulation software, COMPRO, the commercial structural analysis software, A B A Q U S , and the inverse reliability software, IRELAN, are used here to determine the process parameters that must be controlled to Chapter 1: Introduction meet structural performance standards with target reliabilities. In the following section, a brief review of the two software packages developed at U B C , COMPRO and IRELAN, and used in this study is given. 1.2 SOFTWARE TOOLS USED IN THIS STUDY 1.2.1 COMPOSITE PROCESS MODELLING SOFTWARE - C O M P R O In this section, a brief review of the software used in this study to simulate the autoclave processing of a composite structure is given. A more detailed description of the software can be found in Johnston (1997). A modular approach, first introduced by Loos and Springer (1983), is adopted for process modelling. As illustrated in1 Figure 1.1, using this approach, the problem is divided into a series of 'sub-models', or 'modules', each responsible for performing a related set of tasks such as prediction of temperature and resin degree of cure etc. Coupling between these modules is maintained by solving each in sequence as the solution marches forward in time. In this way, the predictions of one module can be used by another. A brief description of each module is given below: 1. Autoclave Simulation Module: At the beginning of each time step, this module simulates the automated autoclave control and response and updates all autoclave environmental variables including autoclave air temperature and pressure, and vacuum bag pressure. 2. Thermochemical Module: This module is responsible for calculation of temperature in the structure of interest and the modelled tooling as well as the degree of cure in composite component. Chapter 1: Introduction 3. Flow Module: This module uses a percolation model to predict resin flow and fibre bed compaction in composite components (see Hubert, 1996). 4. Stress Module: This module is responsible for modelling the internal stresses and stress-induced deformations developed within a component during the processing. Through integration of this module with analyses of component temperature and resin degree of cure (thermochemical module) and resin flow and fibre bed compaction (flow module), this model can examine all five major sources of process-induced stress and deformation identified in the literature (Johnston, 1997). These are: i . Thermal expansion i i . Resin cure shrinkage i i i . Gradients in component temperature and resin degree of cure iv. Resin pressure gradients (resulting in resin flow) v. Mechanical constraints due to tooling 5. Input Module: controls acquisition of all input parameters required for a program run. This may include translation of the problem finite element description from an external preprocessor to a form used internally by COMPRO. There are 8 input files that a user has to specify: Batch file: Contains a list of all projects to be run in a batch. Project file: contains a list of the input files containing the project description. Finite element description file: contains the finite element description (node, element, region, boundary etc.) of the problem. Chapter I: Introduction Boundary and initial condition specification file: contains specification of all model boundaries (thermochemical, stress and flow) and initial conditions (initial temperature, initial resin degree of cure). Program control data files: contain all information used by COMPRO to control its internal operation including modules to be executed, time step control data, output time intervals, etc. Process cycle definition files: contain all information related to specification of the autoclave process cycle including controller specifications, virtual instruments, etc. Lay-up definition files: contain ply lay-up specification information for all composite material regions. Material identification files: contain the identification names of the materials in each model region. 6. Output Module: controls output of all parameters to external devices, for example, to the screen and hard drive. A l l three major solution modules (thermochemical, flow and stress/deformation) employ the finite element method to solve for the component state variables during processing. The component finite element description can be created using virtually any finite element pre-processor which can generate an output text file. Element Type The previous version of COMPRO employed four-noded isoparametric quadrilateral elements as shown in Figure 1.2a. The new version of COMPRO can employ both four-noded and eight-noded isoparametric quadrilateral elements as shown in Figure 1.2b. For the eight-noded Chapter 1: Introduction analysis, the user only needs to define the four-noded finite element mesh. COMPRO internally converts the four-noded elements to eight-noded elements (Osooly, 2001). Region Definition In COMPRO, a region is defined as a collection of elements containing a single type of material. In addition to defining element material type, regions are used in combination with reference boundaries (next section) to define ply lay-up and element local axes orientation. Regions are also used to define which elements represent process tooling, usually removed at the end of the processing cycle, and which are part of the composite structure (whether specifically defined as composite materials or not). To minimize the number of elements required to discretize a given problem, multiple plies of any orientation may be contained within each composite material element as shown in Figure 1.3. It is not required that a whole number of plies be contained within an element (e.g., a ply can contain 2.6 plies), but the ply plane must be parallel to the element local coordinate x'-axis. Ply orientation is defined according to the right-hand rule with the angle measured from the positive x' axis. Ply lay-up in a region is defined with respect to the reference boundary for that region (i.e. ply 1 starts at the reference boundary). Boundary Definition In COMPRO, a boundary is defined by node groups (never by element 'faces' as in some finite element codes). Boundaries are used for three main purposes: to specify 'boundary conditions' (e.g. heat flux, vacuum pressure), for definition of ply lay-up, and for definition of the element local axes orientation as shown in Figure 1.4. Chapter 1: Introduction Material Property Definition Each time a module is run, the relevant material properties for each element are recalculated from state variable values (e.g., temperature, degree of cure and volume fraction) at the element centroid. Orthotropic material behaviour is assumed in all cases. Orthotropic composite properties may be either specified directly or calculated from micromechanics models which use as input fibre volume fraction, isotropic resin properties and transversely isotropic fibre properties. A l l materials are assumed to be 'instantaneously linear elastic'. This means that, although the modulus may riot be constant during processing, at any given instant (i.e., during any time step), they exhibit linear elastic behaviour. Virtual Autoclave Simulation Virtual autoclave simulation is one of the main special features of COMPRO. The autoclave simulation includes the simulation of the behaviour of the autoclave during processing including all environmental control hardware (i.e. heating, cooling, pressurizing etc.). Tool Part Interaction The amount of stress transferred between the tool and the part is controlled by including a thin elastic shear-layer into the finite element mesh. Depending on the value of the elastic modulus assigned to the shear layer, a range of tool part interface conditions can be simulated. At the one extreme, i f the properties of the shear layer are the same as the tooling, high shear stresses arise between the tool and the part (perfect bonding). At the other extreme, by giving the shear layer low values for En and G13 (i.e. longitudinal Young's modulus and transverse shear modulus), Chapter I: Introduction relatively little stress is transferred between the tool and the part. The shear layer is considered to be part of the tooling and it is removed during the final tool removal step. Tool Removal Simulation At the end of the process cycle, a simulation of the tool removal from the component is performed to determine post-processing component shape. The tool removal simulation is performed in a single elastic step and the calculated changes in displacements are added to the pre-tool-displacements to give the final component shape. 1.2.2 INVERSE RELIABILITY ANALYSIS SOFTWARE - I R E L A N In this section, a brief review of the theory behind the inverse reliability analysis software I R E L A N is given. The software was developed in the Department of Civil Engineering, at UBC. A more detailed description of the software can be found in Foschi and L i (1999b). In combination with an analytical model for the system, safety assessment requires forward reliability methods for estimating the reliability or the probability that it will perform as intended. On the other hand, engineering design requires an inverse reliability approach of determining the design parameters to achieve pre-specified target reliabilities. In general, an inverse problem involves finding either a single design parameter to achieve a given single reliability constraint or multiple design parameters to meet specified multiple reliability constraints. Forward Reliability Method The implementation of the reliability analysis is based on a description of the limit state of interest by a performance function G(x) as follows: Chapter 1: Introduction G(x) = G (x/, X2, »,.xn (1.1) where the x = (xj, X2, >, x„)T is an n-dimensional vector of intervening random variables. Some of these may affect the demand on the system, denoted by D, while the others may influence the system capacity C to withstand the demand. By convention, the performance function (Equation 1.1) may be written as Thus, the system will not perform as intended (failure) i f the combination of the intervening random variables results in G < 0. The corresponding probability of such event (Prob(G < 0)) is then called the probability of failure. Conversely, the combination of the intervening random variables resulting in G > 0 will make the system perform as required (survival) and the corresponding probability (Prob(G > 0)) is termed the reliability. The situation G = 0 is a limit-state between failure and survival and the corresponding surface, G(x) = 0, is called the limit-state failure surface. Since, in general, there could be a number of random variables involved in G, the exact calculation requires the joint probability density function of all random variables and integration over the failure region G < 0, However, this exact approach can hardly be applied since the joint probability function / is unknown and very difficult to find. An alternative method is the straightforward, standard computer simulation (Monte Carlo method) which is simple to implement and can converge to the exact solution. However, it could be very computationally demanding, especially when dealing with a low probability of failure. To overcome this, different methods have been G = C-D (1.2) (1.3) G<0 -10-Chapter 1: Introduction proposed to increase the simulation efficiency such as importance sampling, adaptive important sampling, etc. Both the importance sampling technique and the adaptive importance sampling technique concentrate on sampling points in the region which mainly contributes to the failure probability, instead of spreading them out evenly across the whole range of possible values of each random variable. A second alternative is the use of approximate methods that have been developed during the last three decades, such as the FORM/SORM procedures (First Order or Second Order Reliability Methods), and which are based on the calculation of the reliability index j5 (ft = MG/crG where jUq and <rG are the mean and standard deviation of the G function). From this index, by use of the Standard Normal probability distribution function 0 ( ), the probability of failure P/ and the reliability Pr can be estimated approximately A forward reliability software, called R E L A N was developed to carry out a variety of reliability calculations (Foschi, et al., 1999a). FORM/SORM Procedure These methods are the most complex, both mathematically and conceptually, among all probabilistic analysis methods and will therefore not be described in detail. More details of these methods can be found in L i (1999). pf = o(-d) (1.4) and (1.5) Chapter 1: Introduction The approach of these methods is to transform the integral given in Equation 1.3 to an approximately equal integral that can be efficiently evaluated. This is done in the following step-by-step procedure: 1. Transform the design variable into standard normal distributions. That is, transform G(x) = 0 into G(u) = 0 where u is a vector of independent standard normal variables. 2. Identify the most probable point or design point. For a given limit state function, the main contribution to failure probability comes from the region where G is closest to the origin in the transformed design variable space (u-space). The design point is the closest point to the origin in the transformed space (Figure 1.5). 3. Develop a polynomial approximation to the performance function (G) around the design point (linear function in F O R M and quadratic function in SORM). Compute the probability of failure using the newly defined G - function. Response Surface Method For complex structures, the performance function is not available as an explicit function of the random design variables. The performance of the structure can only be evaluated numerically through a structural analysis procedure such as the finite element method. The goal of the response surface methodology (RSM) is to find a predictive equation (such as Equation 1.6) relating a response to a number of input variables. (1.6) -12-Chapter I: Introduction where x„ (i = 1, ,n) are the basic variables and the parameters a, bt, ch (i = 1, ,n) are to be found. Then this equation can be used to determine the response to given values of input variables, instead of having to repeatedly run the time-consuming deterministic structural analysis. Inverse Reliability Method The objective of the inverse procedure is to find, directly, the design vector, d, so that the target reliability level corresponding to each specified limit state would be satisfied. In general, these are also regarded as random variables, with corresponding mean values and standard deviations. The problem may include the determination of the mean alone when the coefficients of variation are specified, the standard deviations when the means are given, or both means and standard deviations. The inverse reliability problem can be solved by trial and error by repeating forward reliability analyses such as FORM/SORJVI, interpolating the design parameters being sought. This is inefficient and time consuming, especially when multiple design parameters and multiple reliability constraints are involved. To approach the inverse reliability problem directly, an efficient method based on FORM/SORM methods is proposed by Foschi and L i (1999b). For cases requiring an intensive computational effort or having a non-smooth limit state function, a response surface method for inverse analysis is also proposed and implemented to speed up the calculations. Performance-Based Probabilistic Design This is a newly added feature to the I R E L A N software based on the neural network method. A full detail of this feature is given in Foschi and L i (2000). A brief review is given in this section. -13-Chapter 1: Introduction Performance-based design of a structural system implies the definition of a set of limit states or performance requirements, with corresponding target reliabilities for each state. With these requirements, the designer must consider the uncertainty in the intervening variables, x, and calculate values for specific design parameters, d, to achieve the desired reliabilities, B. For multiple performance criteria, the problem may not have an exact solution, since the target reliabilities for each limit state may not be all achieved as prescribed. In this case, an optimization procedure is used to minimize the difference between the achieved reliabilities and the corresponding targets reliabilities (Equation 1.8). where *F is the objective function to be minimized, /3j is the target reliability corresponding to the j t h limit state, and Bj(d) is the reliability index obtained with the design parameter d. This approach requires, for the given performance criteria and associated target reliability levels, the determination of a set of design parameters. This can be achieved using any of the inverse methods discussed in Foschi and L i (1999b). When the response is highly non-linear, the convergence and the robustness of procedure like F O R M will be adversely affected. So, this software uses simulation technique (importance sampling technique or adaptive sampling technique) following the calculation of a design point using FORM on an approximate quadratic response surface. The required evaluation of the G„(xhd) = 0 with target reliability B, (1.7) (1.8) -14-Chapter 1: Introduction performance function is carried out with a local interpolation using a response database obtained independently through deterministic analysis. 1.3 RESEARCH OBJECTIVE AND SCOPE OF THESIS The objective of this research is twofold. 1) To present the effect of process-induced stresses and deformations on structural performance. 2) To present a computational framework for analysing the complete sequence of events from the manufacturing stage to the in-service structural performance of a composite structure, while taking into account the stochastic nature of the parameters involved. In this case, we combine a deterministic model for the manufacturing process and structural analysis with a probability-based reliability model to determine the process parameters that must be controlled to meet structural performance standards with target reliabilities. The deterministic process simulation software, COMPRO, the commercial structural analysis software, A B A Q U S , and the inverse reliability software, IRELAN, are used for the analysis. The research presented in this thesis is organized as follows: 1) Chapter 2 - A case study is presented to illustrate the processing effects on structural performance and the probability approach to determine the manner in which processing variables can be controlled to meet structural performance standards with target reliabilities. 2) Chapter 3 - The effect of various parameters on spring-in of curved composite parts are predicted numerically using COMPRO and, where possible, compared with experimental -15-Chapter 1: Introduction data. The case study considered, that of L - and C-shaped composite parts, serve to demonstrate the range of process parameters that have to be taken into account for more complex shaped composite structures. It also highlights the advantages of having a numerical model that can capture the effects of many of these parameters in a virtual manufacturing environment. Such detailed quantitative information on a processed part can then easily be transported to a structural analysis program as input, thereby accounting for processing effects in assessment of structural performance. 3) Chapter 4 - The residual stress development during the curing process both at the ply level and the micro-level due to CTE mismatch and cure shrinkage are quantified. The ability to compute these residual stresses at all level is vital for subsequent failure analysis of composite structures. 4) Chapter 5 - Conclusion and future work. -16-Chapter 1 Introduction 1.4 FIGURES Figure 1.1: Schematic of COMPRO structure and program flow (Johnston, 1997). (a) (b) Figure 1.2: Four-noded and eight-noded isoparametric quadrilateral elements employed in COMPRO. -17-Chapter I Introduction Figure 1.4: Curved region, showing element local coordinates, local axes orientation and the region reference boundary (Johnston, 1997). -18-Chapter 1 Introduction Failure Region G<0 Figure 1.5: Geometric representation of reliability calculation for FORM (Li, 1999). -19-C H A P T E R 2: PROBABILISTIC DESIGN OF A COMPOSITE P L A T E AGAINST B U C K L I N G 2.1 INTRODUCTION The basic structure of an aircraft is made of three types of elements: shells, plates, and beams (stiffeners). The fuselage of the aircraft is made of stiffened shells and the top and bottom surface of the wings are made of stiffened plates as shown in Figure 2.1. The major design goal is to produce a light structure which is strong enough to withstand the aerodynamic loads. Since the skins on the wings and the fuselage of an aircraft are fairly thin due to the weight requirement, their likely mode of failure under the aerodynamic loading is buckling. The traditional way of designing a thin-walled metallic structure against buckling is to determine the buckling load by a simple deterministic analysis and to reduce this buckling load by a knockdown factor. This knockdown factor is to account for the difference between the predicted buckling load and the actual buckling load measured from a test (Starnes et al, 1999). A linear bifurcation analysis is often used for this deterministic analysis on a geometrically perfect shell structure with the nominal values of material properties and idealized boundary conditions. A simply supported boundary condition is often used to simplify the analysis and because of its relatively low structural stiffness it is also the worst case scenario. The fabrication process generally introduces initial imperfections in metallic shell structures and accounts for the difference between the predicted and actual buckling load. The effect is -20-Chapter 2 Probabilistic Design of Composite Plate Against Buckling particularly severe when the imperfection shape coincides with the buckling mode of the structure (Chryssanthopoulos et al, 1995). The current practice of incorporating the initial imperfection into the design is to include the measured imperfection data in the analysis if these data are available or to include some approximate imperfection shape function if these data are not available (Starnes et al, 1999). For a prototype, the imperfection can be measured experimentally and then can be incorporated into the theoretical analysis to predict the buckling load. But in normal production runs, this approach is impractical due to variability in manufacturing methods, material properties etc. The realization that any further advances towards more accurate buckling load predictions of thin shells depend on the availability of extensive information about realistic imperfections data and their correlation with manufacturing processes has led to the establishment of an International Imperfection Data Bank for metallic structures (Arbocz, 1991). In this data bank, the characteristic initial imperfection distribution that a given fabrication process is likely to produce is established and then this information is combined with some kind of statistical analysis of both the initial imperfection and the corresponding critical buckling load. In recent past, metallic thin-walled structures have been replaced by composite structures, especially in high-performance aerospace structural applications, due to their superior specific strength and stiffness. Compared to its metallic counterparts, the analysis of a composite shell structure becomes more complex due to the complex nature of its manufacturing process and the high variability of their material properties. Due to a variety of uncertainties that are involved in fabrication of composite panels, it is impossible to determine the amplitude and the shape of imperfection. The traditional knockdown factor design approach will not be conservative if the composite shell structural behaviour is not fully understood and the traditional sources of design -21-Chapter 2 Probabilistic Design of Composite Plate Against Buckling knockdown factors for predicting shell buckling loads do not include data for composite shell structures (Starnes et al., 1999). The effects of traditional initial geometric shell-wall imperfections on the buckling response of unstiffened thin-walled graphite-epoxy cylindrical shells is analysed by Starnes et al. (1999) in which the initial geometrical imperfections and the thickness variations are experimentally measured and incorporated in non-linear structural analysis. As mentioned earlier, the above approach is impractical in normal production runs. The best way to handle the above problem is to model the manufacturing process of the composite structures in order to predict the processing imperfection accurately and to combine this process model with structural analysis codes. To this end, a two-dimensional finite element code, COMPRO, has been developed by the Composites Group at the University of British Columbia for modelling the autoclave processing of complex composite structures (Johnston et al., 2001). In addition to the complexity of this manufacturing process, variability of composite structures at both the material and the processing level are important and need to be considered. Handling these uncertainties in the traditional knockdown factor method normally leads to an overly conservative or unconservative design. The probabilistic design approach is an alternative to the traditional knockdown factor approach to account for the uncertainties. In the literature, many probabilistic based design approaches are proposed to account for the uncertainties involved in material and manufacturing levels (Shiao et al., 1999; Elseifi et al, 1999; Cederbaum et al, 1996). A l l these probabilistic methods are used to analyse the reliability of the composite structures for the given uncertainties at the material and manufacturing levels. But in actual design practice, we need to solve the inverse problem i.e. for the given reliability targets, the material properties and the processing condition should be found -22-Chapter 2 Probabilistic Design of Composite Plate Against Buckling with certain statistical distribution. For this purpose, Foschi and L i (1999b) have developed a probabilistic design software, IRELAN, which can handle the inverse problem. The acceptable probability of failure is the criterion that is used to determine whether the design is acceptable. Legal, technical and socio-economic considerations are involved in deciding this target reliability level. It should not be decided by the engineer performing the probabilistic analysis rather it should be specified by the agency certifying the structure. The target reliability level for aircraft structure is generally very high and this issue remains unresolved at the present time (Long and Narciso, 1999). In performance-based design of a structural system, the design parameters should be calculated to meet one or more performance requirements with specified target reliabilities. For this multiple performance criteria, the problem may not have an exact solution, since the target reliabilities for each limit state may not be achieved as prescribed. An optimization procedure is proposed by Foschi and L i (2000) in which the difference between the achieved reliabilities and the corresponding targets are minimized. In this case study, a combined deterministic-stochastic simulation approach is presented to demonstrate the manner in which manufacturing-induced variability in composite structures can be controlled to achieve targeted reliability levels in structural performance. The processing finite element code, COMPRO, is integrated with the non-linear structural analysis and reliability analysis software to compute the statistics of parameters that must be controlled at the manufacturing level in order to achieve an optimum structural performance during service. In the literature, thin-walled composite panels are identified as highly sensitive to geometric imperfections. For illustrative purposes, the buckling behaviour of a composite plate in the -23-Chapter 2 Probabilistic Design of Composite Plate Against Buckling presence of process-induced residual stresses and deformations is examined. Since process modelling is not an established source and research is still being carried out in this area, we have chosen this relatively simple problem for which a large body of experimental results have already been generated at UBC. These results serve as a validation of the process modelling undertaken here and therefore ensure confidence in using the output of the process model as input for structural modelling and reliability analysis. Buckling is a potential failure mechanism for thin-walled structural components commonly employed in aircraft structures. 2.2 PROBLEM DEFINITION A unidirectional, carbon fibre-reinforced plastic (CFRP) composite laminate with a span of 1200 mm, a width of 100 mm and a thickness of 1.6 mm (consisting of 8 plies each of thickness 0.2 mm) is considered for this study. This is the most slender of many geometries of flat plates studied experimentally by Twigg (2001). The plate is manufactured in an autoclave under the temperature and pressure cycle shown in Figure 2.2 on an aluminum tool. In the ideal case, the plate is assumed to be straight and stress free and therefore there is no deformation prior to buckling. But in reality, due to processing, various stresses and deformations are induced in the structure. The resulting initial imperfection affects the lateral response of the plate when subjected to a certain axial compressive load (Figure 2.3). The amount of the initial imperfection may depend on many processing parameters such as curing temperature, pressure, heating rate, cooling rate, holding period, fibre volume fraction etc. Some of these parameters are more influential than others. Practically, it is impossible and costly to keep all these processing parameters constant. In this design problem, first, the processing parameters that the imperfections are most sensitive to are identified and then the variability of -24-Chapter 2 Probabilistic Design of Composite Plate Against Buckling these processing parameters are found so that the plate withstands the applied axial load without exceeding the given deformation limit. 2.3 PROCESS MODELLING The part and the tool in COMPRO are modelled using 4200 four-noded isoparametric quadrilateral elements under plane strain condition as shown in Figure 2.4. Because of the symmetry, only half of the geometry is modelled. As discussed in Chapter 1, an elastic shear layer is incorporated into the finite element mesh in order to simulate the tool-part interaction. The shear layer properties are adjusted to match the experimental warpage results of Twigg, (2001). In the experiment, two layers of FEP sheets were placed between the tool and the part to reduce the tool-part interaction. In this case, using a shear layer modulus of 5.52 kPa resulted in an accurate prediction of the deflection profile as shown in Figure 2.5. 2.4 STRUCTURAL MODELLING The predicted deformation, stresses and material properties of the composite at the end of the process modelling are transferred as input to a 3D non-linear structural analysis software, A B A Q U S (ABAQUS 5.8). The material properties transferred over to A B A Q U S were fully cured properties of the laminate computed using the micromechanical models in COMPRO (Johnston et al., 2001). These properties for the nominal case are listed in Table 2.1. Being a 2D plane strain model, COMPRO can only provide information on a 2D plane strain section of the composite part. Therefore, to obtain the 3D spatial distribution of residual stresses and deformations necessary for the structural analyses, the output from COMPRO have to be mapped onto the 3D elements of A B A Q U S . This is accomplished through a superposition method in which the deformations and stresses from COMPRO are assigned as initial conditions to -25-Chapter 2 Probabilistic Design of Composite Plate Against Buckling structural elements in A B A Q U S . This method is explained in detail in the next section (Section 2.5). The virtually manufactured structure is considered to be simply supported and subjected to progressively increasing compressive load P as shown in Figure 2.6. First the critical Euler buckling load, Pcr, of a perfect plate, free of imperfection and with nominal material properties, is calculated using the eigenvalue analysis option in A B A Q U S and this value is found to be 29.5 N . Then this critical load is applied to the plate with initial imperfections. The non-linear relation between the axial load P and the out-of-plane displacement w is obtained by increasing the load P incrementally up to the limiting buckling load Pcr and computing the corresponding displacement. In general design practice, a structure is designed to its maximum capacity. A load larger than the capacity of the structure results in failure and a load lower than the capacity of the structure results in over-conservative design. In this case, the application of the critical buckling load leads to the buckling failure of the structure and a lower load leads to an over-conservative design. Hence, 75 percent of the critical buckling load, which gives a moderate deformation, is considered as the applied load and the deflection at this load level is considered for the subsequent reliability study. 2.5 PREDICTION OF THE 3D STATE OF STRESS This section will discuss the procedure for predicting the three-dimensional state of stress in an element using the knowledge of stresses obtained from a plane strain type analysis such as COMPRO. -26-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 2.5.1 PLANE STRAIN In plane strain problems, it is assumed that the out-of-plane strain is zero. The assumption is valid under the following conditions: 1. There is no loading in the out-of-plane direction and the in-plane loads are independent of the out-of-plane coordinate. 2. The out-of-plane dimension is very large compared to the in-plane dimensions. Under thermal loading conditions, a plane strain element, by virtue of its constraints, will not allow any thermal expansion in the out-of-plane direction thus resulting in higher out-of-plane stresses and more deformation in the plane of the element. 2.5.2 THREE DIMENSIONAL PREDICTION Since the complete state of stress (i.e. all six components of stress) is required in our case studies, an approach based on superposition principle is proposed here to predict the 3D states of deformation and stress using results obtained from 2D plane strain computations. It is well known that the principle of superposition is only valid for linear elastic materials. Since the deformation of the part is constrained by the tool and most of the residual stresses develop during the cooling down part of the process cycle when the evolution of the elastic modulus of the resin has terminated, the application of the superposition method is justifiable. A l l the following calculations are carried out using A B A Q U S . -27-Chapter 2 Probabilistic Design of Composite Plate Against Buckling (a) For applied mechanical load First consider a cube of an isotropic material (the material properties used in A B A Q U S are E = 2 GPa, v = 0.25, a = 25 um/m/°C, however the results are expressed here symbolically) under an applied mechanical load, P, in all three direction as shown in Figure 2.7(a). In a three dimensional problem, the cube is allowed to deform freely in all three directions. From Hooke's law and the equilibrium equation, the relation between the loads and the displacements can be written as P_ I2 P_ J P_ 7 E(l-v) (l + v) ( l -2v) 1 v (1-v) (1-v) V . V 1 (1-v) (1-v) (1-v) (1-v) 1 (2.1a) which results in u = v = w = 8 = P(l - 2v) IE (2.1b) It will be shown that the same results can be obtained by applying the loads separately in two steps and superposing the results. Stepl On the cube, only the loads in 1 and 2 directions are applied and the cube is allowed to deform freely in 1 and 2 directions but it is constrained in the 3 r d direction (similar to plane strain condition) as shown in Figure 2.7(c). The displacements can then be obtained from the following equation (Equation 2.2a). -28-Chapter 2 Probabilistic Design of Composite Plate Against Buckling P_ I2 P_ I2 E(l-v) (l + v)( l -2v) 1 v (1-v) (1-v) 1 (1-v) (1-v) 1 (1-v) (1-v) w = 0 (2.2a) solving for S] results in J , = ( l + v / ( 1 2v)=(\ + v)S IE (2.2b) It can be noted that the displacements are (1+v) times higher than the actual three-dimensional displacements. The constraint in the 3 r d direction results in an internal stress, cx3 = P3II1, in that direction, such that P3 = 2vP (2.2c) Step 2 In the 2 n d step, the load in the 3 r d direction, P, with the stress calculated in the previous step, considered as an initial stress, is applied to the cube as shown in Figure 2.7(d) and the cube is allowed to deform freely in all three directions. The load-displacements relation can then be written as I2 = 0 I2 E(l-v) (l + v) ( l -2v) 1 (1-v) (1-v) 1 (1-v) (1-v) 1 (1-v) (1-v) / I w 7 (2.3a) resulting in the following displacements -29-Chapter 2 Probabilistic Design of Composite Plate Against Buckling W = IE (2.3b) pa - 2v) -v— = -vS IE Now, the actual three-dimensional results (Equation 2.1b) can be obtained by combining (superposing) the results of both step 1 and step 2 (Figure 2.7(b)). (b) For thermal load The superposition principle for thermal loads is similar to mechanical loads, but instead of the applied load, the thermal expansions are superposed. First consider the same cube as in the previous case (Figure 2.8(a)) under thermal load (temperature increment AT = 100 °C). If the cube is allowed to expand freely (Figure 2.8(b)), there will be no stresses in the cube and the deformation will be given by But in the plane strain case (Figure 2.8(c)), since the expansion in the out-of-plane direction (axis 3) is prevented, the out-of-plane stress will be nonzero and the in-plane deformations are (1+v) times greater than the actual deformations, i.e. 8 = laAT (2.4) <j\ =cr2 =0 (2.5a) a\ = -EaAT Therefore, the deformations are given by 8x =S2 =Q + v)laAT S,=0 (2.5b) -30-Chapter 2 Probabilistic Design of Composite Plate Against Buckling Now, i f these stresses are introduced as initial stresses in a cube (Figure 2.9(a)) and the cube is allowed to expand freely, it will expand in the out-of-plane direction and contract in the in-plane direction such that (Figure 2.9(b)): S, = = -via AT (2.6) 83 = la AT The combination of Equations 2.5b and 2.6 (Figure 2.8(c) and Figure 2.9(b)) will give the actual deformation state. The above method can be extended to different numbers of elements by simply dividing the cube into several elements (Figure 2.9(c) and (d)). The same is true for shell elements although there is no plane strain shell element available. This statement can be explained through an example. Consider a plate made of the same material as before, simply supported at all four corners and subjected to thermal loading. The dimensions of the plate and the number of shell elements used to model this plate are given in Figure 2.10. Due to the symmetry, only a quarter of the problem is considered for analysis. First only a row of elements along the x direction is considered for the analysis as shown in Figure 2.11(a). To make it similar to a plane strain solid element, the displacement in the y direction and the rotations about the x and z axes are constrained to be zero. The stresses and deformations of these elements are found for a temperature increment of 10 °C. Then similar columns of elements with initial stresses are combined to form the full structure and allowed to equilibrate under the specified boundary conditions (Figure 2.11 (b)). -31-Chapter 2 Probabilistic Design of Composite Plate Against Buckling The results of the above method are compared with the full 3D analysis of the plate in Figures 2.12 - 2.13. An excellent agreement between the superposition method and the full 3D analysis can be observed. The above method is a 2D shell to 3D shell mapping. But, i f this problem is to be modeled in COMPRO, it should be modeled using plane strain elements as shown in Figure 2.14 (a). From this 2D plane strain analysis, the three dimensional results can be obtained in one of two ways: 1. Combine similar rows of 3D solid elements with initial stresses (obtained from plane strain analysis) in the y direction to form the full plate and allow it to equilibrate under these stresses (all 6 components) to obtain the 3D states of deformation and stress (i.e. 2D solid to 3D solid mapping). 2. Map the plane strain element results to a row of 3D composite shell elements. This mapping can be done by considering each row of the solid elements as a layer of the composite shell elements and averaging the stresses to the mid-point of the shell elements as shown in Figure 2.14 (b). (Each shell element covers 20 solid elements. The stresses in these 20 solid elements are averaged and transferred to the mid-point). Then similar rows of these shell elements with initial stresses are combined to form the full plate and allowed to equilibrate under these stresses. Since the first method requires a huge amount of 3D brick elements, the second method is preferred here. The results from this modified method and the actual 3D results are compared in Figure 2.15-2.16. The differences in the results are due to the averaging of the stresses. Also when the stresses of the solid elements are mapped on to the shell elements, it is assumed that the normal stress in the z direction (out-of-plane stress) and the shear stresses in the x-z and y-z -32-Chapter 2 Probabilistic Design of Composite Plate Against Buckling planes are zero. This assumption is valid for thin parts and it is expected that the results would be less accurate as the part thickness increases. 2.6 PARAMETRIC STUDY The parameters that influence the buckling characteristic of a composite structure can be grouped into three categories: parameters related to material properties and geometry of the structure, parameters related to processing of the structure (these parameters also influence the material properties and geometry of the composite structure), and the parameters related to loading and boundary condition. In the following reliability study, the parameters related to processing of the composite structure are considered as random variables and the other two categories of parameters are considered as deterministic. The following parametric study was performed to identify the most important process variables that have an effect on the warpage of the plate. The process parameters considered in this study are maximum autoclave temperature (7), fibre volume fraction (Vj), ply misalignment (0), autoclave pressure ip), heating rate (HR), hold time (HT) and cooling rate (CR) (here only a limited number of processing parameters are considered for demonstration purposes of the method even though there may be some other parameters which affect the warpage and load capacity of the part). The sensitivity factor, which is defined as the ratio of the percentage increase in warpage to percentage increase in the parameter, is calculated and plotted for each parameter as shown in Figure 2.17. From this figure, it is identified that the part warpage is most sensitive to temperature, fibre volume fraction and fibre misalignment. From now on it is assumed that the part warpage depends on only these three processing parameters. -33-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 2.7 RELIABILITY-BASED DETERMINATION OF PROCESS CONTROL PARAMETERS In the following reliability study, the three processing parameters which are identified as the most influential on the warpage of the composite structure in the parametric study are considered as random variables and the rest of the processing parameters are considered as deterministic. Among the chosen three random variables, the distorted shape of the plate after processing, w, is a function of T, Vf, and 9, while the critical load, PCI is only a function of 0 and Vf, which influence the material properties. Since the variables T, Vf, and 0 are subject to uncertainty due to possible lack of tight process control, the outcomes w and Pct are also random. For the purposes of this study the performance standards for the plate were defined as follows: 1. The deflection w should be within a range, from 20mm to 40mm, with specified reliability. These two limits are the design requirements for the serviceability and ultimate conditions. The lower limit (span/600) will ensure that buckling will not result in a sudden departure from the flat configuration (i.e. there would be some warning before collapse). The upper limit (span/300) will ensure that the plate will not be excessively deformed. These two conditions can be described by a set of two "performance functions", G\ and G2 as follows G] = w(0,Vf,T)-2O (2.5) G2=4O-w(0,Vf,T) 2. The applied load P should result in a sufficient reliability against buckling. This condition is expressed as a performance function G3, -34-Chapter 2 Probabilistic Design of Composite Plate Against Buckling G3=Pcr(0,Vf)-P (2.6) Given that the variables are random, the failure or non-performance probabilities correspond to the events G\ < 0, G2 < 0 and G3 < 0. Provided that the statistics of the variables are defined, these probabilities can be readily evaluated by forward reliability procedures. Although the formulation of the problem is simple, the evaluation of the G function for each combination of random variables is time consuming and generally makes standard simulation procedure quite inefficient. To make the simulation procedure efficient, the actual performance function, G, can be replaced by a response surface and the failure probabilities can be obtained by any forward reliability procedures using the software R E L A N (Foschi et al., 1999a). But it is difficult to find a simple response surface to represent G over the entire domain of interest. In this case, calculation is carried out by local interpolation of the deterministic database (Foschi and L i , 2000). The deterministic database is constructed for all possible ranges of the random variables as shown in Table 2.2. For each combination of the random variables, the stresses and deformations obtained from COMPRO are fed to A B A Q U S to obtain the required results (zb.5%, Pcr)- The main advantage of this method is that the database can be expanded in future by either increasing the number of random variables or increasing the range of the random variables. The probability of failure in each case is obtained by importance sampling simulation around an anchor point, which is determined by using an approximate response surface and the First Order Reliability Method (FORM). The method is explained in detail in Chapter 1. In each case, the probability of failure can be expressed by the associated reliability index fi. -35-Chapter 2 Probabilistic Design of Composite Plate Against Buckling Example 1 For the inverse problem, let us assume that the task is to find the design parameters Vf, T and P , mean values of the corresponding variables, assuming that the remaining relevant statistics are as shown in Table 2.3. That is, knowing the accuracy with which the angle 6 and the variables Vf and T can be controlled, their mean values are to be found (except for 6 = 0.0°). Along with these, the mean of the load P that can be applied is to be determined in such a way that the plate response has the following target reliabilities for each of the three performance criteria: 1. For G\, a failure probability of 0.05, or a target reliability index /3\ = \.645 2. For G2, a failure probability of 0.05, or a target reliability index fii = 1.645 3. For G3, a failure probability of 6.2x10" , or a target reliability index /% = 2.500 The first two conditions imply a 90% confidence in the interval 20mm < w <40mm. The inverse reliability computer program, IRELAN, searches for the optimum solution, calculating the reliability indices for different combinations of the design parameters and minimizing the differences between the calculated and target reliabilities, /? and J3j (/'= 1,2,3), respectively: In each case, the reliability indices were calculated by the forward procedure previously described, using the deterministic database for w and Pcr. Of course, the results are not an exact (2.7) -36-Chapter 2 Probabilistic Design of Composite Plate Against Buckling achievement of the target reliabilities, but rather an optimized solution, and are shown in Table 2.4. Example 2 In a second design example, we assume that the load P is given (with corresponding statistics), and that the mean of the variables Vf and T are the design parameters that have to be calculated to meet the same performance criteria. That is, process parameters have to be found to satisfy performance under given load demands in the structural application. Let us assume that P has an Extreme Type I distribution, with a mean of 15 N and a coefficient of variation 0.25. Statistics for #are as shown in Table 2.3. The results of this second optimization are shown in Table 2.5. As shown, the lower load P, with the same coefficient of variation and target reliabilities, allows a slightly lower Vf and mean temperature T . 2.8 OTHER SOURCES OF VARIABILITY In the above analysis, the influence of variations in part thickness on the buckling is not considered. During processing of a composite part, the thickness of the part changes due to the following sources: inaccuracy in lay-up, drawing vacuum, cure shrinkage, resin flow and the change in resin modulus. From the above COMPRO can predict the thickness changes due to only the last three sources. The effect of the thickness change due to these three sources is explained below. Cure shrinkage of the resin causes permanent thickness reduction of the lamina. Since the weight ratio of the fibre and the resin remains the same, cure shrinkage will not alter the material properties of the lamina. According to the expression for the Euler critical buckling load for a -37-Chapter 2 Probabilistic Design of Composite Plate Against Buckling simply supported column, x percentage of reduction in the thickness causes 1 0 0 ( l - ( l - j ^ ) 3 ) percentage reduction in the buckling load. For the composite material used in the above analysis, cure shrinkage causes 3% reduction in thickness that in turn causes around 9% reduction in the buckling load. During processing the resin modulus changes from a very low to relatively very high value (roughly 3 orders of magnitude increase in modulus). Due to the application of the autoclave pressure, this change in the resin modulus induces permanent deformation in the lamina. If the initial and final transverse moduli of the lamina are Eri and Erf respectively and the autoclave pressure is p, by a simple strength of materials calculation, the permanent deformation of the lamina per unit thickness is given by P(— —) . In this case, p is 586 kPa Eri and Erf values Eri Erf are 183.8 MPa and 10.16 GPa respectively. This results in only 0.3% reduction in the thickness. The third source of thickness reduction is the flow of resin. The excess resin in the lamina is removed by using a bleeder on top/side of the part (Hubert, 1996). The resin flow reduces the thickness of the part and in the meantime the fibre volume fraction of the lamina also increases. By assuming that the longitudinal modulus of the part is V/Eif, an x percentage reduction in the thickness due to resin flow causes 100(1 - (1 - j ^ ) 2 ) percentage reduction in the buckling load. For the flat part considered in the above analysis, top bleed of the resin causes around 7% reduction in thickness that in turn causes 14% reduction in the buckling load capacity. The reduction in thickness due to the three sources mentioned above is the average values. These values can be easily calculated and accounted for in the design. But, the variability of the change -38-Chapter 2 Probabilistic Design of Composite Plate Against Buckling in thickness due to the above three sources with other processing parameters such as temperature, pressure, V/, etc., should be investigated. 2.9 SUMMARY AND CONCLUSIONS The methodology shown here is an effective means of linking manufacturing process control, structural analysis and evaluation of performance reliability for composite structures. Although a flat plate was used here as an application, the method is completely general and could be used with any structural form. The approach adopted here can potentially lead to a new design methodology that uses the as-manufactured rather than as-designed conditions as the basis for structural analysis leading to the reduction or elimination of uncertainties and costs associated with the use of empirical safety factors currently employed in design of composite structures. -39-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 2.10 TABLES Table 2.1: Laminate material properties Longitudinal elastic modulus, E\ \ 126 GPa Transverse elastic modulus, E22 - £33 10.9 GPa Major Poisson's ratio in-plane, V12 0.264 Longitudinal shear modulus, Gn = G13 5.04 GPa Transverse shear modulus, G23 3.28 GPa Table 2.2: Database e Ply misalignment (degree) v f Fiber volume fraction T Temperature (F) A75% Out of plane displacement at P/P c r = 75% (mm) Per Critical buckling load (N) e v f T A75% Peri 0 0.589 350 31.89 29.4772 0 0.589 385 36.33 29.4772 0 0.589 367.5 34.11 29.4772 0 0.589 332.5 29.76 29.4772 0 0.589 315 27.63 29.4772 0 0.648 350 23.05 32.2845 0 0.648 385 26.3 32.2845 0 0.648 367.5 24.67 32.2845 0 0.648 332.5 21.46 32.2845 0 0.648 315 19.9 32.2845 0 0.6185 350 26.81 30.8808 0 0.6185 385 30.57 30.8808 0 . 0.6185 367.5 28.69 30.8808 0 0.6185 332.5 24.98 30.8808 0 0.6185 315 23.18 30.8808 0 0.56 350 38.75 28.0735 0 0.56 385 44.12 28.0735 0 0.56 367.5 41.47 28.0735 0 0.56 332.5 36.19 28.0735 0 0.56 315 33.63 28.0735 0 0.53 350 50.3 26.6698 -40-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 0 0.53 385 57.14 26.6698 0 0.53 367.5 53.79 26.6698 0 0.53 332.5 47.06 26.6698 0 0.53 315 43.82 26.6698 4 0.589 350 50.44 26.7868 4 0.589 385 56.56 26.7868 4 0.589 367.5 53.39 26.7868 4 0.589 332.5 47.25 26.7868 4 0.589 315 43.9 26.7868 4 0.648 350 33.65 29.5006 4 0.648 385 37.77 29.5006 4 0.648 367.5 35.62 29.5006 4 0.648 332.5 31.35 29.5006 4 0.648 315 28.96 29.5006 4 0.6185 350 40.49 28.1437 4 0.6185 385 45.4 28.1437 4 0.6185 367.5 42.84 28.1437 4 0.6185 332.5 37.82 28.1437 4 0.6185 315 35.04 28.1437 4 0.56 350 67.66 25.4767 4 0.56 385 75.86 25.4767 4 0.56 367.5 71.74 25.4767 4 0.56 332.5 63.55 25.4767 4 0.56 315 59.32 25.4767 4 0.53 350 99.31 24.1432 4 0.53 385 116.6 24.1432 4 0.53 367.5 105.3 24.1432 4 0.53 332.5 92.92 24.1432 4 0.53 315 87.33 24.1432 2 0.589 350 37.58 28.7519 2 0.589 385 42.07 28.7519 2 0.589 367.5 39.7 28.7519 2 0.589 332.5 35.12 28.7519 2 0.589 315 32.5 28.7519 2 0.648 350 27.07 31.5359 2 0.648 385 30.41 31.5359 2 0.648 367.5 28.64 31.5359 2 0.648 332.5 25.16 31.5359 2 0.648 315 23.17 31.5359 2 0.6185 350 31.53 30.1556 2 0.6185 385 35.34 30.1556 2 0.6185 367.5 33.32 30.1556 2 0.6185 332.5 29.37 30.1556 2 0.6185 315 27.11 30.1556 2 0.56 350 46.05 27.3717 2 0.56 385 51.56 27.3717 2 0.56 367.5 48.72 27.3717 2 0.56 332.5 43.12 27.3717 2 0.56 315 40.03 27.3717 2 0.53 350 60.53 25.9914 2 0.53 385 67.81 25.9914 2 0.53 367.5 64.13 25.9914 2 0.53 332.5 56.83 25.9914 -41-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 2 0.53 315 51.71 25.9914 -2 0.589 350 37.58 28.7519 -2 0.589 385 42.07 28.7519 -2 0.589 367.5 39.7 28.7519 -2 0.589 332.5 35.12 28.7519 -2 0.589 315 32.5 28.7519 -2 0.648 350 27.07 31.5359 -2 0.648 385 30.41 31.5359 -2 0.648 367.5 28.64 31.5359 -2 0.648 332.5 25.16 31.5359 -2 0.648 315 23.17 31.5359 -2 0.6185 350 31.53 30.1556 -2 0.6185 385 35.34 30.1556 -2 0.6185 367.5 33.32 30.1556 -2 0.6185 332.5 29.37 30.1556 -2 0.6185 315 27.11 30.1556 -2 0.56 350 46.05 27.3717 -2 0.56 385 51.56 27.3717 -2 0.56 367.5 48.72 27.3717 -2 0.56 332.5 43.12 27.3717 -2 0.56 315 40.03 27.3717 -2 0.53 350 60.53 25.9914 -2 0.53 385 67.81 25.9914 -2 0.53 367.5 64.13 25.9914 -2 0.53 332.5 56.83 25.9914 -2 0.53 315 51.71 25.9914 -4 0.589 350 50.44 26.7868 -4 0.589 385 56.56 26.7868 -4 0.589 367.5 53.39 26.7868 -4 0.589 332.5 47.25 26.7868 -4 0.589 315 43.9 26.7868 -4 0.648 350 33.65 29.5006 -4 0.648 385 37.77 29.5006 -4 0.648 367.5 35.62 29.5006 -4 0.648 332.5 31.35 29.5006 -4 0.648 315 28.96 29.5006 -4 0.6185 350 40.49 28.1437 -4 0.6185 385 45.4 28.1437 -4 0.6185 367.5 42.84 28.1437 -4 0.6185 332.5 37.82 28.1437 -4 0.6185 315 35.04 28.1437 -4 0.56 350 67.66 25.4767 -4 0.56 385 75.86 25.4767 -4 0.56 367.5 71.74 25.4767 -4 0.56 332.5 63.55 25.4767 -4 0.56 315 59.32 25.4767 -4 0.53 350 99.31 24.1432 -4 0.53 385 116.6 24.1432 -4 0.53 367.5 105.3 24.1432 -4 0.53 332.5 92.92 24.1432 -4 0.53 315 87.33 24.1432 -42-Chapter 2 Probabilistic Design of Composite Plate Against Buckling Table 2.3: Statistics of random variables Variable Mean Coefficient of variation Distribution e 0.00° 0.50° (Standard Deviation) Normal (To be found) 0.20 Lognormal T (To be found) 0.10 Lognormal P (To be found) 0.25 Extreme Type I Table 2.4: Inverse reliability and optimum design parameters, EXAMPLE 1 Calculated Optimum Design Parameters Performance Requirement Target Reliability Index P1 Achieved Reliability Index p Vf =0.791 f =356.89 °F P = 16.39 N G\ = w - 20mm 1.645 1.652 G2 = 40mm - w 1.645 1.633 G3 = Pa-P 2.500 2.497 Table 2.5: Inverse reliability and optimum design parameters, EXAMPLE 2 Calculated Optimum Design Parameters Performance Requirement Target Reliability Index p1 Achieved Reliability Index p Vf =0.786 T = 350.50 °F G\=w- 20mm 1.645 1.632 Gj = 40mm - w 1.645 1.646 Gi = Pct-P 2.500 2.783 -43-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 2.11 FIGURES Figure 2.1: Basic structural elements (shell, plate and stiffener) of an aircraft. Temperature Autoclave Pressure Bag Pressure 800 700 600 500 A 400 9 300 V S-3 200 t/> at „^ 100 a. 0 -100 -200 50 100 200 250 300 150 Time (min) Figure 2.2: The air temperature and pressure cycle as applied in the autoclave during the processing of the composite plate (Twigg, 2001). -44-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 1 0 100 200 300 400 W,=L/2 (mm) Figure 2.3: Non-dimensional axial load versus transverse deflection curve of the plate with various initial imperfections. P~ is the critical Euler buckling load of the plate. u x Regions 1 -Aluminum Tool - Shear Layer - Composite Part 1.6 mm (1) 0.4 mm (2) 6 mm (3) 600 mm Figure 2.4: Finite element mesh of the plate and the tool used in COMPRO to virtually manufacture the part. Shear layer is used to model the tool - part interaction. -45-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 10 64 <u DC 08 & 4 • • • Experimental COMPRO • -<> 1 i r 1 1 t 200 400 600 800 Distance along part length (mm) 1000 1200 Figure 2.5: The warped profile of the 1200 mm long, 8 plies thick unidirectional composite part as it is removed from the aluminum tool. The prediction of COMPRO is compared with experimental results (Twigg, 2001). Fibre Orientation Figure 2.6: Schematic of the composite plate used for the case study. The applied loading and boundary conditions are also shown. -46-Chapter 2 Probabilistic Design of Composite Plate Against Buckling Figure 2.7: The superposition method to obtain the 3D state of stress from plane strain results for applied mechanical load (a) load P is applied in all three directions on a cube (b) 3D element: free to expand in all three directions (c) 2D plane strain element: load is applied in 1 and 2 directions, free to expand only in 1 and 2 direction, displacement in the 3 r d direction is constrained to be zero (d) load is applied only in the 3 r d direction, free to expand in all three directions. Combination of (c) and (d) gives (b). -47-Chapter 2 Probabilistic Design of Composite Plate Against Buckling (a) to + CO + 8, = ( l+v) l a A T CT, = r j 2 = 0 o \ = - E a A T Figure 2.8: The difference between 3D element and plane strain element (a) temperature increment is applied to a cube (b) 3D element: free to expand in all three directions, no stresses developed (c) 2D plane strain element: free to expand only in 1 and 2 directions, displacement in the 3 r d direction is constrained to be zero. -48-Chapter 2 Probabilistic Design of Composite Plate Against Buckling °"2 (c) (d) Figure 2.9: Superposing the plane strain stresses to obtain the 3D results (a) stresses from the plane strain analysis are applied as initial stresses to the cube (b) the cube is allowed to equilibrate under these initial stresses with free boundary conditions (c) the cube problem is modeled with a number of plane strain elements (d) the plane strain block in (c) is combined to form the full cube and allowed to equilibrate under the stresses. -49-Chapter 2 Probabilistic Design of Composite Plate Against Buckling Figure 2.10: Plate simply supported at the corners under thermal load (finite element mesh is shown for one quarter of the plate). (b) Figure 2.11: Developing the full 3D results using the results of a row of shell elements (a) a row of shell elements (b) full 3D plate. -50-Chapter 2 Probabilistic Design of Composite Plate Against Buckling b -100 + | -150 -200 Actual at y=0.025 Modified at y=0.025 Actual at y=0.975 Modified at y=0.975 0.2 0.4 0.6 Distance along x axis (m) (a) 0.8 50 0 a -so <§-100 I-150 -200 Actual at y=0.025 • Modified at y=0.025 - A c t u a l at y=0.975 Modified at y=fl.975 0.2 0.4 0.6 Distance along x axis (m) (b) 0.8 50 h-100 1-150 -200 Actual at y=0.025 »-Modified at y=0.025 -^Actual at y=0.975 Modified at y=0.975 0.2 0.4 0.6 Distance along x axis (m) (C) 0.8 Figure 2.12: Comparison of stresses from the actual 3D analysis and the modified superposition method (2D shell - 3D shell mapping) (a) stress o x (b) stress G y (c) stress Txy. -51-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 35 0 -4 1 1 1 ' 1 0 0.2 0.4 0.6 0.8 1 Distance along x axis (m) Figure 2.13: Comparison of vertical displacement from the actual 3D analysis and the modified superposition method (2D shell - 3D shell mapping). 1.0 m (500 Elements) (a) Figure 2.14: Mapping the stresses from plane strain solid element to shell elements, (a) The plate problem is modeled with plane strain solid elements and (b) the stresses from the plane strain analysis are mapped to the shell elements. -52-Chapter 2 Probabilistic Design of Composite Plate Against Buckling 2 0 (a) 1 0 0 -10 4 -20 -f -40 -50 *— Actual at y=0.025 •—Modified at y=0.025 Actual at y=0.975 Modified at y=0.975 0.2 0.4 0.6 Distance along x axis (m) (b) 0.8 Hlr- ^- 4h •—Actual aty=0.025 •^Modified at y=0.025 a— Actual at y=0.975 ^Modified at y=0.975 0.2 0.4 0.6 Distance along x axis (m) (c) 0.8 Figure 2.15: Comparison of stresses from the actual 3D analysis and the modified superposition method (2D solid - 2D shell - 3D shell mapping) (a) stress o x (b) stress rjy (c) stress Tjy. -53-Chapter 2 Probabilistic Design of Composite Plate Against Buckling Figure 2.16: Comparison of vertical displacement from the actual 3D analysis and the modified superposition method (2D solid - 2D shell - 3D shell mapping). T -Temperature Vf - Fibre volume fraction 6 - Ply misalignment HR- Heating rate P - Pressure HT - Hold time CR - Cooling rate Vf 0 HR HT CR Figure 2.17: Sensitivity study of the processing parameters that mostly affect the warpage of the part. Sensitivity is defined as the ratio of the percentage increase in warpage to the percentage increase in the parameter. -54-C H A P T E R 3: SPRING-IN OF A N G L E D COMPOSITE L A M I N A T E 3.1 INTRODUCTION In the previous chapter, it was demonstrated using a simple example, how the process induced deformations can affect the structural response of composite structures. In that study, only a limited number of processing parameters were considered for simplicity. In this chapter, we will consider a more complex shape laminate, an angle section, where an attempt will be made to identify and separate the various phenomena that give rise to process-induced residual stresses and deformations. Stresses develop at two different levels during the processing of composite structures: micro-level and macro-level. Macro-mechanical stresses originate at the ply- and laminate-level. These stresses are the source for large-scale deformation of the processed part and may lead to premature failure of the structure during its service life. In most processing models, the macro-mechanical stresses are calculated by treating each lamina as a homogeneous but anisotropic layer. The effective material properties of each lamina are determined either using micro-mechanics equations or micro-mechanical finite element models based on the properties of the resin and fibre and their volume fraction at a given temperature and degree of cure of resin. A number of sources, such as tool-part interaction, anisotropic thermal strain, and cure shrinkage of resin, contribute to these stresses during the autoclave processing. The next chapter deals with the micro-level stresses in more detail. -55-Chapter 3 Spring-in of Angled Composite Laminate The following are identified as the important sources that contribute to the process-induced deformations and stresses (Johnston, 1997): 1. Anisotropic thermal expansion ""i > Anisotropy 2. Chemical shrinkage of resin J 3. Tool-part interaction 4. Lay-up 5. Fibre volume fraction gradient through thickness The spring-in of a curved composite part can be separated into two main components: corner component and warpage component. Among the sources stated above, anisotropy mainly contributes to the corner component while tool-part interaction contributes to the warpage component (Albert, 2001). Many parameters are reported to affect the spring-in. These parameters can be separated into two groups: intrinsic and extrinsic. Intrinsic parameters are related to part geometry and material properties such as initial part angle, corner radius, thickness, lay-up, staking sequence and part shape and extrinsic parameters are related to tooling and processing such as cure cycle, tool surface condition, tool material, cure temperature, and pressure. Over the past three decades, significant research has been conducted to either experimentally measure or model the effect of these parameters on spring-in. On the modelling front, a two-dimensional finite element computer code, COMPRO, was developed as mentioned earlier. Although the predictive capability of COMPRO has been investigated earlier, not every feature in the code has been validated against a carefully designed set of experiments. Recently Albert -56-Chapter 3 Spring-in of Angled Composite Laminate (2001) carried out an extensive experimental program to study the effect of some of the above-mentioned parameters on dimensional stability of various shape composite laminates. In this study, an attempt is made to capture numerically some of the process-induced dimensional changes measured by Albert (2001) using COMPRO. 3.2 REVIEW OF EXPERIMENTS In this section, a brief review of the experimental work carried out by Albert will be given. This review includes only the experiments that are relevant to our studies. A more detailed description of the experiments and the results can be found in Albert (2001). One of the tool geometries used in the experiments is shown in Figure 3.1. The tool surfaces were cleaned with acetone to remove any trace of oil or dirt and two or three coats of release agent (Multishield) were applied over the entire surface to prevent sticking of the part to the tool. In some cases, to reduce the friction between the tool and the part, an FEP (Fluorinated Ethylene Propylene) sheet was placed over the release agent before the part was laid up on the tool. The parts were C- or L-shaped with flange length of 57 mm and thickness consisting of either 8 plies (1.6 mm) or 16 plies (3.2 mm). The web length of the C-shaped part was 102 mm. The parts made were either unidirectional laminates [0]n, where all the fibres follow the curved contours, or quasi-isotropic laminates [0/45/-45/90]ns with 0° and 90° directions shown in Figure 3.2. Once the lay-up was complete, the tool-part assemblies were placed inside an autoclave and cured under a prescribed temperature and pressure cycle (both one-hold and two-hold cure cycles were used). After the part was removed from the tool, approximately 1 cm was trimmed off the -57-Chapter 3 Spring-in of Angled Composite Laminate edge to eliminate the edge thinning effect. The surfaces were polished and scanned. The scanned images were then digitally analyzed. 3.3 PROCESS MODELLING USING C O M P R O In this case study, the stresses and the deformations of an L-shaped composite component are numerically analysed and compared with experimental results of Albert (2001). The composite material used in the modelling is different from the material used in the actual experiments since the complete material data for the latter was not available. Hence, in this study, only the trends, rather than the absolute values, are compared with the experiments. The modelled composite component is made of unidirectional layers of AS4/8552 material. The thermal and mechanical properties of this material are given in Appendix A (Johnston, 1997). The plies are laid up on a tool and virtually processed using COMPRO. The part and the tool dimensions are shown in Figure 3.3. The part and the tool are modeled in COMPRO using 8-noded isoparametric elements (Osooly, 2001) as shown in Figure 3.4. The interaction between the part and the tool is modeled by adding a shear layer in between as indicated in Chapter 1. The mechanical and thermal boundary conditions are also shown in Figure 3.4. The 16-ply unidirectional laminate of AS4/8552 material on aluminum tooling is chosen as the nominal case for this study, although other lay-ups and thicknesses are also examined. The part is cured under the cure cycles (both one and two holds) shown in Figure 3.5. After the simulation of tool removal, the spring-in angle is calculated. The spring-in angle is calculated using a technique analogous to that used in the experiment. The total spring-in is calculated by calculating the included angle of the lines connecting the two end points (one at the corner and one at the end) of the flanges as shown in Figure 3.6. The corner component of this total spring--58-Chapter 3 Spring-in of Angled Composite Laminate in is obtained by calculating the included angle between the tangent lines at the corner (see Figure 3.6) and the warpage component is given by the gradient of the warpage profile of the flange at the corner as shown in Figure 3.7. 3.4 SOURCES OF RESIDUAL STRESSES 3.4.1 ANISOTROPY When an L-shaped composite part is cured in the autoclave without a tool, the part will spring in due to two main sources: anisotropy in the coefficient of thermal expansion (CTE) and cure shrinkage. The carbon fibre has a very low longitudinal CTE compared to its transverse CTE and the CTE of the matrix (Appendix A). This results in the laminate having a higher out-of-plane thermal strain than the in-plane thermal strain during curing. Similarly, the chemical shrinkage of the matrix is constrained by stiff fibres in the longitudinal direction, leading to higher chemical shrinkage strains in the out-of-plane direction than in the in-plane direction. This difference between the in-plane and out-of-plane strains leads to the spring in/up of a curved composite part. This is purely driven by the circular geometry of the corner and mismatch of the in-plane and out-of-plane CTE and chemical shrinkage properties of the part. It is independent of the thickness and the mechanical properties of the part. A simple equation has been proposed for predicting the spring in/up of a curved laminate due to anisotropy (Nelson et al., 1989 and Radford et al, 1993): A0 = A0CRE+A0CS =6 (CTEL-CTET)AT)+0(</>L-^ \+CTER AT (3.1) where 9= initial part angle (°) A0= spring-in angle (°) -59-Chapter 3 Spring-in of Angled Composite Laminate AOCTE = component of spring-in due to thermal expansion anisotropy (°) AOcs = component of spring-in due to cure shrinkage anisotropy (°) AT = change in temperature (°C) </>L - longitudinal chemical shrinkage strain (-) <jtr = through-thickness chemical shrinkage strain (-) CTEi = longitudinal coefficient of thermal expansion (m/m/°C) CTE, = through-thickness coefficient of thermal expansion (m/m/°C) The volumetric cure shrinkage of Hexcel 8552 epoxy resin is 0.099 (Johnston, 1997). Converting the volumetric shrinkage strain, Vf, to linear strain, £-;s, using the following equation leads to an isotropic strain of 3.2%. The lamina longitudinal and transverse chemical shrinkage strains are calculated using the micromechanics equations (Equations B . l l and B. l2 in Appendix B) proposed by Hashin and Rosen (1964) by setting the fibre strains equal to zero. The lamina longitudinal and transverse shrinkage strains are computed to be 5.217E-04 and 1.858E-02, respectively. During the curing process, according to Equation 3.1, the anisotropic coefficient of thermal expansion does not contribute to the spring-in and the cure shrinkage results in a spring-in of 1.6°. In COMPRO, the 16-ply part without the tool is virtually processed and the predicted spring-in angle is 2.61°. The same spring-in is predicted for 8-ply and 24-ply parts as well. The higher value is due to the plane strain assumption in COMPRO. When the same run is repeated without (3.2) -60-Chapter 3 Spring-in of Angled Composite Laminate the Poisson's effect (i.e. when Poisson's ratios of the fibre and the matrix are set to zero), the prediction is 1.44°, which is very close to the analytical solution (Equation 3.1). The error due to the plane strain assumption is also investigated using the commercial software A B A Q U S for the 16-layer part by comparing the spring-in prediction by both plane strain analysis and three-dimensional analysis. In the analysis, only the spring-in due to CTE-anisotropy during the cooling down portion of the cure cycle is considered by assuming the material is fully cured and the temperature is uniform in the part. In the three-dimensional case, only one element of length 2 mm is modelled in the out-of-plane direction. From the results shown in Figure 3.8, the plane strain assumption over-predicts the spring-in for the unidirectional and [0/90] lay-ups and predicts accurately the spring-in for the quasi-isotropic lay-up. This trend can be explained by calculating the equivalent mechanical properties of the above lay-ups. Since the classical laminate plate theory does not account for the through thickness properties, a simple numerical model outlined below is used to predict the equivalent properties. Numerical model to predict the equivalent properties of lay-ups Consider a composite cube of 16 layers as shown in Figure 3.9(a). The cube is modelled in A B A Q U S with a single 20-noded isoparametric solid composite element. The mechanical properties of a lamina are calculated from the fully cured properties of the resin and fibre using the micromechanics equations B . l - B.12 in Appendix B. The stacking sequence of the lay-up starts from face-3a in the 3 r d direction. The boundary conditions applied are (Figure 3.9(b)): Face-la = Symmetric Face-lb = Coupled -61-Chapter 3 Spring-in of Angled Composite Laminate Face-2a = Symmetric Face-2b = Coupled Face-3a = Symmetric Face-3b = Coupled The symmetric boundary condition implies that the deformations of the nodes along that particular face in the direction normal to the face are constrained to be zero. The coupled boundary condition implies that all the nodes of that particular face are coupled to deform by the same amount in the direction normal to the face (this is achieved using the EQUATION function in ABAQUS). To calculate the Young's modulus in the 1 direction (£/), a constant displacement (0.1 units) is applied in the 1 direction to the face-lb and the resultant force in the same direction at the face-la, the resultant deformations of the face-2b and face-3b are calculated. By dividing the resultant force by the applied displacement, we obtain E[. i.e. Ei = (Resultant force in direction l)/(Applied displacement in direction 1) The Poisson's ratios are given by V12 = -(Resultant displacement of face-2b)/(Applied displacement in direction 1) vj3 = -(Resultant displacement of face-Sb)/(Applied displacement in direction 1) The other mechanical properties are also calculated in a similar manner. The coefficients of thermal expansions are calculated by applying a known temperature (100 °C) and calculating the resultant deformations of faces lb, 2b, and 3b in the direction normal to the respective faces. CTE/ = (Displacement offace-lb)/(Applied temperature) -62-Chapter 3 Spring-in of Angled Composite Laminate The equivalent properties so calculated are listed in Table 3.1. The in-plane properties are seen to agree very well with those obtained from CLT (the values obtained from CLT are listed inside a square bracket next to the corresponding values from the above method in Table 3.1). Difference between 3D and plane strain In the 3D case, under a temperature change, the curved composite part can expand freely in the out-of plane direction (direction 2). But in the plane strain case, the out-of plane strain (ei) is fully constrained. This is equivalent to applying the negative of this out-of plane strain to a 3D part. Hence the plane strain condition induces the following additional in-plane strains due to the Poisson's effect. £3 ~ V23£2 (3.3) If the material is isotropic then V21 = V23 and the plane strain condition does not change the included angle of the part. However, for orthotropic materials V21 £ V23 and the plane strain condition changes the included angle according to Equation 3.1. Afl>"-V*k (3.4) {\ + V23£2) Hence, the magnitude of the difference in spring-in between 3D and plane strain depends on the difference between the two Poisson's ratios, namely, (V21 - V23), and 82. For the unidirectional lay-up, both quantities are high and therefore the difference between the 3D and plane strain prediction of spring-in is high. For the [0/90]n lay-up, though (v?/ - V23) is high, 82 is small and therefore the difference is smaller. For the quasi-isotropic lay-up, both quantities are very small (see table 3.1) and so is the difference between 3D and plane strain predictions. -63-Chapter 3 Spring-in of Angled Composite Laminate 3.4.2 TOOL-PART INTERACTION The presence of the tool affects the distortion of the part in two ways: thermally and mechanically. The thermal effect is due to the large thermal mass of the tool that induces a temperature gradient through the part thickness (Figure 3.10). During processing, the tool mechanically interacts with the curved part in two modes: constrains the part from springing in and introduces a shearing force at the tool-part interface due to the CTE mismatch between the part and the tool. Experimental study of tool-part interaction The tool-part interaction of a flat plate is discussed in detail by Twigg (2001). In his work, Twigg (2001) proposed an empirical equation (Equation 3.5) to predict the variation of part warpage with part thickness and length for the parts processed under similar processing conditions. Maximum part warpage, w m a x oc — (3.5) W m a x = C ^ where L is the half-length of part, t is the laminate thickness, and C is a proportionality constant. By assuming first ply slip and constant shear stress along the part length, he proposed that this proportionality constant can be expressed as C = (3.6) E pari where rnel = (Tmterface - Tinterpfy) as shown in Figure 3.11 and Epar, is the elastic modulus of the part. From the average of his experimental results, he deduced a value of 6.2 kPa for T n e t - Based -64-Chapter 3 Spring-in of Angled Composite Laminate on this value and the smaller lengths of the parts considered, the warpage component of the parts studied by Albert (2001) should be negligible (both Albert and Twigg used similar processing conditions and materials). However, Albert observed considerable warpage in her experiments. This may be attributed to the tool surface condition and the size of the part in both experiments. As reported by Twigg (2001), two conditions may prevail at the tool-part interface: sticking and sliding as shown in Figure 3.12. For a smaller part, the sticking condition may be dominant. This provides a rough explanation for the increased warpage in Albert's experiment. More experimental studies have to be done to verify this. Modelling the tool-part interaction As mentioned earlier, the tool-part interaction in the current version of COMPRO is simulated by introducing a shear layer. The introduction of an isotropic shear layer with very low Young's (and shear) modulus removes both the normal constraint and the tool-part interface shear effect. However, the introduction of an orthotropic shear layer with very low Young's (and shear) modulus only in the direction parallel to the tool-part interface removes the tool-part interface shear effect but maintains the constraining effect of the tool normal to the interface. The values of En and Gi$ of the shear layer are lowered by identical orders of magnitude. In the following discussion, the shear layer properties are quoted by the value of En. The order of En dictates the order of G/j. The usual practice is to select the shear layer properties to match the results of an experiment and then the same shear layer properties are used to predict the trend of the variation with other geometric parameters such as part thickness, part length etc. To satisfy Equation 3.5, the shear stress should be constant with these geometric parameters. To verify this, the shear stress in the shear layer between a flat part and an aluminum tool, similar to the model of Twigg (2001), is investigated at a certain time (200 min) during the cure cycle (Figure 3.13-3.14). -65-Chapter 3 Spring-in of Angled Composite Laminate From these figures, it can be seen that the average interfacial shear stress remains constant with part thickness, which complies with Equation 3.5. But the average shear stress does not remain constant with part length. These graphs confirm the finding of shear lag theory that the gradient of the interfacial shear stress curve depends on the shear layer properties. Hence, by scaling the shear layer property by the ratio of the part length, the average shear stress can be maintained constant with part length as shown in Figure 3.15. In this figure, the shear layer properties of the 600 mm part are kept constant and the shear layer properties of the 300 mm and 1200 mm parts are modified by the length ratio. The variation of the maximum warpage of a flat plate with part thickness and part length is shown in Figure 3.16 - 3.17. In these figures, prediction of the COMPRO is compared with the analytical solution using Equation 3.5. In the analytical solution, the warpage of the 8-layer, 600 mm part is kept constant and the variation of the warpage with the part thickness and part length is calculated using Equation 3.5. By introducing an isotropic shear layer with very low Young's and shear moduli, the mechanical interaction of the tool with the part is essentially eliminated while keeping the thermal mass effect of the tool, which introduces a temperature gradient through the thickness. The spring-in from this analysis was compared with the predicted spring-in value of the part without the tool. Only a small variation was observed due to the difference in the through thickness temperature gradient between the two cases. And the parts with three different thicknesses (1.6 mm, 3.2 mm and 4.8 mm) resulted in almost the same value of spring-in (-2.61°). This shows that the temperature gradient does not have a strong influence on the spring-in due to anisotropy, for this case study. It may have some influence on a thicker part. During processing, when the part is constrained from springing in by the tool, stresses build up in the part that cause it to spring in once the tool is removed. The predicted spring-in of 8-ply, 16--66-Chapter 3 Spring-in of Angled Composite Laminate ply and 24-ply parts when the parts are constrained by the tool are given in Table 3.2 and Figure 3.18. According to this result, the contribution of the cure shrinkage to the total spring-in is small (-20%) compared to the contribution from CTE mismatch (-80%). However, when the part is not constrained, cure shrinkage anisotropy is the major source of springing in. This can be explained by observing the stress development through the cure cycle. Stress starts to develop significantly only after the vitrification point (mostly on the cooling down part). But much of the resin shrinkage occurs while the resin is in the rubbery state and therefore without any significant stress development. The shearing effect of the tool is investigated by changing the shear layer properties from no bonding to full bonding (shear layer = tool). The labels used for different shear layers and their material properties are listed in Table 3.3. In the following discussion, the shear layer will be referred to by the labels. The spring-in of the part with different shear layer properties is shown in Figure 3.19. From the figure, it can be observed that only the changes in warpage component with shearing layer properties is significant while the change in corner component is very small. To compare the modelling results with experiment, first an appropriate shear layer property should be selected to match one of the experimental results. Since in the experiment two types of surface conditions were used (Release agent (3 coats), Release agent (3 coats) + FEP), two different shear layers are selected to represent these two surface conditions. A comparison between the predicted and measured results for the spring-in of an 8-layer part under the 2-hold cure cycle, as shown in Figure 3.20, reveals that the selected shear layer properties, namely TS-1 for RA+FEP and TS-2 for R A only, are reasonable choices. The properties of the selected shear layers are listed in Table 3.4. In the following sections, when COMPRO results are compared to the experimental results, the same similarity will be maintained (i.e. RA+FEP = TS-1, R A = TS--67-Chapter 3 Spring-in of Angled Composite Laminate 2). The results for both 8 layers and 16 layers, for the two tool surface conditions, are compared in Figure 3.21 and the parts details are listed in Table 3.5. The question that comes to mind from this result is whether it is possible to conclude that anisotropy (CTE mismatch and cure shrinkage) contributes to the corner component of the spring-in while tool-part interaction (shear) contributes to the warpage component of the spring-in. If anisotropy is the sole contributor to the corner component, then the corner component should not vary with thickness, radius, tool surface condition etc. according to Equation 3.1. But according to Figure 3.21, the corner component does change with tool surface condition and part thickness. This trend can be explained by having a close look at the method of separating the warpage and the corner component. In our measurement method, the separation of warpage and corner components is done at the corner-flange interface. Hence, the warpage component only includes the warpage of the flange. It does not include the warpage of the corner as shown in Figure 3.22. Hence, depending on the curvature, the corner component may consist of the spring-in due to anisotropy and the spring-in due to warpage of the corner length. This warpage portion of the corner component may be the reason the corner component of the spring-in changes with part thickness, tool surface, etc. 3.5 EFFECT OF INTRINSIC AND EXTRINSIC PARAMETERS In this section, the effect of the intrinsic and the extrinsic parameters are numerically predicted by modelling and compared with experiments. Since the material used in the modelling is different from the material used in the experiments and the material data are incomplete, no attempt is made to compare the absolute values. Instead only the trend of variation of the spring-in with the parameters is compared. In comparing the effect of a given parameter between two -68-Chapter 3 Spring-in of Angled Composite Laminate parts, only that parameter is varied, while the rest of the parameters are kept identical. For example, i f the purpose of the analysis is to compare the effect of part thickness, the trend of spring-in variation between two L-shaped parts predicted by the model may be compared with the trend of spring-in variation between two C-shaped parts from the experiment. 3.5.1 EFFECT OF PART THICKNESS The predicted variation of corner and warpage components of the spring-in with part thickness is shown in Figure 3.23. According to Equation 3.1, spring-in due to anisotropy does not change with part thickness. But according to Equation 3.5, the spring-in due to tool-part interaction decreases with part thickness. However, both corner and warpage components are expected to decrease with increasing part thickness as shown in the Figure 3.23. A similar trend was also observed in the experiments as shown in Figure 3.24. In the experiments, while the tool is aluminum, the parts are C-shaped with both 8 and 16 unidirectional layers processed under a 2-hold cure cycle (Table 3.6). 3.5.2 EFFECT OF TOOL MATERIAL In this section, the variations of spring-in with two different tool materials (aluminum and steel) are analysed. Similar to the part thickness, the tool material also does not change the spring-in due to anisotropy. Since the thermal properties of the tool materials are different, the temperature gradient through the part thickness may vary. But it has already been shown that the temperature gradient does not have significant affect the spring-in (Section 3.4.2). However, the warpage component may vary with the tool materials due to the difference in the CTE values. Since aluminum has a higher CTE value than steel, it causes a higher CTE mismatch between the tool -69-Chapter 3 Spring-in of Angled Composite Laminate and the part. Therefore, it is expected that the warpage component for the case of aluminum tool are higher. The comparison of COMRO predictions with the experimental results is shown in Figure 3.25 for 8 layers and Figure 3.26 for 16 layers. In the experiment C-shaped parts of unidirectional layers were processed under a 2-hold cure cycle (Table 3.7 - Table 3.8). In the modelling, selecting the shear layer properties for both tooling materials to match the respective tooling surface condition is not clear. For example SL-2 is selected to represent the R A tool surface for the aluminum tooling. But it is not verified whether the same shear layer can be used to represent the R A tool surface for the steel tooling. Hence, in this case, the mechanical properties of the shear layer are kept constant. But the thermal and mass properties are made to match the properties of the tooling material. In the Figures 3.25 and 3.26, both the experimental and modelling results show the same trend for the corner component. But the experiments show an opposite trend for the warpage component except for the 8-layer part with R A tool surface. 3.5.3 EFFECT OF LAY-UP In this section, the effect of two different lay-ups (unidirectional and quasi-isotropic) on the spring-in are investigated and compared with experimental results. As discussed in Section 3.4.1, the lay-up changes the spring-in due to anisotropy by altering the in-plane and the through thickness properties. According to Equation 3.1, the quasi-isotropic lay-up gives more spring-in due to anisotropy than the unidirectional lay-up (Figure 3.8). But it has already been shown that COMPRO may predict the opposite trend due to the plane strain assumption (Section 3.4.1). -70-Chapter 3 Spring-in of Angled Composite Laminate Similar to anisotropy, it is reported in the literature, that the quasi-isotropic lay-up gives more spring-in due to tool-part interaction than the unidirectional lay-up (Albert, 2001). The reason behind this is given below (the mathematical derivation to support this argument is given in Appendix D). During the heat up portion of the cure cycle, the resin has very low modulus. Hence, the stress due to tool-part interaction mostly occurs at the first ply level and it reduces drastically with the distance away from the first ply. According to this, both unidirectional lay-up and the quasi-isotropic lay-up with 0-ply at the bottom have similar stress distribution during heat up. On the other hand, during the cool down portion of the cure cycle, the resin is fully cured and therefore has a higher modulus. Hence, the whole thickness of the part is involved in the interaction with the tool and the stress due to the tool-part interaction is almost constant through the thickness. The magnitude of this stress depends on the CTE and the axial stiffness of the part and it is different for both lay-ups. On the removal of the tool, the resultant of the stresses at these two stages (heat up and cool down) causes the part to warp. The curvature of the part is given by d2w M M —T = — 0 0 — (3-7) dx2 EI E Therefore, the warpage component of the spring-in should be in the following order ^[0/45/-45/90]2.s- > ^[0]4.v > ^[90/45/-45/0]4s (3-8) The prediction by COMPRO of the warpage of a 600 mm long and 16 plies thick flat part processed on a aluminum tool under the fully bonded tool surface condition shows a different trend, W[0]4.s > W[0/45/-45/90]2.v > W[90/45/-45/0]4.v (3-9) -71-Chapter 3 Spring-in of Angled Composite Laminate The above argument is totally based on the first ply slip, which, in turn, is very sensitive to the initial resin modulus. Hence, warpage of the same flat part is analysed for different initial resin moduli as listed in Table 3.9. It shows that the expected trend is predicted correctly at low resin modulus. This leads to the same conclusion made by Twigg (2001) that the initial resin modulus should be obtained experimentally rather than selected arbitrarily as is done in currently (Johnston, 1997). The effect of the lay-up on the spring-in is shown in Figure 3.27. From this figure, it can be observed that at low shear layer moduli, the corner component of the unidirectional lay-up is higher than the quasi-isotropic lay-up, but the warpage component is consistently lower. The comparison between the experimental results and the COMPRO prediction is shown in Figure 3.28. In the experiments, the parts are L-shaped and made of 8 layers. The flange length of the parts is 89 mm (greater than the model) and the included angle is 45° (less than the model). They are processed on an aluminum tool with R A surface condition under a 2-hold cure cycle (Table 3.10). Good agreement is observed between the experimental and modelling results when the trends in behaviour are compared. 3.5.4 EFFECT OF PART SHAPE In this section, the effect of part shape or geometry on the spring-in is examined. In particular, L and C shapes are compared. The C-shaped part geometry is shown in Figure 3.29. The web length is twice the flange length similar to the experiment. The variation of spring-in between the two shapes predicted by the model is shown in Figure 3.30. The corner components of both the L and C shapes are comparable. But the predicted warpage component of the C-shaped part is much larger than that of the L-shaped part with a higher proportion of the warpage occurring at -72-Chapter 3 Spring-in of Angled Composite Laminate the web. In contrast to the above prediction, the experiment shows no difference between a C and an L shaped part. Note that Albert (2001) only reported the-total warpage component. When the experimental warpage component is actually separated into the warpages of the flange and the web, it is observed that most of the warpage occurs in the web section as shown in Figure 3.31 (the flange warpages are nearly zero). In a previous combined experimental and numerical work (Fernlund et al., 2001), the C-shaped part had a greater spring-in than the L-shaped part by about 20%. For the numerical study, a different version of COMPRO with updated geometry capability was used. The COMPRO version used in this current study does not have the updated geometry option and the prediction for the spring-in of the C-shaped part is much higher than the L-shaped part even if the tool part interaction is switched off by using a very low modulus shear layer. This indicates that the corner of the C-shaped part is prevented from sliding over the tool (corner locking) and, because of that, when the tool expands it stretches the web with it as shown in Figure 3.32. This suggests that, when a tool with an effectively zero CTE, such as Invar, is used, then both the L and C shapes should result in identical spring-in values. In fact, the analysis using COMPRO shows that the spring-in of the C shaped parts on an Invar tool is only 6% higher than the spring-in of the L-shaped part on the same tool. Currently research work is underway to replace the shear layer in COMPRO with a more realistic and robust contact interface model with updated geometry capability. In the meantime, a crude method is introduced in this section to overcome the corner-locking problem. In this method, the shear layer is divided into three regions (web, corner and flange) as shown in Figure 3.33. The web and the flange shear layers have the usual shear layer properties (example En = 6.95E5 Pa). To eliminate the corner-locking, a gap is introduced in the corner (only for the C--73-Chapter 3 Spring-in of Angled Composite Laminate shaped part) by using an isotropic shear layer with very low Young's (and shear) modulus. The results of this modified method is compared with the earlier results for the C-shaped part in Figure 3.34. It can be seen that with this modification both the C-shaped and L-shaped parts show almost the same amount of spring-in. The result is also compared with the experimental result in Figure 3.35 and Figure 3.36. In the experiment, the parts are 8 layers thick and they are processed on both aluminum and steel tool with RA(3) surface condition under a 2-hold cure cycle (Table 3.11 - Table 3.12). A good agreement between the experimental results and the predicted results is observed. Hence, this crude method may serve as a temporary solution to overcome the corner-locking problem in COMPRO until a more scientifically sound method such as the contact surface method is developed. 3.6 SUMMARY AND CONCLUSIONS The accuracy of the reliability design method introduced in Chapter 1 mostly depends on the prediction capability of COMPRO. The key issue in the modelling is to model the tool-part interaction correctly. Extensive research in both experimental and numerical areas has been conducted to model the tool-part interaction with reasonable accuracy in COMPRO. Currently, the tool-part interaction is modeled in COMPRO by introducing a thin shear layer. In this work, the capability of COMPRO to predict the variation of spring-in of an angled part with intrinsic and extrinsic parameters is verified. For this purpose, COMPRO results are compared with the experimental results of Albert (2001). Unfortunately, due to the unavailability of the material data for the composite material used in the experiment, only the trend of variation was examined using the data of a similar material. Although the shear layer method is not a physically accurate representation of the tool-part interaction, it can be used to capture some -74-Chapter 3 Spring-in of Angled Composite Laminate important phenomena during processing. The method proposed here to overcome the corner-locking in a C-shaped part by introducing a gap in the corner should be verified with more experimental data. -75-Chapter 3 Spring-in of Angled Composite Laminate 3.7 TABLES Table 3.1: Evaluated laminate equivalent properties (the values in [ ] are those obtained using C L T ) . Property [0]l6 [0/90]4s [0745/-45/90]2s E, (GPa) 122.800 [122.8] 66.664 [66.64] 48.492 [48.47] E2 (Gpa) 9.920 [9.920] 66.664 [66.64] 48.492 [48.47] E3 (GPa) 9.920 12.092 12.092 Vl2 0.2680 [0.2680] 0.0397 [0.0401] 0.3015 [0.3018] V 1 3 0.2680 0.4207 0.3060 V23 0.4710 0.4207 0.3060 V21 0.0216 0.0397 0.3015 V31 0.0216 0.0763 0.0763 V32 0.4710 0.0763 0.0763 CTE, (um/mm/C) 0.600 [0.600] 3.1515 [3.151] 3.1515 [3.151] C T E 2 (um/mm/C) 28.600 [28.600] 3.1515 [3.151] 3.1515 [3.151] C T E 3 (um/mm/C) 28.600 39.793 39.793 Table 3.2: Spring-in of the L-shaped composite parts when constrained by the tool. Part Lay-up Total Spring-in (deg) C T E Mismatch (%) Cure Shrinkage (%) 8 0.750 72.17 27.83 16 0.659 86.34 13.66 24 0.634 90.88 9.12 Table 3.3: Shear layer properties and their identification labels. Shear Layer Labels Longitudinal Modulus E n (Pa) Shear Modulus G , 3 (Pa) En , s h e a r L a y e r / E n , t o o l SL-0 6.9*10"J 2.6*10"J 10"u SL-1 6.9*103 2.6*103 IO"3 SL-2 6.9*10b 2.6* 10b IO"4 SL-3 6.9*10' 2.6*10' 10"j SL-4 6.9*10* 2.6*10s 10"2 SL-5 6.9*101U 2.6*101U 1 -76-Chapter 3 Spring-in of Angled Composite Laminate Table 3.4: Selected shear layer properties to represent the experimental tool surface conditions. Tool Surface COMPRO Experiment Shear Layer Properties Warpage Component of Spring-in (deg) Tool Surface Condition Warpage Component of Spring-in (deg) TS-1 En =6.9E4 Pa 0.135 Release Agent (3 coats) + FEP 0.14 G i 3 = 2.6E4 Pa TS-2 En =7.6E6 Pa 0.895 Release Agent (3 coats) 0.92 G 1 3 = 2.9E6 Pa Table 3.5: Spring-in variation with tool surface conditions. COMPRO predictions are compared with experimental results. Part Shape Lay-up Tool Material Tool Surface Cure Cycle Spring-in Corner Warpage COMPRO L [0]8 Aluminum TS-1 2 0.79 0.14 L [0]8 Aluminum TS-2 2 0.95 0.90 L [0]l6 Aluminum TS-1 2 0.73 0.08 L [0]l6 Aluminum TS-2 2 0.75 0.33 Experiment C [0]8 Aluminum RA(3)+FEP 2 0.78 0.14 C rois Aluminum RA(3) 2 0.97 0.92 c [0],6 Aluminum RA(3)+FEP 2 0.72 0.07 c [01,6 Aluminum RA(3) 2 0.83 0.28 -77-Chapter 3 Spring-in of Angled Composite Laminate Table 3.6: Spring-in variation with part thickness. COMPRO predictions are compared with experimental results. Part Shape Lay-up Tool Material Tool Surface Cure Cycle Spring-in Corner Warpage COMPRO L rois Aluminum TS-1 2 0.79 0.14 L [0],6 Aluminum TS-1 2 0.73 0.08 L [0]g Aluminum TS-2 2 0.95 0.90 L [01,6 Aluminum TS-2 2 0.75 0.33 Experiment C [0]8 Aluminum RA(3)+FEP 2 0.78 0.14 C [0],6 Aluminum RA(3)+FEP 2 0.72 0.07 c [0]8 Aluminum RA(3) 2 0.97 0.92 c [0],6 Aluminum RA(3) 2 0.83 0.28 Table 3.7: Spring-in variation with tool material for an 8-layer unidirectional part. COMPRO predictions are compared with experimental results. Part Shape Lay-up Tool Material Tool Surface Cure Cycle Spring-in Corner Warpage COMPRO L [0]8 Aluminum TS-1 2 0.79 0.14 L [0]8 Steel TS-1 2 0.79 0.13 L [0]8 Aluminum TS-2 2 0.95 0.90 L [0]8 Steel TS-2 2 0.87 0.50 Experiment C [0]8 Aluminum RA(3)+FEP 2 0.78 0.14 C [0]8 Steel RA(3)+FEP 2 0.65 0.19 c [0]8 Aluminum RA(3) 2 0.97 0.92 c [0]8 Steel RA(3) 2 0.81 0.55 -78-Chapter 3 Spring-in of Angled Composite Laminate Table 3.8: Spring-in variation with tool material for a 16-layer unidirectional part. COMPRO predictions are compared with experimental results. Part Shape Lay-up Tool Material Tool Surface Cure Cycle Spring-in Corner Warpage COMPRO L [0]l6 Aluminum TS-1 2 0.73 0.08 L [0]l6 Steel TS-1 2 0.73 0.07 L [0]l6 Aluminum TS-2 2 0.75 0.33 L [0]l6 Steel TS-2 2 0.75 0.20 Experiment C [0]l6 Aluminum RA(3)+FEP 2 0.72 0.07 C [0]l6 Steel RA(3)+FEP 2 0.60 0.18 c [0]l6 Aluminum RA(3) 2 0.83 0.28 c [0]l6 Steel RA(3) 2 0.77 0.38 Table 3.9: Variation of maximum warpage of a flat part (600mm long, 16 plies thick) with lay-up. Initial Resin Modulus (Pa) Maximum Warpage (mm) [01,6 [0,45,-45,90]2S [90,45,-45,0]2S 4.67E6 38.3 37.2 15.8 4.67E5 28.8 45.4 22.1 4.67E4 17.9 27.5 9.8 Table 3.10: Spring-in variation with part lay-up. COMPRO predictions are compared with experimental results. Part Shape Lay-up Tool Material Tool Surface Cure Cycle Spring-in Corner Warpage COMPRO L [01,6 Aluminum TS-2 2 0.95 0.90 L [0,45,-45,90]s Aluminum TS-2 2 1.01 1.23 L [0]l6 Steel TS-1 2 0.73 0.08 L [0,45,-45,90]2S Steel TS-1 2 0.69 0.10 Experiment L [0]g Aluminum RA(3) 2 0.81 0.71 L [0,45,-45,90]s Aluminum RA(3) 2 0.76 1.22 C [0]l6 Steel RA(3)+FEP 2 081 0.04 C [0,45,-45,90]2S Steel RA(2)+FEP 2 0.95 0.15 -79-Chapter 3 Spring-in of Angled Composite Laminate Table 3.11: Spring-in variation with part shape for aluminum tool. COMPRO predictions are compared with experimental results. Part Shape Lay-up Tool Material Tool Surface Cure Cycle Spring-in Corner Warpage Web Flange COMPRO L [0]8 Aluminum TS-2 2 0.95 0.45 0.45 C [0]s Aluminum TS-2 2 1.09 2.57 0.48 C-modified [0]8 Aluminum TS-2 2 0.83 0.75 0.08 Experiment L [0]s Aluminum RA(3) 2 0.96 0.38 0.42 C [0]8 Aluminum RA(3) 2 0.97 0.83 0.09 Table 3.12: Spring-in variation with part shape for steel tool. COMPRO predictions are compared with experimental results. Part Shape Lay-up Tool Material Tool Surface Cure Cycle Spring-in Corner Warpage Web Flange COMPRO L [0]8 Steel TS-2 2 0.87 0.25 0.25 C [0]s Steel TS-2 2 0.95 1.37 0.28 C-modified [0]8 Steel TS-2 2 0.80 0.46 0.07 Experiment L [0]8 Steel RA(3) 2 0.74 0.29 0.34 C [0]8 Steel RA(3) 2 0.81 0.45 0.10 -80-Chapter 3 Spring-in of Angled Composite Laminate Figure 3.1: Tool geometry used in the experiment (Albert, 2001). web length Z 90° 0° (a) — a / y V i d t h ^flange lengthy 90' — > n° (b) flange angle Figure 3.2: Part geometries used in the experiment (a) C-shaped (b) L-shaped (Albert, 2001). -81-Chapter 3 Spring-in of Angled Composite Laminate £ E 50 mm Fibre Direction Thickness, t »| R Rj = inner radius = 5 mm R 0 = outer radius t = 1.6,3.2,4.8 mm Tool -Aluminum, Steel Part Tool Figure 3.3: The geometry of the L-shaped composite part and the tool used in modelling. c A i r Fixed Convective Heat Sliding Transfer Tool re moval BC Figure 3.4: Finite element representation of the angle laminate and the boundary conditions applied. -82-Chapter 3 Spring-in of Angled Composite Laminate 200 0 Temp - 1 hold Temp-2 holds Vacuum 100 200 Time (min) Pressure "1 hold Pressure 2 holds ~~ 800 600 400 1 200 im 3 w SB -~ -200 300 400 Figure 3.5: Target cure cycles used in autoclave processing of the L-shaped part. Both one and two hold cure cycles (similar to the ones used by Albert, 2001) are used. Original Shape Figure 3.6: Measurement of spring-in of the L-shaped part. -83-Chapter 3 Spring-in of Angled Composite Laminate 0.12 0 10 20 30 40 50 Distance, x (mm) Figure 3.7: Predicted flange warpage profile of a 16-layer L-shaped part. 0.7 0.6 ^ 0 . 5 •-OD 3 0.4 6*0.3 e % 0.2 0.1 • Plane Strain • 3-Dimensional [0]16 [0745/-45/90]2S [0/90]4s Figure 3.8: Comparison of spring-in prediction between plane strain and three-dimensional analysis. -84-Chapter 3 Spring-in of Angled Composite Laminate Figure 3.9: Numerical model to predict the equivalent properties of a composite laminate using a 3D elastic analysis, (a) Schematic diagram of a 16-layer composite element (b) & (c) identification of faces for boundary conditions. -85-Chapter 3 Spring-in of Angled Composite Laminate X=0 < > A z A Plies 23 -n 1 > X 1 P 'y i V / tpiyN/ 'Laminate Figure 3.11: Part geometry for the analytical model. Sliding friction occurs at both the tool-part interface and at the interply region (Twigg, 2001). -86-Chapter 3 Spring-in of Angled Composite Laminate A LDcbond Vi TSliding - h Sticking Condition Sliding Condition > Debond Front Migration Tool Centre Tool End X Coordinate Figure 3.12: Illustration of interfacial shear stress distribution as the debond front travels from the tool end towards the center (Twigg, 2001). 10 50 100 150 200 250 300 Distance from the Centre line (mm) Figure 3.13: The variation of interfacial shear stress distribution with part thickness (600 mm long part with unidirectional layers) according to COMPRO prediction using a shear layer of modulus 6.9 kPa. -87-Chapter 3 Spring-in of Angled Composite Laminate Figure 3.14: The variation of interfacial shear stress distribution with part length (8 unidirectional layer part) according to COMPRO prediction using a shear layer of modulus 6.9 kPa. 10 0 100 200 300 400 500 600 Distance from the Centre line (mm) Figure 3.15: The variation of interfacial shear stress distribution with part length. By modifying the shear layer property with length ratio, the peak shear stresses are made equal (average shear stresses are also equal). -88-Chapter 3 Spring-in of Angled Composite Laminate Figure 3.16: The variation of the warpage of a flat part with part thickness. The prediction of COMPRO is compared with the analytical solution proposed by Twigg (2001). Figure 3.17: The variation of the warpage of a flat part with part length. The prediction of COMPRO, with the shear layer modified by part length ratio, is compared with the analytical solution proposed by Twigg (2001). -89-Chapter 3 Spring-in of Angled Composite Laminate Figure 3.18: Contribution by the CTE anisotropy and cure shrinkage to the total spring-in when the part is constrained by the tool. • Warpage Component • Corner Component SL-0-6.9E-3 Pa SL-1 - 6.9E+5 SL-2 - 6.9E+6 SL-3 - 6.9E+7 SL-4 - 6.9E+8 SL-5 - 6.9E+10 (Tool) o — CO CO CO 05 T i co CO 05 05 - J CO CO CO >—1 co CO CO CO CO 8 Plies 16 Plies 24 Plies Figure 3.19: Spring-in of the parts with different shear layers. Only the Young's modulus is shown in here although the shear modulus is also changed by the same order. -90-Chapter 3 Spring-in of Angled Composite Laminate 1.2 -OX) <u •V OX) s a !/3 0.8 0.4 • -COMPRO • - Experiment TS-1 RA(3)+FLP TS-2 RA(3) Figure 3.20: Comparison between the predicted and measured warpage components of an 8-layer part to select the appropriate shear layer property to represent the tool surface condition used in the experiment. ^1 .5 <u S l a OX) *> -a w .5 1 OX) = Cm 0.5 -COMPRO Experiment + 2 8 Layers • Warpage Component • Comer Component COMPRO 83 ->• 16 Layers Experiment + 2 2S Figure 3.21: Spring-in variation with tool surface condition. COMPRO predictions are compared with experimental results for both 8-layer and 16-layer unidirectional parts. -91-Chapter 3 Spring-in of Angled Composite Laminate Figure 3.23: Variation of spring-in with part thickness predicted by COMPRO. The results are shown for two different tool surfaces. -92-Chapter 3 Spring-in of Angled Composite Laminate -ST 1.5 5 •~ w> u J 1 I wo e • M l U a «» 0.5 • Warpage Component • Corner Component COMPRO 16 TS-1 Expe rime nt 8 16 RA(3)+FEP COMPRO Expe rime nt 16 8 16 TS-2 RA(3) Figure 3.24: Spring-in variation with part thickness. COMPRO results are compared with experimental results for two different tool surface conditions. u lm OX) •a OX) s a 1 0.5 • Warpage Component • Corner Component COMPRO Al St Experiment Al St COMPRO Al St Experiment Al St TS-1 RA(3)+FEP TS-2 RA(3) Figure 3.25: Spring-in variation with tool material for an 8-layer unidirectional part. COMPRO results are compared with experimental results for two different tool surface conditions. -93-Chapter 3 Spring-in of Angled Composite Laminate 1.5 l 4 l a OX) OX) .5 a 0.5 • Warpage Component • Corner Component COMPRO 0 Al St TS-1 Expe rime nt Al St RA(3>fFEP COMPRO Al St TS-2 Expe rime nt Al St RA(3) Figure 3.26: Spring-in variation with tool material for a 16-layer unidirectional part. COMPRO results are compared with experimental results for two different tool surface conditions. No-Shear • Warpage Component • Corner Component No-Shear TS-1 —^16 Layers u Q u Q u TS-2 Figure 3.27: The predicted variation of spring-in with part lay-up for two different thickness of the part. The results for three different shear layer properties are also shown. No-Shear indicates that there is minimal tool part interaction. -94-Chapter 3 Spring-in of Angled Composite Laminate 2.5 u V s--S 1.5 _= '•Z in l 4 0.5 COMPRO U Experiment U 8 Layers -4-• Warpage Component • Corner Component U=[0]n Q = [0/45/-45/90]n COMPRO Experiment U U Q -> 16 Layers Figure 3.28: Spring-in variation with part lay - up. COMPRO results are compared with experimental results for two different thickness of the part. Web-50 mm S s Tool Line of Symmetry Figure 3.29: The geometry of half of the C-shaped composite part and the tool. -95-Chapter 3 Spring-in of Angled Composite Laminate 4.5 91 -oxi <u •a WD S 11.5 C73 • Warpage - Flange • Warpage - Web • Corner No-Shear TS-1 TS-2 8 - Layers <4~ No-Shear TS-1 — • 16 - Layers TS-2 Figure 3.30: Variation of spring-in with part shapes predicted by COMPRO. The results are shown for three different tool surface conditions. No-Shear indicates that there is minimal tool part interaction. -OX) •o OX) s-a. 0.5 • Warpage - Flange • Warpage - Web • Corner RA(3) RA(3)+FEP 8 Layers -4~ RA(3) 16 Layers RA(3)+FEP Figure 3.31: The experimentally measured spring-in of C-shaped parts. In all cases, the warpage of the flange is negligible (Albert, 2001). -96-Chapter 3 Spring-in of Angled Composite Laminate A Figure 3.33: Illustration of the three zones of the shear layer used in the modified method to overcome the corner-locking. -97-Chapter 3 Spring-in of Angled Composite Laminate 4.5 OD WD a !Z1 1.5 4 • Warpage - Flange • Warpage - Web • Corner C m TS-1 TS-2 8 - Layers -4~ TS-1 - • 16 - Layers TS-2 Figure 3.34: Spring-in results by the modified method to overcome the corner-locking problem. C is the original results and C,„ is the modified results. • Warpage - Flange • Warpage - Web • Corner COMPRO Experiment Figure 3.35: Effect of part shape on spring-in. COMPRO predictions by both the original method (C) and the modified method (Cm) are compared with experimental results for aluminum tool. -98-Chapter 3 Spring-in of Angled Composite Laminate • Warpage - Flange • Warpage - Web • Corner i U m i C L L C COMPRO Experiment Figure 3.36: Effect of part shape on spring-in. COMPRO predictions by both the original method (C) and the modified method (Cm) are compared with experimental results for steel tool. -99-C H A P T E R 4: PROCESS-INDUCED RESIDUAL STRESSES 4.1 INTRODUCTION In the previous section the intrinsic and the extrinsic parametres that contribute to the residual stresses in a composite part were investigated. Once the part is removed from the tool, these residual stresses are released and cause the L-shaped part to spring in and warp. But due to the changes in the material properties during the curing process, some of the residual stresses will be locked in after the part is released from the tool. Even though this portion of the residual stresses is small compared to the residual stresses that cause the part to spring in, it can potentially elevate the stresses in the composite structure during its service life and lead to its failure sooner than expected. In addition to the macro-level stresses that develop at the laminate-level, there will be micro-level stresses at the ply-level. Micro-level stresses are defined as those stresses that develop between the fibre and matrix as a result of thermal expansion coefficient mismatch and the resin shrinkage during curing. These stresses are complicated due to the complex arrangements of the fibres in the matrix in a ply. But these stresses are generally self-equilibrating and do not contribute to the macro-level deformations. Although micro-level stresses do not contribute to the macro-level deformation, they may lead to reduction of strength and initiation of cracking. Since the experimental methods are complicated and expensive, the simplest way to evaluate the process-induced micro-level residual stresses is to use numerical methods. -100-Chapter 4 Process-Induced Residual Stresses Wagner and Nairn (1997) proposed a simple concentric cylindrical model to predict the thermal stresses in a composite material. This model is similar to the model used by Rosen and Hashin (1964) to predict the mechanical properties of composite materials. Several researchers used the finite element approach to analyse the thermal stresses in a composite material (Bowles et al., 1991; Hyer, 1997). In all these analyses, only the stresses due to thermal mismatch are considered and the resin properties are assumed to be either constant or temperature dependent. Chen et al. (2001) used a viscoelastic micromechanical model to calculate the residual stresses induced during the curing process. In this work, only the stresses due to thermal mismatch were considered. In this section the residual stresses due to thermal mismatch, cure shrinkage and mechanical load at the micro-level in a specific lamina (Figure 4.1) of the previously studied L-shaped part are investigated. 4.2 MICROMECHANICAL FINITE ELEMENT MODELLING In a composite lamina, the actual fibre arrangement is random. But for simplicity, the fibres are assumed to be arranged in a periodic manner so that a representative volume element (RVE) can be identified and analysed separately. The commonly used periodic arrangements of fibres are square-packed and hexagonal-packed arrays. The names of the arrays are derived from the shape of the polygons that describe the repeating fibre-packing pattern as shown in Figure 4.2. The hexagonal-packed array is the preferred model since it closely represents the actual fibre arrangement as shown in Figure 4.3. If a load is applied in the fibre direction, each fibre will respond in the same way as its neighbouring fibres, since the fibres are embedded in a vast array. So the analysis can be focused only on the R V E that represents the whole fibre-matrix system. -101-Chapter 4 Process-Induced Residual Stresses The R V E has the same fibre volume fraction and material properties as the lamina. In this study, a unidirectional ply with continuous fibres is considered and the following assumptions are made: The fibres are circular in cross section and continuous in the longitudinal direction The fibres are arranged in a periodic square or hexagonal array The fibre properties are isotropic and constant (do not change with temperature) The matrix is cure hardening instantaneously linear elastic (CHILE) as it is in COMPRO The displacements are continuous across the fibre/matrix interface (perfect bonding) The temperature is uniform throughout the body The applied load is uniform over the area of application The R V E for square-packed and hexagonal-packed arrays are shown in Figure 4.4 and the corresponding finite element mesh is shown in Figure 4.5. The dimension of the unit cell in the axial direction (1) is not important, since the stresses and strains are independent of the axial coordinate. So the axial dimension is chosen to keep the element aspect ratio below 2 (here the axial dimension is taken to be 0.25 units). Due to the symmetry, only a quarter of the R V E is considered. The analysis is performed using the commercial finite element software A B A Q U S . 20-noded isoparametric solid elements are used for the analysis. Mesh sensitivity and convergence of the finite element solution is examined by dividing each element into four as shown in Figure 4.5(b). -102-Chapter 4 Process-Induced Residual Stresses 4.3 BOUNDARY CONDITIONS Because of the periodicity, the straight lines outlining the unit cells in Figure 4.6 remain straight when the composite is subjected to normal loading (tensile, compressive) and temperature change. Symmetry boundary condition implies that the corresponding displacement component at that boundary and in the direction normal to that boundary is zero and the couple boundary condition (all the nodes along this boundary are coupled with each other to deform by the same amount in a specified direction) implies that the corresponding displacement component at that boundary and in the direction normal to that boundary is constant. A summary of the boundary conditions is given below (x/, X 2 and xj are any arbitrary distances from the origin in 1, 2 and 3 directions respectively): w(0,X2,X3) = 0 - Symmetry - 1 «(0.025, X 2 y X j ) = constant - Couple - 1 v(x/,0, X 3 ) = 0 - Symmetry - 2 v(x/,l, xi) — constant - Couple - 2 w(x], X2,0) = 0 - Symmetry - 3 w(xj, X 2 , l ) = constant - Couple - 3 The couple boundary is applied using the EQUATION option in A B A Q U S . The boundary conditions are shown schematically in Figure 4.6. 4.4 MATERIAL MODELLING IN A B A Q U S Three types of material properties are necessary for the calculation of stresses: elastic constants, thermal expansion coefficients and cure shrinkage coefficients. A l l materials are assumed to be instantaneously linear elastic i.e. although the material properties vary throughout the processing, -103-Chapter 4 Process-Induced Residual Stresses at any given instant it behaves as a linear elastic material. The material properties of the fibre and the matrix are given in Appendix A (Full details of the material can be found in Johnston Elastic constants of the transversely isotropic fibre are assumed to vary linearly with temperature. For AS4 fibre material, the properties remain constant with the temperature. The isotropic matrix resin is modelled as a cure hardening instantaneously linear elastic material. This means the resin modulus is assumed to be linear elastic at any given time and the modulus increases with the degree of cure. As implemented in COMPRO, in this case too, the heat transfer and the cure kinetics equations are solved independently (Johnston, 1997). Since consideration is given to a single ply in the L-shaped part, the temperature-time history of this ply is imported from COMPRO and a uniform temperature is assumed throughout the fibre and the matrix. The temperature history is shown in Figure 4.7. In COMPRO, the cure kinetics equation (Equation 4.1) is solved by Backward Euler method. But in this case, the equation is solved by Forward Euler method to avoid iteration (Equation 4.2). (1997)). da dt Kam(\-a)" (4.1a) I + e C{a-(aC0+<XcTT)} K = Ae -AEIRT (4.1b) where a Resin degree of cure T Resin temperature Total resin heat of reaction (a= 0 to 1) A Pre-exponential factor AE Activation energy m, n Exponent parameters -104-Chapter 4 Process-Induced Residual Stresses R Gas constant C, aco, OCCT Diffusion constants (aM-a,)_ * ,< ( l -g , . r At 1 + exp{C(«, - aco - a^T,)} (4.2a) K,=Ae-&E'ST-aM = At KrfQ-a.y l + exp{C(a, - a c o -aCTT,)} + a, (4.2b) The accuracy of the results depends on the time step taken: the smaller the time step the more accurate the results. The time step used in COMPRO vary from 1 second to 100 seconds (user can define the minimum and the maximum time step). COMPRO increases and reduces the time step according to the change in a. The degree of cure calculation using the forward Euler method with the time steps of 10 min, 5 min and 2 min are compared with the results of COMPRO in Figure 4.8. A zlr of 2 min is used in the subsequent calculations. The resin modulus development during curing is calculated according to Equation 4.3 (Johnston, 1997). Er=K[\ + aER{T-T0)] K = EI+[(r- - rcx )I(T'C2 - T„ )] (E; - K) rcx < r < rcl (4.3) K = E; r > rC2 T*=(T°+arga)-T T* -T* +T* T 1Cl lC\a ^ 1C\b l where a The resin degree of cure ER0 The resin elastic modulus at very low degrees of cure -105-Chapter 4 Process-Induced Residual Stresses E™ The resin elastic modulus at To and a - 1.0 T Resin temperature T* Difference between resin temperature and resin instantaneous Tg Tga Glass transition temperature at a = 0 Tgb Factor expressing the degree of cure dependence of resin glass transition temperature Tcia T above which resin modulus begins to increase at T = 0 K Tcib Factor expressing the temperature dependence of the T* above which resin modulus begins to increase Tc2 T* above which resin modulus has reached its full value To Resin modulus development reference temperature ctEr Factor expressing the temperature dependence of resin modulus The resin modulus development model is implemented using the UFIELD subroutine in A B A Q U S . The A B A Q U S calculation is verified by a simple spreadsheet calculation throughout the cure cycle. 4.5 STRESS ANALYSIS Since the material is assumed to be linear elastic at a given time, the increment in stress for that time step is given by Hooke's law: ACT = E AS (4.4) Where E is the instantaneous modulus and the total stress is given by a M = cr, + ACT (4.5) -106-Chapter 4 Process-Induced Residual Stresses However, in A B A Q U S , i f we use the UFIELD option to calculate the instantaneous material properties, the stress values are also instantaneous (Equation 4.4). Thus the accumulating part (Equation 4.5) should be done separately. Although the cumulative stresses can be calculated using the U M A T subroutine in A B A Q U S , for simplicity a spreadsheet is used for this calculation. The stresses due to CTE mismatch and cure shrinkage are calculated separately and added together at each time step to find the total stress at any instant during the curing process. 4.5.1 THERMAL STRESSES Equation 4.6 gives the instantaneous thermal stress at any given time. cjM=ECTE(Ti+]-T0) (4.6) Where E is the instantaneous Young's modulus, To is initial temperature and Ti+i is the current temperature.' The above stress is due to the temperature change of (Ti+i - T0). But for the current time step, the required stress should be due to a temperature change of (Ti+i - Ti). So the incremental stress due to thermal mismatch can be calculated by scaling the instantaneous stress by the temperature increment for the current time step, i.e. PM-T,) '+' (TM-T0) ACT, = <rM j : (4.7) where To is the initial temperature, 07+/ is the instantaneous stress due to thermal mismatch at temperature Ti+j and (Ti+i - T,) is the temperature increment at this time step. The cumulative stress at this time is given by Equation 4.5. -107-Chapter 4 Process-Induced Residual Stresses 4.5.2 CURE SHRINKAGE STRESSES The stress due to resin cure shrinkage is calculated by converting the linear resin cure shrinkage strain to an equivalent thermal expansion (shrinkage) coefficient and setting the thermal expansion coefficient of fibre equal to zero. The volumetric cure shrinkage strain is given by Equation 4.5, Vf =0.0 a< a a Vf =Aas+ (Vf -A)a] aci <a<aC2 (4.8) ys _ys« a >ac '•a as =• « - « c i aC2 — acl where Vr Resin volumetric cure shrinkage VrSa° Total volumetric resin shrinkage from a = 0 to 1 as The 'cure shrinkage' degree of cure aci Degree of cure after which the resin shrinkage begins ac2 Degree of cure after which the resin shrinkage stops A Linear cure shrinkage coefficient This volumetric strain is converted to a linear strain by es=[\ + Vf)n-\ (4.9) where ss is the linear shrinkage strain and Vf is the volumetric shrinkage. This is the strain that a lamina undergoes at a given time. The incremental linear cure shrinkage strain for the time step will be (4.10) -108-Chapter 4 Process-Induced Residual Stresses where ef+x and sf+l are the linear cure shrinkage strains at the current time and the previous time respectively. From this value, the incremental cure shrinkage stress and the total cure shrinkage stress are calculated according to Equations 4.4 and 4.5. Applying the incremental strain directly to the resin will be difficult. Instead, this incremental strain is applied as thermal load by converting the incremental strain to an equivalent thermal expansion coefficient given by CTE, = ' + 1 (4.11) where TM is the instantaneous temperature and T0 is the initial temperature. Hence, the incremental cure shrinkage stress will be given by Aas=ECTEs(TM-T0) (4.12) and the total stress can be calculated according to Equation 4.5. The axial stress development in the fibre during the curing process due to thermal mismatch and cure shrinkage is given in Figure 4.9. According to this figure, the contribution of thermal mismatch to the total stress is around 85% and the contribution from cure shrinkage around 15%. Although the cure shrinkage strain is very high compared to the thermal strain, most of the cure shrinkage occurs at low resin modulus. Hence, the contribution from cure shrinkage is low compared to the thermal mismatch. The stress components to be discussed are shown in Figure 4.10. Since the stress values at the fibre-matrix interface are the critical ones, here we consider only those values. The -109-Chapter 4 Process-Induced Residual Stresses circumferential locations around the fibre-matrix will be identified by the angle 0 which will be considered positive clockwise from the axis 3. The normal and shear stresses acting on the interface between the fibre and the matrix, a r and xre, have the same value on the matrix side of the interface as on the fibre side due to the equilibrium condition. These stress components are responsible for the failure of the interface. The normal and shear stresses at the interface between the fibre and the matrix for the AS4/8552 material due to thermal mismatch and cure shrinkage are shown in Figure 4.11. The stresses due to cure shrinkage are lower compared to the stresses due to CTE mismatch. The normal stress, ar, depends strongly on 9, and it is compressive at all circumferential locations. This strong dependence on 0 is due to the interaction between fibres. The compressive stresses are higher at locations where fibres are close to each other i.e. when 6 = 30 and 9 = 90. Compared to normal stresses, shear stresses are smaller in magnitude. Since the normal stresses are compressive, there is no possibility for crack initiation in the matrix but the shear stresses may cause interface failure. The circumferential stress components in the fibre and the matrix do not have the same value (Figure 4.12). The tensile circumferential stress in the matrix may contribute to crack formation in the matrix. The crack in this case will be emanating radially from the interface into the matrix. Note that the ultimate tensile strength of the 8552 resin is around 120 MPa (obtained from the H E X C E L Cooperation product data sheet for 8552 Epoxy Matrix). The longitudinal stress components in the fibre and the matrix are shown in Figure 4.13. The longitudinal stresses in the fibre are compressive while the longitudinal stresses in the matrix are tensile as expected (the CTE of the resin is higher than the CTE of the fibre). The tensile stresses -110-Chapter 4 Process-Induced Residual Stresses in the matrix may contribute to matrix cracking. The longitudinal stresses are quite uniform in the fibre and the matrix. 4.6 CONCENTRIC CYLINDRICAL MODEL The prediction of the micro-level thermal residual stresses using the finite element method is time consuming. Can the residual stresses be predicted by simple calculations? In this section, the thermal residual stresses calculated by both the finite element method (hexagonal array) and the concentric cylinders model proposed by Wagner et al. (1997) are compared. The stresses in the fibre and the matrix by the concentric cylinder model are given in Equation 4.13. ajr = a{ = Af Stresses in the fibre: (4.13a) cr{ = CJ where AJ and CJ are constants. <y"; =Am + Bm(Rl/r)2 Stresses in the matrix: < = Am - Bm (Rx I rf (4.13b) Cr{ = Cm where Am ,Bm and C" are constants and i?, is the radius of the fibre. ar, ae, and cr, are radial, hoop and axial stresses as shown in Figure 4.10. Due to the axisymmetry assumption, all the stresses are independent of 6. Full detail of the Equation 4.13 is given in Appendix C. Since most of the thermal residual stresses develop at the cool down part of the processing cycle, in this case only this part of the cure cycle is considered for comparison of the methods. The material properties of the AS4/8552 material are assumed to remain constant during cool down according to the data available (Johnston, 1996). The thermal residual stresses that develop when the composite is cooled down from 180 °C to room temperature (20 °C), predicted by both finite -111-Chapter 4 Process-Induced Residual Stresses element analysis and elasticity approach, are shown in Figures 4.14 to 4.16. In these figures, the various stress components at the fibre-matrix interface are plotted as a function of circumferential location around the interface (Figure 4.10). The above results shows that the axisymmetric assumption is not good enough in predicting the stresses, particularly the normal and shear stresses. This is due to high fibre volume fraction of the material. It is shown by Hyer (1997) that for low fibre volume fractions (up to 40%), the axisymmetric assumption is quite accurate and the accuracy decreases with increasing fibre volume fraction. Even though the axisymmetric condition is not valid for materials with high fibre volume fraction, the elasticity approach provides a quick and reasonably good estimation of the residual stresses compared to the more time consuming finite element analysis. 4.7 TRANSFERRING MACRO-LEVEL STRESSES TO MICRO-LEVEL High residual stresses coupled with the low ductility of the matrix may lead to the fracture of the fibre-matrix interface. In order to quantify the residual stresses at the fibre-matrix interface, the macro-level stresses should be transferred to the micro-level. These macro-level stresses are applied as mechanical loads on the R V E . The application of a unit mechanical stress (any component) on the R V E produces stresses (all six component) at the interface. The ratio of the applied stress and the resultant interfacial stress is called the multiplication factor or magnification ratio, i.e. corresponding to an applied stress component, there are six magnification factors for the six stress components at the interface. If the actual mechanical stress is known, then the interface stresses can be calculated by multiplying the mechanical stress with the multiplication factor corresponding to each stress component at the interface. Since the properties of the fibre and the matrix vary during the curing process, these magnification factors -112-Chapter 4 Process-Induced Residual Stresses are also not constant during the curing process. For AS4/8552 material, only the resin modulus varies during the curing process. So the multiplication factor also changes according to the changes in resin modulus as shown in Figure 4.17. Hence, the transfer of the mechanical stresses should be done at suitable intervals along the cure cycle. These intervals should be selected according to the composite materials used in the processing. For example, for the material AS4/8552 used here, the magnification factor changes only during a short period of time as shown in Figure 4.17. The time step that was used to calculate the degree of cure (2 min) is used during this period for transferring the macro stresses to micro-level. 4.8 PREDICTION OF ENGINEERING PROPERTIES USING MICROMECHANICS The interactions of the elastic properties of the fibre and the matrix produce the elastic properties of the composite material. Several micromechanical equations have been proposed for the prediction of composite properties from those of the constituent fibre and matrix using variational principles. A survey of the existing micromechanical models was carried out by Hashin(1983). In COMPRO, the engineering properties are calculated by using the micromechanics equations found in the paper by Bogetti and Gillespie (1992) that was originally proposed by Rosen and Hashin (1964). These equations are based on classical elasticity solution by assuming that the volume of fibres and matrix in a composite can be filled with an assemblage of cylindrical fibres and surrounding matrix material, with the fibres being of various sizes to the degree that the fibre-matrix combinations of cylinders completely filled the volume of the composite. This -113-Chapter 4 Process-Induced Residual Stresses notion is called the composite cylinders model. The ratio of the diameter of the fibre to the diameter of the surrounding matrix represents the fibre volume ratio. In this section, the engineering properties calculated by both the elasticity approach and numerical approaches are compared. Since the resin elastic properties change during the autoclave processing, the results are compared for different values of resin modulus. First we consider an isotropic case. In this case, both the fibre and the matrix are isotropic. The material properties of fibre and matrix are selected arbitrarily and given in Table 4.1. The effective material properties from different approaches are given in Table 4.2. A l l these methods quite accurately predicted the longitudinal engineering properties (E\, vn). Even the simple rule of mixture predicts the same results. But there are discrepancies in the prediction of transverse engineering properties as shown in Figure 4.18. In this figure, the transverse modulus predicted by all three approaches is compared for different matrix modulus. At high matrix modulus, both approaches (Hexagonal and Square) agree well with the elasticity approach. But at lower matrix modulus, the difference between the square and elasticity approaches is much greater (around 30%) compared to the difference between hexagonal and elasticity approaches (around 5%). This can be justified by comparing the unit cells of all these approaches. The unit cells for all three approaches (elasticity, square, hexagonal) are shown in Figure 4.19. The cross-sectional areas of all three unit-cells are equal and the fibre volume fraction is 0.573. A l l these three unit-cells are superposed on top of each other and are shown in Figure 4.20. By observing this figure, we can see that the unit cell considered in the elasticity approach is somewhat closer to the hexagonal than the square array unit cells. -114-Chapter 4 Process-Induced Residual Stresses The predicted composite engineering properties for the AS4/8552 material by hexagonal array and elasticity approach are listed in Table 4.3. 4.9 SUMMARY AND CONCLUSIONS From the above study, it is clear that apart from the residual stresses due to the thermal mismatch between the fibre and the resin, the cure shrinkage also contributes to residual stress formation. Although the cure shrinkage strains are higher than the thermal strains, for this resin data set, most of the cure shrinkage occurs when the resin modulus is low thus resulting in low residual stresses. For another material the contribution of the cure shrinkage may be greater. Significant amount of thermal residual stresses develop only during the cooling down part of the curing process. Therefore, the thermal residual stresses can be calculated with reasonable accuracy by only considering the cooling down part of the processing cycle. The micro-level residual stresses can be predicted using the simple concentric cylinders model even though the accuracy of the model reduces with increasing fibre volume fraction. -115-Chapter 4 Process-Induced Residual Stresses 4.10 TABLES Table 4.1: Material properties of isotropic fibre and matrix Property Isotropic Fibre Isotropic Matrix Young Modulus (GPa) 100 1-100 Poisson's Ratio 0.25 0.25 Table 4.2: Material properties of composite for fibre volume fraction of 0.573. Both fibre and matrix are isotropic Resin Modulus (GPa) Elasticity Approach Numerical Approach Square Array Hexagonal Array E i (GPa) E 2 (GPa) Vl2 V23 E , (GPa) E 2 (GPa) Vl2 V 23 Ei (GPa) E 2 (GPa) V l 2 V 23 1 57.73 3.13 0.250 0.329 57.73 4.01 0.250 0.173 57.73 3.26 0.250 0.301 10 61.57 25.44 0.250 0.299 61.57 29.30 0.250 0.208 61.57 26.03 0.250 0.283 20 65.84 42.42 0.250 0.280 65.84 46.03 0.250 0.226 65.84 42.93 0.250 0.271 30 70.11 54.88 0.250 0.268 70.11 57.62 0.250 0.235 70.11 55.23 0.250 0.264 40 74.38 64.65 0.250 0.261 74.38 66.54 0.250 0.241 74.38 64.86 0.250 0.259 50 78.65 72.67 0.250 0.256 78.65 73.89 0.250 0.245 78.65 72.78 0.250 0.255 60 82.92 79.48 0.250 0.253 82.92 80.21 0.250 0.247 82.92 79.54 0.250 0.253 70 87.19 85.44 0.250 0.252 87.19 85.81 0.250 0.248 87.19 85.46 0.250 0.252 80 91.46 90.75 0.250 0.251 91.46 90.90 0.250 0.249 91.46 90.75 0.250 0.251 90 95.73 95.57 0.250 0.250 95.73 95.60 0.250 0.250 95.73 95.57 0.250 0.250 100 100.00 100.00 0.250 0.250 100.00 100.00 0.250 0.250 100.00 100.00 0.250 0.250 Table 4.3: Material properties for an AS4/8552 lamina. Material properties are calculated at both initial (uncured) and final (cured) resin modulus. Resin Modulus (GPa) Elasticity Approach Hexagonal Array E, (GPa) E 2 (GPa) Vl2 V23 CTE1 Hm/m/C CTE3 |^m/m/C E, (GPa) E 2 (GPa) Vl2 V23 CTE1 (xm/m/C CTE3 (im/m/C 4.67E-3 120.3 0.017 0.263 0.525 0.031 28.7 120.3 0.018 0.263 0.479 0.031 25.9 4.67 122.4 9.900 0.268 0.471 0.714 28.5 122.4 9.862 0.269 0.431 0.811 27.3 -116-Chapter 4 Process-Induced Residual Stresses ooeeo ooooo (a) (b) Figure 4.2: The cross section of a unidirectional ply with commonly used periodic fibre arrangements: (a) square array (b) hexagonal array. -117-Chapter 4 Process-Induced Residual Stresses Figure 4.3: Cross section of graphite-reinforced composite material (Hyer, 1997). _ 3 t oooooi ooooo ooopo (b) Figure 4.4: The representative volume element (RVE) concept for (a) square array (b) hexagonal array. -118-Chapter 4 Process-Induced Residual Stresses (a) (b) (c) Figure 4.5: Finite element mesh used for the micro stress analysis: (a) square array (b) square array (fine mesh) (c) hexagonal array. Figure 4.6: The applied boundary condition for the RVE. Due to the symmetry only a quarter of the RVE is considered. "Sym" stands for symmetric boundary and "Coup" stands for the couple boundary. -119-Chapter 4 Process-Induced Residual Stresses 200 0 -I 1 1 1 — — 1 1 1 1 0 50 100 150 200 250 300 350 Time (min) Figure 4.7: The temperature-time history of the element shown in Figure 4.1 as obtained from COMPRO. Figure 4.8: The development of degree of cure of resin as it cures during the autoclave processing of the L-shaped composite part: the ABAQUS results for three different time steps (10 min, 5 min and 2 min) are compared with COMPRO results. The circled portion of the graph is magnified and shown. -120-Chapter 4 Process-Induced Residual Stresses 200 150 100 50 Temperature — Thermal Mismatch Cure Shrinkage 50 100 150 200 Time (min) 250 300 + -7.5 + -15 Xfi V) -22.5 £ t/3 -30 350 Figure 4.9: The development of stresses due to CTE mismatch between the fibre and the matrix and cure shrinkage of resin as the resin cures during processing of an L-shaped composite part. 34 Figure 4.10: The axial (1), normal (r), and the circumferential (0) stress components at the fibre matrix interface. The maximum stress occurs at the interface. -121-Chapter 4 Process-Induced Residual Stresses 5 Figure 4.11: The distribution of normal (or) and shear (T^) stress components along the fibre-matrix interface due to the CTE mismatch and resin cure shrinkage. • CTE Mismatch - Fibre CTE Mismatch - Matrix Cure Shrinkage - Fibre Figure 4.12: The distribution of circumferential (oe) stress components along the fibre-matrix interface due to the CTE mismatch and resin cure shrinkage. -122-Chapter 4 Process-Induced Residual Stresses 50 25 a 4> b -25 •50 0 CTE Mismatch - Fibre • CTE Mismatch - Matrix - A — Cure Shrinkage - Fibre -*— Cure Shrinkage - Matrix 30 6(deg) 60 90 Figure 4.13: The distribution of the axial stress components, oi, along the fibre-matrix interface due to the CTE mismatch and resin cure shrinkage. 45 -^30 93 <B 1 C 55 o •15 Finite Element - Fibre • Elasticity - Fibre - A — Finite Element - Matrix Elasticity - Matrix 30 6 (deg) 60 90 Figure 4.14: Comparison of circumferential stresses between concentric cylinders model (elasticity approach) and finite element solution for a temperature change of -160 °C. -123-Chapter 4 Process-Induced Residual Stresses 0 Finite Element-ar Elasticity-a r Finite element-t^ Elasticity - T ^ 30 6(deg) 60 90 Figure 4.15: Comparison of normal and shear stresses between concentric cylinders model (elasticity approach) and finite element solution for a temperature change of -160 °C. 50 _ 25 08 04 b ^-25 -50 0 Finite Element - Fibre Elasticity - Fibre Finite Element - Matrix Elasticity - Matrix 30 0 (deg) 60 90 Figure 4.16: Comparison of axial stresses between concentric cylinders model (elasticity approach) and finite element solution for a temperature change of-160 °C. -124-Chapter 4 Process-Induced Residual Stresses 85 o s "3 •o o 3.75 2.5 4 1.25 4 o 4 o -*— Resin Modulus • Multiplication Factor 50 100 150 200 Time (min) 250 300 0.4 0.3 | U. s 0-2 1 o.i I 350 Figure 4.17: The variation of the stress multiplication factor for axial stress at the interface during the cure cycle. The multiplication factor is for the fibre at 9 = 0 (Figure 4.10) for an applied axial mechanical stress. 120 8> ioo fad c/l "s •o o 0> c 80 60 40 20 •— Elasticity Approach • Square Array * Hexagonal Array 20 40 60 80 Matrix Modulus, Er (GPa) 100 120 Figure 4.18: Comparison of transverse modulus, E 2 by three different approaches for different matrix modulus. Both fibre and matrix are isotropic. -125-Chapter 4 Process-Induced Residual Stresses (a) (c) Figure 4.19: Unit cells considered in: (a) elasticity approach (b) hexagonal array and (c) square array. Fibre Elasticity Figure 4.20: All unit cells superposed on top of each other to reveal the difference between different approaches. -126-CHAPTER 5: CONCLUSIONS AND FUTURE W O R K The objective of this thesis was to investigate the processing effect on the performance of composite structures. Since the experimental investigation is costly and time consuming, in this work, investigation was done using a combination of available software. The major accomplishments of the current work include: • A probabilistic design tool was proposed to determine the variability in the processing parameters to achieve the target performance reliabilities. A computational framework using commercially available and in-house developed numerical software was proposed for analysing the complete sequence of events from the manufacturing stage to the in-service structural performance of a composite structure. The in-house developed processing software, COMPRO, was used to predict the process-induced deformation and stresses on 2D cross section of the composite structure. An approximate superposition method was used to develop the 3D effects based on the 2D plane strain results obtained from COMPRO. The commercial finite element software, A B A Q U S , was used to analyse deterministically the process effect on the structural behaviour. A parametric study was carried out to determine the most influential process parameters and the inverse reliability analysis software, IRELAN, was used to determine the process parameters that must be controlled to meet structural performance standards with target reliabilities. • The parameters and mechanisms that affect the spring-in of angled thermoset composite parts during processing were investigated. -127-Chapter 5.Conclusions and Future Work In the design example presented earlier, only a limited number of process parameters that affect the structural performance of a simple composite structure, flat plate, was considered for simplicity. In this case, the intrinsic and extrinsic parameters that affect the dimensional stability of a more complex composite structure, angle laminate, were investigated using COMPRO. The prediction capability of COMPRO was examined by comparing the predicted results with available experimental results. Although, using a shear layer to model the tool-part interaction is not a physically reasonable method, it was shown that the shear layer could be used to predict most of the trends observed in the experimental study. A crude method of introducing a gap in the corner was proposed to overcome the corner-locking problem in the C-shaped parts. • The residual stresses developed during processing at micro-level were investigated. The process-induced residual stresses may lead to the premature failure of a composite structure. To apply any failure criteria available in the literature, the stresses at the material point should be known. Apart from the stresses that develop at the laminate level, the mismatch between the fibre and the matrix also produces additional stresses. Stresses due to the CTE mismatch and cure shrinkage were calculated. It is shown that for this particular resin material, the stresses due to cure shrinkage were not significant compared to the stresses due to CTE mismatch. A simple concentric cylinder model available in the literature was suggested to determine the stresses due to CTE mismatch with reasonable accuracy. A magnification factor method was proposed to transfer the laminate level stresses to micro-level during the curing process. 5.1 FUTURE WORK The proposed design tool is only an initial attempt to incorporate the processing effect explicitly into the design as a more scientific method to replace knockdown factors that attempt to do this -128-Chapter 5Conclusions and Future Work in an empirical fashion. Knockdown factors are typically established through an extensive experimental study that often may not be totally relevant to the current design problem. The following are some points to be considered in the future work: • The difficulty found during the course of this study was to find a realistic design problem where the processing effects were actually observed. To develop the confidence level in the proposed design approach, a realistic design problem should be analysed and compared by both the proposed probabilistic design approach and the traditional knockdown factor approach. • The spatial variation of the variabilities in the material and dimensional parameters should be incorporated through a more sophisticated stochastic finite element approach. • The variability of all possible parameters from the raw material level to structural level should be taken into account. For this, a robust inverse reliability system should be established to handle the huge amount of random variables. -129-REFERENCES A B A Q U S , (5.8), Hibbitt, Karlsson & Sorensen, Inc., 1998. Albert, C.I., "Spring-In of Angled Thermoset Composite Laminates", M.A.Sc. Thesis, The University of British Columbia, 2001. Arbocz, J. and Hoi, M . A . M . , "Collapse of Axially Compressed Cylindrical Shells with Random Imperfections", A I A A Journal, Vol . 29, No. 12, 1991, pp. 2247-2256. Bogetti, T.A. and Gillespie, J.W., Jr., "Process-Induced Stress and Deformation in Thick-Section Thermoset Composite Laminates", Journal of Composite Materials, Vol . 26, No. 5, 1992, pp. 626-659. Bowles, D.E. and Griffin, O.H., Jr., "Micromechanics Analysis of Space Simulated Thermal Stresses in Composites. Part I: Theory and Unidirectional Laminates", Journal of Reinforced Plastics and Composites, Vol . 10, 1991, pp. 504-521. Cederbaum, G. and Arbocz, J., "On the Reliability of Imperfection-Sensitive Long Isotropic Cylindrical Shells", Structural Safety, Vol . 18, No. 1, 1996, pp. 1-9. Chen, P.C. and Ramkumar, R.L., "RAMPC - An Integrated Three-Dimensional Design Tool for Processing Composites", 33 r d International SAMPE Symposium, March 7-10, 1988, pp. 1697-1708. Chen, Y. , Xia, Z. and Ellyin, F., "Evaluation of Residual Stresses Induced during Curing Processing Using a Viscoelastic Micromechanical Model", Journal of Composite Materials, Vol . 35, No. 6, 2001, pp. 522-542. Chryssanthopoulos, M.K. , Giavotto, V . and Poggi, C , "Characterization of Manufacturing Effects for Buckling-Sensitive Composite Cylinders", Composites Manufacturing, Vol . 6, No. 2, 1995, pp. 93-101. Elseifi, M.A. and Giirdal, Z., Nikolaidis, E., "Convex/Probabilistic Models of Uncertainties in Geometric Imperfections of Stiffened Composite Panels", A I A A Journal, Vol . 37, No. 4, April, 1999, pp. 468-474. Fernlund, G., Rahman, N . , Courdji, R., et al., "Experimental and Numerical Study of the Effect of Cure Cycle, Tool Surface, Geometry and Lay-up on the Dimensional Fidelity of Autoclave-Processed Composite Parts", Submitted to Composites Part A : Manufacturing. Foschi, R.O., Folz, B., Yao, F., L i , H . and Baldwin, J., " R E L A N : Reliability Analysis Software", Department of Civil Engineering, The University of British Columbia, 1999a. -130-References Foschi, R.O. and L i , H. , " IRELAN: Inverse Reliability Analysis Software", Department of Civil Engineering, The University of British Columbia, 1999b. Foschi, R.O. and L i , H. , "Reliability and Performance-Based Design in Earthquake Engineering", Proceedings, Int. Conf. for Structural Safety and Reliability, ICOSSAR, Newport Beach, CA, 2000. Hashin, Z., "Analysis of Properties of Fibre Composites with Anisotropic Constituents", A S M E Journal of Applied Mechanics, Vol . 46, 1979, pp. 543-550. Hashin, Z., "Analysis of Composites Materials - A Survey", A S M E Journal of Applied Mechanics, Vol . 50, 1983, pp. 481-505. Hayer, M.W., "Stress Analysis of Fiber - Reinforced Composite Materials", WCB/McGraw -Hil l , 1997. Hexcel Corporation Product Data Sheet, "8552 Epoxy Matrix", Hexcel Corporation, Retrieved 16 February 2002 <http://www.hexcelcomposites.com/downloads/techdata/PRE8552D.PDF>. Hubert, P. "Aspects of Flow and Compaction of Laminated Composite Shapes During Cure", Ph.D. Thesis, University of British Columbia, 1996. Hubert, P., and Poursartip, A. , "Study of the Interaction Between Percolation and Shear Flow in Processing of Thermoset Matrix Composites", American Society for Composites, College Station, Texas, September 2000. Johnston, A. , "An Integrated Model of the Development of Process Induced Deformation in Autoclave Processing of Composite Structures", Ph.D. Thesis, The University of British Columbia, 1997. Johnston, A. , Vaziri, R., and Poursartip, A. , "A Plane Strain Model for Process-Induced Deformation of Laminated Composite Structures", Journal of Composite Materials, Vol . 35, No.16, 2001, pp. 1435-1469. Kaw, A . K., "Promal", (1.0), CRC Inc., Boca Raton, FL, 1997. L i , H. , "An Inverse Reliability Method and its Applications in Engineering Design", Ph.D. Thesis, The University of British Columbia, 1999. L i , Y.W., Elishakoff, I. and Starnes, J.H., "Axial Buckling of Composite Cylindrical Shells with Periodic Thickness Variation", Computers and Structures, Vol.56, N o . l , 1995, pp. 65-74. Long, M.W. and Narciso, J.D., "Probabilistic Design Methodology for Composite Aircraft Structures", Final Report, U.S. Department of Transportation, Federal Aviation Administration, June 1999. Loos A C , Springer GS., "Curing of Epoxy Matrix Composites", Journal of Composites Materials, Vol.17, No. 2, 1983, pp.58-75. -131-References Nelson R H and Cairns DS., "Prediction of dimensional changes in composite laminates during cure", International SAMPE Symposium and Exhibition Covina, Ca, Book 2, 1989, pp.2397-2410. Osooly-Talesh, A. , "Implementation of an 8-noded 2D element in COMPRO", Internal Report, Composite Group, The university of British Columbia, November 2001. Radford DW and Diefendorf RJ., "Shape instabilities in composites resulting from laminate anisotropy", Journal of Reinforced Plastics and Composites, Vol.12, 1993, pp.58-75. Radford, D.W. and Rennick, T., "Separating Sources of Manufacturing Distortion in Laminated Composites", Journal of Reinforced Plastics and Composites, Vol . 19, No. 8, 2000, pp. 621-641. Rosen, B.W. and Hashin, Z., "Elastic Moduli of Fibre-Reinforced Materials", A S M E Journal of Applied Mechanics, Vol . 31, 1964, pp. 223-232. Rosen, B.W. and Hashin, Z., "Effective Thermal Expansion Coefficients and Specific Heats of Composite Materials", International Journal of Engineering Science, Vol . 8, 1970, pp. 157-173. Schapery, R.A., "Thermal Expansion Coefficients of Composite Materials Based on Energy Principles", Journal of Composite Materials, Vol . 2, 1968, pp. 380-404. Shiao M . C. and Chamis C.C., "Probabilistic Evaluation of Fuselage-Type Composite Structures", Probabilistic Engineering Mechanics, Vol . 14, 1999, pp. 179-187. Starnes, J.H., Jr., Hilburger, M.W., and Nemeth, M.P., "The Effects of Initial Imperfections on the Buckling of Composite Shells", Paper presented at the A S T M Symposium on Composite Structures: Theory and Practice, Seattle, WA, May 1 7 - 18, 1999. Twigg, GA., "Tool-part interaction in composites processing", M.A.Sc. Thesis, The University of British Columbia, 2001. Wagner, H.D. and Nairn, J.A., "Residual Thermal stresses in Three Concentric Transversely Isotropic Cylinders: Application to Thermoplastic - Matrix Composites Containing a Transcrystalline Interphase", Composites Science and Technology, Vol . 57, 1997, pp. 1289-1301. Zhu, T.L., "A Reliability-Based Safety Factor for Aircraft Composite Structures", Computers & Structures, Vol . 48, No. 4, 1993, pp. 745-748. -132-A P P E N D I X A : MATERIAL PROPERTIES OF HERCULES AS4/8552 UNIDIRECTIONAL PREPREG These data are for a fibre volume fraction of 0.573 and are obtained from Johnston (1997). A . l FIBRE ELASTIC PROPERTIES The fibre is a transversely isotropic material and the isotropic axes are the 2-3 axes. These properties are independent of temperature. Description Units Value Longitudinal Modulus, E l GPa 210 Transverse Modulus, E2 GPa 17.24 Major Poisson Ratio, v l2 0.2 Minor Poisson Ratio, v23 0.25 Longitudinal Shear Modulus, G13 GPa 27.6 A.2 RESIN CURE KINETICS MODEL da _ Kam{\-a)n dt ' 1 + e C M « C 0 + « C T / ' ) } Variable Description Units Values a Resin degree of cure. -T Resin temperature K HR Total resin heat of reaction (a= 0 to 1) J/kg 5.40E5 Ai Pre-exponential factor. /s 1.528E5 AEi Activation energy. J/mol 6.65E4 m Equation superscript. - 0.8129 n Equation superscript. - 2.736 -133-(A.1) (A.2) Appendix A Material Properties of Hercules AS4/8552 Unidirectional Prepreg R Gas constant J/(mol K) 8.3143 C Diffusion constant - 43.09 aco Diffusion constant - -1.684 OCT Diffusion constant / K 5.475E-3 A.3 RESIN CURE SHRINKAGE MODEL Vrs =0.0 a < aa V* = Aas + (vrs°> - A)oc2s acx<a< aC2 V*=Vr a>aC2 as= aC2 ac\ Variable Description Units Values v r > Resin volumetric cure shrinkage -v r S c o Total volumetric resin shrinkage from a = 0 to 1. - 0.099 «s The 'cure shrinkage' degree of cure. -aci Degree of cure after which the resin shrinkage begins. - 0.055 OtC2 Degree of cure after which the resin shrinkage stops. - 0.670 A Linear cure shrinkage coefficient - 0.173 A.4 RESIN MODULUS DEVELOPMENT MODEL E,=E: T <TC] r -T: K=K + 1 - ; c : (E;-K) rcl <r <rC2 E;=E; r>rC2 (A.4> ER=E'r[l + aEr(T-T0)] where = (Tga + Tgb * a) - 7 % = T*CXa + * T -134-Appendix A Material Properties of Hercules AS4/8552 Unidirectional Prepreg Variable Description Units Values a The resin degree of cure. -The resin elastic modulus at very low degrees of cure. Pa 4.67E6 Er°° The resin elastic modulus at T 0 and a = 1.0 Pa 4.67E9 T Resin temperature. °C T* Difference between resin temperature and resin instantaneous T g °C Tga Glass transition temperature at a = 0 K 2.68E2 Tgb Factor expressing the degree of cure dependence of resin glass transition temperature. K 2.20E2 Tela T above which resin modulus begins to increase at T = OK. K -45.7 Tcib Factor expressing the temperature dependence of the T above which resin modulus begins to increase. - 0.00 Tc2 T above which resin modulus has reached its full value. K -12 To Resin modulus development reference temperature. °C 20 a£r Factor expressing the temperature dependence of resin modulus. Pa/°C 0.00 A.5 RESIN POISSON'S RATIO DEVELOPMENT MODEL K = 2 \ Er(\-2v~) £00 (A.5) Variable Description Units Values v r Resin Poisson's ratio -E r Resin elastic modulus. Pa E 0 0 Resin elastic modulus at a = 1.0. Pa 4.67E9 Resin Poisson's ratio at a = 1.0. 0.37 -135-Appendix A Material Properties of Hercules AS4/8552 Unidirectional Prepreg A.6 THERMAL EXPANSION COEFFICIENT Lumped Ply Coefficients CTE, = 0.60 pl°C CTE2= CTE, = 28.6 pI'C -136-A P P E N D I X B : PLY MlCROMECHANICAL EQUATIONS The following equations are extracted from Johnston (1997). However, the original sources of these equations are the papers by Bogetti and Gillespie (1992) and Rosen and Hashin (1964). B. 1 CALCULATION OF ELASTIC CONSTANTS Given as inputs the transversely isotropic mechanical properties of the fibres (E\\f, E^f, G\y, v\y, and V23J), the properties of the isotropic resin (Er, vr) and the fibre volume fraction, Vf, ply mechanical properties are calculated as follows: B.1.1 IN-PLANE MODULI The longitudinal modulus is given by1 Eu=EufVf+Er(l-Vf) + 4(vr-vUf)2kfkrGr(l-Vf)Vf (kf+Gr)kR+(kf-kr)GrVf (B.l) where ' 2(1+ v , ) (B.2) and the transverse modulus is given by E22 — £33 1 (\/4kT) + (\/4G2i) + (v22/Eu) (B.2) There is a typographical error in the equation shown in the paper by Bogettie and Gillespie (1992): qyr-v?3f)kfkrGrQ-vf)Vf Eu =EufVf+EA\-V,) + l U f r 1 [(kf+Gr)kR+(kf-kr)GrVf original source, Rosen and Hashin (1964). . The correct equation is obtained from the •137-Appendix B Ply Micromechanical Equations B.1.2 SHEAR MODULI G12 = G U = G r (Gi3f+Gr) + (Gi3f-Gr)Vf (Guf+Gr)-(GUf-Gr)Vf (B.3) G2i = Gr [k, (GR + G23f) + 2 G 2 3 / G, + kr ( G 2 3 / - Gr )Vf ] kr(G23f +Gr) + 2G23fGr -(kr+2Gr)(G23f -Gr)Vf (B.4) where: G23f = 2(l + v 2 3 / ) (B.5) B.1.3 POISSON'S RATIOS 1^2 = v n = VHfVf+vrQ-Vf) + (vr-vuf)(kr-kf)Gr(\-Vf)Vf (kf+Gr)kr+(kf-kr)GrVf (B.6) 2EukT -EUE22 -Av\3krE 22 2Eukr (B.7) In the above, k is the so-called (by Bogetti and Gillespie) 'isotropic plane strain bulk modulus' defined by: k = 2(\-v-2v2) (B.8) and kr is the effective plane strain bulk modulus of the composite defined by: k Jkf+Gr)kr+{kf-kr)GrVf 1 (kf+Gr)-(kf-kr)Vf -138-Appendix B Ply Micromechanical Equations B.2 CALCULATION OF THERMAL AND CURE SHRINKAGE STRAINS Ply strains in the material principal directions arising from strains of the constituent resin and fibre are calculated from the fibre and resin mechanical properties and strains as follows: e,fEufVf+erEr(\-Vf) EufVf+Er(l-Vf) (B.10) £2 =e3 = (e2f+v,3fe,f)Vf +(£r +*/,*r)(1-*7) (B.H) -[vKfvf+vr{\-vfy\ slfElfVf +eREr(\-Vf) E,Vf+Er(\-Vf) To calculate thermal expansion coefficients (CTE,), replace resin and fibre strains by their respective thermal expansion coefficients. -139-A P P E N D I X C : RESIDUAL THERMAL STRESSES IN TWO CONCENTRIC TRANSVERSELY ISOTROPIC CYLINDERS The following equations are obtained from Wagner (1997). In this case, the fibre is transversely isotropic and the matrix is isotropic. The strain - stress relationship for transversely isotropic materials has the following form in cylindrical coordinates, 1 Er Ex 1 Er % Ex Vu 1 £, E E ag > AT (C.l) where A7" = (T -T0) and T , T0 are the current and initial temperatures. The stresses in the fibre are given by erf =Cf (C.2) where Af and A1 are constants and the stresses on the matrix are given by a"r' = Am + Bm(Rl lr)2 o-; =Am -Bm(RJr)2 (C.3) (J, = c -140-Appendix C Residual Thermal Stresses in Two Concentric Transversely Isotropic Cylinders where A"' ,Bmand C 'are constants, R\ is the fibre radius and r is any arbitrary radius on the matrix. The force balance in the longitudinal direction gives the relationship between these constants: Bm = — A" AJ V (C.4) c f _ -c V m m Now the two unknowns ,4mand C ' can be found by applying the no-slip condition at the interface and the strain/stress relationship given in Equation A. 18: Am=_KAMx-K2M2 ( t _ T q ) K^KA — K2K3 K.M2-K,M, ^ — — — J — L ( r - r 0 ) (C.5) where: K, =2 K2=-\ E\fVfJ f X 1 V. A — + E]f Vf j K3=-\ V 1 - V Em VfJ (C.6) E. M , = (a]f -or J M2 =(arf -aJ V/and Vm are fibre and matrix volume fractions respectively. -141-A P P E N D I X D : T O O L - P A R T INTERACTION D . l BI-METALLIC STRIP Consider a bi-metallic strip under a temperature change of AT as shown in Figure D . l . Due to the mismatch of the thermal expansion coefficient of both materials, the strip will expand and bend. But in the autoclave processing of a flat part on a tool, the bending of the tool-part combination is prevented by the application of high autoclave pressure. Hence, in this case only the expansion of the bi-metallic strip is considered. The expansion, 8, can be written as s=(t]Elg] +t2E2g2)AT (txE,+t2E2) where E is the Young's modulus, a is the coefficient of thermal expansion and t is the thickness. The normal stresses in strip 1 and 2, namely 07 and 02 respectively, are given by _ t2E]E2 (a2 - a] )AT °x (ttE]+t2E2) (D.2) _ r, E] E2 («, - a2 )AT (t^Ex+t2E2) D.2 TOOL-PART INTERACTION DURING HEATING UP Consider a composite part with 4 plies on an aluminum tool as shown in Figure D.2. During the heating up portion of the cure cycle, the degree of cure of the resin is very low. Hence the resin has a very low modulus. Due to the very low modulus of the resin, most of the stress due to tool--142-Appendix D: Tool-Part Interaction part interaction occurs at the first ply level and it drastically reduces with distance away from the first ply. The normal stress distribution across the part thickness may be approximated by a triangular distribution with maximum stress at the bottom ply as shown in Figure D.3. This maximum stress can be calculated by assuming first ply slip as shown in Figure D.2. According to Equation D.2, by assuming a unit thickness for the tool, the maximum stress in the part will be given by E E,(a, - a )A7/ CT, = P' V P—— (D.3) V ply EpJ + E/ ) where EpJ initial longitudinal modulus of the part Et modulus of the tool apJ initial longitudinal coefficient of thermal expansion of the part a, coefficient of thermal expansion of the tool tpiy thickness of a ply D.3 TOOL-PART INTERACTION DURING COOLING DOWN During the cooling down portion of the cure cycle, the resin is fully cured. Hence the resin modulus is high. Due to the high resin modulus, there is no slip between the plies and the whole part is involved in the interaction with the tool as shown in Figure D.4. Hence, the normal stress due to the tool part interaction is constant through the part thickness (Figure D.5). The normal stress in the part will be given by E fE,(a, -a f)AT = P-J ' v PJJ.— m4) (t F + F ~) where Epj final longitudinal modulus of the part apj final longitudinal coefficient of thermal expansion of the part -143-Appendix D: Tool-Part Interaction tiam thickness of the part D.4 TOOL REMOVAL Upon removal of the tool, the resultant normal stress in the part can be calculated by adding a; and (72 as shown in figure D.6. The part warps due to this resultant stress. The moment that causes the part to warp is given by M = F(^--Z) (D.5) where F is the resultant force and Z is the lever arm and given by F = t-^(2a2-cr]) Z = ^ ( 3 < X 2 - 2 < T , ) (D.6) The warpage of the part can be calculated from the moment-curvature relationship d2w _ M dx2 ~ EI (D.7) -144-Appendix D: Tool-Part Interaction D.5 FIGURES t k > f E„a, < iff E 2 , a2 > 1 Figure D. 1 Bi-metallic strip under temperature change, AT. Figure D. 2 First ply slip during the initial stage of curing process (heating up). The resin is not cured and it has very low modulus. -145-Appendix D: Tool-Part Interaction Normal Stress Figure D. 3 The normal stress distribution through part thickness during heating up. > > f • Part 1 i ] 1 1 Tool 1 > f Figure D. 4 Tool-part interaction at the final stage of the curing process (cooling down). The resin is fully cured and it has high modulus. The entire part is moved with the tool. -146-Appendix D: Tool-Part Interaction Normal Stress Figure D. 5 The normal stress distribution through part thickness during cooling down. Figure D. 6 Resultant normal stress distribution through part thickness. -147-
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Effect of processing on the behaviour of laminated...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Effect of processing on the behaviour of laminated composite structures : a numerical and probabilistic… Arafath, Abdul Rahim Ahamed 2002
pdf
Page Metadata
Item Metadata
Title | Effect of processing on the behaviour of laminated composite structures : a numerical and probabilistic approach |
Creator |
Arafath, Abdul Rahim Ahamed |
Date Issued | 2002 |
Description | The manufacturing of composite structures made of fibre-reinforced thermoset matrix materials generally results in residual stress build-up during the curing process. As a consequence, the precise geometry and material properties of the final structure after tool (mould) removal are often difficult to control. To account for such process-induced variations, engineers either have employed empirical knockdown factors, which are determined experimentally for a range of distinct structures, materials and manufacturing processes or have incorporated the measured variations as input into sophisticated finite element analyses. The first approach may lead to an unconservative design while the second approach is not suitable for normal design practice. In this study, a combined deterministic-stochastic simulation approach is presented to demonstrate the manner in which manufacturing-induced variabilities in composite structures can be controlled to achieve targeted reliability levels in structural performance. A finite element code, COMPRO, which deterministically simulates the various physical phenomena during manufacturing of composite structures, is integrated with non-linear structural analysis and reliability analysis to compute the statistics of the parameters that must be controlled at the manufacturing level in order to result in an optimum or reliable structural performance during service. The methodology is demonstrated through a case study that examines the buckling behaviour of a composite plate in the presence of manufacturing-induced imperfections. The intrinsic and extrinsic process parameters that affect the dimensional stability of more complex composite structures such as L- and C-shaped angular laminate sections, are numerically investigated using COMPRO. The predictive capability of COMPRO is examined by comparing the simulation results with available experimental measurements. Finally, to account for all sources of residual stress build-up required for accurate assessment of failure, the micro-level stresses arising from the mismatch between the coefficient of thermal expansions of the fibre and the matrix and the matrix cure shrinkage are computed through micromechanical finite element analysis of a representative volume of the fibre and the matrix. The ply- and laminate-level residual stresses computed by COMPRO during the curing process are then transferred to the micro-level using suitably derived stress magnification factors. |
Extent | 14617682 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2009-08-06 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0063921 |
URI | http://hdl.handle.net/2429/11974 |
Degree |
Master of Applied Science - MASc |
Program |
Civil Engineering |
Affiliation |
Applied Science, Faculty of Civil Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2002-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_2002-0007.pdf [ 13.94MB ]
- Metadata
- JSON: 831-1.0063921.json
- JSON-LD: 831-1.0063921-ld.json
- RDF/XML (Pretty): 831-1.0063921-rdf.xml
- RDF/JSON: 831-1.0063921-rdf.json
- Turtle: 831-1.0063921-turtle.txt
- N-Triples: 831-1.0063921-rdf-ntriples.txt
- Original Record: 831-1.0063921-source.json
- Full Text
- 831-1.0063921-fulltext.txt
- Citation
- 831-1.0063921.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0063921/manifest