A density-matching approach for optimization under uncertainty

Modern computers enable methods for design optimization that account for uncertainty in the system---so-called optimization under uncertainty. We propose a metric for OUU that measures the distance between a designer-specified probability density function of the system response the target and system response's density function at a given design. We study an OUU formulation that minimizes this distance metric over all designs. We discretize the objective function with numerical quadrature and approximate the response density function with a Gaussian kernel density estimate. We offer heuristics for addressing issues that arise in this formulation, and we apply the approach to a CFD-based airfoil shape optimization problem. We qualitatively compare the density-matching approach to a multi-objective robust design optimization to gain insight into the method.


Introduction
Optimization algorithms are critical components of modern computational design strategies in a broad range of applications. In this iterative process, design variables that characterize the design-space are systematically varied to minimize a set of objectives. Given the availability and power of today's computers, simulations that model the design and its physical environment are routinely employed within design optimization routines in industry. This trend is ubiquitous in the aircraft, engine, automotive and shipping industries Abbas-Bayoumi and Becker (2011); Shahpar et al. (2014); Chandra et al. (2011);Papanikolaou (2010).
One common assumption in design optimization is that the physical environment in which a device operates is static and invariant to uncertainties. In other words, uncertainties in the design's physical environment, uncertainties in its manufacturing or even uncertainties inherent within its underlying simulations are not modeled.
A more complete perspective on design optimization takes these uncertainties into account. This necessitates the use of statistical metrics within the design optimization process. The goal is not only to attain ameliorated objective values but also have designs that are robust or de-sensitized to uncertainties. This philosophy has given rise to a plethora of approaches all under the umbrella of optimization under uncertainty.

Optimization under uncertainty
Broadly stated, optimization under uncertainty aids design decisions when the engineered system or its surroundings are subject to uncertainties. These problems can be classified either as robust design optimization (RDO) or reliability-based design optimization (RBDO) problems.

Robust design optimization
Robust design methodologies seek to optimize the mean performance metric while minimizing its variance Messac and Ismail-Yahaya (2002). Typically RDO lends itself to a multi-objective optimization strategy where there is a trade-off between these two objectives. This approach yields a Pareto set of feasible design solutions from which the designer will select a single design that offers a reasonable trade-off between the two objectives.
By and large, multi-objective genetic algorithms are the optimizer of choice for these problems. An alternative single-objective strategy for obtaining the Pareto front is scalarization, which minimizes a weighted sum of the two objectives. Typically multiple instances of the scalarization are required, each with a different set of weights, to obtain points on the Pareto front Allen and Maute (2004). One of the weaknesses of this approach is its inability to obtain parts of the Pareto set that are non-convex Das and Dennis (1997). Another major shortcoming is the difficulty associated in an a priori selection of the weights depending on the shape of the Pareto front Rangavajhala and Mahadevan (2011). An alternative strategy presented in Li et al. (2002) is to employ multi-point optimization with sensitivity analysis, where a common descent direction is sought for all uncertainties. However, this method presents difficulties as the number of uncertainties increases Allen and Maute (2004).
Of all the strategies for tackling RDO problems, the multi-objective optimization strategy has received the most use. Applications range from the design of Formula One brake ducts Axerio-Cilles (2012), compressor blades Seshadri et al. (2014b) to entire compression systems Ghisu et al. (2011), airfoils Tachikawa et al. (2012) and civil engineering structures Doltsinis and Kang (2004). Most of the examples cited use low-fidelity computational simulations, where individual design simulations are obtained within minutes. This approach may not be viable for industrial design cases, where complex multi-physics simulations may take several hours or days.

Reliability-based design optimization
RBDO methods seek to minimize a cost function for an engineered system, while making allowances for a specific risk and target reliability under various sources of uncertainty Frangopol and Maute (2003). This involves the evaluation of reliability metrics that need to be computed within the optimization framework .
RBDO can be formulated as a single-or multi-objective problem with reliability constraints that require tail statistics. Practical applications of RBDO include the optimization of a transonic compressor Lian and Kim (2005), an airfoil with non-linear aeroelasticity Missoum et al. (2010), structures Allen and Maute (2004) and vehicle crash worthiness .

Uncertainty quantification
The main challenge in RBDO is the evaluation of the probabilistic constraints, i.e., the tail statistics. For RDO on the other hand, the key lies in accurate computation of the first few statistical moments. We outline some of the techniques used to obtain these metrics below.

Tail statistics
One goal of RBDO is to determine the probability that an engineered system will fail during service. This behavior is expressed as the multivariate integral of the joint pdf for the input random variables. It is the area under the pdf curve, below or above a certain threshold value (known as the limit state). As the joint pdf does not typically have an explicit analytical expression, computing the integral is not trivial; however, numerous approximations exist. These include most-probable point (MPP)-based methods with first-and second-order approximations (FORM and SORM), the mean value method (also called FOSM) and sampling methods, such as direct Monte Carlo and adaptive importance sampling Bichon et al. (2008).

Moment computations
Moment methods can also be used to obtain the mean and variance of a function through a Taylor series expansion about its nominal value. These approximations are exact only when the uncertainties are Gaussian and the objective function can be expressed as a Taylor series Ghisu et al. (2011). For cases where these assumptions do not hold, direct sampling methods, such as Monte Carlo and Latin hypercube, can be used. Statistical accuracy in the mean and variance will typically require far fewer samples than an estimation of the tail probabilities. If the number of uncertain variables is sufficiently small and the response sufficiently smooth, then numerical quadrature can be used in place sampling methods for higher accuracy.

Aggressive design
The RBO and RBDO methods emphasize the roles of low-order moments and tail probabilities, respectively, in the design optimization process when uncertainty is present in operating conditions. In this paper, we present an alternative single-objective RDO / RDBO strategy in which no moment definitions and no scalarization is required. The approach takes into account the full probability density function of the response with respect to the uncertainties.
We consider the situation where a designer has expressed the desired performance of the system as a pdf, which we call the target pdf. We then seek a point in the design space whose response pdf matches the given target as closely as possible. Using kernel density estimates of the pdf, we formulate this pdf matching as a continuous optimization problem suited for existing optimization software. We refer to this approach to design under uncertainty as aggressive design. The continuous formulation is an improvement over previous work Seshadri et al. (2014a), where we discretized the pdfs with histograms; the histogram approach yields a more difficult optimization problem.
Section 2 details the aggressive design framework and poses the optimization problem. In Section 3 we elaborate on essential heuristics within our methodology. We demonstrate the efficacy of aggressive design on two numerical examples in Section 4.

Aggressive design
Simply stated, aggressive design is pdf matching. It seeks to minimize the distance between the pdf of a quantity of interest under uncertainty and a given target. The objective of aggressive design is to find the design whose pdf most closely resembles the target pdf.

Formulation
Consider a function f = f (s, ω) that represents the response of a physical model with design variables s ∈ S ⊆ R n and uncertain variables ω ∈ Ω. The space S denotes the constraints on the design variables. For simplicity, we assume that f is scalar-valued, f ∈ F ⊆ R, though this can be generalized. We also assume that f is continuous in both s and ω. The uncertain inputs ω are a vector of random variables (or in some cases a single variable) defined on a complete probability space (Ω, Σ, P ) with a pdf p = p(ω) corresponding to the measure P . For a fixed s ∈ S, let q s : F → R + be a probability density function of f (s, ω). The shape of q s will be different for different values of s.
We assume that we are given a pdf that describes the desired behavior of the model, including random variations. This pdf is the target pdf, and it is provided by the designer. For example, we may be given the mean and the standard deviation of a normal distribution, which respectively describe the desired average response and the permissible fluctuations. Denote the target pdf by t : R → R + .
To find the values of the design variables s that bring the model's behavior as close as possible to the designer's target, we pose the following optimization problem: where d(·, ·) is a distance metric between two comparable probability density functions. The values s * correspond to the optimal design under uncertainty. A few comments on this optimization problem are in order. By virtue of being a distance metric, d ≥ 0. However, d(t, q s ) is not generally a convex function of s. Therefore, s * may not be unique, and the optimization problem may need a regularization term to make it well-posed (e.g., Tihkonov regularization).
The minimum value of the objective function d(t, q s * ) measures how well the optimal design meets the designer's specifications. A non-zero value at the minimum means that the model can be improved, e.g., by incorporating more controls or otherwise modifying the relationship between the design variables and the system behavior. If the minimum is deemed to be too far from the ideal, then the designer may request a radical redesign of the system that allows the model to get closer to their specifications.

Computational methods
Next we turn to the computational aspects of solving the optimization problem in Eq. (1). There are several distance metrics for probability density functions. To use scalable, gradient-based optimization methods, we would like to choose a distance metric that ensures that the objective is a differentiable function of s. The most familiar is the (squared) L 2 -norm, To use this norm, we assume that t and q s are square-integrable. Moreover, let us assume that f (s, ω) is bounded for all s and ω, so that Since t is independent of s, we can ignore the last two terms in the optimization. We choose an N-point numerical integration rule on the interval [f ℓ , f u ] with pointsf i and weights w i with i = 1, . . . , N. The discretized objective function becomes For a fixed s, the density q s is, in general, not a known function of f and must be estimated. To ensure that the objective function remains differentiable with respect to s, we can use a kernel density estimate of q s . In particular, we choose a set of M points ω j ∈ Ω and define the functions For a bandwidth parameter h and a kernel K = K h , we approximate q s by Then the vector q s can be approximated where , and e is an M-vector of ones. For computation, we replace q s byq s in the approximate objective functiond in Eq. (5).
We can directly compute the gradient of the approximate objectived with respect to the design variables s. For the kth component of s, denoted s k , Note that ∂q s where K ′ is the derivative of the kernel with respect to its argument. Define Then we can concisely write the derivative of the objectived from Eq. (5) with respect to the kth component of s as The elements of K, K ′ , and f ′ k all depend on s. Define the M × n matrix F ′ by We can write the gradient of the objective function -oriented as a row vector -as This expression puts two restrictions on the types of problems we can consider. First, we need to use a differentiable kernel K in the kernel density estimator of q s in Eq. (8). Second, we must consider functions f (s, ω) that are differentiable with respect to the parameters s. In terms of interfaces from the simulation code for the physical model f (s, ω), we need to evaluate: (i) f given s and ω, and (ii) the gradient ∇ s f given s and ω -just like a standard optimization method without uncertainty. If gradients are not available, then one may consider finite difference approximations. The rest of the operations involve evaluating the kernel K, its derivative, matrix-vector products, and matrix-matrix products.

Discussion and heuristics
In this section we discuss a few important choices in aggressive design: selecting the kernel and its bandwidth in the density estimate, using response surfaces for expensive simulations, and choosing the initial design for the aggressive design optimization.

Kernel and bandwidth parameter selection
Two key components of kernel density estimates are the underlying kernel function and its bandwidth parameter. To ensure that the kernel density estimate in Eq. (9) is differentiable, two conditions must be met. First, the underlying kernel function must be differentiable; we choose a Gaussian kernel, which is infinitely differentiable. Second, the bandwidth parameter must not be too small, as such a choice will lead to spurious oscillations in the kernel density estimate Scott (1992).
Various techniques exist for selecting the bandwidth parameter. In principle, an optimal choice seeks to minimize the mean-squared error of the kernel density estimate and the true distribution. For a Gaussian kernel approximating a Gaussian pdf, this yields: A robust estimate of σ can be made using the formula whereμ denotes the median of the sample and where x i is the i th sample. Each iteration of the optimizer computes a new bandwidth parameter for the kernel density estimate. The above formulation is used in the numerical examples in this paper. It is also the default option in MATLAB's ksdensity function. Two other techniques for bandwidth parameter selection are crossvalidation and plug-in bandwidths Bowman and Azzalini (1997). Both pursue different approaches of minimizing the mean-squared error. By opting for the Gaussian-optimal bandwidth parameter criterion for all distributions, we will encounter some errors in our estimation of the moments of other distribution families. We examine this in some simple cases; Table 1 compares the true mean µ, variance σ 2 and skewness γ for a series of distributions with those obtained from a kernel density estimate of the same distribution. Here kernel density estimates are indicated by a tilde. We used 10 7 random samples of uniform, Gaussian and beta distributions. The kernel density estimate is computed using a Gaussian kernel function with 2000 grid points. For the Gaussian distribution (see the first row of Table 1), this was over the support [−5, 5], while for the other distributions this was over [0, 1]. Bandwidth parameter h values are also shown. As expected, nearly exact results are obtained for the Gaussian distribution. Mean and variance for the other distributions are generally well captured. Skewness values for the uniform and beta distributions are, however, marginally over-predicted. Histograms of the four distributions do reveal that the pdfs are generally well represented (see Fig. 1), indicating that the apparent discrepancies in skewness manifest as a relatively small error.

Optimization Strategy
To leverage the smoothness of the objective function and the ability to compute gradients, the sequential quadratic program (SQP) optimizer is used in this work. At every iteration, SQP chooses a search direction based on the constrained optimization of a second-order least-squares response surface (the quadratic program). This quadratic subproblem is solved using an appropriate projection method. SQP computes the gradients using forward finite differences and approximates the Hessian using the Broyden-Fletcher-Goldfarb-Shanno (BFGS) method to determine a feasible search direction. It then performs a line search in this direction to find a point that satisfies first and second order optimality conditions. At each major iteration the response surface is updated with the BFGS update rule. Convergence is declared when improvement of the objective, or norms of the gradients and Hessian fall below a tolerance.

Aggressive design with multiple objectives
One advantage of kernel density estimates is that they extend to more than one dimension when there are multiple quantities of interest. Multivariate kernel density estimates are products of their individual univariate kernel density estimates, where the bandwidth parameter for the K d needs to be defined in each dimension. For example, one may be interested in finding a design whose joint pdf of two responses of interest matches a given bivariate target. In this case, we can re-cast Eq. (4) as where tensor product extensions of a univariate quadrature rule can be used to approximate the integral.

Surrogate response surfaces
When simulations are expensive and the number of uncertain variables is sufficiently small, it may be more efficient to fit a response surface to a few well chosen runs for use as a cheaper surrogate when approximating the pdfs. Response surfaces for approximating pdfs are ubiquitous in uncertainty quantification, including well-known polynomial approximations (e.g., polynomial chaos). In the example in Section 4.2, we employ a polynomial surrogate in place of the relatively expensive computational fluid dynamics (CFD) solver.
In the context of polynomial surrogates, the vector q s is obtained by sampling a polynomial that interpolates a few samples of q s (f (s, w j )) for particular values of w j . The large set of samples obtained from the surrogate is then used to generate a kernel density estimate of the pdf.
Surrogates can also be used to approximate the gradients of the objective function of with respect to the design variables. In particular, if the gradient is a smooth function of the design variables, then the columns of F ′ in Eq. (15) can be interpolated. This involves obtaining n polynomial interpolants -one for each design parameter sensitivity -and then sampling each polynomial surrogate.

Choosing an initial design
The objective function (Eq. (2)) is not necessarily a convex function of the design variables, so we expect the optimizer's performance to depend strongly on the initial design point. To choose an initial design s 0 such that the corresponding pdf and the target are close, we use the following heuristic: is the mean of the target pdf and is the mean of the nominal design pdf. A similar idea is common in Bayesian inverse problems Martin et al. (2012), where the maximum a posteriori point of the Bayesian posterior is first discovered through a deterministic optimization and then used as the starting point for Markov Chain Monte Carlo.

Numerical experiments
In the following numerical experiments, we apply aggressive design to a simple linear model and an airfoil design problem.

Linear model
Consider a simple linear model given by f (s, ω) = sω + 3.5.
Here ω is a Gaussian random variable representing uncertainty. The design parameter in the above model is s. The objective of this problem is to determine the optimal value of s that minimizes the distance between the pdf of f (s, w) and a given target pdf For a fixed value of s, the output of a linear model with a Gaussian input is also a Gaussian. Thus the pdf of f (s, ω) can be written as a Gaussian distribution with mean 3.5 and standard deviation s: 13 The optimization problem seeks to minimize the distance between q s (f ) and t(f ), i.e.,  Table 2 shows the effect of varying the number of quadrature points used in the integral computation for this problem. Above 1000 points the error in the choice of s is in the fourth decimal place. A nearly identical result for the equivalent number of quadrature points is also observed in Table 3, where the analytical definitions in Eq. (26) are replaced with their kernel density estimates. These values give us some confidence in the use of kernel densities for optimization. It should be noted that for this case, if ω is uniform, then the target is matched exactly.

Computational airfoil design
In this example, we apply aggressive design to the computational design of an airfoil under uncertainty. The airfoil used is a NACA0012 at a Reynolds number of 10 6 and an angle of attack of 5 • . The uncertainty is in the inlet  Mach number which is characterized by a β(2, 2) distribution between Mach numbers of 0.66 and 0.69. Flow computations for this airfoil are carried out by solving the compressible Euler equations using Stanford University's SU 2 flow solver Aerospace Design Lab (2011). For the design problem, the airfoil is parameterized with a total of 16 Hicks-Henne bump functions: 8 on the upper surface and 8 on the lower surface. The design space is the amplitude of each bump, which is blended smoothly with the rest of the airfoil. The amplitude ranges and locations for these bumps are shown in Table 4.
Once a vector of these design parameters is selected, the airfoil mesh is deformed using a torsional spring analogy to match the new airfoil's geometry. This is done using a mesh deformation code, also available in SU 2 . Next, the flow solver is run on the new mesh, where, upon convergence, the lift-to-drag ratio (L/D) values are stored. For this problem our design objective includes pdfs of L/D.
To determine the pdf under the chosen Mach number uncertainty, we build a polynomial interpolant as in stochastic collocation Xiu (2010) of L/D as a function of Mach number. Flow computations are carried out at Mach numbers of 0.66, 0.67, 0.68 and 0.69. The resulting L/D values are then interpolated via a Lagrange polynomial. We then sample from the input beta distribution, and for each sample we evaluate the polynomial response surface as a surrogate for L/D. This produces a set of 10 7 samples of L/D that we use to approximate the density function. Mean, variance and higher moments can also be estimated from these samples. For speed-up, the four CFD simulations at different Mach numbers on the same geometry were run in parallel.

Robust design optimization
We begin our investigations into the computational airfoil design problem by carrying out a multi-objective RDO. The aim of this study is two fold; first to obtain the set of Pareto optimal designs, and second to determine the cost of RDO for design under uncertainty. The optimization objectives are to minimize the inverse of the mean L/D ratio while minimizing the variance in L/D. Using integral notation, we write these as where f is the L/D ratio, ω represents the Mach number uncertainty, and s the vector of Hicks-Henne bump function amplitudes that parameterize the airfoil. The optimization algorithm chosen for the study is NSGA- II Deb et al. (2002); with the default algorithm parameters shown in Table 5. A population size of 100 is selected with 35 generations, yielding a total of 3500 function calls made by the optimizer. Each function call in turn requires 4 CFD computations (as described earlier). Thus for the RDO study, a total of 3500 × 4 = 14, 000 CFD computations are carried out. Fig. 3 plots the results of this optimization (all the solutions evaluated in the optimization). For clarity the mean is plotted on the x-axis and  We now highlight some of the main observations of this study. RDO with its emphasis on low-order moments, cannot control potentially important characteristics of the optimal design pdf -such as skewness and tail probabilities. For the chosen design space, designs with similar mean and variance values happen to yield similar skewness values. However, this may not always be the case. In contrast to RDO, aggressive design attempts to control the entire pdf -instead of the first two moments.
Another key observation is the cost; the entire optimization took approximately a day on a desktop machine. As most GAs do not use explicit gradient information they cannot leverage adjoint methods for RDO. This is puts these methods at a significant disadvantage especially for problems with a large number of design variables, and where relatively cheap adjoints are available. Finally, once a Pareto set of optimal solutions is obtained (as in Fig. 3), the designer is still tasked with inspecting the pdfs of designs with suitable mean and variance values. In other words, a more careful inspection of the actual pdf is necessary -especially in cases where a certain tail probability must be satisfied (e.g. RDBO). Now we investigate the same problem with aggressive design.

Aggressive design
Aggressive design approaches the optimization under uncertainty problem by incorporating the designer's ideal design pdf. In this section, four different target distributions are selected. For each target, two optimizations are carried out: an initial design selection optimization, and an aggressive design optimization to yield the aggressive design. As mentioned earlier, the initial design selection is carried out by solving the following: where the value µ T arget is the mean of the target pdf, and c is a scaling constant for numerical stability. Here, flow computations for each airfoil design are carried out at a Mach number of 0.66, which corresponds to the nominal Mach number. The initial design selected for this problem was the NACA0012 airfoil. The optimization was carried out using MATLAB's sequential quadratic programming (SQP) implementation within the fmincon function. The optimizations were terminated when a tolerance in the design parameters of 10 −5 was met.
The optimized design from the initial design optimization was used as the initial design for the aggressive design optimization. The rationale for this two-step approach was discussed previously. Recall the discretized objective function in aggressive design is where theq s denotes the kernel density estimate of the design pdf evaluated at the quadrature points, t is the target pdf evaluated at the quadrature points, and c is a scaling constant.
Here we use N = 2000 quadrature points on the interval [0, 90]. In general, the value of N should be high enough to ensure integration accuracy and visibly smooth kernel density estimates. We also tested higher values of N, but they did not yield any visible improvement in the optimized design. A trapezoidal quadrature rule is used, though tests were also carried out using Clenshaw-Curtis and Gauss-Legendre rules with no significant difference in the output. Kernel density estimates were computed using samples of the polynomial response surface. A Gaussian kernel function was used with the Gaussian-optimal bandwidth parameter. The optimization was carried out using SQP.
The gradient of the objective function was given in Eq. (15) as a product of the kernel K, its derivative K ′ and design parameter sensitivity matrix F ′ . The steps followed to compute entries of the F ′ matrix for the airfoil design problem are elaborated upon below.
For a given airfoil design, two adjoint computations for the surface sensitivity of the lift and drag coefficients, L and D, are carried out. A gradient projection algorithm is used to project the sensitivities on the perturbed mesh to determine the sensitivity with respect to the design parameters. Once the 16 gradient values of L and D were obtained, a quotient rule was applied to obtain the sensitivity of L/D. These computations were carried out at each of the four Mach numbers. The resulting values were then interpolated using a third-order Lagrange polynomial. Polynomial surrogates for all 16 adjoint sensitivities are shown in Fig. 4 for the NACA0012. From these plots we conclude that a third-order polynomial surrogate is sufficient, and it can be sampled to produce F ′ , which has the form CFD evaluations are required; and to obtain ∇ sd an additional 8 CFD calls are needed (4 for the lift-adjoint and 4 for the drag-adjoint).
We now examine the computational cost, optimized solution and optimization trajectory for four different target distributions; see Table 6. The first two targets; A and B, are selected from the RDO Pareto front and are thus feasible -i.e. there exists a design vector which will match the target exactly. Targets C and D on the other hand are selected to lie just outside the Pareto set of optimal solutions and are thus infeasible -i.e. there exists no design vector that will match the target exactly.
For each case, an initial design optimization is carried out followed by the aggressive design optimization. Specific comments for the individual cases are as follows:  Fig. 3, and in principle would be a worthy robust design candidate. The location of the target and initial design with respect to the RDO set of solutions is shown in Fig. 5(a). The initial design optimization required a total of 172 CFD function evaluations with a wall clock time of 1871 seconds. The convergence history for this optimization is shown in Fig. 5(b). This optimization yielded an initial design for the aggressive design problem with a mean of 44.135 and a variance of 12.427. Fig. 5(a) plots the optimization trajectory for the aggressive design problem. This optimization took 12 SQP iterations with 91 calls made by the optimizer (see Fig. 5(c-d)) with a wall clock time of 7628 seconds. Thus, combining wall clock times for both the initial design and aggressive design optimizations we have a total of 9499 seconds with 1315 CFD evaluations -roughly an order of magnitude less than required by RDO. The obtained aggressive design nearly matches the mean, variance and skewness of the target, as shown in Fig. 5(e). Contours of the C p profiles over the airfoil chord further confirm the similarity in designs. The minor differences observed may be attributed to the errors arising from adjoint sensitivities and response surface approximations. 2. Case B: This target was selected because of its low variance. In the same spirit as Case A, plots for Case B are shown in Fig. 8. The initial design optimization took 882 seconds with 4 SQP iterations and 105 CFD evaluations. The aggressive design optimization took 7244 seconds with 7 SQP iterations and 44 function calls. The aggressive design yielded an airfoil with mean, variance and skewness of 29.224, 1.865 and 0.380 respectively. These values are extremely close to the target values of 29.126, 1.437 and 0.341 respectively. 3. Case C and D: These cases represent two targets that have the same mean and variance values with opposing skewnesses: C has a positive skewness of 1.05 while D has a negative skewness of the same magnitude. The results for these two studies are shown in Figure 7. As both designs have extremely close mean values the initial design obtained form the initial design optimization routine was used for both aggressive design optimizations. There are some interesting points to note in the aggressive design results. First, both designs are Pareto optimal -i.e. they lie on the Pareto front as shown in Figure 7(a,b). Second, the aggressive designs for case C and D are different -i.e. the different skewness incorporated in the target affects the outcome of aggressive design optimization. Third, is the difference in cost. While case C took 2262 seconds (6 iterations with 24 function calls), case D took thrice as long at 7780 seconds (8 iterations with 83 function calls). We infer that this difference in cost is attributed to the fact that for the positively skewed target, the optimal value of the normalized distance is 0.1312; comapred to 0.029 for the negatively skewed target.
In general we find aggressive design to yield solutions that either match or closely approximate the selected target. This matching is done not solely on the mean and variance but on higher-order moments as well. Furthermore, the computational cost of aggressive design is on average an order of magnitude less than RDO.

Conclusions
We present a new methodology for design under uncertainty: aggressive design. Aggressive design seeks to minimize the distance between the pdf of the response and a target pdf representing a designer's ideal performance under uncertainty. In aggressive design, the distance between the two pdfs is expressed as a smooth and differentiable objective function, which enables the use of gradient-based optimization procedures.
We have applied this novel approach to a simple linear model to illustrate the concept, and to the computational design of an airfoil to demonstrate its effectiveness.By applying aggressive design to practical design problems, we are advocating a new design under uncertainty philosophy that accounts for the entire pdf of the design instead of only low-order moments or tail probabilities.