UFDC Home  myUFDC Home  Help 



Full Text  
PAGE 1 1 IDENTIFICATION OF STRUCTURESPECIFIC DAMAGE GROWTH PROPERTIES AND ITS IMPACT ON IMPROVED PROGNOSIS By ALEXANDRA COPPE 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 2010 PAGE 2 2 2010 Alexandra Coppe PAGE 3 3 To my mom Gisela Coppe PAGE 4 4 ACKNOWLEDGMENTS I thank my advisors Profs Nam Ho Kim and Raphael Haftka for their help and great guidance fro m the moment I applied to the position. They have shown such great wisdom, support and availability for both personal and professional matters I could not have wished for better advisors. I thank my family, in particular my parents, Gisela and Gervais Coppe, as well as my sister, Adiza Coppe, for their support and understanding of my choice to move to the United States. For the moral, emotional and financial support in that transition part of my life. I thank my friends, in particular Philippe Bergeron and Joel Mingam for their psychological and emotional support and their all time great advising and mentoring in some fairly tumultuous times. For having been able to stay in touch despite the distance and help me feel home in that new environment despite not being by my side. I thank Matthew Pais for his help as a fellow grad student and labmate and for the collaboration work that resulted in the Chapter 5 of that dissertation that would not have been possible without him. I thank my labmates for the comments they gave me along the way that made my work what it is now, the great discussions, their friendship and the support system I felt I could always rely on, they because a substitute family in that new environment Last but not least I thank Rick Ross, Thiagaraja Krishnamurthy and Tzikang Chen from NASA for their comments and question all along this project. I also thank the Air Force Off ice of Scientific Research Grant FA9550071 0018 and NASA Grant NNX08AC334 for funding this research. PAGE 5 5 TABLE OF CONTENTS page ACKNOWLEDGMENTS .................................................................................................. 4 LIST OF TABLES ............................................................................................................ 7 LIST OF FIGURES .......................................................................................................... 8 ABSTRACT ................................................................................................................... 12 CHAPTER 1 INTRODUCTION AND LITERATURE REVIEW ..................................................... 14 1.1 What is Structural Health Monitoring (SHM)? ................................................. 14 1.2 Advantages/Disadvantages of SHM Compared to Conventional Inspections ..................................................................................................... 16 1.3 Technical Issues in Preventing SHM f rom Being Popular .............................. 17 1.4 Different Type of Prognosis Methods (DataDriven vs. Model Based) ........... 18 1.5 Damage Growth Models ................................................................................. 20 1.6 Sources of Uncertainty ................................................................................... 22 1.7 Noise and Bias ............................................................................................... 23 1.8 Statistical Pr ediction of Remaining Useful Life ............................................... 24 1.9 Prognosis Using Least Square Method .......................................................... 25 1.10 Prognosis Using Bayesian Inference.............................................................. 25 1.11 Prognostic Metrics .......................................................................................... 27 1.12 Objective of Research .................................................................................... 28 1.13 How to Achieve th e Objective ........................................................................ 29 2 UNCERTAINTY REDUCTION OF DAMAGE GROWTH PROPERTIES USING BAYESIAN INFERENCE ........................................................................................ 31 2.1 Introduction .................................................................................................... 31 2.2 Damage Growth Model .................................................................................. 33 2.3 Statistical Characterization of Damage Growth Properties Using Bayesian Inference ........................................................................................................ 36 2.4 Numerical Applications ................................................................................... 45 2.5 Updating Damage Growth Exponent m .......................................................... 48 2.6 Updating Damage Growth Parameter C ......................................................... 54 2.7 Conclusions .................................................................................................... 58 3 UNCERTAINTY IDENTIFICATION OF DAMAGE GROWTH PARAMETERS USING HEALTH MONIT ORING DATA AND NONLINEAR REGRESSION .......... 60 3.1 Introduction .................................................................................................... 60 3.2 Uncertainty Quantification in NonLinear Least Square .................................. 62 PAGE 6 6 3.2.1 Uncertainty in the Linear Least Square Method ................................... 62 3.2.2 Uncertainty in the NonLinear Least Square Method ........................... 65 3.3 Identification of a Single Parameter ............................................................... 66 3.4 Identification of Multiple Parameters .............................................................. 70 3.5 Conclusions .................................................................................................... 73 4 LEAST SQUARESFILTERED BAYESIAN INFERENCE TO REDUCE UNCERTAINTY IN DAMAGE GROWTH PROPERTIES ........................................ 75 4.1 Introduction .................................................................................................... 75 4.2 Characterization of Damage Growth Properties using Bayesian Inference .... 76 4.3 Characterization of Damage Properti es Using Least Square Fit .................... 80 4.4 Least SquareFiltered Bayesian (LSFB) Method for Estimating RUL ............. 84 4.5 Conclusions .................................................................................................... 90 5 IDENTIFICATION OF EQUIVALENT DAMAGE GROWTH PARAMETERS WHEN USING WRONG RANGE OF STRESS INTENSITY FACTOR ................... 91 5.1 Introduction .................................................................................................... 91 5.2 Crack Growth Models ..................................................................................... 93 5.2.1 Damage Growth Model Used for LSFB and RUL Estimation .............. 93 5.2.2 Damage Growth Model Used for Measurement Data Generation ....... 94 5.3 Least Square Filtered Bayesian (LSFB) Method ............................................ 98 5.4 Results ......................................................................................................... 100 5.4.1 Center Crack in a Finite Plate ............................................................ 101 5.4.2 Edge Crack in a Finite Plate .............................................................. 105 5.4.3 Center Crack in a Plate with Holes .................................................... 108 5.5 Conclusions .................................................................................................. 112 6 CONCLUSIONS ................................................................................................... 114 APPENDIX A BAYESIAN INFERENCE USING DAMAGE GROWTH INFORMATION .............. 117 B LEAST SQUARE FILTERED BAYESIAN TO UPDATE JOINT PDF .................... 128 LIST OF REFERENCES ............................................................................................. 133 BIOGRAPHICAL SKETCH .......................................................................................... 140 PAGE 7 7 LIST OF TABLES Table page 2 1 Geometry, loading and damage growth parameters of 7075T651 aluminum alloy .................................................................................................................... 46 2 2 Statistical characteristics of final PDF of m with different combinations bias/noise. Illustration with one simulated set of measurements. ....................... 51 2 3 Statistical characteristics of updated PDF of C with different bias/noise ............ 56 PAGE 8 8 LIST OF FIGURES Figure page 2 1 Illustration of fuselage panel with a throughthe thickness crack ........................ 33 2 2 Illustration of the estimation of Paris model parameters using a log log plot of crack growth rate ................................................................................................ 37 2 3 Flowchart of likelihood calculation in Bayesian inference ................................... 42 2 4 Cumulative distribution function (CDF) of the RUL. ............................................ 45 2 5 Updated probability density functions of m ( mtrue = 3.8, Ctrue = 1.5E 10, b = 0 mm, V = 1 mm). Illustration with one synthetic set of measurements. ................ 49 2 6 Effect of bias on updated PDF of m ( mtrue = 3.8, Ctrue = 1.5E 10, V = 1 mm) Illustration with one simulated set of measurements. ......................................... 50 2 7 Effect of noise on the updated PDF of m ( mtrue = 3.8, Ctrue = 1.5E 10, b = 0 mm). Illustration with one simulated set of measurements. ................................ 51 2 8 Distribut ion (onesigma intervals) of 5 percentile (95% conservative) RUL obtained using 100 sets of measurements compared to the true RUL ............... 53 2 9 Distribution (onesigma intervals) of error between the t rue RUL and the maximum likelihood of the estimated RUL distribution ....................................... 54 2 10 Updated PDF of C ( mtrue = 3.8, Ctrue = 1.5E 10, b = 0 mm, V = 1 mm) .............. 55 2 11 Effect of bias on the updated PDF of C ( mtrue = 3.8, Ctrue = 1.5E 10, V = 1 mm) .................................................................................................................... 55 2 12 Effect of noise on final PDF of C ( mtrue = 3.8, Ctrue = 1.5E 10, b = 0 mm) ........... 56 2 13 Distribution (onesigma intervals) of 95% conservative RUL obtained using 100 sets of measurements compared to the true RUL ....................................... 57 2 14 Distribution (onesigma intervals) of error between the true RUL and the mean of the estimated RUL distribution .............................................................. 58 3 1 Comparison of the derived standard error with the simulated standard error ..... 69 3 2 Identified value of m ........................................................................................... 69 3 3 Comparison of the SE obtained using least square fit and Bayesian inference .. 70 3 4 Uncertainty in a0, b and m using least square fit. ............................................... 72 PAGE 9 9 3 5 Identified value of a0, b and m using least square fit. ........................................ 73 4 1 Distribution (mean one standard deviation) of the 5th percentile of RUL for b = +2mm and V = 1mm, using Bayesian inference .............................................. 78 4 2 Distribution (onesigm a intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for b = +2mm and V = 1mm, using Bayesian inference .......................................................................... 79 4 3 Distribution (mean one standard deviation obtained using 1,000 MCS simulations) of fitted results of a0, m and b /2. ..................................................... 82 4 4 5th percentile of RUL for b = 2mm and V = 1mm, using least square fit .............. 83 4 5 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for b = +2mm and V = 1mm, using least square fit, positive error corresponding to unconser vative estimates ............................................................................................................ 84 4 6 Comparison of the true, measured and fitted damage sizes for b = 2mm and V = 1mm. Illustration with one simulated set of measurements. ......................... 85 4 7 Fitted (dots) and measured (stars) damage at early and late stage in damage growth compared to the actual damage size (dotted line). ................................. 86 4 8 Distribution (mean ones standard deviation) of the 5th percentile of RUL for b = +2mm and V = 1mm, using the LSFB method .............................................. 87 4 9 Distribution (onesigma intervals) of error between the true RUL and the maximum li kelihood of the estimated RUL distribution for b = +2mm and V = 1mm, using LSFB, positive error corresponding to unconservative estimates ... 87 4 10 Comparison of the distribution (mean ones standard deviation) of the 5th percentile of RUL using the three methods ........................................................ 88 4 11 Comparison of the distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL using the three methods, positive error corresponding to unconservative estimates .................. 89 4 12 Comparison of the distribution (onesigma intervals) of m for one set of measurements usi ng the three methods. ........................................................... 89 5 1 Comparison of stress intensity response for some correction factors and crack sizes. Plate width is 200 mm. .................................................................... 95 5 2 Comparison of theoretical and XFEM crack growth curves using different plates heights, h for XFEM in order to validate the XFEM code and assess the validity of the theoretical model for the chosen plate geometry .................... 97 PAGE 10 10 5 3 Theoretical and XFEM prediction of f ( ) ............................................................. 98 5 4 A center crack in a finite plate. .......................................................................... 101 5 5 Correction factor for center crack ..................................................................... 102 5 6 Comparison of XFEM crack growth data with crack growth predicted from LSFB analysis. .................................................................................................. 102 5 7 Updated distribution of mLSFB using one set of data for a center crack in a finite plate. ........................................................................................................ 103 5 8 Distribution (mean one standard deviation) of 5th percentile of RUL for a center crack in a finite plate compared to the true RUL (black line) and the RUL obtained using the true m and K for a centre crack in an infinite plate (dark gray line). ................................................................................................. 104 5 9 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated R UL distribution for a center crack in a finite plate ......................................................................................................... 104 5 10 Edge crack in a finite plate ............................................................................... 105 5 11 Correction factor for edge crack ...................................................................... 105 5 12 Comparison of XFEM crack growth data with crack growth predicted from LSFB analysis. .................................................................................................. 106 5 13 Updated distribution of mLS FB using one set of data for an edge crack in a finite plate. ........................................................................................................ 107 5 14 Distribution (mean one standard deviation) of 5th percentile of RUL for an edge crack in a finite plate compared to the t rue RUL (black line) and the RUL obtained using the true m and K for a centre crack in an infinite plate (dark gray line). ................................................................................................. 107 5 15 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for an edge crack in a finite plate ......................................................................................................... 108 5 16 Center crack in a finite plate with holes. ........................................................... 108 5 17 Correction factor for plate with holes. ............................................................... 109 5 18 Comparison of XFEM crack growth data with crack growth predicted from LSFB analysis. .................................................................................................. 110 5 19 Updated distribution of mLSFB using one set of data for an center crack in a finite plate with holes. ....................................................................................... 110 PAGE 11 11 5 20 Distribution (mean one standard deviation) of 5th percentile of RUL for an edge crack in a finite plate compared to the true RUL (black line) and the RUL obtained using the true m and K for a centre crack in an infinite plate (dark gray line). ................................................................................................. 111 5 21 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for an center crack in a finite plate with holes ........................................................................................ 112 PAGE 12 12 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 IDENTIFICATION OF STRUCTURESPECIFIC DAMAGE GROWTH PROPERTIES AND ITS IMPACT ON IMPROVED PROGNOSIS By Alexandra Coppe December 2010 Chair: Nam Ho Kim Co chair: Raphael Tuvia Haftka Major: Mechanical Engineering Structural health moni toring (SHM) employs sensor data to monitor fatigueinduced damage growth in service. Damage growth information can be used to improve the characterization of the material properties that govern damage propagation for the structure being monitored, turning aircrafts into flying fatigue laboratories. Initially, t hese properties are often widely distributed between nominally identical structures because of differences in manufacturing processes and aging effects. SHM data, in particular, measured crack growth are used to narrow the distribution of damage growth parameters using the Bayesian inference technique. The improved accuracy in damage growth parameter s allows a more accurate prediction of the remaining useful life ( RUL ) of the monitored structural component. It can also help in predicting damage growth o f other similar components. In the absence of actual SHM data, we had to simulate measured data by applying error model s to damage sizes obtained using Paris law as a damage growth model. We consider t hat there are two kinds of errors, random noises resulting from the PAGE 13 13 measurement environment and a deterministic bias resulting from the sensor s calibration or modeling error. The Bayesian inference method is used for progressively reducing the uncertainty in structure specific damage growth parameters in spite of noise and bias in sensor measurements. However, the Bayesian inference method is computationally intensive due to uncertainty propagation in likelihood calculation, which results in a difficulty i n handling multi parameter identification, and may not be feasible to use with an extremely large number of measurement data. On the other hand, l east square fitting of damage growth parameters is efficient but it does not provide good statistical informa tion on the uncertainty in their estimates and in RUL estimates. It is in particular efficient to identify deterministic variables such as bias In the proposed research, we combined the two approaches by using the least square approach to filter the data and then perfor m ing Bayesian inference. The proposed approach is applied to crack growth in fuselage panels due to cycles of pressurization. It is shown that the proposed method rapidly converges to accurate damage parameters with much smaller uncertainties. Fairly accurate damage growth parameters were also obtained with measurement errors of 5mm. Using the identified damage parameters, it is shown that the 95% conservative RUL converges to the true RUL from the conservative side. PAGE 14 14 CHAPTER 1 INTRODUCTION AND LITERATURE REVIEW 1.1 What is S tructural H ealth M onitoring (SHM)? Structural health monitoring (SHM) is the process of implementing damage detection and characterization s ystems for structures. In SHM, d amage is defined as changes in the material and/or geometric properties of the structur e which includes cracks, corrosion degradation of material properties, etc. These changes may be dangerous because they affect the systems performance and safety In the case of fuselage panels in an air plane, for example, the SHM system can detect cracks in the panel s due to repeated pressurizations and provide warnings so that the panel s can be replaced before the cracks become unstable. There are many SHM technologies available [1 3] among which we have: e lectrical strain gauges e lectrical crack wires a coustic emission a coustoultrasonic l aser vibrometry c omparative vacuum monitoring (CVM) o p tical fib er Bragg grating (FBG), ME MS e lectromagnetic (eddy current sensors) foils s ensitive coating e nvironmental degradation monitoring sensors micro wave sensors and i maging ultrasonics These methods are applicable to different cases and with different levels of accuracy However, we will not discuss th em in detail here as SHM technologies are not the main focus of this work. The work presented here is not applied to a specific SHM method, but rather provides means to use the overall properties of SHM such as the low cost, frequency of measurements, and error model in order to improve our knowledge on the behavior of damage. We intend the current method to be applicable to any SHM method that allows us to estimate the size of damage. PAGE 15 15 The SHM process involves observations of a system over time using periodically sampled measurements from an array of sensors, extraction of damagesensitive features from these measurements, and statistical analysis of these features to determine the current health state of the system. For a long term, t he output of this process can be used to periodically update information regarding the ability of the structure to perform its intended function in light of the inevitable aging and degradation resulting from operational environments. After extreme events, such as earthquakes or blast loading, SHM can be used on civil structures such as buildings or bridges for rapid condition screening and aims to provide, in near real time, reliable information regarding the integrity of the structure. SHM is a technology under development which is implemented to prevent engineering structures to undergo catastrophic failure. It was first and foremost applied to civil structures such as buildings or bridges [4] using materials such as reinforced composites or reinforced concretes [5] Recently, SHM has received a lot of attention in aircraft structures, which is the focus of this research, and has a lot of promising results [1] Research related to SHM can be divided into two main categories, diagnosis and prognosis [6] Diagnos is deals with identifying the damage properties, such as damage location and size. Prognosis focuses on predicting the future behavior of the damage and the structure. The proposed research focuses on prognosis, especially when the measured data include uncertainties. Continual online SHM is based on dynamic processes through the diagnosis of early damage detection, then prognosis of health PAGE 16 16 status and remaining life this should allow a longer usage of a structure with the same level of reliability [7] 1.2 Advantages/ D isadvantages of SHM C ompared to C onventional I nspections The traditi onal inspection method is called preventive maintenance (PM) or manual inspection. I t involves the care and servicing by personnel for the purpose of maintaining equipment and facilities in satisfactory operating condition. It provides for systematic inspe ction, detection, and correction of incipient failures either before they occur or before they develop into major defects [8, 9] In this research, any kind of inspections that involve human intervention in the measurements is considered manual inspection. It also includes hybrid methods in which a human inspector performs the inspection using electronic sensors. While manual inspection is generally considered to be worthwhile, there are risks such as equipment failure or human error involved when performing PM, just as in any maintenance operations [10] On the other hand, SHM based maintenance (SHMM) removes or diminishes the human factor involved in inspection, and therefore, reduc es the risks and uncertainties involved [11] Even though SHM is not yet accepted by aircraft manufacturer s for current airplane structures and so far only laboratory level studies have been performed on SHM, some clear advantages have been shown for SHM over traditional manual inspectio n [12, 13] Among the se advantages are the low cost of damage assessment that leads to the possibility of performing damage assessment much more frequently I n fact the interval of damage assessment can be decreased from every 5,000 cycles in manual inspection to possibly every cycle. This capability of frequent assessment can compensate for the loss in accuracy compared to manual inspection [14] S ince SHM PAGE 17 17 can afford to perform assessment more frequently, it is possible to let the damage grow to a larger size before having to intervene and perform maintenance. P erforming frequent damage assessment will give us a large amount of data on the health state of the structure. Another advantage of SHM is that the error in damage size estimation is expected to be lower as a result of a reduction in human and environmental effects. This aspect compensates for the fact that SHM devices have a larger minimum detectable crack than preventive maintenance, whose order s of magnitude are of 5 mm and 0.3 mm, respectivel y. Another aspect to consider is cost SHM based maintenance requires additional installation cost as well as fuel cost related to the increase in weight resulting from the presence of SHM systems onboard, but on the other hand it drastically reduces the c ost of damage assessment We believe that so far the se advantages have not been used to their full extent, and we want here to get one step closer to that. 1.3 Technical I ssues in P reventing SHM from B eing P opular SHM is not yet a fully accepted and devel oped technology I t is still under research and need to make its proof. One of the main issues a ppears to be that engineers do not trust it yet, because it ha s its own issues that we are going to discuss below [15, 16] The first issue is the increase in weight resulting from the addition of the SHM system to the structure, including sensors, wires, and equipment. T his implies redesigning structures or modif ying existing structure s. Another issue is the reliability of the sensors and thei r various characteristics I n particular an important issue is related to aging of sensors T heir performance deteriorates over time s [17] PAGE 18 18 Health monitoring could be applied at many levels on an aircraft such as electronics, engine [18] and structure [19] and we are focused on the last one. But the application has limitations so far when it comes to actual structures because it is mostly focused on flat panels in the current laboratory studies However, a irplanes are composed of many curved panels with various connecting joints [20] 1.4 Different T ype of P rognosis M ethods ( D ataD riven vs. M odel B ased) P rognosis uses i nformation from collected data to predict future status of the system Prognosis techniques can be categorized based on the usage of information: (1) physics based, (2) datadriven, and (3) hybrid methods. The physics based methods, or model based methods [21] assume that the s tructures behavior can be predicted based on physical models of structure. Dynamic stochastic equations, lumpedparameter model s [22] and functional models [23] correspond to this category. In the case of SHM, crack growth model s [22, 24, 25] are often used for damage growth. The datadriven methods [26] include least square regression [27, 28] Gaussian process regression [29, 30] neural network, and relevance vector machine [29, 31] The d atadriven methods ha ve advantages when the system s are so complex that no simple physical model is available. The hybrid methods [32] use the advantages from both methods I t uses physical model s and identifies model parameter s using measured data, and includes such methods as partic le filters [33] and Bay esian techniques [34, 35] It is generally accepted that uncertainty is the most challenging part for prognosis [34, 36] Sources of uncertainty are initial state estimation, current state estimation, fai lure threshold, measurement, future load, future environment, and models. In order to PAGE 19 19 address the u ncertainty problem various methods have been proposed, such as confi dence intervals [37] relevance vec tor machine [29] Gaussian process regression [29, 30] an d particle fil t ers [33, 38] Although both physics based and datadriven methods use damage information measured by SHM systems, the major difference comes from the utilization of the physical model s that causes damage. For example, in t he case of monitoring crack in a panel, SHM systems measure its growth as a function of time. The datadriven methods use this information to predict the future behavior of the crack though extrapolation. On the other hand, the physics based model s utilize both loading and measured damage information to predict the future behavior of the crack. The physical model s relate the loading and crack growth. Since loading is the main source of uncertainty, the physics based methods can reduce uncertainty in prognos is significantly compared to datadriven methods. In addition, the physics based methods can predict the behavior of cracks in longer range of time. The datadriven methods are more appropriate for the case when detecting damage/fault requires a quick corr ecting action, such as an engine failure. In the case of crack growth in a panel, however, there is a significant time delay between the initial detection and replacing/repairing the crack. In such a case, the physics based model s are more appropriate because a more accurate maintenance period can be predicted. In this research, we will focus on the physics based methods in identifying parameters that control the physical model and predicting the service life before maintenance actions. PAGE 20 20 1.5 Damage G rowth M o dels The theory of f atigue crack growth has developed rapidly, starting in the early 1960s with the proposition of the initial Paris law [39] Paris applied Irvins crack stress analysis [40] to fatigue crack growth and opened a new way of modeling this important phenomen on. Paris method is nowadays widely accepted due to i ts simplicity. Since then knowledge about crack behavior has been constantly improved by the development of new models and the increased knowledge about factors involved. Although many developments in the theory have been presented during last four decades there are still remaining issues in fatigue crack growth. For example, most theories are based on constant cyclic loading conditions. When variable loadings are applied, it is still unclear theoretically and practically how the theory can incorporate wit h them. In the same context, when abnormally large loads are applied during the service life, the crack growth rate will change. Thus, the load h istory dependency of crack growth is an important and difficult issue in this field. Another important issue is the uncertain nature of crack growth physics. Since the cracks start from the micro scale level, it is extremely difficult to predict the crack growth pattern. M icro level crack growth is not yet fully understood and crack initiation time is another probl em because it is hard to identify Even in laboratory, two nominally identical specimens usually end up different crack propagation pattern. Thus, addressing uncertainty in crack propagation is an important issue. There have been no significant advances for the abovementioned issues, because multiple scales are involved, starting from atomistic scale to macro scale. There are recent developments in multi scale modeling for crack propagation, but it will take several years for them to be practical. PAGE 21 21 In additi on to the above mentioned issues, there have been a lot of research papers improving the original Paris model. Many researchers observed that the original model works only for limited range of load ratio, R which is the ratio between minimum and maximum s tress intensity factors. Although the effect of crack tip stress intensity factor s and their reversing plastic zone were already very well known but for some reason the effect of crack closure was ignored in the original development. This was initially obs erved by Elber [41, 42] in the late 1960s. Af ter that various service load damage models appeared in the 1970s. Among them, a finite element based model developed by Newman [43] was the most promisi ng approach. It is now generally accepted that this discrepancy is caused by crack closure during the cycle of applied loads. The most recent work using a partial crack closure model has been developed by Paris, Tada and Donald [39] The se issues have been addressed by introducing effective range of stress intensity factor, which accounts for plastic deformation at the crack tip as well as crack closure effect. In addition to the crack closure, different materials have different crack propagation slopes, which need to be obtained from experiments. The issue is that the unit of the crack growth equation becomes inconsistent when noninteger slopes are used. One of the recent developments from the original contributor, P. C. Paris, is to use a universal model that normalizes the range of stress intensity factor by elastic modulus [44] Since we are interested in illustrating how damage size estimation resulting from SHM inspection can be used to identify damage gr owth properties we choo se to use the original Paris law [39] E ven if it is far from being flawless it remains the simplest PAGE 22 22 approach. Even if more compl ex models, such as finite element based NASGRO software, can be used with more number of parameters a similar approach can be used. 1.6 Sources of U ncertainty Uncertainties are a key part of damage growth identification (measurement as well as prediction) because they make prognosis difficult T here are different types of uncertainties which are not all considered at this point. Some of those uncertainties are related to the loading history like pressurization amplitude or range of stress intensity factor [45] The se uncertaint ies are not considered in the current work because pressurization in the fus elage is relatively regular. However, they need to be addres sed in the case of crack growth in wing panels in which loading conditions are complicated. Another kind of uncertaint y is related to the damage itself, such as uncertainties in the initial damage size D ue to the variability in the manufacturing process every fuselage panel is assumed to have a small crack to start with, as specified in the Department of Defense Joint Service Specification Guide for aluminum alloys In this research, however, the uncertainties in initial crack size are not considered bec ause we are looking at a detected damage size from SHM This leads us to a related uncertainty which is the uncertainty in damage size measurement, which is a key part because it is directly related to the detection device [14] The uncertainty in measured damage size can actually be divided into two categories: a deterministic error that stands for the unknown bias that can exist in SHM systems, and a random noise that represents the variability in m easurements that can result from variability in the inspection environment. PAGE 23 23 We are also addressing uncertainties related to the parameters directing damage growth [46] which is the center of this work. These parameters depend on the damage growth model T he uncertainties of these parameters are in general large because they are meant to cover any possible case of manufacturing and aging conditions. However, for a specific panel, it is expected that the uncertainties of these parameters are very small, or even possibly they may have deterministic values. The proposed research focuses on finding these narrow distributions of panel specific parameters. Another goal of this work is to show that using damage growth model as an extrapolation tool allows one to use a simpler model despite t he fact that it might be erroneous. 1.7 Noise and B ias As described in Section 1. 6 the current technology of diagnosis and prognosis using sensor based SHM has issu es associated with uncertainties in sensor data, damage growth models, loading conditions and material and geometric properties. The first is related to identifying the current health status, while the others are related to predicting the health status in the future. Uncertainties in sensor data can be classified in two categories: systematic departure due to bias and random variability due to noise [47, 48] The former is caused by calibration error, sensor location and device error, while the latter is caused by measurement environment. Note that bias may tend to vary as the crack grows due to the nature of the error ; for example, the placement of the sensor with respect to the growth direction of the crack. However, we assume the bias to be constant over the entire life of the structure [47] Bias and noise are key variables in this work because of the absence of actual data that forces us to simulate data reflecting the reality of experimental data as close PAGE 24 24 as possible The problem of noise in data is already being addressed [25] but bias is often neglected in the modeling of measurement error. It is a different variable to address because it implies a dependence from one measurement to another at different times. Noise and bias can be dealt with in different ways; bias has the advantage that it can be identified as an additional variable in the model [49] The noise on the other hand cannot be identified due to its randomness but its effect can be reduced if there are enough data points [50] 1.8 Statistical P rediction of R emaining U seful L i fe Identifying damage parameters will allow us to improve our knowledge on a given panel or structure. The direct application of the improved knowledge is estimation of the remaining useful life (RUL), in other words, the amount of time left before the str ucture fails. Having better knowledge of damage growth parameters of a specific structure should lead to more accurate prognosis. The RUL is directly derived from the damage propagation law chosen, in this case Paris law [39] and in our case the failure c riterion is the size of damage [51] Various methods are used when it comes to statistical prognosis, such a maximum entropy approach [52] as well as Bayesian inference, Kalman filter and particle filter. In this work the idea is to estimate the distribution of RUL. For a given noise and bias model, a set of synthetic data can be c onstructed from the damage growth model. Since this set of data is not from real measurements, this process is repeated multiple times to estimate statistical properties. Since the predicted RUL becomes a distribution for a given set of data, the repeated simulation will yield distribution of RUL distribution. In order to present this distribution of distributions, the distribution of a specific percentile of the RUL distribution as a result of the presence of measurement errors will be plotted PAGE 25 25 That is, we calculate the 5th percentile of the RUL distribution which is the 95% confidence estimator of RUL in order to have a conservative estimate as a result of the uncertainties presented in Section 1.6. We end up with an estimation of the distribution of the 5th percentile of the distribution of RUL as a result of the noise in the data. The reason we choose that percentile rather than the mean is the result of the various uncertainties that lead to the need to have a conservative estimate that is still very close to the actual RUL. Having an estimate that is too conservative leads to not using the aircraft to its full capacity. 1.9 Prognosis U sing L east S quare M ethod When it comes to identifying parameters in a model or a system using experimental measurements the most straightforward and commonly used method is the least square fit [53] which is fairly simple and computationally efficient One of its clear advantages is that it allows one to increase the number of parameters to identify without increasing the computational cost significantly, unlike other methods such as Bayesian inference described below. Another advantage is that it deals very well with a large sample of data and can be used to reduce the effect of noise in a set of data by smoothing it out; it also can identify the bias if included in the model. Another way to filter out noise in a set of data is Kal man filter [54] or particle filters [38, 55] They are not investigated in the current work but are becoming more and more popular when it comes to the model parameter s identification. 1.10 Prognosis U sing Bayesian I nference Bayesian methods have become very popular over the past year s in various domains such as finance and engineering. Bayesian methods are useful in various PAGE 26 26 applications. These methods rest on the fundamental idea developed in Bayes theorem in the late 18th century The idea is to improve the knowledge on a theory by integrating new information obtained experimentally to the past/assumed knowledge [56] Below is a simple expression of Bayes theorem. ()  ) ) ) ( ( ( PBAP A A P P B B ( 1 1 ) In th e above equation A is the variable we are trying to identify and B is the observed information. P ( A  B ) is what we will call the posterior which i s the conditional probability of A knowing B P ( A ) is the prior (or assumed) probability or marginal probability of A. It is "prior" in the sense that it does not take into account any information about B P ( B ) is the prior or marginal probability of B and acts as a normalizing constant. And the last but not the least part of the equation is P ( B  A ) which is the conditional probability of B given A I t is also called likelihood function, and it is the key notion in Bayesian inference. Bayesian inference uses an estimation of the likelihood of a hypothesized value for the variable we are looking to identify using information observed experimentally. This process is repeated when additional information is obtained. A key notion in Bayesian inference is the likelihood function which represents the impact that the information has on the belief in the hypothesis. If it is likely that the information would be observed when the hypothesis under consideration is true, then that evidence will have a large effect on the likelihood function. Multiplying the prior probability of the hypothesis by the likelin ess would result in a larger posterior probability of the hypothesis given the evidence. Conversely, if it is unlikely that the information would be observed if the hypothesis under consideration is true; i.e. the prior probability that the information wo uld be PAGE 27 27 observed is high then the factor would reduce the posterior probability for the hypothesis Bayes' theorem therefore measures how much new evidence should alter a belief in a hypothesis. As expres sed earlier the key notion in Bayesian inference is the likelihood function and there are many interpretations for it Most of the time a hypothesis is made on the likelihood as being a given distribution [57] A frequent assumption is that it is distributed following a Gaussian distribution [58] Bayesian has many applications when it comes to damage growth, such as damage location and extent [57] material properties [59] crack initiation time [60] and probability of fail ure [61] Baye sian methods have the advantage of giv ing a rigorous formal approach to modeling uncertainties, but they are computationally expensive and require conditional independence in the variables of the problem Furthermore, there can be some issues with data acquisition [48] The uniqueness of the proposed research is that the likelihood is not assumed to be a particular distribution, such as Gaussian. Rather it is calculated numerically by propagating uncertainty using Monte Carlo simulations (MCS). MCS is used to propagate the uncertainties in crack growth estimation due to errors in measurements and uncertainty in material consta n ts Thi s allows us to do make one less assumption and as a result be closer from the actual behavior. 1.11 Prognostic M etrics With the development of SHM and prognostics in past few years various metrics have been developed to estimate the goodness of presented prognos tic methods, this resulted in a need for a standardized approach to compare the different approaches [6264] PAGE 28 28 Metrics can be sorted based on the type of prognostic such as cost benefit estimation, future behavior prediction or in our case endof life prediction. They can also be sorted based on the end user, researcher, designer, policy marker or operator. What we are interested here is endof life prediction without a specific end user [65] We are interested in assessing the accuracy of our method which means that we are looking to quantify the prec ision of our method [66, 67] Since we are dealing with distributed RUL estimation we made the decision to use the standard deviation and the mean of the error between the maximum likelihood of the distribution of RUL and the true RUL. 1.12 Objective of R esearch Compared to manual inspections, the accuracy of SHM is still poor. The minimum size of detectable damages of SHM is much larger than that of manual inspection methods. In addition, the measured data have the abovementioned noise and bias. Thus, the major challenge in SHM based prognosis is how to accurately predict the damage growth when the measured data include both noise and bias. However unlike manual inspection, SHM may provide frequent measurement s of damage, allowing us to follow damage growth. This in turn, should allow us to reduce the large uncertainty in the material properties that govern damage growth using SHM data. The objecti ve of the research is to reduce uncertainties in damage growth properties using noisy/biased sensor data in order to improve prognosis technology. The research will focus on identifying damage growth parameters of the monitored panel, followed by statistic ally predicting the conservative remaining useful life. PAGE 29 29 An application will be that by identifying the damage parameters and using the damage growth model as an extrapolation device, one can use a simple model that might be erroneous and still be able to predict RUL fairly accurately. 1.13 How to A chieve the O bjective The uncertainty in the damage growth parameters is normally large because of variability in manufacturing and ageing of the monitored structure d The idea is that SHM will give a large amount of data on damage sizes by monitoring the damage growth closely. Methods such as Bayesian inference and least square fit can then be used to identify the damage growth parameters using the se data. The approach is illustrated for a throughthethickness ce nter crack in an aircraft fuselage panel which grows through cycles of pressurization. A simple damage growth model, Paris model, with two damage parameters is utilized. However, more advanced damage growth models can also be used, which usually comes wit h more parameters. Using this simple model we aim to demonstrate that noisy SHM data can be used to identify the damage growth parameters in Paris law of a particular panel. This process can be viewed as turning every aircraft into a flying fatigue laboratory. Reducing uncertainty in damage growth parameters can reduce in turn the uncertainty in predicting remaining useful life (RUL); i.e., in prognosis. Since no actual inspection data are available, we have to simulate them which means that we have to define a model for measurement error and apply it to the data obtained using the model. A probabilistic approach using Bayesian statistics is first employed to progressively improve the accuracy of predicting damage parameters under noise and bias of sensor measurements. It is then compared to the least square fit identification method and the final goal is to combine the strengths of both methods in the least  PAGE 30 30 squarefiltered Bayesian (LSFB) method. The final goal being to identify the damage parameter distri bution with enough accuracy to improve the estimat ion of the RUL of the structure compared to the same estimate we would obtain using the handbook distribution of that parameter For an infinite plate with a center crack problem, we use analytical equation for the stress intensity factor. For more complex geometries, such as center or edge crack in a finite plate with/without holes, we use the extended finite element method to estimate the stress intensity factor and damage growth. Instead of developing mo re complex damage growth model with more number of parameters, we will show that simple damage growth model can still be used to predicting RUL by identifying effective damage growth parameters. The se equivalent parameters compensating for the modeling err ors should allow us to have a fairly accurate and conservative estimate of the RUL. PAGE 31 31 CHAPTER 2 UNCERTAINTY REDUCTIO N OF DAMAGE GROWTH PROPERTIES USING BAYESIAN INFERENCE 2.1 I ntroduction Bayesian methods have become popular in the past years to identify variables in a model For example, it can be used for fatigue crack by identifying parameters in damage growth models. As discussed in Chapter 1, Bayesian inference has mainly been used to extrapolate the crack behavior by updating the crack size distribution rather than the material properties that govern crack growth. Although the crack size distribution is important to diagnose the current health status, the crack growth properties of the material are important for the prognosis purpose. These properties are used to predict how fast a crack will grow in the future. Initially these properties are widely distributed due to variability in manufacturing process and to cover all possible cases of loading history A s a result they lead to overly conservative estimate of remaining useful life (RUL). The objective of this chapter is to charac terize the crack growth properties using Bayesian inference as an intermediate step toward predicting the RUL of the structure, which will allow us to improve our knowledge on the entire panel rather than the specific damage that is being monitored. Compar ed to manual inspections, the accuracy of structural health monitoring ( SHM ) is relatively poor. The minimum detectable size of damage using SHM is much larger than that of manual inspection methods. In addit ion, the measured data have noise and bias as di scussed in Chapter 1. Thus, the major challenge in SHM based prognosis is how to accurately predict the damage growth properties when the measured data include both noise and bias. Although noise is comm only discussed, bias is often ignored in the literature However, unlike manual inspection, SHM may PAGE 32 32 provide frequent measurements of damage sizes, allowing us to track damage growth. This in turn, may allow us to reduce the uncertainty in the material properties that govern damage growth. The uncertainty in these properties is normally large because of variability in manufacturing and aging of the structure. The main objective of this chapter is to demonstrate the reduction in uncertainty in these parameters using the abundance in SHM data, in spite of noise and bias by using Bayesian inference In other words we want to use the numerous data obtained from SHM in order to know more about damage growth behaviors of a specific panel A statistical approach using Bayesian inference is employed to progressively improve the accuracy of predicting damage growth parameters under noise and bias of sensor measurements. The proposed approach is demonstrated using a throughthethickness crack in an aircraft fuselage panel which grows through cycles of pressurization. A simple damage growth model by Paris [39] with two damage growth parameters is utilized. However, more advanced damage growth models can also be used, which usually come with more parameters. Using this simple model we aim to demonstrate that noisy SHM da ta can be used to identify the damage growth parameters of the monitored panel. This process can be viewed as turning every aircraft into a flying fatigue laboratory because accurate material parameters will be identified during the operation of aircraft Reducing uncertainty in damage growth parameters can reduce, in turn, the uncertainty in predicting RUL; i.e. prognosis. Improved knowledge on RUL can have practical consequences such as increased time between visual inspections, or a reduction in hardwar e testing when SHM is combined with manual inspection. PAGE 33 33 The chapter is organized as follows. In Section 2. 2, a simple damage growth model based on Paris model is presented. In S ection 2. 3 the measurement model used in this paper is introduced which shows how error in measurements due to SHM is added to the model presented in S ection 2. 2 It also presents how Bayesian inference is used to identify damage growth parameters. In S ection 2. 4 the model used to illustrate the method presented in this research is introduced. In Section 2. 5, the updating of damage parameter m is presented as well as the prognosis results resulting from it. In Section 2. 6, results similar to the one presented in S ection 2. 5 but obtained by updating the other damage parameter, C is presented, followed by conclusions in Section 2. 7. 2.2 Damage G rowth M odel Damage in a structure starts at the microstructure level, such as dislocations and gradually grows to the level of detectable macrocracks through nucleation and growth. Damage in t he microstructure level grows slowly, is difficult to detect, and is not dangerous for structural safety. Thus, SHM often focuses on macrocracks, which grow relatively quickly due to fatigue loadings. Figure 21 Illustratio n of fuselage panel with a t hroughthe thickness crack PAGE 34 34 In this work we consider a fatigue crack growth in a fuselage panel as illustrated in Figure 21 with initial half crack size ai subjected to fatigue loads with constant amplitude due to repeated pressurization of fuselage. The hoop stress varies between a maximum value of to a minimum value of zero in one flight ; o ne cycle of fatigue loading represents one flight. Like many other researchers (e.g., [68, 69] ), we use the damage growth model by Paris [39] as mda CK dN ( 2 1 ) where a is the half crack size in meters N the number of cycles (flights), da/dN the crack growth rate in meters/cycle K the range of stress intensity factor in MPa meter The above model has two damage growth parameters, C and m If the crack growth rate is plotted against the range of stress intensity factor in the loglog scale m is the slope of the curve and C is the y intercept at K = 1. K of stress intensity factor for a center cracked panel is calculated as a function of the stress range and the half crack length a in Eq. (2.2) and the hoop stress due to the pressure differential p is given by Eq. (2.3) Ka ( 2 2 ) pr t ( 2 3 ) where r is the fuselage radius, and t is the panel thickness. Equation (2.2) does not include a geometric correction factor due to finite size of the panel, and Eq. (2.3) does not include corrections due to the complexi ty of the fuselage construction, so that they are both approximate. PAGE 35 35 The number of cycles N of fatigue loading that makes a crack to grow from the initial half crack size ai to the final half crack aN can be obtained by integrating Eq. (2.1) 1/21/21/2N imm a Ni m aaa da N Cm Ca ( 2 4 ) Alternatively, the half crack size aN after N cycles of fatigue loading can be obtained by solving Eq. (2.4) for aN. 2 2 1/21 2m m m Nim aNC a ( 2 5 ) The panel will fail when the crack reaches a critical half crack size, aC. Here we assume that this critical crack size is when the stress intensity factor exceeds the fracture toughness KIC. This leads to the following expression for the critical crack size (again neglecting finite panel effects) 2 IC CK a ( 2 6 ) In the model presented above we have various sources of uncertainties. Uncertainty in fracture toughness results from testing variability. Another uncertainty is in the pressure differential, which results from issues related to its measurement during a flight as well as issues related to the unpredictability of that variable in the future flights. The uncer tainty in fracture toughness is not considered in this work because airline companies normally have a threshold crack size, smaller than the critical crack size, to repair the crack. Therefore, the end of life or RUL is defined as the number of cycles for which the crack grows up to the predetermined, critical crack size. The uncertainty in the pressure differential is not considered either as its level is very small PAGE 36 36 and has a minimal effect on crack growth. Another uncertainty is in damage size measurement Finally the last and most important uncertainties are related to the damage parameters, m and C which is the main focus of this research. 2.3 Statistical C haracterization of D amage G rowth P roperties U sing Bayesian I nference Damage growth parameters, C a nd m are critical factors to determine the growth of damage. These parameters are normally estimated by fitting fatigue test data under controlled laboratory environment. However, uncertainty in these parameters is normally large not only at a material level because of variability in manufacturing and aging of the specific panel, but also at a specimen level because of variability related to the testing process. However, a specific panel in an airplane may have a much narrower distribution of damage parameters, or even have deterministic values. We therefore use Bayesian inference to identify these panel specific parameters. As can be seen in Figure 22 the exponent m is the slope of the fatigue crack curve in the loglog scale, while the parameter C corresponds to y K = 1, of the fatigue curve. In order to make the presentation simple we assume that the parameter C has a known deterministic value, and thus, uncertainty is only in m However, the uncertainty in C ca n also be considered using the same concept. Fr o m the scattered data, the upper and lower bounds of m can be estimated using loglog plots of crack growth rate illustrated in Figure 22 Since the prior knowledge is limited, we a ssume that m is uniformly distributed between these two bounds. Then, the goal is to narrow the distribution of the exponent using Bayesian inference with measured damage sizes. PAGE 37 37 Figure 22 Illustration of the estimation of Par is model parameter s using a log log plot of crack growth rate Bayesian inference is based on the Bayes theorem of conditional probability. It is used to obtain the updated (also called posterior) probability of a random variable by using new information available for the variable. In this research, since the probability a is of interest, we use the following form of Bayes theorem [70] :  ini updt inilamfm fm lamfmdm ( 2 7 ) where fini is the initial (or prior) probability density function (PDF) of m fupdt is the updated (or posterior) PDF of m and l ( a m ) is called the likelihood function, which is defined as the likelihood of obtaining the measured damage size a for a given value of m The denominator in Eq. (2. 7) can be considered as a normalizing constant that makes fupdt satisfy the property of PDF; i.e., the area of the PDF should be one. Since fini is given or assumed, the most important step in Bayesian inference is to calculate the likelihood function, w hich determines the uncertainty structure of the posterior PAGE 38 38 distribution. It represents uncertainty in test data, which includes modeling error and measurement variability. In the literature, the likelihood function is often assumed to be Gaussian or to hav e another analytical expression [58] This assumption is made not because of physics but because of convenience. Since the posterior distribution strongly depends on the likelihood function, any assumptions on the lik elihood function may lead to errors in the posterior distribution. The main contribution of this chapter is to rigorously show the process of calculating the likelihood function by propagating uncertainties through the physical model. The likelihood funct ion is designed to integrate the information obtained from SHM measurement to the knowledge about the distribution of m. The physical interpretation of the likelihood is the PDF value of the true crack size at measured crack size for given m Although the true crack size would be a single value, it is considered to be randomly distributed in the viewpoint of measured crack size due to various uncertainties in the process. Thus, it is important to estimate the distribution of true crack size. In general, the measured crack size includes the effect of bias and noise of the sensor measurement as well as uncertainty in input loads In order to simulate the effect of noise and bias, a simple model, which is called synthetic data, is introduced. Let aN be the true half crack size, b the bias, and vN the noise at the current cycle N The measured half crack size, meas Na is then given as 22meas NNNaabv ( 2 8 ) The measurement bias b reflects a deterministic bias, such as calibration bias, while the noise vN reflects random noise. For subsequent simulated measurements, the PAGE 39 39 bias b rem ains constant, while we assume that the noise vN is uniformly distributed with in the range of [ V + V ] At a given SHM measurement the measured half crack size in Eq. (2.8) has the same distribution type with the noise. Since there is no information regarding the distribution of noise, in this work it is assumed to be uniformly distributed with mean at zero. However, it is also possible to model the noise as a normal distribution, which is close to the white noise. Thus, the measured crack size is also uniformly distributed and can be defined as : ~Uniform/2/2 ; /2/2meas NNNaabVabV ( 2 9 ) The quantity defined above only involves measurement error In general, however, the crack growth model may have modeling error, which is related to numerical simulation. I n order to calculate the likelihood function, we introduce a simulated half crack size, sim Na that involves a modeling error sim Ne, as sim sim NNNamaem ( 2 10) We use the superscript sim for crack size that includes the modeling error because it includes propagated uncertainty through numerical simulation. The simulated crack size depends on Paris parameters, m and C as well as the initial crack size. Since we only consider uncertainty in m and the initial damage size a0, Eq. (2.10) depends only on them The idea of calculating likelihood is to identify the damage growth parameter m by comparing the measured crack size meas Na with the simulated crack size sim Na with g iven m The difference between these two sizes can be defined as PAGE 40 40 sim meas NNdmama ( 2 11) If the PDFs of meas Na and sim Na are available, then the PDF of d can also be calculated. The likelihood l ( a m ) is then defined as the value of this PDF at d ( m ) = 0 Since the probability of that event is theoretic ally zero if we use MonteCarlo simulation ( MCS) to estimate the likelihood, we will get zero for the likelihood. Since MCS is a discrete process, it is not trivial to calculate the PDF value directly. Instead, we will use the probability of d with being a small constant as a definition of likelihood:  lamPd ( 2 12) Note that if the right hand side is divided by 2 and if approaches zero, then the likelihood becomes the value of PDF at d ( m ) = 0 In view of Eq. (2.7) since the posterior distribution will be normalized, the above definition works for likelihood although it is given in the form of probability. If we calculate l ( a m ) purely by sampling meas Na and sim Na th en the tolerance needs to be large enough to include enough samples to reduce sampling errors. On the other hand if is too large, we will incur errors due to nonlinearity in the likelihood function. In general, since the measurement error that controls meas Na is independent of the modeling error that controls sim Na separable sampling can be performed, and samples of d in Eq. (2.11) can be calculated by comparing all possible combinations from the two sets of samples [71] This can significantly improve computational efficiency since the analytical PDF of meas Na is available from Eq. (2.9) The PDF of sim Na is not available analytically, because it is obtained by propagating uncertainti es through the crack growth model. PAGE 41 41 The definition of likelihood in Eq. (2.12) can be expanded by  00 lamPdPdPd ( 2 13) Using conditional expectation on the second term on the right hand side we obtain 00 0 0sim N sim Nsimmeas NN simmeas simsim NN simNN a sim simsim measN simNN aPdPaa Paafada Fafada ( 2 14) w here f sim is the PDF of sim Na and Fmeas is the CDF of meas Na The last relation is obtained from the definition of CDF; i.e., by considering meas Na as the only random variable, sim meas N sim N measNPaaFa S imilarly the first term can be written as 00 0sim N sim Nsimmeas simsim NN simNN a sim simsim measN simNN aPd Paafada Fafada ( 2 15) Thus, by combining Eqs. (2.14) and (2.15) the likelihood can be written as 2sim N sim Nsim sim simsim measN measN simNN a sim simsim measNsimNN alamFaFafada fafada ( 2 16) where the central finite difference approximation is used in the second relation, which becomes accurate as 0 As explained before, since t he posterior PDF will be normalized, the coefficient 2 can be ignored. The above expression is in particular convenient for separable MCS because the analytical expression of fmeas is known, and fsim can be evaluated by propagating uncertainty through numerical simulation. Let M be the number of samples in M CS, the likelihood can then be calculated by PAGE 42 42 1 1sim Nsim simsim measNsimNN a M sim measNi ilamfafada fa M ( 2 17) Input data: 0 measa meas Na Discretize m For every m i : M samples of: 00 simmeas jaav w ith ~Uniform,jv VV 2 2 1 2 0.1 2i i im m m sim sim i Nm NC aa 1 ,1sim i mea i sNi Ma lamf M Figure 23 Flowchart of likelihood calculation in Bayesian inference As shown in Figure 21 first, the range of m is divided by 100 grids, and the likelihood is going to be calculated at each grid point of m Second, input random sampl es, such as initial crack size, noise and pressure, are generated according to their distribution types. These input random samples are propagated through the Paris model with a given value of m to produce M samples of crack size sim Na Third, the values of PDF sim measNfa are evaluated for all samples, whose average is used as the likelihood. The numerical experiments showed that M = 2,000 is enough to obtain a smooth distribution of the likelihood function. Note that li kelihood calculation is computationally intensive because Eq. (2.17) needs to be evaluated for every m in the range of Eq. (2.7) In addition, the Bayesian inference in Eq. (2.7) is repeated at every inspection cycle. PAGE 43 43 However, this process removes the assumption of likelihood distribution. In fact, due to nonlinear relation of the Paris model, the distribution of sim Na does not have any analytical distribution type. It is noted that in the above pr ocess of likelihood calculation, the bias, which is unknown, is ignored. This can cause inaccuracy in likelihood calculation. Theoretically, it is possible to identify both m and bias simultaneously. However, this can cause significant increase in computat ional cost because twodimensional grid (for m and bias) is required. In this chapter, the effect of bias is ignored in likelihood calculation, even if the measured data includes it. Once the posterior distribution of m is obtained from Bayesian inference, it can be used to estimate the RUL, which is the expected life from the current cycle to failure. In this research, failure is defined when the crack size reaches the critical crack size aC in Eq. (2.6) From Eq. (2.5) the RUL can be estimated by 1/2 1/21/2 m mtrue CN faa N Cm ( 2 18) Note that the RUL, Nf is also randomly distributed. Thus, it only makes sense to estimate the RUL in a probabilistic sense. The dis tributions of m and are given from Bayesian inference and Eq. (2.3), respectively. Although the true crack size, true Na should be a deterministic value, i t has to be estimated from the measured crack size, meas Na Thus, it needs to be considered as a random variable. For a given noise and bias model, the true crack size can be estimated by ~Uniform/2/2 ; /2/2true meas meas NNNaabVabV ( 2 19) PAGE 44 44 The distribution of RUL is calculated at every inspection cycle using MCS with 50,000 samples We generate 50,000 true damage sizes based on the measured damage size use Eq. (2.19) with b = 0 F or each of those sample we generate an m based on the updated distribution and then estimate the RUL using Eq. (2.18) Since predicting RUL is an extrapolation process, the input uncertainties are normally amplified in predicting uncertainty in RUL. In order to predict the RUL safely, we choose the 5th percentile as a conservative estimate of RUL. Figure 24 shows an example of predicted cumulative distribution function (CDF) of RUL of a panel under cyclic constant stress with initial half crack size of 10mm after 500 and 1,500 cycles. The total life of the structure is about 2,500 cycles if mean values of all input parameters are used. Since the useful life is consumed at every cycle, i t can be observed that the estimated mean values of RUL are about 2,000 and 1,000 after 500 and 1,500 cycles respectively However, these mean values would provide about 50% chance of over predicting the RUL for given uncertainties in input parameters and SHM measurements. In order to have 95% conservative estimate we choose the CDF value of 0.05, which corresponds to 1,000 and 81 0 cycles, respectively (two vertical lines in Figure 24 ) At 500 cycles the conservative RUL estimate has ratio of 0.5 (= 1,000/2,000) compared to the mean RUL, while at 1,500 cycles the ratio becomes 0.81 (= 810/1,000). This happens because the knowledge on damage growth parameter m is improved through Bayesian inference; i.e., the uncertainty in m is reduced. PAGE 45 45 A B Figure 24 Cumulative distribution function (CDF) of the RUL. A) CDF of RUL after updating the distribution of damage growth parameter m using Bayesian inference the horizontal line indicating the 95% conservative estimate of RUL. B ) detail of the CDF focused on the 95% conservative estimate of RUL 2.4 Numerical A pplications In this research, synthetic SHM measurements are utilized to demonstrate the process of Bayesian inference and predicting RUL. Depending on manufacturing and assembl y processes, the actual damage growth parameters for individual aircraft may be different. For a specific panel, we assume that there exists a true value of deterministic damage growth parameters ( m = 3.8 and C = 1.5E 10). In the following numerical simulation, the true damage will grow according to the true values of damage growth parameters. Then, the measured damage sizes are obtained by adding bias and noise from the true damage size. In order to simplify the presentation, we consider the distributions of m and C separately, which means that when we consider one variable as being uncertain we consider the other one as being known with the true value. Typical material properties for 7075T651 aluminum alloy are presented in Table 2 1 The applied fuselage pressure differential is 0.06 MPa obtained from Niu [45] and the range of stress is given by Eq. (2.3) Paris model parameters m and C are obtained using a crack growth rate plot published by Newman [43] Note that due to scatter of the 0 800 3,000 5,000 7,000 9,000 0 0.05 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Remaining useful lifeCDF of RUL RUL after 500 cycles RUL after 1,500 cycles 0 200 400 600 800 1000 0 0.05 0.1 Remaining useful lifeCDF of RUL RUL after 500 cycles RUL after 1,500 cycles PAGE 46 46 crack growth rate, the exponent m and log( C ) are assumed to be uniformly distributed between the lower and upper bounds. Table 21 Geometry, loading and damage growth parameters of 7075T651 a luminum alloy Property Distribution type and value Radius of fuselage, r (meter) Deterministic 3.25 Thickness of panel, t (meter) Determin istic 0.00248 Pressure differential, p (MPa) Lognormal (0.06,0.003) Fracture toughness, KIC (MPa meter ) Deterministic 30 True damage growth parameter, m true 3.8 True damage growth parameter, C true 1.5E 10 Initial distribution of m Uniform (3 .3, 4.3) Initial distribution of log(C) Uniform (log(5E 11), log(5E 10)) Noise, v (mm) Uniform ( V, +V), V = 1.0, or 3.0 Bias, b (mm) Deterministic, 2.0, 0.0, or 2.0 From the preliminary damage growth analysis, it was found that the distribution of pressure differential p has negligible effect on damage growth because the effect of randomness is averaged out. Thus, in the following analysis, the applied pressure differential is assumed to be deterministic, 0.06MPa, the mean of the distribution obtained from Niu [45] On the other hand, the uncertainty in the applied pressure di fferential is considered in the calculation of likelihood function. In general, the minimum size of detectable damage using SHM is much larger than that of manual inspection. Although different SHM techniques may have different minimum detectable size, we chose an initial half crack size of ai = 10 mm which is Lognormal (mean, standard deviation), modeled as constant in simulations Uniform (lower bound, upper bound) PAGE 47 47 large enough to be detected by most SHM methods. In addition, this size of damage will provide significant crack growth data between two consecutive inspections There are discussions of estimating i nitial crack size in the literature [72] but this is irrelevant to this research. The initial crack size in this research is when the crack is detected first time using SHM, which is much larger than the initial crack size of pristine material. In addition, the operating cycle N in this research starts from the cycle at which the crack is initially detected. In the following sections, two cases are considered. The first is updating parameter m and the second is updating parameter C Although it is wellknown that these two parameters are statistically correlated, we consider them separately in this chapter. The correlation can be identified if two parameters are updated simultaneously, which will be discussed in Appendix B First, one set of measured crack sizes is generated at every SHM measurement interval N by adding noise and bias to the true crack sizes that are calculated from the Paris model (see Eqs. (2.5) and (2.8) ). Then, at every measurement interval, Bayesian inference is used to update the PDF of m Once the PDF of m is available, Eq. (2.18) is used to estimate the distribution of RUL of which the 5th percentile is used as a conservative estimate of RUL. Since we used synthetic data by adding random noise, the result may vary with different sets of samples. Therefore, the above process is repeated with 100 sets of measurements and the mean plus and minus one standard deviation intervals are plotted. Below is a flow chart of Bayesian inference and predicting RUL. Figure 23 explains the detailed algorithm of the Bayesian procedure at the Nth inspection. The first step i s to generate M samples of true initial damage sizes, 0 sima PAGE 48 48 based on the initial measured damage size and Eq. (2.19) with b = 0 as for the RUL estimation, for each of the se simulated damage sizes we calculate the corresponding damage size after N cycles using Eq. (2.5) and we use it then in the likelihood calculation using Eq. (2.17) Note that in this case the effect of bias is ignored in the likelihood calculation. It can be added into the model as an additional variable. It has been ignored here, because its effect on the RUL estimation is negligible compared to the additional computational expense it would generate. In the flow chart in Figure 24, the measured crack size information is used in calculating the likelihood function. However, considering the fact that these parameters are related to damage growth, it is possible to calculate the likelihood function using the damage growth data, not the damage size date. Appendix A summarizes Bayesian inference procedure when damage growth information is used in likelihood calculation, along with updated distributions of damage growth parameters. In general, this approach is more sensitive to the noise in the measured data. 2.5 Updating D amage G rowth Exponent m In this section, we present how the distribution of the damage growth parameters for a fuselage panel can be narrowed using SHM measurement data and Bayesian inference. It is assumed that initially the panel has a 20mm throughthe thickness crack that is monitored by SHM system. The true crack grows according to damage growth parameters mtrue and Ctrue in Table 21 W e assume that SHM measurements are performed at every 100 cycles (i.e., N = 100). S ince the crack grows slowly and the noise and bias of measurements are in general large, too frequent measurements may not provide significant information about PAGE 49 49 the crack growth. The synthetic measured crack size data are generated by adding ra ndom noise and deterministic bias to the true crack size data. In Bayesian inference, only the measured crack size data are used. As a first example, we consider to update parameter m only, while the true value of parameter Ctrue is assumed to be known. St arting from the initial uniform distribution, the PDF of m is progressively updated using Bayesian inference with measured damage sizes. The range of m [3. 3 4 .3 ] is divided by 100 grids and PDF value is calculated at each grid point. The noise in crack detection is first assumed to be Uniform(1, +1) and the bias to be zero. SHM measurements are conducted until the crack reaches its half critical size aC defined in Eq. (2.6) that has a value about 42.7 mm. Figure 25 plots the updated PDFs of m at every 1,000 cycles for one set of measurements N ote that another set of measurements might lead to slightly different updated distribution but not significantly different It is clear that as the crack grows, the PDF of m becomes narrower and it converges to the true value of mtrue = 3.8. It is also noted that the convergence becomes faster as the crack size increases because the crack growth is faster for a larger crack. Figure 25 Updated probability density functions of m ( mtrue = 3.8, Ctrue = 1.5E 10, b = 0 mm, V = 1 mm). Illustration with one synthetic set of measurements. 3.6 3.7 3.8 3.9 4 0 5 10 15 20 25 30 35 40 mPDF of m a2,400 = 41.4 mm a0 = 10 mm a1,000 = 14.9 mm a2,000 = 27.8 mm PAGE 50 50 Figure 26 shows the effect of bias on the final updated PDF of m The noise in crack detection is assumed to be Uniform(1, +1) and two different biases are used: b = 2 mm and b = +2 mm It is clear that bias shifts the maximum likelihood point (the peak of PDF) from that of the true value; the negative bias overestimates the PDF of m while the positive bias underestimates it This can be expected because, for example, the positive bias yields larger crack sizes that grow in the rate of in smaller sized cracks. Since larger cracks should grow faster than smaller cracks, the identification process estimates that the crack grows slowly with a positive bias error. Figure 26 Effect of bias on updated PDF of m ( mtrue = 3.8, Ctrue = 1.5E 10, V = 1 mm) Illustration with one simulated set of measurements. Figure 27 shows the effect of noise on the PDF of m when b = 0 mm It is obvious that noise has an effect on the standard deviation but it does not shift the dis tribution as the bias does. The smaller the noise is, the narrower the final PDF of m In addition, it is clear in the case of large noise that the distribution of m is not symmetric. It is closer to the lognormal distribution. This may attribute to the fa ct that the Paris model is nonlinear process of crack growth. 3.6 3.7 3.8 3.9 4 0 5 10 15 20 25 30 35 40 45 mPDF of m b = +2 mm, V = 1 mm Initial distribution mtrue = 3.8 b = 2 mm, V = 1 mm PAGE 51 51 Figure 27 Effect of noise on the updated PDF of m ( mtrue = 3.8, Ctrue = 1.5E 10, b = 0 mm) Illustration with one simulated set of measurements. Table 21 shows statistical characteristics, such as the maximum likelihood, mean and standard deviation of PDF of m corresponding to Figure 26 and Figure 2 7 Note that in the figures each of these sets of characteristics are obtained using a single set of measurements It can be observed that the mean and maximum likelihood values are affected by the bias, while the standard deviation increases with a larger noise. As expected, a positive bias (true crack size is smaller than measured one) leads to underestimation of m Table 22 Statistical characteristics of final PDF of m with different combinations bias/noise. Illustration with one simulated set of measurements. Effect of noise Effect of bias Bias, noise (mm) b = 0, V = 1 b = 0, V = 3 b = 2, V = 1 b = +2, V = 1 Max. likelihood 3.79 3.80 3.85 3.76 Mean 3.80 3.81 3.85 3.76 Standard deviation 0.01 0.02 0.01 0.01 Once the PDF of m is obtained, it can be used to predict the RUL of the monitored panel. Since the PDF is updated at every SHM measurement, the predicted RUL will vary at every measurement interval N In predicting RUL, 50,000 samples of m true Na and are generated, and Eq. (2.18) is used to calculate samples of Nf. In order to have a safe prediction of RUL the 5th percentile of Nf samples is used as a conservative 3.6 3.7 3.8 3.9 4 0 5 10 15 20 25 30 35 40 45 mPDF of m Initial distribution b = 0 mm, V = 3 mm mtrue = 3.8 b = 0 mm, V = 1 mm PAGE 52 52 estimate of RUL. I n order to estimate the accuracy of the method we also calculate the error between the maximum likelihood of the estimated distribution of RUL and true RUL. Since we used synthetic data by adding random noise, the result may vary with different sets of dat a. Thus, the above process is repeated with 100 sets of measurement data and mean one standard deviation intervals are plotted for the 5th percentile as well as the error. Figure 28 shows the conservative intervals of RUL with two different combinations of noise and bias. These combinations correspond to extreme cases; the most and least conservative estimates of RUL. In order to compare the predicted RUL with the true one, the true RUL is also plotted in the figure. Since the u seful life is consumed at each cycle, the true RUL is a diagonal line. If the predicted RUL is above this diagonal line, it means the predicted RUL is longer than the true one, which is unconservative. Therefore, a conservative prediction should stay below the line. On the other hand, accuracy is measured by how close the prediction to the true line. Note that initially the difference between the true and predicted RULs is significant because uncertainty is large at early stage. However, the predicted RUL converges to the true one as more updates are performed. In addition, the variability of estimated RUL is also gradually reduced. Thus, it can be concluded that the proposed Bayesian inference can estimate panel specific damage growth parameters as well as the RUL. It is also observed that the positive bias yield more unconservative prediction because it underestimate the damage growth parameter m The actual crack grows faster than the estimated one. PAGE 53 53 Figure 28 Distribution (o nesigma intervals) of 5 percentile ( 95% conservative) RUL obtained using 100 sets of measurements compared to the true RUL Figure 29 shows the distribution of the error between the maximum likelihood of the estim ated RUL distribution and the true RUL with two different combinations of noise and bias as defined below maxtrue ffENN ( 2 20) w here max fN is the maximum likelihood of the estimated distribution of RUL and true fN is the true RUL. Note that positive values of the error correspond to unconservative es timates. A single value of E is calculated at each set of SHM data. After repeating this calculation with 100 sets of SHM data, the onesigma confidence interval of E is plotted in Figure 29. It can be observed that at first the m aximum likelihood estimat e is unconservative, this can be results from the fact that the large values of m are least likely to generate the damage size observed. The distribution converges from upper l imit of the initial distribution, m = 4.3, as can be observed in Figure 25 Towards the end of life of the structure, the maximum likelihood estimate converges to the true value and it confirms the conclusion drawn from the previous figure that the distribution of PAGE 54 54 RUL becomes narrower as m is better characterized. It can also be concluded that positive bias lead to more unconservative results. Figure 29 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution 2 .6 Updating D amage G rowth P arameter C In this section, as a second example, the damage growth parameter C is updated, while m = mtrue is used, starting from the initial distribution given in Table 21 Since the Paris model is lin ear in loglog scale and since is the y intercept at K = 1 log( C ) is updated instead of C The updating process is the same as updating the distribution of m described in Section 5 with the same type of likelihood function and the same noise and bias. Starting from the initial uniform distribution, the PDF of log( C ) is progressively updated using Bayesian inference with measured damage sizes. The noise in crack detection is assumed to be Uniform(1, +1) and the bias to be zero. SHM measurements are conducted until the crack reaches its half critical size, aC, defined in Eq. (2.6) that has a value about 42.7 mm. Figure 210 plots th e updated PDFs of C at every 1,000 cycles for a single set of measurements, note that as for m other sets of measurements might lead to slightly but not drastically different updated distributions It PAGE 55 55 is clear that as the crack grows, the PDF of C becomes narrower and it converges to the true value of Ctrue = 1.5E 10. It is noted that the convergence becomes faster as the crack size increases because the crack growth is faster for a larger crack. Figure 210 Updated PDF of C ( mtrue = 3.8, Ctrue = 1.5E 10, b = 0 mm, V = 1 mm) The effects of noise and bias turn out to be similar to the case of updating m Figure 211 shows the effect of bias on the final updated PDF of C The noise in crack detection i s assumed to be Uniform( 1, +1) and two different biases are used: b = 2 mm and b = +2 mm As for m bias appears to shift the maximum likelihood point from that of the true value. T he negative bias overestimates the PDF of C while the positive bias underestimates it. Figure 211 Effect of bias on the updated PDF of C ( mtrue = 3.8, Ctrue = 1.5E 10, V = 1 mm) 1E10 1.5E10 2E10 3E10 0 1 2 3 4 5 6 7 8 9 10 x 1010 CPDF of C a0 = 10 mm a1,000 = 14.9 mm a2,000 = 27.8 mm a2,400 = 42.4 mm 1E10 1.5E10 2E10 3E10 0 1 2 3 4 5 6 7 8 9 10 x 1010 CPDF of C Initial distribution b = 2 mm, V = 1 mm Ctrue = 1.5E10 b = +2 mm, V = 1 mm PAGE 56 56 Figure 212 shows the effect of noise on the PDF of C when b = 0 mm It is obvious that noi se increases the standard deviation but it does not shift the distribution as the bias does. The smaller the noise is, the narrower the final PDF of C Figure 212 Effect of noise on final PDF of C ( mtrue = 3.8, Ctrue = 1.5E 1 0, b = 0 mm) Table 23 shows statistical characteristics, such as the maximum likelihood, mean and standard deviation of PDF of C corresponding to Figure 211 and Figure 212. N ote that as previously, each of these sets of characteristics are obtained using a single set of measurements. It can be observed that the mean is affected by the bias, while the standard deviation is large with a large noise. As expected, a positive bias (true c rack size is smaller than measured) leads to underestimation of C From Figure 211 and Figure 212 it can be concluded that the effects of noise and bias are similar to both updating m and C Table 2 3 Statistical characteristics of updated PDF of C with different bias/noise Effect of noise Effect of bias Bias, noise (mm) b = 0, V = 1 b = 0, V = 3 b = 2, V = 1 b = +2, V = 1 Max. likelihood 1.5E 10 1.5E 10 1.7E 10 1.4 E 10 Mean 1.5E 10 1.6E 10 1.7E 10 1.4E 10 Standard deviation 3.8E 12 1.3E 11 5.1E 12 4.5E 12 1E10 1.5E10 2E10 3E10 0 1 2 3 4 5 6 7 8 9 10 x 1010 CPDF of C Initial distribution Ctrue = 1.5E10 b = 0 mm, V = 1 mm b = 0 mm, V = 3 mm PAGE 57 57 Similar to the case of updating m the distribution of the 5th percentile of predicted RULs from the updated C at every measurement interval are plotted in Figure 213. Again, it is clear that the conservative estimate of RUL converges to the true RUL from the safe side. In addition, the variability of estimated RUL is also gradually reduced. Figure 213 D istribution (onesigma intervals) of 95% conservative RUL obtained using 100 sets of measurements compared to the true RUL Figure 214 shows the distribution of the error between the true RUL and the mean of the estimated distribution of RUL. As previously the behavior observed can be explained by the convergence of the distribution of C from the upper limit of the initial distribution. It can also be observed that the distribution of RUL narrows as the distribution of C is identif ied and it converges to the true RUL. Thus it can be concluded that the proposed Bayesian inference can estimate panel specific damage growth parameters as well as can predict the RUL. PAGE 58 58 Figure 214 Distribution (onesigma int ervals) of error between the true RUL and the mean of the estimated RUL distribution 2.7 Conclusions In this chapter, the Bayesian inference method is employed to identify panel specific damage growth parameters using damage sizes measured from SHM sensors The actual measurement environment is modeled by introducing deterministic bias and random noise. The likelihood function is calculated by comparing measured crack size with simulated crack size, which requires uncertainty propagation though the physics model that governs the crack growth. Due to many uncertainties involved, the RUL is predicted statistically and used 95% conservative estimation. Through numerical examples, it is shown that the probability distributions of the two damage growth parameters m and C are effectively narrowed. However, the effect of bias and bias remains, and they affect the identification of true damage growth parameters. The identification is sensitive to the error in the data, two sets of measurements will lead to two esti mates of RUL. PAGE 59 59 The identified distributions of parameters are used to estimate the RUL with 95% confidence. In all combined cases with noise and bias, the proposed method converges to the true RUL from the conservative side. In the more general approach, it is possible to update both m and C using their joint PDF. In addition, the unknown bias can also be considered as uncertain variable and can be updated together. However, as the number of variables increases, the computational cost increases significantly because the proposed method is based on MCS for uncertainty propagation. PAGE 60 60 CHAPTER 3 UNCERTAINTY IDENTIFI CATION OF DAMAGE GRO WTH PARAMETERS USING HEALTH MONITORING DATA AND NONLINEAR REGRESSION 3.1 Introduction In the previous chapter, it is shown that S HM can not only provide damage diagnosis but also predict the remaining useful life (RUL) by identifying damage growth parameters. Bayesian inference has been used to reduce uncertainty in the damage growth parameters using measured damage size information. Bayesian inference is a powerful method of quantifying uncertainty in the model parameters with experimental data. It can take into account the prior knowledge on the unknown parameters and improve it using experimental observations. However, in the cas e of SHM, the advantage of the prior information can be overpowered by the amount of data available. That is, the effect of prior information becomes insignificant when numerous SHM data are used in Bayesian inference. In addition, when many parameters are updated simultaneously, Bayesian inference becomes computationally expensive due to multi dimensional integration. On the other hand, the traditional linear least square method [73] can be used to identify deterministic parameters when the model is a linear function of the parameters and uncertainty is Gaussian. This method is in particular powerful when many data are available, which is the case for SHM data. Unlike Bayesian inference, this method cannot include prior information. By assuming that the noise in the experimental data is Gaussian, it is possible to estimate the uncertainty in the identified parameters, which is also Gaussian. When the physical model is a nonlinear function of model parameters and uncertainty is not Gaussian, it is not straightforward to apply for uncertainty PAGE 61 61 quantification in the least square method. As will be shown in the numerical example, the damage growth in aircraft structures is governed by a nonlinear equation whose parameters need to be identified. In this chapter, we proposed a linear perturbation concept to quantify uncertainty in the nonlinear least square method. First, the nonlinear mathematical programming problem is solved to find the model parameters that minimize the square error between the model and experiment. Then, the nonlinear model is linearized about the identified parameter values, and then the uncertainty quantification for the linear least square method is applied. This approach can introduce two errors into the estimate of the uncertainty: (1) linearization error and (2) error associated with assumption of Gaussian noise. In addition, it is assumed that noises at different experiments are uncorrelated. This chapter presents a nonlinear least square method to identify the same damage growth parameters as in the previous chapter using a throughthethickness crack in an aircraft fuselage panel which grows through cycles of pressurization. We present here two least square fit probl ems based on the same model, one identifies a single variable while the other one identifies three variables simultaneously The objective is to examine the accuracy of uncertainty quantification using nonlinear least squares. The uncertainty is derived analytically and compared to a Monte Carlo estimate to examine its accuracy for both problems presented. For the onevariable problem it is then compared to the uncertainty obtained using Bayesian inference. The chapter is organized as follows. In Section 3 2, the derivation of the uncertainty quantification for the general least square methods, linear and non linear is presented. In Section 3. 3 the one variable identification problem is discussed, showing PAGE 62 62 the derivation of the analytical uncertainty quantif ication, as well as results comparing the calculated uncertainty with the one obtained using Monte Carlo simulations. In Section 3. 4, the threevariable problem is discussed and used to introduce the issues related to the uncertainty estimation resulting f rom potentially correlated variable and highly nonlinear problems. Concluding remarks are presented in Section 3. 5. 3.2 Uncertainty Quantification in NonLinear Least Square The least square method is commonly used for identifying unknown parameters of a physical model using experimental observations, which normally include random noise. Thus, if the experiment is repeated, it is likely that different values of the parameters may be identified. In this section, a method of calculating the distribution of t he identified parameters in the nonlinear least square method will be presented. In order to make the presentation easy to understand, estimation of parameter uncertainty in the linear least square method is discussed first, followed by that of the nonlinear least square method. 3.2.1 Uncertainty in the Linear Least Square Method In re gression, the response function () yt is approximated by an analytical function (,) yt with vector of parameters whose dimension is di ) m( n : yt yt ( 3 1 ) where is the approximation error. The objective of regression is to estimate the parameters b ased on ny data, whi ch are given in the form of ,, 1,...,ii ytyin th at may contain measurement noise. In regression the parameters are estimated by minimizing the sum of the squares of the discrepancies between the measurements and PAGE 63 63 (,) yt The regression model is called linear when the approximate function is a linear function of as 1 ,n ii it yt ( 3 2 ) where it are basis functions. In general, the exact values of can only be found when the number of experimental data is infinite and the noise has zero mean (unbiased measurements). With finite ny, the values are only estimate, which will be denoted by b in this chapter. The vector of errors (discrepancies) can be written as 1121 1 11 1 1222 2 22 2 12yy y yy yn n nn n nnnnttt ey b ttt ey b ey b ttt ( 3 3 ) Or, symbolically eyX.b ( 3 4 ) The parameters b are estimated by minimizing the root meansquare error defined as 1T RMS yn e ee ( 3 5 ) After substituting Eq. ( 3.4) into Eq. ( 3.5 ) and minimizing the root meansquare error, the following linear regression equation is obtained: TTXX.bXy ( 3 6 ) which can be solved for the estimate b of parameters. PAGE 64 64 Because the experimental data includes random noise, the estimated parameters will be different for different sets of experimental data. The objective is to estimate the uncertainty in the estimated parameters due to the random noise in the experimental data. The uncertainty in the parameters can be found assuming that the random noise has a Gaussian distribution with standard deviation (STD) of ; i.e., 2~0, N and that the noise at different measurements is uncorrelated. The unbiased estimate of the standard deviation can be obtained from 2T yn n ee ( 3 7 ) The sensitivity of estimated parameters with respect to small differences in data can be calculated using the covariance matrix of b defined as TEE bbbbb ( 3 8 ) where E b is the expected value of b Using Eq. ( 3.6) the covariance matrix can be obtained as 1 2 T bXX ( 3 9 ) The diagonal components of b is the square of standard deviation of b which represents a measure of the sensitivity of estimated parameters with respect to the noise. Since the standard deviation of the noise is unknown in advance, its estimate in Eq. ( 3.7) can be used. Thus, the standard error (SE) of parameter bi can be obtained by 1T i iiSEb X X ( 3 10) The above standard error is indeed the estimate of the standard deviation of bi. PAGE 65 65 3.2.2 Uncertainty in the NonLinear Least Square Method Physical models cannot often be represented as a linear combination of unknown parameters as in Eq. ( 3.2) Thus, instead of solving a linear regression Eq. ( 3.6) a nonlinear optimization problem needs to be solved to minimize the root meansquare error in Eq. ( 3.5) In this work, Matlab lsqnonlin funct ion is used to solve the nonlinear regression problems. The identified parameters are denoted by b* It is noted that a tight convergence criterion should be used in lsqnonlin function because it is possible that some parameters are insensitive to the error. In the following, the nonlinear equation will be linearized at the optimum point with respect to the identified parameters in order to estimate the uncertainty in the parameters. To linearize the problem, the first order Taylor series expansion method c an be used at bb b ,,,i i iy ytyt tb b bb*b* ( 3 11) By moving yt b* to the left hand side, the equation for the residual can be defined as ,,i i iy yyt tb b r b* b* ( 3 12) Equation ( 3.12) can be considered as a linear least square problem with unknown parameters b and the gradients ib y becomes the basis vector in Eq. (3.2) Thus, the uncertainty in parameters can be calculated using the same procedure described in the previous subsection. For that pur pose, the vector of basis functions can be written as PAGE 66 66 2 22 2 111 1 2 1 2 1 yyyn n nnn nyyy yyy yy bbb bbb b y bb X ( 3 13) Then, Eq. ( 3.10) can be used to estimate the standard error of b which can also be considered as the standard error of b if the problem is linear. Due to the nonlinearity, the st andard error of b will be different from that of b However, if the nonlinearity is small, or if the uncertainty in b is small, then the difference between them will be small. In order to verify the proposed uncertainty quantification method of nonlinear least square method, Monte Carlo simulation can be used to estimate the uncertainty in the identified parameters. In this approach, it is assumed that the experiments are repeated many times, and the parameters are identified for each experiment, from which the distribution of identified parameters can be estimated. In the following two sections we respectively illustrate how the uncertainty resulting from least square fit compares to the one obtained using Bayesian inference using a one variable problem and how the uncertainty quantification behaves in the case of a multi variable, highly nonlinear problem. 3.3 Identification of a Single Parameter The problem here is to identify the damage growth parameter using measured damage sizes. We fit a damage gr owth law, in this case Paris law in Eq. ( 2.1) to a set of damage size measurement s and identify the Paris law exponent m PAGE 67 67 In th e damage growth model, several unknown par ameters are involved. First, the damage growth parameters, C and m need to be identified. In addition, the initial damage size, a0, is often unknown, and needs to be identified too. Bayesian inference is used in Chapter 2 to identify the unknown parameter s with measured damage size. However, due to computational challenge in Bayesian inference, we only identified either unknown damage growth parameter, m or C separately by assuming all other parameters are known. In this section, uncertainty quantification of nonlinear least square method is used to identify the uncertainty in m and compare it with that of Bayesian inference. If we define ameas as the measured data and amodel as the damage size corresponding to aN in the damage growth model in Eq. (2.5) The least square fit problem can then be stated as minimizing the L2norm of error r with meas iiyN ra m ( 3 14) with model iya ( 3 1 5 ) The least square fit problem can then be defined as 2minm i ir ( 3 16) Let the identified Paris exponent from the above nonlinear least square problem be m *, the standard error ( SE) in m can then be obtained (see Eqs. ( 3.7) to ( 3.10) ) as 1 11TXX SEm ( 3 17) With the Nx1 matrix X defined as PAGE 68 68 y m X m ( 3 18) with 2 2 1 1yn i ie N ( 3 19) and *sim ii ida rm dm e ( 3 20) As defined in Section 2.2 the measured data, ameas, are actually simulated by apply ing an error model defined in Eq. ( 2.8) to the modeled damage size, amodel. It has to be noted here that when fitting only m we assume the bias to be zero and that the initial damage size is known accurately. Although this assumption is not practical, this can make us to compare the uncertainty in parameter m with that of Bayesian inference. Since the uncertainty in identifying m results from the noise in the data, we can also quantify it using Monte Carlo Simulation (MCS) by simulating 1,000 sets of data, performing a fit for each and then calculating the standard deviation for those data. We can afford to do this here because we are simulating the measured data using the error model. Figure 31 shows the comparison between the standard error and the standard deviation from the MCS. PAGE 69 69 Figure 31 Comparison of the derived standard error with the simulated standard error Figure 32 shows the identified value of m corresponding to the uncertainty illustrated in Figure 31 it can be observed that like for Bayesian inference, neglecting the presence of bias leads to a slight underestimation of the damage growth parameter m Figure 32 Identified value of m It can be observed that the derived standard error (solid line) fits very well the estimated standard deviation (dashed line); the plots can hardly be distinguished. The next step is then to compare the standard error of m to the standard deviation from Bayesian inference. In order to do that we use the Bayesian inference method 0 500 1000 1500 2000 2500 105 104 103 102 101 Number of cycles at inspectionUncertainty in m Standard error Standard deviation 0 500 1000 1500 2000 2500 3.73 3.74 3.75 3.76 3.77 3.78 3.79 3.8 3.81 Number of cycles at inspectionIdentified value of m Identified value of m True m PAGE 70 70 developed and presented in Chapter 2. Fi gure 33 illustrates the comparison between the standard deviation resulting from least square fit identification and the one resulting from Bayesian inference. Figure 33 Comparison of the SE obtained using least square fit and Bayesian inference It can be observed that least square fit is able to identify m much more accurately than Bayesian inference; at least 3 times more accurately. This can be explained by the fact that Bayesian inference is performed every 100 cycles, w hile least square fit uses data at every cycle. I f the same intervals were used with least square fit then a similar standard deviation would be achieved. 3.4 Identification of M ultiple P arameters As mentioned before, Bayesian inference becomes computationally expensive when multiple parameters are identified simultaneously However, the proposed method is straightforward for quantifying uncertainty of multiple parameters. Unlike the idealized situation of the previous section, we now assume that the initi al crack a0 and the bias b also need to be identified simultaneously Then the threevariable problem can then be defined as 0 500 1000 1500 2000 2500 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 Number of cycles at inspectionUncertainty in m Least square fit Bayesian inference PAGE 71 71 0 with ,,meas iiyN Ra mab ( 3 21) with model ib ya ( 3 22) The least square fit problem can then be defined as 0,, 2minmab i iR ( 3 23) The derivation of the SE is the same as previously 1 11TXX SEm ( 3 24) 2 1 2 0TXX SEa ( 3 25) 3 1 3TXX SEb ( 3 26) W ith the Nx3 matrix X defined as 0 *** X m yyy ab ( 3 27) And 2 2 3i ie N ( 3 28) As described in the previous section, we compare the analytical estimate of SE to the simulated STD obtained using MCS. This can be found in Figure 34 It can be observed that in this case the derived standard error does not match the simulated standard deviation very well for the first 1,000 cycles. The least square fit method overestimates the uncertainty in parameters. However, it converges to the correct standard deviation beyond 1,500 cycles. There are several explanations for this PAGE 72 72 discrepancy. First, the least square method predicts a larger standard error because not many data are available at the early stage. The nonlinearity of the nonlinear least square problem can also contribute to the discrepancy. Another aspect is the close relationship between a0 and b where these two parameters compensate for each other and lead to an ill conditioned XTX matrix. As the damage grows, it becomes easier to distinguish between the effects of a0 and b which is why the standard errors converge to the standard deviation from MCS. A B C Figure 34 Uncertainty in a0, b and m using least square fit A) Uncertainty in a0. B) Uncertainty in b C) Uncertainty in m. Figure 35 shows the identified value of the parameters a0, m and b at each inspection. It can be observed that the identification converges to true value of the 0 500 1000 1500 2000 2500 105 104 103 102 101 100 Number of cycles at inspection Uncertainty in a0 Standard error Standard deviation 0 500 1000 1500 2000 2500 105 104 103 102 101 100 Number f cycles at inspectionUncertainty in b Standard error Standard deviation 0 500 1000 1500 2000 2500 103 102 101 100 101 102 Number of cycles at inspectionUncertainty in m Standard error Standard deviation PAGE 73 73 parameter a little after 1,500 cycles which corresponds to the point where the S E and the STD start to agree in Figure 34 It can also be observed that by indentifying more parameters we increase the accuracy in identifying the damage growth parameter m not the we used the same set of measurements as for Figure 32 We are here able to identify mtrue instead of underestimating it. A B C Figure 35 Identified value of a0, b and m using least square fit. A) Identified value of a0. B) Identified value of b C) Identified value of m. 3.5 Conclusions This chapter presents an easy way of estimating the uncertainty of damage growth parameters using linearization of the nonlinear least square fitting process It has been shown that the method yields very good results when compared to MCS estimation for a 0 500 1000 1500 2000 2500 8 8.5 9 9.5 10 10.5 11 x 103 Number of cycles at inspectionIdentified value of a0 Identified value of a0 True value of a0 0 500 1000 1500 2000 2500 1 1.5 2 2.5 3 x 103 Number of cycles at inspectionIdentified value of b/2 Identified value of b/2 True value of b/2 0 500 1000 1500 2000 2500 3.78 3.8 3.82 3.84 3.86 3.88 3.9 3.92 3.94 3.96 3.98 Number of cycles at inspection Identified value of m Identified value of m True value of m PAGE 74 74 single variable case. In addition, the uncertainty of identified parameter is compared with that of Bayesian inference. For the multiple variable case, the interaction between variables can cause difficulties when not enough data are a vailable. Since Bayesian inference is computationally expensive for multiple variable identification, MCS were used to compare the quantif ied uncertainty. It is shown that the proposed method successfully identify the uncertainty of all parameters when enough data are present. PAGE 75 75 CHAPTER 4 LEAST SQUARESFILTERED BAYESIAN INFERENCE TO REDUCE UNCERTAINT Y IN DAMAGE GROWTH PRO PERTIES 4.1 Introduction The goal of this work is to improve prognosis (i.e., reducing uncertainty in damage growth parameters and uncertai nty in the remaining useful life) using noisy/biased structural health monitoring (SHM) data. In Chapter 2, we discussed how Bayesian inference can be used to identify damage growth parameters and to predict the remaining useful life (RUL) using SHM data. Although the method leads to a good estimate of the distribution of a conservative estimate of RUL, we did not identify the bias in the data because of computational cost in multidimensional identification. On the other hand the least square fit in Chapt er 3 is particularly good at indentifying deterministic parameters such as the bias. Another interesting characteristic of least square fit is that it smoothes out the noise when applied to a large enough sample of data. T he main objective of this chapter is to demonstrate the further reduction in uncertainty that may be achieved using Bayesian inference in combination with least square fit identification. In Chapter 2, a probabilistic approach using Bayesian inference was employed to progressively improve the accuracy of predicting damage parameters under noise and bias of sensor measurements It was discussed in more detail in Chapter 2 but we neglected the effect of bias in the likelihood calculation. It is discussed that including bias is possible but i t is computationally expensive. In Chapter 3, t he result obtained using Bayesian inference is then compared to that of the least square fit identification method. T he goal of this chapter is to combine the strength of both methods in order to further reduc e the uncertainty in damage growth parameters The PAGE 76 76 basic idea is to use the least square fit as a preprocessor to reduce variability in the measured data and to identify the deterministic bias, followed by Bayesian inference to reduce uncertainties in damage growth parameters. Therefore, we name the proposed method as the least squarefiltered Bayesian (LSFB) method. The final goal is to identify the damage growth parameter s distribution with enough accuracy to improve the estimat ion of the RUL of the structure. The approach is demonstrated for the model defined in Section 2 of Chapter 2, a throughthethickness crack in an aircraft fuselage panel that grows through cycles of pressurization. A simple damage growth model, Paris model, with two damage parameters is utilized. However, more advanced damage growth models can also be used, which usually comes with more parameters. Using this simple model we aim to demonstrate that SHM can be used to identify the damage parameters of a particular panel. This proc ess can be viewed as turning every aircraft into a flying fatigue laboratory. Reducing uncertainty in damage growth parameters can reduce in turn the uncertainty in predic ting RUL; i.e. in prognosis. The chapter is organized as follows. In S ection 4. 2 th e results are given for Bayesian inference, on the identification of the damage growth parameters and the estimate on remaining useful life. In Section 4. 3 the corresponding results are given for least square fit method. In Section 4. 4, comparable results are presented obtained using the LSFB method. 4.2 Characterization of D amage Growth P roperties using Bayesian I nference Depending on manufacturing and assembly processes, the actual damage parameters for individual aircraft can vary For a specific panel we assume that there exists a true value of the deterministic damage parameters. In the following numerical PAGE 77 77 simulation, the true damage will grow according to the true value of the damage parameters. On the other hand, the measured damage size will have bias and random noise. As a first step, we consider the distributions of m and C separately, which means that when we consider one variable as being uncertain we consider the other one as being known. The damage growth parameter m is a critical factor to determine the growth of damage. This parameter is normally measured using fatigue tests under controlled laboratory tests. However, uncertainty in this parameter is normally large at the material level because of variability in manufacturing and aging of the specific panel and also at a specimen level because of variability related to testing process It is possible to curve fit the data and obtain estimates of this parameter for the individual panel. However, curve fits do not take into account prior information on the distribution of these parameters, nor statistical information on the measurement uncertainty. We therefore use Bayesian statistics to identify these parameters. T he exponent m is the slope of the crack growth rate vs. cycle curve in the loglog scale [74] As a first step in developing a prognosis methodology, we assumed that an accurate value of C is known, while that of m is uncertain. Since the range of the exponent m is generally known from literature or material handbooks, we assume that the exponent is uniformly distributed between the lower and upper bounds. Then, the goal is to narrow the distribution of the exponent using the Bayesian statistics with measured damage sizes. The details of the Bayesian inference procedure are discu ssed in the Chapter 2, some similar results are presented here in order to be PAGE 78 78 compared to the other two methods discussed in this chapter. The distribution s shown below are obtained by performing MCS with multiple sets of measurements. We will discuss the application to the calculation of the RUL at every inspection interval In order to show the value of our method we compare the RUL calculated using the actual value of m, mtrue, and the distribution (mean one standard deviation) of the 5th percentile of the distribution of RUL obtained using the updated distribution of m at each inspection, for the case of negative bias, + 2mm and a noise of amplitude 1mm Another metric we use is the error between the maximum likelihood of the estimated distribution of RUL and the true RUL. Th e distribution of the 5th percentile of RUL is shown in Figure 41 It converges to the true RUL, but it is sensitive to errors in the data as it has been observed in Chapter 2 and as a result some estimates are unconservative. Figure 41 Distribution (mean one standard deviation) of the 5th percentile of RUL for b = +2mm and V = 1mm, using Bayesian inference Figure 42 shows the error between the maximum likelihood of the estimated distribution of RUL and the true RUL as defined in Eq. (2.20) note that positive values of the error corres pond to unconservative estimates. As observed previously Bayesian PAGE 79 79 inference method appears to be very sensitive to error in the data but in average it appears to lead to an accurate estimate of RUL. Figure 42 Distribution (on e sigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for b = +2mm and V = 1mm, using Bayesian inference One of the major advantages of SHM is that measurements can be performed frequently. Thus, the update in Eq. (2.7) can be performed as frequently as needed. However, since the damage grows slowly and the bias and noise of measurements are in general large, too frequent measurements may not help to narrow down the distribution of damage parameters using Bayesian inference, because it does not deal well with large samples of data. In addition, the discreti zation of m may result in an updated distribution equal to zero if the distribution is narrowed down too much. As a result of this characteristic we end up not using the amount of data to its full extent; even if the measurement can be done at every cycle, we use the measurement data at every 100 cycles. Another issue with Bayesian inference is its computational cost; it is very expensive to identify multiple parameters. As a result of this, bias is ignored in the calculation of the likelihood. Theoretically it is possible to identify the bias by consider ing it as an additional unknown parameter, but it would be computationally expensive as it PAGE 80 80 requires higher order integration. This is where least square fit comes into picture because of its ability to identify deterministic variables very accurately and r educe noise in the data. In this section, Bayesian inference is used to identify a single damage growth parameter, either m or C and assumed the other parameter is already known. In practice, when both parameters are unknown, Bayesian inference needs to update the joint PDF of both parameters. In general, this can be achieved by dividing the ranges of uncertain parameters into a grid and calculate the joint PDF value at each grid point. If the range is divided by 100 x 100 grid, the updating process includes 10,000 times calculation of likelihood, which requires uncertainty propagation for given parameter values. Thus, the updating process easily becomes computationally impractical. In addition, the bias and initial true crack size are uncertain because their values are unknown. They can also be considered as unknown parameters similar to damage growth parameters and to be identified through Bayesian inference. Then, a four dimensional joint PDF needs to be updated, which makes the Bayesian inference further computationally impractical. Before we present a new method of addressing this computational issue, we will present a conventional method of identifying unknown parameters in this section, along with its advantages and disadvantages. 4.3 Characterization of D amage P roperties U sing L east S quare F it Least square fit is the easiest and most commonly used way of identifying model parameters by minimizing the difference between measured data and predicted data from the physics model. In our application, the Paris model is used with the following unknown parameters: initial crack size, ai, damage growth parameter, m and bias, b PAGE 81 81 The parameter is still assumed to be known in order to compare with the results in the previous section. The least square fit problem can be formulated as 2 2 2 1/2 ,,min with 1 2im m meas m jjj i amb jm abaaNC a ( 4 1 ) The minimization problem is unconstrained except for lower and upper bounds on the identified parameters. In order to compare with the Bayesian inference, the same bounds on m are used with the Bayesian inference method. The bounds of bias are 3 and 3mm as w e assume that this is reasonable lower and upper limit s. Since the initial crack size is within the combined ranges of bias and noise from the measured initial crack size, its lower and upper bounds are selected accordingly. The least square fit is performed using the lsqnonlin function in MATLAB. At measurement cycle N all previous measurement data are used in the least square fit. Thus, more and more data are used in least square fit as the number of cycles increases. Although Bayesian inference is per formed at every N = 100 interval, the least squarefit is performed at every cycle because the fit is more accurate with more data. Similar to Bayesian inference, the identified parameters depend on synthetic measurement data. Thus, 1,000 sets of measurem ent data are produced by adding deterministic bias and random noise to the true crack sizes. Consequently, at every measurement cycle, we obtain 1,000 sets of identified parameters. In order to show how these parameters are distributed, Figure 43 shows intervals of mean plus and minus one standard deviation for two different combinations of noise and bias: (1) b = 2 mm and V = 3 mm and (2) b = 2 mm and V = 1 mm. These combinations are chosen because they represent extreme cases of noise for a given bias, and the cases with a PAGE 82 82 negative bias will give similar results. It can be observed that the parameters from the least square fit converge to the true values of the parameters as the number of cycles increases. Different from Bayesian inference, the effect of bias and noise is insignificant because the effect of noise is reduced by the amount of data used and the bias is identified using that method. The convergence is slow when the crack size is small, and it is fast when the crack si ze is large. This happens because a larger crack has faster crack growth, which leads to similar behavior observed in Bayesian inference. A B C Figure 43 Distributio n (mean one standard deviation obtained using 1,000 MCS simulations ) of fitt ed results of a0, m and b /2. A) Distribution of fitted a0. B) Distribution of fitted m C ) Distribution of fitted b /2 The main interest of prognosis is to predict the RUL at inspection. In order to show the value of our method, we compare RUL calculated using the actual value of m mtrue, PAGE 83 83 and the distribution of the 5th percentile of the estimated distribution of RUL using the identified m Figure 44 shows the mean one standard deviation of the 5th percentile of the distribution of remaining useful life for a bias of 2mm and a noise of maximum amplitude 1mm. It can be observed that least square fit leads to a very good identification of the parameters and extrapolation of RUL, but it does not give a very good estimation of a conservative estimate of RUL. This leads to slightly unconservative results even for the 5th percentile. Figure 44 5th percentile of RUL for b = 2mm a nd V = 1mm, using least square fit Figure 45 shows the error between the maximum likelihood of the estimated distribution of the RUL and the true RUL. Note that positive values of the error correspond to unconservative estimates. As observed previously, least square fit gives a good estimate of the value of the RUL. However, it appears unable to estimate the standard deviation of the RUL accurately, which results to mostly unconservative estimate of the RUL. PAGE 84 84 Figure 45 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for b = +2mm and V = 1mm, using least square fit, positive error corresponding to unconservative estimates 4.4 Least S quareF iltered Bayesian ( LSFB ) Method for Estimating RUL As described in the previous sections, both Bayesian inference and least square fit present advantages and limitations, but they appear to be complementary. The ability of least square fit to identify the bias and reduce the noise by including every measurement makes it a useful tool to preprocess the data. Then in order to identify the distribution of RUL, Bayesian inference can use the processed data. The idea is to measure the crack size i nformation at small intervals, and send these crack size data to the least square fit at every hundred cycles to reduce noise and to estimate initial crack size and bias. The fit to the measurement data has much less noise, and the effect of bias is reduced because it is identified in the model Figure 46 and Figure 47 illustrate the noise reduction and bias compensation of the least squares fit. The fit then provides crack sizes that can then be used in Bayesian inference in order to narrow down the distribution of m and obtain a more accurate prognosis. PAGE 85 85 Figure 46 illustrates one case showing the true crack size as the dashed line compared to the measurement every 100 cycles ( stars) obtained with a bias of 2mm and a noise of amplitude 1mm, and the damage sizes obtained using least square fit (dots). Note that least square fit uses measurements at every cycle. The dots represent the data that are used with Bayesian inference to update the distribution of m and estimate the distribution of RUL. Figure 46 Comparison of the true, measured and fitted damage sizes for b = 2mm and V = 1mm Illustration with one simulated set of measurements. Figure 46 illustrates the behavior of the crack over the entire life of the panel. Due to bias the general trend of the crack size is shifted down from the true crack size, and due to noise the crack growth is not consistent. In some cases, the measured crack size is reduced from the previous measurement, which is not possible physically. Since the bias is larger than the noise the measurements never come close to the true values. Crack sizes are always underestimated by at least 1mm. On the other hand, the estimated crack sizes using parameters from the least square fit are close to the true crack sizes and show much more consistent behavior. As can be observed in Figure 43 the identification of parameters is so close to the true value because with a noise of amplitude 1mm we can have a very good identification of the variables. If we had used PAGE 86 86 a larger noise for this illustration case we would have had a poorer approximation early on. Although SHM measurements are performed every cycle, data are shown at every 100th cycle in order consistently compare with that of Bayesian inference. Figure 47 on the other hand shows the behavior over two inspection intervals at the beginning (a) and toward the end of life (b). It shows how the least square fit improves as the crack grows, because we have more information. But in both cases, it can be observed that the least square fitted data are closer to the true values than the measured ones. This confirm s the fact that least square fit reduces the effect of bias and noise. A B Figure 47 Fitted (dots) and measured (stars) damage at early and late stage in damage growth compared to the actual damage size (dotted line) A)Fit at the early stage. B) Fit at the late stage Figure 48 shows the 5th percentile of remaining useful life obtained using LSFB for data with a bias of +2 mm and a noise of amplitude 1 mm shown previously for the two other methods. The figure shows the mean of the percentile plus/minus one standard deviation. The estimation converges but remains on the conservative side, and is narrower, hence less sensitive to the errors in measurements. 200 250 300 350 400 7.5 8 8.5 9 9.5 10 10.5 11 11.5 12 Number of cyclesDamage size 1900 1950 2000 2050 2100 22 23 24 25 26 27 28 29 30 31 Number of cyclesDamage size PAGE 87 87 Figure 48 Distr ibution (mean ones standard deviation) of the 5th percentile of RUL for b = +2mm and V = 1mm, using the LSFB method Figure 49 shows the error between the maximum likelihood of the estimated distribution of the RUL and the true RUL N ote that positive values of the error correspond to unconservative estimates. As observed previously LSFB gives a good and conservative estimate of the value of the RUL. Therefore, combining least square and Bayesian inference gives better results t han either of them separately. Figure 49 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for b = +2mm and V = 1mm, using LSFB, positive error corr esponding to unconservative estimates Figure 410 compares the three methods presented in this chapter. It can be observed that the LSFB method is a good compromise between least square fit and PAGE 88 88 Bayesian inference. It is much less sensitive to the noise in the data, the variability in the distribution is much smaller and the estimated 5th percentile is conservative. Figure 410 Comparison of the distribution (mean ones standard deviation) of the 5th percentile of RUL using the three methods Figure 411 compares the error between the maximum likelihood of the estimated distribution of the RUL and the true RUL for the three methods. As observed previously LSFB gives the most c onservative and most accurate estimate of the three methods. This can be explained by the fact that it combines the strength of both methods : the accuracy of least square fit with uncertainty quantification of Bayesian. T his results in a method that is bot h accurate and conservative. Both Bayesian and least square fit were fairly accurate but neither was conservative; the first one as a result of its sensitivity to error in the data, the second because of its poor ability to estimate the uncertainty in RUL. PAGE 89 89 Figure 411 Comparison of the distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL using the three methods, positive error corresponding to unconservative estimates Figure 412 shows the updated/identified distribution of m at every inspection for the three methods discussed in this chapter using one set of measurements. It can be observed that LSFB leads to similar results as Bayesian inference but with a slightly better convergence that results from the fact that we are not ignoring the effect of the bias. Note that that plot does not show the effect of the noise though since it is only one set of measurements, this can only be observed in t he RUL estimation show in the previous two figures. Figure 412 Comparison of the distribution (onesigma intervals) of m for one set of measurements using the three methods. PAGE 90 90 4.5 Conclusions We presented prognosis results for three methods: Bayesian inference, least square fit, and a combination of both, least squarefiltered Bayesian. The results show that the first two methods are good at identifying the damage growth parameters or estimating the RUL. They also show limitations that can lead to unconservative results Bayesian as a result of its sensitivity to the error in the measurements and least square fit as a result of its inability to estimate the uncertainty in RUL. However, the two methods can be combined to take advantage of the strengths of both methods. The combination, LSFB, led to estimation of RUL that converges to the actual RUL while remaining on the conservative side. Another advantage of LSFB is that it is less sensitive to the errors in measurements. Note t hat even though the results and conclusions presented in this chapter are for a specific case of bias/noise combination similar results and conclusions were obtained for other cases. PAGE 91 91 CHAPTER 5 IDENTIFICATION OF EQ UIVALENT DAMAGE GROW TH PARAMETERS WHEN USING WRONG RANGE OF STRESS INTENSITY FACTOR 5.1 Introduction In this chapter equivalent model parameters which are different from the true values are identified which result the same prediction in the model based prognosis model. Once these model parameters are identified, they can be used to predict the future behavior of the system. M any physical models are limited to simple conditions. For example, the Paris model [39] describes the rate of crack growth in terms of material properties and the stress intensi ty factor. The simplest available expression for the stress intensity factor is the infinite plate with a throughthethickness center crack under mode I loading. In reality the stress intensity factor is a complicated function of applied loading, boundary conditions, crack position, geometry, and material properties. Although there are many correction factors for taking into account for finite plate size or edge cracks [75] still they are l imited in representing complex engineering systems. In such case where no algebraic relationship is know n analytical methods such as finite element methods can be used to find the stress intensity factor for a given crack geometry. The objective of this c hapter is to demonstrate that in model based identification, one can use simple models to predict the remaining useful life (RUL) even if they do not model accurately the actual behavior of the monitored damage. This is accomplished through the identificat ion of an equivalent damage growth parameter that compensates for the difference between the model and the true stress intensity factor. In this chapter, a square plate is chosen as a geometry for explanation. The addition of cracks and holes to the plate causes the crack tip state of stress to PAGE 92 92 experience finite plate effects in both the horizontal and vertical directions as well as stress concentrations caused by the addition of holes to the plate. As no handbook solution is known which considers all these effects, the damage growth is simulated using the extended finite element method (XFEM) for calculating numerical stress intensity factors and Paris law is used to grow the crack. This numerical stress intensity factor is considered true compared to the analytical stress intensity factor based on infinite plate with center crack. XFEM [76] allows for discontinuities to be modeled independently of the finite element mesh, which avoids costly remes hing as the crack grows. The stress intensity factors which are the driving force for crack growth are calculated using the domain form of the contour integrals [77] In practice, the actual damage sizes would be measured using structural health monitoring (SHM) systems to detect damage location and size. In this chapter, like in the previous ones, instead of using actual measurement data, synthetic data are generated to demonstrate the insensitivity of RUL to errors in the stress intensity model. First, the true values of Paris model parameters are assumed. Then, the true crack will grow according to the given parameters and prescribed operating and loading conditions. Thus, the true crack size at every measurement time is known. With the true cra ck size, the RUL is defined when the crack size reaches the critical crack size, which is a function of material, operating, and loading conditions. It is assumed that the measurement instruments may have a deterministic bias and random noise. These bias a nd noise are added to the true crack sizes, to generate synthetic measured crack sizes. Then, these data are used to identify the damage growth parameters and thus the RUL. In this way, it is possible to evaluate the accuracy of prognosis method. PAGE 93 93 Of the ma ny methods presented before for parameter identification, the least squarefiltered Bayesian method introduced in Chapter 4 is used to identify damage growth parameters using the synthetic data. This method applies nonlinear least square method to the meas urement data, so that the magnitude of noise can be reduced, followed by Bayesian inference, [35] to identify a probability distribution for model parameters. The identified distribution of damage growth parameters can then be used to estimate the distribution of RUL. An important question that is explored in this chapter is whether or not a simple stress intensity factor model can be used for general crack geometries for the purpose of prognosis. The key concept is that the Paris model can be consi dered as an extrapolation tool. Thus, even if the simplified stress intensity factor expression is grossly inaccurate, LSFB will identify equivalent damage growth parameters, different from the true ones, such that the model accurately predicts future damage growth behavior. The chapter is organized into the following sections. In Section 5. 2, the crack growth model is introduced. In Section 5. 3, the least squarefiltered Bayesian method is summarized. Results are presented in Section 5. 4, three problems wi th increasingly complicated geometry, in the sense that the center crack in an infinite plate model is an increasingly erroneous predictor of the actual state of stress at the crack tip. Concluding remarks and future work are presented in Section 5. 5. 5.2 Crack G rowth M odel s 5.2.1 Damage Growth Model Used for LSFB and RUL Estimation The Paris model [39] gives the fatigue crack growth rate as a function of material properties C and m and the stress intensity factor range K as PAGE 94 94 mda CK dN ( 5 1 ) This model is created from experimental observation. For a center crack in an infinite plate in Mode I loading, the stress intensity factor range K is given as Ka ( 5 2 ) where is the applied stress range and a is the c haracteristic crack size. The characteristic crack length ai at the ith cycle derived from Eqs. ( 5.1 ) and ( 5.2) is given as 2 2 1 2 01 2m mm iim a aNC a ( 5 3 ) where a0 is the initial crack size and Ni the number of cycles at the ith measurement. Similarly, the number of cycles to failure for a center crack in an infinite plate can be derived by integrating Eq. ( 5.1) as 11 222 1mm Ci faa N C m ( 5 4 ) where aC is the critical crack size. Note that Nf is uncertain because the initial crack size and damage growth parameters are uncertain. Although the critical crack size can be uncertain, it can be specified by airliner as a criterion to fix the damage, in this case we use this as a deterministic value. 5.2.2 Damage Growth Model Used for Measurement Data Generation In general, the accuracy of Eq. ( 5.2 ) depends on geometrical effects, boundary conditions, crack shape, and crac k location. A more general expression [75] is Kfa ( 5 5 ) PAGE 95 95 where f ( is the correction factor, given as the ratio of the true stress intensity factor to the value predicted by Eq. ( 5.2) The value of is given in terms of the geometry and characteristic crack size and is problem dependent. An example of the effect that the correction factor f ( can have on the stress intensity factor curve for a range of crack sizes is shown in Figure 1 for a center crack in an infinite plate, center crack in a finite plate, and an edge crack in a finite plate [75] For this case the assumed plate width for the finite models was 0.2 m. Figure 51 Comparison of stress intensity response for some correction factors and crack sizes. Plate width is 200 mm. For a complex geometry, analytical expressions given in Eqs. ( 5.4 ) and ( 5.5) do not exist. In such a case, numerical methods can be used to calculate the stress intensity factor K and Eq. ( 5.1) can be numerically integrated to calculate the crack size as a function of the number of cycles using XFEM as follows Modeling crack growth in a traditional finite element framework is a challenging engineering task. The fini te element framework is not well suited for modeling crack growth because the domain of interest is defined by the mesh. At each increment of crack growth, at least the domain surrounding the crack tip must be remeshed such that PAGE 96 96 the updated crack geometry is accurately represented. If a large number of cycles are to be considered, this repeated remeshing can consume a large amount of the computational time for the analysis which was avoided here through the use of XFEM XFEM allows discontinuities to be represented independently of the finite element mesh [76] Arbitrarily oriented discontinuities can be modeled by enriching all elements cut by a discontinuity using enrichment functions satisfying t he discontinuous behavior and additional nodal degrees of freedom. For the case of a domain containing a crack and voids [78] the approximation is: h III I IuxVNuHab ( 5 6 ) where NI are the finite element shape function, V is the void enrichment function, H is the Heaviside enrichment function, are the crack tip enrichment functions, and uI, aI, and bI are the classical and enriched degrees of freedom (DOF). To decrease the computational time for the repeat ed solutions, a reanalysis algorithm (Pais, 2010) is used which takes advantage of the large constant portion of the global stiffness matrix. The mixedmode stress intensity factors KI and KII for the given cracked geometry were calculated using the domain form of the interaction integrals [77] The direction of crack growth was calculated using the maximum circumferential stress criterion [77] The effective stress intensity factor [79] given as 44 48effIIIK KK ( 5 7 ) was used to convert the mixedmode stress intensity factors into a single value for used in Paris law. The crack growth at a given cycle is given as PAGE 97 97 m effaCK ( 5 8 ) The implementation of XFEM used here was verified using the center crack in a finite pl ate problem For this problem the theoretical finite correction factor based on the equations of elasticity for a center crack in a finite plate [75] is given as 243 240 s 0 c 5 e1 f ( 5 9 ) where = a / w and a and w are the half crack length and half plate width. This model assumes that the plate is finite in the x direction and infinite in the y direction. A comparison of the crack lengths as a function of the number of cycles was first performed to ensure the accuracy of the XFEM dat a provided to the identification routine. A comparison of the results is shown in Figure 52 Figure 52 Comparison of t heoretical and XFEM crack growth curves using different plates heights, h for X FEM in order to validate the XFEM code and assess the validity of the theoretical model for the chosen plate geometry It was noticed that for the chosen plate geometry square plate of width and height 2h of 0.2m the handbook solution and XFEM models predicted different crack growth curves. Increasing the height of the plate leads to good agreement with the theoretical 0 500 1000 1500 2000 2500 10 15 20 25 30 35 40 Number of Cycles [N]Crack Length [mm] Handbook Solution XFEM, h = 0.125 m XFEM, h = 0.100 m PAGE 98 98 values indicating that the chosen crack configuration has a finite effect from both the vertical and horizontal directions The resulting difference in f ( ) caused by the vertical finite effect as a function of the number of cycles is shown in Figure 53 It is clear that as the height (h) of the plate increases, the correction factor obtained from XFEM agrees well with the theoretical value. Figure 53 Theoretical and XFEM prediction of f ( ) 5.3 Least S quare F iltered Bayesian (LSFB) M ethod Bayesian inference and least square fit are often used for identifying unknown model parameters and present advantages and limitations, but as discussed in Chapter 5, they appear to be complementary. Least square fits ability to identify the bias and reduce the noise makes it a useful tool to process the data in order to identify the distribution of R UL using Bayesian inference Note that we chose to update Paris exponent m here but similar results could be obtained by updating C or both parameters together. The LSFB method processes information collected at every cycle by least square fit in order to reduce the noise, and identify the bias, b as well as the initial crack size a0. The least square problem is expressed as 0 500 1000 1500 2000 2500 10 15 20 25 30 35 40 Number of Cycles [N]Crack Length [mm] Theoretical XFEM, h = 0.125 m XFEM, h = 0.100 m PAGE 99 99 0,, 2minmeas amb i iiaba ( 5 10) where meas ia are the synthetic measured crack sizes with noise model to simulate measurement data, and ai is the damage size at ith cycle. The LSFB method in this chapter uses th e K given by Eq. ( 5.2) and an effective value of m is identified resulting in the same solution to Eq. ( 5.1) The identified values of a0, m and b are then used to generate a new estimate of the damage size at the ith cycle using Eq. ( 5.3) they are referred to as filtered data. The se data are then used in Bayesian inference in order to narrow down the distribution of m and obtain a more accurate prognosis. The identified a0 and b are considered as deterministic and only uncertainty in m is considered in Bayesian inference. Bayesian inference is based on the Bayes theorem on conditional probability. It is used to obtain the updated (also called posterior) probability of a random variable by using new information. In this chapter, since the probability distribution of m given a is of interest, the following form of Bayes theorem is used [70]  ini updt inilamfm fm lamfmdm ( 5 11) where fini the assumed (or prior) probability density function (PDF) of m fupdt the updated (or posterior) PDF of m and l ( a m ) is the likelihood function, which is t he probability of obtaining the characteristic crack length a for a given value of m T he derivation of the likelihood function can be found in Chapter 2. The likelihood function is designed to integrate the information obtained from SHM measurement to the knowledge about the distribution of m Instead of assuming an analytical form of the likelihood function, uncertainty in measured crack sizes is PAGE 100 100 propagated and estimated using the Monte Carlo simulation (MCS). Although this process is computationally expensive, it will provide accurate information for the posterior distribution. Once the distribution of m has been identified at cycle Ni, it can be used to predict the RUL. The distribution of RUL is calculated at every SHM measurement cycle Ni using MCS a nd the RUL is estimated using Eq. ( 5.4) derived from Paris law. This allows us to estimate the distribution and from there obtain the 5th percentile. The 5th percentile of samples of Nf is used as a conservative estimate of RUL in order to have a safe prediction. Since random noise is added to the synthetic data, the result may vary with different sets of data. Thus, the above process is repeated with 100 sets of measurement data and mean plus and minus one standard deviation intervals are plotted. In order to show the value of the LSFB method, the RUL calculated using the distribution of mLSFB and the distribution (mean one standard deviation) of the 5th percentile of the distribution of RUL obtained using the updated distribution of m at each inspection are compared. 5.4 Results For each example an aluminum 7075 square plate with edge lengths of 0.2 m and thickness 2.48 mm and an initial crack size of 0.01 m is used. Aluminum 7075 has Youngs modulus E of 71.7 GPa, Poissons ratio of 0.33, critical mode I stress intensity factor KIC of 30 MPa m Paris Law constant C of 1.5E 10, and an assumed, deterministic Paris Law exponent m of 3.8. The plate is assumed to be an aircraft panel with radius 3.25 m, which undergoes pressurization cycles of amplitude 0.06 MPa. The PAGE 101 101 relatively large initial crack size is chosen because many SHM sensors cannot detect small cracks. In addition, there is no significant crack growth when the size is small. However, this size is still too small to threaten the safety of an airplane. True crack growth data was calculated using the extended finite element method using stress calculated from the pressurization model. XFEM simulations were performed on a structured mesh of square linear quadrilateral elements with characteristic length of 1 mm. Each cycle of fatigue crack propagation was modeled until the equivalent mode I stress intensity factor exceeded KIC. The characteristic crack length at each iteration was then used in the identification of an equivalent Paris Law exponent through the use of the least squarefiltered Bayesian method with the simplified stress intensity formula, Eq. ( 5.2) 5.4.1 Center Crack in a Finite Plate The first example considered is that of a center crack in a finite plate as shown in Figure 54 Only the right half of the plate was modeled with XFEM through the use of symmetry. Figure 54 A center crack in a finite plate. PAGE 102 102 The corresponding curve of the correction factor f ( ) which this edge crack represented is given in Figure 55 For this case, it was found that failure occurred at 2 070 cycles with a corresponding crack length of 37.5 mm. If we consider an arbitrary critical crack size of 24 mm the total RUL is then reduced to 1,686 cycles. Figure 55 Correction factor for center crack As the LSFB analysis results in a final distribution of m the predicted crack lengths for this distribution are plotted and compared directly to the XFEM data in Figure 56 The XFEM data fall within the bounds of the LSFB identification. Figure 56 Comparison of XFEM crack growth data with crack growth predicted from LSFB analysis PAGE 103 103 Figure 57 shows the updated distribution of m using LSFB. As expected it compensates for the error in K by over estimating m This is expected because the damage grows faster than it would if it was a center crack in an infini te plate. Figure 57 Updated distribution of mLSFB using one set of data for a center crack in a finite plate. Figure 58 shows in grey the distribution (mean one standard deviation obtained from 10 0 sets of different measurements) of 5th percentile of RUL discussed in Section 2 for that geometry, compared to the actual remaining useful life for an arbitrarily chosen deterministic critical damage size aC of 24 mm. It can be observed that the estimate of RUL converges to the true RUL from the conservative side. As expected the slight over estimation in the correction factor leads to a slight over estimation in the RUL. 3.7 3.75 3.8 3.85 3.9 3.95 4 0 5 10 15 20 25 30 35 40 45 50 mPDF of m PAGE 104 104 Figure 58 Distribution (mean one standard deviation ) of 5th percentile of RUL for a center crack in a finite plate compared to the true RUL ( black line) and the RUL obtained using the true m and K for a centre crack in an infinite plate (dark gray line) Figure 59 shows the error between the maximum likelihood estimation of the estimated RUL distribution and the true RUL as described in Eq. ( 2.20) it can be observed that the estimate is very conservative at the beginning and it becomes then unconservative but with a smaller amplitude. Figure 59 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for a center crack in a finite plate PAGE 105 105 5.4.2 Edge Crack in a Finite Plate Next, an edge crack in a finite plate was considered as shown in Figure 510. For this case the boundary conditions were fixing the lower left hand corner and allowing the top left corner to only move in the vertical direction. Figure 510 Edge crack in a finite plate The correction factor corresponding to the finite edge effect which this edge crack represented is given in Figure 511. For this case, it was found that failure occurred at 1 018 cycles with a corresponding crack length of 27.2 mm. If we consider an arbitrary critical crack size of 24 mm the total RUL is then reduced to 955 cycles. Figure 511 Correction factor for edge crack. PAGE 106 106 As the LSFB analysis results in a final distribution of m the predicted crack lengths for this distribution are plotted and compared directly to the XFEM data in Figure 512. The XFEM data fall within the bounds of the LSFB identification. Figure 512 Comparison of XFEM crack growth data with crack growth predicted from LSFB analysis. Figure 513 shows the updated distribution of m using LSFB. As expected it compensates for the error in K by over estimating m It can be observed that the larger the error in K the more m is overestimated to compensate for it. Compared to the case of a center crack in an infinite plate the range of the identified distribution of m is wider, which is likely caused by the increased difference between the actual and assum ed models for the stress intensity factor. PAGE 107 107 Figure 513 Updated distribution of mLSFB using one set of data for an edge crack in a finite plate. Figure 514 shows the distribution of 5th percentile of RUL discussed in Section 2 for that geometry, compared to the actual remaining useful life. As for the previous geometry it can be observed that the estimate of RUL converges to the actual value from the conservative side. In this case the error in the cor rection factor is up to 30% but it leads to an overestimation in the RUL of almost 100%. Figure 514 Distribution (mean one standard deviation) of 5th percentile of RUL for an edge crack in a finite plate compared to the true RUL (black line) and the RUL obtained using the true m and K for a centre crack in an infinite plate (dark gray line) Figure 515 shows the error between the maximum likelihood of the estimated distribution of the RUL and the true RUL. As observed previously LSFB leads to a 3.8 3.9 4 4.1 4.2 0 5 10 15 20 25 30 mPDF of m PAGE 108 108 somewhat unc onservative estimate of the RUL if you consider the maximum likelihood but it converges to the true value fairly accurately. Figure 515 Distribution (onesigma intervals) of error between the true RUL and the maximum likeliho od of the estimated RUL distribution for an edge crack in a finite plate 5.4.3 Center Crack in a Plate with Holes The final example considers differences between the actual and predicted model that may be caused by localized stress concentrations in struct ures. Four holes are inserted into the plate as shown in Figure 516. Only the right half of the plate was modeled with XFEM through the use of symmetry. Figure 516 Center crack in a finite plate wit h holes. PAGE 109 109 Unlike the other examples presented, the authors are unaware of an approximation to f ( ). Th erefore, the correction factor is obtained from XFEM as shown in Figure 517. For this case, it was found that failure occurred at 625 cycles with a corresponding crack length of 24.2 mm. If we consider an arbitrary critical crack size of 24 mm the total RUL is then reduced to 605 cycles. Figure 517 Correction factor for plate with holes. As the LSFB analysis results in a final distribution of m the predicted crack lengths for this distribution are plotted and compared directly to the XFEM data in Figure 518. The XFEM data fall within the bounds of the LSFB identification. The identified crack size distribution is wider than in the previous cases, it corresponds to the model being increasingly far away from reality. PAGE 110 110 Figure 518 Comparison of XFEM crack growth data with crack growth predicted from LSFB analysis. Figure 513 shows the updated distribution of m using LSFB. As expected it compensates for the error in K by over estimating m The same conclusion can be drawn as previously, the larger the error in K the more m is overestimated to compensate for it. Figure 519 Updated distribution of mLSFB using one set of data for an center crack in a finite plate with holes. Figure 520 shows the distribution of 5th percentile of RUL discussed in Section 2 for that geometry, compared to the actual remaining useful life for a critical damage size of 24 mm. As for the previous geometries it can be observed that the estimate of RUL 4 4.05 4.1 4.15 4.2 4.25 4.3 0 2 4 6 8 10 12 14 16 18 mPDF of m PAGE 111 111 converges to the actual value from the conservative side. It has to be observed that the estimation is not as accurate but this can be explained by the fact that the geometry is very differ ent from the one assumed in the model and the number of cycles to failure is much smaller. The same conclusion that was drawn for the edge crack can be drawn in this case, the up to 40% error in the correction factor leads to a 200% error in the RUL estima tion Figure 520 Distribution (mean one standard deviation) of 5th percentile of RUL for an edge crack in a finite plate compared to the true RUL (black line) and the RUL obtained using the true m and K for a centre crack in an infinite plate (dark gray line) The error between the maximum likelihood of the estimated distribution of the RUL and the true RUL presented in Figure 521 shows that despite the fact that we are using simplistic model in which the range of stress intensity factor does not account for the complexity of the geometry we are able to estimate the RUL not only with accuracy but also fairly conservatively. PAGE 112 112 Figure 521 Distributio n (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for an center crack in a finite plate with holes 5.5 Conclusions Effective damage growth parameters were identified using the LSFB method for cases of finite plates with complex geometric effects. For prognosis purpose, t he stress intensity factor relationship was assumed to follow the center crack in an infinite plate and the Paris Law exponent m was identified which compensates for the incorr ect stress intensity factor relationship. Damage growth was simulated at each loading cycle though the use of the extended finite element method. This represents the versatility of the proposed method in that it does not require a priori know ledge of the c orrection factor f ( ) The m aximum likelihood value of the updated distribution of m and the RUL curves show good agreement with the simulated results. It is especially encouraging that the RUL converges from the conservative side. Although t he method is demonstrated here updating only one parameter, m of Paris law, the same idea can be applied to the parameters m and C together. This should allow for even more accurate results because it would allow for more flexibility in PAGE 113 113 fitting the equivalent model. The feasibility of using XFEM in the calculation of the likelihood function will also be explored which may identify the true m and C PAGE 114 114 CHAPTER 6 C ONCLUSION S In this work we first developed and employed the Bayesian inference method to identify panel spec ific damage growth parameters using damage sizes measured from SHM sensors. The actual measurement environment is modeled by introducing deterministic bias and random noise. The likelihood function is calculated by comparing measured crack growth with simulated crack growth, which requires uncertainty propagation though the physics model that governs the crack growth. Due to many uncertainties involved, the RUL is predicted statistically and used 95% conservative estimation. Through numerical examples, it i s shown that t he probability distributions of the two damage growth parameters, m and C are effectively narrowed and converged to the true values. Note that since m and C are strongly correlated if wrong information i s used in the model for the determinis tic assumed to be known variable, our method would compensate for that error I n fact it identifies the crack growth behavior rather than material parameter s. The large number of SHM data diminishe s the effect of noise, and thus, the identified damage growth parameters are relatively insensitive to them However, it does not completely remove the effect of bias and noise remains and they affect the identification of true damage growth parameters. The identified distributions of parameters are then used to estimate the RUL with 95% confidence. In all combined cases with noise and bias, the proposed method converges to the true RUL from the conservative side. In a more general approach, it is possible to update both m and C using their joint PDF (see Appendi x B) In addition, the unknown bias can also be considered as an PAGE 115 115 uncertain variable and can be updated together. However, as the number of variables increases, the computational cost increases significantly because the proposed method is based on MCS for uncertainty propagation. In the future, we will investigate the possibility of reducing computational cost by utilizing surrogate modeling techniques. The abovementioned computational cost related limitations are what led us to then present and compare similar prognosis results for two other methods: least square fit, and LSFB method. Those results show that even though the Bayesian and least square fit are very good at identifying the damage parameter or estimating the remaining useful life, they both have limitations, they cannot identify the damage with the extreme accuracy and estimate the distribution of RUL at the same time. But they are complementary, such that they can be combined and by using the advantages of both methods we come up to a third oneLSFB The LSFB giv es an estimation of RUL that converges to the actual RUL faster than Bayesian inference while remaining on the conservative side. Another advantage of LSFB is that it is less sensitive to the errors in measurements. After developing the L SFB method, we applied it to more complex damage geometry, using extended finite element method to generate the damage grow th that are used to generate the measure damage sizes. Effective damage growth parameters are identified using stress intensity facto r for a center crack in an infinite plate which is obviously making wrong assumptions. It is observed that the equivalent parameters compensate for the modeling error this results in a fairly accurate estimate of the RUL. This shows that using a damage gr owth model as an extrapolation device allows us to PAGE 116 116 use a simpl e model despite the fact that it is making wrong assumption and still have an accurate estimate of RUL. PAGE 117 117 APPENDIX A B AYESIAN INFERENCE USING DAMAGE GROWTH INFORMATION I n Chapter 2, the measured crack growth data from diagnosis is used to characterize the damage growth parameters. However, since these parameters are related to damage growth rate, not damage size, it would be more appropriate if damage growth information is used to identify these p arameters. This chapter summarizes likelihood calculation procedures using damage growth information. Numerical results shows that the results are similar to that of Chapter 2. N and N is the current cycle. The half crack growth between two measurements can be defined as measmeasmeas NNNNaaa (A. 1 ) As explained in Chapter 2 Bayesian inference is based on on Bayes theorem on conditional probability It is used to obtain the updated (also called posterior) probability of a random variable by using new information available for the variable. In this work since the probability distribution of m a is of interest, we use the following form of Bayes theorem [70] :  ini updt inilfm fm lfmd am m am (A. 2 ) where fini is the assumed (or prior) probability density function (PDF) of m fupdt is the updated (or posterior) PDF of m and l a m ) the likelihood function, which is here the a for a given value of m The step that differentiates this version of Bayesian inference from the on presented in Chapter 2 is the likelihood definiti on and calculat ion. PAGE 118 118 The likelihood function is designed to integrate the information obtained from structural health monitoring ( SHM ) measurement to the knowledge about the distribution of m The physical interpretation of the likelihood is the PDF value of the true crack growth at measured crack growth for given m. Although the true crack growth would be a single value, it is considered to be randomly distributed in the viewpoint of measured crack growth due to various uncertainties in the process. Thus, it is important to estimate the distribution of true crack growth. In general, the measured crack growth includes the effect of bias and noise of the sensor measurement as well as uncertainty in input loads. Let aN be the true half crack size, b the bias, and vN the noise at the current cycle N The measured half crack size, meas Na is then given as 2 2meas NNNa abv (A. 3 ) For subsequent simulated measurements, the bias b remains constant, while we assume that the noise vN is uniformly distributed within the range of [ V + V ] The above expression can be used to define the crack growth between consecutive two SHM measurements as follows measmeasmeas NNNNNNaaaav (A. 4 ) aN is the true crack growth and vN = vN vN N. Although vN and vN N have the same range of [ V, +V], they are independent. At a given SHM measurement, the measured half crack size in Eq. (A. 3 ) has the same distribution type with the noise, while the measured crack growth in Eq. (A. 4 ) has the same distribution type with vN. When the distribution of noise is Gaussian, both of them will also be Gaussian. Since there is no information regarding the distribution of PAGE 119 119 noise, in this work it is assumed to be uniformly distributed with mean at zero. Thus, the measured crack size is also uniformly distributed. On the other hand, it can be easily shown using Monte Carlo simulation (MCS) that vN is triangularly distributed with mean aN. Thus, the respective distributions can be defined as: ~ Uniform/2/2 ; /2/2 ~ Triangular ; ; meas NNN meas N NNNbVbV aVa a V aaa a (A. 5 ) The quantities defined above only involve measurement error. In general, however, the crack growth model may have modeling error, which is related to numerical simulation. In order to calculate the likelihood function, we introduce a simulated half crack size, sim Na that involves a modeling error, sim Ne, as sim sim NNNamaem (A. 6 ) In the above equation, the superscript sim is used for modeling error because it also includes propagated uncertainty through numerical simulation. The simulated crack size depends on Pari s parameters, m and C as well as the initial crack size. Since we only consider uncertainty in m Eq. (A. 6 ) only depends on it. Similarly, the simulated crack growth can be written as sim sim NNNamaem (A. 7 ) Different from measurement errors, the uncertainty in sim Ne is not well ch aracterized; it often requires MCS through the physics model that governs the crack growth. PAGE 120 120 The idea of calculating likelihood is to identify the damage growth parameter m by comparing the measured crack growth, meas Na with the simulat ed crack growth, sim Na with given m The difference between these two growths can be defined as sim meas NNdmm aa (A. 8 ) If the PDFs of meas Na and sim Na are available, then the PDF of d can also be calculated. The likelihood l a m ) is then defined as the value of this PDF at d ( m )=0. Since this rarely happens, we will use MCS to calculate the likelihood. Since MCS is a discrete process, it is not trivial to calculate the PDF directly. Instead, we will use the probability of  d  with being a small constant as a definition of likelihood:   lPd am (A. 9 ) Note that if the right hand side is divide d by 2 and if approaches zero, then the likelihood becomes the value of PDF at d ( m )=0. In the viewpoint of Eq. (A. 2 ) since the posterior distribution will be normalized, the above definition works for likelihood although it is given in the form of probability. If we calculate l a m ) purely by sampling meas Na and sim Na then the tolerance needs to be large enough to include enough samples to reduce sampling errors. On the other hand if is too large, we will incur errors due to nonlinearity in the likelihood function. In general, since the measurement error that controls meas Na is independent of the modeling error that controls sim Na separable sampling can be performed, and samples of d in Eq. (A. 8 ) can be calculated by comparing all possible combinations of the two sets of samples [71] In addition computational efficiency can significantly be improved PAGE 121 121 since the analytical PDF of meas Na is available from Eq. (11). The PDF of sim Na is not available analytically, because it is obtained by propagating uncertainties through the crack growth model. The definition of likelihood in Eq. (A.9) can be expanded by 00   lPdPdPd am (A. 10) Using conditional expectation on the second term on the right hand side we obtain 00 0 0sim N sim Nsimmeas NN simmeas simsim NN simNN a sim simsim measN simNN aPdP Pf aa aaaa aaa d Ffd (A. 11) where fsim( x) is the PDF of sim Na and Fmeas( x) is the cumulative density function ( CDF ) of meas Na The last relation is obtained from the definition of CDF; i.e., by considering meas Na as the only random variable, si meas N N measN m simPF aaa Similarly, the first term can be written as 00sim N sim Nsimmeas simsim NN simNN a sim simsim measN simNN aaaaa a Pd P fd Ffd aa (A. 12) Thus, by combining Eqs. (A.11) and (A.12) the likelihood can be written as 2 sim N sim Nsim sim simsim measN measN simNN a sim simsim measNsimNN aam a a aa lFFfd ff aa d a (A. 13) where the central finite difference approximation is used in the second relation, which becomes exact when normalized, the coefficient 2 can be ignored. The above expression is in particular PAGE 122 1 22 convenient for separable MCS because the analytical expression of fmeas( x) is known, and fsim( x) can be evaluated by propagating uncertainty through numerical simulation. Let M be the number of samples in MCS, the likelihood can then be calculated by 1 1sim Nsim simsim measNsimNN a M sim measNi ilf a fd maaa fa M (A. 14) First, input random samples, such as noise and pressure, are generated according to their distribution types. These input random samples are propagated through the Paris model to produce samples of crack growth sim Na Second, the values of PDF sim measNf a are evaluated for all samples, whose average is the likelihood. The numerical experiments showed that M = 2,000 is enough to obtain a smooth distribution of the likelihood function. Note th at likelihood calculation is computationally intensive because Eq. (A.14) needs to be evaluated for every m in the range of Eq. (A.2) In addition, the Bayesian inference in Eq. (A.2) is repeated at every inspection cy cle. Similar results as in Chapter 2 are shown here, Figure A 1 shows the updated PDFs of m at every 1,000 cycles for one set of measurements, note that another set of measurements might lead to slightly different updated distribution but not significantly different. It is clear that as the crack grows, the PDF of m becomes narrower and it converges to the true value of mtrue = 3.8. It is noted that the convergence becomes faster as the crack size increases because the crack growth is faster for a larger crack. PAGE 123 123 Figure A 1 Updated probability density functions of m (mtrue 3.8, Ctrue 1.5E 10, b = 0 mm, V = 1 mm), illustration with one simulated set of measurements. Figure A 2 shows the effect of bias on the final updated PDF of m The noise in crack detection is assumed to be Uniform(1, +1) and two different biases are used: b = 2 mm and b = +2 mm It is clear that bias shifts the maximum likelihood point (the peak of PDF) from that of the true value; the negative bias overestimates the PDF of m while the positive bias underestimates it Figure A 2 Effect of bias on updated PDF of m (mtrue 3.8, Ctrue 1.5E 10, V = 1 mm), i llustrat ion with one simulated set of measurements. 3.4 3.6 3.8 4 4.2 0 5 10 15 20 25 30 35 mPDF of m a0 = 10 mm a1,000 = 14.9 mm a2,000 = 27.8 mm a2,400 = 41.4 mm 3.4 3.6 3.8 4 4.2 0 5 10 15 20 25 30 35 40 mPDF of m b = +2mm; V = 1mm Initial distribution of m b = 2mm; V = 1mm mtrue = 3.8 PAGE 124 124 Figure A 3 shows the effect of noise on the PDF of m when b = 0 mm It is obvious that noise has an effect on the standard deviation but it does not shift the distribution as the bias does. The smaller the noise is, the narrower the final PDF of m Figure A 3 Effect of noise on the updated PDF of m (mtrue 3.8, Ctrue 1.5E 10, b = 0 mm), i llustration with one simulated set of measurements Table A 1 shows statistical characteristics, such as the maximum likelihood, mean and standard deviation of PDF of m corresponding to Figure A 2 and Figure A 3 Note that as those figures each of those sets of characteristics are obtained using a single set of measurements It can be observed that the mean and maximum likelihood values are minimally affected by the bias and noise. However, the standard deviation increases with a large noise. As expec ted, a positive bias (true crack size is smaller than measured one) leads to underestimation of m Table A 1 Statistical characteristics of final PDF of m with different combinations bias/noise, illustration with one simulated set of measurements. Effect of noise Effect of bias Bias, noise (mm) b = 0, V = 1 b = 0, V = 3 b = 2, V = 1 b = +2, V = 1 Max. likelihood 3.80 3.80 3.82 3.78 Mean 3.80 3.80 3.82 3.78 Standard deviation 0.01 0.04 0.01 0.01 3.4 3.6 3.8 4 4.2 0 5 10 15 20 25 30 mPDF of m b = 0mm; V = 1mm mtrue = 3.8 b = 0mm; V = 3mm Initial distribution of m PAGE 125 125 Once the PDF of m is obtained, it can be used to predict the RUL of the monitored panel. Since the PDF is updated at every SHM measurement, the predicted RUL will vary at every measurement interval N AS previously the RUL is calculated using MCS, 50,000 samples of m true Na and are generated, and Eq. (2.18) is used to calculate samples of Nf. In order to have a safe prediction of RUL the 5th percentile of Nf sam ples is used as a conservative estimate of RUL. Since we used synthetic data by adding random noise, the result may vary with different sets of data. Thus, the above process is repeated with 100 sets of measurement data and mean one standard deviation in tervals are plotted. Figure A 4 shows these conservative intervals of RUL with two different combinations of noise and bias. These combinations correspond to extreme cases; the most and least conservative estimates of RUL. In order to compare the predicted RUL with the true one, the true RUL is also plotted in the figure. Note that initially the difference between the true and predicted RULs is significant because uncertainty is large at early stage. However, the predicted RUL conv erges to the true one from the safe side as more updates are performed. In addition, the variability of estimated RUL is also gradually reduced. Thus, it can be concluded that the proposed Bayesian inference can estimate panel specific d amage growth parameters as well as can predict the RUL while maintaining conservativeness. PAGE 126 126 Figure A 4 Distribution (onesigma intervals) of 5 percentile (95% conservative) RUL obtained using 100 sets of measurements compared to the true RUL Figure A 5 shows the distribution of the error between the true RUL and the mean of the estimated distribution of RUL. As previously the behavior observed can be explained by the convergence of the distribution of C from the upper limit of the initial distribution. It can also be observed that the distribution of RUL narrows as the distribution of C is identified and it converges to the true RUL. Thus it can be concluded that the proposed Bayesian inference can es timate panel specific damage growth parameters as well as can predict the RUL. PAGE 127 127 Figure A 5 Distribution (onesigma intervals) of error between the true RUL and the m aximum likelihood of the estimated RUL distribution, positive error corresponding to unconservative estimate Comparing the results obtained using damage growth to those using damage size shown in Chapter 2 it can be observed that the updated distribution is wider and as a result the estimated RUL is m ore conservative. PAGE 128 128 APPENDIX B LEAST SQUARE FILTERE D BAYESIAN TO UPDATE JOINT PDF In the work presented up to that point, least square filtered Bayesian is used to identify a single damage growth parameter, m and assumed the other parameter to be already known. In practi ce, when both parameters are unknown, Bayesian inference needs to update the joint PDF of both parameters. In general, this can be achieved by dividing the ranges of uncertain parameters into a grid and calculate the joint PDF value at each grid point. If the range is divided by 100100 grid, the updating process includes 10,000 times calculation of likelihood, requires uncertainty propagation for given parameter values. As can be seen in Figure 22 the exponent m is the slope of the curve in the loglog scale and C is value of log(da /d N K ) = 1. As a first step in developing a prognosis methodology, we assumed that the accurate value of C is known, while that of m is uncertain now we assume both variables as being distri buted. Since the range of the variables m and C is generally known from literature or material handbooks, we assume that they are uniformly distributed between the lower and upper bounds. Then, the goal is to narrow the joint distribution using the Bayesi an statistics with measured damage sizes. The data used for the updating of material properties are crack size measurements The Bayesian inference equation used is the same as for m note that we chose to update log( C ) in the joint distribution instead of C because of the difference in order of magnitude: ,log,log ,log ,log,loglogini updt inilamCfmC fmC lamCfmCdmdC (B. 15) PAGE 129 129 The likelihood function is designed to integrate the information obtained from structural health monitoring ( SHM ) measurement to the knowledge we have about the distribution of the variables. The calculation of the likelihood function is very similar to the one that can be found in Chapter 2, the main difference is that the estimation is done for every couple ( mi, log( Cj). Figure B 1 shows the final update of the joint distribution, it can be observed that the correlation between m and C is very clear from that distribution even though no correlation was initially assumed. Figure B 1 Joint probability density function, illustration with one simulated set of measurements Once the joint distribution has been identified at cycle N, it can be used to predict th e remaining useful life (RUL). The distribution of RUL is calculated at every SHM measurement cycle N using MCS as well but with a larger sample than the one used to calculate the likelihood function, 50,000 samples of true crack sizes, m and C are generat ed using the updated distribution and the RUL is estimated in order to estimate the distribution of RUL. This allows us to estimate the distribution and from there obtain the 5th percentile. m log(C) 3.4 3.6 3.8 4 4.2 23.5 23 22.5 22 21.5 PAGE 130 130 Since we used synthetic data by adding random noise, the result m ay vary with different sets of data. Thus, the above process is repeated with 100 sets of measurement data and mean one standard deviation intervals are plotted. In order to show the value of our method we compare RUL calculated using the actual value: of m mtrue, and C Ctrue with the distribution (mean one standard deviation) of the 5th percentile of the distribution of RUL obtained using the updated joint distribution of m and C and the distribution (mean one standard deviation) of the 5th percentile of the distribution of RUL obtained using the updated distribution of m at each inspection, for the case of negative bias, 2mm and a noise of amplitude 1mm, this is shown in Figure B 2 It can be observed that LSFB method to update the joint PDF allows us to estimate the RUL converging to the true RUL from the conservative side but it does not improve much compared to the results obtained updating only m Figure B 2 Distribution (mean one standard deviation) of the 5th percentile of RUL for b = 2mm and V = 1mm, using LSFB to update the joint PDF Figure B 1 shows the error between the maximum likelihood of the estimated distribution of the RUL and the true RUL, note that positive values of the error PAGE 131 131 correspond to unconservative estimates. It can be observed that using LSFB method to update the joint PDF allows us to estimate the RUL converging to the true RUL from the conservative side, it is more conservative than when only m is updated but it does not converge as fast, the error is also larger Figure B 3 Distribution (onesigma intervals) of error between the true RUL and the maximum likelihood of the estimated RUL distribution for b = 2m m and V = 1mm, using LSFB to update the joint PDF, positive error corresponding to unconservative estimates Updating the joint distribution made the correlation between m and C very clear and one of the main conclusion that can be drown from the updated distribution is that there is not a single combination of parameters that will lead to the observed damage growth behavior. That strong correlation is also why updating one parameter while updating the other one leads to very similar results. We presented here a comparison between updating the two damage growth parameters in Paris law and updating one while assuming the other one as being constant. It can be concluded that we do not gain accuracy in the estimation of remaining useful by updating both parameters but it is computationally more expensive. PAGE 132 132 This can be explained by the strong correlation between the parameters. More results will be added to support those remarks, among which the effect of assuming the wrong value for C when updating m. PAGE 133 133 LIST OF REFERENCES 1. Boller C, Meyendorf N. Stateof the art in structural health monitoring for aeronautics. International Symposium on NDT in Aerospace. Furth/Bavaria, Germany; 2008. 2. Speckmann H, Henrich R. Structural health monitoring ( SHM) Overview on technologies under development. 16th World Conference on NDT, Montreal, Canada, 2004. 3. Staszewski W, Boller C, Tomlinson G. Health monitoring of aerospace Structures Wiley: New York, NY, 2004. 4. Chang P, Flatau A, Liu S. Review paper: health monitoring of civil infrastructure. Journal of Structual Health Monitoring 2003; 3(3):257267. 5. Wang CS, Wu F, Chang FK. Structural health monitoring from fiber reinforced composites to steel reinforced concrete. Smart Materials and Structure s 2001; 10:548552. 6. Farrar C, Lieven N. Damage prognosis: the future of structural health monitoring. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 2007; 365(1851):623632. 7. Boller C, Buderath M. Fatigue in aero structures where structural health monitoring can contribute to complex subject. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 2007; 365(1851):561587. 8. Barlow R, Hunter L. Optimum pre ventive maintenance policies. Operations Research 1960; 8:90100. 9. Nakagawa T. Sequential imperfect preventive maintenance policies. Journal of Applied Probability 1986; 23(2):536542. 10. Megaw E. Factors affecting visual inspection accuracy. Applied Ergonomics 1979; 10(1):2732. 11. Mital A, Govindaraju M, Subramani B. A comparison between manual and hybrid methods in parts inspection. Integrated Manufacturing Systems 1998; 9(6):344349. 12. Beral B, Speckmann H. Structural health monitoring (SHM) f or aircraft structures: a challenge for system developers and aircraft manufactures. In: Chang FK, editor. 4th Int Workshop on Structural Health Monitoring Stanford University, Stanford, CA: DEStech Publications, Inc; 2003. p. 1229. PAGE 134 134 13. Telgkamp J, Schmidt H. Benefits by the application of structural health monitoring (SHM) systems on civil transport aircraft. In: Chang FK, editor. 4th Int Workshop on Structural Health Monitoring. Stanford University, Stanford, CA: DEStech Publications, Inc; 2003. p. 285. 14. Giurgiutiu V, Zagrai A, Bao J. Piezoelectric wafer embedded active sensors for aging aircraft structural health monitoring. Structural Health Monitoring 2002; 1(1):4161. 15. Doebling S, Farrar C, Prime M, Shevitz D. Damage identification and health monitoring of structural and mechanical systems from changes in their vibration characteristics: a literature review. LA 13070MS, Los Alamos National Lab., NM (United States); 1996. 16. Farrar C, Worden K. An introduction to structural health monitoring. Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 2007; 365(1851):303315. 17. Gorinevsky D, Kim S, Boyd S, Gordon G, Beard S, Chang F. Optimal estimation of accumulating damage trend from a series of SHM images. 6th International Workshop on Structural Health Monitoring; 2007. 18. Tumer I, Bajwa A. A survey of aircraft engine health monitoring systems. AIAA Joint Propulsion Conference. Los Angeles, CA: Citeseer; 1999. 19. Giurgiutiu V, Cuc A. Embedded NDE for structural health monitoring, damage detection, and failure prevention. Shock and Vibration Reviews 2005; 37(2):83105. 20. Lee J, Takatsubo J, Toyama N, Kang D. Health monitoring of complex curved structures using an ultrasonic wavefield propagation imaging system. Measurement Science and Technology 2007; 18(12):38163824. 21. Luo J, Namburu M, Pattipati K, Qiao L, Kawamoto M, Chigusa S. Model based prognostic techniques [maintenance applications]. IEEE Systems Readiness Technology Conference; 2003. p. 330340. 22. Li C, Lee H. Gear fatigue crack prognosis using embedded model, gear dynamic model and fracture mechanics. Mechanical systems and signal processing 2005; 19(4):836846. 23. Berruet P, Toguyeni A, Craye E. Structural and functional approach for dependability in FMS. IEEE International Conference on Systems, Man, and Cybernetics; 1999. 24. Ray A, Patankar R. A stochastic model of fatigue crack propagation under variableamplitude loading. Engineering Fracture Mechanics 1999; 62(45): 477493. PAGE 135 135 25. Ray A, Tangirala S. A stochastic modeling of fatigue crack dynamics for online failure prognostics. IEEE Transactions on Control Systems Technology 1996; 4(4):443451. 26. Schwabacher M. A survey of datadriven prognostics. AIAA Infotech@Aer ospace Conference. Reston, VA; 2005. 27. Montanari G. Aging and life models for insulation systems based on PD detection. IEEE Transactions on Dielectrics and Electrical Insulation 1995; 2(4):667675. 28. Xue Y, McDowell D, Horstemeyer M, Dale M, Jordon J. Microstructure based multistage fatigue modeling of aluminum alloy 7075T651. Engineering Fracture Mechanics 2007; 74(17):28102823. 29. Goebel K, Saha B, Saxena A, Mct N, Riacs N. A comparison of three datadriven techniques for prognostics. 62nd Meet ing of the Society For Machinery Failure Prevention Technology (MFPT) Virginia Beach, VA, 2008. 30. Srivastava A, Das S. Detection and prognostics on low dimensional systems. IEEE Transactions on Systems, Man, and Cybernetics, Part C: Applications and Rev iews 2009; 39(1):4454. 31. Tipping M. The relevance vector machine in advances in neural information processing systems. MIT Press, Cambridge; 2000. 32. Yan J, Lee J. A hybrid method for online performance assessment and life prediction in drilling operations. IEEE International Conference on Automation and Logistics; 2007. p. 25002505. 33. Orchard M, Kacprzynski G, Goebel K, Saha B, Vachtsevanos G. Advances in uncertainty representation and management for particle filtering applied to prognostics. In ternational Conference on Prognostics and Health Management, Denver, CO, 2008. 34. Saha B, Goebel K. Uncertainty management for diagnostics and prognostics of batteries using bayesian techniques. IEEE Aerospace Conference, 2008. 35. Sheppard J, Kaufman M, Inc A, Annapolis M. Bayesian diagnosis and prognosis using instrument uncertainty. IEEE Autotestcon, 2005 2005:417423. 36. Engel S, Gilmartin B, Bongort K, Hess A. Prognostics, the real issues involved with predicting life remaining. IEEE Aerospace Co nference, Big Sky, MT, 2000. 37. Gu J, Barker D, Pecht M. Uncertainty assessment of prognostics of electronics subject to random vibration. AAAI Fall Symposium on Artificial Intelligence for Prognostics 2007. PAGE 136 136 38. Orchard M, Wu B, Vachtsevanos G. A partic le filter framework for failure prognosis. World Tribology Congress III. Washington, DC; 2005. 39. Paris P, Tada H, Donald J. Service load fatigue damage a historical perspective. International Journal of fatigue 1999; 21:3546. 40. Irwin G. The stresses and strains near the tip of a crack traversing a plate. Journal of applied Mechanics, ASME 1957. 41. Elber W. The significance of crack closure. Damage tolerance in aircraft structures, ASTM STP 1971; 486:230242. 42. Newman J, Elber W. Mechanics of Fatigue Crack Closure. ASTM International, 1988. 43. Newman J. Crack growth under variable amplitude and spectrum loading. TMS Symposium Proceedings to Honor P Paris 1997. 44. Paris P, Lados D, Tada H. Reflections on identifying the real DKeffective in the threshold region and beyond. Engineering Fracture Mechanics 2008; 75(34):299305. 45. Niu M. Airframe structural design. In: Conmilit Press LTD., editor. Fatigue, Damage Tolerance and Fail Safe Design. Hong Kong, 1990. p. 538570. 46. Zhao T, Jiang Y. Fatigue of 7075T651 aluminum alloy. International Journal of Fatigue 2008; 30(5):834849. 47. Center NR. Accuracy, Error, Precision, and Uncertainty available at http ://www.ndt ed.org/GeneralResources/ErrorAnalysis/UncertaintyTerms.htm 48. Organization RT. Crack size errors Nature and effect on POD estimation. available at http://ftp.rta.nato.int/public//PubFullText/RTO/TR/RTOTR AVT 051///TR AVT 051ANNG.pdf 49. Neath A, Cavanaugh J. Bayesian estimation of linear statistical bias. International journal of pure and applied mathematics 2006; 32(2):255. 50. Maeder M, Zuberbuhler A. Nonlinear least squares fitting of multivariate absorption data. Analytical Chemistry 1990; 62(20):22202224. 51. Collins J. Failure of Materials in Mechanical Design: Analysis, Prediction, Prevention. Wiley interscience, 1993. PAGE 137 137 52. G uan X, Liu Y, Saxena A, Celaya J, Goebel K. Entropy based probabilistic fatigue damage prognosis and algorithmic performance comparison. Annual Conference of the Prognostics and Health Management Society 2009, San Diego, CA, 2009. 53. Bjorck A. Numerical Methods for Least Squares Problems Society for Industrial Mathematics, 1996. 54. Welch G, Bishop G. An introduction to the Kalman filter. University of North Carolina at Chapel Hill, Chapel Hill, NC 1995. 55. Pfeiffer M. A Brief Introduction to Particle Filters. 2007, available at http://www.igi.tugraz.at/pfeiffer/documents 56. Gelman A, Carlin J, Stern H, Rubin D. Bayesian Data Analysis CRC Press, 2004. 57. Lam H, Lee Y, Sun H, Cheng G, Guo X Application of the spatial wavelet transform and Bayesian approach to the crack detection of a partially obstructed beam. ThinWalled Structures 2005; 43(1):1 21. 58. Li G, Yuan F, Haftka R, Kim N. Bayesian segmentation for damage image using MRF prior. Sensors and Smart Structures Technologies for Civil, Mechanical, and Aerospace Systems, SPIE Conference. San Diego, CA; 2009. 59. Gogu C, Haftka R, Leriche R, Molimard J, Vautrin A, Sankar B. Bayesian statistical identification of orthotropic elastic constants accounting for measurement and modeling errors. 11th AIAA Non Deterministic Approaches Conference, Palm Springs, CA, 2009. 60. Heller R, Stevens G. Bayesian estimation of crack initiation times from service data. Journal OF Aircraft 1978; 15:794798. 61. Zhang R, Mahadevan S. Model uncertainty and Bayesian updating in reliability based inspection. Structural Safety 2000; 22(2):154160. 62. Saxena A, Celaya J, Balaban E, Goebel K, Saha B, Saha S, Schwabacher M. Metrics for evaluating performance of prognostic techniques. International Conference on Prognostics and Health Management (PHM) Denver, CO, 2008. 63. Saxena A, Celaya J, Saha B, Saha S, Goebel K. On Applying the Prognostic Performance Metrics. Annual Conference of the Prognostics and Healt h Management Society, 2009, San Diego, CA, 2009. 64. Saxena A, Celaya J, Saha B, Saha S, Goebel K. Metrics for offline evaluation of prognostic performance. International Journal of Prognosis and Health Management 2010. PAGE 138 138 65. Vachtsevanos G, Wang P. Fault prognosis using dynamic wavelet neural networks. IEEE Systems Readiness Technology Conference; 2001. p. 857870. 66. Hoaglin D, Mosteller F, Tukey J. Understanding Robust and Exploratory Data Analysis. In: Sons JW, editor. Probability and Mathematical Stat istics, 1983. 67. Vachtsevanos G, Lewis F, M R, Hess A, Wu B. Intelligent Fault Diagnosis and Prognosis for Engineering Systems John Wiley & Sons: Hoboken, NJ, 2006. 68. Harkness H. Computational methods for fracture mechanics and probabilistic fatigue. Ph.D. Thesis Northwestern University 1994. 69. Lin K, Rusk D, Du J. Equivalent level of safety approach to damagetolerant aircraft structural design. Journal of Aircraft 2002; 39(1):167174. 70. An J, Acar E, Haftka R, Kim N, Ifju P, Johnson T. Being conservative with a limited number of test results. Journal OF Aircraft 2008; 45(6):19691975. 71. Smarslok B, Alexander D, Haftka R, Carraro L, Ginsbourger D. Separable Monte Carlo applied to laminated composite plates reliability. 49th AIAA/ASME/ASCE/A HS/ASC Structures, Structural Dynamics, and Materials Conference, Schaumburg, Ill, 2008. 72. Sankararaman S, Ling Y, Shantz C, Mahadevan S. Inference of equivalent initial flaw size under multiple sources of uncertainty International Journal of Fatigue 2011; 33(2):7589. 73. Lawson C, Hanson R. Solving Least Squares Problems Society for Industrial Mathematics, 1995. 74. DuQuesnay J, Underhill P. Fatigue life scatter in 7xxx series aluminum alloys. International Journal of Fatigue 2010; 32(2):398402. 7 5. Mukamai Y. Stress Intensity Factors Handbook Pergamon Press: New York, NY, 1987. 76. Moes N, Dolbow J, Belytschko T. A finite element method for crack growth without remeshing. International Journal of Numerical Methods in Engineering 1999; 46:131150. 77. Shih C, Asaro R. Elastic plastic analysis of crack on bimaterial interface, Part I: small scale yield. Journal of applied Mechanics, ASME 1988; 55:299316. 78. Daux C, Moes N, Dolbow J, Sukumar N, Belytschko T. Arbitrary branched and intersecting cracks with the extended finite element method. International Journal for Numerical Methods in Engineering 2000; 48:17411760. PAGE 139 139 79. Tanaka K. Fatigue crack propagation from a crack inclined to the cyclic tension axis. Engineering Fracture Mechanics 1974; 6( 3):493498. PAGE 140 140 BIOGRAPHICAL SKETCH Alexandra Coppe was born in 1984 in Cahors, France. She went to college at the Univers ity Paul Sabatier in Toulouse, France, where she graduated with a Bachelor of S cience in applied mathematics in 2005 and an M aster of S cience in applied mathematics in 2007. She moved to the United States in April 2007 to join the M ultidisciplinary O ptimization group in the D epartment of M echanical and A erospace E ngineering at the University of Florida as a visiting student in order to finish her m asters thesis. There, s he started her Ph.D. in Fall 2007 working Prof essor Raphael Haftka and Professor Nam Ho Kim. 