Citation |

- Permanent Link:
- https://ufdc.ufl.edu/UFE0017266/00001
## Material Information- Title:
- Reliability-Based Design and Load Tolerance Evaluation Using Stochastic Response Surface and Probabilistic Sensitivities
- Copyright Date:
- 2008
## Subjects- Subjects / Keywords:
- Design analysis ( jstor )
Design engineering ( jstor ) Design optimization ( jstor ) Fatigue ( jstor ) Fatigue life ( jstor ) Mathematical variables ( jstor ) Polynomials ( jstor ) Random variables ( jstor ) Sensitivity analysis ( jstor ) Stress cycles ( jstor )
## Record Information- Source Institution:
- University of Florida
- Holding Location:
- University of Florida
- Rights Management:
- All applicable rights reserved by the source institution and holding location.
- Embargo Date:
- 3/1/2007
## UFDC Membership |

Downloads |

## This item has the following downloads: |

Full Text |

RELIABILITY-BASED DESIGN AND LOAD TOLERANCE EVALUATION USING STOCHASTIC RESPONSE SURFACE AND PROBABILISTIC SENSITIVITIES By HAOYU WANG A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2006 Copyright 2006 by Haoyu Wang To my family ACKNOWLEDGMENTS I would like to express my appreciation to my advisor, Professor Nam-Ho Kim, for his endless encouragement and continuous support during my Ph.D. research. Without his guidance, inspiration, experience and willingness of sharing his knowledge, this work would have never been possible. Dr. Kim made a tremendous contribution to this dissertation, as well as my professional and personal life. I must express my gratitude to the members of my supervisory committee, Professor Raphael T. Haftka, Professor Stanislav Uryasev, Professor Nagaraj K. Arakere, and Professor Ashok V. Kumar, for their willingness to review my Ph.D. research and provide constructive comments to help me complete this dissertation. Special thanks go to Professor Raphael T. Haftka for not only his guidance with several technical issues during my study, but also comments and suggestions during group meetings which were extremely helpful for improving my work. Special thanks are also given to Professor Nestor V. Queipo from the University of Zulia at Venezuela, for his interaction in my research and collaboration in publishing papers during his visiting at the University of Florida. My colleagues in the Structural and Multidisciplinary Optimization Lab at the University of Florida also deserve my gratitude. In particular, I thank Dr. Xueyong Qu, Dr. Amit A. Kale, Dr. Erdem Acar, Tushar Goel, Long Ge, Saad Mukrus for their encourage and help. My parents deserve my deepest appreciation for their constant love and support. Lastly, I would like to thank my beautiful and lovely wife, Zhilan, for her love, patience and support during my study. TABLE OF CONTENTS A C K N O W L E D G M E N T S ................................................................................................. iv LIST OF TABLES ..................................... .. .......... ...................................... ix LIST OF FIGURES ........................................ .............. xi ABSTRACT ............................................ ............. ................. xiii CHAPTER 1 IN T R O D U C T IO N .................................................................. .. ... .... ............... 1 M motivation ..................................................................................................... ...............1 Obj active ............................................. ............................... 2 Scope .................................................. .......................3 O u tlin e ........................................................................... .................................... . 5 2 LITERA TURE SU RVEY ..........................................................................................7 Uncertainty and Reliability Analysis of Structural Applications ..............................7 R eliability-B ased D esign O ptim ization........................ ...........................................11 Sensitivity in R liability A analysis ...........................................................................12 D im ension R education Strategy ................................................................................13 R obust D design ............................................................................... ...................... 13 F atigu e L ife P reduction .......................................................................... ...............14 3 UNCERTAINTY ANALYSIS USING STOCHASTIC RESPONSE SURFACE .... 19 In tro d u ctio n ............... .................................................... ..................... 1 9 Description of Uncertainty M odel ...........................................................................20 Stochastic Response Surface Method (SRSM)...................................................22 Polynomial Chaos Expansion (PCE) in Gaussian Space .................................22 Numerical Example of Stochastic Response Surface.......................................26 Improving Efficiency of SRS Using Local Sensitivity Information ..........................30 Continuum-Based Design Sensitivity Analysis ................................................31 Constructing SRS Using Local Sensitivity.......................................................34 Numerical Example Torque Arm M odel.......................................................36 Sum m ary ......................................................................................... ............. ....... 37 4 RELIABILITY-BASED DESIGN OPTIMIZATION...........................................39 G general R B D O M odel .................. ........ ... .............. ............. ........................... 39 Reliability Index Approach (RIA) and Performance Measure Approach (PMA)......41 Probability Sensitivity A analysis (PSA ) ................................................. ................ 43 Probability Sensitivity Analysis in FORM ..................................... ................ 44 Probability Sensitivity Analysis Using SRSM...............................................47 Reliability-Based Design Optimization Using SRSM...........................................48 R B D O w ith R IA .................................................................................................. 4 9 R B D O w ith Inverse M measure ......................................................... ................ 52 S u m m a ry .................................................................................................................. .. 5 3 5 GLOBAL SENSITIVITY ANALYSIS FOR EFFICIENT RBDO......................... 54 In tro d u ctio n ................................................................................................................. 5 4 Sensitivity A analysis ............................................................. .. .... ............... ........ ..... 55 Variance-Based Global Sensitivity Analysis (GSA).............................................56 Global Sensitivity Analysis Using Polynomial Chaos Expansion ..........................58 Adaptive Reduction of Random Design Space Using GSA in RBDO.................... 59 S u m m a ry ..................................................................................................................... 6 5 6 FATIGUE RELIABILITY-BASED LOAD TOLERANCE DESIGN................... 66 In tro d u ctio n ................................................................................................................. 6 6 F atigu e L ife P reduction ............................................... .......................................... 67 Crack Initiation Fatigue Life Prediction ............... .............. ..................... 68 Variable Amplitude Loading and Cumulative Damage .................................70 Model Preparation for Fatigue Reliability Analysis ........................ ..................... 71 Finite Elem ent M odel .. .. ................. ............................................... 71 D ynam ic L oad H history ................................................................... ................ 73 Uncertainty in Material Properties and S-N Curve Interpolation .....................74 Uncertainty Modeling of Dynamic Loadings........................................................75 Linear E stim ation of Load Tolerance.................................................... ................ 76 Variability of Dynamic Load Amplitude ............... ...................................77 V ariability of M ean of Dynam ic Load ........................................... ................ 82 Safety Envelope Concept for Load Tolerance Design ..........................................84 Numerical Path Following Algorithm................ ...................................85 Example for Multi-Dimensional Load Envelope ...........................................88 C conservative D distribution Type............................................................. ................ 90 S u m m a ry .................................................................................................................. ... 9 2 7 ROBUST DESIGN USING STOCHASTIC RESPONSE SURFACE................... 94 In tro d u ctio n .. ..................... ... .... ............................................... ..................... 9 4 Performance Variance Calculation Using SRS .............. .....................................96 Variance Sensitivity ................................ ......... ......................97 Robust Design Two-layer Beam...... .......... ........ ..................... 102 Dynamic Response of Two-Layer Beam ..................................... ................ 102 Robust D esign for Tw o-Layer Beam ....... .......... ....................................... 103 Global Sensitivity A analysis ...... ............. ............ ..................... 106 Robust Design by Tolerance Control ....... ... .... ...................... 107 S u m m ary ................................................................................................ .... ........... 1 1 1 8 SUMMARY AND RECOMMENDATIONS .....................................................113 APPENDIX A SAMPLING-BASED PROBABILITY SENSITIVITY ANALYSIS FOR DIFFERENT DISTRIBUTION TYPE...... .... .........................115 Normal Distribution X, N (,/ o-2) ................................1... 15 C ase 1: 0 = u .............. ............................................... ................. .......... 115 C ase 2 : 0 = cr .............. ................................................................ ........ 116 U uniform D distribution .............. .................. ................................................ 117 L og-N orm al D distribution ...................................... ......................... ............... 119 B NATURAL FREQUENCY OF CANTILEVER COMPOSITE BEAM............... 122 B en d in g M o m en t ..................................................................................................... 12 2 Geom etric Properties of Composite Beam ....... .......... ...................................... 122 Effective Compliance for Composite Beam ...... .... ..................................... 123 Effective Mass for Composite Beam...... .... ...... ..................... 123 LIST O F REFEREN CE S .. .................................................................... ............... 125 BIOGRAPH ICAL SKETCH .................. ............................................................. 134 LIST OF TABLES Table page 3-1. The type of polynomials and corresponding random variables for different Askey- Chaos (N >0 denotes a finite integer)................................................... ................ 22 3-2. Root mean square error of PDF compared with the exact PDF of performance fu n ctio n y = ex.............. ....................................................................... 2 7 3-3. Comparison of probability of G>520MPa obtained from different uncertainty analysis methods (Full sampling without using local sensitivity).........................30 3-4. Comparison of probability of G>520MPa obtained from different uncertainty analysis methods (reduced sampling using local sensitivity) ; .....................37 3-5. Comparison of probability of G>520MPa obtained with/without local sensitivity (7/27 sam pling points) using 2nd order SRS........................................ ................ 37 4-1. Probability sensitivity with respect to random parameters (unit: centimeter)............46 4-2. Computational efficiency of analytical method for probability sensitivity .............47 4-3. Definition of random design variables and their bounds. The values of design variables at optimum design are shown in the 5th column (unit: centimeter). ........50 4-4. Reliability Index of active constraint at optimal design........................................52 5-1. Variances of the Hermite bases up to the second order.........................................58 5-2. Global sensitivity indices considering only main factors for the torque arm model at the initial design. Only three random variables (u2, u6, and us) are preserved w hen a threshold value of 1.0% is in place ........................................ ................ 63 5-3. Comparison of the number of random variables in each design cycle. The threshold of 1.0% is used. The first constraint is listed ................. ..................... 64 6-1. Q quality of response surface ........................................ ........................ ................ 78 6-2. T -statistic of the coefficients ....................................... ....................... ................ 79 7-1. Random variables for cantilevered beam structure ............................... ................ 99 7-2. Variance estimation of linear performance (strength)................... .................. 100 7-3. Variance estimation of nonlinear performance (deflection)............................... 101 7-4. Sensitivity of variance for linear performance (strength)...................................102 7-5. Sensitivity of variance for nonlinear performance (deflection)............................. 102 7-6. Random parameters for the composite beam structure ................. ...................104 7-7. Sensitivities of objective functions at the initial design (ts = 6tm, tp = 0.2[tm, L = 1000[tm ) ........................................................................................................ 105 7-8. Total sensitivity indices for the composite beam structure (ts = 6|tm, tp = 0.2[tm, L = 10 0 0[tm ) ............................................................................. . ... ............... 10 7 7-9. Sensitivity of variance for linear performance (strength)...................................109 7-10. Sensitivity of variance for nonlinear performance (deflection)........................... 109 7-11. Random variables and cost-tolerance functions....... .................. .................. 110 7-12. Random variables and cost-tolerance functions...................................................111 A-1: Accuracy of proposed probability sensitivity method for normal distribution using 200,000 sam pling M CS ....... ........... ............ ..................... 116 A-2: Accuracy of proposed probability sensitivity method for uniform distribution using 200,000 sam pling M CS ....... ........... ............ ..................... 119 A-3: Accuracy of proposed probability sensitivity method for Log-normal distribution using 200,000 sam pling M CS ....... ........... ............ ..................... 121 LIST OF FIGURES Figure page 3-1. Limit state function divides the safe region from the failure region .......................21 3-2. PDF of performance function y(x) = ex ..........................................27 3-3. Shape design param eters for the torque arm ......................................... ................ 28 3-4. PDF of performance function G(x) -torque arm model......................................29 3-5. Variation of a structural domain according to the design velocity field V(x)............32 3-6. PDF of performance function G(x) Torque model at initial design (SRS with se n sitiv ity ) ................................................................................................................ 3 6 4-1. Flow chart for reliability-based design optimization.............................................43 4-2. Optimum design and stress distribution of the torque arm model with 8 random v a ria b le s ............................................................................................................. .. 5 1 4-3. Optimization history of cost function (mass) for the torque arm model with 8 ran d om v ariab les. ..................................................................................................... 5 1 4-4. PDF of the performance function at the optimum for the torque-arm problem .........52 5-1. Global sensitivity indices for torque arm model at initial design.............................. 59 5-2. Adaptive reduction of unessential random design variables using global sensitivity indices in RBDO. Low-order SRS is used for global sensitivity analysis, while a high-order SRS is used to evaluate the reliability of the system. .61 5-3. Optimum designs for the full SRS (solid line) and adaptively reduced SRS (dotted line). Because some variables are fixed, the interior cutout of the design from the adaptively reduced SRS is larger than that from the full SRS. ........................64 6-1. Flow chart for fatigue life prediction..................................................... ................ 67 6-2. R ain-fl ow and hysteresis .......................................... .......................... ................ 70 6-3. Front loader frame of CAT 994D wheel loader (subject to 26 channels of dynamic lo a d in g ) .................................................................................................................. ... 7 2 6-4. Finite elem ent m odel for front fram e .................................................... ................ 73 6-5. M material S-N curve w ith uncertainty...................................................... ................ 74 6-6. Illustration of one channel of dynam ic loads......................................... ................ 75 6-7. Reliability index /Twith respect to random parameter Iy ................ ............... 81 6-8. Probability of failure Pf with respect to random parameter [y .............................. 81 6-9. Reliability index / with respect to random parameter a ................. ..................... 84 6-10. Probability of failure Pf with respect to random parameter [a ................................84 6-11. Safety envelope for tw o variables ....................................................... ................ 85 6-12. Predictor-corrector algorithm ..................................... ...................... ................ 86 6-13. C construction of load envelope............................................................. ................ 89 6-14. Safety envelop for fatigue reliability of CAT 994D front loader frame................90 6-15. Reliability index / with respect to random parameter a ....................................... 91 6-16. Probability of failure Pf with respect to random parameter a ............................... 91 6-17. 2-D safety envelope for different distribution type with same random parameters .92 7-1. Cantilever beam subject to two direction loads..................................... ................ 99 7-2. Piezoelectric cantilevered composite beam.......... ...................................... 103 7-3. Pareto optimal front for the robust design of the composite beam........................ 106 B-1: Free body diagram of two-layer beam............... ..........................122 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy RELIABILITY-BASED DESIGN AND LOAD TOLERANCE EVALUATION USING STOCHASTIC RESPONSE SURFACE AND PROBABILISTIC SENSITIVITIES By Haoyu Wang December 2006 Chair: Nam-Ho Kim Major: Mechanical Engineering Uncertainty is inevitable in structural design. This research presents an efficient uncertainty analysis technique based on stochastic response surfaces (SRS). The focus is on calculating uncertainty propagation using fewer number of function evaluations. Due to sensitivity analysis, the gradient information of the performance is efficiently calculated and used in constructing SRS. Based on SRS, reliability-based design optimization (RBDO) is studied intensively in this research. Probability sensitivity analysis using the sampling technique is also proposed. Since the computational cost of RBDO increases significantly proportional to the increasing number of random variables, global sensitivity analysis is introduced to adaptively reduce unessential random variables. It has been shown that the global sensitivity indices can be calculated analytically because the SRS employs the Hermite polynomials as bases. Traditional structural design focuses on designing a reliable structure under well characterized random factors (dimensions, shape, material properties, etc). Variations of these parameters are relatively small and well characterized. However, everyday engineering life tends to use the existing structural part in a different applications instead of designing a completely new part. In this research, a reliability-based safety envelope concept for load tolerance is introduced. This shows the capacity of the current design as a future reference for design upgrade, maintenance and control. The safety envelope is applied to estimate the load tolerance of a structural part with respect to the reliability of fatigue life. Stochastic response surface is also applied on robust design in this research. It is shown that the polynomial chaos expansion with appropriate bases provides an accurate and efficient tool in evaluating the performance variance. In addition, the sensitivity of the output variance, which is critical in the mathematical programming method, is calculated by consistently differentiating the polynomial chaos expansion with respect to the design variables. A reliability-based robust design method that can reduce the variance of the output performance as well as the deviation of the mean value is proposed using SRS and efficient sensitivity analysis. Numerical examples are shown to verify accuracy of the sensitivity information and the convergence of the robust design problem. CHAPTER 1 INTRODUCTION Motivation A typical mechanical design procedure includes two steps: first, a design space is defined and a mathematical model is established, which includes the objective function and required constraints. Second, a proper optimization algorithm is selected properly based on this mathematical model to solve the design problem. In engineering design, the deterministic optimization model has been studied intensively to reduce the objective function by pushing design to the limits of system failure boundaries. However, everything in the real world involves uncertainties, and so does the design of mechanical components. After realizing deterministic design leaves very little or no room for tolerances of the imperfections in design, manufacturing and variety of service conditions, design engineers incorporate a safety factor into the structural design to leave safety margins. Without considering uncertainties and probabilistic quantification, deterministic design with a safety factor may be either unsafe or too conservative. Motivated by overcoming the bottleneck of the deterministic design, the reliability- based design optimization (RBDO) model has become popular in past two decades since uncertainties exist everywhere in every phase of the structure system, from design and manufacturing to service and maintenance. If elements in the mathematical model are considered to be probabilistic with certain types of random distribution, the design problem becomes a typical RBDO problem. The probabilistic elements can be design variables, material properties, applied loads, etc. One of the most important issues in RBDO is a good model of uncertainty propagation in mathematical models. Besides RBDO, which only considers the failure mode as a constraint in the probabilistic point of view, robust design will also be considered in this research in order to design a structure 'less sensitive' to the existing uncertainty factors. In optimization point of view, that means minimization of performance variance. For a certain design, it is also important to consider the service capability of the system subject to applied loads since engineers tend to use the same design in different applications instead of a completely new design. Another motivation of this research is the load tolerance design. A good estimation of load tolerance shows the capacity of the current design, future reference for design upgrade, maintenance and control. Since static or quasi-static loading is rarely observed in modern engineering practice, the majority of engineering design projects involves machine parts subjected to fluctuating or cyclic loads. Such loads induce fluctuating or cyclic stresses that often result in failure by fatigue. In addition, because service loads are subjective, which means the load characteristic of one operator may be completely different from that of the other, it is necessary to consider the uncertainties while estimating the load tolerance of dynamic systems. Objective Uncertainty in the design parameters makes structural optimization a computationally expensive task due to the significant number of structural analyses required by traditional methods. Critical issues for overcoming these difficulties are those related to uncertainty characterization, uncertainty propagation, ranking of design variables, and efficient optimization algorithms. Conventional approaches for these tasks often fail to meet constraints (computational resources, cost, time, etc.) typically present in industrial environments. The first objective of this research is to develop a computationally efficient method for uncertainty propagation. Local and global sensitivities can then be used to improve the efficiency of estimating uncertainty propagation. Besides efficiency, the accuracy and applicability of the methods to a wide range of applications need to be addressed. The second objective of this research is to develop a computationally efficient RBDO and robust design framework based on proposed uncertainty analysis. In the gradient-based algorithm, the sensitivity information is required during the optimization procedure. The computational cost can be significantly saved if the gradient can be obtained analytically, instead of using the finite difference method. The probabilistic sensitivity analysis is utilized to calculate the gradient of the reliability constraints. In the framework of robust design, sensitivity analysis of performance variance is also studied. Traditional structural design usually makes assumption on randomness of factors involved in modeling a structural system such as design variables, material properties, etc. However, it is also important to consider the capacity of the system subject to uncertain loads. The final objective of this research is to present a reliability-based load design method, which provides the safety envelope, for a structure subject to fatigue failure. Scope In the standard framework of RBDO, constraints are provided in terms of the probabilities of failure. The uncertainties involved in the system are modeled by assuming random input variables with a certain type of probabilistic distribution. RBDO achieves the design goal by meeting the requirement of structural reliability constraints. The RBDO involving a computationally demanding model has been limited by the relatively high number of required analyses for uncertainty propagation during the design process. The scope of this research is to present an efficient uncertainty propagation technique based on stochastic response surfaces (SRS) constructed using model outputs at heuristically selected collocation points. The efficiency of the uncertainty propagation approach is critical since the response surface needs to be reconstructed at each design cycle. In order to improve the efficiency, the performance gradient, calculated from local sensitivity analysis, is used. Even if the local sensitivity information can reduce the number of required simulations, the dimension of the SRS is still increased according to the number of random variables. If the contribution of a random variable is relatively small to the variance of the model output, it is possible to consider the random variable as a deterministic one. In this research, the global sensitivity index is used for that purpose. The role of the global sensitivity is to quantify the model input's contributions to the output variability, hence establishing which factors influence the model prediction the most so that i) resources can be focused to reduce or account for uncertainty where it is most appropriate, or ii) unessential variables can be fixed without significantly affecting the output variability. Reliability constraint in RBDO requires probability sensitivity analysis for gradient-based algorithms. In this research, both FORM-based and sampling-based reliability sensitivity analysis are investigated. The analytical expression for probability sensitivity based on SRS is derived and used for RBDO. Variations in dynamic loads are usually too complicated to be predicted. A simplified uncertainty modeling technique based on the mean and amplitude of the load history is proposed. Using the uncertainty in the load history, a reliability-based safety envelope is constructed that can provide load tolerance of the current design. In addition, the effect of different distribution types is investigated so that the design engineers can choose the conservative distribution type. This research involves uncertainty modeling and quantification, design sensitivity analysis, fatigue life prediction, reliability-based design optimization (RBDO) and robust design. Methodologies investigated or applied in reliability analysis include moment- based methods such as first- and second-order reliability method (FORMISORM), approximation methods such as Monte Carlo Simulation (MCS) with stochastic response surface method (SRSM). Furthermore, sensitivity analysis for reliability constraints of RBDO is investigated to improve the computational cost involved in reliability analysis and design. Performance variance and sensitivity are calculated based on SRSM for robust design. Computationally affordable reliability-based optimization and robust design method, and safety envelope for load tolerance are presented in this work. Outline A literature survey is presented in chapter 2, which includes all aspects involved in this research such as reliability analysis, reliability-based design optimization, robust design, sensitivity analysis, dimension reduction strategy and fatigue analysis. Chapter 3 describes the uncertainty modeling and widely used reliability analysis methods. A stochastic response surface method (SRSM) coupled with the sensitivity analysis of performance measure is introduced. It is shown that the local sensitivity information improves computational efficiency significantly by reducing required number of samples. Convergence and accuracy of the proposed SRSM scheme are also discussed in this chapter. In Chapter 4, the mathematical model is defined for RBDO. RBDO using either direct probability measure or inverse measure is investigated and compared. The difference of numerical procedures between RBDO and deterministic optimization are also compared. As required by RBDO, probability sensitivity analysis is studied in this chapter. In Chapter 5, a dimension reduction strategy is proposed by introducing the concept of variance-based global sensitivity analysis, which saves the computational resources further by fixing the unessential design variables. Chapter 6 demonstrates a fatigue reliability-based load tolerance design by using reliability sensitivity information. A reliability-based safety envelope is constructed by path following continuation method. Chapter 7 proposes an optimization model for robust design where SRS is used to calculate the performance variance and its sensitivity. Chapter 8 concludes this research followed by recommendations for future research work. CHAPTER 2 LITERATURE SURVEY Uncertainty and Reliability Analysis of Structural Applications Reliability-based design optimization (RBDO) provides tools for making decision within a feasible domain of design variables with consideration of uncertainties underneath. In the past decades, tremendous amount of work has been carried out in this area and it is still moving forward. Compared to deterministic optimization, design variables included in RBDO are random and usually modeled with specific distribution types, so do the random parameters such material properties as Young's modulus and Poisson's ratio. Usually random parameters do not change during optimization, but their effects to the probability propagation must be counted due to its uncertainty. Reliability, which is defined as the probability that a system response does not exceed the limit threshold, is often considered as constraints. The system response is a function of design variables and random parameters, which is called a performance function in this research. Performance function is usually implicit and nonlinear prediction of random variables, making probabilistic description of a system response difficult. Several approximation methods for reliability analysis have been developed in the literatures. Among them, Monte Carlo Simulation (MCS) (Metropolis and Ulam 1949; Rubinstein 1981) has been widely used due to its simplicity and dependability. However, the large sample size required in MCS in order to reduce the noise and to reach a certain level of accuracy makes it practically formidable in computationally intensive engineering applications, such as Finite Element Analysis (FEA). Even improved version of MCS are developed, such as importance sampling, Latin Hypercube Sampling (Wyss and Jorgensen 1998), Stratified Sampling, etc, they are still expensive in structural reliability analysis. Moment-based methods (Breitung 1984; Haldar and Mahadevan 2000; Hasofer and Lind 1974) have been developed to provide less expansive calculation of the probability of failure compared to MCS. However, they are limited to a single failure mode. As the most widely used moment-based methods, the development of the theory of First- and Second-Order Reliability Method (FORMISORM) is claimed to be finished and only technical work left to do (Rackwitz 2000). FORM/SORM are based on the linear/quadratic approximation of the limit state function around most probable point(MPP), which is defined in standard normal space as the closest point from the origin on the response surface. For highly nonlinear problems, predictions of reliability from FORMISORM are not accurate enough because they approximate the response using a linear or quadratic function. The response surface method (Khuri and Cornell 1996; Myers and Montgomery 1995) is proposed to resolve this difficulty. This method typically employs polynomials bases to approximate the system performance and facilitate reliability analysis with little extra computational cost by combining with MCS. Since the accuracy of MCS with fixed sample size relies on the seeking level of probability of failure which sometimes is extremely low in structural design, the probability calculated by MCS near optima is too rough to represent the true value of failure probability. Reliability analysis using safety factor (Wu et al. 2001) or probability sufficiency factor (PSF) (Qu and Haftka 2004) is proposed to ameliorate this effect. With the PSF as the constraints in RBDO, the variation of magnitude of constraints is usually several orders of magnitude lower than that of the probability of failure, and so is the magnitude of the numerical noise caused by MCS. One of the significant advantages of the moment-based approach is that the sensitivity of the system reliability or probability of failure can be obtained with very little extra computation (Yu et al. 1998). However, moment based approach such as FORM/SORM still has limitations when the performance function is highly nonlinear (Ghanem and Ghiocel 1996). The evaluation of the probabilistic constraints may have large errors in this case. Mahadevan and Shi (Mahadevan and Shi 2001) presented a multipoint linearization method (MPLM) for the reliability analysis of nonlinear limit states, which determines the multiple linearization points through the secant method. It increases the complexity of the problem with limited accuracy improvement. The response surface method can approximate the system response and with little extra computation for MCS, the probability of failure can easy to be obtained. Compared to the conventional deterministic design response surface, Stochastic Response Surface (SRS) (Isukapalli et al. 1998) has the advantage that it only approximates the function around most probability region which highly improved accuracy. Another advantage of SRS is the choice of basis function. The monomial bases(Qu et al. 2000) are widely used due to its simplicity. Other polynomial bases are also being studied intensively such as radial basis function (RBF) (Krishmamurthy 2003), orthogonal polynomials (Xiu et al. 2002),etc. Since Ghanem and Spanos proposed the spectral approach of stochastic finite element method (Ghanem and Spanos 1991), the homogeneous Polynomial Chaos Expansion (PCE) has been widely utilized to represent the uncertainties due to the nature of stochastic process. To make better approximation with less model analyses, sampling methods are studied intensively. Different sampling methods were studied and brought in different applications recent years, such as Latin Hypercube Sampling (LHS) (Choi et al. 2003; Qu et al. 2000) and collocation sampling method (Webster et al. 1996). In the collocation method, Webster and Tatang derived a set of polynomials from the probability density function of each input parameter such that the roots of each polynomial are spread out over the high probability region of the parameter by deriving orthogonal polynomials. Because the uncertainty is usually evaluated by transforming all the random variables and parameters into the Gaussian space, the corresponding orthogonal polynomials are Hermite polynomials. To obtain additional accuracy of SRS, moving least square (MLS) method (Youn 2001, Dec; Youn and Choi 2004) is proposed by introducing weight functions. The number of simulations can be reduced if the sensitivity information is available. Isukapalli (Isukapalli et al. 2000) used an automatic differentiation program to obtain the sensitivity and utilized it in constructing the response surface. However, the computational cost for automatic differentiation is usually very high(Van Keulen et al. 2004), which reduces the significance of the method. Design sensitivity analysis can provide analytical sensitivity information of response with little extra computation (Kim et al. 2000). Thus, coupling the regression based stochastic response surface method (SRSM) with sensitivity can save large amount of computational cost, especially when the required number of design variables is large (Kim et al. 2004b). Several methods(Lauridensen et al. 2001; Malkov and Toropov 1991; Rijpkema et al. 2000; Van Keulen et al. 2000) have been proposed to use sensitivity information in constructing response surface. Vervenne (Vervenne 2005) proposed a gradient-enhanced response surface method based on above mentioned methods. He developed a two-step approach is proposed: first, different response surfaces using function values and derivatives are constructed separately; Second, these response surfaces are combined together to form a single response surface which fits as good as possible for both function value and response surfaces. In his study, several types of response surface and different combination scheme have been compared. Reliability-Based Design Optimization As mentioned in the previous section, FORM/SORM performs reliability analysis through linear/quadratic approximation of the performance function at MPP. Thus, searching MPP is the main task for moment-based RBDO. However, most advanced MPP search methods such as two point adaptive nonlinear approximation method (TPA) (Grandhi and L.P. 1998; Wang and Grandhi 1995; Xu and Grandhi 1998) or hybrid mean value (HMV) method(Youn 2001, Dec; Youn et al. 2003) can not make significant improvement of efficiency in the computational cost(Du and Chen 2002b). In conventional RBDO, the probability constraint is described by the reliability index, which in FORM is the shortest distance from the origin to the limit state in standard normal space. This approach is called reliability index approach (RIA). By modifying the formulation of probabilistic constraints, Tu proposed an inverse measure approach, called Performance Measure Approach (PMA) (Tu 1999; Tu and Choi 1997; Tu et al. 1999; Tu 2001) which is proved to be consistent with the RIA but is inherently robust and more efficient if the probabilistic constraint is inactive. Both RIA and PMA employ double loop strategy with analysis loop (inner loop for reliability analysis) nested within the synthesis loop (outer loop for design optimization). Due to the nature of double loop optimization, the computational cost is usually high. A couple of new strategies were proposed to improve the efficiency (Yang and Gu 2004). Sequential Optimization and Reliability Assessment (SORA) method (Du and Chen 2002b) decouples optimization loop from the reliability analysis loop and each deterministic optimization loop followed by a series of MPP searches. This method shifts the boundaries of violated constraints to the feasible direction based on the reliability results obtained in the previous cycle. Thus it improves design quickly from cycle to cycle and ameliorates the computational efficiency. Other single loop methods(Chen et al. 1997; Kwak and Lee 1987; Liang et al. 2004; Wang and Kodiyalam 2002) are also developed to provide efficient RBDO. In this method, the relationship between random variables and its mean is found through the Karush-Kuhn-Tucher (KKT) optimality condition. The double loop RBDO formulation is transformed to a single loop deterministic optimization problem and expensive MPP search is avoided. However, there is no guarantee that an active reliability constraint converges to its own MPP, and the required reliability may not be satisfied. Sensitivity in Reliability Analysis When RBDO problems are solved using gradient-based optimization algorithms, sensitivities of reliability or probability of failure with respect to the design parameters are required. Probability sensitivity can be used to identify those insignificant random variables during the design stage. In the moment-based approaches such as FORM, the sensitivity can be obtained accompanied by the reliability analysis without extra function evaluation once MPP is located (Karamchandani and Cornell 1992; Yu et al. 1997). Wu(Wu 1994) proposed an adaptive importance sampling(AIS) method to calculate reliability and AIS-based reliability sensitivity coefficients. Liu et al(Liu et al. 2004) compare four widely-used probability sensitivity analysis(PSA) methods, which include Sobol' indices, Wu's sensitivity coefficients, the MPP based sensitivity coefficients and the Kullback-Leibler entropy based method. The merits behind each method are reviewed in details. Dimension Reduction Strategy In reliability analysis, the computational cost of multidimensional integration is high. Xu and Rahman(Rahman and Xu 2004; Xu and Rahman 2004) use series expansions to decompose the multidimensional problem to lower dimensional integration, such as univariate and bivariate integration. Compared to multidimensional integration, the total computation of univarate integration is much lower. Recent development in statistics introduces global sensitivity analysis (GSA)(Saltelli et al. 2000; Saltelli et al. 1999; Sobol 1993; Sobol 2001), which studies how the variance in the output of a computational model can be apportioned, qualitatively and quantitatively, to different sources of variation. Considering the contribution of the variance of design variables to performance variances are not of same importance, Kim et al proposed an adaptive reduction method using total sensitivity indices to reduce the problem dimensions (Kim et al. 2004a). Robust Design Robust design, known as Taguchi parameter design (Taguchi 1986; Taguchi 1987), is to design a product in such a way that the performance variance is insensitive to variation of design variables which is beyond the control of designer. Wang & Kodiyalam(Wang and Kodiyalam 2002) formulated robust design as an optimization problem by minimizing the variation of system response. Since the material cost has to be considered as well as manufacturing cost, Chen and Du's formulation compromises cost reduction with performance variance control(Du and Chen 2002a). A robust design can also be achieved by using traditional optimization techniques to minimize the performance sensitivities. Chen & Choi formulated the robust design by minimizing a total cost function and sum of squares of magnitudes of first-order design sensitivities(Chen and Choi 1996), which requires the evaluation of second-order sensitivity analysis. This is a different philosophy compared to the variance based approach. It is more focus on the local behavior of the system performance and can achieve local robustness. The final design by minimizing local sensitivity cannot guarantee the robustness of system globally if the input variances are considerable. By summarizing approaches popularly applied in robust design, Park et al. (Park et al. 2006) define robust design methodologies into two different category: Taguchi method and robust optimization. Under the context of multi-scale and multi-disciplinary applications, Allen et al.(Allen et al. 2006) reviewed robust design methods and categorizes robust design into four different types based on the sources of variability. Fatigue Life Prediction In 1829, Albert found that a metal subjected to a repeated load will fail at a stress level lower than that required to cause failure on a single application of the load. Then the question comes out: how parts fail under time-varying or non-static conditions? Such phenomenon is called fatigue. The first approach developed to carry out fatigue analysis is the nominal stress method, which is still widely used in applications where the applied stress varies with constant amplitude within the elastic range of the material and the number of cycles to failure is large. The nominal stress method works well in high cycle fatigue analysis but does not fit for the low cycle fatigue analysis where the material has a significant part in the plastic region. August Wohler(Wohler 1860) carried out experiments to obtain a plot of cyclic stress level versus the logarithm of life in mid-19th century, which is well known as S-N curve. Basquin proposed a stress-life(S-N) relationship(Basquin 1910) which can be plotted as a straight line using log scales. S-N approach is applicable to situations where cyclic loading is essentially elastic, so the S-N curve should be confined on the life axis to numbers greater than about 105 cycles in order to ensure no significant plasticity occurs. Most basic fatigue data are collected in the laboratory by testing procedures which employ fully reversed loading. However, most realistic service loads involve non-zero mean stresses. Therefore, the influence of mean stress on fatigue life should be considered so that the fully reversed laboratory data can be used in the evaluation of real service life. Since the tests required to determine the influence of mean stress are quite expensive, several empirical relationships(Gerber 1874; Goodman 1899; Soderberg 1939) which related alternating stress amplitude to mean stress have been developed. Among the proposed relationships, two are widely used, which are based on Goodman(Goodman 1899) and Gerber(Gerber 1874). S-N approach works well when the cyclic loading is essentially elastic, which means in high cycle fatigue life evaluation. While using this method, it assumes that most of the life is consumed by nucleating cracks (around 0.01 mm) and nominal stresses and material strength control fatigue life. Accurate determinations of miscellaneous effects factor Kf for each geometry and material are also required. The advantage of S-N approach is apparent since changes in material and geometry can easily be evaluated and large empirical database for steel with standard notch shape is available. However, the limitation should also be accounted. This method does not consider the effects of plasticity, and mean stress effect evaluation is often in error. As the matter of fact, the requirement of empirical Kf for good results is also a kind of disadvantage. As mentioned above, when the cyclic loads are relatively large and have a significant amount of plastic deformation, the components will suffer relatively short lives. This type of fatigue behavior is called low-cycle fatigue or strain-controlled fatigue. The analytical procedure in dealing with strain-controlled fatigue is called the strain-life, local stress-strain or critical location approach. In 1950's, Coffin and Manson(Coffin 1954; Manson 1954) suggested that the plastic strength component of a fatigue cycle may also be considered in fatigue life prediction by a simple power law. In order to account for the mean stress effects, two correction methods are proposed by Morrow and Smith, Topper & Watson (STW)(Smith et al. 1970), respectively. Local Strain-Life (e-N) method assumes that the local stresses and strains control the fatigue behavior. In this method, the plastic effects and mean stress effects are considered well. The limitation is that it also needs the empirical Kf. In the local strain- life approach, the most of the life is consumed by micro-crack growth (0.1-1mm). To account for macro-crack growth(>lmm), the fracture mechanics-based crack propagation method is proposed(Hoeppner and Krupp 1974; Paris 1964; Paris and Erdogan 1963). In this method, major assumption is that nominal stress and crack size control the fatigue life and the initial crack size is determined accurately. It is the only method to directly deal with cracks. However, the complex sequence effects and accurate initial crack size are difficult to be determined. Linear elastic fracture mechanics (LEFM) is a new branch of engineering. The earliest work was done by Inglis(Inglis 1913) but the major developments were carried out by Griffth(Griffith 1921) at Royal Aircraft Factory(RAF,UK) in 1921 and Irwin(Irwin 1956) in the USA in 1956. In LEFM theory, the driving force for a crack to extend is not the stress or strain but the stress intensity factor, known as K. The stress intensity factor uniquely describes the crack tip stress field independent of global geometry by embodying both the stress and the crack size. The relationship of the crack growth in the sense that the rate of crack growth, da/dN, with respect to the cyclic range of the stress intensity factor, AK, was derived by Paul C. Paris(Paris et al. 1961) in 1961, known as Paris Law. In reality, mechanical component are seldom subjected to purely constant amplitude loading history. The irregular stress history must be counted as a series of constant amplitude stresses. In addition, it is difficult to define a cycle in an irregular stress history. Since the reverse of stress curve can be easily found according to the sign change of the stress history curve, cycle counting techniques such as rain-flow counting method(Matsuishi and Endo 1968) are developed to combine reversals to form cycles. After that, cumulative damages can be calculated by Miner's law(Miner 1945). For most realistic structures or components, stress or strain fields are multi-axial. Fatigue life prediction methods for multi-axial loading also have been developed(Bannatine et al. 1990; Fuchs and Stephens 1980; Miller et al. 1966). In addition to some traditional method such as maximum principle stress/strain method, maximum shear stress/strain method and Von Mises' effective stress/strain method; Miller and Brown formulized the critical plane approach(Brown and Miller 1973) from the observation that the stress and strain normal to the plane with maximum shear has been recognized to strongly influence the development of fatigue crack. No consensus has been reached on the methods of multi-axial fatigue life prediction. All these methods have their own advantages in the specific application. So far, fatigue life analysis has been separated into two categories, (a) crack initiation, including S-N and e-N method, and (b) crack propagation. The criteria for the fatigue life of a component in engineering design depend on material properties or work conditions. In general, the automotive industry usually applies crack initiation criteria because of the nature of the product and use. On the other hand, the aircraft industry mainly uses crack propagation criteria by periodic inspection and fatigue crack monitoring to achieve and maintain structural safety. CHAPTER 3 UNCERTAINTY ANALYSIS USING STOCHASTIC RESPONSE SURFACE Introduction Uncertainty modeling and reliability analysis are the key issues in the reliability based design process. Uncertainty modeling can be decomposed into three fundamental steps: i) uncertainty characterization of model inputs, ii) propagation of uncertainty, and iii) uncertainty management/decision making. The uncertainty in model inputs can be represented in terms of standardized normal random variables (srv) with mean zero and variance equal to one. The selection is supported by the fact that they are widely used and well-behaved. For other types of random variables, an appropriate transformation must be employed. It is assumed that the model inputs are independent so each one is expressed directly as a function of a srv through a proper transformation. Devroye(Devroye 1986) presents the required transformation techniques and approximations for a variety of probability distributions. More arbitrary probability distributions can be approximated using algebraic manipulations or by series expansions. For uncertainty propagation, Monte Carlo Simulation (MCS) may be the most common choice because of the accuracy and robustness, but the dilemma of MCS is that the required large number of samples that makes it impractical for computationally demanding models. There are several remedies to reduce the number of samples in MCS, such as importance sampling(Melchers 2001) and separable MCS(Smarslok and Haftka 2006). However, they require special knowledge of the problem or special form of the response. Several computationally efficient methods were proposed in last two decades with reasonable accuracy in many structural problems, such as first- and second- reliability method (FORM/SORM), and response surface method (RSM). The stochastic response surface (SRS) can be viewed as an extension of classical deterministic response surfaces for model outputs constructed using uncertain inputs and performance data collected at heuristically selected collocation points. The polynomial expansion uses Hermite polynomial bases for the space of square-integrable probability density function (PDF) and provides a closed form solution of model outputs from a significant lower number of model simulations than those required by conventional methods such as modified Monte Carlo methods and Latin hypercube sampling. In this chapter, a surrogate-based uncertainty model using stochastic response surface (SRS) is introduced. Reliability analysis using Monte Carlo simulation on this surrogate model shows promising results in terms of accuracy and efficiency. The proposed method is compared with the first-order reliability method (FORM) and MCS. Description of Uncertainty Model When the inputs of a system are uncertain or described as random variables/parameters, the output or response from this system will have a stochastic behavior as well. Let us assume that these random inputs are given in an n-dimensional vector X with continuous joint distribution function fx(x). As shown in Figure 3-1, the system state can have a Boolean description such that the system fails when the limit state G(X) < 0 The probability of failure Pfcan then be defined as a cumulative distribution function (CDF) over the failure region, as P,= fx(x)dx. (3.1) G(X)<0 X2 G(X)<0 Failure region Limit state G(x)=0 G(X)>0 Safe region x1 Figure 3-1. Limit state function divides the safe region from the failure region Equation (3.1) is called the reliability integral. Since the integral domain defined by limit state function G(X) is complex in the multi-dimensional random space, the reliability integral is difficult to calculate. As introduced in the previous section, by transforming random variables from the original random space to the standard normal space, the limit state function can be expressed as a function of a set of srvs {uf} Then, Pfcan be expressed in standard Gaussian space as Pf = p, (u)du (3.2) G(U)<0 where p(*) is the standard normal PDF and U is the vector of standard random variables. The transformation between X and U is denoted as U=T(X). In FORM, the probability level of a system is usually represented by the reliability index or safety index /Y. For instance, if ((*) is the CDF of the standard random variable, the failure probability can often be represented by the reliability index 8/ = -( '(P ). Stochastic Response Surface Method (SRSM) Polynomial Chaos Expansion (PCE) in Gaussian Space Orthogonal polynomials have many useful properties in the solution of mathematical and physical problems. Just as Fourier series provide a convenient method of expanding a periodic function in a series of linearly independent terms, orthogonal polynomials provide a natural way to solve, expand, and interpret solutions to many types of important differential equations. Orthogonal polynomials associated with the generalized polynomial chaos (Askey- Chaos) are different according to different weight functions. The type of polynomials is decided by the match between the specific weight function and the standard probability density function (PDF). The corresponding type of polynomials and their associated random variables are listed in Table 3-1. Table 3-1. The type of polynomials and corresponding random variables for different Askey-Chaos (N>0 denotes a finite integer) Random variable Orthogonal Support range polynomials Gaussian Hermite (-00,00) Continuous Gamma Laguerre [0, co) Beta Jacobi [a,b] Uniform Legendre [a,b] Poisson Charlier {0,1,2,... } Discrete Binomial Krawtchouk {0,1,2,...,N} Negative Binomial Meixner {0,1,2,... } Hypergeometric Hahn {0,1,2,...,N} For example, in Table 3-1, Hermite polynomial chaos expansion requires the weight functions to be Gaussian probability density function, and it satisfies the following orthogonal relation: fJk fk(xk)Fk(xk)F (xk)dk =C, Vi, j (3.3) where fk(xk) is Gaussian PDF for variable xk, Fk(xk) is the Hermite polynomial basis, and upper indices i,j denote for two different polynomials. In this research, the uncertainty propagation is based on stochastic response surfaces (polynomial chaos expansion). The SRS(Isukapalli et al. 1998; Webster et al. 1996) can be view as an extension of classical deterministic response surfaces(Khuri and Cornell 1996; Myers and Montgomery 1995) for model outputs constructed using uncertain inputs and performance data collected at heuristically selected collocation points. The polynomial expansion uses Hermite polynomial bases for the space of square- integrable probability density function (PDF) and provides a closed form solution of model outputs from a significant lower number of model simulations than those required by conventional methods such as modified Monte Carlo methods and Latin hypercube sampling. Let n be the number of random variables and p the order of polynomial. The model output can then be expressed in terms of the srv {u,} as: G= al+ --afI,(u,)+ 1 a2 +I, I)- I a 3(, ,k)+- (3.4) 1=1 1 j= 1 i-1 ]=1 k=1 where GP is the model output, the ap, as,... are deterministic coefficients to be estimated, and the Fp(ui,...,Up) are multidimensional Hermite polynomials of degree: P (ut_) ( 1=)Pe 12UU _9P e -1/2UU (3.5) ( ..Ou t (-- ) u O& 9 ...Ou/p where u is a vector ofp independent and identically distributed normal random variables selected among the n random variables that represent the model input uncertainties. Equation (3.4) is also called polynomial chaos expansion. The Hermite polynomials Frp(u,..., up) are set of orthogonal polynomials with weighting function e /2, which has the same form with the PDF of standard random variables. In this research, a modified version ofHermite polynomial(Isukapalli et al. 1998) is used. The first four terms are u, u2 1, u3 3u, and u4 6u2 + 3, when a single random variable is involved. The use of the Hermite polynomials has two purposes: (1) they are used to determine the sampling points, and (2) they are used as bases for polynomial approximation. In general, the approximation accuracy increases with the order of the polynomial, which should be selected reflecting accuracy needs and computational constraints. The expressions for the 2nd- and 3rd-order polynomials in n dimensions (later used in the numerical examples) are: 2nd-order: n n n-1 n G2(u) (2) + 2) + (2) (u2 )+) 2) (3.6) 1=1 =1 1=1 j>1 3rd-order: G(3(u) =a(3) + 3)u + a3 (u 21) =1 1=1 n n-1 n + a) (u3 -3u)+ ( 3) (3.7) =1 z=1 y>1 n n n-2 n-1 n 1=1 j=1,ji iZ=1 j>i k>j The number of unknown coefficients is determined by dimension of the design space n. For 2nd and 3rd order expansion, if the number of unknowns is denoted by 2), M3), respectively: N(2) = 1+ 2n + n(n 1) (3.8) 2 N(3) = 1+3n+ 3n(n 1) n(n 1)(n 2) (3.9) 2 6 For n = 2, 4, and 8, for example, M2) = 6, 15, and 45; and N3) = 10, 35, and 165, respectively. The coefficients in the polynomial chaos expansion are calculated using the least square method, considering samples of input/output pairs. Since all inputs are represented using srv, more accurate estimates for the coefficients can be expected, in the sense of statistics, if the probability distribution of the u,'s is considered. The idea of Gaussian Quadrature of numerical integral can be borrowed to generate collocation points(Webster et al. 1996). In Gaussian Quadrature, the function arguments are given by the roots of the next higher polynomial. Similarly, the roots of the next higher order polynomial are used as the points at which the approximation being solved, which is proposed as the orthogonal collocation method by Villadsen and Michelsen (Villadsen and Michelsen 1978). For example, to solve for a three dimensional second order polynomial chaos expansion, the roots of the third order Hermite polynomial, -/r 0 and ,,r3 are used, thus the possible collocation points are (0,0,0),( -3 -x/3 ),(-/3 ,0, f3 ),etc.. There are 27 possible collocation points in this case. However, in equation(3.9), there are only 10 unknown coefficients. Similarly, for higher dimensional systems and higher order approximations, the number of available collocation points is always greater than the number of unknowns, which introduces a problem of selecting the appropriate collocation points. For a good approximation in polynomial chaos expansion, the choice of collocation points is critical. Hence, a set of points near the high probability region is heuristically selected among the roots of the one-order higher polynomial under restrictions of symmetry and closeness to the mean. Since the origin always corresponds to the highest probability in standard Gaussian space, the exclusion of the origin as a collocation point could potentially lead to a poor estimation. Thus, when the roots of high-order polynomial do not include zero, it is added in addition to the standard orthogonal collocation method. The Hermite polynomials orthogonall with respect to the Gaussian PDF) provide several attractive features, namely, more robust estimates of the coefficients with respect to those obtained using non-orthogonal polynomials(Gautschi 1996); it converges to any process with finite second order moments(Cameron and Martin 1947); and the convergence is optimal (exponential) for Gaussian processes(Xiu et al. 2002). In addition, the selected SRS approach includes a sampling scheme collocationn method) designed to provide a good approximation of the model output (inspired by the Gaussian Quadrature approach) in the higher probability region with limited observations. Once the coefficients are calculated, statistical properties of the prediction, such as mean and variance can be analytically obtained, and sensitivity analyses can be readily conducted. Numerical Example of Stochastic Response Surface As an illustration of the efficiency and convergence properties of the SRS approach, consider the construction of the PDF associated with a simple analytical function represented by: y(x)= ex (3.10) with x being a normally distributed random variable, as N(0,0.42). Note that in this case the analytical expression of the PDF is known. The SRS for 2nd- and 3rd-order polynomials are shown in Eqs.(3.11) and (3.12), respectively. y2) = 1.0833+0.4328u + 0.0833(u2 -1) (3.11) (3 = 1.0843 + 0.4333u + 0.0863(u2 -1) + 0.0112(u3 3u) (3.12) f0.8 00.6 0.4 0.2 0 1 2 3 4 5 y=ex x-N(0,0.42) Figure 3-2. PDF of performance function y(x) = ex In this particular example, the accuracy of the proposed SRS is compared with the analytical solution. Figure 3-2 shows the PDF obtained from MCS applied to the SRS and from the exact solution. A good agreement is observed in the PDF approximation, and the root mean square errors decreases with higher order polynomials, showing the convergence of the proposed SRS (Table 3-2). Table 3-2. Root mean square error of PDF compared with the exact PDF of performance function y=ex Polynomial order Errors 2 0.03835 3 0.00969 To illustrate the effectiveness of the SRS in the application to the structural problem, consider a torque-arm model in Figure 3-3 that is often used in shape optimization(Kim et al. 2003). The locations of boundary curves have uncertainties due to manufacturing tolerances, modeled as probabilistic distributions. Thus, the relative locations of corner points of the boundary curves are defined as random variables with x,~N(d,, 0.12). The mean values d, of these random variables are chosen as design parameters, while the standard deviation remains constant during the design process. X6 5 A X87 2 -2789N C ------------5066N Symmetric Design Figure 2-3: Shape design parameters for the torque Figure 3-3. Shape design parameters for the torque arm The torque arm model consists of eight random variables. In order to show how the SRS is constructed and the PDF of the model output is calculated, we choose the three random parameters (x2, x6, and x8) that contribute most significantly to the stress performance at points A and B in Figure 3-3. In the deterministic analysis with mean value, the maximum stress of A = 319MPa occurs at location A. The stress limit is established to be Gmax = 800MPa. In the reliability analysis the performance function is defined such that G < 0 is considered a failure. Thus, the performance function can be defined as G(x) = omax- A(X). The number of unknown coefficients is a function of the dimension n of the random variables. For 2nd- and 3rd-order expansion, the numbers of coefficients, denoted by N2 and N3, are 10 and 20, respectively. There are 27 possible collocation points and 10 unknown coefficients in the case of 2nd-order expansion. For robust estimation, the number of collocation points in general should be at least twice the number of coefficients. In this particular example, all possible collocation points are used. After coefficients are obtained, MCS with 100,000 samples is used to obtain the PDF. Figure 3-4 shows the PDF associated with G(x) when different orders of polynomial approximations are used. The PDF obtained from the direct MCS with 100,000 sample points is also plotted. It is clear that the PDF from the 3rd-order is much closer than that of the 2nd-order to the PDF from the MCS. 0.035 1 SRS p=2(27pts) 03- SRS p=3(125pts) 03----- MCS 0 025 0.02 I 0 0015 0.01 400 420 440 460 480 500 520 540 560 580 600 Response(Mpa) Figure 3-4. PDF of performance function G(x) torque arm model In order to compare the accuracy of the probability estimation through proposed SRS, let us check the probability of response being larger than 520MPa. In Table 3-3, the probability obtained from MCS is regarded as the reference. The relative error (e) of failure probability from MCS estimation with sample size of N can be calculated using the following equation: s = k -- p(3.13) where k denotes the confidence level, for confidence level of 95%, k=1.96, which can be verified from standard normal table. Thus, in Table 3-3, number of MCS sample is 100,000, the error in Pfwill be less than 5% with 95% confidence. As shown in Table 3-3, it is clear that the SRS provides a convergent probability result as the order increases. With third order SRS, the accuracy of reliability analysis is significantly improved, compared to FORM. Table 3-3. Comparison of probability of G>520MPa obtained from different uncertainty analysis methods (Full sampling without using local sensitivity) Method FORM 2nd order 3rd order SRS MCS SRS Prob. of G>520MPa 1.875% 2.061% 1.682% 1.566% Relative error* 19.732% 31.609% 7.407% - *Relative error: prob(approx.) prob(MCS) x 100% prob (MCS) Improving Efficiency of SRS Using Local Sensitivity Information In the proposed SRS, the number of sampling points depends on the number of unknown coefficients. Although the proposed method is accurate and robust, we have to address the curse of dimensionality: as the number of random variables increases, the number of coefficients rapidly increases, as can be seen in Eqs. (3.8) and (3.9). In addition to the efficient collocation method, the number of simulations can be reduced even further when local sensitivity is available. Recently, Isukapalli et al.(Isukapalli et al. 2000) used an automatic differentiation program to calculate the local sensitivity of the model output with respect to random variables and used them to construct the SRS. Their results showed that local sensitivity can significantly reduce the number of sampling points as more information is available. The computational cost of the automatic differentiation, however, is often higher than that of direct analysis(Van Keulen et al. 2004). However, in the application to the structural analysis, local sensitivity can be obtained at a reasonable computational cost. At each sampling point, the local sensitivity is a partial derivative of the limit state with respect to random variables. Hence, if local sensitivity information is available, then n+1 data at each sampling point can be used for constructing the proposed SRS, which significantly reduces the required number of sampling points. Continuum-Based Design Sensitivity Analysis In this research, the continuum-based design sensitivity analysis(Choi and Kim 2004a) is utilized to calculate the gradient of the performance function with respect to random variables. Even if the idea can be used in a broader context, only structural problems are considered in this research. Let z be the displacement and ibe the displacement variation that belongs to the space Z of kinematically admissible displacements. For given body force f and surface traction force t, the variational equation in the continuum domain Q is formulated as a(z, ) = I(T), (3.14) for all i e Z. In Eq. (3.14), the structural bilinear and load linear forms are defined, respectively, as a(z,' Y) = ,, j(z)E, (z)dQ (3.15) l() = ffdQ + t, dF (3.16) where e,, are components of the engineering strain tensor, and a,, are components of the stress tensor. In linear elastic materials, the constitutive relation can be given as 0 (z) = C"k -, (Z) (3.17) where the constitutive tensor cjki is constant. The summation rule is used for the repeated indices. In order to solve Eq.(3.14) numerically, the finite-element-based method or the meshfree method can be employed, which ends up solving the following form of matrix equation: [K]{D} = {F} (3.18) where [K] is the stiffness matrix, {F} the discrete force vector, and {D} the vector of nodal displacements. The major computational cost in solving Eq.(3.18) is related to L-U factorization of the coefficient matrix. As will be shown later, the efficiency of sensitivity calculation comes from the fact that sensitivity analysis uses the same coefficient matrix that is already factorized when Eq.(3.18) is solved. In design sensitivity analysis, the variational Eq.(3.14) is differentiated with respect to design variables. In shape design, the design variable does not appear explicitly in the governing equation. Rather, the shape of the domain that a structural component occupies is treated as a design variable. Thus, a formal procedure is required to obtain the design sensitivity expression. As shown in Figure 3-5, suppose that the initial structural domain Q is changed into the perturbed domain QO in which the parameter T controls the shape perturbation amount. By defining the design changing direction to be V(x), the material point at the perturbed design can be denoted as x, = x + TV(x). The solution z,(x,) of the structural problem is assumed a differentiable function with respect to shape design. The sensitivity ofz,(x,) at X, is defined as S= lim z(x + rV(x))- z(x) (3.19) TryO z- Initial domain Perturbed domain /r V(x) rT Figure 3-5. Variation of a structural domain according to the design velocity field V(x) The design sensitivity equation is obtained by taking the material derivative of the variational equation(3.14) The derivative of the structural energy form then becomes d d- aQ (z ,T) 0 = a0(z, ) +aa(z,i) (3.20) The first term on the right-hand side represents an implicit dependence on the design through the state variable, while the second term, the structural fictitious load, denotes an explicit dependence on the design velocity V(x), defined as a (z, z) = [J<(1)g (z) + (z)cy (z) + ()7,(z)divV]dQ (3.21) where S- Z) OZVk + jV(3.22) If the applied load is independent of displacement, i.e., conservative, then +(I)V ffQ [; 0 -L, V +zf/ ]dQ S J[( + Vtx ] dF (3.23) is the external fictitious load form, where Vn is the normal component of the design velocity on the boundary, and K is the curvature of the boundary. The design sensitivity equation is obtained from Eq. (3.20) to (3.23) as a(z,z) = ()- a(z,7) (3.24) for all e Z. Note that by substituting z into z, the left-hand side of the design sensitivity equation (3.24) takes the same form as that of the response analysis in Eq.(3.14). Thus, the same stiffness matrix [K] can be used for sensitivity analysis and response analysis, with a different right-hand side. Once the sensitivity z of the field vector is calculated, the sensitivity of the performance function with respect to design variable u, can be calculated using the chain rule of differentiation, as dy(z;x) _y(z;x) y(z;x) (3 duA au z When finite element analysis is used, the sensitivity equation can be solved inexpensively because the coefficient matrix is already factorized when solving Eq.(3.14) and the sensitivity equation uses the same coefficient matrix. The computational cost of sensitivity analysis is usually less than 20% of the original analysis cost. The computational efficiency of the uncertainty propagation approach is critical to RBDO since as previously stated at each design cycle the updated version of the PDF for the constraint function (related to model outputs) is required. Constructing SRS Using Local Sensitivity In SRS, the system response can be approximated as polynomial expansion when k sampling data are available, the linear regression equation can be written as Y1 1 ui u 2 ** ao 2r a, {y} = [A]{a} (3.26) yk U 3 23 The above equation is the standard form for linear regression to solve for unknown coefficients {a}. When the sensitivity information is available, additional information at each sampling point can be used in calculating the coefficients. By differentiating Eq.(3.4) with respect to random variable ui and by applying the same regression process in Eq.(3.26), we have {dy = A {a} (3.27) du, au, Equations (3.26) and (3.27) can be combined to obtain the following regression equations: y dy dul dy du ,T{ Let {Y}= y d A OA u. {a} (3.28) OA AA aA = A A T, Eq. (3.28) can be 9u, ou" written as {Y} = [B]{a} (3.29) Thus, the coefficients of SRS can be obtained using least square regression, such {a} = ([B]' [B]) [Bl] {Y} (3.30) Note that the sensitivity can be calculated using the transformation of au, random variables, as -Y (3.31) As introduced in the previous section, the local sensitivity c9y / 9x can be obtained implicitly through Eq. (3.25), where design variable is represented by u, instead of x, since notation x has been used as space coordinate. Since the transformation between srv and variables with other types of distribution are also mathematically well developed, cx, /au, can be obtained explicitly. Therefore, Eq.(3.30) provides an explicit solution for obtaining coefficients of SRS. Numerical Example Torque Arm Model In order to show the effectiveness of the proposed SRS with local sensitivity, the same torque arm problem with previous example is used. All conditions are the same as before. By taking advantage of using sensitivity information to build stochastic response surface, the number of collocation points is reduced significantly. Here for the second- order polynomial chaos expansion, 7 points have been used, while 31 points for the third- order case. At each sampling point, the function value and sensitivity information are used to construct the SRS. The PDF obtained from the SRS with sensitivity is plotted in Figure 3-6 along with that from MCS with 100,000 samples. In the case of 2nd-order, the SRS with sensitivity is more accurate than the SRS without sensitivity (Figures 3-4 & 3-6). In order to calculate the accuracy, the probability of G > 520MPa is calculated using FORM, second- and third-order SRS (Table 3-4). Since no analytical solution is available, MCS with 100,000 samples is used as a reference. Both SRS are more accurate than FORM. 0.035 SRS-SEN p=2(7 pts) -- SRS-SEN p=3(31 pts) 03 ----- MCS(100,000 pts) 0.025 S0.02 - & 0015 \ S001 0.005 \ 400 420 440 460 480 500 520 540 560 580 600 Response (MPa) Figure 3-6. PDF of performance function G(x) Torque model at initial design (SRS with sensitivity) Table 3-4. Comparison of probability of G>520MPa obtained from different uncertainty analysis methods (reduced sampling using local sensitivity); Method FORM 2nd order SRS 3rd order SRS MCS Prob. of G>520MPa 1.875% 1.520% 1.545% 1.566% Relative error* 19.732% 2.937% 1.341% - *Relative error: prob(approx.)- prob(MCS) x 100% prob (MCS) Table 3-5 compares the probability of G>520 MPa of second order SRS with/without using local sensitivity with that of MCS, which is regarded as the reference of exact value. With local sensitivity and seven sampling points, SRS provides more accurate probabilistic result than that without utilizing local sensitivity and twenty-seven sampling points. The accuracy is improved by using local sensitivity while computational cost is reduced. Table 3-5. Comparison of probability of G>520MPa obtained with/without local sensitivity (7/27 sampling points) using 2nd order SRS 2nd order SRS 2rd order SRS using 27 sampling using 7 sampling MCS Method points without points with (100,000 samples) sensitivity sensitivity Prob. of G>520MPa 0.2061% 1.520% 1.566% Relative error* 31.6091% 2.937% - *Relative error: prob(approx.)- prob(MCS) x 100% prob(MCS) Summary In this chapter, a stochastic response surface method (SRSM) using polynomial chaos expansion is used in calculating structural reliability. Compared with FORM, which is based on the linear approximation at the most probability point, it provides more accurate result for nonlinear responses. In addition, orthogonal polynomial basis provide 38 a convergent behavior as the order of polynomial is increased. A nonlinear function has been used as numerical example to show its accuracy and convergence. Since continuum based sensitivity results were obtained during structure analysis, the computational cost is further reduced by providing gradient information while fitting response surface. SRSM has been applied on a structural problem to show its effectiveness. When sensitivity information is provided, numerical results show that even lower number of sampling point can provide more accurate result. CHAPTER 4 RELIABILITY-BASED DESIGN OPTIMIZATION Although statistical methods of uncertainties quantification have been studied intensively for decades, traditional deterministic design optimization still takes no advantage in these scientific advances and compensates uncertainties based on experience; for example, the safety factor. Such an optimization scheme usually yields either unsafe or too conservative design due to the lack of uncertainty quantification. In order to impose existing knowledge of uncertainty to engineering design process, reliability-based design optimization (RBDO) methodologies have been proposed and developed(Chandu and Grandi 1995; Chen et al. 1997; Du and Chen 2002b; Enevoldsen and Sorensen 1994; Grandhi and L.P. 1998; Kim et al. 2004b; Kwak and Lee 1987; Liang et al. 2004; Tu 1999; Tu and Choi 1997; Youn et al. 2003), where the system reliability or probability of failure is used to evaluate the system performance. Compared to the deterministic optimization, RBDO provides margins on reliability by quantifying the uncertainty in the response of structural system due to input uncertainty. General RBDO Model Design optimization has been introduced to structural engineering for decades(Arora 2004; Haftka and Gurdal 1991; Vanderplaats 2001). Its methodologies have been well developed mathematically, and applications in product development are flourishing. The underlying design philosophy is to reduce the cost by pushing the design to its performance margin. In traditional deterministic design, an optimization problem is formulated as minimize Cost(d) subject to G, (d) < Gaowable, j = 1, 2,..., np (4.1) dL _ corresponding maximum constraint allowable; and d denotes the vector of the deterministic design variables. The objective is to minimize the cost while meeting the system constraints. A system design depends on the system design variables. In deterministic optimization, both objective and constraints only depend on the design vector d which contains all deterministic design variables d,. In reliability-based design, design is based on a randomly distributed system vector, e.g., denoted by X, which contains random variable Xi. In RBDO, the mean value ti or the standard deviation ci of the system variable Xi can be used as the design variable. In some cases, uncontrollable random variables may contribute to the uncertainty of the performance. Instead of directly setting constraints on deterministic performance, the RBDO problem(Chandu and Grandi 1995; Enevoldsen and Sorensen 1994; Grandhi and L.P. 1998; Wu and Wang 1996) can generally be defined by setting constraints to be uncertainty measures, such as probability of failure. It is then formulated as minimize Cost(d) subject to P(GJ(x) f, is the prescribed failure probability limit for thejth constraint. Reliability Index Approach (RIA) and Performance Measure Approach (PMA) In the RBDO formulation described in the previous section, each prescribed failure probability limit Pf is often represented by the reliability target index as, 8 = >(73) . Hence, any probabilistic constraint in Eq. (4.2) can be rewritten using equation as FG (0) _< (-/6) (4.3) where FG(O)=P(G<0) is the cumulative distribution function(CDF) of G at the failed region. Equation(4.3) can also be expressed in another way through inverse transformations S=-- (FG (0)) > f (4.4) where fls is traditionally called the reliability index. The expression of probability constraint in Eq. (4.4) leads to the so called reliability index approach (RIA)(Tu and Choi 1997; Tu 2001; Youn 2001). The two forms of constraint described in equations (4.2)and (4.4) are basically the same. In FORM/SORM based RBDO, an inner loop optimization is used to find the most probability point (MPP) in the standard Gaussian space. RIA may cause singularity because /, approaches infinity or negative infinity when the failure probability is zero or one. In that case, inner loop optimizer may fail to find the MPP. There is an alternative way to avoid singularity(Tu and Choi 1997) based on a different concept of reliability measure. For any given target probability, a certain level of performance can be reached to meet the reliability requirement. Tu(Tu and Choi 1997) proposed an inverse measure approach called performance measure approach (PMA) based on FORM by transforming Eq.(4.4) to g* = FGI ((-fl)) > 0 (4.5) where g* is named the target probabilistic performance measure. In PMA, Eq.(4.5) is used as probabilistic constraint of RBDO. PMA has been proved to be consistent with RIA in prescribing the probabilistic constraint, but their differences in probabilistic constraint evaluation can be significant (Tu 1999). PMA is more robust in FORM/SORM than RIA based on the fact that RIA may yield singularity; that is, 8f, approaches infinity or negative infinity. In addition, for an inactive probabilistic constraint, PMA is more efficient than RIA. Known as an inverse measure approach, PMA can also be implemented on the sampling based uncertainty estimation method. For example, in MCS, a performance measure that meets reliability requirement can be obtained from the order statistics of sampled performance values. Figure 4-1 shows the general numerical procedure of RBDO. The effect and efficiency of inverse measure approach has been investigated for FORM/SORM(Tu 1999; Tu and Choi 1997). In this research, RIA and PMA as two different philosophies for probability constraint evaluation are also addressed for SRS-based RBDO. FORM/SORM: inner loop optimization for MPP search Reliability Analysis: RIA I PMA SRSM: DOE I Sensitivity Analysis of Cost and Constraints Update Design in Optimizer No Converge? Yes Figure 4-1. Flow chart for reliability-based design optimization Probability Sensitivity Analysis (PSA) Similar to the traditional design sensitivity, where sensitivity quantifies the effect of deterministic design variable to the structure response, probability sensitivity provides the quantitative estimation of the changing of failure probability or reliability with respect to the changes of random parameters, such as means or standard deviations of random design variables. In RBDO, the gradient based optimizer needs sensitivity information to carry out optimization. Automatic differentiation using finite differentiation leads to a significantly extra computational cost, especially when there are many design variables. In RBDO, if constraints are set with the probability of failure being less than a certain threshold, the gradient of probability with respect to the random input is required. In this research, probability sensitivity analysis is utilized to calculate the gradient information. First, the probability sensitivity calculation in FORM is introduced by taking advantage of structural sensitivity analysis. It can be shown that one can obtain accurate probability sensitivity without extra simulation cost. Since SRSM shows an advantage for nonlinear response, sampling based probability sensitivity is also introduced. For inverse measure approach, sensitivity for both FORM and sampling based RBDO can be obtained. Probability Sensitivity Analysis in FORM In first order reliability method (FORM), reliability index (3) can be obtained by following equation /7 (UTTU*)12 (4.6) where U* is the vector of MPP. The derivative of failure probability with respect to the design variables in FORM can then be written as 0r7 0)7 0,8 0r7 0t7 where (p(*) is the PDF of the standard random variable. Thus, the sensitivity of the failure probability is directly related to that of the reliability index, which can be obtained by =. (VuT -I2 -_ (4.8) 0o7 0Qa /7 0a For a random variable 7 = 0, /p 1 .T FT(X*,0) aT(X*,0) aX' -, -U8 + X 0 = T(x ) I ((4.9) 1 U.T OT(X',0) Since the reliability index and the most probable point are available from the reliability analysis, the sensitivity can be easily obtained. If the computationally expensive structure analysis code does not come with sensitivity analysis, a finite difference method is widely used to provide gradient information for searching the most probability point (MPP). The computational cost of finite difference method is proportional to the number of design variables. Using design sensitivity analysis, we can avoid the finite difference calculation and provide more accurate gradient information to the line search for MPP. In the finite difference method, the gradient of the limit state in the standard normal space is defined as Vg(U) = lim g(U + AU)-g(U) (4.10) AU-O AU Every iteration in line search needs to perturb each design variable to evaluate the gradient. If the sensitivity information can be obtained from a structural analysis code, there is a more efficient way to obtain the gradient information for MPP search. The gradient Vg(U) can be computed as Vg(U) = Vg(X) T (U) (4.11) where T: X->U. The transformation T from original random design space to the standard Gaussian space can usually be obtained explicitly, and the gradient Vg(X*) is provided by design sensitivity analysis. In this section, the torque arm model described in Chapter 3 is used to evaluate the accuracy of the probability sensitivity analysis using FORM. At the initial design, the probabilistic parameters of eight random variables are considered as design variables. Each random variable is assumed to be normally distributed with a mean of zero and a standard deviation 0.1. The sensitivity of reliability index is calculated based on Eq.(4.9). Since the transformation T is an explicit function of probabilistic parameters, the sensitivity can easily be calculated with reliability analysis. Table 4-1 shows the sensitivity results with respect to mean (0f/ /p, ) and standard deviation (6,8/ /cr ). The accuracy of the sensitivity is compared with that of the finite difference method with 1% perturbation size. Table 4-1. Probability sensitivity with respect to random parameters (unit: centimeter) design af A,8 A,/8 Ap, x 100% A8 Af Au, x 100% A9u, A/, C9,8 C9, Au, C9 , xi 0.376 0.377 100.26 -0.030 -0.030 100.00 x2 5.243 5.243 100.00 -5.773 -5.775 100.03 x3 0.034 0.034 100.00 -0.000 -0.000 100.00 x4 0.106 0.106 100.00 -0.002 -0.002 100.00 x5 0.055 0.055 100.00 -0.001 -0.001 100.00 x6 -7.244 -7.244 100.00 -11.022 -11.011 99.90 x7 -0.140 -0.140 100.00 -0.004 -0.004 100.00 x8 -4.457 -4.457 100.00 -4.171 -4.173 100.05 In Table 4-1, the first column represents eight random variables that have normal distributions. All random variables are assumed to be independent. Since the mean value and the standard deviation are considered as probabilistic parameters, there are 16 cases in the sensitivity calculations. The second and fifth columns represent the sensitivity results obtained from the analytical derivative, while the third and sixth columns are sensitivity results from the finite difference method. A very good agreement between the two methods is observed. Table 4-2 shows the computational efficiency of the proposed analytical sensitivity calculation. The gradient information is provided from design sensitivity analysis in MPP search in the standard HL-RF method(Hasofer and Lind 1974; Liu and Kiereghian 1991). The computational savings are about 90% compared to the case when only the function values are provided. Once the reliability analysis is finished, the sensitivity of reliability index requires additional 17 function evaluations for the finite difference method, while only a single analysis is enough for the proposed method because the analytical expression in Eq.(4.9) and (4.11) is used. Table 4-2. Computational efficiency of analytical method for probability sensitivity Finite differential method Analytical method number of analyses in MPP search 90 10 number of analyses in sensitivity calculation 17 1 Total number of analysis 107 11 Probability Sensitivity Analysis Using SRSM In RBDO, the probability of failure can be formulated as P, = f (x)dx (4.12) G(X)<0 where G(X) < 0 is the failure region and f(*) is the joint probability density function(PDF). By introducing an indication function I(G(X) < 0) such that I=1 if G(X) < 0 and I=0 otherwise, Eq.(4.12) can be rewritten as Pf = f/(G(x) < 0)f (x)dx (4.13) where QX denotes the entire random design space. Since Eq.(4.12) is used as a constraint in RBDO, the sensitivity of Pfis required. The derivative of failure probability can be written as P= J I(G(x) < 0) x = J I(G(x) < 0) f(X) f (x)dx x 0 x (4.14) f I(G(u) <0). cf(x) p(u)du f(x)aO x=T 1(u) where Qu denotes the entire standard normal space. Explicit expression of Eq.(4.14) for different distribution types and numerical examples are derived in Appendix A. The accuracy of the sensitivity results are also presented in Appendix A for the case of various distribution types. Reliability-Based Design Optimization Using SRSM Although RIA and PMA are theoretically consistent in prescribing the probability constraint, there are still significant differences in probabilistic constraint evaluation. The RBDO based on RIA and PMA may lead to either different convergence or efficiency. In this section, an RBDO problem is formulated for the same torque-arm model in Chapter 3 using the concepts of RIA and PMA. A 3rd-order SRS is constructed for uncertainty analysis for both RIA and PMA. When the reliability index is used as a constraint in RBDO, it sometimes experiences numerical difficulty because it can have a value of infinity for very safe design. When SRSM is used in evaluating the probabilistic constraint in RBDO, the problem of singularity can be avoided naturally since the value of failure probability can always be obtained from MCS. The accuracy and convergence of SRSM have been illustrated in the previous chapter. Although SRSM usually requires more performance evaluation compared to FORM, it is still an affordable and applicable approach to obtain more accurate results for the highly nonlinear system. RBDO with RIA In this section, the RBDO problem of the torque arm model is solved using RIA. Stochastic response surface is used in uncertainty analysis to evaluate probability constraints. RBDO formulation of Eq. (4.2) can be used straightforwardly to solve the problem. For the torque-arm problem, the objective is to minimize the weight while meeting the requirement of reliability constraint. If we define that the structure fails when stresses in this structure reach yield stress, such that G,(x)= a,(x)-c <0 (4.15) where x is the random input variables, cr (x) is stress response for ith constraint, ar, denotes yield stress. The RBDO problem is then defined as Minimize Mass(d) subject to P(G,(x) where ft is the target reliability index and 1( ) is the cumulative density function of srv. During the optimization, a fit = 3 is used, which corresponds to 99.87% reliability. Since the maximum stress location can shift, the probabilities of failure at eight different locations are chosen as constraints in Eq. (4.16). The constraints can be evaluated using the SRS. The SRS needs to be reconstructed at each design cycle. Table 4-3 shows the properties of the random variables and the lower and upper bounds of their mean values (design variables). Note that the design variables are the relative change of the corner points and the initial values of all design variables are zero. The lower and upper bounds are chosen such that the topology of the boundary is maintained throughout the whole design process. Table 4-3. Definition of random design variables and their bounds. The values of design variables at optimum design are shown in the 5th column (unit: centimeter). Random dL d dU dopt Standard Distribution Variables (Initial) (optimum) Deviation type di -3.0 0.0 1.0 -2.5226 0.1 Normal d2 -0.5 0.0 1.0 -0.4583 0.1 Normal d3 -1.0 0.0 1.0 -0.9978 0.1 Normal d4 -2.7 0.0 1.0 -2.4663 0.1 Normal d5 -5.5 0.0 1.0 -2.1598 0.1 Normal d6 -0.5 0.0 2.0 1.9274 0.1 Normal d7 -1.0 0.0 7.0 2.3701 0.1 Normal d8 -0.5 0.0 1.0 -0.0619 0.1 Normal The design optimization problem is solved using the sequential linear programming technique. The optimization process converges after 14 design cycles and 27 performance evaluations where the relative convergence criterion has been met in two consecutive designs. The optimum values of the design variables are shown in the 5th column in Table 4-3. Figure 4-2 shows the optimum design and analysis results at the mean values. The major changes occur at design parameters d4, d/5, d and cd. Even if the maximum stress constraint is set to 800MPa, the active constraint at optimum design converges to a lower value so that the variance of the input parameters can be accounted for. In Figure 4-2 the maximum value shows 739MPa, which is the extrapolated nodal stress. The actual element stress at the active constraint is about 618MPa, which is much lower than the extrapolate stress show on the figure. Figure 4-3 provides the optimization history of the cost function. As a result of the optimization process, the mass of the structure is reduced from 0.878kg to 0.497kg (a reduction of about 43.4 %). 7 39+00 Nt$'rUO 211 700 211,UOO 1 bR+DOE 1 06+00 5 04-003 Figure 4-2. Optimum design and stress distribution of the torque arm model with 8 random variables. 09 0.85 - 08 - 0.75 - 07 - S0.65 06 - 0.55 - 05 - 0.45 - 04 0 2 4 6 B 10 12 14 iterations Figure 4-3. Optimization history of cost function (mass) for the torque arm model with 8 random variables. In order to observe the impact of the accuracy of the uncertainty propagation procedure at the optimum design, a 3rd-order SRS is considered at the optimum design (Figure 4-4). Table 4-4 shows the values of the reliability indices for the active constraint at the optimum design obtained from MCS withl00,000 samples, the proposed SRS approach, and the FORM. The MCS result is used as a reference mark to compare the other two methods, which has about 1.5% error in estimating reliability index with confidence level of 95%. The proposed SRS approach exhibits a lower error than the corresponding to the FORM and compares very well with the exact result (namely, 3) and that obtained using MCS (0.6% error) 10 iPU {P ot vespoPsP .t Optimum Figure 4-4. PDF of the performance function at the optimum for the torque-arm problem Table 4-4. Reliability Index of active constraint at optimal design Reliability Index Error (%) MCS 3.0307 - SRS 3.0115 0.633 FORM 2.9532 2.556 RBDO with Inverse Measure As we discussed in section 3.2, enlightened by PMA, an inverse measure approach can also be applied in SRS based RBDO. In this section, the inverse measure approach is applied on reliability based design optimization for the torque arm model. As discussed in section 3.2, the design problem can be formulated as Minimize Mass(d) subject to G,(u) > 0, i= 1,...,NC (4.17) dL < d < is N and allowed maximum probability of failure is Pf, then G, can be found by ordering samples and selecting pth smallest sample. p=NxPf (4.18) Thus, the evaluation of reliability constraints is transformed to find the pth order statistic of sampling. One of the advantages of this approach is that the sensitivity of performance based constraint measure can be obtained directly through structure sensitivity analysis: 8G(u) a8G(x*) 8G(x*) QT- (u*;d) (4.19) 00 Q0 ax, Q0 For example, if 0 = /,, and X ~ Normal(u, cr2) OG(u) _G(x*) _G(x*) QaT-(u*;d) (4.20) G(x*) Summary In this chapter, reliability based design optimization using stochastic response surface is discussed. Procedures for both RIA and PMA are investigated and formulated. A torque arm problem shown in Chapter 3 has been used to demonstrate the feasibility of RBDO using SRS. Since accurate, sensitivity calculation is important to the convergence of gradient based optimizer, probability sensitivity using FORM and SRS is presented. It is shown that probability sensitivity in SRS based sampling approach can be also obtained with minimal increase of computational cost. If the SRS is accurate enough, the accuracy of sensitivities obtained also have convergent with respect to the increasing of the sampling size. CHAPTER 5 GLOBAL SENSITIVITY ANALYSIS FOR EFFICIENT RBDO Introduction In industrial applications, a system usually involves considerable number of random variables. As stated in the previous two chapters, the increasing number of random variables boosts the computational cost of reliability analysis significantly. Structural reliability analysis which involves a computationally demanding model is limited by the relatively high number of required function analyses for uncertainty. Even if the local sensitivity information described in chapter 3 can reduce the number of required simulations, the dimension of the SRS will still increase according to the number of random variables. Design engineers want to reduce the number of variables based on their contribution to the output performance. However, it is challenging to identify the importance of a random variable in the process of uncertainty propagation. Those random variables with extremely low contribution to the performance variance can be filtered out to reduce the computational cost of uncertainty propagation. Recent development in statistics introduces global sensitivity analysis (GSA)(Saltelli et al. 2000; Saltelli et al. 1999; Sobol 1993; Sobol 2001), which studies how the variance in the output of a computational model can be apportioned, qualitatively and quantitatively, to different sources of variations. In this chapter, global variance- based sensitivity analysis has been applied on structural model to illustrate different roles of random variables in uncertainty propagation. Effort has been made to reduce the dimension of random space in the RBDO process. To reduce the number of simulations required to construct the SRS even further, unessential random variables are fixed during the construction of the SRS. A random variable is considered unessential (and hence it is fixed) if its contribution to the variance of the model output is below a given threshold. Global sensitivity indices considering only main factors are calculated to quantify the model input contributions to the output variability hence establishing which factors influence the model prediction the most so that: i) resources can be focused to reduce or account for uncertainty where it is most appropriate, or ii) unessential variables can be fixed without significantly affecting the output variability. The latter application is the one of interest in the context of this chapter. The RBDO problem in the previous chapter was solved with all random variables. However, some random variables did not significantly contribute to the variance of the stress function. Thus, a large amount of computational cost can be saved if the random variables whose contribution to the variance of the output is small are considered as deterministic variables at their mean values. This section describes how the global sensitivity indices considering only main factors can be used for deciding unessential random variables during the construction of stochastic response surfaces. Sensitivity Analysis As defined by Saltelli(Saltelli et al. 2000), sensitivity analysis studies the relationships between information flowing in and out of the model. In engineering design application, sensitivity usually refers to the derivative of performance measure with respect the input design variable. This derivative is used directly in deterministic design as sensitivity information is required by the gradient-based optimizer. In this research, this derivative is called local sensitivity, which has been applied in constructing SRS in the previous chapters. Local sensitivity analysis concentrates on the local impact of the design variables. It is carried out by computing partial derivatives of the output performance with respect to the input variable at the current design point. In structural optimization, substantial research has been done on the local sensitivity(Choi and Kim 2004a; Choi and Kim 2004b). Another choice of sensitivity measure explores what happens to the performance variance if all design variables are allowed a finite variation. Global sensitivity techniques apportion the output uncertainty to the uncertainty in the input factors. A couple of techniques have been developed in recent two decades. In this research, Sobol's sensitivity indices(Sobol 1993; Sobol 2001), which is based on the decomposition of function into summands of increasing dimension, will be discussed and applied to the constraint evaluation of RBDO. Variance-Based Global Sensitivity Analysis (GSA) Variance-based methods are the most rigorous and theoretically sound approaches (Chen et al. 2005; Saltelli et al. 2000; Saltelli et al. 1999; Sobol 1993; Sobol 2001)for global sensitivity calculations. This section describes the fundamentals of the variance- based approach and illustrates how the polynomial chaos expansions are particularly suited for this task. The variance based methods: (i) decompose the model output variance as the sum of partial variances, and then, (ii) establish the relative contribution of each random variable (global sensitivity indexes) to the model output variance. In order to accomplish step (i), the model output is decompose as a linear combination of functions of increasing dimensionality as described by the following expression: n n n '1 lj>1 subject to the restriction that the integral of the weighted product of any two different functions is zero. Formally, f .. f p(x)flx,..., ...,x,)fJl,"(_xJl ( ...,xJ)dx = 0, for il,...,i, j i,..., j (5.2) where p(x) is the joint PDF of input random variable x. If, for example, the weighting function is the uniform distribution for the random variables or the Gaussian probability distribution, the functions of interest can be shown to be Legendre and Hermite orthogonal polynomials, respectively. The model output variance can now be calculated using a well-known result in statistics. The result establishes that the variance of the linear combination of random variables (U,) can be expressed as: V bo0 + -bU1 = Eb2V(U,) + 2E-ECOV(U7,UJ) (5.3) Hence, the model output variance can be shown to be: V(f)= -aV(f) a -jaV(fj) +... a2...(2...) (5.4) where the terms represent partial variances and each V( ) may be found by definition as: V(f) [ f () E(f (x))]2p(x)dx (5.5) In the above formula, f(x) represents the function under consideration and the symbol E( ) denotes expected value. There are no covariance terms in Eq.(5.4) because of the orthogonality property shown in Eq.(5.2). The global sensitivity index S, that considers only main factor is called main sensitivity index, which associated with each of the random variables which is represented by Eq. (5.4): A sensitivity index that considers the interaction of two or more factors is called interaction sensitivity index. Thus, as denoted by Chan and Saltelli(Saltelli et al. 2000), the summation of all sensitivity indices, involving both main and interaction effect of i-th random variable, is called total sensitivity index. Sobol(Sobol 2001) suggested to use total sensitivity indices to fix unessential variables. If total sensitivity index for certainty variable extremely small compare to 1, that means the contribution of the variable is neglectable and the variable can be considered as a deterministic one. Global Sensitivity Analysis Using Polynomial Chaos Expansion The polynomial chaos expansion is particularly suited for computing global sensitivity indices because: (1) The model output is already decomposed as a sum of functions of increasing dimensionality; (2) the functions are orthogonal with respect to the Gaussian measure (Hermite polynomials); and (3) the variance of the bases are analytically available. For example, the variances of the functions associated with a two dimensional chaos expansion of order 2 are shown in Table 5-1. Table 5-1. Variances of the Hermite bases up to the second order Function f V(f) Xi 1 X2 1 2 Xi2 1 2 X1X2 1 2 X22 1 2 In addition, given the polynomial chaos expansion (i.e., the coefficients of the linear combination of Hermite polynomials), the model output variance and global sensitivity indices can be easily computed using Eqs. (5.4) and (5.6), respectively. In both equations, the variances V( -), are readily available for polynomial chaos expansions of arbitrary order and number of variables. In order to show the effect of the global sensitivity, let us consider the torque arm model presented in chapter 3. The eight random design variables are normally distributed with standard deviation of 0.1. A cubic stochastic response surface with eight variables in standard normal space has been constructed to approximate the stress response. Global sensitivity indices for main factor have been calculated to quantify the contribution of each random variable to performance variability. Figure 5-1 shows that three design variables(x2,X6,x8) makes most contribution to the total variance, influences from other variables are extremely small. That means, at the initial design, that the randomness of those variables with low global sensitivity indices can be ignorable. If only three variables are considered in constructing SRS instead of eight, the computational cost will be saved significantly, while maintaining the same level of variability. 1.595% 0.05% 0.001%. 0.021% 15.784% 15% EX1i 0.000% EX2 OX3 OX4 EX5 EX6 EX7 82.332% 1 X8 Figure 5-1. Global sensitivity indices for torque arm model at initial design Adaptive Reduction of Random Design Space Using GSA in RBDO The idea of adaptive reduction of random variables is based on the main factor of each random variable. If it is smaller than a threshold, it is fixed in constructing SRS. For that purpose, a linear polynomial chaos expansion is enough to obtain the main factors of GSI. The current algorithm with linear polynomial chaos expansion for the reduction of random variables can be modified to use a sensitivity index that accounts for interactions. These interactions will only appear in higher order polynomial chaos expansions. The choice of a non-linear polynomial chaos expansion would reduce the computational efficiency of the proposed approach with unclear significant advantages. As stated at the beginning of this section, once the global sensitivity indices are calculated, variables that have the least influence on the model prediction (unessential variables) can be identified and eventually held fixed without significantly affecting the output variability. The procedure is adaptive because the global sensitivity indices are calculated at each design iteration and as a result different sets of random variables may be fixed throughout the RBDO process. A flow chart of RBDO using this strategy is shown in Figure 5-2. Figure 5-2. Adaptive reduction of unessential random design variables using global sensitivity indices in RBDO. Low-order SRS is used for global sensitivity analysis, while a high-order SRS is used to evaluate the reliability of the system. At the initial design stage, a lower-order stochastic response surface is constructed using all random variables. In this particular example the first-order SRS is constructed using 17 sampling points. At the initial design, the first-order SRS with eight random variables can be expressed as, G1 a0 + alul + a22 +' + a44 + a5u5 a6u6 a7u7 + a8u8 = 4.95 + 0.0063u, + 0.117u2 + 0.Onnn( ,. 0.0019u4 (5.7) +0.01121 ., 0.052u6 0.0002u7 0.016us One useful aspect of the polynomial chaos expansion is that the coefficients in Eq. (5.7) are a measure of the contribution of the corresponding random variable to the variation of the output, and these coefficients will not change significantly in higher- order SRS. On the other hand, typically the global sensitivity index associated with a particular variable is responsible for most of its contribution to the output variance. Thus, evaluating the global sensitivity indices using the first-order SRS can be justified. Since all random variables are transformed into standard normal random variables, the variance of G1 can be evaluated analytically. Using Eq. (5.6) and (5.7) and assuming the design variables are independent, the global sensitivity index can be calculated as: a2 S= (5.8) j=1 Note that the denominator in Eq. (5.8) is the total variance of G1 using the first- order approximation. Thus, the global sensitivity index, S,, is the ratio of the contribution of i-th random variable to the total variance. If the global sensitivity index of a specific variable is less than a threshold value, the variable is considered as deterministic and fixed at its mean value. In order to show the advantage of fixing unessential random variables, the global sensitivity indices of the torque-arm model are calculated. Table 5-2 shows the global sensitivity indices of the torque-arm model using the first-order SRS at the initial design. The total variance of stress function is 1.670x 10-2. Based on the global sensitivity indices, there are only three random variables whose contribution is greater than 1.0%; i.e., u2, U6, and u8. Thus in the reliability analysis, only these three random variables are used in constructing the third-order SRS, which now requires only 19 sampling points for 10 unknown coefficients. All other random variables are considered as deterministic variables at their mean values. If the total number of sampling points for both lower- (17) and higher-order (19) polynomial expansions are compared with the higher-order SRS using all random variables (89), a significant reduction of the number of sampling points was achieved. Table 5-2. Global sensitivity indices considering only main factors for the torque arm model at the initial design. Only three random variables (u2, u6, and u8) are preserved when a threshold value of 1.0% is in place. SRV Variance GSI (%) ui 3.916 x 105 0.235 u2 1.369x 102 82.0 u3 6.403 x 10-9 0.00003834 U4 3.667 x 10-6 0.02197 us 6.864x10-6 0.04109 u6 2.702x10-3 16.179 u7 4.818x10-8 0.0002885 u8 2.538x 10-4 1.519 The RBDO problem, defined in Eq.(5.2) in chapter 4 is now solved using the proposed adaptive reduction of random variables. The optimization algorithm converges after the 17-th iteration. As shown in Figure 5-3, the optimum design using the adaptively reduced SRS is slightly different from that with all random variables in chapter 4. The former has a longer interior cutout than the latter. This can be explained from the fact that the model with reduced random variables has less variability than the full model. Furthermore, the optimum value achieved using the adaptively reduced SRS converges to a lower value than the one without adaptive reduction). The total mass of the torque arm is reduced by 54.8%. The difference between the two approaches is approximately 1.8%. The number of active random variables associated with the modeling of the first constraint during the design iterations are listed in Table 5-3. On average, four random variables were preserved, which implies that only 29 sampling points were required for constructing the SRS. This is three times less than the SRS approach without adaptive reduction (89 sampling points). Figure 5-3. Optimum designs for the full SRS (solid line) and adaptively reduced SRS (dotted line). Because some variables are fixed, the interior cutout of the design from the adaptively reduced SRS is larger than that from the full SRS. Table 5-3. Comparison of the number of random variables in each design cycle. The threshold of 1.0% is used. The first constraint is listed. Iter Full SRS Reduced SRS 1 8 3 2 8 3 3 8 3 4 8 3 5 8 4 6 8 4 7 8 5 8 8 4 17 8 4 Summary In this chapter, a dimension reduction technique using global sensitivity indices is introduced. Since the variances of the Hermite polynomial bases are analytically available, the SRS is suitable to compute global sensitivity indices. Its application to RBDO is also presented. In the RBDO procedure, the global sensitivity indices that are calculated using the lower-order SRS are used to fix unessential random variables, a higher-order SRS with reduced dimension is then used in evaluating the probability constraint. Fixing the unessential random variables accelerates the design optimization process. The RBDO result obtained in this way is compared to that from previous chapter, which shows little difference because of the loss of variability in fixing random variables. CHAPTER 6 FATIGUE RELIABILITY-BASED LOAD TOLERANCE DESIGN Introduction Traditional reliability-based structural design usually makes assumptions on randomness of factors involved in modeling a structural system such as design variables, material properties, etc. These parameters are relatively well controlled so that the variability is usually small. However, it is also important to consider the capacity of the system subject to working conditions, e.g., uncertain loadings, because the uncertainty in load or force is much larger than that of others. The variability of the load is often ignored in the design stage and is difficulty to quantify it. Without knowing the accurate uncertainty characteristics of input load, it is hard to rely on the reliability of the output. In this chapter, a different approach from the traditional RBDO is taken by asking how much load a system can support. The amount of load, which a structural system can support, becomes an important information for evaluating a design. As an illustration, the fatigue reliability-based load tolerance of the front loader frame of CAT 994D wheel loader is studied. Besides the uncertainty in the material properties, which can be incorporated in S-N curve(Ayyub et al. 2002; Chopra and Shack), uncertainties are also investigated on both mean and amplitude of a given dynamic load. Either the variation of load amplitude or mean may affect the fatigue life of the structural system. This research presents a reliability based load design method, which provides the load envelope for a structure subject to fatigue failure mode. Both one dimensional and multidimensional problems are addressed. Since service loads are subjective such that the load characteristic of one operator may completely different from that of others. In order to perform reliability analysis, it is necessary to know uncertainty characteristics of inputs. However, distribution type and parameters of loads are often unknown. In this chapter, instead of modeling variability in parameters by assuming specific type of random distribution, the effect of different distribution types on the system response is investigated by introducing the concept of conservative distribution type, which provides a safer way to model uncertainties. Fatigue Life Prediction Recent developments in the computer-aided analysis provide a reasonable simulation for fatigue life prediction at early design stage for components under complex dynamic loads. For most automotive components, fatigue analysis means to find crack initiation fatigue life. Figure 6-1 illustrates the procedure to the crack initiation fatigue life prediction. FE Model Quasi-xtic Stress FE analysis .. Inference Indexes < Dynamic / Dynamic / Loading )-- Superposition -- Stress History Peak & Value Rain-flow Multi-axial Life Editing Counting Prediction Figure 6-1. Flow chart for fatigue life prediction Crack Initiation Fatigue Life Prediction Two major crack initiation fatigue life prediction methods are stress-based and strain-based methods. The stress-life (S- N) approach employs relationship between the stress amplitude and the fatigue life. This method is based primarily on linear elastic stress analysis. The advantage of stress-life approach is apparent since changes in material and geometry can easily be evaluated and large empirical database for steel with standard notch shape is available. However, the effects of plasticity are not considered in this method. The local strain-life( N ) method assumes that the local strains control the fatigue behavior. The plastic effects are considered well in this method. It is similar to the stress-life approach in that it uses E N curve instead of S N curve, but differs in that the strain is the variable related to the life, and also in that plastic deformation effects are specifically considered. Machine parts are usually required to be durable and able to undertake high numbers of life cycles. The front loader frame of CAT 994D wheel loader is one of such case. The critical position of fatigue failure is usually at welding joints. Because the stress-life method is works well for the brittle material and provides a reasonable approximation for a high cycle fatigue crack initiation life, by taking advantage of the availability of a large amount of available uniaxial fatigue data, stress-life method is employed in this chapter. Since the crack is usually initiated along the component surface, for saving unnecessary computation, FE based fatigue analysis chooses element along surface to calculate the fatigue life. For multi-axial application, the principal stress method has been applied using the planes perpendicular to the surface. Fatigue lives are calculated on eighteen planes spaced at 10 degree increment. On each plane the principal stresses are used to calculate the time history of the stress normal to the plane. It has been shown that this method should only be used for fatigue analysis of 'brittle' metals like cast iron and very high strength steels, as it provides nonconservative results for most ductile metals. Based on the factor that material in our application is cast iron and the interested region is welding joints, the principal stress algorithm can offer the fatigue life calculation with reasonable accuracy. Using superposition of dynamic loadings and the quasi-static FE analysis, the dynamic stresses in the component are used to analyze multi-axial fatigue, based on principal stress using conventional S-N curve(Fe-safe 2004). Most basic fatigue data are collected in the laboratory by testing procedures which employ fully reversed loading. However, realistic service loading usually involves nonzero mean stresses. Therefore, the influence of mean stress on fatigue life should be considered so that the fully reversed laboratory data can be used in the evaluation of real service life. Since the tests required to determine the influence of mean stress are quite expensive, several empirical relationships which related alternating stress amplitude to mean stress have been developed. Among the proposed relationships, two are widely used, which are Goodman and Gerber models. Goodman: (S/ S)+ + (SmS) = 1 Gerber: (Sa/Se)+(Sm/S)2 =1 where S,: Alternating stress amplitude; S,: Endurance stress limit Sm: Mean stress; S,: Ultimate strength 70 Experience shows that test data tend to fall between the Goodman and Gerber curves. In the application of fatigue life prediction of front loader frame of CAT 994D, Goodman relation is applied to address the mean stress effect. Variable Amplitude Loading and Cumulative Damage In real application, components are usually subject to complex dynamic loading which has variable amplitude. It requires identifying cycles and assessing fatigue life for each cycle. The rain-flow counting method(Matsuishi and Endo 1968) is the most commonly used cycle counting technique. This method defined cycles as closed stress- strain hysteretic loops as shown in the figure below: Stress D 0 1- T im e Figure 6-2. Rain-flow and hysteresis An algorithm of rain-flow counting can be developed based on ASTM standard description. Although the rain-flow counting method is not based on an exact physical concept to account for fatigue damage accumulation, it is expected to provide a more realistic representation of the loading history. Cumulative damage of each cycle can be obtained by the Palmgren-Miner hypothesis, which is referred to as the linear damage rule: D= (6.1) in i-th cycle. In Eq. (6.1), D, is the fraction of the damage; n, is the counted number of cycles forj-th stress range; N, is the cycles to fail; and n is the total number of stress ranges counted from rain flow. Failure is predicted to occur if N D >1 (6.2) where N is the number of cycles. Thus, the fatigue life can be calculated as the number of applied load cycle until the cumulative damage reaches 1: Life cycles = N (6.3) YD, 1=1 Model Preparation for Fatigue Reliability Analysis Finite Element Model Figure 6-3 shows the component of a front loader frame of CAT 994D wheel loader, which is subjected to 26 channels of dynamic loading. As show in Figure 6-4, the finite element model consists of 49,313 grid points and 172,533 elements (24 beam, 280 gap, 952 hexagon, 1016 pentagon, 226 quadrilateral, 160,688 tetrahedron, 9,144 triangular, 203 rigid body elements). In order to apply for the displacement boundary 72 conditions and loads, pins are modeled using beam and gap elements. The existence of gap elements makes the problem nonlinear. However, if the gap status does not change during analysis, we can still consider the problem to be linear. 1 APin Lf X 2 APin Lf Y 3 APin Lf Z Pos Upper Hitch 4 APin Lf Z Neg 5 APin Rt X 6 APin Rt Y 7 APin Rt Z Pos 8 APin Rt Z Neg 9 GPin Lf X 10 GPin Lf Y Steer 11 GPin Rt X ylinder Pin 12 GPin Rt Y Lower Hitch 13 YPin Lf X 14 YPin Lf Y Figure 6-3. Front loader frame of CAT 994D wheel loader (subject to 26 channels of dynamic loading) G-Pin Bore- A-Pin Bore Axle Pad Casi,... Y-Pin Bore . 15 YPin Rt X 16 YPin Rt Y 17 UHitch X 18 UHitch Z 19 LHitch X 20 LHitch Z 21 Steer Rt X 22 Steer Rt Z 23 Steer Lf X 24 Steer Lf Z 25 LHitch Y Pos 26 LHitch_Y Neg Clamped Clamped Figure 6-4. Finite element model for front frame Dynamic Load History In Figure 6-3, a total of 26 degrees-of-freedom are chosen for the application of dynamic loads. All loads are located in the center of the pins and the hitches. Even if the dynamic loadf(t) is applied to the system, it is assumed that the inertia is relatively small and the method of superposition can be applied. Thus, only a linear static analysis is enough with the unit load applied to each degree-of-freedom. The stress value from the unit load is called the stress influence coefficient. The dynamic stress can be obtained by multiplying these stress influence coefficients with the dynamic load history. Measured data of 26 channels are used for the dynamic loads with the duration of 46 minutes. This duration is defined as a working cycle. The dynamic load is sampled such that 9,383 data points are available for each channel. Uncertainty in Material Properties and S-N Curve Interpolation Based on available material properties and the component's working conditions, principle stress analysis using the Goodman model is used as the algorithm for fatigue life prediction. From superposition of quasi-static linear finite element analysis and dynamic loading, the stress data are obtained for each element. These stresses can be regard as 'true stress', which means S-N curve can be applied directly on principle stress life method without considering the stress concentrate factor. The S-N curve can be interpolated from nominal stress-life data. Considering the uncertainties of material properties, this interpolation will be implemented in a random manner. A lognormal distribution in S-N curve can be assumed to simplify the randomness. Although there is no rigorous statistical evaluation was performed, but this assumption seems reasonable empirically(Ayyub et al. 2002; Chopra and Shack). Figure 6-5 shows the S-N curve obtained from stress- life data for cast iron used in the front frame. Solid line is the nominal S-N curve and two dashed lines represent the variation in S-N curve. 104 S-N Curve(Logarithmic scale) I 6MP MPa UP -a Figure 6-5. Material S-N curve with uncertainty Uncertainty Modeling of Dynamic Loadings Dynamic loadings are usually very complicated and may involve a lot of uncertainties. Figure 6-6 shows one channel of the dynamic load. The mean and amplitude of dynamic loading usually plays the most important role in fatigue life estimation. Therefore, for the purpose of illustration and simplification, uncertainties can be model based on these two quantities. Figure 6-6. Illustration of one channel of dynamic loads By combining the effects of the randomness of mean and amplitude of the loads, two load capacity coefficients a and y are defined for the mean and amplitude, respectively. The dynamic load can be parameterized as f(t) = mean + y(f,(t)- fmean) (6.4) where a y are random parameters to describe the uncertainties of the loads called load capacity coefficient (LCC). In Eq. (6.4), f0(t) is the original dynamic loads and /fa is the mean value of the initial loads. Due to the random parameters, the dynamic load f(t) shows probabilistic behavior. Equation (6.4) provides a simple two dimensional model of uncertainty in dynamic load history. Note that when both a and y equal to one and fixed, the original deterministic loading history f (t) is recovered. In the following section, one-dimensional problem will first be investigated by fixing one of them. For example, if we fix a at 1 and treat y as a random variable, then uncertainty is modeled for the amplitude of the load. Linear Estimation of Load Tolerance The major challenge of the research is to estimate the load tolerance with respect to the reliability of fatigue life performance, which depends on the load history and uncertainty characterization. Identifying the load distribution is one of the most difficult tasks in the uncertainty analysis because different operating conditions will yield completely different distribution types. At this point, it is assumed that the load has a specific uncertainty characteristic (distribution type and corresponding parameters). When the variance of the load is fixed, for example, it is possible to construct the safety envelope by gradually changing the mean value of the applied load, which requires a large number of reliability analyses. When nonlinearity of the system is small, it is possible to estimate the safety envelope using the sensitivity information at the current load without requiring further reliability analyses. This estimation is based on the first- order Taylor series expansion method. For illustration, one-dimensional models (only considering single random variable) for the effect of amplitude and mean are separately investigated to meet the reliability requirement of fatigue life under uncertainty. In these one dimensional cases, the variation in S-N curve is ignored. Variability of Dynamic Load Amplitude In order to consider the variability of the dynamic load, the mean value of the load is first assumed to remain constant, while the amplitude of the load is varied randomly. From Eq. (6.4), the uncertainty caused by the amplitude can be represented using the following decomposition of the dynamic load: f(t) = fn + (fO(t) -/mean) (6.5) When /= 1, the original load history is recovered. When y= 0, the dynamic load becomes a static load with the mean value. In this definition, cannot take a negative value. Since yis a random variable, it is necessary to assume the distribution type and distribution parameters for y First we assume that y is normally distributed with the mean of one and the standard deviation of 0.25 (COV=0.25). The standard deviation is estimated from the initial dynamic load history. Since the first-order reliability analysis is performed using the standard random variable, we convert into the standard random variable u by (6.6) =1+0.25u where u ~ N(0,12), y N(1,0.252), / = mean, u= standard deviation. For any given sample point u corresponding ycan be obtained from Eq. (6.6), and using ya new dynamic load history can be obtained from Eq. (6.5). By applying this dynamic load history, we can calculate the fatigue life of the system. Since this procedure includes multiple steps, we can construct a (stochastic) response surface for the fatigue life and then perform the reliability analysis using the response surface. Since the fatigue life changes in several orders of magnitudes, it would be better to construct the response surface for the logarithmic fatigue life. In one dimensional case, five collocation points are available from the DOE introduced in chapter 3. Using these five sampling points, a cubic stochastic response surface is constructed as a surrogate model for the logarithmic fatigue life as L(y) Loglo(Life) = 5.7075 -0.7223u- 0.0581(u2 -1) + 0.0756(u3- 3u) (6.7) Note that in one dimensional case, the five collocation points available from previous DOE scheme are sometimes not enough to construct a high fidelity SRS, a couple of complementary sampling points can be chosen to construct a new SRS, which spread evenly between the original collocation points, i.e., four more point in the middle of intervals of the original five points have been chosen. The corresponding SRS logarithmic fatigue life becomes L(y) Logl(Life)= 5.6976- 0.6826u-0.0541(u2 -1)+0.0617(u3- 3u) (6.8) Various quantities for estimating the quality of SRS are shown in Table 6-1 for both five and nine sampling points scheme. The table shows nine points DOE schemes fit the data better based on significant improvement of PRESS (prediction error sum of squares). Table 6-2 lists the t-statistics for the evaluation of each coefficient in the above response surface. Although using more sampling points can increase the fidelity of estimation, it also increases the computational cost. In our specific problem, considering the saving of computation, the five sampling DOE scheme is sufficient. Table 6-1. Quality of response surface Error RMSE SSe R2 R2adj PRESS statistics 5 sampling 0.0406 0.0082 0.9980 0.9921 8.5671 DOE 9 sampling 0.0511 0.0235 0.9965 0.9943 0.1503 DOE Table 6-2. T-statistic of the coefficients coefficient 1 2 3 4 t-statistics of 5 122.7590 15.9279 3.5759 4.0828 sampling DOE t-statistics of 9 229.1431 33.9519 4.9313 6.4894 sampling DOE The response surface in Eq. (6.7) shows that the mean of logarithmic fatigue life is 5.6976 and the standard deviation is about 0.6826. It also shows that the contribution of the higher-order terms is relatively small, compared with the constant and linear terms. Thus, we can conclude that the performance is mildly nonlinear with respect to the random variable. Since the required life of the working component is 60,000 hours and each cycle corresponding to 46 minutes, the target of the fatigue life can be written in the logarithmic scale by Ltarget = Loglo (60,000 hours) = Logo, (78,261 cycles). (6.9) &4.9 The system is considered to be failure when the predicted life from Eq. (6.7) is less than the target life in Eq. (6.9). Accordingly, we can define the probability of failure as Pf,[L(ca)-Ltarget <0]< Parget, (6.10) where Ptarget is the target probability of failure. For example, when Ptarget = 0.1, the probability of failure should be less than 10%. Even though the interpretation of Eq. (6.10) is clear, it is often inconvenient because the probability changes in several orders of magnitudes. In reliability analysis, it is common to use the reliability index, which uses the notion of the standard random variable. Equation (6.10) can be rewritten in terms of the reliability index as Pf (-t) < target target), (6.11) where /fis called the reliability index and D is the cumulative distribution function of the standard random variable. When Ptarget = 0.1, /target 1.3. The advantage of using the reliability index will be clear in the following numerical results. With the response surface in Eq. (6.7), reliability analysis is carried out using the first-order reliability method (FORM) at /,= 1. The results of reliability analysis are as follows: P =17.81% = 0.922456 (6.12) =_ -3.972 where 08//80r is the sensitivity of the reliability index with respect to the mean value of y. Since Ptarget = 0.1 and targett = 1.3, the current system does not satisfy the reliability requirement. From the mild nonlinear property of the response, we can estimate the mean value of y that can satisfy the required reliability. The linear approximation of the mean value can be obtained from /es-mte = 1 -( target)/ = 0.9049, (6.13) which means that the mean value of y needs to be decreased about 10% from the original load amplitude in order to satisfy the required reliability. In order to verify the accuracy of the estimated result, several sampling points are taken and reliability analyses are performed. Figure 6-7 shows the reliability index with respect to u,, while Figure 6-8 shows the probability of failure Pfwith respect to /r. The solid line is linearly approximated reliability using sensitivity information. The reliability 81 index is almost linear and the estimation using sensitivity is close to the actual reliability index. When the target probability of failure is 0.1 and y has the distribution of N(/u,0.252), the safety envelope can be defined as 0_< y_< 0.9049. (6.14) Effect of Load Amplitude 25 1 7 Linear predition 2 ------ ---- --- ------- ------- A Exact value -1 CI I I 15-- --- ------- 0.7 OB 0.9 1 1.1 12 1.3 1.4 1.5 1.6 Figure 6-7. Reliability index /?with respect to random parameter wr 100 Effect of load amplitude P -- Linear prediction -------- ---------- ----- --- -- A Exact value -------- 105 -- --- :::::::::: ::::::: 10ge 07 08 0.9 1 11 12 13 14 1.5 16 Figure 6-8. Probability of failure Pfwith respect to random parameter [r The result means that current design, considering 25% standard deviation in the load amplitude, is not enough to achieve 90% reliability. The structure should work under milder working condition, which means either lower the mean of the load amplitude by about 10% or provide more accurate estimation of the initial load to reduce the variance. Variability of Mean of Dynamic Load Since both mean and amplitude are used to describe the dynamic load history, both of their effects under uncertainty are studied separately. When the mean value of the load is assumed to be varied randomly and the load amplitude remains as the initial load, the uncertainty in load can be modeled as: f(t) = afmean +(fo (t) fmean) (6.15) Same as the case of load amplitude, when a= 1, the applied load is identical to the original load history. When a= 0, the applied load has the same amplitude with the original load history but the mean value is zero. In this definition, a can be a negative value. Since a is a random variable, it is necessary to assume the distribution type and distribution parameters for a. First we assume that a is normally distributed with the mean of one and the standard deviation of 0.25 (COV=0.25). Since the first-order reliability analysis is performed using the standard random variable, we convert a into the standard random variable u by ,a =U+ (6.16) = 1+0.25u where u~ N(0,12), a -N(1,0.252), p= mean, a= standard deviation. By following the same procedure with previous section, using nine sampling points, we can construct a cubic stochastic response surface as a surrogate model for the logarithmic fatigue life as L(a) -Loglo(Life) =5.6906 -0.0905u- 0.0013(u2 -1)-0.0003(u3-3u). (6.17) If we compare the response surface in Eq. (6.17) with the case of amplitude change in Eq.(6.7), the mean values of the both cases are close but the standard deviations are quite different. From this result, we can conclude that the variance of the mean value does not contribute significantly to the variance of the fatigue life. Using the response surface in Eq. (6.17), reliability analysis is carried out using the first-order reliability method (FORM) at / = 1. The results of reliability analysis are as follows: Pf, 08 ,8 6.3435 (6.18) 8= -4.012 Since Ptarget = 0.1 and Aarget = 1.3, the current system satisfies the reliability requirement. The linear approximation of the mean value can be obtained from estimate = 1(/f -1 target )/0=1 2.257, (6.19) which means that the system satisfies the reliability requirement even if the mean value of a is increased up to 225% from the original load. This observation is consistent with the conventional notion of fatigue analysis in which the effect of the amplitude is significant while that of the mean is not. In order to verify the accuracy of the estimated result, several sampling points are taken and reliability analyses are performed. Figure 6-9 shows the reliability index with respect to /a, while Figure 6-10 shows the probability of failure Pfwith respect to /e. The solid line is linearly approximated reliability using sensitivity information. The reliability index is almost linear and the estimation using sensitivity is close to the actual 84 reliability index. When the target probability of failure is 0.1 and a has the distribution of N(/,a,0.252), the safety envelope can be defined as 0 < /l < 2.257. (6.20) Effect of mean of load Figure 6-9. Reliability index /?with respect to random parameter [a effect of mean lU Pf 10,5 rI S i target Linear estimation --- ------- ------- A Exact value 10 25 I i I i I i I 0 0.5 1 1.5 2 2.5 3 3.5 4 Figure 6-10. Probability of failure Pf with respect to random parameter ta Safety Envelope Concept for Load Tolerance Design The safety envelope is defined as the magnitudes of the input design variables when the system fails. When design variables are loads, it is called the load envelope. In one dimensional case, this is simply the range of the allowed loads, e.g., the range of mean value of uc or y in the previous section. In multi-dimensional case, the combination of various loads constitutes an envelope, which is convex in linear systems. Such information is very useful as a capacity of the current design, a future reference for design upgrade, maintenance and control. Figure 6-11 shows a schematic illustration of the safety envelope when two variables are involved. In such a complex situation, a systematic way of searching the boundary of the safety envelope needs to be developed. S Safety envelope SSafe N Figure 6-11. Safety envelope for two variables When the relationship between the safety of the system and the applied loads is linear or mildly nonlinear, linear approximation can produce a very effective way of estimating the safety envelope. In context of reliability based safety measure, the target of safety envelope is that failure probability cannot reach over the prescribed value. Numerical Path Following Algorithm According to the reliability based safety envelope concept introduced above, when target reliability has been specified, a safety envelope can be constructed using numerical path following algorithm to search the boundary of the safety envelope(Allgower and Georg 1990). In this research, a systematic way of searching the boundary of the safety envelope is proposed using a predictor-corrector method, which is similar to the Euler- Newton continuation method(Allgower and Georg 1990; Kwak and Kim 2002). When the relationship between the safety of the structure and the applied loads is linear or mildly nonlinear, this approach can produce an efficient way of estimating the safety envelope. In the context of reliability-based safety measure, the boundary of the safety envelope is the location where the probability of failure is equal to the target probability. (k+l1 k+1) B Boundary ,3 = constant Figure 6-12. Predictor-corrector algorithm The predictor-corrector algorithm(Figure 6-12) is explained below with two random variables, a and y First, the distribution type of random variables is assumed. The effect of different distribution types on the safety envelope can be investigated by following the same procedure as in the previous section. It is clear that the two parameters must have non-negative values, which means that the safety envelope only exists in the first quadrant. The capacity of the structure with respect to ([ta, [ty) is interesting. If the required probability of failure is Ptarget (i.e., target t = '(Ptarget)), the following steps can been taken to construct the safety envelope: |

Full Text |

PAGE 1 RELIABILITY-BASED DESIGN AND LO AD TOLERANCE EVALUATION USING STOCHASTIC RESPONSE SURFACE AND PROBABILISTIC SENSITIVITIES By HAOYU WANG A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLOR IDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2006 PAGE 2 Copyright 2006 by Haoyu Wang PAGE 3 To my family PAGE 4 iv ACKNOWLEDGMENTS I would like to express my a ppreciation to my advisor, Professor Nam-Ho Kim, for his endless encouragement and continuous s upport during my Ph.D. research. Without his guidance, inspiration, experience and willingn ess of sharing his knowledge, this work would have never been possible. Dr. Kim made a tremendous contribution to this dissertation, as well as my prof essional and personal life. I must express my gratitude to the members of my supervisory committee, Professor Raphael T. Haftka, Professor Stanis lav Uryasev, Professor Nagaraj K. Arakere, and Professor Ashok V. Kumar, for their will ingness to review my Ph.D. research and provide constructive comments to help me co mplete this disserta tion. Special thanks go to Professor Raphael T. Haftka for not only his guidance w ith several technical issues during my study, but also comments and sugge stions during group meetings which were extremely helpful for improving my work. Special thanks are also given to Professo r Nestor V. Queipo from the University of Zulia at Venezuela, for his interaction in my research and collaboration in publishing papers during his visiting at the University of Florida. My colleagues in the Structural and Mu ltidisciplinary Optimization Lab at the University of Florida also deserve my grat itude. In particular, I thank Dr. Xueyong Qu, Dr. Amit A. Kale, Dr. Erdem Acar, Tushar Goel, Long Ge, Saad Mukrus for their encourage and help. My parents deserve my deepest appreciati on for their constant love and support. PAGE 5 v Lastly, I would like to tha nk my beautiful and lovely wife, Zhilan, for her love, patience and support during my study. PAGE 6 vi TABLE OF CONTENTS page ACKNOWLEDGMENTS.................................................................................................iv LIST OF TABLES.............................................................................................................ix LIST OF FIGURES...........................................................................................................xi ABSTRACT.....................................................................................................................xiii CHAPTER 1 INTRODUCTION........................................................................................................1 Motivation.....................................................................................................................1 Objective...................................................................................................................... .2 Scope.......................................................................................................................... ...3 Outline........................................................................................................................ ..5 2 LITERATURE SURVEY.............................................................................................7 Uncertainty and Reliability Analys is of Structural Applications.................................7 Reliability-Based Design Optimization......................................................................11 Sensitivity in Reliability Analysis..............................................................................12 Dimension Reduction Strategy...................................................................................13 Robust Design.............................................................................................................13 Fatigue Life Prediction...............................................................................................14 3 UNCERTAINTY ANALYSIS USING STOCHASTIC RESPONSE SURFACE....19 Introduction.................................................................................................................19 Description of Uncertainty Model..............................................................................20 Stochastic Response Surface Method (SRSM)...........................................................22 Polynomial Chaos Expansion (PCE) in Gaussian Space....................................22 Numerical Example of Stochastic Response Surface..........................................26 Improving Efficiency of SRS Using Local Sensitivity Information..........................30 Continuum-Based Design Sensitivity Analysis...................................................31 Constructing SRS Using Local Sensitivity..........................................................34 Numerical Example Torque Arm Model..........................................................36 Summary.....................................................................................................................37 PAGE 7 vii 4 RELIABILITY-BASED DESIGN OPTIMIZATION................................................39 General RBDO Model................................................................................................39 Reliability Index Approach (RIA) and Performance Measure Approach (PMA)......41 Probability Sensitivity Analysis (PSA)......................................................................43 Probability Sensitivity Analysis in FORM..........................................................44 Probability Sensitivity Analysis Using SRSM....................................................47 Reliability-Based Design Optimization Using SRSM................................................48 RBDO with RIA..................................................................................................49 RBDO with Inverse Measure..............................................................................52 Summary.....................................................................................................................53 5 GLOBAL SENSITIVITY ANALYSIS FOR EFFICIENT RBDO............................54 Introduction.................................................................................................................54 Sensitivity Analysis....................................................................................................55 Variance-Based Global Sensitivity Analysis (GSA)..................................................56 Global Sensitivity Analysis Using Polynomial Chaos Expansion.............................58 Adaptive Reduction of Random Design Space Using GSA in RBDO.......................59 Summary.....................................................................................................................65 6 FATIGUE RELIABILITY-BASED LOAD TOLERANCE DESIGN......................66 Introduction.................................................................................................................66 Fatigue Life Prediction...............................................................................................67 Crack Initiation Fatigue Life Prediction..............................................................68 Variable Amplitude Loading and Cumulative Damage......................................70 Model Preparation for Fatigue Reliability Analysis...................................................71 Finite Element Model..........................................................................................71 Dynamic Load History........................................................................................73 Uncertainty in Material Properties and S-N Curve Interpolation........................74 Uncertainty Modeling of Dynamic Loadings.............................................................75 Linear Estimation of Load Tolerance.........................................................................76 Variability of Dynamic Load Amplitude............................................................77 Variability of Mean of Dynamic Load................................................................82 Safety Envelope Concept for Load Tolerance Design...............................................84 Numerical Path Following Algorithm.................................................................85 Example for Multi-Dimensional Load Envelope................................................88 Conservative Distribution Type..................................................................................90 Summary.....................................................................................................................92 7 ROBUST DESIGN USING STOCHA STIC RESPONSE SURFACE......................94 Introduction.................................................................................................................94 Performance Variance Calculation Using SRS..........................................................96 Variance Sensitivity....................................................................................................97 Robust Design Two-layer Beam............................................................................102 PAGE 8 viii Dynamic Response of Two-Layer Beam..........................................................102 Robust Design for Two-Layer Beam................................................................103 Global Sensitivity Analysis...............................................................................106 Robust Design by Tolerance Control.......................................................................107 Summary...................................................................................................................111 8 SUMMARY AND RECOMMENDATIONS..........................................................113 APPENDIX A SAMPLING-BASED PROBABILITY SENSITIVITY ANALYSIS FOR DIFFERENT DISTRIBUTION TYPE.....................................................................115 Normal Distribution 2(,)iiiXN .......................................................................115 Case 1: i....................................................................................................115 Case 2: i....................................................................................................116 Uniform Distribution................................................................................................117 Log-Normal Distribution..........................................................................................119 B NATURAL FREQUENCY OF CANT ILEVER COMPOSITE BEAM..................122 Bending Moment......................................................................................................122 Geometric Properties of Composite Beam...............................................................122 Effective Compliance for Composite Beam.............................................................123 Effective Mass for Composite Beam........................................................................123 LIST OF REFERENCES.................................................................................................125 BIOGRAPHICAL SKETCH........................................................................................... 134 PAGE 9 ix LIST OF TABLES Table page 3-1. The type of polynomials and correspond ing random variables for different AskeyChaos (N 0 denotes a finite integer)........................................................................22 3-2. Root mean square error of PDF co mpared with the exact PDF of performance function y = ex.............................................................................................................27 3-3. Comparison of probability of G>520M Pa obtained from different uncertainty analysis methods (Full sampling w ithout using local sensitivity)............................30 3-4. Comparison of probability of G>520M Pa obtained from different uncertainty analysis methods (reduced sampling using local sensitivity) ;................................37 3-5. Comparison of probability of G>520MPa obtained with/without local sensitivity (7/27 sampling points) using 2nd order SRS.............................................................37 4-1. Probability sensitivity with respect to random parameters (unit: centimeter)............46 4-2. Computational efficiency of analy tical method for probability sensitivity................47 4-3. Definition of random design variable s and their bounds. The values of design variables at optimum design are shown in the 5th column (unit: centimeter).........50 4-4. Reliability Index of active constraint at optimal design.............................................52 5-1. Variances of the Hermite bases up to the second order..............................................58 5-2. Global sensitivity indices considering on ly main factors for the torque arm model at the initial design. Only three random variables ( u2, u6, and u8) are preserved when a threshold value of 1.0% is in place..............................................................63 5-3. Comparison of the number of rando m variables in each design cycle. The threshold of 1.0% is used. The first constraint is listed............................................64 6-1. Quality of response surface........................................................................................78 6-2. T-statistic of the coefficients......................................................................................79 7-1. Random variables for can tilevered beam structure....................................................99 PAGE 10 x 7-2. Variance estimation of linear performance (strength)..............................................100 7-3. Variance estimation of non linear performance (deflection).....................................101 7-4. Sensitivity of variance fo r linear performance (strength).........................................102 7-5. Sensitivity of variance for nonlinear performance (deflection)................................102 7-6. Random parameters for the composite beam structure............................................104 7-7. Sensitivities of objective f unctions at the initial design ( ts = 6m, tp = 0.2m, L = 1000m).................................................................................................................105 7-8. Total sensitivity indices for the composite beam structure ( ts = 6m, tp = 0.2m, L = 1000m)..............................................................................................................107 7-9. Sensitivity of variance fo r linear performance (strength).........................................109 7-10. Sensitivity of variance for nonlinear performance (deflection)..............................109 7-11. Random variables and cost-tolerance functions.....................................................110 7-12. Random variables and cost-tolerance functions.....................................................111 A-1: Accuracy of proposed probability sensitivity method for normal distribution using 200,000 sampling MCS................................................................................116 A-2: Accuracy of proposed probability se nsitivity method for uniform distribution using 200,000 sampling MCS................................................................................119 A-3: Accuracy of proposed probability se nsitivity method for Log-normal distribution using 200,000 sampling MCS................................................................................121 PAGE 11 xi LIST OF FIGURES Figure page 3-1. Limit state function divides the safe region from the failure region..........................21 3-2. PDF of performance function () x y xe ....................................................................27 3-3. Shape design parameters for the torque arm..............................................................28 3-4. PDF of performance function G(x) torque arm model............................................29 3-5. Variation of a structural domain according to the design velocity field V(x)............32 3-6. PDF of performance fuction G(x) Torque model at in itial design (SRS with sensitivity)................................................................................................................36 4-1. Flow chart for reliability-based design optimization..................................................43 4-2. Optimum design and stress distributi on of the torque arm model with 8 random variables...................................................................................................................51 4-3. Optimization history of cost function (mass) for the torque arm model with 8 random variables......................................................................................................51 4-4. PDF of the performance function at the optimum for the torque-arm problem.........52 5-1. Global sensitivity indices for to rque arm model at initial design...............................59 5-2. Adaptive reduction of unessential random design variables using global sensitivity indices in RBDO. Low-order SRS is used for global sensitivity analysis, while a high-order SRS is used to evaluate the reliability of the system..61 5-3. Optimum designs for the full SRS (solid line) and adaptively reduced SRS (dotted line). Because some variables are fixed, the interior cutout of the design from the adaptively reduced SRS is larger than that from the full SRS...........................64 6-1. Flow chart for fatigue life prediction..........................................................................67 6-2. Rain-flow and hysteresis............................................................................................70 6-3. Front loader frame of CAT 994D wheel loader (subject to 26 channels of dynamic loading).....................................................................................................................72 PAGE 12 xii 6-4. Finite element model for front frame.........................................................................73 6-5. Material S-N curve with uncertainty...........................................................................74 6-6. Illustration of one channel of dynamic loads..............................................................75 6-7. Reliability index with respect to random parameter ...........................................81 6-8. Probability of failure Pf with respect to random parameter ...................................81 6-9. Reliability index with respect to random parameter ...........................................84 6-10. Probability of failure Pf with respect to random parameter .................................84 6-11. Safety envelope for two variables............................................................................85 6-12. Predictor-corrector algorithm...................................................................................86 6-13. Construction of load envelope..................................................................................89 6-14. Safety envelop for fatigue reliab ility of CAT 994D front loader frame...................90 6-15. Reliability index with respect to random parameter .........................................91 6-16. Probability of failure Pf with respect to random parameter .................................91 6-17. 2-D safety envelope for different dist ribution type with sa me random parameters.92 7-1. Cantilever beam subject to two direction loads..........................................................99 7-2. Piezoelectric cantilevered composite beam..............................................................103 7-3. Pareto optimal front for the robust design of the composite beam...........................106 B-1: Free body diagram of two-layer beam.....................................................................122 PAGE 13 xiii Abstract of Dissertation Pres ented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy RELIABILITY-BASED DESIGN AND LO AD TOLERANCE EVALUATION USING STOCHASTIC RESPONSE SURFACE AND PROBABILISTIC SENSITIVITIES By Haoyu Wang December 2006 Chair: Nam-Ho Kim Major: Mechanic al Engineering Uncertainty is inevitable in structural desi gn. This research presents an efficient uncertainty analysis technique based on stocha stic response surfaces (SRS). The focus is on calculating uncertainty propagation using fewer number of function evaluations. Due to sensitivity analysis, the gradient info rmation of the performance is efficiently calculated and used in constructing SRS. Based on SRS, reliability-based design optim ization (RBDO) is studied intensively in this research. Probability sensitivity an alysis using the sampling technique is also proposed. Since the computational cost of RBDO increases signifi cantly proportional to the increasing number of random variables, gl obal sensitivity analysis is introduced to adaptively reduce unessential random variable s. It has been shown that the global sensitivity indices can be calculated analyt ically because the SR S employs the Hermite polynomials as bases. PAGE 14 xiv Traditional structural desi gn focuses on designing a reli able structure under well characterized random factors (dim ensions, shape, material prop erties, etc). Va riations of these parameters are relatively small a nd well characterized. However, everyday engineering life tends to use the existing structural part in a different applications instead of designing a completely new part. In this research, a reliability-based safety envelope concept for load tolerance is introduced. This shows the capac ity of the current design as a future reference for design upgrade, maintena nce and control. The safety envelope is applied to estimate the load tole rance of a structural part with respect to the reliability of fatigue life. Stochastic response surface is also applied on robust design in this research. It is shown that the polynomial chaos expansion with appropriate bases provides an accurate and efficient tool in evaluati ng the performance variance. In addition, the sensitivity of the output variance, which is critical in the mathematical programming method, is calculated by consistently differentiating the polynomial chaos expansion with respect to the design variables. A reliability-based robust design method that can reduce the variance of the output performance as well as the deviation of the mean value is proposed using SRS and efficient sensitivity analysis Numerical examples are shown to verify accuracy of the sensitivity in formation and the convergence of the robust design problem. PAGE 15 1 CHAPTER 1 INTRODUCTION Motivation A typical mechanical design procedure incl udes two steps: first, a design space is defined and a mathematical model is established, which includes the objective function and required constraints. Second, a proper opt imization algorithm is selected properly based on this mathematical model to solve th e design problem. In engineering design, the deterministic optimization model has been studied intensively to reduce the objective function by pushing design to the limits of system failure boundaries. However, everything in the real world involves uncertain ties, and so does the design of mechanical components. After realizing deterministic design leaves very li ttle or no room for tolerances of the imperfections in desi gn, manufacturing and variety of service conditions, design engineers incorp orate a safety factor into th e structural design to leave safety margins. Without considering uncer tainties and probabilistic quantification, deterministic design with a safety factor may be either unsafe or too conservative. Motivated by overcoming the bottleneck of the deterministic design, the reliabilitybased design optimization (RBDO) model has become popular in past two decades since uncertainties exist everywhere in every phase of the structure system, from design and manufacturing to service and maintenance. If elements in the mathematical model are considered to be probabilistic with certai n types of random dist ribution, the design problem becomes a typical RBDO problem. Th e probabilistic elemen ts can be design variables, material properties, applied loads, etc. One of the most important issues in PAGE 16 2 RBDO is a good model of uncertainty propaga tion in mathematical models. Besides RBDO, which only considers the failure mode as a constraint in the probabilistic point of view, robust design will also be considered in this research in order to design a structure less sensitive to the existing uncertainty factors. In optimization point of view, that means minimization of performance variance. For a certain design, it is also important to consider the service capability of the system subject to applied loads since engine ers tend to use the same design in different applications instead of a completely new desi gn. Another motivation of this research is the load tolerance design. A good estimation of load toleran ce shows the capacity of the current design, future reference for design upgrade, maintenance and control. Since static or quasi-static loading is rare ly observed in modern engineer ing practice, the majority of engineering design projects i nvolves machine parts subjecte d to fluctuating or cyclic loads. Such loads induce fluctuating or cy clic stresses that often result in failure by fatigue. In addition, because service load s are subjective, which means the load characteristic of one operator may be complete ly different from that of the other, it is necessary to consider the uncertainties wh ile estimating the load tolerance of dynamic systems. Objective Uncertainty in the design parameters makes structural optimization a computationally expensive task due to the significant number of structural analyses required by traditional methods. Critical issues for overcoming these difficulties are those related to uncertainty characterization, uncertainty propagation, ranking of design variables, and efficient optim ization algorithms. Conventional approaches for these tasks PAGE 17 3 often fail to meet constraints (computational resources, cost, time, etc.) typically present in industrial environments. The first objective of this re search is to develop a com putationally efficient method for uncertainty propagation. Local and global sensitivities can then be used to improve the efficiency of estimating uncertainty propa gation. Besides efficiency, the accuracy and applicability of the methods to a wide range of applications need to be addressed. The second objective of this research is to develop a computationally efficient RBDO and robust design framework based on proposed uncertainty analysis. In the gradient-based algorithm, the sensitivity info rmation is required during the optimization procedure. The computational cost can be si gnificantly saved if the gradient can be obtained analytically, instead of using the finite difference method. The probabilistic sensitivity analysis is utilized to calculate the gradient of the reliability constraints. In the framework of robust design, sensitivity analysis of performance variance is also studied. Traditional structural de sign usually makes assumpti on on randomness of factors involved in modeling a st ructural system such as design variables, material properties, etc. However, it is also important to cons ider the capacity of th e system subject to uncertain loads. The final objective of this res earch is to present a reliability-based load design method, which provides the safety enve lope, for a structure subject to fatigue failure. Scope In the standard framework of RBDO, c onstraints are provided in terms of the probabilities of failure. The uncertainties involved in the system are modeled by assuming random input variables with a certai n type of probabilistic distribution. RBDO achieves the design goal by meeting the requirement of structural reli ability constraints. PAGE 18 4 The RBDO involving a computationally de manding model has been limited by the relatively high number of required analyses fo r uncertainty propagation during the design process. The scope of this research is to present an efficient uncertainty propagation technique based on stochastic response surf aces (SRS) constructed using model outputs at heuristically selected collocation points. The efficiency of the uncertainty propagation approach is critical since the response surf ace needs to be reconstructed at each design cycle. In order to improve the efficiency, th e performance gradient, calculated from local sensitivity analysis, is used. Even if the local sensitivity inform ation can reduce the number of required simulations, the dimension of the SRS is st ill increased accordi ng to the number of random variables. If the contribution of a ra ndom variable is rela tively small to the variance of the model output, it is possibl e to consider the random variable as a deterministic one. In this research, the global sensitivity index is used for that purpose. The role of the global sensitivity is to quan tify the model inputs contributions to the output variability, hence establishing which factors influence the model prediction the most so that i) resources can be focused to reduce or account for uncertainty where it is most appropriate, or ii) une ssential variables can be fixed without significantly affecting the output variability. Reliability constraint in RBDO require s probability sensitivity analysis for gradient-based algorithms. In this resear ch, both FORM-based and sampling-based reliability sensitivity analysis are investigat ed. The analytical expression for probability sensitivity based on SRS is derived and used for RBDO. PAGE 19 5 Variations in dynamic loads are usually too complicated to be predicted. A simplified uncertainty modeling technique base d on the mean and amplitude of the load history is proposed. Using the unc ertainty in the load history, a reliability-based safety envelope is constructed that ca n provide load tolerance of the current design. In addition, the effect of different distri bution types is investigated so that the design engineers can choose the conservative distribution type. This research involves uncertainty mode ling and quantification, design sensitivity analysis, fatigue life prediction, reliability-b ased design optimizati on (RBDO) and robust design. Methodologies investigated or applie d in reliability analysis include momentbased methods such as firstand seco nd-order reliability method (FORM/SORM), approximation methods such as Monte Carlo Si mulation (MCS) with stochastic response surface method (SRSM). Furthermore, sensitivity analysis for reliability constraints of RBDO is investigated to improve the computati onal cost involved in reliability analysis and design. Performance variance and sensit ivity are calculated based on SRSM for robust design. Computationally affordable reliability-based optimization and robust design method, and safety envelope for load tolerance are presented in this work. Outline A literature survey is presented in chapter 2, which includes all aspects involved in this research such as reliability analysis reliability-based desi gn optimization, robust design, sensitivity analysis, dimension reduction strategy and fatigue analysis. Chapter 3 describes the uncertainty modeli ng and widely used reliability analysis methods. A stochastic respons e surface method (SRSM) coupl ed with the sensitivity analysis of performance measure is introduced It is shown that the local sensitivity information improves computational effici ency significantly by reducing required PAGE 20 6 number of samples. Convergence and accuracy of the proposed SRSM scheme are also discussed in this chapter. In Chapter 4, the mathematical model is defined for RBDO. RBDO using either direct probability measure or inverse m easure is investigated and compared. The difference of numerical pro cedures between RBDO and deterministic optimization are also compared. As required by RBDO, probability sensitivity analysis is studied in this chapter. In Chapter 5, a dimension reduction stra tegy is proposed by introducing the concept of variance-based global sensitivity analys is, which saves the computational resources further by fixing the unesse ntial design variables. Chapter 6 demonstrates a fatigue reliabili ty-based load tolera nce design by using reliability sensitivity information. A reliability -based safety envelope is constructed by path following continuation method. Chapter 7 proposes an optimization model for robust design where SRS is used to calculate the performance variance and its sensitivity. Chapter 8 concludes this research followed by recommendations for future research work. PAGE 21 7 CHAPTER 2 LITERATURE SURVEY Uncertainty and Reliability Analysis of Structural Applications Reliability-based design optimization ( RBDO) provides tools for making decision within a feasible domain of design variab les with consideration of uncertainties underneath. In the past decades, tremendous amo unt of work has been carried out in this area and it is still moving forward. Compared to deterministic optimization, design variables included in RBDO are random and usually modeled with specifi c distribution types, so do the random parameters such material properties as Youngs modulus and Poissons ratio. Usually random parameters do not change during optimiza tion, but their effect s to the probability propagation must be counted due to its uncertainty. Reliability, which is defined as the probability that a system response does not exce ed the limit threshold, is often considered as constraints. The system response is a function of design variables and random parameters, which is called a performance f unction in this research Performance function is usually implicit and nonlin ear prediction of random vari ables, making probabilistic description of a system response difficult. Several approximation methods for reliability analysis ha ve been developed in the literatures. Among them, Monte Carlo Simula tion (MCS) (Metropolis and Ulam 1949; Rubinstein 1981) has been widely used due to its simplicity and dependability. However, the large sample size re quired in MCS in order to reduce the noise and to reach a certain level of accuracy makes it pr actically formidable in computationally intensive PAGE 22 8 engineering applications, such as Finite El ement Analysis (FEA). Even improved version of MCS are developed, such as importanc e sampling, Latin Hypercube Sampling (Wyss and Jorgensen 1998), Stratified Sampling, etc, they are still expe nsive in structural reliability analysis. Moment-based methods (Breitung 1984; Ha ldar and Mahadevan 2000; Hasofer and Lind 1974) have been developed to provide less expansive calc ulation of the probability of failure compared to MCS. However, they are limited to a single failure mode. As the most widely used moment-based methods, the development of the theory of Firstand Second-Order Reliability Method (FORM/SORM ) is claimed to be finished and only technical work left to do (Rackwitz 2000). FORM/SORM are based on the linear/quadratic approximation of the lim it state function around most probable point(MPP), which is defined in standard normal space as the closest point from the origin on the response surface. For highly nonlinear problems, predictions of reliability from FORM/SORM are not accurate enough because they approximate the response using a linear or quadratic function. The response surface method (Khuri and Cornell 1996; Myers and Montgomery 1995) is proposed to resolve this difficult y. This method typically employs polynomials bases to approximate the system performance and facilitate reliability analysis with little extra computational cost by combining with MC S. Since the accuracy of MCS with fixed sample size relies on the seeking level of probability of failure which sometimes is extremely low in structural design, the probability calculated by MCS near optima is too rough to represent the true valu e of failure probability. Reliab ility analysis using safety factor (Wu et al. 2001) or probability suffici ency factor (PSF) (Qu and Haftka 2004) is PAGE 23 9 proposed to ameliorate this effect. With the PSF as the constraints in RBDO, the variation of magnitude of constraints is usually several orders of magnitude lower than that of the probability of failure, and so is the magnitude of the numerical noise caused by MCS. One of the significant advantages of th e moment-based approach is that the sensitivity of the system reliability or proba bility of failure can be obtained with very little extra computation (Yu et al. 1998). However, moment based approach such as FORM/SORM still has limitations when the performance function is highly nonlinear (Ghanem and Ghiocel 1996). The evaluation of the probabilistic constraints may have large errors in this case. Mahadevan a nd Shi (Mahadevan and Shi 2001) presented a multipoint linearization method (MPLM) for th e reliability analysis of nonlinear limit states, which determines the multiple lineari zation points through the secant method. It increases the complexity of the probl em with limited accuracy improvement. The response surface method can approximate the system response and with little extra computation for MCS, the probability of failure can easy to be obtained. Compared to the conventional deterministic design re sponse surface, Stochastic Response Surface (SRS) (Isukapalli et al. 1998) has the advantage that it only approximates the function around most probability region which highly improved accuracy. Another advantage of SRS is the choice of basis function. The monomial bases(Qu et al. 2000) are widely used due to its simplicity. Other pol ynomial bases are also being studied intensively such as radial basis function (RBF) (Krishmamurthy 2003), orthogonal polynomials (Xiu et al. 2002),etc. Since Ghanem and Spanos proposed the spectral approach of stochastic finite element method (Ghanem and Spanos 1991), the homogeneous Polynomial Chaos Expansion (PCE ) has been widely utilized to represent PAGE 24 10 the uncertainties due to the nature of stocha stic process. To make better approximation with less model analyses, sampling methods ar e studied intensively. Different sampling methods were studied and brought in different applications recent years, such as Latin Hypercube Sampling (LHS) (Choi et al. 2003; Qu et al. 2000) and collocation sampling method (Webster et al. 1996). In the collo cation method, Webster and Tatang derived a set of polynomials from the probability density function of each input parameter such that the roots of each polynomial are spread out over the high probability region of the parameter by deriving orthogonal polynomials Because the uncertainty is usually evaluated by transforming all the random variables and parameters into the Gaussian space, the corresponding orthogonal polynomials are Hermite polynomials. To obtain additional accuracy of SRS, moving leas t square (MLS) method (Youn 2001, Dec; Youn and Choi 2004) is proposed by in troducing weight functions. The number of simulations can be reduced if the sensitivity information is available. Isukapalli (Isukapall i et al. 2000) used an automa tic differentiation program to obtain the sensitivity and utilized it in cons tructing the response surface. However, the computational cost for automatic differentiati on is usually very hi gh(Van Keulen et al. 2004), which reduces the signifi cance of the method. Design sensitivity analysis can provide analytical sensitivity information of response with little extra computation (Kim et al. 2000). Thus, coupling the regression based stochastic re sponse surface method (SRSM) with sensitivity can save large amount of computational co st, especially when the required number of design variab les is large (Kim et al. 2004b). Several methods(Lauridensen et al. 2001; Malkov and Toropov 1991; Rijpkema et al. 2000; Van Keulen et al. 2000) have been pr oposed to use sensitivity information in PAGE 25 11 constructing response surface. Vervenne (Ver venne 2005) proposed a gradient-enhanced response surface method based on above men tioned methods. He developed a two-step approach is proposed: first, different re sponse surfaces using function values and derivatives are constructed separately; Sec ond, these response surfaces are combined together to form a single response surface wh ich fits as good as possible for both function value and response surfaces. In his study, several types of response surface and different combination scheme have been compared. Reliability-Based Design Optimization As mentioned in the previ ous section, FORM/SORM perf orms reliability analysis through linear/quadratic a pproximation of the performance function at MPP. Thus, searching MPP is the main task for moment-based RBDO. However, most advanced MPP search methods such as two point ad aptive nonlinear approximation method (TPA) (Grandhi and L.P. 1998; Wang and Grandhi 199 5; Xu and Grandhi 1998) or hybrid mean value (HMV) method(Youn 2001, Dec; Youn et al. 2003) can not make significant improvement of efficiency in the computational cost(Du and Chen 2002b). In conventional RBDO, the probability cons traint is describe d by the reliability index, which in FORM is the shortest distance from the origin to the limit state in standard normal space. This approach is called reliability index approach (RIA). By modifying the formulation of probabilistic constraints, Tu proposed an inverse measure approach, called Performance Measure Appr oach (PMA) (Tu 1999; Tu and Choi 1997; Tu et al. 1999; Tu 2001) which is proved to be consistent with the RIA but is inherently robust and more efficient if the probabilisti c constraint is inactiv e. Both RIA and PMA employ double loop strategy with analysis loop (inner loop for re liability anal ysis) nested within the synthesis loop (outer loop for design optimization). PAGE 26 12 Due to the nature of double loop optimiza tion, the computational cost is usually high. A couple of new strategies were proposed to improve the efficiency (Yang and Gu 2004). Sequential Optimization and Reliabili ty Assessment (SORA) method (Du and Chen 2002b) decouples optimization loop from the reliability analysis loop and each deterministic optimization loop followed by a se ries of MPP searches. This method shifts the boundaries of violated constraints to the feasible direction based on the reliability results obtained in the previous cycle. Thus it improves design quickly from cycle to cycle and ameliorates the computational e fficiency. Other single loop methods(Chen et al. 1997; Kwak and Lee 1987; Liang et al 2004; Wang and Kodiyalam 2002) are also developed to provide efficient RBDO. In this method, the relationship between random variables and its mean is found through the Karush-Kuhn-Tucher (KKT) optimality condition. The double loop RBDO formulati on is transformed to a single loop deterministic optimization problem and expensive MPP search is avoided. However, there is no guarantee that an active reliability constraint co nverges to its own MPP, and the required reliability may not be satisfied. Sensitivity in Reliability Analysis When RBDO problems are solved using gradient-based optimization algorithms, sensitivities of reliability or probability of failure with resp ect to the design parameters are required. Probability sensit ivity can be used to identify those insignificant random variables during the design stag e. In the moment-based approaches such as FORM, the sensitivity can be obtained accompanied by the reliability analysis without extra function evaluation once MPP is located (Karamcha ndani and Cornell 1992; Yu et al. 1997). Wu(Wu 1994) proposed an adaptive importa nce sampling(AIS) method to calculate reliability and AIS-based reliability sensitiv ity coefficients. Liu et al(Liu et al. 2004) PAGE 27 13 compare four widely-used probability sensit ivity analysis(PSA) methods, which include Sobol indices, Wus sensitivity coefficients, the MPP based sensitivity coefficients and the Kullback-Leibler entropy based method. Th e merits behind each method are reviewed in details. Dimension Reduction Strategy In reliability analysis, the computational cost of multidimensional integration is high. Xu and Rahman(Rahman and Xu 2004; Xu and Rahman 2004) use series expansions to decompose the multidimensional problem to lower dimensional integration, such as univariate and bivariate integrations. Compared to multidimensional integration, the total computation of univa rate integrations is much lower. Recent development in statistics intr oduces global sensitivity analysis (GSA)(Saltelli et al. 2000; Saltelli et al. 1999; Sobol 1993; Sobol 2001), which studies how the variance in the output of a computational model can be appor tioned, qualitatively a nd quantitatively, to different sources of variation. Considering the contribution of the variance of design variables to performance vari ances are not of same importance, Kim et al proposed an adaptive reduction method using total sens itivity indices to reduce the problem dimensions (Kim et al. 2004a). Robust Design Robust design, known as Taguchi parame ter design (Taguchi 1986; Taguchi 1987), is to design a product in such a way that the performance variance is insensitive to variation of design variab les which is beyond the control of designer. Wang & Kodiyalam(Wang and Kodiyalam 2002) formul ated robust design as an optimization problem by minimizing the variation of system re sponse. Since the material cost has to be considered as well as manufacturing cost, Ch en and Dus formulation compromises cost PAGE 28 14 reduction with performance variance control( Du and Chen 2002a). A robust design can also be achieved by using traditional optimization techniques to minimize the performance sensitivities. Chen & Choi formulated the robust design by minimizing a total cost function and sum of squares of magnitudes of first-order design sensitivities(Chen and Choi 1996), which requires the evaluation of second-order sensitivity analysis. This is a different philosophy compared to the variance based approach. It is more focus on the local be havior of the system performance and can achieve local robustness. The final design by minimizing local sensitivity cannot guarantee the robustness of system globally if the input variances are considerable. By summarizing approaches popularly applied in robust design, Park et al. (Park et al. 2006) define robust design methodologies into two different category: Taguchi method and robust optimization. Under the c ontext of multi-scale and multi-disciplinary applications, Allen et al.(Allen et al 2006) reviewed robust design methods and categorizes robust design into four differen t types based on the sources of variability. Fatigue Life Prediction In 1829, Albert found that a metal subjected to a repeated load will fail at a stress level lower than that required to cause failure on a single application of the load. Then the question comes out: how parts fail under time -varying or non-static conditions? Such phenomenon is called fatigue. The first approach developed to carry out fatigue analysis is the nominal stress method, which is still wi dely used in applicat ions where the applied stress varies with constant amplitude within the elastic range of the material and the number of cycles to failure is large. The nominal stress method works well in high cycle fatigue analysis but does not fit for the low cy cle fatigue analysis where the material has a significant part in the plastic region. PAGE 29 15 August Whler(Wohler 1860) carried out expe riments to obtain a plot of cyclic stress level versus the loga rithm of life in mid-19th cen tury, which is well known as S-N curve. Basquin proposed a stress-life( S-N ) relationship(Basquin 1910) which can be plotted as a straight line using log scales. S-N approach is applicable to situations where cyclic loading is essentially elastic, so the S-N curve should be confined on the life axis to numbers greater than about 105 cycles in order to ensure no significant plasticity occurs. Most basic fatigue data are collected in the laboratory by testing procedures which employ fully reversed loading. However, mo st realistic service loads involve non-zero mean stresses. Therefore, the influence of mean stress on fatigue life should be considered so that the fully reversed laborator y data can be used in the evaluation of real service life. Since the tests required to dete rmine the influence of mean stress are quite expensive, several empirical relationsh ips(Gerber 1874; Goodman 1899; Soderberg 1939) which related alternating stress amplitude to mean stress have been developed. Among the proposed relationships, two are widely used, which are based on Goodman(Goodman 1899) and Gerber(Gerber 1874). S-N approach works well when the cyclic loading is essentially elastic, which means in high cycle fatigue life evaluation. Wh ile using this method, it assumes that most of the life is consumed by nucleating cr acks (around 0.01 mm) and nominal stresses and material strength control fatigue life. Accura te determinations of miscellaneous effects factor Kf for each geometry and mate rial are also required. The advantage of S-N approach is apparent since changes in material and geometry can easily be evaluated and large empirical data base for steel with standard notch shape is available. However, the limitation should also be accounted. This method does not PAGE 30 16 consider the effects of plastic ity, and mean stress effect eval uation is often in error. As the matter of fact, the requirement of empirical Kf for good results is also a kind of disadvantage. As mentioned above, when the cyclic lo ads are relatively large and have a significant amount of plastic deformation, the components will suffer relatively short lives. This type of fatigue behavior is called low-cycle fatigue or strain-controlled fatigue. The analytical procedure in dealing with strain-contro lled fatigue is called the strain-life, local stress-strain or criti cal location approach. In 1950s, Coffin and Manson(Coffin 1954; Manson 1954) suggested that the plastic streng th component of a fatigue cycle may also be considered in fatigue life prediction by a simple power law. In order to account for the mean stress effects, two correction methods are proposed by Morrow and Smith, Topper & Watson (STW)(Smith et al. 1970), respectively. Local Strain-Life (-N ) method assumes that the local stresses and strains control the fatigue behavior. In this method, the plastic effects and mean stress effects are considered well. The limitation is that it also needs the empirical Kf. In the local strainlife approach, the most of the life is cons umed by micro-crack growth (0.1-1mm). To account for macro-crack growth(>1mm), the fracture mechanics-based crack propagation method is proposed(Hoeppner and Krupp 1974; Pari s 1964; Paris and Erdogan 1963). In this method, major assumpti on is that nominal stress and crack size control the fatigue life and the initial crack size is determined accurately. It is the only method to directly deal with cracks. However, the complex sequence effects and accurate initial crack size are difficult to be determined. PAGE 31 17 Linear elastic fracture mechanics (LEFM) is a new branch of engineering. The earliest work was done by Inglis(Inglis 1913) but the major developments were carried out by Griffth(Griffith 1921) at Royal Aircraft Factory(RAF,UK) in 1921 and Irwin(Irwin 1956) in the USA in 1956. In LEFM theory, the driving force for a crack to extend is not the stress or strain but the stress intensity factor, known as K. The stress intensity factor uniquely describes the crack tip stress field independent of global geomet ry by embodying both the stress and the crack size. The relationship of the crack growth in the sense th at the rate of crack growth, da/dN, with respect to the cyc lic range of the stress intensity factor, K, was derived by Paul C. Paris(Paris et al 1961) in 1961, known as Paris Law. In reality, mechanical component are seldom subjected to purely constant amplitude loading history. The irregular stress history must be counted as a series of constant amplitude stresses. In addition, it is difficult to define a cycle in an irregular stress history. Since the reverse of stre ss curve can be easily f ound according to the sign change of the stress history curve, cycle c ounting techniques such as rain-flow counting method(Matsuishi and Endo 1968) are developed to combine reversals to form cycles. After that, cumulative damages can be calculated by Miners law(Miner 1945). For most realistic structures or components, stress or strain fi elds are multi-axial. Fatigue life prediction methods for mu lti-axial loading also have been developed(Bannatine et al. 1990; Fuchs and Stephens 1980; Miller et al. 1966). In addition to some traditional method such as maximum principle stress/strain method, maximum shear stress/strain method and Von Mises effective st ress/strain method; Miller and Brown formulized the critical plane approach(Brown and Miller 1973) from PAGE 32 18 the observation that the stress and strain normal to the plane with maximum shear has been recognized to strongly influence the development of fatigue crack. No consensus has been reached on the methods of multi-axia l fatigue life prediction. All these methods have their own advantages in the specific application. So far, fatigue life analysis has been separated into two categories, (a) crack initiation, including S-N and -N method, and (b) crack propaga tion. The criteria for the fatigue life of a component in engineering de sign depend on material properties or work conditions. In general, the automotive industr y usually applies crack initiation criteria because of the nature of the product and use. On the other hand, the aircraft industry mainly uses crack propagation criteria by periodic inspection and fatigue crack monitoring to achieve and maintain structural safety. PAGE 33 19 CHAPTER 3 UNCERTAINTY ANALYSIS USING STOCHASTIC RESPONSE SURFACE Introduction Uncertainty modeling and reliability analysis are the key issues in the reliability based design process. Uncertainty modeling can be decomposed into three fundamental steps: i) uncertainty characterization of mode l inputs, ii) propagati on of uncertainty, and iii) uncertainty management/decision making. The uncertainty in model inputs can be represented in terms of standa rdized normal random variables ( srv ) with mean zero and variance equal to one. The selec tion is supported by the fact that they are widely used and well-behaved. For other types of random variab les, an appropriate tr ansformation must be employed. It is assumed that the model inputs are independent so each one is expressed directly as a function of a srv through a proper transforma tion. Devroye(Devroye 1986) presents the required transformation techni ques and approximations for a variety of probability distributions. More arbitrary proba bility distributions can be approximated using algebraic manipulations or by series expansions. For uncertainty propagation, Monte Carl o Simulation (MCS) may be the most common choice because of the accuracy and robu stness, but the dilemma of MCS is that the required large number of samples that makes it impractical for computationally demanding models. There are several remedies to reduce the number of samples in MCS, such as importance sampling(Melchers 2001) and separable MCS(Smarslok and Haftka 2006). However, they require special knowledge of the problem or special form of the response. Several computationally efficient methods were proposed in last two decades PAGE 34 20 with reasonable accuracy in many structural problems, such as firstand secondreliability method (FORM/SORM), and response surface method (RSM). The stochastic response surface (SRS) can be viewed as an ex tension of classical de terministic response surfaces for model outputs constructed usi ng uncertain inputs and performance data collected at heuristically selected collo cation points. The polynomial expansion uses Hermite polynomial bases for the space of s quare-integrable probability density function (PDF) and provides a closed form solution of model outputs from a significant lower number of model simulations than those required by conventional methods such as modified Monte Carlo methods and Latin hypercube sampling. In this chapter, a surrogate-based uncer tainty model using stochastic response surface (SRS) is introduced. Reliability anal ysis using Monte Carlo simulation on this surrogate model shows promis ing results in terms of accu racy and efficiency. The proposed method is compared with the firstorder reliability met hod (FORM) and MCS. Description of Uncertainty Model When the inputs of a system are uncertain or described as random variables/parameters, the output or response from this system will have a stochastic behavior as well. Let us assume that thes e random inputs are given in an n-dimensional vector X with continuous joint distribution function() fXx. As shown in Figure 3-1, the system state can have a Boolean description such that the system fails when the limit state()0 G X. The probability of failure Pf can then be defined as a cumulative distribution function (CDF) ove r the failure region, as ()0()f GPfdX Xxx (3.1) PAGE 35 21 Figure 3-1. Limit state function divides the safe region from the failure region Equation (3.1) is called the reliability integr al. Since the integral domain defined by limit state function G ( X ) is complex in the multi-dimensional random space, the reliability integral is difficult to calculate. As introduced in the previous section, by transforming random variables from the original random space to the standard nor mal space, the limit state function can be expressed as a function of a set of srvs { ui} Then, Pf can be expressed in standard Gaussian space as ()0()f GPdU Uuu (3.2) where () is the standard normal PDF and U is the vector of standard random variables. The transformation between X and U is denoted as U =T( X ). In FORM, the probability level of a system is usually represented by the reliability index or safety index For instance, if () is the CDF of the st andard random variable, the failure probability can often be represented by the reliability index 1() f P. Limit state G ( x ) = 0 G ( X )<0 Failure region G ( X )>0 Safe region x1 x2 PAGE 36 22 Stochastic Response Surface Method (SRSM) Polynomial Chaos Expansion (PCE) in Gaussian Space Orthogonal polynomials have many useful properties in the solution of mathematical and physical problems. Just as Fourier series provide a convenient method of expanding a periodic functi on in a series of linearly independent terms, orthogonal polynomials provide a natural way to solve, ex pand, and interpret solutions to many types of important differential equations. Orthogonal polynomials associated with the generalized polynomial chaos (AskeyChaos) are different according to different we ight functions. The type of polynomials is decided by the match between the specific weig ht function and the standard probability density function (PDF). The corresponding type of polynomi als and their associated random variables are listed in Table 3-1. Table 3-1. The type of polynomials and co rresponding random variables for different Askey-Chaos (N 0 denotes a finite integer) Random variable Orthogonal polynomials Support range Gaussian Hermite (, ) Gamma Laguerre [0, ) Beta Jacobi [a,b] Continuous Uniform Legendre [a,b] Poisson Charlier {0,1,2,} Binomial Krawtchouk {0,1,2,,N} Negative Binomial Meixner {0,1,2,} Discrete Hypergeometric Hahn {0,1,2,,N} For example, in Table 3-1, Hermite polynomial chaos expansion requires the weight functions to be Gaussian probabil ity density function, and it satisfies the following orthogonal relation: ()()(), ,kij kkkkkkkij x f xxxdxCij (3.3) PAGE 37 23 where ()kk f x is Gaussian PDF for variable k x ()i kk x is the Hermite polynomial basis, and upper indices i,j denote for two different polynomials. In this research, the uncertainty propag ation is based on stochastic response surfaces (polynomial chaos expansion). The SR S(Isukapalli et al. 1998; Webster et al. 1996) can be view as an extension of cla ssical deterministic response surfaces(Khuri and Cornell 1996; Myers and Montgomery 1995) for model outputs constructed using uncertain inputs and performance data collec ted at heuristically selected collocation points. The polynomial expansion uses Hermite polynomial bases for the space of squareintegrable probability density function (P DF) and provides a closed form solution of model outputs from a significant lower number of model simulations than those required by conventional methods such as modified Monte Carlo methods and Latin hypercube sampling. Let n be the number of random variables and p the order of polynomial. The model output can then be expressed in terms of the srv { ui} as: 0123 111111()(,)(,,)j nnini ppppp iiijijijkijk iijijkGaauauuauuu (3.4) where Gp is the model output, the ,,...pp iijaa are deterministic coefficients to be estimated, and the p(u1,,up) are multidimensional Hermite polynomials of degree p : 1/21/2(,,)(1)TTp p pip ipuuee uu = uuuu (3.5) where u is a vector of p independent and identically di stributed normal random variables selected among the n random variables that represent the model input uncertainties. Equation (3.4) is also called polynomial chao s expansion. The Hermite polynomials (,,)pipuu are set of orthogonal polynomial s with weighting function 2/2 ue, which PAGE 38 24 has the same form with the PDF of standard random variables. In this research, a modified version of Hermite polynomial(Isukapa lli et al. 1998) is used. The first four terms are u, u2 1, u3 3u, and u4 6u2 + 3, when a single rando m variable is involved. The use of the Hermite polynomials has two pur poses: (1) they are us ed to determine the sampling points, and (2) they are used as ba ses for polynomial approximation. In general, the approximation accuracy increases with the order of the polynomial, which should be selected reflecting accuracy needs and computational constraints. The expressions for the 2ndand 3rd-order polynomials in n dimensions (later used in the numerical examples) are: 2nd-order: 1 (2)(2)(2)(2)2(2) 0 111()(1)nnnn iiiiiijij iiijiGaauauauu u (3.6) 3rd-order: (3)(3)(3)(3)2 0 11 1 (3)3(3) 11 21 (3)2(3) 11,1()(1) (3) ()nn iiiii ii nnn iiiiiijij iiji nnnnn ijjijiijkijk ijjiijikjGaauau auuauu auuuauuu u (3.7) The number of unknown coefficients is de termined by dimension of the design space n For 2nd and 3rd order expansion, if the number of unknowns is denoted by N(2), N(3), respectively: (2)(1) 12 2 nn Nn (3.8) (3)3(1)(1)(2) 13 26 nnnnn Nn (3.9) For n = 2, 4, and 8, for example, N(2) = 6, 15, and 45; and N(3) = 10, 35, and 165, respectively. PAGE 39 25 The coefficients in the polynomial chaos expansion are calculated using the least square method, considering samples of input/ output pairs. Since all inputs are represented using srv more accurate estimates for the coeffici ents can be expected, in the sense of statistics, if the probabi lity distribution of the uis is considered. The idea of Gaussian Quadrature of numerical inte gral can be borrowed to gene rate collocation points(Webster et al. 1996). In Gaussian Quadrature, the f unction arguments are given by the roots of the next higher polynomial. Similarly, the roots of the next higher order polynomial are used as the points at which the approximation being solved, which is proposed as the orthogonal collocation method by Villadsen a nd Michelsen (Villadsen and Michelsen 1978). For example, to solve for a three di mensional second order polynomial chaos expansion, the roots of the third order Hermite polynomial, 3 ,0 and 3are used, thus the possible colloca tion points are (0,0,0),( 3 3 3 ),( 3,0, 3),etc.. There are 27 possible collocation points in this case. However, in equation(3.9), there are only 10 unknown coefficients. Similarly, for higher dimensional systems and higher order approximations, the number of available collocation points is always greater than the number of unknowns, which introduces a problem of selecting the appropriate collocation points. For a good approximation in polynomial chaos expansion, the choice of collocation points is critical. Hence, a set of points near the hi gh probability region is heuristically selected among the roots of the one-order higher polynomial under restrictions of symmetry and closeness to the mean. Since the origin always corresponds to the highest probability in standard Gaussian space, the exclusion of the origin as a collocat ion point could potentially lead to a poor PAGE 40 26 estimation. Thus, when the roots of high-orde r polynomial do not include zero, it is added in addition to the standard orthogonal collocation method. The Hermite polynomials (orthogonal with re spect to the Gaussian PDF) provide several attractive features, namely, more robust estimates of the coefficients with respect to those obtained using n on-orthogonal polynomials(Gautsch i 1996); it converges to any process with finite second order mo ments(Cameron and Martin 1947); and the convergence is optimal (exponential) for Gaussian processes(Xiu et al. 2002). In addition, the selected SRS approach in cludes a sampling scheme (col location method) designed to provide a good approximation of the model outpu t (inspired by the Gaussian Quadrature approach) in the higher probability region with limited observations. Once the coefficients are calculated, st atistical properties of the pr ediction, such as mean and variance can be analytic ally obtained, and sens itivity analyses can be readily conducted. Numerical Example of Stoc hastic Response Surface As an illustration of the efficiency and convergence properties of the SRS approach, consider the constr uction of the PDF associated with a simple analytical function represented by: () x y xe (3.10) with x being a normally distributed random variable, as N(0,0.42). Note that in this case the analytical expression of the PDF is known. The SRS for 2ndand 3rd-order polynomials are shown in Eqs.(3.11) and (3.12), respectively. (2)21.08330.43280.0833(1) yuu (3.11) (3)231.08430.43330.0863(1)0.0112(3) y uuuu (3.12) PAGE 41 27 Figure 3-2. PDF of performance function () x y xe In this particular example, the accuracy of the proposed SRS is compared with the analytical solution. Figure 3-2 shows the PDF obtained from MCS applied to the SRS and from the exact solution. A good agreement is observed in the PDF a pproximation, and the root mean square errors decreases with higher order polynomia ls, showing the convergence of the proposed SRS (Table 3-2). Table 3-2. Root mean square error of PDF compared with the exact PDF of performance function y = ex Polynomial order Errors 2 0.03835 3 0.00969 To illustrate the effectiveness of the SR S in the application to the structural problem, consider a torque-arm model in Fi gure 3-3 that is often used in shape optimization(Kim et al. 2003). The locations of boundary curves have uncertainties due PAGE 42 28 to manufacturing tolerances, modeled as probabilistic distributions. Thus, the relative locations of corner points of the boundary cu rves are defined as random variables with xi~N( di, 0.12). The mean values di of these random variables are chosen as design parameters, while the standard deviation re mains constant during the design process. Figure 3-3. Shape design parameters for the torque arm The torque arm model consists of eight ra ndom variables. In or der to show how the SRS is constructed and the PDF of the model output is calculated, we choose the three random parameters ( x2, x6, and x8) that contribute most si gnificantly to the stress performance at points A and B in Figure 3-3. In the deterministic analysis with mean value, the maximum stress of A = 319MPa occurs at loca tion A. The stress limit is established to be max = 800MPa. In the reliability anal ysis the performance function is defined such that G 0 is considered a failure. Thus, the performance function can be defined as G( x ) = maxA( x ). The number of unknown coeffi cients is a function of the dimension n of the random variables. For 2ndand 3rd-order expans ion, the numbers of coefficients, denoted by N2 and N3, are 10 and 20, respectivel y. There are 27 possible collocation points and 10 unknown coefficients in the case of 2nd-order expansion. For robust estimation, the number of collocation points in genera l should be at least twice the number of coefficients. In this particular example, all possible collocation points are Figure 2-3: Shape design parameters for the torque 5066N 2789N x1 x2 x5 x6 x7 x8 x3 x4 Symmetric Design A B PAGE 43 29 used. After coefficients are obtained, MCS with 100,000 samples is used to obtain the PDF. Figure 3-4 shows the PDF associated with G( x ) when different orders of polynomial approximations are used. The PDF obtained from the direct MCS with 100,000 sample points is also plotted. It is clear that the PDF from the 3rd-order is much closer than that of the 2nd-order to the PDF from the MCS. Figure 3-4. PDF of performance function G( x ) torque arm model In order to compare the accuracy of the probability estimation through proposed SRS, let us check the probability of response being larger than 520MPa. In Table 3-3, the probability obtained from MCS is regarded as the reference. The relative error () of failure probability from MCS estimation with sample size of N can be calculated using the following equation: 1 f f P k NP (3.13) where k denotes the confidence level, for confid ence level of 95%, k=1.96, which can be verified from standard normal table. Thus in Table 3-3, number of MCS sample is 100,000, the error in Pf will be less than 5% with 95% confidence. PAGE 44 30 As shown in Table 3-3, it is clear that the SRS provides a convergent probability result as the order increases. With third order SRS, the accuracy of reliability analysis is significantly improved, compared to FORM. Table 3-3. Comparison of probability of G> 520MPa obtained from different uncertainty analysis methods (Full sampling w ithout using local sensitivity) Method FORM 2nd order SRS 3rd order SRS MCS Prob. of G>520MPa 1.875% 2.061% 1.682% 1.566% Relative error* 19.732% 31.609% 7.407% *Relative error: (.)() 100% () probapproxprobMCS probMCS Improving Efficiency of SRS Usin g Local Sensitivity Information In the proposed SRS, the number of sampling points depends on the number of unknown coefficients. Although the proposed met hod is accurate and robust, we have to address the curse of dimensionality: as the number of random variables increases, the number of coefficients rapidly incr eases, as can be seen in Eqs. (3.8) and (3.9). In addition to the efficient collocation method, the number of simulations can be reduced even further when local sensitivity is available. Recently, Isukapalli et al.(Isukapalli et al. 2000) used an automatic differentiation program to calculate the local sensitivity of the model output with respec t to random variables and used them to construct the SRS. Their results showed that local sensitivity can significantly reduce the number of sampling points as more informati on is available. The computational cost of the automatic differentiation, however, is ofte n higher than that of direct analysis(Van Keulen et al. 2004). However, in the applic ation to the structur al analysis, local sensitivity can be obtained at a reasonable computational cost. At each sampling point, the local sensitivity is a part ial derivative of the limit st ate with respect to random PAGE 45 31 variables. Hence, if local sensitivity inform ation is available, then n+1 data at each sampling point can be used for construc ting the proposed SRS, which significantly reduces the required number of sampling points. Continuum-Based Design Sensitivity Analysis In this research, the continuum-based design sensitivity analysis(Choi and Kim 2004a) is utilized to calculate the gradient of the performa nce function with respect to random variables. Even if the idea can be us ed in a broader context, only structural problems are considered in this research. Let z be the displacement and zbe the displacement variation that belongs to the space of kinematically admissible displacements. For given body force f and surface traction force t, the variational equation in the continuum domain is formulated as (,)() al zzz, (3.14) for all z In Eq. (3.14), the structural bilinear and load linear forms are defined, respectively, as (,)()()ijijad zzzz (3.15) ()Tiiiilfzdtzd z (3.16) where ij are components of the engineering strain tensor, and ij are components of the stress tensor. In linear elastic materials, the constitutive relation can be given as ()()ijijklklc zz (3.17) where the constitutive tensor cijkl is constant. The summation rule is used for the repeated indices. In order to solve Eq.(3.14) numerically, the finite-element-based method or the meshfree method can be employed, which ends up solving the following form of matrix equation: PAGE 46 32 []{}{} KDF (3.18) where [K] is the stiffness matrix, {F} the discrete force vector, and {D} the vector of nodal displacements. The major com putational cost in solving Eq.(3.18) is related to L-U factorization of the coefficient matrix. As will be shown later, the efficiency of sensitivity calculation comes from the fact that sensitivit y analysis uses the same coefficient matrix that is already factorized when Eq.(3.18) is solved. In design sensitivity analysis, the variational Eq.(3.14) is differentiated with respect to design variables. In shape design, the design variable does not appear explicitly in the governing equation. Rather, the shape of the dom ain that a structural component occupies is treated as a design variable. Thus, a formal procedure is required to obtain the design sensitivity expression. As shown in Figure 3-5, suppose that the initial structural domain is changed into the perturbed domain in which the parameter controls the shape perturbation amount. By defining the design changing direction to be V(x), the material point at the perturbed design can be denoted as x = x + V(x). The solution z(x) of the structural problem is assumed a differentiable function wi th respect to shape design. The sensitivity of z(x) at x is defined as 0(())() lim zxVxzx z (3.19) Figure 3-5. Variation of a st ructural domain according to the design velocity field V(x) x x Initial domain Perturbed domain V(x) PAGE 47 33 The design sensitivity equation is obtained by taking the material derivative of the variational equation(3.14) The derivative of the structural energy form then becomes 0(,)(,)(,) d aaa dVzzzzzz (3.20) The first term on the right-hand side repres ents an implicit dependence on the design through the state variable, while the second term the structural fictitious load, denotes an explicit dependence on the design velocity V(x), defined as '(,)[()()()()()()]VV VijijijijklklijijacdivdzzzzzzzzV (3.21) where 1 () 2j V ikk ij kjkiz zVV x xxx z (3.22) If the applied load is independent of displacement, i.e., conservative, then '()[] []j i Vijii jj i ijiin jV f lzVzfd xx t zVztVd x z (3.23) is the external fictitious load form, where Vn is the normal component of the design velocity on the boundary, and is the curvature of the boundary. The design sensitivity equation is obtained from Eq. (3.20) to (3.23) as ''(,)()(,)VValazzzzz (3.24) for all z. Note that by substituting into z, the left-hand side of the design sensitivity equation (3.24) takes the same form as that of the response analysis in Eq.(3.14). Thus, the same stiffness matrix [K] can be used for sensitivity analysis and response analysis, with a different right-hand side. PAGE 48 34 Once the sensitivity of the field vector is calculated, the sensitivity of the performance function with respect to design variable ui can be calculated using the chain rule of differentiation, as (;)(;)(;)iidyyy duu zxzxzx z z (3.25) When finite element analysis is used, the sensitivity equa tion can be solved inexpensively because the coefficient matrix is already factorized when solving Eq.(3.14) and the sensitivity equation uses the same coef ficient matrix. The computational cost of sensitivity analysis is usually less than 20% of the original analysis cost. The computational efficiency of the uncertainty propagation approach is critical to RBDO since as previously stated at each de sign cycle the updated version of the PDF for the constraint function (related to model outputs) is required. Constructing SRS Using Local Sensitivity In SRS, the system response can be a pproximated as polynomial expansion when k sampling data are available, the linear regression equation can be written as 2 10 11 2 21 22 2 3311 11 11kNya uu ya uu ya uu y Aa (3.26) The above equation is the standard form for linear regression to solve for unknown coefficients { a }. When the sensitivity information is available, additional information at each sampling point can be used in calcul ating the coefficients. By differentiating Eq.(3.4) with respect to random variable ui and by applying the same regression process in Eq.(3.26), we have iid duu yA a (3.27) PAGE 49 35 Equations (3.26) and (3.27) can be combined to ob tain the following regression equations: 11nnd duu d duu yA yA a yA (3.28) Let 1T ndd dudu yy Yy, T inuu AA BA, Eq. (3.28) can be written as YBa (3.29) Thus, the coefficients of SRS can be obtained using least square regression, such that 1 TT aBBBY (3.30) Note that the sensitivity iu y can be calculated usi ng the transformation of random variables, as i iii x y uxu y (3.31) As introduced in the previous section, the local sensitivity /i x y can be obtained implicitly through Eq. (3.25), where design vari able is represented by ui instead of xi since notation x has been used as space coordina te. Since the transformation between srv and variables with other types of distribution are also math ematically well developed, /ii x u can be obtained explicitly. Therefore, Eq.(3.30) provides an explicit solution for obtaining coefficients of SRS. PAGE 50 36 Numerical Example Torque Arm Model In order to show the effec tiveness of the proposed SRS with local se nsitivity, the same torque arm problem with previous exam ple is used. All conditions are the same as before. By taking advantage of using sensitivit y information to build stochastic response surface, the number of collocation points is reduced significantly. Here for the secondorder polynomial chaos expansion, 7 points ha ve been used, while 31 points for the thirdorder case. At each sampling point, the func tion value and sensitivity information are used to construct the SRS. The PDF obtained from the SRS with sensit ivity is plotted in Figure 3-6 along with that from MCS with 100,000 samples. In the case of 2nd-order, the SRS with sensitivity is more accurate than the SRS without sensitivity (Figures 3-4 & 3-6). In order to calculate the accuracy, the probability of 520 GMPa is calculated using FORM, secondand third-order SRS (Table 3-4). Since no analyt ical solution is available, MCS with 100,000 samples is used as a reference. Both SRS are more accurate than FORM. Figure 3-6. PDF of performance fuction G(x) Torque model at in itial design (SRS with sensitivity) PAGE 51 37 Table 3-4. Comparison of probability of G> 520MPa obtained from different uncertainty analysis methods (reduced sampling using local sensitivity) ; Method FORM 2nd order SRS3rd order SRS MCS Prob. of G>520MPa 1.875% 1.520% 1.545% 1.566% Relative error* 19.732% 2.937% 1.341% *Relative error: (.)() 100% () probapproxprobMCS probMCS Table 3-5 compares the probability of G>520 MPa of second order SRS with/without using local sensitivity with that of MCS, which is regarded as the reference of exact value. With local sensitivity a nd seven sampling points, SRS provides more accurate probabilistic result than that without utilizing loca l sensitivity and twenty-seven sampling points. The accuracy is improved by using local sensitivity while computational cost is reduced. Table 3-5. Comparison of probability of G>520MPa obtained with/without local sensitivity (7/27 sampling points) using 2nd order SRS Method 2nd order SRS using 27 sampling points without sensitivity 2rd order SRS using 7 sampling points with sensitivity MCS (100,000 samples) Prob. of G>520MPa 0.2061% 1.520% 1.566% Relative error* 31.6091% 2.937% *Relative error: (.)() 100% () probapproxprobMCS probMCS Summary In this chapter, a stochastic respons e surface method (SRSM) using polynomial chaos expansion is used in calculating structural reliability. Compared with FORM, which is based on the linear approximation at the most probability point, it provides more accurate result for nonlinear responses. In addition, orthogonal polynomial basis provide PAGE 52 38 a convergent behavior as the order of pol ynomial is increased. A nonlinear function has been used as numerical example to show its accuracy and convergence. Since continuum based sensitivity results were obtained during structure analysis, the computational cost is fu rther reduced by providing gradie nt information while fitting response surface. SRSM has been applied on a structural problem to show its effectiveness. When sensitivity information is provided, numerical results show that even lower number of sampling point can provide more accurate result. PAGE 53 39 CHAPTER 4 RELIABILITY-BASED DESIGN OPTIMIZATION Although statistical methods of uncertainties quantific ation have been studied intensively for decades, traditional determ inistic design optimization still takes no advantage in these scientific advances and compensates uncertainties based on experience; for example, the safety factor. Such an optimization scheme usually yields either unsafe or too conservative design due to the lack of uncertainty quantification. In order to impose existing know ledge of uncertainty to e ngineering design process, reliability-based design optimization (RB DO) methodologies have been proposed and developed(Chandu and Grandi 1995; Chen et al. 1997; Du and Chen 2002b; Enevoldsen and Sorensen 1994; Grandhi and L.P. 1998; Kim et al. 2004b; Kwa k and Lee 1987; Liang et al. 2004; Tu 1999; Tu and Choi 1997; Youn et al. 2003), where the system reliability or probability of failure is used to evalua te the system performance. Compared to the deterministic optimization, RBDO provides margins on reliability by quantifying the uncertainty in the response of structural system due to input uncertainty. General RBDO Model Design optimization has been introdu ced to structural engineering for decades(Arora 2004; Haftka and Gurdal 1991; Vanderplaats 2001). Its methodologies have been well developed mathematically, a nd applications in product development are flourishing. The underlying desi gn philosophy is to reduce the cost by pushing the design to its performance margin. In traditional deterministic design, an optimization problem is formulated as PAGE 54 40 minimize Cost() subjectto (),1,2,, jjallowable LUGGjnp d d ddd (4.1) where ()jG d is the constraint functi on, for example, stress; jallowableG is the corresponding maximum cons traint allowable; and ddenotes the vector of the deterministic design variables. The objective is to minimize the cost while meeting the system constraints. A system design depends on the system design variables. In deterministic optimization, both objective and constr aints only depend on the design vector d which contains all deterministic design variables di. In reliability-based design, design is based on a randomly distributed syst em vector, e.g., denoted by X, which contains random variable Xi. In RBDO, the mean value i or the standard deviation i of the system variable Xi can be used as the design variable In some cases, uncontrollable random variables may contribute to the unc ertainty of the performance. Instead of directly setting constraints on deterministic performance, the RBDO problem(Chandu and Grandi 1995; Enevoldsen and Sorensen 1994; Grandhi and L.P. 1998; Wu and Wang 1996) can generally be defined by setting constraints to be uncertainty measures, such as probability of failure. It is then formulated as ,minimize Cost() subjectto(()0),1,2,, jfj LUPGPjnp d x ddd (4.2) where the cost can be any f unction of the design variable d = [ di]T, ( i =1, 2, n ) and fjP is the prescribed failure probability limit for the jth constraint. PAGE 55 41 Reliability Index Approach (RIA) and Performance Measure Approach (PMA) In the RBDO formulation described in the pr evious section, each prescribed failure probability limit f P is often represented by the reliability target index as 1()tfP. Hence, any probabilistic constraint in Eq. (4.2) can be rewritten using equation as (0)()GtF (4.3) where FG(0)=P(G<0) is the cumulative dist ribution function(CDF) of G at the failed region. Equation(4.3) can also be expressed in another way through inverse transformations 1((0)) s GtF (4.4) where s is traditionally called the reliability index. The expression of probability constraint in Eq. (4.4) leads to the so called reliabili ty index approach (RIA)(Tu and Choi 1997; Tu 2001; Youn 2001). The two forms of constraint described in equations (4.2)and (4.4) are basically the same. In FORM/SORM based RBDO, an inner loop op timization is used to find the most probability point (MPP) in the standard Gaussian space. RIA may cause singularity because s approaches infinity or negative infinity when the failure probability is zero or one. In that case, inner loop optimizer may fail to find the MPP. There is an alternative way to avoid singularity(Tu and Choi 1997) based on a different concept of reliability measure. Fo r any given target probability, a certain level of performance can be reached to meet the reliability requirement. Tu(Tu and Choi 1997) proposed an inverse measure approach cal led performance measure approach (PMA) based on FORM by transforming Eq.(4.4) to *1(())0GtgF (4.5) PAGE 56 42 where g* is named the target probabilistic performance measure. In PMA, Eq.(4.5) is used as probabilistic constraint of RBDO. PM A has been proved to be consistent with RIA in prescribing the probabi listic constraint, but their differences in probabilistic constraint evaluation can be significant (Tu 1999). PM A is more robust in FORM/SORM than RIA based on the fact that RIA may yield singularity; that is, s approaches infinity or negative infinity. In addition, for an inactive probabilistic constraint, PMA is more efficient than RIA. Known as an inverse measure approach PMA can also be implemented on the sampling based uncertainty estimation method. For example, in MCS, a performance measure that meets reliability requirement can be obtained from the order statistics of sampled performance values. Figure 4-1 shows the general numerical procedure of RBDO. The effect and efficiency of inverse measure approach has been investigated for FORM/SORM(Tu 1999; Tu and Choi 1997). In th is research, RIA and PMA as two different philosophies for probability constraint evaluation are also addre ssed for SRS-based RBDO. PAGE 57 43 Figure 4-1. Flow chart for reliab ility-based design optimization Probability Sensitivity Analysis (PSA) Similar to the traditional design sensitivit y, where sensitivity quantifies the effect of deterministic design variable to the stru cture response, probability sensitivity provides the quantitative estimation of the changing of fa ilure probability or re liability with respect to the changes of random parameters, such as means or standard deviations of random design variables. In RBDO, the gradient based optimizer n eeds sensitivity information to carry out optimization. Automatic differentiation using fi nite differentiation l eads to a significantly extra computational cost, especially wh en there are many design variables. In RBDO, if constraints are set with th e probability of failure being less than a certain threshold, the gradient of probability with respect to the ra ndom input is required. PAGE 58 44 In this research, probability sensitivity anal ysis is utilized to calculate the gradient information. First, the probability sensitivity calcu lation in FORM is introduced by taking advantage of structural sensitivity analysis. It can be shown that one can obtain accurate probability sensitivity without extra simulati on cost. Since SRSM shows an advantage for nonlinear response, sampling based probability se nsitivity is also introduced. For inverse measure approach, sensitivity for both FORM and sampling based RBDO can be obtained. Probability Sensitivity Analysis in FORM In first order reliability met hod (FORM), reliability index () can be obtained by following equation 1/2()T**UU (4.6) where U* is the vector of MPP. The derivative of failure probability with respect to the design variables in FORM can then be written as ()() ()fP (4.7) where () is the PDF of the standard random vari able. Thus, the sensitivity of the failure probability is directly related to that of th e reliability index, which can be obtained by 1/2()1 T T*** *UUU U (4.8) For a random variable i PAGE 59 45 1 11 1iii ii iT T TT T*** *** ** *T(X, )T(X, )X U X T(X, )(X, )X UU X T(X, ) U (4.9) Since the reliability index and the most probable point are available from the reliability analysis, the sensitivity can be easily obtained. If the computationally expensive structure analysis code does not come with sensitivity analysis, a finite difference method is widely used to provide gradient information for searching the most probability point (MPP). The computational cost of finite difference method is proportional to the number of design variables. Using design sensitivity analysis, we can avoid the finite difference calculation and pr ovide more accurate gradient information to the line search for MPP. In the finite difference method, the gradient of the limit state in the standard normal space is defined as 0()() ()lim gg g UUUU U U (4.10) Every iteration in line search needs to pert urb each design variab le to evaluate the gradient. If the sensitivity information can be obtained from a struct ural analysis code, there is a more efficient way to obtain th e gradient information for MPP search. The gradient () g U can be computed as 1() ()() gg TU UX U (4.11) where T : XU. PAGE 60 46 The transformation T from original random design sp ace to the standard Gaussian space can usually be obtained explicitly, and the gradient () g *X is provided by design sensitivity analysis. In this section, the torque arm model descri bed in Chapter 3 is used to evaluate the accuracy of the probability sensitivity analys is using FORM. At the initial design, the probabilistic parameters of ei ght random variables are cons idered as design variables. Each random variable is assumed to be norma lly distributed with a mean of zero and a standard deviation 0.1. The sensitivity of reliability index is ca lculated based on Eq.(4.9). Since the transformation T is an explicit function of probabilistic parameters, the sensitivity can easily be calculate d with reliability analysis. Table 4-1 shows the sensitivity results with respect to mean (/i ) and standard deviation (/i ). The accuracy of the sensitivity is compared with that of the finite difference method with 1% perturbation size. Table 4-1. Probability sensitivity with resp ect to random parameters (unit: centimeter) design i i / 100% /i i i i / 100% /i i x1 0.376 0.377 100.26 -0.030 -0.030 100.00 x2 5.243 5.243 100.00 -5.773 -5.775 100.03 x3 0.034 0.034 100.00 -0.000 -0.000 100.00 x4 0.106 0.106 100.00 -0.002 -0.002 100.00 x5 0.055 0.055 100.00 -0.001 -0.001 100.00 x6 -7.244 -7.244 100.00 -11.022 -11.011 99.90 x7 -0.140 -0.140 100.00 -0.004 -0.004 100.00 x8 -4.457 -4.457 100.00 -4.171 -4.173 100.05 In Table 4-1, the first column represents eight random variables that have normal distributions. All random variables are assumed to be independent. Since the mean value and the standard deviation are considered as probabilistic parameters, there are 16 cases PAGE 61 47 in the sensitivity calculations. The second and fifth columns represent the sensitivity results obtained from the analytical deriva tive, while the third and sixth columns are sensitivity results from the finite differ ence method. A very good agreement between the two methods is observed. Table 4-2 shows the computational efficien cy of the proposed analytical sensitivity calculation. The gradient information is provi ded from design sensitivity analysis in MPP search in the standard HL-RF method(Haso fer and Lind 1974; Liu and Kiereghian 1991). The computational savings are about 90% co mpared to the case when only the function values are provided. Once the reliability analysis is finished, the sensitivity of reliability index requires additional 17 f unction evaluations for the fin ite difference method, while only a single analysis is enough for the proposed method because the analytical expression in Eq.(4.9) and (4.11) is used. Table 4-2. Computational efficiency of an alytical method for pr obability sensitivity Finite differential method Analytical method number of analyses in MPP search 90 10 number of analyses in sensitivity calculation 17 1 Total number of analysis 107 11 Probability Sensitivity Analysis Using SRSM In RBDO, the probability of fa ilure can be formulated as ()0()f GPfdXxx (4.12) where ()0 G X is the failure region and () f is the joint probability density function(PDF). By introducing an indication function (()0) IG Xsuch that I=1 if ()0 G X and I=0 otherwise, Eq.(4.12) can be rewritten as PAGE 62 48 (()0)()fPIGfdxxxx (4.13) where X denotes the entire random design space. Since Eq.(4.12) is used as a constraint in RBDO, the sensitivity of Pf is required. The derivative of failure probability can be written as 1()() (()0)(()0)() () () (()0)() ()fP ff I GdIGfd f f IGd f xx ux=T(u)xx xxxxx x x uuu x (4.14) where udenotes the entire standard normal space. Explicit expression of Eq.(4.14) for different distribution types and numerical examples are derived in Appendix A. The accuracy of the sensitivity re sults are also presented in Appendix A for the case of various distribution types. Reliability-Based Design Optimization Using SRSM Although RIA and PMA are theoretically cons istent in prescribing the probability constraint, there are still signi ficant differences in probabilisti c constraint evaluation. The RBDO based on RIA and PMA may lead to eith er different converg ence or efficiency. In this section, an RBDO problem is formulated for the same torque-arm model in Chapter 3 using the concepts of RIA and PMA. A 3rd-order SRS is constructed for uncertainty analysis for both RIA and PMA. When the reliability index is used as a constraint in RBDO, it sometimes experiences numerical difficulty because it can have a value of infinity for very safe design. When SRSM is used in evaluating the probabilistic constraint in RBDO, the problem of singularity can be avoided naturally since the value of failure probability can always be obtained from MCS. The accuracy and convergence of SRSM have been PAGE 63 49 illustrated in the previous chapter. Alt hough SRSM usually requires more performance evaluation compared to FORM, it is still an a ffordable and applicable approach to obtain more accurate results for the highly nonlinear system. RBDO with RIA In this section, the RBDO problem of the torque arm model is solved using RIA. Stochastic response surface is used in uncer tainty analysis to evaluate probability constraints. RBDO fo rmulation of Eq. (4.2) can be used strai ghtforwardly to solve the problem. For the torque-arm problem, the obj ective is to minimize the weight while meeting the requirement of reliability constraint If we define that th e structure fails when stresses in this structure r each yield stress, such that ()()0iiyG xx (4.15) where x is the random input variables, ()i x is stress response for i th constraint, y denotes yield stress. The RBDO problem is then defined as Minimize() subjectto(()0)(),1,,iit LUMass PGiNC d x ddd (4.16) where t is the target reliability index and ( ) is the cumulative density function of srv. During the optimization, a = 3 is used, which corresponds to 99.87% reliability. Since the maximum stress location can shift, the probabilities of failure at eight different locations are chosen as constraints in Eq. (4.16). The constraints can be evaluated using the SRS. The SRS needs to be rec onstructed at each design cycle. Table 4-3 shows the properties of the random variables and the lower and upper bounds of their mean values (design variable s). Note that the de sign variables are the PAGE 64 50 relative change of the corner points and the initial values of all design variables are zero. The lower and upper bounds are chosen such that the topology of the boundary is maintained throughout the whole design process. Table 4-3. Definition of random design variab les and their bounds. The values of design variables at optimum design are shown in the 5th column (unit: centimeter). Random Variables dL d (Initial) dU dopt (optimum) Standard Deviation Distribution type d1 .0 0.0 1.0 -2.5226 0.1 Normal d2 .5 0.0 1.0 -0.4583 0.1 Normal d3 .0 0.0 1.0 -0.9978 0.1 Normal d4 .7 0.0 1.0 -2.4663 0.1 Normal d5 .5 0.0 1.0 -2.1598 0.1 Normal d6 .5 0.0 2.0 1.9274 0.1 Normal d7 .0 0.0 7.0 2.3701 0.1 Normal d8 .5 0.0 1.0 -0.0619 0.1 Normal The design optimization problem is solved using the sequential linear programming technique. The optimization process converges after 14 design cycl es and 27 performance evaluations where the relative convergence cr iterion has been met in two consecutive designs. The optimum values of the design va riables are shown in the 5th column in Table 4-3. Figure 4-2 shows the optimum design and analysis results at the mean values. The major changes occur at design parameters 4, 5,6 and 7. Even if the maximum stress constraint is set to 800MPa, the active constraint at optimum design converges to a lower value so that the variance of the input parameters can be accounted for. In Figure 4-2 the maximum value shows 739MPa, which is the extrapolated nodal stress. The actual element stress at the active constraint is about 618MPa, which is much lower than the extrapolate stress show on the figure. Fi gure 4-3 provides the optimization history of the cost function. As a result of the optimi zation process, the mass of the structure is reduced from 0.878kg to 0.497kg (a reduction of about 43.4 %). PAGE 65 51 Figure 4-2. Optimum design and stress distri bution of the torque arm model with 8 random variables. Figure 4-3. Optimization history of cost function (mass) for the torque arm model with 8 random variables. In order to observe the impact of th e accuracy of the uncertainty propagation procedure at the optimum design, a 3rd-orde r SRS is considered at the optimum design (Figure 4-4). Table 4-4 shows the values of th e reliability indices for the active constraint at the optimum design obtained from MCS with100,000 samples, the proposed SRS approach, and the FORM. The MCS result is used as a reference mark to compare the other two methods, which has about 1.5% e rror in estimating reliability index with confidence level of 95%. The proposed SRS a pproach exhibits a lower error than the PAGE 66 52 corresponding to the FORM and compares very well with the exact resu lt (namely, 3) and that obtained using MCS (0.6% error) Figure 4-4. PDF of the performance function at the optimum for the torque-arm problem Table 4-4. Reliability Index of ac tive constraint at optimal design Reliability Index Error (%) MCS 3.0307 SRS 3.0115 0.633 FORM 2.9532 2.556 RBDO with Inverse Measure As we discussed in section 3.2, enlighten ed by PMA, an inverse measure approach can also be applied in SRS based RBDO. In th is section, the inverse measure approach is applied on reliability based design optimization for the torque arm model. As discussed in section 3.2, the design problem can be formulated as *Minimize() subjectto()0,1,,i LU M ass GuiNC d ddd (4.17) where Gi is the i -th constraint. If the total sample size of Monte-Carlo simulation for SRS is N and allowed maximum proba bility of failure is Pf, then Gi can be found by ordering samples and selecting p th smallest sample. PAGE 67 53 f p NP (4.18) Thus, the evaluation of reliability co nstraints is transformed to find the p th order statistic of sampling. One of the advantages of this approach is that the sensitivity of performance based constraint measure can be obtained directly through structure sensitivity analysis: ***1*()()()(;)iGuGxGxTud x (4.19) For example, if i and 2(,)iiXNormalu ***1* *()()()(;) ()iiii iGuGxGxTud x Gx x (4.20) Summary In this chapter, reliability based desi gn optimization using stochastic response surface is discussed. Procedures for both RI A and PMA are investigated and formulated. A torque arm problem shown in Chapter 3 has be en used to demonstrate the feasibility of RBDO using SRS. Since accurate, sensitivity calculation is im portant to the convergence of gradient based optimizer, probability sensitivity using FORM and SRS is presented. It is shown that probability sensitivity in SRS based sa mpling approach can be also obtained with minimal increase of computational cost. If the SRS is accurate enough, the accuracy of sensitivities obtained also have convergent with respect to the increasing of the sampling size. PAGE 68 54 CHAPTER 5 GLOBAL SENSITIVITY ANALYSIS FOR EFFICIENT RBDO Introduction In industrial applications, a system us ually involves considerable number of random variables. As stated in the previ ous two chapters, the increasing number of random variables boosts the computational cost of reliability anal ysis significantly. Structural reliability analys is which involves a computa tionally demanding model is limited by the relatively high number of require d function analyses for uncertainty. Even if the local sensitivity information descri bed in chapter 3 can reduce the number of required simulations, the dimension of the SRS will still increase according to the number of random variables. Design engineers want to reduce the number of variables based on their contribution to the output performance. However, it is challenging to identify the importance of a random variable in the proce ss of uncertainty propa gation. Those random variables with extremely low contribution to the performance variance can be filtered out to reduce the computational cost of uncertainty propagation. Recent development in statistics intr oduces global sensitivity analysis (GSA)(Saltelli et al. 2000; Sa ltelli et al. 1999; Sobol 1993; Sobol 2001), which studies how the variance in the output of a computational model can be apportioned, qualitatively and quantitatively, to different sources of va riations. In this chapter, global variancebased sensitivity analysis has been applied on st ructural model to illustrate different roles of random variables in uncertainty propagati on. Effort has been made to reduce the dimension of random space in the RBDO process. PAGE 69 55 To reduce the number of simulations require d to construct the SRS even further, unessential random variables are fixed duri ng the construction of the SRS. A random variable is considered unessent ial (and hence it is fixed) if its contribution to the variance of the model output is below a given thresh old. Global sensitivity indices considering only main factors are calculated to quantify the model input contribu tions to the output variability hence establishing which factors influence the m odel prediction the most so that: i) resources can be focused to reduce or account for uncertainty where it is most appropriate, or ii) unessential variables can be fixed without signi ficantly affecting the output variability. The latter application is th e one of interest in the context of this chapter. The RBDO problem in the previous chapter was solved with all random variables. However, some random variables did not sign ificantly contribute to the variance of the stress function. Thus, a large amount of com putational cost can be saved if the random variables whose contribution to the variance of the output is small are considered as deterministic variables at their mean valu es. This section describes how the global sensitivity indices considering only main f actors can be used for deciding unessential random variables during the constructi on of stochastic response surfaces. Sensitivity Analysis As defined by Saltelli(Salte lli et al. 2000), sensitivity analysis studies the relationships between information flowing in and out of the model. In engineering design application, sensitivity usually refers to the derivative of performance measure with respect the input design variable This derivative is used dire ctly in deterministic design as sensitivity information is required by the gradient-based optimizer. In this research, this derivative is called local sensitivity, wh ich has been applied in constructing SRS in PAGE 70 56 the previous chapters. Local sensitivity analys is concentrates on the local impact of the design variables. It is carried out by co mputing partial derivatives of the output performance with respect to th e input variable at the curren t design point. In structural optimization, substantial research has been done on the local sensitivity(Choi and Kim 2004a; Choi and Kim 2004b). Another choice of sensitivity measure expl ores what happens to the performance variance if all design variables are allowe d a finite variation. Global sensitivity techniques apportion the output uncertainty to the uncertainty in the input factors. A couple of techniques have been developed in re cent two decades. In this research, Sobols sensitivity indices(Sobol 1993; Sobol 2001), which is based on the decomposition of function into summands of increasing dimension, will be discussed and applied to the constraint eval uation of RBDO. Variance-Based Global Sens itivity Analysis (GSA) Variance-based methods are the most rigo rous and theoretically sound approaches (Chen et al. 2005; Saltelli et al. 2000; Saltelli et al. 1999; Sobol 1993; Sobol 2001)for global sensitivity calculations. This section describes the fundamentals of the variancebased approach and illustrates how the pol ynomial chaos expansions are particularly suited for this task. The variance based methods: (i) decompos e the model output variance as the sum of partial variances, and then (ii) establish the relative contribution of each random variable (global sensitivity indexes) to th e model output variance. In order to accomplish step (i), the model output is decompose as a linear combination of functions of increasing dimensionality as described by the following expression: PAGE 71 57 01212...12 11()()(,)(,,,)nnn iiiijijijnnn iijifaafxafxxafxxx==>=++++x(5.1) subject to the restriction that the integral of the weighted product of any two different functions is zero. Formally, 1111,,,,11()(,,)(,,)0,for,,,,ssssiiiijjjjsspfxxfxxdiijj = xx (5.2) where p(x) is the joint PDF of input random variable x. If, for example, the weighting function is the uniform distribut ion for the random variables or the Gaussian probability distribution, the functions of interest can be shown to be Legendre and Hermite orthogonal polynomials, respectively. The model output variance can now be calculated using a well-known result in statistics. The result establishes that the variance of the linear combination of random variables ( Ui) can be expressed as: 2 0 111()2(,)nnnn iiiiij iiijiVbbUbVUCOVUU===> +=+ (5.3) Hence, the model output variance can be shown to be: 222 1212 11()()()()nnn iiijijnn iijiVfaVfaVfaVf==>=+++ (5.4) where the terms represent par tial variances and each V( ) may be found by definition as: () []2()()()() VffxEfxpxdx= (5.5) In the above formula, f ( x ) represents the function under consideration and the symbol E( ) denotes expected value. There are no covariance terms in Eq.(5.4) because of the orthogonality property shown in Eq.(5.2). The global sensitivity index Si that considers only main factor is called main sensitivity index, which associated with each of the random variables which is represented by Eq. (5.4): PAGE 72 58 2(()) ()iii iaVfx S Vf = (5.6) A sensitivity index that considers the inter action of two or more factors is called interaction sensitivity index. Thus, as denoted by Chan and Saltelli(Saltelli et al. 2000), the summation of all sensitivity indices, involvi ng both main and interaction effect of i-th random variable, is called to tal sensitivity index. Sobol (Sobol 2001) suggested to use total sensitivity indices to fix unessential variab les. If total sensitivity index for certainty variable extremely small compare to 1, that means the contribution of the variable is neglectable and the variable can be considered as a deterministic one. Global Sensitivity Analysis Usin g Polynomial Chaos Expansion The polynomial chaos expansion is part icularly suited for computing global sensitivity indices because: (1) The model output is already decomposed as a sum of functions of increasing dimensionality; (2) the functions are ort hogonal with respect to the Gaussian measure (Hermite polynomials) ; and (3) the variance of the bases are analytically available. For example, the vari ances of the functions associated with a two dimensional chaos expansion of order 2 are shown in Table 5-1. Table 5-1. Variances of the Hermite bases up to the second order In addition, given the polynomial chaos e xpansion (i.e., the coefficients of the linear combination of Hermite polynomial s), the model output variance and global sensitivity indices can be easily computed using Eqs. (5.4) and (5.6), respectively. In both Function f V(f) X1 1 X2 1 X1 2 1 2 X1X2 1 X2 2 1 2 PAGE 73 59 equations, the variances V( ), are readily available for pol ynomial chaos expansions of arbitrary order and number of variables. In order to show the effect of the global sensitivity, let us consider the torque arm model presented in chapter 3. The eight ra ndom design variables are normally distributed with standard deviation of 0.1. A cubic stoc hastic response surface with eight variables in standard normal space has been constructed to approximate the stress response. Global sensitivity indices for main factor have been calculated to quantif y the contribution of each random variable to performance variab ility. Figure 5-1 shows that three design variables(x2,x6,x8) makes most contribution to the tota l variance, influences from other variables are extremely small. That means, at the initial design, that the randomness of those variables with low global sensitivity indices can be ignorable. If only three variables are considered in constructing SRS instead of eight, the computational cost will be saved significantly, while maintaining the same level of variability. Figure 5-1. Global sensitivity indices fo r torque arm model at initial design Adaptive Reduction of Random De sign Space Using GSA in RBDO The idea of adaptive reduction of random variables is based on the main factor of each random variable. If it is smaller than a thre shold, it is fixed in constructing SRS. For PAGE 74 60 that purpose, a linear polynomial chaos expans ion is enough to obtain the main factors of GSI. The current algorithm with linear pol ynomial chaos expansion for the reduction of random variables can be modified to use a sens itivity index that accounts for interactions. These interactions will only appear in higher order polynomial chaos expansions. The choice of a non-linear polynomial chaos e xpansion would reduce the computational efficiency of the proposed approach with unclear significant advantages. As stated at the beginning of this sect ion, once the global sensitivity indices are calculated, variables that have the least influence on the model prediction (unessential variables) can be identified and eventually held fixed without signi ficantly affecting the output variability. The procedure is adaptive because the global sensitivity indices are calculated at each design iteration and as a re sult different sets of random variables may be fixed throughout the RBDO process. A flow chart of RBDO using this strategy is shown in Figure 5-2. PAGE 75 61 Figure 5-2. Adaptive reduction of unessent ial random design variables using global sensitivity indices in RBDO. Low-order SRS is used for global sensitivity analysis, while a high-order SRS is used to evaluate the reliability of the system. Low-rder Stochastic Response Surface Reliability-based Design Optimization Min. s.t.(()0)ifi LUCost PGP X XXXX1 X2 Xn Random Design Variables Optimized? No Stop Yes Global Sensitivity X1 X2 X3 X4 Xn Fixing Unessential Variables Update Design Higher-Order Stochastic Response Surface PAGE 76 62 At the initial design stage, a lower-order stochastic response surface is constructed using all random variables. In this particul ar example the first-or der SRS is constructed using 17 sampling points. At the initial desi gn, the first-order SRS with eight random variables can be expressed as, 1 01122334455667788 1234 56784.950.00630.1170.000080.0019 +0.00260.0520.00020.016 Gaauauauauauauauau uuuu uuuu =++++++++ =+++ (5.7) One useful aspect of the polynomial chaos e xpansion is that the coefficients in Eq. (5.7) are a measure of the contribution of the corresponding random variable to the variation of the output, and these coefficients will not change significantly in higherorder SRS. On the other hand, typically the global sensitivity inde x associated with a particular variable is responsible for most of its contri bution to the outpu t variance. Thus, evaluating the global sensitivity indices using th e first-order SRS can be justified. Since all random variables are transformed into st andard normal random variables, the variance of G1 can be evaluated analytically. Using Eq. (5.6) and (5.7) and assuming the design variables are independent, the global se nsitivity index can be calculated as: 2 2 1 i i n j ja S a (5.8) Note that the denominator in Eq. (5.8) is the total variance of G1 using the firstorder approximation. Thus, the global sensitivity index, Si, is the ratio of the contribution of i-th random variable to the total variance. If the global sensitivity index of a specific variable is less than a threshold value, the variable is considered as deterministic and fixed at its mean value. PAGE 77 63 In order to show the advantage of fixi ng unessential random variables, the global sensitivity indices of the torque-arm model are calculated. Table 5-2 shows the global sensitivity indices of the torque-arm model usin g the first-order SRS at the initial design. The total variance of stress function is 1.670 2. Based on the global sensitivity indices, there are only thr ee random variables whose contri bution is greater than 1.0%; i.e., u2, u6, and u8. Thus in the reliability analysis only these three random variables are used in constructing the third-order SRS, which now requires only 19 sampling points for 10 unknown coefficients. All other random vari ables are considered as deterministic variables at their mean values. If the total number of sampling points for both lower(17) and higher-order (19) polynomial expans ions are compared with the higher-order SRS using all random variables (89), a signi ficant reduction of the number of sampling points was achieved. Table 5-2. Global sensitivity indices consid ering only main factors for the torque arm model at the initial design. Only three random variables ( u2, u6, and u8) are preserved when a threshold value of 1.0% is in place. SRV Variance GSI (%) u1 3.916 5 0.235 u2 1.369 2 82.0 u3 6.403 9 0.00003834 u4 3.667 6 0.02197 u5 6.864 6 0.04109 u6 2.702 3 16.179 u7 4.818 8 0.0002885 u8 2.538 4 1.519 The RBDO problem, defined in Eq.(5.2) in chapter 4 is now solved using the proposed adaptive reduction of random variables. The optim ization algorithm converges after the 17-th iteration. As shown in Figure 5-3, the optimum design using the adaptively reduced SRS is slightly different from that with all random variable s in chapter 4. The PAGE 78 64 former has a longer interior cutout than the latt er. This can be explained from the fact that the model with reduced random variables ha s less variability th an the full model. Furthermore, the optimum value achieved using the adaptively reduced SRS converges to a lower value than the one w ithout adaptive reduction). The total mass of the torque arm is reduced by 54.8%. The difference between th e two approaches is approximately 1.8%. The number of active random variables asso ciated with the modeling of the first constraint during the design iter ations are listed in Table 53. On average, four random variables were preserved, which implies th at only 29 sampling points were required for constructing the SRS. This is three times le ss than the SRS approach without adaptive reduction (89 sampling points). Figure 5-3. Optimum designs for the full SR S (solid line) and adaptively reduced SRS (dotted line). Because some variables ar e fixed, the interior cutout of the design from the adaptively reduced SRS is larger than that from the full SRS. Table 5-3. Comparison of the number of ra ndom variables in each design cycle. The threshold of 1.0% is used. Th e first constraint is listed. Iter Full SRS Reduced SRS 1 8 3 2 8 3 3 8 3 4 8 3 5 8 4 6 8 4 7 8 5 8 8 4 17 8 4 PAGE 79 65 Summary In this chapter, a dimension reduction tec hnique using global sens itivity indices is introduced. Since the varian ces of the Hermite polynomia l bases are analytically available, the SRS is suitable to compute global sensitivity indices Its application to RBDO is also presented. In the RBDO procedur e, the global sensitivity indices that are calculated using the lower-order SRS are us ed to fix unessential random variables, a higher-order SRS with reduced dimension is then used in evalua ting the probability constraint. Fixing the unessential random variables accelerates the design optimization process. The RBDO result obtained in this wa y is compared to that from previous chapter, which shows little difference because of the loss of variability in fixing random variables. PAGE 80 66 CHAPTER 6 FATIGUE RELIABILITY-BASED LOAD TOLERANCE DESIGN Introduction Traditional reliability-based structur al design usually makes assumptions on randomness of factors involved in modeling a stru ctural system such as design variables, material properties, etc. These parameters are relatively well controlled so that the variability is usually small. However, it is al so important to consider the capacity of the system subject to working conditions, e.g., un certain loadings, because the uncertainty in load or force is much larger than that of others. The variability of the load is often ignored in the design stage and is difficulty to quantify it. Without knowing the accurate uncertainty characteristics of input load, it is hard to rely on the reliability of the output. In this chapter, a different approach fr om the traditional RBDO is taken by asking how much load a system can support. The amount of load, which a structural system can support, becomes an important info rmation for evaluating a design. As an illustration, the fatigue reliability-ba sed load tolerance of the front loader frame of CAT 994D wheel loader is studied. Besides the uncertainty in the material properties, which can be incorporated in S-N curve(Ayyub et al. 2002; Chopra and Shack), uncertainties are also investigated on both mean and amplitude of a given dynamic load. Either the variation of load amp litude or mean may affect the fatigue life of the structural system. This research pres ents a reliability based load design method, which provides the load envelope for a structur e subject to fatigue failure mode. Both one dimensional and multidimensional problems are addressed. PAGE 81 67 Since service loads are subjective such that the load characteristic of one operator may completely different from that of others. In order to perform reliability analysis, it is necessary to know uncertainty characteristics of inputs. Ho wever, distribution type and parameters of loads are often unknown. In this chapter, instead of m odeling variability in parameters by assuming specific type of ra ndom distribution, the effect of different distribution types on the system response is investigated by introducing the concept of conservative distribution type which provides a safer wa y to model uncertainties. Fatigue Life Prediction Recent developments in the computer-a ided analysis provide a reasonable simulation for fatigue life prediction at ea rly design stage for components under complex dynamic loads. For most automotive components, fatigue analysis means to find crack initiation fatigue life. Figure 6-1 illustrates the procedure to the crack initiation fatigue life prediction. Figure 6-1. Flow chart fo r fatigue life prediction PAGE 82 68 Crack Initiation Fatigue Life Prediction Two major crack initiation fatigue life prediction methods are stress-based and strain-based methods The stress-life (SN ) approach employs relationship between the stress amplitude and the fatigue life. This method is based primarily on linear elastic stress analysis. The advantage of stress-life approach is apparent since changes in material and geometry can easily be evaluated and large empirical database for steel with standard notch shape is available. However, th e effects of plasticity are not considered in this method. The local strain-life(N ) method assumes that the local strains control the fatigue behavior. The plastic eff ects are considered we ll in this method. It is similar to the stress-life approach in that it uses N curve instead of SN curve, but differs in that the strain is the variable relate d to the life, and also in that plastic deformation effects are specifically considered. Machine parts are usually required to be durable and able to undertake high numbers of life cycles. The front loader fram e of CAT 994D wheel loader is one of such case. The critical position of fatigue failure is usually at welding joints. Because the stress-life method is works well for the br ittle material and provides a reasonable approximation for a high cycle fatigue crack in itiation life, by taking advantage of the availability of a large amount of available uniaxial fatigue data, stress-life method is employed in this chapter. Since the crack is usually initiated al ong the component surface, for saving unnecessary computation, FE based fatigue analysis chooses element along surface to calculate the fatigue life. For multi-axial app lication, the principal stress method has been applied using the planes perpendicular to the surface. Fatigue lives are calculated on eighteen planes spaced at 10 degree incremen t. On each plane the principal stresses are PAGE 83 69 used to calculate the time history of the stre ss normal to the plane. It has been shown that this method should only be used for fatigue an alysis of brittle me tals like cast iron and very high strength steels, as it provides nonconser vative results for most ductile metals. Based on the factor that material in our application is cast iron and the interested region is welding joints, the principal stress algorithm can offer the fatigue life calculation with reasonable accuracy. Using superposition of dynamic loadings and the quasi-static FE analysis, the dynamic stresses in the component are used to analyze multi-axial fatigue, based on principal stress using conventional S-N curve(Fe-safe 2004). Most basic fatigue data are collected in the laboratory by testing procedures which employ fully reversed loading. However, realistic service load ing usually involves nonzero mean stresses. Therefore, the influe nce of mean stress on fatigue life should be considered so that the fully reversed laborator y data can be used in the evaluation of real service life. Since the tests required to determine the influence of mean stress are quite expensive, several empirical relationships wh ich related alternating stress amplitude to mean stress have been developed. Among th e proposed relationships, two are widely used, which are Goodman and Gerber models. Goodman: (/)(/)1aemuSSSS Gerber: 2(/)(/)1aemuSSSS where aS : Alternating stress amplitude; eS : Endurance stress limit mS : Mean stress; uS : Ultimate strength PAGE 84 70 Experience shows that test data tend to fall between the Goodman and Gerber curves. In the application of fatigue life pr ediction of front loader frame of CAT 994D, Goodman relation is applied to a ddress the mean stress effect. Variable Amplitude Loading and Cumulative Damage In real application, components are us ually subject to complex dynamic loading which has variable amplitude. It requires iden tifying cycles and asse ssing fatigue life for each cycle. The rain-flow counting met hod(Matsuishi and Endo 1968) is the most commonly used cycle counting technique. This method defined cycles as closed stressstrain hysteretic loops as shown in the figure below: Stress 0 A B C D E F G H I(=A) Time B D F H C A=I E G Figure 6-2. Rain-flow and hysteresis PAGE 85 71 An algorithm of rain-flow counting can be developed based on ASTM standard description. Although the rain -flow counting method is not based on an exact physical concept to account for fatigue damage accumu lation, it is expected to provide a more realistic representation of the loading history. Cumulative damage of each cycle ca n be obtained by the Palmgren-Miner hypothesis, which is referred to as the linear damage rule: 1 n j i j jn D N (6.1) in i-th cycle. In Eq. (6.1), Di is the fraction of the damage; nj is the counted number of cycles for j-th stress range; Nj is the cycles to fail; and n is the total number of stress ranges counted from rain flow. Fa ilure is predicted to occur if 11N i iD (6.2) where N is the number of cycles. Thus, the fatigue life can be calculated as the number of applied load cycle until the cumulative damage reaches 1: 11 N i iLifecycles D (6.3) Model Preparation for Fatigue Reliability Analysis Finite Element Model Figure 6-3 shows the component of a fr ont loader frame of CAT 994D wheel loader, which is subjected to 26 channels of dynamic loading. As show in Figure 6-4, the finite element model consists of 49,313 gr id points and 172,533 elements (24 beam, 280 gap, 952 hexagon, 1016 pentagon, 226 quadr ilateral, 160,688 tetrahedron, 9,144 triangular, 203 rigid body elements). In or der to apply for the displacement boundary PAGE 86 72 conditions and loads, pins are modeled usi ng beam and gap elements. The existence of gap elements makes the problem nonlinear. However, if the gap status does not change during analysis, we can still consider the problem to be linear. Figure 6-3. Front loader frame of CAT 994D wheel loader (subject to 26 channels of dynamic loading) 1. APin_Lf_X 2. APin_Lf_Y 3. APin_Lf_Z_Pos 4. APin_Lf_Z_Neg 5. APin_Rt_X 6. APin_Rt_Y 7. APin_Rt_Z_Pos 8. APin_Rt_Z_Neg 9. GPin_Lf_X 10. GPin_Lf_Y 11. GPin_Rt_X 12. GPin_Rt_Y 13. YPin_Lf_X 14. YPin_Lf_Y 15. YPin_Rt_X 16. YPin_Rt_Y 17. UHitch_X 18. UHitch_Z 19. LHitch_X 20. LHitch_Z 21. Steer_Rt_X 22. Steer_Rt_Z 23. Steer_Lf_X 24. Steer_Lf_Z 25. LHitch_Y_Pos 26. LHitch_Y_Neg Steer Cylinder Pin PAGE 87 73 Figure 6-4. Finite elemen t model for front frame Dynamic Load History In Figure 6-3, a total of 26 degrees-of-fr eedom are chosen for the application of dynamic loads. All loads are located in the center of the pins and the hitches. Even if the dynamic load f ( t ) is applied to the system, it is assume d that the inertia is relatively small and the method of superposition can be applie d. Thus, only a linear static analysis is enough with the unit load applied to each degr ee-of-freedom. The stress value from the unit load is called the stress influence coe fficient. The dynamic st ress can be obtained by multiplying these stress influence coefficients with the dynamic load history. Measured data of 26 channels are used for the dynamic loads with the durati on of 46 minutes. This duration is defined as a working cycle. The dynamic load is sampled such that 9,383 data points are available for each channel. Clamped Clamped PAGE 88 74 Uncertainty in Material Properties and S-N Curve Interpolation Based on available material properties and the components working conditions, principle stress analysis usi ng the Goodman model is used as the algorithm for fatigue life prediction. From superpos ition of quasi-static linear finite element analysis and dynamic loading, the stress data are obtaine d for each element. These stresses can be regard as true stress, which means S-N curve can be applied di rectly on principle stress life method without considering the stress concentrate factor. The S-N curve can be interpolated from nominal stress-life data. Considering the uncertainties of material properties, this interpolation will be implemented in a random manner. A lognormal distribution in S-N curve can be assumed to simplify the randomness. Although there is no rigorous statistical evalua tion was performed, but this assumption seems reasonable empirically(Ayyub et al. 2002; Chopra and Shack). Figure 6-5 shows the S-N curve obtained from stresslif e data for cast iron used in the front frame. Solid line is the nominal S-N curve and two dashed lines represent the variation in S-N curve. Figure 6-5. Material S-N curve with uncertainty PAGE 89 75 Uncertainty Modeling of Dynamic Loadings Dynamic loadings are usually very co mplicated and may involve a lot of uncertainties. Figure 6-6 shows one cha nnel of the dynamic load. The mean and amplitude of dynamic loading usually plays the most important ro le in fatigue life estimation. Therefore, for the purpose of illustration and simplification, uncertainties can be model based on these two quantities. Figure 6-6. Illustration of one channel of dynamic loads By combining the effects of the randomness of mean and amplitude of the loads, two load capacity coefficients and are defined for the mean and amplitude, respectively. The dynamic load can be parameterized as 0()(())meanmeanftfftf (6.4) where are random parameters to describe the uncertainties of the loads called load capacity coefficient (LCC). In Eq. (6.4), 0() f t is the original dynamic loads and mean f is the mean value of the initial loads. Due to the random parameters, the dynamic load () f t shows probabilistic behavior. Equation (6.4) provides a simple two dimensional model of PAGE 90 76 uncertainty in dynamic load hi story. Note that when both and equal to one and fixed, the original deterministic loading history0() f t is recovered. In the following section, one-dimensional problem will first be investigated by fixing one of them. For example, if we fix at 1 and treat as a random variable, then uncertainty is modeled for the amplitude of the load. Linear Estimation of Load Tolerance The major challenge of the research is to es timate the load tolerance with respect to the reliability of fatigue life performan ce, which depends on the load history and uncertainty characterization. Id entifying the load distribution is one of the most difficult tasks in the uncertainty analysis becaus e different operating conditions will yield completely different distribution types. At this point, it is assumed that the load has a specific uncertainty character istic (distribution type a nd corresponding parameters). When the variance of the load is fixed, for ex ample, it is possible to construct the safety envelope by gradually changing the mean va lue of the applied load, which requires a large number of reliability analyses. When nonlinearity of the system is small, it is possible to estimate the safety envelope usi ng the sensitivity information at the current load without requiring further reliability analyses. This estimation is based on the firstorder Taylor series expansion method. Fo r illustration, one-dimensional models (only considering single random variab le) for the effect of amplit ude and mean are separately investigated to meet the reliability requirem ent of fatigue life under uncertainty. In these one dimensional cases, the varia tion in S-N curve is ignored. PAGE 91 77 Variability of Dynamic Load Amplitude In order to consider the va riability of the dynamic load, the mean value of the load is first assumed to remain constant, while th e amplitude of the load is varied randomly. From Eq. (6.4), the uncertainty caused by the amplitude can be represented using the following decomposition of the dynamic load: 0mean()(())meanftfftf (6.5) When = 1, the original load history is recovered. When = 0, the dynamic load becomes a static load with the mean value. In this definition, cannot take a negative value. Since is a random variable, it is necessary to assume the distribution type and distribution parameters for First we assume that is normally distributed with the mean of one and the standard deviation of 0.25 (COV=0.25). The standard deviation is estimated from the initial dynamic load history. Since the first-order re liability analysis is performed using the standard random variable, we convert into the standard random variable u by 10.25 u u (6.6) where u ~ N(0,12), ~ N(1,0.252), = mean, = standard deviation. For any given sample point u corresponding can be obtained from Eq. (6.6), and using a new dynamic load history can be obtained from Eq. (6.5). By applying this dynamic load history, we can calculate the fatig ue life of the system. Since this procedure includes multiple steps, we can construct a (stochastic) response surface for the fatigue life and then perform the reliability an alysis using the response surface. PAGE 92 78 Since the fatigue life changes in several or ders of magnitudes, it would be better to construct the response surface for the logarith mic fatigue life. In one dimensional case, five collocation points are available from th e DOE introduced in chapter 3. Using these five sampling points, a cubic stochastic re sponse surface is constr ucted as a surrogate model for the logarithmic fatigue life as 23 10()()5.70750.72230.0581(1)0.0756(3) LLogLifeuuuu (6.7) Note that in one dimensional case, the five collocation points available from previous DOE scheme are sometimes not enough to construct a high fidelity SRS, a couple of complementary sampling points can be chosen to construct a new SRS, which spread evenly between the original collocation points, i.e., four more point in the middle of intervals of the original five point s have been chosen. The corresponding SRS logarithmic fatigue life becomes 23 10()()5.69760.68260.0541(1)0.0617(3) LLogLifeuuuu (6.8) Various quantities for estimating the quali ty of SRS are shown in Table 6-1 for both five and nine sampling points scheme. Th e table shows nine points DOE schemes fit the data better based on significant improveme nt of PRESS (prediction error sum of squares). Table 6-2 lists the tstatistics for the evaluation of each coefficient in the above response surface. Although using more samp ling points can increase the fidelity of estimation, it also increases the computationa l cost. In our specific problem, considering the saving of computation, the five sampling DOE scheme is sufficient. Table 6-1. Quality of response surface Error statistics RMSE SSe R2 R2 adj PRESS 5 sampling DOE 0.0406 0.0082 0.9980 0.9921 8.5671 9 sampling DOE 0.0511 0.0235 0.9965 0.9943 0.1503 PAGE 93 79 Table 6-2. T-statistic of the coefficients coefficient 1 2 3 4 t-statistics of 5 sampling DOE 122.7590 15.9279 3.5759 4.0828 t-statistics of 9 sampling DOE 229.1431 33.9519 4.9313 6.4894 The response surface in Eq. (6.7) shows that the mean of logarithmic fatigue life is 5.6976 and the standard deviation is about 0.6826. It also sh ows that the c ontribution of the higher-order terms is relatively small, co mpared with the constant and linear terms. Thus, we can conclude that the performan ce is mildly nonlinear with respect to the random variable. Since the required life of the working component is 60,000 hours and each cycle corresponding to 46 minutes, th e target of the fatigue life can be written in the logarithmic scale by target10 10(60,000hours) (78,261cycles) 4.9 LLog Log (6.9) The system is considered to be failure when the predicted life from Eq. (6.7) is less than the target life in Eq. (6.9). Accordingly, we can define the probability of failure as targettarget[()0]fPLLP (6.10) where Ptarget is the target probability of failure. For example, when Ptarget = 0.1, the probability of failure should be less than 10%. Even though the interpretation of Eq. (6.10) is clear, it is often inconvenient becau se the probability changes in several orders of magnitudes. In reliability an alysis, it is common to use the reliability index, which uses the notion of the standard random variable. Equation (6.10) can be rewritten in terms of the reliability index as PAGE 94 80 targettarget()()fPP (6.11) where is called the reliability index and is the cumulative distribution function of the standard random variable. When Ptarget = 0.1, target 1.3. The advantage of using the reliability index will be clear in the following numerical results. With the response surface in Eq. (6.7), reliability analysis is carried out using the first-order reliability method (FORM) at = 1. The results of reliability analysis are as follows: 17.81% 0.922456 3.972fP (6.12) where / is the sensitivity of the reliability in dex with respect to the mean value of Since Ptarget = 0.1 and target = 1.3, the current system does not satisfy the reliability requirement. From the mild nonlinear property of the re sponse, we can estimate the mean value of that can satisfy the required reliability. The linear approximation of the mean value can be obtained from 1 1arg1()/0.9049estimatetet (6.13) which means that the mean value of needs to be decreased about 10% from the original load amplitude in order to satisfy the required reliability. In order to verify the accuracy of the estimated result, several sampling points are taken and reliability analyses are performed. Figure 6-7 show s the reliability index with respect to while Figure 6-8 shows the probability of failure Pf with respect to The solid line is linearly approximated reliability using sensitivity information. The reliability PAGE 95 81 index is almost linear and the estimation using sensitivity is close to the actual reliability index. When the target probabi lity of failure is 0.1 and has the dist ribution of N(,0.252), the safety envelope can be defined as 00.9049 (6.14) Figure 6-7. Reliability index with respect to random parameter Figure 6-8. Probability of failure Pf with respect to random parameter The result means that current design, c onsidering 25% standard deviation in the load amplitude, is not enough to achieve 90% reliability. The structure should work under Pf Ptarget target PAGE 96 82 milder working condition, which means either lower the mean of the load amplitude by about 10% or provide more accurate estimation of the initial load to reduce the variance. Variability of Mean of Dynamic Load Since both mean and amplitude are used to describe the dynamic load history, both of their effects under uncertainty are studied se parately. When the mean value of the load is assumed to be varied randomly and the load amplitude remains as the initial load, the uncertainty in load can be modeled as: 0()(())meanmeanftfftf (6.15) Same as the case of load amplitude, when = 1, the applied load is identical to the original load history. When = 0, the applied load has the same amplitude with the original load history bu t the mean value is zero. In this definition, can be a negative value. Since is a random variable, it is necessary to assume the distribution type and distribution parameters for First we assume that is normally distributed with the mean of one and the standard deviation of 0.25 (COV=0.25). Since the first-order reliability analysis is performed using the standard random variable, we convert into the standard random variable u by 10.25 u u (6.16) where u ~ N(0,12), ~ N(1,0.252), = mean, = standard deviation. By following the same procedure with pr evious section, using nine sampling points, we can construct a cubic stochastic response surface as a surrogate model for the logarithmic fatigue life as 23 10()()5.69060.09050.0013(1)0.0003(3) LLogLifeuuuu. (6.17) PAGE 97 83 If we compare the response surface in Eq. (6.17) with the case of amplitude change in Eq.(6.7), the mean values of the both cases are close but the standard deviations are quite different. From this result, we can conclude that the variance of the mean value does not contribute significantly to the variance of the fatigue life. Using the response surface in Eq. (6.17), reliability analysis is carried out using the first-order reliability method (FORM) at = 1. The results of reliability analysis are as follows: 810 6.3435 4.012fP (6.18) Since Ptarget = 0.1 and target = 1.3, the current system satisfies the reliability requirement. The linear approximation of the mean value can be obtained from 1 1target1()/2.257estimate (6.19) which means that the system satisfies the reli ability requirement even if the mean value of is increased up to 225% from the original load. This observation is consistent with the conventional notion of fatigue analysis in which the effect of the amplitude is significant while that of the mean is not. In order to verify the accuracy of the estimated result, several sampling points are taken and reliability analyses are performed. Figure 6-9 show s the reliability index with respect to while Figure 6-10 shows th e probability of failure Pf with respect to The solid line is linearly approximated reli ability using sensitivity information. The reliability index is almost linear and the esti mation using sensitivity is close to the actual PAGE 98 84 reliability index. When the target probability of failure is 0.1 and has the dist ribution of N(,0.252), the safety envelope can be defined as 02.257 (6.20) Figure 6-9. Reliability index with respect to random parameter Figure 6-10. Probability of failure Pf with respect to random parameter Safety Envelope Concept for Load Tolerance Design The safety envelope is defined as the magnitudes of the input design variables when the system fails. When design variables ar e loads, it is called the load envelope. In one dimensional case, this is simply the range of the allowed loads, e.g., the range of Pf Ptarget target PAGE 99 85 mean value of or in the previous section. In mu lti-dimensional case, the combination of various loads constitutes an envelope, which is convex in linear systems. Such information is very useful as a capacity of the current design, a future reference for design upgrade, maintenance and control. Figure 6-11 shows a schematic illustration of the safety envelope when two variables are involved. In such a complex situation, a systematic way of searching the boundary of th e safety envelope needs to be developed. Figure 6-11. Safety envelope for two variables When the relationship between the safety of the system and the applied loads is linear or mildly nonlinear, linear approxima tion can produce a very effective way of estimating the safety envelope. In context of reliability based safety measure, the target of safety envelope is that failure probability cannot reach over the prescribed value. Numerical Path Following Algorithm According to the reliability based safety envelope concept intr oduced above, when target reliability has been specified, a safety envelope can be constructed using numerical path following algorithm to search the boundary of the safety envelope(Allgower and Safety envelope Safe PAGE 100 86 Georg 1990). In this research, a systematic way of searching the boundary of the safety envelope is proposed using a predictorcorrect or method, which is similar to the Euler Newton continuation method(Allgower a nd Georg 1990; Kwak and Kim 2002). When the relationship between the safety of the st ructure and the applied loads is linear or mildly nonlinear, this approach can produce an efficient way of estimating the safety envelope. In the context of reliabilitybased safety measur e, the boundary of the safety envelope is the location where the probability of failure is equal to the target probability. Figure 6-12. Predictor-corrector algorithm The predictorcorrector algorithm(Figure 6-12) is explained below with two random variables, and First, the distribution type of random variables is assumed. The effect of different distri bution types on the safety envelo pe can be investigated by following the same procedure as in the prev ious section. It is clear that the two parameters must have non-ne gative values, which means th at the safety envelope only exists in the first quadrant. The capac ity of the structure with respect to (, ) is interesting. If the required probability of failure is Ptarget (i.e., target (Ptarget)), the following steps can been taken to construct the safety envelope: PAGE 101 87 Step 1: Set k =1. Set the move limit l and a small parameter Initialize (, )=( ). Step 2: Find an initial state ( 1) such that ( 1)=target. Step 3: Determine the tr ial state (predictor). The trial state can be obtained by movi ng in the tangent direction of the boundary of the envelope. From the firs t-order Taylor series expansion of ( k, k)=constant and from the move limit of l the following two conditions can be used to determine the trial increments: 22()()trtrl (6.21) 0kktrtr (6.22) Of the two possible direc tions, the one that provides a clockwise (or counter clockwise) direction is selected. Then the trial state can be obtained by trktr trktr (6.23) According to the convex property of th e envelope, the trial state in Eq.(6.23) can be either inside or outside the envelope. The reliability index at the trial state is tr=( tr, tr). Step 4: Return to the boundary of the envelope (corrector). Since the trial state is not on the bounda ry, it needs to be returned to the boundary of the envelope. The correcting di rection is perpendicular to the trial direction. targettrktrcrcr (6.24) PAGE 102 88 0crtrcrtr (6.25) Then, the new state on the boundary of the envelope can be obtained by 1 1 kprcr kprcr (6.26) Step 5: Stop if 1100,,kk Step 6: Otherwise, set k=k +1 and go to step 3. As schematically explained in Figure 6-12, the limit of the envelope is first found in one parameter while is fixed (Point A ). The reliability result and sensitivity information are calculated at this point, fr om which the new search direction is found using sensitivity information and linear Taylor series expansion. The trial state can be obtained by moving the parameters by l in the search direction. From the trial state, the location of the boundary can be recovered by mo ving in the perpendicu lar direction to the search direction. Using li near search, a new position B on the envelope can be found. This sequence can be repeated to create a closed safety envelope. As expected, the accuracy of the safety envelope can be improved by using a smaller size of the move limit. Example for Multi-Dimensional Load Envelope In the uncertainty model of the dynamic loading formulated in Eq. (6.4) for the front loader frame of CAT 994D, suppose both mean() and amplitude() of the dynamical are normally distributed with standard deviation of 0.25. As discussed in previous section, uncertainty in S-N curve wi ll be modeled as lognormal distributed with a constant coefficient of variance of 0.1 in this problem. A reliability based load tolerance PAGE 103 89 design method can be applied to construct the safety envel ope for the dynamic loads with respect to the mean values of both ra ndom parameters related to the load. Figure 6-13. Construction of load envelope Figure 6-13 provides the proce dures for the construction of fatigue reliability based load envelope. In the front loader frame problem, it is obvious that the mean of and are positive value. Load tolerance can only be de fined as nonnegative values. If the required probability of failure is 10% (arg1.3tet ), following steps can be taken to construct the load envelope: Step 1: Find the initial point (0 ) on the envelope ( (,)1.3 ); Step 2: Find the next solution on load e nvelope using Euler-Newton continuation to meet the constraint1.3 ; Step 3: Since only,0 is meaningful, continue step 2 until the curve end in this region; Figure 6-14 shows the two-dimensional safety envelope for the loader frame while LCCs are both normally distributed. It is clea r from the figure that the system has much PAGE 104 90 more safety margin in the average value than that of the amplitude which is consistent with the one-dimensional case. Figure 6-14. Safety envelop for fatigue re liability of CAT 994D front loader frame Conservative Distribution Type In the previous section, the load ha s been assumed a specific uncertainty characteristic (distribution type and corresponding parame ters). Identifying the load distribution, however, is one of the most difficult tasks in the uncertainty analysis because different operating conditions will yiel d completely different distribution types. Thus, design engineers often look for a conserva tive distribution type. As can be seen in Figure 6-15 and Figure 6-16, for example, l og-normal distribution is more important when the amplitude of the applied load is small, whereas normal distribution is more important when the amplitude is large. Us ing the sensitivity information and linear approximation, it would be possible to predic t which distribution t ype has a significant effect on the load tolerance. Once dominant di stribution type is selected, the detailed load tolerance can be constructed. PAGE 105 91 Figure 6-15. Reliability index with respect to random parameter Figure 6-16. Probability of failure Pf with respect to random parameter In two-dimensional case stated in the prev ious section, based on the same average value and standard deviation of LCCs but different distribution types, e.g., lognormal distribution, the effect of normal and lognormal distribution on safety envelope is compared in Figure 6-17. Pf PAGE 106 92 Figure 6-17. 2-D safety enve lope for different distribu tion type with same random parameters It is obvious that the st ructure working unde r dynamic load modeled with normal distribution is in more severe situation than that of lognor mal distribution, which means uncertainties of dynamic load modeled by norma l distribution is more conservative than the lognormal distribution while both of them have the same random parameters, i.e., mean and standard deviation. Following th e same procedure, the safety boundary for different distribution type can be found. Thus, a conservative safety envelope is constructed by considering all the possible distribution type associated with different possibilities of working conditions. The comp leted safety envelope provides higher confidence of operation and offers a reference for future design upgrade. Summary In this chapter, fatigue reliability based load tole rance design for industrial equipment has been studied. A concept of safe ty envelope has been introduced for fatigue reliability based load tolerance design. The sy stematic road map of safety envelope has PAGE 107 93 been presented. FE-based fatigue evaluati on, SRS-based reliability and sensitivity analysis, path following algorithm are integrat ed to construct a design reference for a structure. By considering the difficulties to obtain the uncertainty characteristics, conservative distribution type is considered to offer safer design of load without complete knowledge of uncertainty properties. The pro cedure of safety enve lope construction has been presented for a two dimensional load model with considering the material uncertainty. PAGE 108 94 CHAPTER 7 ROBUST DESIGN USING STOCHASTIC RESPONSE SURFACE Introduction In the previous chapter, the objective of RBDO is to minimize the cost of product while meeting reliability level of performance. In quality engineering, it has been realized that the deviation from target value of performance due to the uncontrollable input variances/noises results in quality loss, whic h is a measure of dissatisfaction from the customers point of view. Thus, robust desi gn, which targets on making the performance of the product insensitive (robus t) to the noise factors, has been gaining increasing attention in recent re search activities. Traditionally, the performance variance is evaluated either using the Monte-Carlo simulation (MCS) or linear approximation. The computational cost of MCS and the lack of accuracy of the linear approxi mation have been issues in robust design. In this chapter, an efficient and accurate method of evalua ting the performance variance is proposed using the polynomial chaos expansion (Ghane m and Spanos 1991; Isukapalli et al. 1998). The proposed method has comparable accuracy with MCS, while requiring much less computational cost. By selecting appropria te bases of the surrogate model, the performance variance is calculated analytical ly. In addition, the derivatives of the performance variance with respect to design variables and input random parameters are calculated consistently with the variance calcul ation method, which is critical information for design optimization algorithm. PAGE 109 95 In general, the robust design problem s hould not be formulated to reduce the variance alone. Even if robus tness is a requirement from quality point of view, a good design should also satisfy the requirement of the performance. In most of cases, quality and performance requirement are two comp eting design objectives. Thus, the robust design problem becomes a multi-objective optimization problem. In multi-objective optimization, there are multiple optimum designs in a sense that one objective function cannot be reduced further wit hout increasing other objective f unctions. The optimal set is referred to as the Pareto optimal set and yi elds a set of possible answers from which the engineer may choose the desired va lues of the design variables. In this chapter, a numerical example of a cantilever beam with both linear and nonlinear performance is used to show the advantage of SRS-based variance calculations. The variance calculated from the proposed method is compared with that from traditional approximation using first-order Tyler series expansion. Robust design for the natural frequency of a micro-scale cantilever composite beam is also presented. Since the objective is to design a structure whose lo west natural frequenc y reaches resonance frequency with lower variability, the robus t design is modeled as a multi-objective optimization problem with two competitive targets, one for performance mean and another for the standard deviation of the pe rformance. A Pareto optimal front is then obtained by setting one objective as a constraint a nd by gradually changing its constraint limit. It turns out that contro lling design variables makes le ss change of variance of the performance than performance mean. It is mo re important to contro l the input variance itself rather than the design variable in this problem. Global sensitivity is then introduced PAGE 110 96 at the optimum design to address the importan ce of different input variance, which means that random variables should be paid more attention to reduce total performance variance. Performance Variance Calculation Using SRS One important issue for robust design is to evaluate the performance variance. Traditionally, a linear approximation using Tayl or series expansion is often employed for that purpose (Koch et al. 2004). However, th e error of approximation increases according to the nonlinearity of th e performance. In addition, the c oupled effect of input variances cannot be counted in the linear model. As we introduced in chapter 5, perfor mance variance can be obtained by SRS at each design stage. The advantage of the pol ynomial chaos expansion becomes clear in evaluating the variance. In general, the pol ynomial chaos expansion in a surrogate model provides an analytical solu tion for the variance, since Hermite basis functions are orthogonal with respect to the Gaussian m easure and the varian ce of the bases are analytically available. It is convenient to obtain output va riance though the sum of partial variances based on the coefficients of each term. If the polynomial bases are generally defined as ()i with being the vector of standard random variable, the SRS in Eq. (2.4) can be re-written as 1()()N ii iga== (7.1) where g is the approximated system performance and N is the number of coefficients in SRS. Since the above expression is linear with respect to the unknown coefficients, the performance variance can be written as 2 1Var()Var[()]N ii iga== (7.2) PAGE 111 97 Thus, the analytical expression of the perf ormance variation can be obtained if the variations of the polynomial bases are ava ilable. When input va riables are SRV, the analytical variations of Hermite bases can be found in Ghanem and Spanos (Ghanem and Spanos 1991). Variance Sensitivity The robust design problem in this paper is formulated as an optimization problem that minimizes the performance variation in Eq.(7.2). In gradient -based optimization algorithms, calculation of sensitivity information is a critical issue for saving the computational cost and making the algorith m converge. The finite difference method requires a complete recalculation of the perf ormance variation (Haf tka and Gurdal 1991). The goal is to calculate the gradient in formation without carrying out a complete recalculation of the pe rformance variance. From the fact that the SRV in the polynomial chaos remains constant while the design cha nges, the regression coefficients only depend on design variables. Thus, in the proposed pol ynomial chaos expansion, the gradient of the performance variance with respect to jth design parameter, dj, can be written as 1Var() 2Var[()]N i ii jj iga a dd= = (7.3) It is clear that the deriva tives of the regression coeffi cients are enough to calculate the derivative of performance variation. In the linear regression method, the co efficients of SRS are obtained from 1()TT =aXXXg (7.4) where g = [ g1, g2, gM]T is the vector of performance functions at sampling points, and X is the matrix of bases at sampling points, defined as PAGE 112 98 Y11211 12222 12()()...() ()()...() [] ()()...()N N MMNM MN == iX() (7.5) In the above equation, M is number of sampling points, and N is the number of bases. Then, the derivatives of the coefficients can be obtained from 1(()TT jjdd = ag XXX (7.6) The last term, /jd g, is the derivative of performa nce function at sampling points, which can be calculated using design sensit ivity analysis (see C hoi and Kim (Choi and Kim 2004a; Choi and Kim 2004b)). By substituting Eq. (7.6)into Eq.(7.3), the derivative of performance variation can be obtained. This procedure of calculati ng sensitivity of the performance variation is much more efficient than the traditional finite difference method because most information, such as a and X, is already available from the performance variation calculation. The only term required for sensitivity analysis is /jd g. When finite element analysis is used as a computational tool for calculating the performance function, sensitivity analysis provides an efficient tool for calculating the performance derivative. In the context of structural analysis, for example, the discrete system is often represented using a matrix equation of the form [K]{D} = {F}. The performance function g in Eq. (7.4) can be expressed as a function of the nodal solution {D}. Thus, the sensitivity of th e performance can be easily calculated if that of the nodal solution is available. When the design variab les are defined, the matrix equation can be differentiated with respect to them to obtain []{}jjjddd = DFK KD (7.7) PAGE 113 99 Equation (7.7) can be solved inexpensively because the matrix [K] is already factorized. The computational cost of sensitiv ity analysis is usually less than 20% of the original analysis cost, so, local sensitivity can in fact be obtained efficiently. As an illustrative example, a cantilevered beam (Figure 7-1) is taken from literature (Qu and Haftka 2004; Wu et al. 2001). Two failu re modes are considered in this example: (1) the maximum stress of the beam should be less than the streng th of the material [Eq.(7.8)], and the tip deflection should be less than the allowable displacement [Eq.(7.9) ]. These two constraints can be expressed by 1 22600600 ()0 gRYX wtwt = + (7.8) () () 22 3 20 224 0LYX gD Ewt tw= + (7.9) where R is the yield strength, E is the elastic modulus, X and Y are the independent horizontal and vertical loads shown in Figure 7-1. D0 is the allowable tip displacement which is given as 2.25 in. Figure 7-1. Cantilever beam subject to two direction loads Two cross-sectional dimensions, w and t, are considered as controllable design variables. Five uncontroll able random variables are defined in Table 7-1. Table 7-1. Random variables for cantilevered beam structure Random variable X Y R E Distributio n type Normal (500,1002) lb Normal (1000,1002) lb Normal (40000,20002)psi Normal (29E6,(1.45e6)2)psi t w Y X L = 100" PAGE 114 100 It is obvious that the strength constraint defined in Eq.(7.8) is a linear function of the random inputs. For linear performance, th e variance can be anal ytically obtained as 22 1 22600600 Var()Var()()Var()()Var() gRYX wtwt =++ (7.10) Using this property, the accuracy of the proposed variance estimation in Eq.(7.10) can be verified. Table 7-2 shows the comp arison between the variance from the SRSbased method and that from analytical appro ach. The variance is calculated at the deterministic optimal design (w = 1.9574", t = 3.9149"). Table 7-2. Variance estimation of linear performance (strength) Analytical variance Variance from SRS (3rd-order) 2.4e7 2.4e7 In the case of strength constraint, it is po ssible to find the anal ytical expression of the variance. However, in the case of non linear performance, such as deflection constraint in Eq.(7.9), there is no easy way of calcu lating this analytical expression except for the first-order approximation. Th e linear approximation of the deflection constraint becomes () () ()222 222 2linearVar()Var()Var()Var()ggg gXYE XYE =++ (7.11) Due to the error involved in the linear approximation, MC S is the only method that can verify the accuracy of variance calculati on. Since MCS is a sampling-based method, the estimated variance always has variability. Let 2 be the variance of a random variable and let s2 be the unbiased estimator of 2. When n number of samples are used, the variance of the MCS-estimated performance variance can be predicted by (Ang and Tang 1975) PAGE 115 101 4 4 2 43 Var()() 1 n s nn = (7.12) where 4 4() EX= is the fourth central mo ment of random variable X 4 4/ is called kurtosis. For the nonlinear performance in Eq.(7.9), the third-order SRS is used to approximate the deflection and the expression in Eq. (7.2) is used to evaluate the performance variance. Table 7-3 compares the variance obtained from these three methods. As expected, the linear approximation has about 1% error compared with MCS, while SRS-based variance is within the confid ence range of MCS. The error in the linear approximation will increase proportionally to the nonlinearity of the function. Table 7-3. Variance estimation of no nlinear performance (deflection) 2MCSVar() g (500,000 samples) 2linearVar() g 2SRSVar() g 0.1947 (standard deviation = 3.9248E 4) 0.1966 0.1948 Based on the accuracy of the proposed method in calculating performance variance, the variance sensitivity in Eq. (7.3) is also tested using th e cantilevered beam model. Table 7-4 and Table 7-5 show the variance sensitivities obtained from the proposed method compared with those fr om the central finite differe nce method (FDM). In FDM, the design variables are perturbed by 2% and the variance is recalculated using the SRS. When the performance is linear with resp ect to random variables, the analytical sensitivity can be obtained, fo r example, by differentiating Eq. (7.10) with respect to the design variables. In Table 7-4, the sensitivity obtained from SRS agrees well with that of the analytical sensitivity. The finite differe nce sensitivity shows a small error because the variance is still a nonlinea r function with respect to the design variable. PAGE 116 102 Table 7-4. Sensitivity of variance for linear performa nce (strength) Var/w (SRS) Var/w (FDM) Var/w (Analytical)Var/t (SRS) Var/t (FDM ) Var/t (Analytical) 3.6784E7 3.6801E7 3.6785E7 1.2261E7 1.2265E7 1.2261E7 Table 7-5. Sensitivity of variance fo r nonlinear performa nce (deflection) Var/w (SRS) Var/w (FDM ) Var/t (SRS) Var/t (FDM ) 0.6538 0.6544 0.0712 0.0712 Robust Design Two-Layer Beam Dynamic Response of Two-layer Beam The robust design problem formulation is demonstrated using a cantilevered, composite beam, shown in Figure 7-2. When an electric field is applied to the piezoelectric (PZT) part, it will generate a bending moment and deform the beam. On the other hand, when the base is oscillating with a specific frequency, the deformation of the beam will induce an electric field through PZT, which can be used as an energy reclamation device. System dynamic response of the composite beam is highly coupled and the closed-form solution is difficult to obtain (De Rosa 1994; Jang and Bert 1989a; Jang and Bert 1989b). In this chapter, a lump ed element modeling technique (LEM, (Li et al. 2006)) is used to obtain the approximate solution for the system. Under the quasistatic assumption, the LEM can estimate the first fundamental natu ral frequency with accuracy. First, the effectiv e mechanical compliance ( Ce) and the effective mass ( Me) can be calculated by lumping the total strain en ergy and kinetic energy, respectively. The detailed procedure is summarized in the Appe ndix B. The first natural frequency is then calculated using the followi ng expression (Li et al. 2006). PAGE 117 103 11 2n eef CM= (7.13) ts tp b L PZT Shim Figure 7-2. Piezoelectric can tilevered composite beam When the composite beam is used as an energy reclamation device, the maximum efficiency can be obtained when the excita tion frequency and the natural frequency are resonant. Thus, the design goal is to find desi gn variables such that the natural frequency is as close as possible to the excitation freque ncy. However, due to the uncertainty of the material properties, the performance function (natural frequency) in Eq. (7.13)is not a deterministic quantity. Thus, the additional de sign goal is to minimize the variance of the natural frequency due to the input random vari able inputs. Robust Design for Two-Layer Beam When a robust design problem is formulated in such a way that only the variance of the output is minimized, the optimization problem may find an inappropriate design without considering the mean value of the pe rformance. Thus, it would be appropriate to consider both the variance and the mean value simultaneously. In this section, the robust design problem is formulated as a two-objec tive optimization: one for the variance and the other for the mean value. When two objec tives are competing with each other, there will be no single optimum design. Instead, a Pareto optimal front can be constructed, which represents the best combinations be tween the competing objective functions. Due to the uncertainty in inputs, all constraint s are modeled as relia bility constraints. PAGE 118 104 In the two-layer composite beam problem, the goal is to design a structure with natural frequency close to the prescribed value. Considerin g the uncertainties involved in input variables however, the natural frequency at any design will have certain variations, which should also be minimized. In addition, th e reliabilities for the stress and deflection constraints should be considered. In reliab ility-based robust design, the reliability constraints are imposed by pushing the mean valu e to certain levels of standard deviation in the conservation direction. Thus, the robust design problem is formulated as 102 00MinimizeandVar() s.t.()Var()0 ()Var()0f wgfgf RkR DkwD = = + + (7.14) where f is the mean of the first natural frequency; f0 is the excitation frequency; is the maximum stress; R is the material strength, which is assumed as 11,743Pa; w is the tip deflection and D0 is the allowable maximum tip defl ection, which is 7.138 nanometer; and k is the user-defined constant for specific target reliability level. It is assumed that uncertainties only exist in the material prope rties such as elastic modulus and material density. Table 7-6 lists the random parameters of these quantities. All random variables are assumed to be normally distributed and th e standard deviation for the elastic modulus is 10% of the mean value and that of the density is 5%. Table 7-6. Random parameters for the composite beam structure Random variable Mean Standard deviation Youngs modus of shim (Es) 169 GPa 16.9 GPa Density of shim (s) 2330 kg/m3 116 kg/m3 Youngs modus of PZT (Ep) 60 GPa 6 GPa Density of PZT (p) 7500 kg/m3 375 kg/m3 PAGE 119 105 In the composite beam problem, three de sign variables are defined: beam length L shim thickness ts, and PZT layer thickness tp. The robust design problem involves three deterministic design variables and four random parameters. For given design variables, the SRS for the performance functions, such as natural frequency, stress, and tip deflection, are constructed according to Eq.( 2.4). Then, the performance variances are calculated from Eq. (7.2) and variance sensitivity from Eq.(7.3). The values and sensitivities of the two objective functions g1 and g2 with respect to the three design variables are summarized in Table 7-7 at the initial design ( ts = 6m, tp = 0.2m, L = 1000m). Table 7-7 shows that for given design, the mean of the frequency will change at least 15 times more than the frequency vari ance. Thus, it is easier to change the mean values than to change the frequency vari ance. This observation leads to the idea of controlling the input variances directly rather than controlling the design variables in the following section. Table 7-7. Sensitivities of objectiv e functions at the initial design ( ts = 6m, tp = 0.2m, L = 1000m) g1(Hz) g2(Hz) g1/ts (Hz/m) g1/tp (Hz/m) g1/L (Hz/m) g2/ts (Hz/m) g2/tp (Hz/m) g2/L (Hz/m) 834.88 144.54 491.02 562.05 5.335 27.15 36.55 0.29 Since two objective functions are compe ting with each other, there will be no single optimum design. In such a case, the va lue of one objective function is fixed and then the minimum value of the other objectiv e function can be found. By repeating this procedure for different values, a Pareto optimal front can be constructed. Figure 7-3 shows the Pareto optimal front of the tw o-objective optimization problem in Eq.(7.14). PAGE 120 106 All points in the Pareto front are optimum de sign in a sense that one objective function cannot be reduced further without incr easing the other ob jective function. Figure 7-3. Pareto optimal front for th e robust design of the composite beam Global Sensitivity Analysis In Figure 7-3, the change in the mean value (abscissa) is more significant than that in the standard deviation (ordinate), which is consistent with the observation in Table 7-7. This result indicates that when the design variables are deterministic, it is relatively easier to change the mean value rather than the performance variance. The performance variance can be changed more effectively by controlling input variance. However, controlling input variance accompanies manuf acturing cost or large numbers of coupon tests. Thus, in practice, it is important to find the contribu tion of random variables to the performance variance and then spend more re sources in controlling the most significant random variable. In Table 7-8, the contribution of rand om input variables to the performance variance is summarized in terms of total sens itivity indices in Eq.( 4.6). It can be found PAGE 121 107 that the contributions of s and Es are more than 99% of the performance variance. Thus, it will be meaningful to reduce the variance of the shim rather than that of the PZT. Table 7-8. Total sensitivity indices for the composite beam structure ( ts = 6m, tp = 0.2m, L = 1000m) totalpES totalpS totalSES totalSS 0.00% 0.96% 85.87% 13.17% Robust Design by Tolerance Control In the previous section, the input vari ances were considered as uncontrollable variables and only deterministic design variables were consid ered. However, in tolerance design, the design variables are fixed, wh ile the variances of random variables are changed to further reduce the output variances. However, in such a problem, the optimum design will reduce all input variances to zero. Thus, the robust design will turn out to be zero variance. In practice, reducing input variance requires cost. Different costs are anticipated in reducing the variance of different inputs. The co st of controlling indi vidual input variance can be represented by a cost-tolerance mode l (Chase and Greenwood 1988). Thus, a more appropriate robust design problem will be: fo r a given investment, how much individual variance should be reduced in order to mi nimize the performance variance. Based on the total budget of controll ing input variability, the robust design problem can be written as []12 total 1MinimizeVar(,,...,) s.t.()n n ii ig CC = (7.15) where i is the standard deviation of ith random variable, Ci is the cost function of controlling ith standard deviation, and Ctotal is the total investment. PAGE 122 108 Similar to the robust design problem with deterministic desi gn variables, the optimization problem in Eq. (7.15) requires the derivative of the performance variance. The only difference now is that the derivative is taken with respect to the input variance. By substituting jth design variable dj in Eq. (7.3) to jth random parameter j the gradient of output variance in Eq. (7.15) with respect to jth random parameter can be written as 1Var() 2Var[()]N i ii jj iga a= = (7.16) Similarly, the derivatives of the coefficients can be obtained from 1()TT jj = ag XXX (7.17) Since all random variable are assumed to be independent, we have the following chain rule of differentiation: 1 1() ()jj jj jjT T = gg (7.18) where Tj is the transformation of jth random variable from original random space to standard normal space: ()jjjTx= (7.19) Therefore, the sensitivity of performa nce variance with respect to random parameters can be obtained by combining Eqs.(7.16),(7.17) and (7.18) if the derivative 1/()/jjjTx = gg is available. As an illustration of the effectiveness and convergence properties of the proposed approach, the cantile vered beam model (Figure 7-1) is used. Based on the accurately estimated performance variance, the variance se nsitivities with respec t to input variances are calculated using the proposed method in Eq.(7.16). Table 7-9 and Table 7-10 show the sensitivities obtained fr om the proposed method along with those from the finite PAGE 123 109 difference method. Since the analytical sensitiv ity is available for th e linear performance, Table 7-9 also lists the analytical sensitivity It turns out that the proposed SRS-based sensitivity calculation method provides accu rate sensitivity information. Since the proposed method only requires th e calculation of performance sensitivity at sampling points [Eq.(7.18)], the computational cost will be much less than that of the finite difference method. Table 7-9. Sensitivity of variance for linear performa nce (strength) Var/ X (SRS) Var/ X (FDM ) Var/ X (Analytic ) Var/ Y (SRS) Var /Y (FDM ) Var /Y (Analytic ) Var /R (SRS) Var /R (FDM ) Var /R (Analytic) 3.2e5 3.2e5 3.2e5 8.0e4 8.0e4 8.0e4 4000 4000 4000 Table 7-10. Sensitivity of variance fo r nonlinear performance (deflection) Var/X (SRS) Var/X (FDM) Var/Y (SRS) Var/Y (FDM) Var/E (SRS) Var/E (FDM) 3.41e-3 3.41e-3 5.77e-5 5.77e-5 2.66e-8 2.66e-8 With performance variance and its sensitiv ity, the robust design problem with input variance control in Eq. (7.16) can be solved efficien tly. Consider the robust design problem that minimizes the variance of natu ral frequency with strength and deflection constraints, as 11 22 1MinimizeVar() s.t.()()0 ()()0 ()n iitot iEgkg Egkg CC = (7.20) PAGE 124 110 where g1 and g2 are strength and deflect ion constraints in Eq. (7.8) and (7.9), respectively. In this optimization problem, the deterministic design variables, w and t are pre-determined ( w = 2.73, t = 3.50) from the previous optimization. Now, the optimization is performed by changing the standa rd deviations of random input variables. In Eq.(7.20), is the first natural frequency of the beam defined as 2 2 4() 23 EItE L AL == (7.21) and E () and () represent the expect value and st andard deviation of random output, respectively, and Ci(i) is the cost-tolerance function for the ith random variable. For a specific boundary condition, the term, is constant. Thus, th e objective function to control the variance of natural frequency is modified to 2Minimize VarVar 23 tE = (7.22) Table 7-11 lists random variables and cost-tolerance functions (Chase and Greenwood 1988) for the random variables. Table 7-11. Random variables a nd cost-tolerance functions Variables X Y R E Distributio n N(500,12)lb N(1000,22)l b N(40000,32)psi N(29E6,42) psi N(0.28, ,52) Costtolerance 1.5+200/1 1.5+200/2 1.5+1.6* 107/ 32 200Exp(4 10-6) 18Exp(1005) i 25 1 20050 2 4001000 3 4000106 4 3*1 060.01 5 0.05 To demonstrate the robust design, total cost of controlling variance at the initial design has been chosen as cost constraint Thus, the design goal is to minimize the performance variance, while maintaining the same cost with initial variance control. Table 7-12 shows that the standard deviation of natural frequency reduced from 452.5 Hz PAGE 125 111 to 325 Hz by redistributing the input vari ances. Since the natu ral frequency is independent of the applied loads and the two constraints are not active, the final design increased the variances of the first th ree random variables. The optimum design maintains the variance of the el astic modulus and halves the density, which is more cost effective than reducing the va riance of the el astic modulus. Table 7-12. Random variables a nd cost-tolerance functions Design Variables Initial design Optimal design 1 100 200 2 100 400 3 2000 4000 4 1.45*106 1.45*106 5 0.02 0.010762 Objective 452.5 325.0 Ctotal 17.8274 17.8274 Summary In this chapter, SRS-based variance calcu lation is proposed to facilitate robust design application. Accurate variance sensitiv ity analysis is presented for the gradientbased optimizer. A simple cantilevered beam with two failure modes, one as linear and another as nonlinear, is used to illustra te the accuracy and robustness of variance calculation. Robust design for the natural frequency of a cantilevered, com posite beam showed that because controlling deterministic de sign variables makes less change of the performance variance than that of the performance mean, we found it is more important to control the input variance itse lf rather than the design vari able in our specific problem. Global sensitivity is then in troduced to address which random variables should be paid more attention to reduce total performance variance. PAGE 126 112 Finally, a cost model based robust design is proposed to control the input variance, an alternative way of tolerance design. Design sensitivity analysis of performance variance with respect to input variance has been proposed in mathematical programming. The cantilever beam model is used to illustrate the effectiveness of tolerance design. PAGE 127 113 CHAPTER 8 SUMMARY AND RECOMMENDATIONS Although reliability-based design optimiza tion has been studied intensively for decades, industrial applications in structur al design are still limited by the significant amount of computational cost as well as accuracy in reliab ility analysis. RBDO techniques are under development and have large room for improvement. To improve the efficiency of RBDO, there are two key factors: one is by introducing new RBDO strategies; and the othe r is by developing efficient and accurate uncertainty analysis methodologies. In the firs t category, two strategi es are investigated in this research. The conventional RBDO strategy sets the required reliability as a constraint, while the inverse measure approach uses the performance measure at required reliability as a constraint. The mathematical equivalence of these two strategies has been discussed in section 3.2. However, since constraints are estab lished using different measures, the convergence and optima may result in different ways. In the second category, the improvement of efficiency of uncertainty analysis has always been a main concern. In this research, SRSM using local sensitivity information is implemented in uncertainty analysis. The convergence and accuracy of SRSM are investigated. During the implementation of SRSM to RBDO model, SRS-based probability sensitivity analysis is devel oped and tested in order to improve the convergence and efficiency of RBDO. The e fficiency of SRSM is further improved significantly by utilizing variance-based gl obal sensitivity analysis to reduce the dimension of random space. PAGE 128 114 In the manufacturing environment, the cost of controlling manufacture variations is often more than making the design insensitive to these variations. Thus, it is necessary to study robust design under the context of uncertain ty. In this research, SRS-based variance calculation is proposed to facilitate the robus t design application. It is shown that the proposed variance calculation is more accurate than the conventiona l first-order Taylor series expansion. The proposed method can in clude the higher-order terms as well as interactions. Accurate variance sensitivity an alysis is further presented for the gradientbased optimizer. Numerical example shows that it is sometimes more important to control the input variance its elf rather than the design va riables. A cost model based robust design is then proposed to control the input variance with the same cost, an alternative way of tolerance design. Design sensitivity analysis of performance variance with respect to input variance has been proposed in mathematical programming. Numerical example is also used to illustra te the effectiveness of tolerance design. As an integrated process, fatigue reliab ility-based load tole rance design involves finite element analysis, fatigue life predicti on, reliability analysis and path following technique. The randomness of dynamic loads is s ubjective and difficult to control. In this research, different possible distribution type s are considered to provide a conservative and safer load design. Implementing more dist ribution types of uncertainty in the load tolerance design and finding a more efficient way to construct the safety envelope are recommended for future research. PAGE 129 115 APPENDIX A SAMPLING-BASED PROBABILITY SENSITIVITY ANALYSIS FOR DIFFERENT DISTRIBUTION TYPE The derivative of failure probability can be written as 1()() (()0)(()0)() () () (()0)() ()fP ff I GdIGfd f f IGd f xx ux=T(u)xx xxxxx x x uuu x (A.1) where is the random parameter, x, udenote entire original design space and standard normal space, respectively; ()0 G X is failure region and () f is joint probability density function (PDF); () is standard normal PDF ; (()0) IG X is an indication function such that I =1 if ()0 G X and I =0 otherwise. In this section, explicit mathematical expressions of sampling-based probability sensitivity analysis for different distributi on types are derived. Nu merical examples of analytical functions have been used to check the accuracy of the derivation. Normal Distribution 2(,)iiiXN Case 1: i 1 1()1 () ()iii iiii x u f f x=T(u) x=T(u)x x (A.2) 11() (()0)() () 11f ii N j ji j iP f I Gd fu Iu N ux=T(u)x uuu x (A.3) PAGE 130 116 Case 2: i 1 12 2()11 ()()(1) ()iiii i iiiiixx f u f x=T(u) x=T(u)x x (A.4) 11() (()0)() () 11 (1)f ii N jj jii j iP f I Gd f Iuu N ux=T(u)x uuu x (A.5) The accuracy of above sensitivity formula tion is evaluated using a simple linear function. The random variable is normally dist ributed. In this particular example, FORM provides an exact solution. Consider a simple one-dimensional linear function represented by 2()1.63 (0,0.4)Gxx xNormal (A.6) A failure criterion for the performance function is set to be G 0. In FORM, probability sensitivity can be obtained through the sensitivity of reliability index through ()fP (A.7) Table A-1 compares the accuracy of Pf and its sensitivities with those from FORM. In the sampling-based probabilistic sensit ivity calculation, 200,000 samples are used. A good agreement between two approaches is observed. Table A-1: Accuracy of proposed probabil ity sensitivity method for normal distribution using 200,000 sampling MCS FORM Sampling based approach(200,000 samples) Difference Pf 0.0915 0.0914 0.11% dPf/d x 0.4100 0.4109 0.22% dPf/d x 0.5469 0.5484 0.27% PAGE 131 117 Uniform Distribution The probability density function for unifo rm distribution can be written as: 1 (), f xaxb ba (A.8) 2 ab 12 ba (A.9) Thus, 3 3 a b (A.10) This distribution can be mode led using a step function, as 1 () 0 x a Hxa x a (A.11) However, it is not possible to calculate the sensitivity of the step function. Thus, an arctangent function is used to a pproximate the step function, as 1 ()[()()] 1 arctan()arctan() () fxHxaHxb ba cxacxb ba (A.12) when c 1 arctan()arctan()() () cxacxbfx ba The sensitivity of failure probability Pf to the mean of random variable xi can be written as ()f iiP f d Xx x (A.13) For an N dimensional system, by assuming all system random variables are independent, the joint probabil ity function is defined as 11 ()arctan()arctan() ()N iiii i ii f cxacxb ba x (A.14) The derivative of this joint proba bility function can be written as PAGE 132 118 ()()() ()()ii iiiii iiab fff ab ff ab xxx xx (A.15) From Eq.(A.13), sensitivity of failure probabi lity with respect to mean value of design variable can be derived as 2222 1 1 2222() (()0) 11 (()0)[] ()1()1() 1 arctan()arctan()... () 11 [] ()1()1() (()0){ 1 (f ii iiiiii N jjjjN j jj ji iiiiii iP f IGd c IG bacxbcxa cxacxbdxdx ba c bacxbcxa IG b X Xx xx x x 12222 2222}() arctan()arctan() ) 11 [] 1()1() (()0)() arctan()arctan() 11 1()1() arctan()arctiiii i iiii iiii jj iiii j ii f d cxacxb a c cxbcxa IGd cxacxb cxbcxa c N cxa X U x=T(u)xx uuu 11an()N j j iicxb x=T(u) (A.16) The accuracy of the above sensitivity formul a is used to calcul ate the sensitivity of the linear function in Eq. (A.6) with uniform distribution, as ()1.63 (,) 2,1 Gxx x Uniformab whereab (A.17) When G 0, the performance function is considered to be fail. PAGE 133 119 For this linear function of uniform distri buted variable x, the exact probability sensitivity is 11 3 ba Table A-2 shows sensitivity re sults obtained from FORM and the sampling based approach. Although FORM cannot provide the exact solution for random variable with non-normal distribution, it is still a good approximation for reliability analysis. The sampling based approach is also a good estimation of both failure probability and probabi lity sensitivity. Table A-2: Accuracy of proposed probability sensitivity method for uniform distribution using 200,000 sampling MCS FORM Sampling based approach (c=10000,200,000 samples) Exact Pf 0.1555 0.1556 0.1556 (7/45) d Pf/d x 0.3334 0.3333 0.3333 (1/3) Log-Normal Distribution When a random variable x has lognormal distribution, the probability density function is defined as 2 2(ln) 21 () 2xfxe x (A.18) where 22ln(1) v 2ln 2 v a Transformation from standard nor mal space can be return as: ()UXUe T (A.19) The sensitivity of failure probability Pf to the mean of random variable xi can be written as PAGE 134 120 1() (()0)() ()f iiP f I Gd f ux=T(u)x uuu x (A.20) For an n dimensional system, by assuming all system random variables are independent, the joint probabil ity function is defined as 2 2[ln()] 2 11 () 2ii ix n i iife x x (A.21) Thus, we can have ()()()ii iiiiifff xxx (A.22) where 2()ln() ()iiii ii iifxx fx 2ln() () ()ii iix f f x x (A.23) 2 21 ()(1)ii iiiiv av (A.24) 2 3()()[ln()] ()iiiiii ii iiifxfxx fx 2 3[ln()] ()() ()ii iiix ff f xx x (A.25) 2 2()(1)i iiv av (A.26) 2 22 2 2 23ln() ()1 ()()(1) [ln()] 1 ()(1)iii iiiii ii iiixv f f av x v av x x (A.27) 1 12 22 2 2 1 23() (()0)() () ln() 1 ()(1) 1 [ln()] 1 ()(1)f ii j iii N iiii j j ii iiiP f IGd f xv av N x v av ux=T(u) x=T(u)x uuu x (A.28) PAGE 135 121 The correctness of the derivation is prove d by comparing the sensitivity result of sampling based approach with that of FO RM. Linear performance function in Eq.(A.6) is used as test function. 2()1.63 (0.5,0.4) Gxx xLogNormal (A.29) Because of the simplicity of the perf ormance function, the exact solution for probability of failure can be found by directly integrating the lognormal PDF from 1.6/3 to infinity(1.6-3x<0). As shown in Tabl e A-3, the sampling-based approach with c=10,000 and 200,000 samples is more accurate than FORM. Table A-3: Accuracy of proposed proba bility sensitivity method for Log-normal distribution using 200,000 sampling MCS FORM Sampling based approach (c=10000,200,000 samples) Exact Pf 0.3287 0.3287 0.3287 dPf/d x 1.1757 1.1763 1.1765 PAGE 136 122 APPENDIX B NATURAL FREQUENCY OF CANT ILEVER COMPOSITE BEAM Bending Moment As indicated in Figure 7-2, the cantilev er composite beam subjects to a bending moment ( M0) at the ends of the piezoceramic. This is caused by induced strain from applied voltage (Cattafesta et al. 2000). Fi gure A-1 replaces the mass of the composite beam as an equivalent uniform load ( q ) due to its weight. R and Mr are the reaction force and bending moment at the clamp. Figure B-1: Free body diagra m of two-layer beam Thus, the bending moment in the com posite beam can be expressed as 2 0(0) 2rx MxxLMRxMq =+ (B.1) where RqL = and 2/2rMqL Geometric Properties of Composite Beam Before we calculate the effective compliance and lumped mass, geometric properties such as location of neutral axis and flexural rigidity of the composite beam are required in static analysis of the beam. M0 R M0 Mr q L PAGE 137 123 If we define c2 as the location of the neutral ax is from the bottom of piezoceramic and ( EI )c as equivalent flexural ri gidity in composite beam, they can be calculated by the following two expressions: () 2 222 ,p s sspp ssppt t EttE c EtEt ++ = + (B.2) ()csscppcEIEIEI =+ (B.3) where Isc and Ipc are the moment of inertia of the sh im and PZT layer with respect to its own neutral axis, respectively. Effective Compliance for Composite Beam To find the effective compliance for the composite beam in Eq. (7.13), we need to use total potential energy in the beam as shown Eq. (B.4): () () 2 2 2 02L cEI dwx PEdx dx = (B.4) where 2 432() 24()6()4()cccqqLqL wxxxx EIEIEI = + (B.5) Equation (B.5) is obtained by conventional Eu ler-Bernoulli beam theory. Thus, by lumping the overall potential strain energy at th e tip, an effective short circuit mechanical compliance for the composite beam will be calculated as 2() 2Ftip ew C PE = (B.6) Effective Mass for Composite Beam In order to calculate the e ffective lumped mass in Eq. (7.13), total kinetic energy in the composite beam [Eq. (B.7)] will be used. ()2 02L LcKEwxdx = (B.7) PAGE 138 124 where Lc is the equivalent mass density of the composite beam and ()wx is the velocity in the beam. For a simple harmonic motion, the veloci ty of the beam are related to the displacement by ()()wxjwx = (B.8) () wx is then expressed as () ()F Ftip tipwx wxw w = (B.9) Effective mass for the composite beam fr om its deflection shape is obtained by lumping the kinetic energy of the beam at its tip: ()2 22 02 ()FFL Lc e tiptipKE Mwxdx ww == (B.10) PAGE 139 125 LIST OF REFERENCES Allen, J. K., Seepersad, C., Choi, H., a nd Mistree, F. (2006). "Robust Design for Multiscale and Multidisciplinary Applications." Journal of Mechanical Design 128(4), 832-843. Allgower, E. L., and Georg, K. (1990). Numerical Continuation Methods: An Introduction Springer-Verlag, Berlin. Ang, A. H.-S., and Tang, W. H. (1975). Probability Concepts in Engineering Planning and Design, Basic Principles Wiley. Arora, J. (2004). Introduction to Optimum Design Academic Press. Ayyub, B. M., Assakkaf, I. A., Kihl, D. P., and Siev, M. W. (2002). "Reliability-Based Design Guidelines for Fatigue of Ship Structures." Naval Engineers Journal 114(2), 113-138. Bannatine, J. A., Comer, J. J., and Handrock, J. L. (1990). Fundamentals of Metal Fatigue Analysis Prentice-Hall, Englewood Cliffs, N. J. Basquin, O. H. (1910). "The Exponential Law of Endurance Tests." Proc Am Soc Test Mater 10, 625-630. Breitung, K. (1984). "Asymptotic Approximations for Multinormal Integrals." Journal of Engineering Mechanics 110(3), 357-366. Brown, M. W., and Miller, K. J. (1973). "A Theory for Fatigue under Multiaxial Stressstrain Conditions." Proceedings of the Instituti on of Mechanical Engineer 187, 745-756. Cameron, R. H., and Martin, W. T. (1947). "The Orthogonal Development of Nonlinear Fuctionals in Series of Fouriers-Hermite Functionals." Ann. Math. 48, 385-392. Cattafesta, L., Mathew, J., and Kurdila, A. (2000). "Modeling and De sign of Piezoelectric Actuators for Fluid Control," San Diego, CA. Chandu, S. V. L., and Grandi, R. V. (1995). "General Purpose Procedure for Reliability Based Structural Optimization unde r Parametric Uncertainties." Advances in Engineering Software 23, 7-14. PAGE 140 126 Chase, K. W., and Greenwood, W. H. (1988). "Design Issues in Mechanical Tolerance Analysis." Manufacturing Review, ASME 1(1), 50-59. Chen, C. J., and Choi, K. K. (1996). "Robus t Design using Second-order Shape Design Sensitivity Analysis." Center for Comput er-Aided Design, University of Iowa, Iowa City, IA. Chen, W., Jin, R., and Sudjianto, A. ( 2005). "Analytical Variance-based Global Sensitivity Analysis in Simulationbased Design under Uncertainty." Journal of Mechanical Design 127, 875-886. Chen, X., Hasselman, T. K., and Neill, D. J. (1997). "Reliability Based Structural Design Optimization for Pract ical Application." 38th AIAA/ASME/ASCE/AHS/ASC Structure,Structural Dynamics,and Materials Conference Kissimmee,Florida. Choi, K. K., and Kim, N. H. (2004a). Structural Sensitivity Analysis and Optimization I: Linear Systems Springer. Choi, K. K., and Kim, N. H. (2004b). Structural Sensitivity Analysis and Optimization II: Nonlinear Systems and Applications Springer. Choi, S. K., Grandhi, R. V., Canfield, R. A ., and Pettit, C. L. (2003). "Polynomial Chaos Expansion with Latin Hypercube Sampli ng for Predicting Response Variability." 44th AIAA/ASME/ASCE/AHS Structures ,Structural Dynamics,and Materials Conference Norfolk,Virginia. Chopra, O. K., and Shack, W. J. (2003). "Rev iew of the Margins for ASME Code Fatigue Design Curve Effects of Surface Roughne ss and Material Variability." Argonne National Laboratory, Argonne, IL. Coffin, L. F., Jr. (1954). "A Study of the Effe cts of Cyclic Thermal Stresses in a Ductile Metal." ASME Transaction 16, 931-950. De Rosa, M. A. (1994). "Free Vibrations of Stepped Beams with Elastic Ends." Journal of Sound and Vibration 173(4), 563-567. Devroye, L. (1986). Noneuniform Random Variate Generation Springer-Verlag, New York. Du, X., and Chen, W. (2002a). "Efficien t Uncertainty Analysis Methods for Multidisciplinary Robust Design." AIAA Journal 4(3), 545-552. Du, X., and Chen, W. (2002). "Sequential Optimization and Reliability Assessment Method for Efficient Pr obabilistic Design." ASME 2002 Design Engineering Technical Conferences and Computers and Information in Engineering Conference Montreal,Canada. PAGE 141 127 Enevoldsen, I., and Sorensen, J. D. (1994). "R eliability-Based Optimization in Structural Engineering." Structural Safety 15, 169-196. Fe-safe. (2004). "Software package, Version 5.1." Safe Technology, Sheffield, England. Fuchs, H. O., and Stephens, R. I. (1980). Metal Fatigue in Engineering John Whiley & Sons, New York. Gautschi, W. (1996). "Orthogonal Polynomials: Applications and Computations." Acta Numerica, A. Iserles, ed., Cambridge University Press, Cambridge, 45-119. Gerber, H. (1874). "Bestimmung der Zulassi gen Spannungen in Eisen-konstruktionen." Z Bayer Arch Ingenieur-Vereins 6, 101-110. Ghanem, R., and Ghiocel, D. (1996). "A Comparative Analysis of FORM/SORM and Polynomial Chaos Expansions for Highly Nonlinear Systems." 11th ASCE Engineering Mechanics Conference Fort Lauderdale,Florida, 535-538. Ghanem, R. G., and Spanos, P. D. (1991). Stochastic Finite Elements: A Spectral Approach Springer-Verlag, New York. Goodman, J. (1899). Mechanics Applied to Engineering Longmans Green, London. Grandhi, R. V., and L.P., W. (1998). "Reliabi lity-Based Structural Optimization using Improved Two-Point Adaptive Nonlinear Approximations." Finite Elements in Analysis And Design 29(1), 35-48. Griffith, A. A. (1921). "The Phenomen a of Rupture and Flow in Solids." Phil. Trans. Roy. Soc. of London A221, 163-197. Haftka, R. T., and Gurdal, Z. (1991). Elements of Structural Optimization (Solid Mechanics and Its Applications) Springer. Haldar, A., and Mahadevan, S. (2000). Reliability Assessment using Stochastic Finite Element Analysis John Wiley & Sons Inc. Hasofer, A. M., and Lind, N. C. (1974). "Exact and Invariant Second-Moment Code Format." Journal of the Engineering Mechanics 100(EM1), 111-121. Hoeppner, D. W., and Krupp, W. E. (1974) "Prediction of Component Life by Application of Fatigue Cr ack Growth Knowledge." Engineering Fracture Mechanics 6, 47-70. Inglis, C. E. (1913). "Stresses in A Plate due to the Presence of Cracks and Sharp Corners." Transaction of the Institute of Naval Architects 55, 219-241. PAGE 142 128 Irwin, G. R. (1956). "Onset of Fast Crack Propagation in High Strength Steel and Aluminum Alloys." Sagamore Research Conference Proceedings 289-305. Isukapalli, S. S., Roy, A., and Georgopoulos, P. G. (1998). "Stochastic Response Surface Methods(SRSMs) for Uncertainty Propagatio n: Application to Environmental and Biological Systems." Risk Analysis 18(3), 351-363. Isukapalli, S. S., Roy, A., and Ge orgopoulos, P. G. (2000). "Efficient Sensitivity/Uncertainty Analysis Using the Combined Stochastic Response Surface Method and Automated Differentiat ion: Application to Environmental and Biological Systems." Risk Analysis 20(5), 591-602. Jang, S. K., and Bert, C. W. (1989a). "Fr ee Vibration of Stepped Beams:exact and numerical solutions." Journal of Sound and Vibration 130(2), 342-346. Jang, S. K., and Bert, C. W. (1989b). "Free Vibration of Stepped Beams:higher mode frequencies and effects of steps on frequency." Journal of Sound and Vibration 132(1), 164-168. Karamchandani, A., and Cornell, C. A. ( 1992). "Sensitivity Estimation within First and Second Order Reliability Methods." Structural Safety 11, 95-107. Khuri, A. I., and Cornell, J. A. (1996). Response Surface: Design and Analysis Marcel Dekker,Inc., New York. Kim, N. H., Choi, K. K., and Botkin, M. E. (2003). "Numerical Method for Shape Optimization Using Meshfree Method." Structure and Mul tidisciplinary Optimization 24, 418-429. Kim, N. H., Choi, K. K., and Chen, J. S. (2000). "Shape Design Sensitivity Analysis and Optimization of Elasto-plasticity with Frictional Contact." AIAA Journal 38(9), 1742-1753. Kim, N. H., Wang, H., and Queipo, N. V. (2004). "Adaptive Reduction of Design Variables using Global Sensitivity Indice s in Reliability-Bas ed Optimization." 10th AIAA/ISSMO Multidisciplinary A nalysis and Optimization Conference Albany, New York. Kim, N. H., Wang, H., and Queipo, N. V. (2004). "Efficient Shape Optimization Technique using Stochastic Response Surfaces and Local Sensitivities." ASCE Joint Specialty Conference on Probabilisti c Mechanics and Structural Reliability Albuquerque, New Mexico. Koch, P. N., Yang, R.-J., and Gu, L. (2004). "Design for Six Sigma through Robust Optimization." Structural and Multidis ciplinary Optimization 26, 235-248. PAGE 143 129 Krishmamurthy, T. (2003). "Response Surface Approximation with Augmented and Compactly Supported Radial Basis Functions." 44th AIAA/ASME/ASCE/AHS Structures,Structu ral Dynamics, and Materials Conference Norfolk,Virginia. Kwak, B. M., and Kim, J. H. (2002). "Concept of Allowable Load Set and Its Application for Evaluation of Structural Integrity." Mechanics of structures and machines 30(2). Kwak, B. M., and Lee, T. W. (1987). "Sensi tivity Analysis for Reliability-Based Optimization Using an AFOSM Method." Computers & Structures 27(3), 399406. Lauridensen, S., Vitali, R., Van Keulen, F ., Haftka, R. T., and Madsen, J. I. (2001). "Response Surface Approximation Using Gradient Information." The Fourth World Congress of Structural and Multidisciplinary Optimization Dalian, China. Li, Y., Papila, M., Nishida, T., Cattafest a, L., and Sheplak, M. (2006). "Modeling and Optimization of A Side-implanted Piezores istive Shear Stress Sensor." San Diego, CA. Liang, J., Mourelatos, Z. P., and Tu, J. (2004). "A Single-Loop Method for ReliabilityBased Design Optimization." Proceedings of DETC'04, ASME 2004 Design Engineering Technical Conferences and Computers and Information in Engineering Conference Salt Lake City, Utah. Liu, H., Chen, W., and Sudjianto, A. (2004). "Probability Sensitivity Analysis Methods for Design under Uncertainty." 10th AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference Albany, New York. Liu, P. L., and Kiereghian, A. D. (1991). "Optimization Algorithms For Structural Reliability." Structural Safety 9, 161-177. Mahadevan, S., and Shi, P. (2001). "Multiple Linearization Method for Nonlinear Reliability Analysis." Journal of Engineering Mechanics 127(11), 1165-1173. Malkov, V. P., and Toropov, V. V. (1991). "Simulation Approach to Structural Optimization Using Design Sensitivity Analysis." Engineering Optimization in Design Process, Proc. Int. Conf ., Lecture Notes in Engineering Karlsruhe,Germany, 225-231. Manson, S. S. (1954). "Behavior of Material s Under Conditions of Thermal Stress." National Advisory Comittee for Aeronautics, Cleveland. Matsuishi, M., and Endo, T. (1968). "Fatigue of Metals Subjected to Varying StressFatigue Lives Under Random Loading." Proceedings of the Kyushu District Meeting 37-40. PAGE 144 130 Melchers, R. E. (2001). Structural Reliability Analysis and Prediction John Wiley & Sons, New York. Metropolis, N., and Ulam, S. ( 1949). "The Monte Carlo Method." Journal of the American Statistical Association 44(247), 335-341. Miller, W. R., Ohji, K., and Marin, J. (1966) "Rotating Principal Stress Axes in High Cycle Fatigue." ASME Paper No. 66-WA/Met 9 New York. Miner. (1945). "Cumulative Damage in Fatigue." Journal of Applied Mechancis 12, A159-A164. Myers, R. H., and Montgomery, D. C. (1995). Response Surface Me thodology: Process and Product Optimization Using Designed Experiments John Wiley & Sons,Inc, New York. Paris, P. C. (1964) "The Fracture Mechanics Approach to Fatigue." Fatigue-An Interdisplinary Approach, Proceedings of the Tenth Sagamore Army Materials Research Conference Syracuse, N.Y., 107-132. Paris, P. C., and Erdogan, F. (1963). "A Cr itical Analysis of Crack Propagation Laws." Journal of Basic Engineering, ASME Transactions, Series D 85(4), 528-534. Paris, P. C., Gomez, M. P., and Anderson, W. P. (1961). "A Rational Analytic Theory of Fatigue." The Trend in Egineering 13, 9-14. Park, G. J., Lee, T. H., Lee, K. H., and Hwang, K. H. (2006). "Robust Design: An Overview." AIAA Journal 44(1), 181-191. Qu, X., and Haftka, R. T. (2004). "Relia bility-based Design Optimization using Probabilistic Sufficiency Factor." Structural and Multidisciplinary Optimization 27(5), 314-325. Qu, X., Venkataraman, S., and Haftka, R. T. (2000). "Response Surface Options for Reliability-Based Optimization of Composite Laminate." 8th ASCE Speciality Conference on Probabilistic Mechan ics and Structural Reliability Notre Dame,Indiana. Rackwitz, R.(2000). "Reliability An alysis-Past,Present and Future." 8th ASCE Specialty Conference on Probabilistic Mec hanics and Struct ure Reliability Notre Dame,Indiana. Rahman, S., and Xu, H. (2004). "A Univar iate Dimension-reduction Method for Multidimensional Intergration in Stochastic Mechanics." Probabilistic Engineering Mechanics Article in Press. PAGE 145 131 Rijpkema, J., Etman, L., and Schoofs, A. (2000). "Metamodel Based Design Optimization Including Design Sensitivities." Proc. SMSMEO Workshop Lyngby,Denmark. Rubinstein, R. Y. (1981). Simulation and the Monte Carlo Method John Wiley & Sons,Inc. Saltelli, A., Chan, K., and Scott, E. M. (2000). Sensitivity Analysis John Wiley & Sons, Ltd. Saltelli, A., Tarantola, S., and Chan, K. (1999). "Quantitative Model-independent Method for Global Sensitivity Analysis of Model Output." Technometrics 41(1), 39-56. Smarslok, B. P., and Haftka, R. T. (2006). "T aking Advantage of Separable Limit States in Sampling Procedures." 47th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, an d Material Conference Newport, RI. Smith, R. N., Watson, P., and Topper, T. H. (1970). "A Stress-Strain Function for the Fatigue of Metal." Journal of Materials 5, 767-778. Sobol, I. M. (1993). "Sensitivity Analysis for Nonlinear Mathematical Models." Mathematical modeling and computational experiment 1, 407-144. Sobol, I. M. (2001). "Global Se nsitivity Indices for Nonlinea r Mathematical Models and Their Monte Carlo estimates." Mathematics and computers in simulation 55, 271-280. Soderberg, C. R. (1939). "Factor of Safety and Working Stress." Transaction of American Society Mechanical Engineers 52, 13-28. Taguchi, G. (1986). Introduction to Quality Engineering UNIPUB, White Plains,New York. Taguchi, G. (1987). System of Experimental Design: Engineering Methods to Optimize Quality and Minimize Cost UNIPUB/Kraus Internati onal, White Plains,New York. Tu, J. (1999). "Design Potential Concept fo r Reliability-Based Design Optimization," PH.D, Univ. of Iowa, Iowa City,Iowa. Tu, J., and Choi, K. K. (1997). "A Performance Measure Approach in Reliability Based Structural Optimization." Techincal Report R97-02 Center for Computer-Aided Design, University of Iowa, Iowa City,Iowa. Tu, J., Choi, K. K., and Park, Y. H. (1999) "A New Study on Reliability-Based Design Optimization." Journal of Mechanical Design 121(4), 557-564. PAGE 146 132 Tu, J., et al. (2001). "Design Potential Me thod for Robust System Parameter Design." AIAA Journal 39(4), 667-677. Van Keulen, F., Haftka, R. T., and Kim, N. H. (2004). "Review of Options for Structural Design Sensitivity Analysis. Part1: Linear Systems." Computer methods in Applied Mechanics and Engineering (to appear). Van Keulen, F., Liu, B., and Haftka, R. T. (2000). "Noise and Discontinuity Issues in Response Surface Based on Functions and Derivatives." Proceeding of the 41st AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference and Exhibit Atlanta, GA,USA. Vanderplaats, G. N. (2001). Numerical Optimization Techniques for Engineering Design Vanderplaats Research and Development, Inc. Vervenne, K. (2005). "Gradient-based A pproximate Design Optimization," Ph.D. Dissertation, Delft Un iversity of Technology. Villadsen, J., and Michelsen, M. L. (1978). Solution of Differentia l Equation Models by Polynomial Approximation Prentice Hall, Englewood Cliffs, New Jersey. Wang, L., and Kodiyalam, S. (2002). "An Effi cient Method for Probabilistic and Robust Design with Non-normal Distributions." 43rd AIAA/ASME/ASCE/AHS/ASC Structures,Structural Dynam ics,and Materials Conference Denver,Colorado. Wang, L. P., and Grandhi, R. V. (1995). "Imp roved 2-Point Function Approximations for Design Optimization." AIAA Journal 33(9), 1720-1727. Webster, M. D., Tatang, M. A., and McRae, G. J. (1996). "Application of the Probabilistic Collocation Method for an Un certainty Analysis of a Simple Ocean Model." 4 MIT Joint Program on the Scien ce and Policy of Global Change. Wohler, A. (1860). "Versuche uber die Fe stigkeit der Eisenbahnwagen-Achsen." Zeitschrift fur Bauwesen Wu, Y.-T. (1994). "Computational Methods fo r Efficient Structural Reliability and Reliability Sensitivity Analysis." AIAA Journal 32(8), 1717-1723. Wu, Y.-T., Shin, Y., Sues, R., and Cesare, M. (2001). "Safety-factor Based Approach for Probability-based Design Optimization." 42nd AIAA/ASME/ASCE/AHS/ASC Structure,Structural Dynamics,and Materials Conference Seattle,Washington. Wu, Y. T., and Wang, W. (1996). "A New Met hod for Efficient Reliability-Based Design Optimization." Probabilistic Mechanics & Structur al Reliability: Proceedings of the 7th Special Conference 274-277. PAGE 147 133 Wyss, G. D., and Jorgensen, K. H. (1998). "A User's Guide to LHS: Sandia's Latin Hypercube Sampling Software." Sandia National Laboratories, Albuquerque,New Mexico. Xiu, D., Lucor, D., Su, C.-H., and Karniadaki s, G. E. (2002). "Stochastic Modeling of Flow-Structure Interactions Using Generalized Polynomial Chaos." Journal of Fluids Engineering 124(1), 51-59. Xu, H., and Rahman, S. (2004). "A Gene ralized Dimension-reduction Method for Multimensional Integration in Stochastic Mechanics." International Journal for Numerical Methods in Engineering Article in Press. Xu, S., and Grandhi, R. V. (1998). "Effective Two-Point Function Approximation for Design Optimization." AIAA Journal 36(12), 2269-2275. Yang, R. J., and Gu, L. (2004). "Experien ce with Approximate Reliability-Based Optimization Methods." Structural and Multidisciplinary Optimization 26(1), 152-159. Youn, B. D. (2001). "Advances in Reliability-Based Design Optimization and Probability Analysis," PH.D, Univ. of Iowa, Iowa City,Iowa. Youn, B. D. (2001). "Advances in Reliability-Based Design Optimization and Probability Analysis," PH.D, Univ. of Iowa, Iowa City,Iowa. Youn, B. D., and Choi, K. K. (2004). "A New Response Surface Methodology for Reliability-Based Design Optimization." Computers & Structures 82, 241-256. Youn, B. D., Choi, K. K., and Park, Y. H. (2003). "Hybrid Analysis Method for Reliability-based Design Optimization." Journal of Mechanical Design 125, 221223. Yu, X., Choi, K. K., and Chang, K. H. (1997) "Reliability and Durability Based Design Sensitivity Analysis and Optimization." Center for Computer-Aided Design, University of Iowa, Iowa City,Iowa. Yu, X. M., Chang, K. H., and Choi, K. K. (1998). "Probabilistic Structural Durability Prediction." AIAA Journal 36(4), 628-637. PAGE 148 134 BIOGRAPHICAL SKETCH Haoyu Wang was born in Jiangyin, Chin a, on May 7th, 1976. He received his Bachelor of Science in mechanical design and manufacture in July 1998 from Nanjing University of Science and Technology, a nd a Master of Science in mechanical manufacture and automation from Southeast Un iversity in April 2002, both in China. His interest in conducting research motivated him to join the Univ ersity of Florida in August 2002, to pursue his Ph.D. degree in mechanical engineering. |