Wayne E King
email: weking@llnl.gov

Athanasios Arsenlis
email: arsenlis@llnl.gov

Charles Tong
email: tong10@llnl.gov

William L Oberkampf

Uncertainties in Predictions of Material Performance using Experimental Data that is Only Distantly Related to the System of Interest

Keywords: uncertainty quantification, model form uncertainty, model uncertainty, nuclear energy, neutron damage, ion damage, irradiation effects scaling, extrapolation

There is a need for predictive material "aging" models in the nuclear energy industry, where applications include life extension of existing reactors, the development of high burnup fuels, and dry cask storage of used nuclear fuel. These problems require extrapolating from the validation domain, where there is available experimental data, to the application domain, where there is little or no experimental data. The need for predictive material aging models will drive the need for associated assessments of the uncertainties in the predictions. Methods to quantify uncertainties in model predictions, using experimental data that is only distantly related to the application domain, are discussed in this paper.

Introduction

There is a growing need to make predictions of material performance in extreme environments and over very long periods of time where there is little or no experimental data. This is particularly the case in the prediction of the effects of "aging" on material performance where desired material lifetimes can exceed by a large margin what is practical for validation under normal application conditions. In large, complex engineering systems, the costs are often too high and/or the times too long to carry out desired validation experiments under actual operating conditions. In the case of aging, the required extrapolations can be orders of magnitude beyond the validation domain.

Further, model development can require reliance on accelerated experiments. Experiments can be accelerated by changing temperature to take advantage of Arrhenius behavior or the rate of application of the experimental forcing function can be increased. Results from accelerated experiments, which are also outside of the application domain, require extrapolation, perhaps over orders of magnitude in rate, to the actual operating conditions, by way of a model.

The need for predictive models is particularly acute in the nuclear energy industry, where applications include (i) the desire for life extension of existing reactors to 80 years, (ii) the development of high burnup fuels, and (iii) the imperative for dry cask storage possibly to hundreds of years.

The nuclear materials community has adopted an experimental strategy involving experiments that are not conducted in commercial power plants. The nuclear industry and regulators require a "sound and defensible case for the relevance of these techniques to actual service conditions" [START_REF] Crabtree | Science for Energy Technology: Strengthening the Link between Basic Research and Industry[END_REF]. This implies the ability (i) to make predictions across diverse irradiation energy spectra, irradiation rates, and irradiating particles including thermal neutrons, fast neutrons, and energetic ions, (ii) to scale the prediction to the relevant reactor conditions, and (iii) then to extrapolate to the application domain of interest.

Predictions are most useful in the presence of quantified uncertainties. In an engineered system, high uncertainty can lead to excessive conservatism and thus necessarily increased margins, which can adversely affect cost, schedule, and system performance. Thus, in addition to the predictive model, estimates of uncertainties in predictions are also required. Uncertainty quantification (UQ) provides a framework within which uncertainties for predictions can be estimated. The case that we are interested in here is uncertainties due to prediction in domains where there is little or no experimental data. This is generally referred to as model extrapolation.

Extrapolations using physics-based models differ from extrapolations using regression curve fits of the system responses quantities of interest. When making an assessment of uncertainties, there are several sources that must be considered, including model inputs and numerical solution approximation. When making an extrapolation using a model, the assumptions associated with the mathematical model itself result in a source of uncertainty usually referred to as model form uncertainty 1 , and it must also be considered [START_REF] Roy | A Comprehensive Framework for Verification, Validation, and Uncertainty Quantification in Scientific Computing[END_REF][START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF]. The inclusion of model form uncertainty represents a specialized field within UQ [START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF]. Approaches to dealing with uncertainties in model inputs are well established and can be implemented. The approach to uncertainties associated with model form is less well established, particularly for cases where models are assessed in a validation domain that does not fully overlap or overlap at all with the intended application domain.

In this paper, we discuss a predictive modeling problem in the nuclear energy area that requires large extrapolation. We also review possible approaches to the extrapolation problem. Throughout, we focus on the special case of predictions of models validated with experimental data that are only distantly related to the system of interest.

Background

2.1 Accelerated experimentation for nuclear energy materials applications [START_REF] King | Report on the Workshop on Accelerated Nuclear Energy Materials Development[END_REF] The nuclear industry needs models that predict the time dependence of microstructural and fission product evolution in structural materials and fuels.

The most challenging extreme environment to study is that of high irradiation dose. Models developed to address this extreme are difficult to validate because of the inability to reach these doses using existing neutron-irradiation facilities in reasonable amounts of time and at modest costs. Furthermore, reactor facilities are problematic experimental venues for combining the various aspects of the extreme environments into a quantitative in situ study of material behavior. Understanding radiation damage using ion irradiation is not a new idea. It has a long history of significant contributions spanning several decades. In fact, much of our understanding of material behavior under irradiation comes from well-controlled ion-irradiation experiments.

However, a key challenge is the scaling, or extension, of ion irradiation experiments and data to actual in-service conditions. Fig. 1 illustrates the particular case of scaling and extrapolation for the damage rate parameter. Scaling refers to use of models to bridge two unconnected validation domains. Extrapolation refers to the use of models to project into an application domain where there is no experimental data. The plot shows the range of damage expected for advanced reactors, GEN IV reactors, GEN III reactors, and light water reactors (LWRs). The plot further shows the damage levels that could be obtained in a 5-year irradiation experiment in a number of test reactors and the damage levels that could be obtained by ion irradiation in 5 days [START_REF] Fluss | Workshop on Science Applications of a Triple Beam Capability for Advanced Nuclear Energy Materials[END_REF]. Consequently, a scientifically defensible argument for the applicability of models developed using accelerated experiments to neutron irradiation environments is critically needed. This should include rate scaling, effects of recoil energy spectra, and the ability to extrapolate to dose regimes not explored by neutrons. This irradiation-effects scaling is identified as a priority research direction in the Science for Energy Technology Workshop Report [START_REF] Crabtree | Science for Energy Technology: Strengthening the Link between Basic Research and Industry[END_REF]. By definition, models developed for materials under neutron irradiation conditions at the extreme of high irradiation dose cannot be validated because little or no neutron irradiation data exists in that domain. Consequently, UQ will be particularly important in dose regimes that have not been explored using neutrons.

Other Relevant Science and Technology Application Areas

"Aging" of materials for nuclear energy applications is not the only area requiring scaling across and extrapolation outside of the validation domains. Fusion energy applications have a similar problem. In fusion machines, materials that can withstand very high levels of radiation damage are required. Without a fusion-relevant neutron source, the fusion materials community has adopted a research strategy similar to the fission community. Fig. 1. Plot of radiation damage, measured in displacements per atom (dpa), as a function of the rate at which the damage is produced [START_REF] Fluss | Workshop on Science Applications of a Triple Beam Capability for Advanced Nuclear Energy Materials[END_REF].

Uncertainty in the calculation of the depletion of nuclear fuels is composed of a number of individual components. Some of these include uncertainty in the cross section at a given neutron energy, uncertainty in core composition and other externally driven parameters such as power level and temperatures, and model approximations made to accommodate computer modeling capabilities.

These uncertainties will affect the prediction of the evolution of the fuel isotopic components in time. The isotopic component distribution at any given time represents the integration of the depletion conditions and uncertainties over all previous time. This effect can compound the effects of uncertainty or approximation. This compounding effect may also be limited by compensating effects. A common figure of merit for depletion is the energy extraction per mass of fuel, which is quantified as giga-watt-days per metric tonne of initial heavy metal (GWD/t). Present commercial fuel is depleted to about 60 GWD/t and the effects of uncertainties are well benchmarked for this range of operation. Advanced fuels may go beyond today's 60 GWd/t burnup. To operate outside the range of experience will require identifying the compounding effects for a given fuel type at a specific higher burnup so that they can be applied to the known uncertainties derived from inside the benchmarked region.

The effects of aging of materials are also of interest in both intermediate and long-term storage of nuclear waste. Analyses of the possible behavior of radioactive waste in a repository at Yucca Mountain (YM), Nevada, were conducted between 1982 and 2008. Early analyses (termed performance assessments or PAs) were for selecting the site and determining the feasibility of the disposal concept. Numerous parameter values for the numerical models were required and were mostly assigned by individual analysts and scientists for the early analyses. However, in 1987, Congress asked the YM Project (YMP) to evaluate the viability of a repository at YM [START_REF] Doe | Viability Assessment of a Repository at Yucca Mountain[END_REF]. For this viability assessment (PA-VA), the YMP formed five panels to examine: (1) groundwater flow in the unsaturated zone, (2) groundwater flow and radionuclide transport in the saturated zone, (3) the near-field effects of heat on the region around the engineered barrier system, (4) waste form degradation, and (5) waste package (WP) degradation. These five panels assigned parameter values by aggregating disparate data available in the literature, prior to completion of project experiments, and estimated the uncertainty present as literature-based information was often for conditions and spatial and temporal scales that differed from those required for the PA-VA. The analysis underlying the license application for the Yucca Mountain repository in 2008 [START_REF] Doe | Yucca Mountain Repository Safety Analysis Report[END_REF] considered a total of 392 uncertain analysis inputs (see Ref. [START_REF]Total System Performance Assessment Model/Analysis for the License Application[END_REF], Table K3-3, for a complete listing of these inputs and additional sources of detailed information). The Waste Isolation Pilot Plant (WIPP), a repository for transuranic radioactive waste in southern New Mexico, also had similar needs [START_REF] Helton | Characterization of Subjective Uncertainty in the 1996 Performance Assessment for the Waste Isolation Pilot Plant[END_REF][START_REF] Rechard | Assignment of Probability Distributions for Parameters in the 1996 Performance Assessment for the Waste Isolation Pilot Plant. Part 1: Description of Process[END_REF][START_REF] Rechard | Assignment of Probability Distributions for Parameters in the 1996 Performance Assessment for the Waste Isolation Pilot Plant. Part 2. Application of Process[END_REF].

Another closely related application is the Qualification Alternatives to the Sandia Pulsed Reactor (QASPR) project at Sandia National Laboratories [START_REF] Myers | The Science of QASPR (Qualification Alternatives to the Sandia Pulsed Reactor) in Sandia's Physical, Chemical, and Nano[END_REF][START_REF] Paulsen | QASPR Silicon Device Prototype Exercise and Results Summary[END_REF]. In this case, pulsed ion beams are being used to understand radiation effects in semiconducting materials. The challenge is to extrapolate those results to the relevant pulsed neutron environment in the absence of a relevant neutron source.

There is a related field known as accelerated testing. Accelerated tests are used to obtain timely information on product-life or performance degradation over time [START_REF] Meeker | Pitfalls of accelerated testing[END_REF]. Ref. [START_REF] Nelson | Accelerated Testing: Statistical Models, Test Plans and Data Analyses[END_REF] provides a comprehensive discussion of useful models and statistical methods for accelerated testing. Ideally, predictions from accelerated tests should be based on models of physics of failure. In practice, however, users of accelerated tests often use a combination of past experience and empirical fitting of data to statistical models. Although these procedures seem to have been adequate in the past, it is generally recognized that the path forward for large extrapolations must be based on physics-based modeling.

The Problem Recast in More General Terms

One of the goals of UQ is to provide a means to evaluate a model's predictive capability. Fig. 2(a) schematically illustrates the synergistic use of modeling, experiments, and UQ to make a prediction [START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF]. In this case, the system response quantity of interest is shown as a function of two system or environmental parameters. The application domain is highlighted in light brown, is the range of parameters #1 and #2 that are of interest in the application. The validation domain, highlighted in burgundy, is the range of parameters #1 and #2 where validation experimentation are carried out. The response surface is also shown with the application domain highlighted. In this example, the validation domain fully contains the application domain and predictions could be obtained through various types of interpolation over the validation domain. Likewise, uncertainties in predictions could be quantified by interpolation of uncertainties over the validation domain, or by direct calculation using the input uncertainties in the model. One of the main reasons why we rely on modeling and simulation is to make predictions in domains where there is little or no experimental or observational data. Fig. 2(b) shows the situation where the application domain is no longer within the validation domain and there is no overlap of the application domain with the validation domain. While there is a clear path to quantify uncertainties in a prediction when the application domain is inside the validation domain, the approach to quantifying uncertainties for the case shown in Fig. 2(b) is less clear.

Although Fig. 2 is illustrative, it fails to capture the complexity inherent in today's multiphysics simulations. Another view of the extrapolation issue is given in Fig. 3, which illustrates the validation hierarchy for a complex engineering system model [START_REF] Oberkampf | Verification and Validation in Computational Fluid Dynamics[END_REF]. At each level of the hierarchy, the problem is broken down into smaller and smaller pieces until it is reduced to the unit problem (or unit mechanism) level. This is a more physics-based, or system-based, perspective than in Fig. 2. Each box at the unit problem level contains individual physics models that can be validated using targeted experiments that may not be in the same conditions as would apply to the next higher level in the hierarchy. At each higher level, the individual effects are brought together forming coupled physics and coupled subsystems and systems interactions.

In the case of small extrapolations, we would expect this physics-based approach to enable extrapolation of uncertainties, including model form uncertainty, to the application domain. However, for large extrapolations, it is possible that new physics could appear at higher levels, or that unexpected coupling could emerge at the higher levels. These possibilities are present in the nuclear materials aging problem, where (i) scaling is required across diverse irradiation energy spectra, irradiation rates, and irradiating particles and (ii) large extrapolation is required to the application domain. Extrapolation of uncertainty includes both (a) extrapolation to the application domains where there is little or no experimental data, Fig. 2, and (b) effects present in higher levels of system complexity in the validation hierarchy, Fig. 3, that are not anticipated in the modeling [START_REF] Cacuci | Sensitivity and Uncertainty Analysis, Data Assimilation, and Predictive Best-Estimate Model Calibration[END_REF]. Aero/Thermal Fig. 3. A hierarchical validation structure for the hypersonic cruise missile (from [START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF]).

Prediction-Coupling of Accelerated Experiments with Physics Models [4]

It takes a long time to develop a new material or investigate the properties of materials that undergo low-dose-rate irradiation, such as the pressure vessel or core internals including the fuel and cladding. Experiments using neutron irradiation can take up to 7 years, including the irradiation time, the radioactive cool-down time, and the post-irradiation examination. Incorporating the effects of high temperature, stresses and a corrosive environment along with irradiation make the problem multidimensional and extremely complicated. Translating that into a program to satisfy a regulatory requirement for a new material or new fuel design can lead to a multi-decadal process. Such a timescale is unacceptable to efficient progress, and yet, it is the present-day norm.

One pathway to accelerate this process is to carry out accelerated experiments either inside the reactor core (in the case of the low-dose-rate regimes for the pressure vessel) or using external radiation sources such as ion beams (in the case of core internals that would see high neutron doses over their lifetimes). Using ion beams, one can investigate a large parameter space (in terms of external forcing functions) for irradiation effects on microstructure and macroscale properties. The phase space includes temperature, ion type, ion dose rate, ion energy, and total dose. It also includes the ability to apply in situ mechanical loading, chemical environments, and coolant fluids.

In some special cases, it is possible to create microstructures and material properties that are very similar to those that would be found in a particular nuclear reactor irradiation experiment. However, microstructures, mechanical properties, or other physical properties that deviate from neutron irradiations significantly, will also be observed. The challenge is to employ all of those observations to develop a science-based understanding of material degradation and performance, specifically, to establish a scientific basis for the key mechanisms of material performance.

The best approach to quantifying such scientific understanding is to build a model that captures all of the relevant physics, which is where modeling and simulation come into play. One seeks to understand the ion-beam forcing function and the material response to that forcing function. With that understanding, if a different boundary condition was applied (in terms of say temperature, ion type, or dose rate), it is reasonable to expect to be able to predict the material response. That is, one could have sufficient confidence that, with this robust model, interpolation and reproduction of an experimental result is possible.

The question then becomes: How does one extrapolate, given a model, to high-dose neutron-irradiation environments? Compared with ion irradiations, neutron dpa rates are much, much lower (∼ 10 2 -10 3 lower). In addition, physical mechanisms, such as transmutation and chemistry changes, occur simultaneously with the neutron bombardment and displacements in a material. With a robust model, the boundary conditions can be altered, while not perturbing any of the model internals, and an extrapolation can be made to project material performance to the neutron-irradiation environment. To make that extrapolation, researchers must also quantify the quality or accuracy of the extrapolation. This forward extrapolation and the qualification of the quality of the predicted extrapolation is where uncertainty quantification becomes important.

Approaches to Extrapolation of Uncertainty

Calibrate Model Parameters Over the Validation Domain and Ignore Model Form Uncertainty

This approach uses the physics model as it is (assuming no model form uncertainties) together with available experimental data for validation and calibration (given a number of model input parameters). The posterior distributions of the calibration parameters are then used in forward uncertainty propagation (through the computational model) to predict the extrapolated configuration and to estimate the corresponding uncertainties. To speed up calibration, response surface or surrogate modeling using Gaussian process, polynomial regression or polynomial chaos is often used. This approach has the inherent assumption that the computational model has captured all essential physics (except that there are uncertainties about some physics parameters that can be estimated using data) and the surrogate models are adequate for the extrapolated regime.

A simple example of this method is shown in Fig. 4, where a prediction and an associated uncertainty are required for the time required for an object to drop

Experimental data:

The time it takes an object to drop from each of 6 floors of a tower is recorded.

There is an uncertainty in the measured drop times of about ±0.2 seconds. Predictions for times are desired for drops from floors 7 through 10, which do not exist yet. A simple UQ example using experimental data and a computational model to predict drop times from new heights [START_REF] Higdon | Uncertainty Quantification and Error Analysis[END_REF].

from yet unconstructed floors of a building [START_REF] Higdon | Uncertainty Quantification and Error Analysis[END_REF]. Measurements are made of the time to drop the object from the first six floors and a model is constructed. Using the available experimental data from the drops, a Bayesian inference methodology is used to update/calibrate the uncertain input parameters in the model. Once the model is calibrated it is then used for prediction outside the calibration range. This is the most commonly used form of extrapolation: the model form uncertainty is assumed to be zero. Continuing with the example in Fig. 4, when comparing the model against experimental data, systematic errors are observed which cannot be resolved satisfactorily by calibration. Model bias (or discrepancy) and its dependence on the drop height are evident. As shown in Fig. 5, an extrapolation would have led to a predicted time and uncertainty that, in fact, would not have predicted what would have been observed had an experiment been carried out. This illustrates that in this case, as in many others, the model form uncertainty can dominate the extrapolated uncertainties. The experiments revealed the presence of physics that was not accounted for in the original model.

Calibrate Parameters and Identify a Discrepancy Function to Characterize Inadequacy

This approach explicitly assumes a functional form for the discrepancy between the simulation and the actual physical process. The full method estimates from data simultaneously both the posteriors for the calibration parameters as well as the parameters in the discrepancy function. Some simplifications were suggested to reduce the complexity of the estimation process (see [START_REF] Kennedy | Bayesian Calibration of Computer Models[END_REF] for details). This approach assumes both the functional forms of the simulation model and the discrepancy stay the same in the extrapolated regime. Kennedy and O'Hagan proposed a Bayesian approach that represents the model form uncertainty using a discrepancy function (in terms of some input parameters) to characterize model inadequacy [START_REF] Kennedy | Bayesian Calibration of Computer Models[END_REF]. This discrepancy function is created using experimental data as well as a selected regression or statistical emulator such as Gaussian process. Predictions are performed by incorporating the discrepancy function evaluated at the extrapolated points, in addition to the uncertainty due in posterior distributions of the calibrated parameters. This approach has the inherent assumption that the discrepancy function essentially captures the misrepresented physics.

In their approach, the system is modeled by

z = ζ(x) + e = ρη(x, θ) + δ(x) + e (1)
where x is a vector of input parameters; θ is a vector of calibration parameters; η(x, θ) denotes the function for the simulation model; z is the observation; e is the observation error (independent normal distribution); ζ(x) is the true value of the process being modeled; ρ is an unknown regression parameter to be determined; and δ(x) is a function to describe model inadequacy, which is independent of η(x, θ).

The full Bayesian calibration of this system is very complicated. Instead, Kennedy and O'Hagan proposed a multi-step approach:

-Build a Gaussian process model for η(x, θ) based on sampling different values of x and θ (that is, to estimate the "hyper-parameters" used to describe a Gaussian process model).

-Use data ({z}) to estimate the regression parameter ρ; the standard deviation of the observation error e; the hyper-parameters in the Gaussian process model of the model inadequacy function δ(x). -Use the data ({z}) and the model (z = ρη(x, θ) + δ(x) + e) for calibration to get the posterior distribution of θ. -Use the model (z = ρη(x, θ)+δ(x)+e) and the posteriors of θ for prediction and uncertainty analysis.

Model inadequacy is defined as "the difference between the true mean value of the real world process and the code output at the true value of the inputs" [START_REF] Kennedy | Bayesian Calibration of Computer Models[END_REF]. There is debate across the field regarding the use of this definition of model inadequacy to describe the difference between a simulation and an experiment (see Sects. 2.2.3, 2.4, and Sects. 12. 1-12.3 in [3].

This method is suitable if there is (a) sufficient experimental data from different input parameter configurations to characterize the discrepancy function, (b) reason believe that the assumption concerning the discrepancy functional form will be valid in the extrapolation regime, and (c) that the discrepancy function is more significant than the observation error.

An alternative approach is where the candidate response surface methods for η(x, θ) in Eq. ( 1) and the discrepancy functions are based on generalized polynomial chaos (or stochastic collocation) methods [START_REF] Eldred | Comparison of Non-Intrusive Polynomial Chaos and Stochastic Collocation Methods for Uncertainty Quantification[END_REF]. The advantages of polynomial chaos UQ methods are their efficiency and their utility for representing and propagating large uncertainties through complex models [START_REF] Najm | Uncertainty Quantification and Polynomial Chaos Techniques in Computational Fluid Dynamics[END_REF]. Both intrusive and non-intrusive applications of the polynomial chaos method are reviewed in [START_REF] Najm | Uncertainty Quantification and Polynomial Chaos Techniques in Computational Fluid Dynamics[END_REF].

This method is suitable if there is sufficient experimental evidence that the discrepancy function has polynomial form in the input parameters, and that this form will also be valid in the extrapolation regime. Moreover, another requirement is that the discrepancy function be more significant than the observation error.

Validation Metric Approach [2,3]

Oberkampf and Roy argue that model form uncertainty should be estimated as part of the process of model validation. They estimate the model form uncertainty in the validation domain using a validation metric which they define as "a mathematical operator that requires two inputs: the experimental measurements of the system response quantities of interest and the prediction of the system response quantities at the conditions used in the experimental measurements" [START_REF] Roy | A Comprehensive Framework for Verification, Validation, and Uncertainty Quantification in Scientific Computing[END_REF]. Oberkampf and Roy in Sects. 13.2, 13.4, and 13.5 of [START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF] and in [START_REF] Ferson | Model Validation and Predictive Capability for the Thermal Challenge Problem[END_REF], describe two validation metrics: the confidence interval approach and the method of comparing cumulative distribution functions (CDFs) from the model and the experiment. In the confidence interval approach, they define the validation metric for model form uncertainty as the difference between the mean of model prediction and the estimated mean of the experimental data. In the CDF method, the validation metric is defined as the area between the experimental and simulation CDFs.

Once a validation metric is estimated over the validation domain, the critical issue is how this error structure should be extrapolated to the application conditions of interest. One simple method for extrapolation is to construct a regression fit of the error structure over the validation domain using a low degree polynomial function [START_REF] Roy | A Comprehensive Framework for Verification, Validation, and Uncertainty Quantification in Scientific Computing[END_REF][START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF][START_REF] Ferson | Model Validation and Predictive Capability for the Thermal Challenge Problem[END_REF]. The regression function is then evaluated at the application conditions, along with the statistical prediction interval at those conditions. The estimate of the model form uncertainty is increased by the prediction interval not only because of the imprecision of the regression function to fully represent the model form uncertainty, but also because of the random measurement uncertainty that is present in each experimental measurement. A level of statistical confidence is chosen for the prediction interval, say 90 or 95%, and the upper bound on the prediction interval is used as the estimate of the model form uncertainty at the application conditions of interest.

This estimated model form uncertainty is considered as an epistemic uncertainty, i.e., an uncertainty whose source is lack of knowledge as opposed to randomness, for the prediction of the system response quantities of interest at the application conditions. It has been found [START_REF] Roy | A Comprehensive Framework for Verification, Validation, and Uncertainty Quantification in Scientific Computing[END_REF][START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF][START_REF] Ferson | Model Validation and Predictive Capability for the Thermal Challenge Problem[END_REF] that even if the model form uncertainty is relatively small over the validation domain, but the magnitude of the extrapolation is large in the multi-dimensional input space over which data are available, the estimate of model form uncertainty is typically quite large at the application conditions of interest. The model form uncertainty is clearly represented to the user of the simulation results, e.g., a designer or decision maker, as a probability-box, or p-box. The p-box is an interval-valued CDF, where the range of possible probabilities of the system response quantity reflects the epistemic uncertainty due to the model form.

Method of Alternate Plausible Models [3]

An approach for assessing uncertainty, both model form uncertainty and parametric uncertainty, is to compare predictions from alternative plausible models. This method is also referred to as the method of competing models. While simple in concept, it is not commonly used because of the time and expense of developing multiple models for a system. Examples of applications of this method include hurricane forecasting (Fig. 6) (e.g., see [START_REF] Vijaya Kumar | Multimodel Superensemble Forecasting of Tropical Cyclones in the Pacific[END_REF][START_REF] Krishnamurti | Improved Weather and Seasonal Climate Forecasts from Multimodel Superensemble[END_REF]), climate prediction, and long-term storage of nuclear waste.

The approach requires multiple models developed by independent groups. This approach does not actually provide an estimate of model form uncertainty; it only provides an indication of the similar or dissimilar nature of each model prediction. Because of the cost and time involved, this approach will likely be limited to application in matters of very high priority.

Since in many applications simulating the full physics is very consuming in terms of time and computational resources, a viable approach is to use simplified models for some of the physics components in the system (for example, see [START_REF] Oden | Control of Modeling Error in Calibration and Validation Processes for Predictive Stochastic Models[END_REF][START_REF] Qian | Building Surrogate Models Based on Detailed and Approximate Simulations[END_REF]). For example, in computational fluid dynamics, a popular simplification Fig. 6. Superensemble track forecast of Hurricane Lenny with predicted tracks of some member models and associated superensemble track shown [START_REF] Vijaya Kumar | Multimodel Superensemble Forecasting of Tropical Cyclones in the Pacific[END_REF].

is the use of the Reynolds-averaged Navier Stokes (RANS) model in place of the more complex large eddy simulation (LES) or even direct numerical simulation (DNS). These simplified models are often benchmarked or calibrated against the more complex counterparts. Outside the benchmark/calibration regime, the predictions from these two models of different physics fidelity can be compared. As a result, the model form uncertainty in the lower fidelity model could be estimated by comparison of the predictions with the higher fidelity model predictions at a limited number of conditions that are similar to the application conditions of interest. If it is concluded that the lower fidelity model accuracy is judged to be inadequate for the application of interest, then one may (a) increase the modeling fidelity of the lower fidelity model so as to attain the needed accuracy, or (b) characterize the model form uncertainty in some appropriate way so that the predictive uncertainty is recognized by the user of the simulation results.

Discussion

The Role of Model Form Uncertainty

For the case of extrapolations to application domains outside of the validation domain or to a higher level of system model, it is likely that model form uncertainty will dominate the extrapolated uncertainty. Extrapolating model form uncertainty is complex because it is extrapolating the error structure of a model, combined with the uncertainty in the experimental data, in a high dimensional space. Uncertainty can be reduced, compared with for example, extrapolating a regression fit of the measured system response quantities, by taking advantage of the physics incorporated into the model [START_REF] Oberkampf | Verification and Validation in Scientific Computing[END_REF]. Take for example Fig. 7, where there is no physical basis to the quadratic or linear models used to describe the data, and then extrapolate to large values of x. The fitted models over the range of the data are indistinguishable but in the range of extrapolation, the quadratic falls far outside of the confidence interval (CI) of the linear fit. A physics based model would help constrain the extrapolated uncertainty.

Use of UQ to Manage Extrapolation Uncertainties

The goal of a simulation is to produce a prediction along with an estimate of the effect of all of the relevant uncertainties on the system response quantities of interest. In addition to estimating the uncertainty of an extrapolation, UQ can serve to reduce that uncertainty through methodologies and mathematical methods for [START_REF] Tong | PSUADE User's Manual (Version 1.2.0)[END_REF]:

tuning (or calibrating) a simulation model to match with experimental results, establishing the integrity of (i.e., validate) a simulation model, assessing the region of validity of a simulation model, characterizing the output uncertainties of a simulation model, identifying the major sources of uncertainties of a model, providing information on which additional experiments are needed to improve the understanding of a model.

These methodologies can be used to understand how the model performs over the validation domain and improve it, if necessary. As in Fig. 5, these improvements serve to improve the quality of the extrapolation and reduce the uncertainty in the extrapolation. While there might be significant uncertainty on the magnitude of extrapolated uncertainties, the fact that the uncertainty has been reduced and by what fraction is useful in itself.

Missing Physics (Unknown-Unknowns)

In distant extrapolations, there is always the potential for missing physics to be present either as new unit mechanisms or as new coupling of mechanisms at higher levels in the validation hierarchy, see Fig. 3. Therefore, we discuss two elements aimed at reducing the unquantified uncertainty arising from missing physics for the case where there is no recognized disagreement between experiment and simulation in the validation domain. (If there is disagreement in the validation domain, as there is in Fig. 5, an unknown-unknown becomes a knownunknown).

The Role of Peer Review to Address Missing Physics. Because extrapolation is more of a physics endeavor than a statistics endeavor, scientific peer review plays a critical role in the extrapolation process. A nearly analogous issue was faced by Theofanous and co-workers in the application of the risk oriented accident analysis methodology (ROAAM) for low probability, high consequence hazards [START_REF] Theofanous | On the Proper Formulation of Safety Goals and Assessment of Safety Margins for Rare and High-consequence Hazards[END_REF]. The basic premise is that once the selected sample of the community of experts in the problem area is convinced that the model reflects to the extent possible all of the relevant physics, the problem may be considered characterized. By this we mean that what is obtained is the best that can be done at the present time with the committed resources. One outcome may be that additional resources need to be committed to the problem or that additional resolution simply cannot be obtained, as there is no known path forward to gain such resolution. This peer review process must be traceable and scrutable [START_REF] Theofanous | On the Proper Formulation of Safety Goals and Assessment of Safety Margins for Rare and High-consequence Hazards[END_REF].

The Importance of Data Assimilation to Mitigate Missing Physics. Unfortunately, as systems become more complex, a point is reached at which whole system models simply cannot be validated, in the sense of comparison with experimental results. The potential for new unit processes or coupling at higher levels of the validation hierarchy must be acknowledged when making an extrapolation [START_REF] Cacuci | Sensitivity and Uncertainty Analysis, Data Assimilation, and Predictive Best-Estimate Model Calibration[END_REF]. This does not diminish the value of the extrapolated uncertainties but requires additional attention through data assimilation.

Take for example the case of hurricane forecasting (see, for example, [START_REF] Cacuci | Sensitivity and Uncertainty Analysis, Data Assimilation, and Predictive Best-Estimate Model Calibration[END_REF][START_REF] Kalnay | Atmospheric Modeling, Data Assimilation, and Predictability[END_REF]). These forecasts save millions of dollars by limiting evacuation areas. But in the early stages of the prediction of the track of a storm, models sometimes predict tracks that do not coincide with the actual track that the storm eventually follows. Hurricane forecasters effectively use data assimilation to constantly update their model predictions and uncertainties. As data is assimilated over time, the uncertainties decrease because the period of prediction (i.e., extrapolation) becomes shorter as a hurricane nears the region of interest (i.e., the location where it will make landfall).

The same will be true for prediction of aging of nuclear reactor materials. As materials are irradiated in a real fission environment, an accompanying data assimilation effort must be in place, the model must be constantly updated, missing physics revealed, and as a result, uncertainties in predictions can be expected to decrease.

Conclusions

The problem posed by the nuclear industry, referred to as irradiation effects scaling [START_REF] Crabtree | Science for Energy Technology: Strengthening the Link between Basic Research and Industry[END_REF], is an ideal example of a high impact uncertainty quantification problem requiring scaling across validation domains and extrapolations to application domains where there is little or no experimental data. The need for predictive material aging models will drive the need for associated uncertainties in the predictions.

The case of extrapolation is more of a physics endeavor than a statistics endeavor. The goal is to produce a prediction with scientifically defensible and acceptable uncertainty. Most extrapolation methods do not deal with missing physics, i.e., they only estimate (extrapolate) the error structure of knownunknowns. Therefore, the process involves using validating experimentation and more detailed physics-based models that capture the essential physics thus enabling the required scaling and extrapolation and reducing uncertainties.

The idea that uncertainty increases when extrapolating outside of the validation domain is clear. However, exactly how the uncertainty increases is not well understood, in addition to being model and situation dependent [START_REF] Cacuci | Sensitivity and Uncertainty Analysis, Data Assimilation, and Predictive Best-Estimate Model Calibration[END_REF]. This is likely due to the fact that the uncertainty increase is very tightly coupled with the physics-basis of the model. Inaccuracies inherent to models that approximate the relevant physics, i.e., model form uncertainty, will likely dominate the uncertainties. Methods to quantify uncertainties in predictions of models, particularly model form uncertainty, are lacking and are topics requiring further fundamental research.

  a p o la ti o n

Fig. 2 .

 2 Fig. 2. Possible relations of the validation domain to the application domain. (a) Complete overlap and (b) no overlap [3].

  Simulated drop times:A (trivial) computational model is developed to predict the drop times as a function of height. The simulated drop times (red line) are systematically too low when compared with the experimental data (symbols). the error bars around the observed drop times show the observation uncertainty.

  Fig. 4.A simple UQ example using experimental data and a computational model to predict drop times from new heights[START_REF] Higdon | Uncertainty Quantification and Error Analysis[END_REF].

  Neglecting model form uncertainty results in an extrapolation of the model and the uncertainty that does not intersect the experimental data.

Fig. 5 .

 5 Fig.5. The impact of ignoring model form uncertainty on the extrapolation and uncertainties[START_REF] Higdon | Uncertainty Quantification and Error Analysis[END_REF].

Fig. 7 .

 7 Fig. 7. Model dependence of the extrapolation with equivalent goodness of fit in the validation domain [27].

Model form uncertainty is referred to as model uncertainty, model bias, or structural uncertainty.

Acknowledgements. The authors are grateful to (i) Jon Helton of Sandia National Laboratories for his contributions to the document in general and to Sect. 2.2 in particular, (ii) Rob Rechard of Sandia National Laboratories for providing input on the Yucca Mountain Project, (iii) David Higdon of Los Alamos National Laboratory for permission to use Fig. 4 and Fig. 5, (iv) T. N. Krishnamurti of The Florida State University for permission to use Fig. 6, (vi) Gary King of Harvard University for permission to use Fig. 7 (vii) Richard Kochendarfer AREVA Federal Services LLC for a contribution to Sect. 2.2, (viii), Len Lorence of Sandia National Laboratories for providing input on the QASPR project, and (ix) Tony O'Hagen of Sheffield University, William Meeker of Iowa State University, J. Tinsley Oden of The University of Texas at Austin, Habib

Najm of Sandia National Laboratories, and Janet Allen of University of Oklahoma for extensive technical input. This work performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344.


Rafi L Muhanna
email: rafi.muhanna@gtsav.gatech.edu

Robert L Mullen

Interval Based Finite Elements for Uncertainty Quantification in Engineering Mechanics

Keywords: interval, finite elements, reliability, aleatory, epistemic, probability bounds

This paper illustrates how interval analysis can be used as a basis for generalized models of uncertainty. When epistemic uncertainty is presented as a range and the aleatory is based on available information, or when random variables are assigned an interval probability, the uncertainty will have a Probability Bound (PB) structure. When Interval Monte Carlo (IMC) is used to sample random variables, interval random values are generated. Interval Finite Element Method (FEM) is used to propagate intervals through the system and sharp interval solutions are obtained. Interval solutions are sorted and PBs of the system response are constructed. All relevant statistics are calculated characterizing both aleatory and epistemic uncertainty. The above mentioned sequence is presented in this work and illustrative examples are solved.

Introduction

The presence of uncertainty in all aspects of life is evident. However, quantifying uncertainty is always advancing. There are various ways in which types of uncertainty might be classified. One is to distinguish between "aleatory" (or stochastic) uncertainty and "epistemic" uncertainty. The first refers to underlying, intrinsic variability of physical quantities, and the latter refers to uncertainty which might be reduced with additional data or information, or better modeling and better parameter estimation [START_REF] Melchers | Structural Reliability Analysis and Prediction[END_REF].

Probability theory is the traditional approach to handling uncertainty. This approach requires sufficient statistical data to justify the assumed statistical distributions. Analysts agree that, given sufficient statistical data, probability theory describes the stochastic uncertainty well. However, probabilistic modeling cannot handle situations with incomplete or little information on which to evaluate a probability, or when that information is nonspecific, ambiguous, or conflicting [START_REF] Ferson | Different Methods Are Needed to Propagate Ignorance and Variability[END_REF], [START_REF] Sentz | Combination of Evidence in Dempster-Shafer Theory[END_REF],and [START_REF] Walley | Statistical Reasoning with Imprecise Probabilities[END_REF]. Many generalized models of uncertainty have been developed to treat such situations, including fuzzy sets and possibility theory [START_REF] Zadeh | Fuzzy Sets as a Basis for a Theory of Possibility[END_REF],

Dempster-Shafer theory of evidence [START_REF] Dempster | Upper and Lower Probabilities Induced by a Multi-Valued Mapping[END_REF], [START_REF] Shafer | A Mathematical Theory of Evidence[END_REF], random sets [START_REF] Kendall | Foundations of a Theory of Random Sets[END_REF], probability bounds [START_REF] Berleant | Automatically Verified Reasoning with Both Intervals and Probability Density Functions[END_REF], [START_REF] Ferson | Different Methods Are Needed to Propagate Ignorance and Variability[END_REF], and [START_REF] Ferson | Constructing Probability Boxes and Dempster-Shafer Structures[END_REF], imprecise probabilities [START_REF] Walley | Statistical Reasoning with Imprecise Probabilities[END_REF], convex models [START_REF] Ben-Haim | Convex Models of Uncertainty in Applied Mechanics[END_REF], and others.

These generalized models of uncertainty have a variety of mathematical descriptions. However, they are all closely connected with interval analysis [START_REF] Moore | Interval Analysis[END_REF]. For example, the mathematical analysis associated with fuzzy set theory can be performed as interval analysis on different α levels [START_REF] Lodwick | Special Issue: Interface between Fuzzy Set Theory and Interval Analysis[END_REF] and [START_REF] Muhanna | Development of Interval Based Methods for Fuzziness in Continuum Mechanics[END_REF]. Fuzzy arithmetic can be performed as interval arithmetic on α cuts. A Dempster-Shafer structure [START_REF] Dempster | Upper and Lower Probabilities Induced by a Multi-Valued Mapping[END_REF] and [START_REF] Shafer | A Mathematical Theory of Evidence[END_REF] with interval focal elements can be viewed as a set of intervals with probability mass assignments, where the computation is carried out using the interval focal sets. Probability bounds analysis [START_REF] Berleant | Automatically Verified Reasoning with Both Intervals and Probability Density Functions[END_REF], [START_REF] Ferson | Different Methods Are Needed to Propagate Ignorance and Variability[END_REF], and [START_REF] Ferson | Constructing Probability Boxes and Dempster-Shafer Structures[END_REF] is a combination of standard interval analysis and probability theory. Uncertain variables are decomposed into a list of pairs of the form (interval, probability). In this sense, interval arithmetic serves as the calculation tool for generalized models of uncertainty. A short description of probability bounds is given in the next section.

Probability Bounds

Probability Bounds (PB) identifies a specific mathematical framework for analysis when precise discrete probabilities (or PDF) are not completely known [START_REF] Walley | Statistical Reasoning with Imprecise Probabilities[END_REF]. Probability bounds are normally associated with epistemic sources of uncertainty where the available knowledge is insufficient to construct precise probabilities. Calculations on PB can be conducted using Monte Carlo methods combined with an interval finite element method [START_REF] Zhang | Interval Monte Carlo Methods for Structural Reliability[END_REF], or by using various discretization methods.

The arithmetic for piecewise constant discretization on PB can be found in publications such as [START_REF] Yager | Arithmetic and Other Operations on Dempster-Shafer Structures[END_REF], [START_REF] Williamson | Probabilistic Arithmetic I: Numerical Methods for Calculating Convolutions and Dependency Bounds[END_REF], [START_REF] Berleant | Bounding the Results of Arithmetic Operations on Random Variables of Unknown Dependency Using Intervals[END_REF], [START_REF] Ferson | Constructing Probability Boxes and Dempster-Shafer Structures[END_REF]. Let the behavior of the system be modeled by a function y = f (x), where x = (x 1 , x 2 , . . . , x n ) is the parameter vector, and each x i is represented by a probability-bounds structure. The CDF of a particular value y * , F (y * ), is to be determined. The objective of probabilitybounds analysis is to obtain an interval to bound the precise probability F (y * ) of interest (in the classical sense). The numerical implementation of probabilitybounds analysis can be done using interval analysis and the Cartesian product method. The general procedure can be found in publications such as [START_REF] Tonon | Reliability Analysis of Rock Mass Response by Means of Random Set Theory[END_REF], [START_REF] Tonon | A Search Algorithm for Calculating Validated Reliability Bounds[END_REF], [START_REF] Joslyn | Convergence Properties of an Interval Probabilistic Approach to System Reliability Estimation[END_REF], and [START_REF] Ferson | Constructing Probability Boxes and Dempster-Shafer Structures[END_REF].

For example, Fig. 1 depicts probability-bounds circumscribing a normal distribution whose mean is sure to lie within the interval [5.6, 6] and whose standard deviation is 1. Such bounds can result from the addition of a normally distributed random variable with a mean of 5.6 and a standard deviation of 1 with a variable bounded between 0 and 0.4.

Probability-bounds structure can also arise by forming probability distributions of intervals. In this context, probability-bounds structure is mathematically analogous to a standard discrete probability distribution except that the probability mass is assigned to an interval rather than to a precise point; thus, the probability-bounds structure can be specified as a list of pairs of the form (interval, probability mass).

The advantage of the probability-bounds approach is that it can capture a wider range of uncertainties than the standard probabilistic approach. On the other hand, the standard probability distribution and interval number are two degenerate cases of probability-bounds structure; thus the probability-bounds approach provides a general framework for handling problems with a mixture of interval-based information and standard probabilistic information.

When and if additional knowledge of a system is obtained, the probability bounds can be refined. Just as interval and scalar calculations can be easily mixed, conventional precise probabilities and probability bounds may be mixed in a single calculation. However, how to accomplish this in a computationally efficient method for finite element analysis is unresolved. In this work interval Monte Carlo (IMC) will be used.

Interval Monte Carlo

When a random variable is described in a probability-bound structure, one natural approach for sampling such a random variable is the use of interval Monte Carlo (IMC). A discretization approach can be used as well [START_REF] Zhang | Interval Monte Carlo Methods for Structural Reliability[END_REF]; however, this approach has not yet been developed for the general case. IMC has been proposed to generate fuzzy numbers in [START_REF] Ferson | Hybrid Arithmetic[END_REF] and for the first time for structural reliability in [START_REF] Zhang | Interval Monte Carlo Methods for Structural Reliability[END_REF].

The first step in the implementation of interval Monte Carlo simulation is the generation of intervals in accordance with the prescribed probability bounds. The inverse transform method is often used to generate random numbers. Consider a random variable x with CDF F (x). If (u 1 , u 2 , . . . , u m ) is a set of values from the standard uniform variate, then the set of values

x i = F -1 x (u i ), i = 1, 2, . . . , m (1)
will have the desired CDF F (x). The inverse transform method can be extended to perform random sampling from a probability bound. Suppose that an imprecise CDF F (x) is bounded by F (x) and F (x), as shown in Fig. 2. For each u i in Eq. ( 1), two random numbers are generated

x i = F -1
x (u i ) and

x i = F -1 x (u i ) (2)
Such a pair of x i and x i form an interval [x i , x i ] which contains all possible simulated numbers from the ensemble of distributions for a particular u i . The method is graphically demonstrated in Fig. 2 for the one-dimensional case. The next step is to solve for the generated interval values using interval finite elements described in the next section.

0 1.2 8 ) (x F X ) (x F X X F X (x) u i i x i x Fig. 2.
Generation of random number from distribution with probability-bound structure.).

Interval Finite Element Methods

One of the main features of interval arithmetic is its capability of providing guaranteed results. However, it has the disadvantage of overestimation if variables have multiple occurrences in the same expression. For example, if x is an interval, evaluation of the function f (x) = x -x using naive rules of interval arithmetic will not return zero but rather an interval that contains zero. Such situations lead to extremely pessimistic results, and have discouraged some researchers of pursuing further developments using interval representations.

The Finite Element Method (FEM) is a numerical method for solving differential and partial differential equations with enormous applications in different fields of the sciences and engineering. Interval Finite Element Methods (IFEM) have been developed for the last 15 years to handle the analysis of systems for uncertain parameters described as intervals. Since the early development of IFEM during the mid-1990s of the last century [START_REF] Koyluoglu | Interval Algebra to Deal with Pattern Loading and Structural Uncertainty[END_REF], [START_REF] Muhanna | Development of Interval Based Methods for Fuzziness in Continuum Mechanics[END_REF], [START_REF] Rao | Fuzzy Finite Element Approach for Analysis of Imprecisely Defined Systems[END_REF], [START_REF] Nakagiri | Finite Element Interval Estimation by Convex Model[END_REF], [START_REF] Rao | Analysis of Uncertain Structural Systems using Interval Analysis[END_REF], and [START_REF] Rao | Numerical Solution of Fuzzy Linear Equations in Engineering Analysis[END_REF], researchers have focused, among other issues, on two major problems: the first is how to obtain solutions with reasonable bounds on the system response that make sense from a practical point of view, or in other words, with the least possible overestimation of their bounding intervals; the second is how to obtain reasonable bounds on the derived quantities that are functions of the system response.

The most successful approaches for overestimation reduction are those which relate the dependency of interval quantities to the physics of the problem being considered; details on these developments can be found in the works of the authors and their collaborators [START_REF] Muhanna | Development of Interval Based Methods for Fuzziness in Continuum Mechanics[END_REF], [START_REF] Muhanna | Uncertainty in Mechanics Problems-Interval-based Approach[END_REF], [START_REF] Zhang | Nondeterministic Linear Static Finite Element Analysis: An Interval Approach[END_REF], and [START_REF] Rao | A New Interval Finite Element Formulation with the Same Accuracy in Primary and Derived Variables[END_REF].

A brief description of IFEM formulation is presented below, details can be found in the authors' work [START_REF] Rao | A New Interval Finite Element Formulation with the Same Accuracy in Primary and Derived Variables[END_REF]. The two major issues in this formulation are:

1. Reducing overestimation in the bounds on the system response due to the coupling and transformation in the conventional FEM formulation as well as due to the nature of used interval linear solvers [START_REF] Muhanna | Uncertainty in Mechanics Problems-Interval-based Approach[END_REF]). 2. Obtaining the secondary (derived) variables such as forces, stresses, and strains of the conventional displacement FEM along with the primary variables (displacements) and with the same accuracy of the primary ones. Previous interval methods calculate secondary variable from interval solutions of displacement which result in a significant overestimation.

Discrete Structural Models

In steady-state analysis, the variational formulation for a discrete structural model within the context of the Finite Element Method (FEM) is given in the following form of the total potential energy functional [START_REF] Bathe | Finite Element Procedures[END_REF] and [START_REF] Gallagher | Finite Element Analysis Fundamentals[END_REF].

Π = 1 2 U T KU -U T P (3)
with the conditions ∂Π ∂U i = 0 for all i

where Π, K, U , and P are total potential energy, stiffness matrix, displacement vector, and load vector respectively. For structural problems the formulation will include both direct and indirect approaches. For the direct approach the strain is selected as a secondary variable of interest, where a constraint can be introduced as C 2 U = . For the indirect approach constraints are introduced on displacements of the form C 1 U = V in such a way that Lagrange multipliers will be equal to the internal forces. C 1 and C 2 are matrices of orders m × n and k × n, respectively, where m is the number of displacements constraints, k is the number of strains, and n is the number of displacements degrees of freedom. We note that V is a constant and is a function of U . We amend the right-hand side of Eq. ( 3) to obtain

Π * = 1 2 U T KU -U T P + λ T 1 (C 1 U -V ) + λ T 2 (C 2 U -) (5)
where λ 1 and λ 2 are vectors of Lagrange multipliers with dimensions m and k, respectively. Invoking the stationarity of Π * , that is δΠ * = 0, we obtain

    K C T 1 C T 2 0 C 1 0 0 0 C 2 0 0 -I 0 0 -I 0         U λ 1 λ 2     =     P V 0 0     (6)
The solution of Eq. ( 6) will provide the values of the dependent variable U and the derived ones λ 1 , λ 2 , and at the same time.

The present interval formulation, which will be introduced in the next section, is an extension of the Element-By-Element (EBE) finite element technique developed in the work of the authors [START_REF] Muhanna | Uncertainty in Mechanics Problems-Interval-based Approach[END_REF].

The main sources of overestimation in the formulation of IFEM are the multiple occurrences of the same interval variable (dependency problem), the width of interval quantities, the problem size, and the problem complexity, in addition to the nature of the used interval solver of the interval linear system of equations. While the present formulation is valid for the FEM models in solid and structural mechanics problems, the truss model will be used here to illustrate the applicability and efficiency of the present formulation without any loss of generality.

The current formulation modifies the displacement constraints used in the previous EBE formulation to yield the element forces as Lagrange Multipliers directly and the system strains. All interval quantities will be denoted by nonitalic boldface font. Following the procedures given in [START_REF] Rao | A New Interval Finite Element Formulation with the Same Accuracy in Primary and Derived Variables[END_REF] we obtain the interval linear system

    K C T 1 C T 2 0 C 1 0 0 0 C 2 0 0 -I 0 0 -I 0         U λ 1 λ 2     =     P 0 0 0     (7)
where K is an interval matrix of dimension (dof × dof), where dof is the sum of element degrees of freedom and the free node degrees of freedom. It consists of the individual elements' local stiffness and zeros at the bottom corresponding to the free nodes' degrees of freedom.

The accuracy of the system solution depends mainly on the structure of Eq. ( 7) and on the nature of the used solver. The solution of the interval system [START_REF] Shafer | A Mathematical Theory of Evidence[END_REF] provides the enclosures of the values of dependent variables which are the interval displacements U, interval element forces λ 1 , the multiplier λ 2 , and the element's interval strains . An iterative solver is discussed in the next section.

Iterative Enclosures

The best known method for obtaining very sharp enclosures of interval linear system of equations that have the structure introduced in Eq. ( 7) is the iterative method developed in the work of Neumaier and Pownuk [START_REF] Neumaier | Linear Systems with Large Uncertainties, with Applications to Truss Structures[END_REF]. The current formulation results in the interval linear system of equations given in [START_REF] Shafer | A Mathematical Theory of Evidence[END_REF] which can be transformed to have the following general form:

(K + BDA)u = a + F b ( 8
)
where D is diagonal. Furthermore, defining

C := (K + BD 0 A) -1 (9)
where D 0 is chosen to ensure invertablility (often D 0 is selected as the midpoint of D), the solution u can be written as

u = (Ca) + (CF )b + (CB)d. ( 10
)
To obtain a solution with tight interval enclosure we define two auxiliary interval quantities,

v = Au (11) d = (D 0 -D)v
which, given an initial estimate for u, we iterate as follows:

v = {(ACa) + (ACF )b + (ACB)d} ∩ v, d = {(D 0 -D)v} ∩ d ( 12
)
until the enclosures converge, from which the desired solution u can be straightforwardly obtained. This formulation allows obtaining the interval displacement U and the accompanied interval derived quantities λ 1 , λ 2 , and with the same accuracy. A number of examples are introduced in the next section.

Examples

Fig. 3 shows a planar truss. The reliability of this structure using single deflection criteria has been presented by Zhang, et al. [START_REF] Zhang | Interval Monte Carlo Methods for Structural Reliability[END_REF]. Two limit states are considered in this work, the serviceability and the strength. The deflection limit at the mid-span is set to 2.4 cm, and the allowable stress in any member is set 200 MPa. Interval linear elastic analyses are performed. The element stress is calculated using conventional interval methods (stresses are calculated from interval displacements) as well as the improved algorithm outlined in this paper. The cross-sectional areas for the 15 members and the three loads are identified as the basic random variables. All the 18 random variables are assumed to be mutually statistically independent. Assume that based on experience, the cross-sectional areas can be modeled by normal distributions, and the loads modeled by lognormal distributions. Suppose the statistics for the random variables were estimated from limited samples of data. Table 1 gives the available sample statistics for the cross-sectional areas and the logarithm of the loads (Ln P). The Youngs modulus is assumed deterministic (200 GPa).

Table 1. Sample statistics for the basic random variables (truss in Fig. 3).

Variables

Sample From the sample size, one can calculate confidence bounds on the mean of the random variables. We will use these bounds to construct the bounding functions defining the probability bounds for the random variables. Two cases are considered (1) the uncertain means for the (logarithm of) loads only and (2) uncertain cross-sectional area as well as loading. The interval Monte Carlo method is used to obtain probability bounds information on the variables associated with the limit states as well as interval estimates for the failure probability. All results are calculated using 10,000 realizations.

Fig. 4 presents the probability bounds for the central deflection of the truss subject to uncertain loading. Using a limit state for deflection of 2.4 cm, the bounds on the probability of survival is given by the intersection of a vertical line at 2.4 cm and the probability bounds. The failure probability bounds are then calculated from the survival bounds and are [1.34, 1.93]%. Fig. 5 presents the probability bounds results for the maximum absolute value of the stress the structure. Two different bounds are presented; the wider bounds are calculated using standard interval methods while the narrower bounds result from the current method. The calculated failure probabilities for the strength-based limit state [0.01, using the current method and [< 0.01, 9.95]% using the previous method. While both results bound the failure probability, the previous method significantly overestimates the bounding values.

Fig. 6 presents the probability bounds for the central deflection of the truss subject to uncertain load and element cross-sectional areas. Using the same 2.4 cm limit state, the bounds on the probability of failure increase to [0.11, 5.39]%. It should be expected that increasing the uncertainty associated with the analysis will increase the width of the bounds on failure probability. Fig. 7 the structure. Again, two different bounds are presented: the new method of this work and previous interval method. Using the same limit state, the failure probabilities bound are: [0.04, 3.71]% using the new method and [< 0.01, 48.88]% using the previous method. The overestimation using previous method of calculations of element stress clearly renders the previous method useless in engineering design.

Conclusion

An interval Monte Carlo finite element method with improved calculation of the bounds of secondary quantities is presented. Using element stress as an example, the improved sharpness of the bounds is illustrated. This work resolves the issue of engineering design with limit states calculated from secondary quantities (stress) that existed in previous interval Monte Carlo finite element analyses.

DISCUSSION

Speaker: Rafi Muhanna

Figure 2

 2 Figure 1. Excel graph

Figure 3

 3 Figure 3 is attached separately.

Fig. 1 .

 1 Fig. 1. Probability-bounds associated with normal distribution with mean = [5.6, 6] and standard deviation = 1.).

Fig. 3 .

 3 Fig. 3. Truss structure.

Figure 4 . 4 .Figure 5 .

 445 Figure 4. Excel graph Fig. 4. Probability bounds for central deflection, uncertain load.

Fig. 5 .Fig. 6 .

 56 Fig. 5. Probability bounds of maximum absolute stress in the structure, uncertain load.

Fig. 7 .

 7 Figure 7. Excel graphFig. 7. Probability bounds of maximum absolute stress in the structure, uncertain load and element cross-sectional area.

Table 2 .

 2 Probability bounds for the basic random variables (truss in Fig.3).

			mean Sample standard deviation No. of samples
	A1 -A6 (cm 2 )	10.32	0.516	30
	A7 -A15 (cm 2 ) 6.45	0.323	30
	Ln P1	3.2122	0.071474	20
	Ln P2	3.9982	0.071474	20
	Ln P3	3.2122	0.071474	20

  Richard Hanson : Does the use of interval arithmetic in finte element methods scale to the size of the structure show in your slide (large!)?Rafi Muhanna : Yes. We have tested that and the results are reported in following paper: Muhanna, R., Mullen, R., and Zhang, H., "Interval Finite Element as a Basis for Generalized Models of Uncertainty in Engineering Mechanics," Reliable Computing Journal, Springer Netherlands, Vol. 13, No. 2, pp. 173-194, April 2007, and it has been found that the computer time scaled cubically to the number of interval variables which is similar to the conventional finite elements. This talk illustrates that the use of interval arithmetic to get bounds not grossly excessive requires close analysis and exploitation of a problem's properties. Here the ideas that work on static analysis by finite elements of elastic structures also works on elliptic partial differential equations. But the methods presented here do not always succeed on dynamical systems, nor on solutions of nonlinear systems of equations with interval coefficients. The right person to ask about interval arithmetic's successes is Professor Ulrich Kulisch in attendence here.Ulrich Kulisch : Yes, for dynamical systems double precision interval and floating-point arithmetic are not very successful. Here long interval arithmetic is the appropriate tool. For the logistic equation x n+1 := 3.75x n (1-x n ), n ≥ 0, and the initial value x 0 = 0.5 double precision floating-point and interval arithmetic totally fail (no correct digit) after 30 iterations. Long interval arithmetic still computes correct digits of a guaranteed enclosure after 2790 iterations.Rafi Muhanna : Note also that nonlinear structural problems have been successfuly solved and reported in: ICASP'11: Applications of Statistics and Probability in Civil Engineering, Faber, Köhler and Nishijima (eds.), Taylor & Francis Group, London (2011), with the title: "Interval finite elements for nonlinear material problems." Linear dynamic problems are addressed in: Modares, M., Mullen, R., L. and Muhanna, R. L. "Natural Frequency of a Structure with Bounded Uncertainty," Journal of Engineering Mechanics, ASCE, Vol. 132, No. 12, pp. 1363-1371, 2006.Ronald Boisvert : Have you written your own interval arithmetic software infrastructure, or are you using tools developed elsewhere? Rafi Muhanna : The results in this presentation are calculated using MATLAB with the INTLAB toolbox.Jeffrey Fong : Log normal or normal distributions are not suitable for the modeling of ultimate tensile strength data because they do not predict a minimum strength. A 3-parameter Weibull distribution with a positive location parameter is a better model for predicting a minimum strength. Do you have experience using Weibull in interval finite elements? Rafi Muhanna : In this work we did not model the material strength. The authors have studied probability bounds using a 3-parameter Weibull distribution. The results are not yet published.

	William Kahan :