Structural Reliability under Uncertainty in Moments:
Distributionally-robust Reliability-based Design Optimization

Yoshihiro Kanno 222Mathematics and Informatics Center, The University of
Tokyo, Hongo 7-3-1, Tokyo 113-8656, Japan. E-mail<EMAIL_ADDRESS>
###### Abstract

This paper considers structural optimization under a reliability constraint,
where the input distribution is only partially known. Specifically, when we
only know that the expected value vector and the variance-covariance matrix of
the input distribution belong to a given convex set, we require that, for any
realization of the input distribution, the failure probability of a structure
should be no greater than a specified target value. We show that this
distributionally-robust reliability constraint can be reduced equivalently to
deterministic constraints. By using this reduction, we can treat a
reliability-based design optimization problem under the distributionally-
robust reliability constraint within the framework of deterministic
optimization, specifically, nonlinear semidefinite programming. Two numerical
examples are solved to show relation between the optimal value and either the
target reliability or the uncertainty magnitude.

> Keywords
>
> Reliability-based design optimization; uncertain input distribution; worst-
> case reliability; robust optimization; semidefinite programming; duality.

##  1 Introduction

Reliability-based design optimization (RBDO) is a crucial tool for structural
design in the presence of uncertainty [45, 2, 49, 36]. It adopts a
probabilistic model of uncertainty, and evaluates the probability that a
structural design satisfies (or, equivalently, fails to satisfy) performance
requirements. An underlying premise is that complete knowledge on statistical
information of uncertain parameters is available. In practice, however, it is
often difficult to obtain statistical information with sufficient accuracy.
This incents recent intensive study of RBDO with incomplete statistical
information [9, 10, 16, 21, 22, 24, 25, 35, 37, 38, 46, 50, 51, 53].

Another methodology dealing with uncertainty in structural design is robust
design optimization [6, 20, 31]. Although there exist several different
concepts in robust design optimization, in this paper we focus attention on
the worst-case optimization methodology, which is called robust optimization
in mathematical optimization community [4]. This methodology adopts a
possibilistic model of uncertainty, i.e., specifies the set of possible values
that the uncertain parameters can take. We call this set an uncertainty set.
Then, the objective value in the worst case is optimized, under the condition
that the constraints are satisfied in the worst cases.

This paper deals with RBDO when the input distribution is only partially
known. Specifically, we assume that the true expected value vector and the
true variance-covariance matrix are unknown (i.e., the true values of the
first two moments of the input distribution are unknown), but they are known
to belong to a given closed convex set. For example, suppose that the input
distribution is a normal distribution, and we only know that each component of
the expected value vector and the variance-covariance matrix belongs to a
given closed interval. Then, for each possible realization of pair of the
expected value vector and the variance-covariance matrix, there exists a
single corresponding normal distribution. The set of all such normal
distributions is considered as an uncertainty set of the input
distribution.111This uncertainty set is dealt with in section 3.  As another
example, suppose that distribution type of the input distribution is also
unknown. Then the uncertainty set is the set of all probability distributions,
the expected value vector and the variance-covariance matrix of which belong
to a given set.222This uncertainty set is dealt with in section 4.1.  Among
probability distributions belonging to a specified uncertainty set defined as
above, the worst-case distribution is the one with which the failure
probability takes the maximum value. Our methodology is that we require a
structure to satisfy the reliability constraint evaluated with the worst-case
distribution. In other words, for any probability distribution belonging to
the uncertainty set, the failure probability should be no greater than a
specified target value. Thus, the methodology guarantees robustness of the
structural reliability against uncertainty in the input distribution.333More
precisely, the uncertainty here means the uncertainty in the expected value
vector and the variance-covariance matrix of the input distribution.  The
major contribution of this paper is to show, under some assumptions, that this
structural requirement is equivalently converted to a form of constraints that
can be treated in conventional deterministic optimization. As a result, a
design optimization problem under this structural requirement can be solved
with a deterministic nonlinear optimization approach.

Recently, RBDO methods with uncertainty in the input distribution have
received considerable attention, because in practice it is often that the
number of available samples of random variables is insufficient. For example,
Gunawan and Papalambros [16] and Youn and Wang [50] proposed Bayesian
approaches to compute the confidence that a structural design satisfies a
target reliability constraint, when both a finite number of samples and
probability distributions of uncertain parameters are available. Noh et al.
[37, 38] proposed Bayesian methods to adjust an input distribution model to
limited data, with a given confidence level. When intervals of input variables
are given as input information, Zaman et al. [52] and Zaman and Mahadevan [51]
use a family of Johnson distributions to represent the uncertainty. Cho et al.
[9] and Moon et al. [35] assume that the input distribution types and
parameters follow probability distributions. The failure probability is
therefore a random variable, the confidence level of a reliability constraint,
i.e., the probability that the failure probability is no greater than a target
value, is specified. To reduce computational cost of this method, Jung et al.
[25] proposed a so-called reliability approach, inspired by the performance
measure approach [33, 34]. Subsequently, to further reduce computational cost,
Wang et al. [46] proposed to use the second-order reliability method for
computation of the failure probability. Ito et al. [21] assume that each of
random variable follows a normal distribution with the mean and the variance
modeled as random variables, and show that RBDO with a confidence level can be
converted to a conventional form of RBDO by altering the target reliability
index value. Zhang et al. [53] proposed to use the distributional probability
box (the distributional p-box) [41] for RBDO with limited data of uncertain
variables. Kanno [29, 30] and Jekel and Haftka [23] proposed RBDO methods
using order statistics. These methods, based on the order statistics, do not
make any assumption on statistical information of uncertain parameters, and
use random samples of uncertain parameters directly to guarantee confidence of
the target reliability.

As reviewed above, most of existing studies on RBDO with uncertainty in the
input distribution [9, 21, 35, 25, 46] consider probabilistic models of input
distribution parameters and/or distribution types. Accordingly, a confidence
level evaluates how the satisfaction of structural reliability is reliable. In
contrast, in this paper we consider a possibilistic model of input
distribution parameters. Hence, what this approach guarantees is a level of
robustness [3] of the satisfaction of structural reliability. A possibilistic
model might be, in general, less information-sensitive, and hence useful when
reliable statistical information of input distribution parameters is
unavailable.

From another perspective referring to Schöbi and Sudret [41], the uncertainty
model treated in this paper can be viewed as follows. Uncertainty in a
structural system is often divided into aleatory uncertainty and epistemic
uncertainty [39]. Aleatory uncertainty, i.e., natural variability, is
reflected by an (uncertain) input distribution. Epistemic uncertainty, i.e.,
state-of-knowledge uncertainty, is reflected by uncertainty in the input
distribution moments. Thus, in our model, aleatory uncertainty is
probabilistic, while epistemic uncertainty is possibilistic. In other words,
state-of-knowledge uncertainty is represented as an uncertainty set of the
input distribution moments.

Throughout the paper, we assume that only design variables possess
uncertainty, and that variation of a performance requirement can be
approximated as a linear function of uncertain perturbations of the design
variables. Also, we do not consider an optimization problem with variation of
structural topology. As for an uncertainty model of moments of the input
distribution, we consider two concrete convex sets. We show that the robust
reliability constraint, i.e., constraint that the structural reliability is no
less than a specified value for any possible realizations of input
distribution moments, can be reduced to a system of nonlinear matrix
inequalities. This reduction essentially follows the idea presented by El
Ghaoui et al. [12] for computing the worst-case value-at-risk in financial
engineering.444Diverse extensions of the methodology in El Ghaoui et al. [12]
can be found in literature on so-called distributionally robust optimization
[11, 15, 47].  We can deal with nonlinear matrix inequality constraints within
the framework of nonlinear semidefinite programming (nonlinear SDP) [48]. In
this manner, we can convert an RBDO problem under uncertainty in the input
distribution moments to a deterministic optimization problem. It is worth
noting that there exist several applications of linear and nonlinear SDPs, as
well as eigenvalue optimization, to robust design optimization of structures
[5, 32, 17, 18, 43, 19, 26, 44, 28].

The paper is organized as follows. In section 2, we consider the reliability
constraint when the input distribution is precisely known, and show some
fundamental properties. Section 3 presents the main result; we consider
uncertainty in the expected value vector and the variance-covariance matrix of
the input distribution, and examine the constraint that, for all possible
realizations of the input distribution, the failure probability is no greater
than a specified value. Section 4 discusses some extensions of the obtained
result. Section 5 presents the results of numerical experiments. Section 6
presents some conclusions.

In our notation, ⊤ denotes the transpose of a vector or matrix. All vectors
are column vectors. We use $I$ to denote the identity matrix. For two matrices
$X=(X_{ij})\in\mathbb{R}^{m\times n}$ and $Y=(Y_{ij})\in\mathbb{R}^{m\times
n}$, we denote by $X\bullet Y$ the inner product of $X$ and $Y$ defined by
$X\bullet
Y=\mathop{\mathrm{tr}}\nolimits(X^{\top}Y)=\sum_{i=1}^{n}\sum_{j=1}^{n}X_{ij}Y_{ij}$.
For a vector $\boldsymbol{x}=(x_{i})\in\mathbb{R}^{n}$, the notation
$\|\boldsymbol{x}\|_{1}$, $\|\boldsymbol{x}\|_{2}$, and
$\|\boldsymbol{x}\|_{\infty}$ designate its $\ell_{1}$-, $\ell_{2}$-, and
$\ell_{\infty}$-norms, respectively, i.e.,

$\displaystyle\|\boldsymbol{x}\|_{1}$
$\displaystyle=|x_{1}|+|x_{2}|+\dots+|x_{n}|,$
$\displaystyle\|\boldsymbol{x}\|_{2}$
$\displaystyle=\sqrt{\boldsymbol{x}^{\top}\boldsymbol{x}},$
$\displaystyle\|\boldsymbol{x}\|_{\infty}$
$\displaystyle=\max\\{|x_{1}|,|x_{2}|,\dots,|x_{n}|\\}.$

For a matrix $X=(X_{ij})\in\mathbb{R}^{m\times n}$, define matrix norms
$\|X\|_{1,1}$, $\|X\|_{\mathrm{F}}$, and $\|X\|_{\infty,\infty}$ by

$\displaystyle\|X\|_{1,1}$
$\displaystyle=\sum_{i=1}^{n}\sum_{j=1}^{n}|X_{ij}|,$
$\displaystyle\|X\|_{\mathrm{F}}$ $\displaystyle=\sqrt{X\bullet X},$
$\displaystyle\|X\|_{\infty,\infty}$ $\displaystyle=\max\\{|X_{ij}|\mid
i=1,\dots,m,\ j=1,\dots,n\\}.$

Let $\mathcal{S}^{n}$ denote the set of $n\times n$ symmetric matrices. We
write $Z\succeq 0$ if $Z\in\mathcal{S}^{n}$ is positive semidefinite. Define
$\mathcal{S}_{+}^{n}$ by $\mathcal{S}_{+}^{n}=\\{Z\in\mathcal{S}^{n}\mid
Z\succeq 0\\}$. For a positive definite matrix $Z\in\mathcal{S}^{n}$, the
notation $Z^{1/2}$ designates its symmetric square root, i.e.,
$Z^{1/2}\in\mathcal{S}^{n}$ satisfying $Z^{1/2}Z^{1/2}=Z$. We use $Z^{-1/2}$
to denote the inverse matrix of $Z^{1/2}$. We use
$\mathsf{N}(\boldsymbol{\mu},\varSigma)$ to denote the multivariate normal
distribution with an expected value vector $\boldsymbol{\mu}$ and a variance-
covariance matrix $\varSigma$. For a random variable $x\in\mathbb{R}$, its
expected value and variance are denoted by $\mathrm{E}[x]$ and
$\mathrm{Var}[x]=\mathrm{E}[(x-\mathrm{E}[x])^{2}]$, respectively.

##  2 Reliability constraint with specified moments

In this section, we assume that the expected value vector and the variance-
covariance matrix of the probability distribution of the design variable
vector are precisely known. We first recall the reliability constraint, and
then derive its alternative expression that will be used in section 3 to
address uncertainty in the probability distribution.

Let $\boldsymbol{x}\in\mathbb{R}^{n}$ denote a design variable vector, where
$n$ is the number of design variables. Assume that performance requirement in
a design optimization problem is written in the form

$\displaystyle g(\boldsymbol{x})\leq 0,$ (1)

where $g:\mathbb{R}^{n}\to\mathbb{R}$ is differentiable. For simplicity,
suppose that the design optimization problem has only one constraint; the case
where more than one constraints exist will be discussed in section 4.

Assume that $\boldsymbol{x}$ is decomposed additively as

$\displaystyle\boldsymbol{x}=\tilde{\boldsymbol{x}}+\boldsymbol{\zeta},$ (2)

where $\boldsymbol{\zeta}$ is a random vector and $\tilde{\boldsymbol{x}}$ is
a constant (i.e., non-random) vector. Therefore, in a design optimization
problem considered in this paper, the decision variable to be optimized is
$\tilde{\boldsymbol{x}}$. We use $\boldsymbol{\mu}\in\mathbb{R}^{n}$ and
$\varSigma\in\mathcal{S}^{n}$ to denote the expected value vector and the
variance-covariance matrix of $\boldsymbol{\zeta}$, respectively, i.e.,

$\displaystyle\boldsymbol{\mu}$
$\displaystyle=\mathrm{E}[\boldsymbol{\zeta}],$ $\displaystyle\varSigma$
$\displaystyle=\mathrm{E}\bigl{[}(\boldsymbol{\zeta}-\mathrm{E}[\boldsymbol{\zeta}])(\boldsymbol{\zeta}-\mathrm{E}[\boldsymbol{\zeta}])^{\top}\bigr{]}.$

It is worth noting that $\varSigma$ is positive definite. Throughout the
paper, we assume that, among parameters in a structural system, only
$\boldsymbol{\zeta}$ possesses uncertainty. Also, we restrict ourselves to
optimization without change of structural topology; i.e., we do not consider
topology optimization.555In topology optimization, it would be proper to
consider the design variables of removed structural elements as non-random
variables. In this paper we do not discuss this issue.

For simplicity and clarity of discussion, we assume
$\boldsymbol{\zeta}\sim\mathsf{N}(\boldsymbol{\mu},\varSigma)$ in section 2
and section 3. In fact, the results established in these sections can be
extended to the case that the type of probability distribution is unknown; we
then require that the reliability constraint should be satisfied for any
probability distribution with moments belonging to a specified set. We defer
this case until section 4.

Since $\boldsymbol{x}$ is a random vector, $g(\boldsymbol{x})$ is a random
variable. Therefore, constraint (1) should be considered in a probabilistic
sense, which yields the reliability constraint

$\displaystyle\mathrm{P}_{\mathsf{N}(\boldsymbol{\mu},\varSigma)}\\{g(\boldsymbol{x})\leq
0\\}\geq 1-\epsilon.$ (3)

Here, $\epsilon\in]0,1]$ is the specified upper bound for the failure
probability. Let $g^{\mathrm{lin}}(\boldsymbol{x})$ denote the first-order
approximation of $g(\boldsymbol{x})$ centered at
$\boldsymbol{x}=\tilde{\boldsymbol{x}}$, i.e.,

$\displaystyle
g^{\mathrm{lin}}(\boldsymbol{x})=g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\,(\simeq
g(\boldsymbol{x})).$

Throughout the paper, we consider an approximation of constraint (3)

$\displaystyle\mathrm{P}_{\mathsf{N}(\boldsymbol{\mu},\varSigma)}\\{g^{\mathrm{lin}}(\boldsymbol{x})\leq
0\\}\geq 1-\epsilon,$ (4)

i.e.,

$\displaystyle\mathrm{P}_{\mathsf{N}(\boldsymbol{\mu},\varSigma)}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0\\}\geq 1-\epsilon.$
(5)

Therefore, the corresponding RBDO problem has the following form:

$\displaystyle\mathop{\mathrm{Minimize}}$ $\displaystyle f(\boldsymbol{x})$
(6a) $\displaystyle\mathop{\mathrm{subject~{}to}}$
$\displaystyle\boldsymbol{x}\in X,$ (6b)
$\displaystyle\mathrm{P}_{\mathsf{N}(\boldsymbol{\mu},\varSigma)}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0\\}\geq 1-\epsilon.$
(6c)

Here, $f:\mathbb{R}^{n}\to\mathbb{R}$ is the objective function,
$X\subseteq\mathbb{R}^{n}$ is a given closed set, and constraint
$\boldsymbol{x}\in X$ corresponds to, e.g., the side constraints on the design
variables.

From the basic property of the normal distribution, we can readily obtain the
following reformulation of the reliability constraint.

######  Theorem 2.1.

Define $\kappa$ by

$\displaystyle\kappa=-\varPhi^{-1}(\epsilon),$

where $\varPhi$ is the (cumulative) distribution function of the standard
normal distribution $\mathsf{N}(0,1)$. Then, $\tilde{\boldsymbol{x}}\in X$
satisfies (5) if and only if it satisfies

$\displaystyle g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}+\kappa\|\varSigma^{1/2}\nabla
g(\tilde{\boldsymbol{x}})\|_{2}\leq 0.$ (7)

###### Proof.

Since $g^{\mathrm{lin}}(\boldsymbol{x})$ follows the normal distribution, it
is standardized by

$\displaystyle
z=\frac{g^{\mathrm{lin}}(\boldsymbol{x})-\mathrm{E}[g^{\mathrm{lin}}(\boldsymbol{x})]}{\sqrt{\mathrm{Var}[g^{\mathrm{lin}}(\boldsymbol{x})]}}\sim\mathsf{N}(0,1).$

By using this relation, we can eliminate $g^{\mathrm{lin}}(\boldsymbol{x})$
from (4) (i.e., (5)) as

$\displaystyle\mathrm{P}_{\mathsf{N}(0,1)}\left\\{z\leq-\frac{\mathrm{E}[g^{\mathrm{lin}}(\boldsymbol{x})]}{\sqrt{\mathrm{Var}[g^{\mathrm{lin}}(\boldsymbol{x})]}}\right\\}\geq
1-\epsilon.$

This inequality is equivalently rewritten by using the distribution function
$\varPhi$ as

$\displaystyle-\frac{\mathrm{E}[g^{\mathrm{lin}}(\boldsymbol{x})]}{\sqrt{\mathrm{Var}[g^{\mathrm{lin}}(\boldsymbol{x})]}}\geq\varPhi^{-1}(1-\epsilon)=-\varPhi^{-1}(\epsilon).$

By direct calculations, we see that the expected value of
$g^{\mathrm{lin}}(\boldsymbol{x})$ is

$\displaystyle\mathrm{E}[g^{\mathrm{lin}}(\boldsymbol{x})]=g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}$

and the variance is

$\displaystyle\mathrm{Var}[g^{\mathrm{lin}}(\boldsymbol{x})]$
$\displaystyle=\mathrm{E}\bigl{[}(g^{\mathrm{lin}}(\boldsymbol{x})-\mathrm{E}[g^{\mathrm{lin}}(\boldsymbol{x})])^{2}\bigr{]}$
$\displaystyle=\mathrm{E}\bigl{[}\bigl{(}(g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta})-(g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu})\bigr{)}^{2}\bigr{]}$
$\displaystyle=\mathrm{E}\bigl{[}\bigl{(}\nabla
g(\tilde{\boldsymbol{x}})^{\top}(\boldsymbol{\zeta}-\boldsymbol{\mu})\bigr{)}^{2}\bigr{]}$
$\displaystyle=\mathrm{E}\bigl{[}\nabla
g(\tilde{\boldsymbol{x}})^{\top}(\boldsymbol{\zeta}-\boldsymbol{\mu})(\boldsymbol{\zeta}-\boldsymbol{\mu})^{\top}\nabla
g(\tilde{\boldsymbol{x}})\bigr{]}$ $\displaystyle=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\mathrm{E}\bigl{[}(\boldsymbol{\zeta}-\boldsymbol{\mu})(\boldsymbol{\zeta}-\boldsymbol{\mu})^{\top}\bigr{]}\nabla
g(\tilde{\boldsymbol{x}})$ $\displaystyle=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\varSigma\nabla g(\tilde{\boldsymbol{x}}),$

which concludes the proof. ∎

In section 3, we deal with the case in which $\boldsymbol{\mu}$ and
$\varSigma$ are known imprecisely. To do this, we reformulate
$\kappa\|\varSigma^{1/2}\nabla g(\tilde{\boldsymbol{x}})\|_{2}$ in (7) into a
form suitable for analysis. The following theorem is obtained in the same
manner as El Ghaoui et al. [12, Theorem 1].

######  Theorem 2.2.

For $\kappa>0$, $\varSigma\in\mathcal{S}_{+}^{n}$, and $\nabla
g(\tilde{\boldsymbol{x}})\in\mathbb{R}^{n}$, we have

$\displaystyle\kappa\|\varSigma^{1/2}\nabla
g(\tilde{\boldsymbol{x}})\|_{2}=\min_{\varLambda\in\mathcal{S}^{n},\,z\in\mathbb{R}}\left\\{\varSigma\bullet\varLambda+\kappa^{2}z\left|\begin{bmatrix}\varLambda&\nabla
g(\tilde{\boldsymbol{x}})/2\\\ \nabla g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\
\end{bmatrix}\succeq 0\right.\right\\}.$

###### Proof.

We first show that the left side of the equation can be reduced to

$\displaystyle\kappa\|\varSigma^{1/2}\nabla g(\tilde{\boldsymbol{x}})\|_{2}$
$\displaystyle=\max_{\boldsymbol{\zeta}\in\mathbb{R}^{n}}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\mid\|\varSigma^{-1/2}\boldsymbol{\zeta}\|_{2}=\kappa\\}.$
(8)

To see this, we apply the Lagrange multiplier method to the equality
constrained maximization problem on the right side of (8). Namely, the
Lagrangian $L_{1}:\mathbb{R}^{n}\times\mathbb{R}\to\mathbb{R}$ is defined by

$\displaystyle L_{1}(\boldsymbol{\zeta};\mu)=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}+\frac{\mu}{2}(\kappa^{2}-\boldsymbol{\zeta}^{\top}\varSigma^{-1}\boldsymbol{\zeta}),$

where $\mu\in\mathbb{R}$ is the Lagrange multiplier. The stationarity
condition of $L_{1}$ is

$\displaystyle\frac{\partial L_{1}}{\partial\boldsymbol{\zeta}}$
$\displaystyle=\nabla
g(\tilde{\boldsymbol{x}})-\mu\varSigma^{-1}\boldsymbol{\zeta}=\boldsymbol{0},$
$\displaystyle\frac{\partial L_{1}}{\partial\mu}$
$\displaystyle=\kappa^{2}-\boldsymbol{\zeta}^{\top}\varSigma^{-1}\boldsymbol{\zeta}=0.$

By solving this stationarity condition, we can find that

$\displaystyle\boldsymbol{\zeta}=\dfrac{1}{\mu}\varSigma\nabla
g(\tilde{\boldsymbol{x}}),\quad\mu=\frac{\|\varSigma^{1/2}\nabla
g(\tilde{\boldsymbol{x}})\|_{2}}{\kappa}$

are optimal. Hence, the optimal value is

$\displaystyle\nabla g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\frac{\kappa\varSigma\nabla
g(\tilde{\boldsymbol{x}})}{\|\varSigma^{1/2}\nabla
g(\tilde{\boldsymbol{x}})\|_{2}},$

which is reduced to the left side of (8).

Next, observe that the right side of (8) is further reduced to

$\displaystyle\max_{\boldsymbol{\zeta}\in\mathbb{R}^{n}}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\mid\|\varSigma^{-1/2}\boldsymbol{\zeta}\|_{2}=\kappa\\}$
$\displaystyle=\max_{\boldsymbol{\zeta}\in\mathbb{R}^{n}}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\mid\|\varSigma^{-1/2}\boldsymbol{\zeta}\|_{2}\leq\kappa\\}$
$\displaystyle=\max_{\boldsymbol{\zeta}\in\mathbb{R}^{n}}\left\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\left|\begin{bmatrix}\varSigma&\boldsymbol{\zeta}\\\
\boldsymbol{\zeta}^{\top}&\kappa^{2}\end{bmatrix}\succeq 0\right.\right\\}.$
(9)

Here, the last equality follows from the fact that the positive semidefinite
constraint is equivalent to the nonnegative constraint on the Schur complement
of $\varSigma$ in the corresponding matrix, i.e.,
$\kappa^{2}-\boldsymbol{\zeta}^{\top}\varSigma^{-1}\boldsymbol{\zeta}\geq 0$;
see [7, appendix A.5.5]. It is worth noting that the last expression in (9) is
an SDP problem.

Finally, we shall show that the right side of the proposition in this theorem
corresponds to the dual problem of the SDP problem in (9). Since this dual
problem is strictly feasible, the proposition follows from the strong duality
of SDP [8, section 11.3]. We can derive the dual problem of (9) as follows.
The Lagrangian is defined by

$\displaystyle
L_{2}(\boldsymbol{\zeta};\varLambda,\boldsymbol{\lambda},z)=\begin{dcases*}\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}+\begin{bmatrix}\varLambda&\boldsymbol{\lambda}\\\
\boldsymbol{\lambda}^{\top}&z\\\
\end{bmatrix}\bullet\begin{bmatrix}\varSigma&\boldsymbol{\zeta}\\\
\boldsymbol{\zeta}^{\top}&\kappa^{2}\end{bmatrix}&if
$\begin{bmatrix}\varLambda&\boldsymbol{\lambda}\\\
\boldsymbol{\lambda}^{\top}&z\\\ \end{bmatrix}\succeq 0$,\\\
+\infty&otherwise,\end{dcases*}$ (10)

where $z\in\mathbb{R}$, $\boldsymbol{\lambda}\in\mathbb{R}^{n}$, and
$\varLambda\in\mathcal{S}^{n}$ are the Lagrange multipliers. Indeed, since the
positive semidefinite cone satisfies [27, Fact 1.3.17]

$\displaystyle\inf_{S\in\mathcal{S}^{n}}\\{S\bullet T\mid S\succeq
0\\}=\begin{dcases*}0&if $T\succeq 0$,\\\ -\infty&otherwise,\end{dcases*}$
(11)

we can confirm that the SDP problem in (9) is equivalent to

$\displaystyle\max_{\boldsymbol{\zeta}}\inf_{\varLambda,\,\boldsymbol{\lambda},\,z}L_{2}(\boldsymbol{\zeta};\varLambda,\boldsymbol{\lambda},z).$

The dual problem is defined by

$\displaystyle\min_{\varLambda,\,\boldsymbol{\lambda},\,z}\sup_{\boldsymbol{\zeta}}L_{2}(\boldsymbol{\zeta};\varLambda,\boldsymbol{\lambda},z).$
(12)

Since (10) can be rewritten as

$\displaystyle
L_{2}(\boldsymbol{\zeta};\varLambda,\boldsymbol{\lambda},z)=\begin{dcases*}(\nabla
g(\tilde{\boldsymbol{x}})+2\boldsymbol{\lambda})^{\top}\boldsymbol{\zeta}+\varSigma\bullet\varLambda+\kappa^{2}z&if
$\begin{bmatrix}\varLambda&\boldsymbol{\lambda}\\\
\boldsymbol{\lambda}^{\top}&z\\\ \end{bmatrix}\succeq 0$,\\\
+\infty&otherwise,\end{dcases*}$

we have

$\displaystyle\sup_{\boldsymbol{\zeta}}L_{2}(\boldsymbol{\zeta};\varLambda,\boldsymbol{\lambda},z)=\begin{dcases*}\varSigma\bullet\varLambda+\kappa^{2}z&if
$\begin{bmatrix}\varLambda&\boldsymbol{\lambda}\\\
\boldsymbol{\lambda}^{\top}&z\\\ \end{bmatrix}\succeq 0$, $\nabla
g(\tilde{\boldsymbol{x}})+2\boldsymbol{\lambda}=\boldsymbol{0}$,\\\
+\infty&otherwise.\end{dcases*}$

Therefore, the dual problem in (12) corresponds to the right side of the
proposition of the theorem. ∎

##  3 Worst-case reliability under uncertainty in moments

In this section, we consider the case that the moments (in this paper, the
expected value vector and the variance-covariance matrix) of the design
variable vector are uncertain, or not perfectly known. Specifically, they are
only known to be in a given set, called the uncertainty set. We require that a
structure satisfies the reliability constraint for any moments in the
uncertainty set. In other words, we require that the failure probability in
the worst case is not larger than a specified value. We show that this
requirement can be converted to a form of conventional constraints in
deterministic optimization.

###  3.1 Convex uncertainty model of moments

Let $U_{\boldsymbol{\mu}}\subset\mathbb{R}^{n}$ and
$U_{\varSigma}\subset\mathcal{S}_{+}^{n}$ denote the uncertainty sets, i.e.,
the sets of all possible realizations, of $\boldsymbol{\mu}$ and $\varSigma$,
respectively. Namely, we only know that $\boldsymbol{\mu}$ and $\varSigma$
satisfy

$\displaystyle\boldsymbol{\mu}$ $\displaystyle\in U_{\boldsymbol{\mu}},$
$\displaystyle\varSigma$ $\displaystyle\in U_{\varSigma}.$

Assume that $U_{\boldsymbol{\mu}}$ and $U_{\varSigma}$ are compact convex
sets. For notational simplicity, we write $(\boldsymbol{\mu},\varSigma)\in U$
if $\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}$ and $\varSigma\in U_{\varSigma}$
hold.

Recall that we are considering the reliability constraint in (5) with a
linearly approximated constraint function. The robust counterpart of (5)
against uncertainty in $\boldsymbol{\mu}$ and $\varSigma$ is formulated as

$\displaystyle\mathrm{P}_{\mathsf{N}(\boldsymbol{\mu},\varSigma)}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0\\}\geq
1-\epsilon,\quad\forall(\boldsymbol{\mu},\varSigma)\in U,$ (13)

i.e., we require that the reliability constraint should be satisfied for any
normal distribution corresponding to possible realizations of
$\boldsymbol{\mu}$ and $\varSigma$. This requirement is equivalently rewritten
as

$\displaystyle\max_{(\boldsymbol{\mu},\varSigma)\in
U}\bigl{\\{}\mathrm{P}_{\mathsf{N}(\boldsymbol{\mu},\varSigma)}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0\\}\bigr{\\}}\geq
1-\epsilon.$ (14)

That is, the reliability constraint should be satisfied in the worst case.

The following theorem presents, with the aid of Theorem 2.1 and Theorem 2.2,
an equivalent reformulation of (14).

######  Theorem 3.1.

$\tilde{\boldsymbol{x}}\in X$ satisfies (14) if and only if there exists a
pair of $z\in\mathbb{R}$ and $\varLambda\in\mathcal{S}^{n}$ satisfying

$\displaystyle g(\tilde{\boldsymbol{x}})+\max\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}\mid\boldsymbol{\mu}\in
U_{\boldsymbol{\mu}}\\}+\max\\{\varSigma\bullet\varLambda\mid\varSigma\in
U_{\varSigma}\\}+\kappa^{2}z\leq 0,$ (15)
$\displaystyle\begin{bmatrix}\varLambda&\nabla g(\tilde{\boldsymbol{x}})/2\\\
\nabla g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\ \end{bmatrix}\succeq 0.$ (16)

###### Proof.

It follows from Theorem 2.1 that (14) is equivalent to

$\displaystyle\max_{(\boldsymbol{\mu},\varSigma)\in
U}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}+\kappa\|\varSigma^{1/2}\nabla
g(\tilde{\boldsymbol{x}})\|\\}\leq 0.$ (17)

Furthermore, application of Theorem 2.2 yields

$\displaystyle\max_{(\boldsymbol{\mu},\varSigma)\in
U}\left\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}+\min_{z,\,\varLambda}\left\\{\varSigma\bullet\varLambda+\kappa^{2}z\left|\begin{bmatrix}\varLambda&\nabla
g(\tilde{\boldsymbol{x}})/2\\\ \nabla g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\
\end{bmatrix}\succeq 0\right.\right\\}\right\\}\leq 0.$ (18)

In the expression above, we see that $U$ is compact and convex, and the
feasible set for the minimization is convex. Also, the objective function is
linear in $\boldsymbol{\mu}$ and $\varSigma$ for fixed $z$ and $\varLambda$,
and is linear in $z$ and $\varLambda$ for fixed $\boldsymbol{\mu}$ and
$\varSigma$. Therefore, the minimax theorem [8, Theorem 8.8] asserts that (18)
is equivalent to

$\displaystyle
g(\tilde{\boldsymbol{x}})+\min_{z,\,\varLambda}\max_{(\boldsymbol{\mu},\varSigma)\in
U}\left\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}+\varSigma\bullet\varLambda+\kappa^{2}z\left|\begin{bmatrix}\varLambda&\nabla
g(\tilde{\boldsymbol{x}})/2\\\ \nabla g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\
\end{bmatrix}\succeq 0\right.\right\\}\leq 0.$

This inequality holds if and only if there exists a feasible pair of
$z\in\mathbb{R}$ and $\varLambda\in\mathcal{S}^{n}$ satisfying

$\displaystyle g(\tilde{\boldsymbol{x}})+\max_{(\boldsymbol{\mu},\varSigma)\in
U}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}+\varSigma\bullet\varLambda+\kappa^{2}z\\}\leq
0,$

which concludes the proof. ∎

The conclusion of Theorem 3.1 is quite abstract in the sense that concrete
forms of $U_{\boldsymbol{\mu}}$ and $U_{\varSigma}$ are not specified. To use
this result into design optimization in practice, we have to reduce
$\max\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}\mid\boldsymbol{\mu}\in
U_{\boldsymbol{\mu}}\\}$ and
$\max\\{\varSigma\bullet\varLambda\mid\varSigma\in U_{\varSigma}\\}$ in (15)
to tractable forms. This is actually performed in section 3.2 and section 3.3,
where we consider two specific models of $U_{\boldsymbol{\mu}}$ and
$U_{\varSigma}$.

###  3.2 Uncertainty model with $\ell_{\infty}$-norm

Let $\tilde{\boldsymbol{\mu}}\in\mathbb{R}^{n}$ and
$\tilde{\varSigma}\in\mathcal{S}^{n}$ denote the best estimates of
$\boldsymbol{\mu}$ and $\varSigma$, respectively, where $\tilde{\varSigma}$ is
positive definite. In this section, we specialize the results of section 3.1
to the case that the uncertainty sets are given as

$\displaystyle U_{\boldsymbol{\mu}}$
$\displaystyle=\\{\tilde{\boldsymbol{\mu}}+A\boldsymbol{z}_{1}\mid\|\boldsymbol{z}_{1}\|_{\infty}\leq\alpha,\
\boldsymbol{z}_{1}\in\mathbb{R}^{m}\\},$ (19) $\displaystyle U_{\varSigma}$
$\displaystyle=\\{\tilde{\varSigma}+BZ_{2}B^{\top}\mid\|Z_{2}\|_{\infty,\infty}\leq\beta,\
Z_{2}\in\mathcal{S}^{k}\\}\cap\mathcal{S}_{+}^{n}.$ (20)

Here, $\boldsymbol{z}_{1}\in\mathbb{R}^{m}$ and $Z_{2}\in\mathcal{S}^{k}$ are
unknown vector and matrix reflecting the uncertainty in $\boldsymbol{\mu}$ and
$\varSigma$, respectively, $A\in\mathbb{R}^{n\times m}$ and
$B\in\mathbb{R}^{n\times k}$ are constant matrices, and $\alpha$ and $\beta$
are nonnegative parameters representing the magnitude of uncertainties.

######  Example 3.1.

A simple example of the uncertainty set in (19) is a box-constrained model.
For example, if we put $\tilde{\boldsymbol{\mu}}=\boldsymbol{0}$ and $A=I$
with $m=n$, (19) is reduced to

$\displaystyle
U_{\boldsymbol{\mu}}=\\{\boldsymbol{z}_{1}\in\mathbb{R}^{n}\mid\|\boldsymbol{z}_{1}\|_{\infty}\leq\alpha\\}.$

This means that the expected value vector $\boldsymbol{\mu}$ belongs to a
hypercube centered at the origin, with edges parallel to the axes and with an
edge length of $2\alpha$. In other words, each component $\mu_{j}$ of
$\boldsymbol{\mu}$ can take any value in $[-\alpha,\alpha]$. Similarly, a
simple example of the uncertainty set in (20) is the one with $B=I$ and $k=n$,
i.e.,

$\displaystyle
U_{\varSigma}=\\{\tilde{\varSigma}+Z_{2}\mid\|Z_{2}\|_{\infty,\infty}\leq\beta,\
Z_{2}^{\top}=Z_{2}\\}\cap\mathcal{S}_{+}^{n}.$

This means that, roughly speaking, the variance-covariance matrix $\varSigma$
has componentwise uncertainty. More precisely, for each $i$, $j=1,\dots,n$ we
have

$\displaystyle\tilde{\varSigma}_{ij}-\beta\leq\varSigma_{ij}\leq\tilde{\varSigma}_{ij}+\beta,\quad\varSigma_{ji}=\varSigma_{ij},$
(21)

and besides $\varSigma$ should be positive semidefinite. It is worth noting
that, even if $\tilde{\varSigma}$ and $\beta$ satisfy
$\tilde{\varSigma}-\beta\boldsymbol{1}\boldsymbol{1}^{\top}\succ 0$ and
$\tilde{\varSigma}+\beta\boldsymbol{1}\boldsymbol{1}^{\top}\succ 0$ (here,
$\boldsymbol{1}$ denotes an all-ones column vector), (21) does not necessarily
imply $\varSigma\succ 0$. Indeed, as for an example with $n=2$, consider

$\displaystyle\tilde{\varSigma}=\begin{bmatrix}3&2\\\ 2&3\\\
\end{bmatrix},\quad\beta=2.$

Then we have

$\displaystyle\tilde{\varSigma}-\beta\boldsymbol{1}\boldsymbol{1}^{\top}=\begin{bmatrix}1&0\\\
0&1\\\ \end{bmatrix}\succ
0,\quad\tilde{\varSigma}+\beta\boldsymbol{1}\boldsymbol{1}^{\top}=\begin{bmatrix}5&4\\\
4&5\\\ \end{bmatrix}\succ 0,$

and, for example, we see that

$\displaystyle\varSigma=\begin{bmatrix}2&3\\\ 3&2\\\ \end{bmatrix}$

satisfies (21) but $\varSigma\not\succ 0$. $\blacksquare$

To derive the main result in this section stated in Theorem 3.2, we need the
two technical lemmas. Lemma 3.1 explicitly computes the value of
$\max\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}\mid\boldsymbol{\mu}\in
U_{\boldsymbol{\mu}}\\}$ in (15). Lemma 3.2 converts
$\max\\{\varSigma\bullet\varLambda\mid\varSigma\in U_{\varSigma}\\}$ to a
tractable form.

######  Lemma 3.1.

For $U_{\boldsymbol{\mu}}$ defined by (19), we have

$\displaystyle\max_{\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}\\}=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|A^{\top}\nabla
g(\tilde{\boldsymbol{x}})\|_{1}.$

###### Proof.

Substitution of (19) into the left side yields

$\displaystyle\max_{\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}\\}$ $\displaystyle=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\max_{\boldsymbol{z}_{1}}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}A\boldsymbol{z}_{1}\mid\|\boldsymbol{z}_{1}\|_{\infty}\leq\alpha\\}$
$\displaystyle=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\max_{\boldsymbol{z}_{1}}\\{(A^{\top}\nabla
g(\tilde{\boldsymbol{x}}))^{\top}\boldsymbol{z}_{1}\mid\|\boldsymbol{z}_{1}\|_{\infty}\leq
1\\}.$

It is known that the dual norm of the $\ell_{\infty}$-norm is the
$\ell_{1}$-norm [7, appendix A.1.6], i.e.,

$\displaystyle\max_{\boldsymbol{t}\in\mathbb{R}^{n}}\\{\boldsymbol{s}^{\top}\boldsymbol{t}\mid\|\boldsymbol{t}\|_{\infty}\leq
1\\}=\|\boldsymbol{s}\|_{1}.$

Therefore, we obtain

$\displaystyle\max_{\boldsymbol{z}_{1}}\\{(A^{\top}\nabla
g(\tilde{\boldsymbol{x}}))^{\top}\boldsymbol{z}_{1}\mid\|\boldsymbol{z}_{1}\|_{\infty}\leq
1\\}=\|A^{\top}\nabla g(\tilde{\boldsymbol{x}})\|_{1},$

which concludes the proof. ∎

######  Lemma 3.2.

For $U_{\varSigma}$ defined by (20), we have

$\displaystyle\max_{\varSigma\in
U_{\varSigma}}\\{\varLambda\bullet\varSigma\\}=\min_{\varOmega\in\mathcal{S}_{+}^{k}}\\{\tilde{\varSigma}\bullet(\varLambda+\varOmega)+\beta\|B^{\top}(\varLambda+\varOmega)B\|_{1,1}\\}.$

###### Proof.

We shall show that the right side corresponds to the dual problem of the SDP
problem on the left side. Therefore, this proposition follows from the strong
duality of SDP [8, section 11.3], because the dual problem is strictly
feasible.

As preliminaries, for a convex cone defined by
$K=\\{(s,S)\in\mathbb{R}\times\mathcal{S}^{k}\mid\|S\|_{1,1}\leq s\\}$,
observe that its dual cone is given by [7, Example 2.25]

$\displaystyle\\{(t,T)\in\mathbb{R}\times\mathcal{S}^{k}\mid st+S\bullet T\geq
0\ (\forall(s,S)\in K)\\}$
$\displaystyle=\\{(t,T)\in\mathbb{R}\times\mathcal{S}^{k}\mid\|T\|_{\infty,\infty}\leq
t\\},$

from which we obtain

$\displaystyle\inf_{s\in\mathbb{R},\,S\in\mathcal{S}^{k}}\\{st+S\bullet
T\mid\|S\|_{1,1}\leq s\\}=\begin{dcases*}0&if $\|T\|_{\infty,\infty}\leq
t$,\\\ -\infty&otherwise.\end{dcases*}$ (22)

By using definition (20) of $U_{\varSigma}$, the left side of the proposition
of this theorem is reduced to

$\displaystyle\max_{\varSigma\in
U_{\varSigma}}\\{\varLambda\bullet\varSigma\\}=\max_{Z_{2}\in\mathcal{S}^{k}}\\{\varLambda\bullet(\tilde{\varSigma}+BZ_{2}B^{\top})\mid\|Z_{2}\|_{\infty,\infty}\leq\beta,\
\tilde{\varSigma}+BZ_{2}B^{\top}\succeq 0\\}.$ (23)

The Lagrangian of this optimization problem is defined by

$\displaystyle
L(Z_{2};v,V,\varOmega)=\begin{dcases*}\varLambda\bullet(\tilde{\varSigma}+BZ_{2}B^{\top})+(\beta
v+Z_{2}\bullet V)\\\
\qquad{}+\varOmega\bullet(\tilde{\varSigma}+BZ_{2}B^{\top})&if
$\|V\|_{1,1}\leq v$, $\varOmega\succeq 0$,\\\ +\infty&otherwise,\end{dcases*}$
(24)

where $v\in\mathbb{R}$, $V\in\mathcal{S}^{k}$, and
$\varOmega\in\mathcal{S}^{n}$ are the Lagrange multipliers. Indeed, by using
(11) and (22), we can confirm that problem (23) is equivalent to

$\displaystyle\max_{Z_{2}}\inf_{v,\,V,\,\varOmega}L(Z_{2};v,V,\varOmega).$

The dual problem is then defined by

$\displaystyle\min_{v,\,V,\,\varOmega}\sup_{Z_{2}}L(Z_{2};v,V,\varOmega).$
(25)

Since (24) can be rewritten as

$\displaystyle
L(Z_{2};v,V,\varOmega)=\begin{dcases*}Z_{2}\bullet(V+B^{\top}(\varLambda+\varOmega)B)\\\
\qquad{}+\tilde{\varSigma}\bullet(\varLambda+\varOmega)+\beta v&if
$\|V\|_{1,1}\leq v$, $\varOmega\succeq 0$,\\\ +\infty&otherwise,\end{dcases*}$

we obtain

$\displaystyle\sup_{Z_{2}}L(Z_{2};v,V,\varOmega)=\begin{dcases*}\tilde{\varSigma}\bullet(\varLambda+\varOmega)+\beta
v&if $\|V\|_{1,1}\leq v$, $\varOmega\succeq 0$,
$V=-B^{\top}(\varLambda+\varOmega)B$,\\\ +\infty&otherwise.\end{dcases*}$

Therefore, the dual problem in (25) is explicitly written as follows:

$\displaystyle\mathop{\mathrm{Minimize}}_{v\in
V,\,\varLambda\in\mathcal{S}^{k},\,\varOmega\in\mathcal{S}^{k}}$
$\displaystyle\tilde{\varSigma}\bullet(\varLambda+\varOmega)+\beta v$
$\displaystyle\mathop{\mathrm{subject~{}to}}$
$\displaystyle\|B^{\top}(\varLambda+\varOmega)B\|_{1,1}\leq v,$
$\displaystyle\varOmega\succeq 0.$

Constraint $\|B^{\top}(\varLambda+\varOmega)B\|_{1,1}\leq v$ becomes active at
an optimal solution, which concludes the proof. ∎

We are now in position to state the main result of this section. By using
Theorem 3.1, Lemma 3.1, and Lemma 3.2, we obtain the following fact.

######  Theorem 3.2.

Let $U_{\boldsymbol{\mu}}$ and $U_{\varSigma}$ be the sets defined by (19) and
(20), respectively. Then, $\tilde{\boldsymbol{x}}\in X$ satisfies (14) if and
only if there exists a pair of $z\in\mathbb{R}$ and $W\in\mathcal{S}^{n}$
satisfying

$\displaystyle g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|A^{\top}\nabla
g(\tilde{\boldsymbol{x}})\|_{1}+\tilde{\varSigma}\bullet
W+\beta\|B^{\top}WB\|_{1,1}+\kappa^{2}z\leq 0,$ (26)
$\displaystyle\begin{bmatrix}W&\nabla g(\tilde{\boldsymbol{x}})/2\\\ \nabla
g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\ \end{bmatrix}\succeq 0.$ (27)

###### Proof.

It follows from Lemma 3.1 and Lemma 3.2 that (15) and (16) in Theorem 3.1 are
equivalently rewritten as

$\displaystyle g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|A^{\top}\nabla
g(\tilde{\boldsymbol{x}})\|_{1}+\tilde{\varSigma}\bullet(\varLambda+\varOmega)+\beta\|B^{\top}(\varLambda+\varOmega)B\|_{1,1}+\kappa^{2}z\leq
0,$ $\displaystyle\varOmega\succeq 0,$
$\displaystyle\begin{bmatrix}\varLambda&\nabla g(\tilde{\boldsymbol{x}})/2\\\
\nabla g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\ \end{bmatrix}\succeq 0.$

Put $W=\varLambda+\varOmega$ to see that this is reduced to

$\displaystyle g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|A^{\top}\nabla
g(\tilde{\boldsymbol{x}})\|_{1}+\tilde{\varSigma}\bullet
W+\beta\|B^{\top}WB\|_{1,1}+\kappa^{2}z\leq 0,$ $\displaystyle
W-\varLambda\succeq 0,$ $\displaystyle\begin{bmatrix}\varLambda&\nabla
g(\tilde{\boldsymbol{x}})/2\\\ \nabla g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\
\end{bmatrix}\succeq 0.$

This is straightforwardly equivalent to (26) and (27). ∎

It should be emphasized that Theorem 3.2 converts the set of infinitely many
reliability constraints in (13) to two deterministic constraints, i.e., (26)
and (27). The latter constraints can be handled within the framework of
conventional (deterministic) optimization.

###  3.3 Uncertainty model with $\ell_{2}$-norm

In this section, we consider the uncertainty sets defined by

$\displaystyle U_{\boldsymbol{\mu}}$
$\displaystyle=\\{\tilde{\boldsymbol{\mu}}+A\boldsymbol{z}_{1}\mid\|\boldsymbol{z}_{1}\|_{2}\leq\alpha,\
\boldsymbol{z}_{1}\in\mathbb{R}^{m}\\},$ $\displaystyle U_{\varSigma}$
$\displaystyle=\\{\tilde{\varSigma}+BZ_{2}B^{\top}\mid\|Z_{2}\|_{\mathrm{F}}\leq\beta,\
Z_{2}\in\mathcal{S}^{k}\\}\cap\mathcal{S}_{+}^{n}.$

######  Example 3.2.

As a simple example, put $\tilde{\boldsymbol{\mu}}=\boldsymbol{0}$ and $A=I$
with $m=n$ to obtain

$\displaystyle
U_{\boldsymbol{\mu}}=\\{\boldsymbol{z}_{1}\in\mathbb{R}^{n}\mid\|\boldsymbol{z}_{1}\|_{2}\leq\alpha\\}.$

This means that the expected value vector $\boldsymbol{\mu}$ belongs to a
hypersphere centered at the origin with radius $\alpha$. Similarly, putting
$B=I$ and $k=n$ we obtain

$\displaystyle
U_{\varSigma}=\\{\tilde{\varSigma}+Z_{2}\mid\|Z_{2}\|_{\mathrm{F}}\leq\beta,\
Z_{2}^{\top}=Z_{2}\\}\cap\mathcal{S}_{+}^{n}.$

This means that the variance-covariance matrix $\varSigma$ satisfies

$\displaystyle\sum_{i=1}^{n}\sum_{j=1}^{n}(\varSigma_{ij}-\tilde{\varSigma}_{ij})^{2}\leq\beta^{2}$

and is symmetric positive semidefinite. $\blacksquare$

In a manner parallel to the proofs of Lemma 3.1 and Lemma 3.2, we can obtain

$\displaystyle\max_{\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}}\\{\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\mu}\\}$ $\displaystyle=\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|A^{\top}\nabla
g(\tilde{\boldsymbol{x}})\|_{2},$ $\displaystyle\max_{\varSigma\in
U_{\varSigma}}\\{\varLambda\bullet\varSigma\\}$
$\displaystyle=\min_{\varOmega\in\mathcal{S}^{k}}\\{\tilde{\varSigma}\bullet(\varLambda+\varOmega)+\beta\|B^{\top}(\varLambda+\varOmega)B\|_{\mathrm{F}}\mid\varOmega\succeq
0\\}.$

Here, the facts

$\displaystyle\max_{\boldsymbol{t}\in\mathbb{R}^{n}}\\{\boldsymbol{s}^{\top}\boldsymbol{t}\mid\|\boldsymbol{t}\|_{2}\leq
1\\}=\|\boldsymbol{s}\|_{2}$

and

$\displaystyle\\{(t,T)\in\mathbb{R}\times\mathcal{S}^{k}\mid st+S\bullet T\geq
0\ (\forall(s,S):\|S\|_{\mathrm{F}}\leq s)\\}$
$\displaystyle=\\{(t,T)\in\mathbb{R}\times\mathcal{S}^{k}\mid\|T\|_{\mathrm{F}}\leq
t\\}$

have been used.

Accordingly, analogous to Theorem 3.2, we obtain the following conclusion:
$\tilde{\boldsymbol{x}}\in X$ satisfies (14) if and only if there exists a
pair of $z\in\mathbb{R}$ and $W\in\mathcal{S}^{n}$ satisfying

$\displaystyle g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|A^{\top}\nabla
g(\tilde{\boldsymbol{x}})\|_{2}+\tilde{\varSigma}\bullet
W+\beta\|B^{\top}WB\|_{\mathrm{F}}+\kappa^{2}z\leq 0,$
$\displaystyle\begin{bmatrix}W&\nabla g(\tilde{\boldsymbol{x}})/2\\\ \nabla
g(\tilde{\boldsymbol{x}})^{\top}/2&z\\\ \end{bmatrix}\succeq 0.$

###  3.4 Truss optimization under compliance constraint

In this section, we present how the results established in the preceding
sections can be employed for a specific RBDO problem. As a simple example, we
consider a reliability constraint on the compliance under a static external
load. We assume linear elasticity and small deformation.

For ease of comprehension, consider design optimization of a truss. In this
context, $x_{j}$ denotes the cross-sectional area of truss member $j$
$(j=1,\dots,n)$, where $n$ is the number of members. We attempt to minimize
the structural volume of the truss, $\boldsymbol{c}^{\top}\boldsymbol{x}$,
under the compliance constraint, where $c_{j}$ denotes the undeformed member
length. Let $\pi(\boldsymbol{x})$ denote the compliance corresponding to a
static external load. The first-order approximation of the compliance
constraint is written as

$\displaystyle\pi(\tilde{\boldsymbol{x}})+\nabla\pi(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq\bar{\pi},$

where $\bar{\pi}$ $(>0)$ is a specified upper bound for the compliance.
Accordingly, the design optimization problem to be solved is formulated as
follows:

$\displaystyle\mathop{\mathrm{Minimize}}_{\tilde{\boldsymbol{x}}}$
$\displaystyle\boldsymbol{c}^{\top}\tilde{\boldsymbol{x}}$ (28a)
$\displaystyle\mathop{\mathrm{subject~{}to}}$
$\displaystyle\tilde{\boldsymbol{x}}\geq\bar{\boldsymbol{x}},$ (28b)
$\displaystyle\mathrm{P}_{\mathsf{N}(\boldsymbol{\mu},\varSigma)}\\{\pi(\tilde{\boldsymbol{x}})+\nabla\pi(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq\bar{\pi}\\}\geq
1-\epsilon,\quad\forall(\boldsymbol{\mu},\varSigma)\in U.$ (28c)

Here, the specified lower bound for the member cross-sectional area, denoted
by $\bar{x}_{j}$ $(j=1,\dots,n)$, is positive, because in this paper we
restrict ourselves to optimization problems without variation of structural
topology.

As for uncertainty sets of the moments, consider, for example,
$U_{\boldsymbol{\mu}}$ and $U_{\varSigma}$ studied in section 3.3. For
simplicity put $A=B=I$ so as to obtain

$\displaystyle U_{\boldsymbol{\mu}}$
$\displaystyle=\\{\tilde{\boldsymbol{\mu}}+\boldsymbol{z}_{1}\mid\|\boldsymbol{z}_{1}\|_{2}\leq\alpha,\
\boldsymbol{z}_{1}\in\mathbb{R}^{n}\\},$ $\displaystyle U_{\varSigma}$
$\displaystyle=\\{\tilde{\varSigma}+Z_{2}\mid\|Z_{2}\|_{\mathrm{F}}\leq\beta,\
Z_{2}\in\mathcal{S}^{n}\\}\cap\mathcal{S}_{+}^{n}.$

From the result in section 3.3, we see that problem (28) is equivalently
rewritten as follows:

$\displaystyle\mathop{\mathrm{Minimize}}_{\tilde{\boldsymbol{x}},\,z,\,W}$
$\displaystyle\boldsymbol{c}^{\top}\tilde{\boldsymbol{x}}$ (29a)
$\displaystyle\mathop{\mathrm{subject~{}to}}$
$\displaystyle\tilde{\boldsymbol{x}}\geq\bar{\boldsymbol{x}},$ (29b)
$\displaystyle\pi(\tilde{\boldsymbol{x}})+\nabla\pi(\tilde{\boldsymbol{x}})^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|\nabla\pi(\tilde{\boldsymbol{x}})\|_{2}+\tilde{\varSigma}\bullet
W+\beta\|W\|_{\mathrm{F}}+\kappa^{2}z\leq\bar{\pi},$ (29c)
$\displaystyle\begin{bmatrix}W&\nabla\pi(\tilde{\boldsymbol{x}})/2\\\
\nabla\pi(\tilde{\boldsymbol{x}})^{\top}/2&z\\\ \end{bmatrix}\succeq 0.$ (29d)

Here, $\tilde{\boldsymbol{x}}\in\mathbb{R}^{n}$, $z\in\mathbb{R}$, and
$W\in\mathcal{S}^{n}$ are variables to be optimized. It is worth noting that
problem (29) is a nonlinear SDP problem.

The remainder of this section is devoted to presenting a method for solving
problem (29) that will be used for the numerical experiments in section 5.

The method sequentially solves SDP problems that approximate problem (29), in
a fashion similar to sequential SDP methods for nonlinear SDP problems [48,
32]. Let $\tilde{\boldsymbol{x}}^{k}$ denote the incumbent solution obtained
at iteration $k-1$. Define $\boldsymbol{h}_{k}\in\mathbb{R}^{n}$ by

$\displaystyle\boldsymbol{h}_{k}=\nabla\pi(\tilde{\boldsymbol{x}}^{k}).$

At iteration $k$, we replace $\nabla\pi(\tilde{\boldsymbol{x}})$ in (29c) and
(29d) with $\boldsymbol{h}_{k}$. Moreover, to deal with
$\pi(\tilde{\boldsymbol{x}})$ in (29c), we use the fact that $s\in\mathbb{R}$
satisfies

$\displaystyle\pi(\boldsymbol{x})\leq s$

if and only if

$\displaystyle\begin{bmatrix}K(\tilde{\boldsymbol{x}})&\boldsymbol{p}\\\
\boldsymbol{p}^{\top}&s\\\ \end{bmatrix}\succeq 0$ (30)

is satisfied [27, section 3.1], where
$K(\tilde{\boldsymbol{x}})\in\mathcal{S}^{d}$ is the stiffness matrix of the
truss, $\boldsymbol{p}\in\mathbb{R}^{d}$ is the external load vector, and $d$
is the number of degrees of freedom of the nodal displacements. It is worth
noting that, for trusses, $K(\tilde{\boldsymbol{x}})$ is linear in
$\tilde{\boldsymbol{x}}$. Therefore, (30) is a linear matrix inequality with
respect to $\tilde{\boldsymbol{x}}$ and $s$, and hence can be handled within
the framework of (linear) SDP. By this means, we obtain the following
subproblem that is solved at iteration $k$ for updating
$\tilde{\boldsymbol{x}}^{k}$ to $\tilde{\boldsymbol{x}}^{k+1}$:

$\displaystyle\mathop{\mathrm{Minimize}}_{\tilde{\boldsymbol{x}},\,z,\,w,\,s}$
$\displaystyle\boldsymbol{c}^{\top}\tilde{\boldsymbol{x}}$ (31a)
$\displaystyle\mathop{\mathrm{s.{\,}t.}}$
$\displaystyle\tilde{\boldsymbol{x}}\geq\bar{\boldsymbol{x}},$ (31b)
$\displaystyle\begin{bmatrix}k(\tilde{\boldsymbol{x}})&\boldsymbol{p}\\\
\boldsymbol{p}^{\top}&s\\\ \end{bmatrix}\succeq 0,$ (31c) $\displaystyle
s+\boldsymbol{h}_{k}^{\top}\tilde{\boldsymbol{\mu}}+\alpha\|\boldsymbol{h}_{k}\|_{2}+\tilde{\varSigma}\bullet
w+\beta\|w\|_{\mathrm{F}}+\kappa^{2}z\leq\bar{\pi},$ (31d)
$\displaystyle\begin{bmatrix}w&\boldsymbol{h}_{k}/2\\\
\boldsymbol{h}_{k}^{\top}/2&z\\\ \end{bmatrix}\succeq 0.$ (31e)

Since this is a linear SDP problem, we can solve this problem efficiently with
a primal-dual interior-point method [1].

##  4 Extensions

This section discusses some extensions of the results obtained in section 3.

###  4.1 Robustness against uncertainty in distribution type

An important extension is that the obtained results can be applied to the case
that, not only the moments, but also the type of probability distribution are
unknown. In this case, we consider any combination of all types of probability
distributions and all possible moments (expected value vectors and variance-
covariance matrices) in the uncertainty set, and require that the failure
probability is no greater than a specified value. This robustness against
uncertainty in distribution type is important as the input distribution in
practice is not necessarily known to be a normal distribution.

Recall that, in section 2 and section 3, we assumed that the design variables,
$\boldsymbol{x}$, follows a normal distribution. Then we consider the robust
reliability constraint in (14). For the sake of clarity, we restate this
problem setting in a slightly different manner. We have assumed that random
vector $\boldsymbol{\zeta}$ can possibly follow any normal distribution
satisfying $\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}$ and $\varSigma\in
U_{\varSigma}$. We use $\mathcal{P}_{\mathsf{N}}$ to denote the set of such
normal distributions, i.e.,

$\displaystyle\mathcal{P}_{\mathsf{N}}=\\{\mathsf{N}(\boldsymbol{\mu},\varSigma)\mid\boldsymbol{\mu}\in
U_{\boldsymbol{\mu}},\ \varSigma\in U_{\varSigma}\\}.$ (32)

In other words, $\mathcal{P}_{\mathsf{N}}$ is the set of all possible
realizations of the input distribution. We write
$p\in\mathcal{P}_{\mathsf{N}}$ if $p$ is one of such realizations. With this
new notation, (14) can be rewritten equivalently as

$\displaystyle\sup_{p\in\mathcal{P}_{\mathsf{N}}}\bigl{\\{}\mathrm{P}_{p}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0\\}\bigr{\\}}\geq
1-\epsilon.$ (33)

For $U_{\boldsymbol{\mu}}$ and $U_{\varSigma}$ defined in section 3.2, Theorem
3.2 shows that (33) is equivalent to (26) and (27).

We are now in position to consider any type of probability distribution. Only
what we assume is that the input distribution satisfies $\boldsymbol{\mu}\in
U_{\boldsymbol{\mu}}$ and $\varSigma\in U_{\varSigma}$, where, for a while, we
consider $U_{\boldsymbol{\mu}}$ and $U_{\varSigma}$ defined in section 3.2. We
use $\mathcal{P}$ to denote the set of such distributions, i.e.,

$\displaystyle\mathcal{P}=\\{p\mid p\text{ is a distribution with
$\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}$ \& $\varSigma\in
U_{\varSigma}$}\\}.$ (34)

Then, instead of (33), we consider the following constraint:

$\displaystyle\sup_{p\in\mathcal{P}}\bigl{\\{}\mathrm{P}_{p}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0\\}\bigr{\\}}\geq
1-\epsilon.$ (35)

That is, we require that the reliability constraint should be satisfied for
any input distribution $p$ satisfying $p\in\mathcal{P}$. A main assertion of
this section is that, by simply setting

$\displaystyle\kappa=\sqrt{\frac{1-\epsilon}{\epsilon}}$ (36)

instead of $\kappa=-\varPhi^{-1}(\epsilon)$, constraint (35) is equivalent to
(26) and (27) in Theorem 3.2. We can show this fact in the following manner.
Let $\mathcal{P}(\boldsymbol{\mu},\varSigma)$ denote the set of probability
distributions, the expected value vector and the variance-covariance matrix of
which are $\boldsymbol{\mu}$ and $\varSigma$, respectively. Observe that, with
$\mathcal{P}(\boldsymbol{\mu},\varSigma)$, (35) can be rewritten equivalently
as

$\displaystyle\sup_{(\boldsymbol{\mu},\varSigma)\in
U}\left\\{\sup_{p\in\mathcal{P}(\boldsymbol{\mu},\varSigma)}\bigl{\\{}\mathrm{P}_{p}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq
0\\}\bigr{\\}}\right\\}\geq 1-\epsilon.$ (37)

With relation to the inner supremum, consider the condition

$\displaystyle\sup_{p\in\mathcal{P}(\boldsymbol{\mu},\varSigma)}\bigl{\\{}\mathrm{P}_{p}\\{g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0\\}\bigr{\\}}\geq
1-\epsilon.$ (38)

El Ghaoui et al. [12, Theorem 1] show that (38) holds if and only if (7) of
Theorem 2.1 holds with $\kappa$ defined by (36). Therefore, all the subsequent
results established in section 2 and section 3 hold by simply replacing the
value of $\kappa$ with the one in (36). Thus, the robust reliability
constraint with unknown distribution type is also reduced to the form in (26)
and (27) of Theorem 3.2.

The result in section 3.3, which are established for the $\ell_{2}$-norm
uncertainty model, is also extended to the case of unknown distribution type
by replacing $\kappa$ with the value in (36).

###  4.2 Multiple constraints

In section 2 and section 3, we have restricted ourselves to the case that the
design optimization problem has a single performance requirement, (1). In this
section, we discuss treatment of multiple constraints.

Suppose that the performance requirement is written as

$\displaystyle g_{i}(\boldsymbol{x})\leq 0,\quad i=1,\dots,m.$

The first-order approximation yields

$\displaystyle g_{i}^{\mathrm{lin}}(\boldsymbol{x})\leq 0,\quad i=1,\dots,m,$

where
$g_{i}^{\mathrm{lin}}(\boldsymbol{x})=g_{i}(\tilde{\boldsymbol{x}})+\nabla
g_{i}(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}$ $(i=1,\dots,m)$.
Suppose that we impose a distributionally-robust reliability constraint for
each $i=1,\dots,m$ independently, i.e.,

$\displaystyle\sup_{p\in\mathcal{P}}\bigl{\\{}\mathrm{P}_{p}\\{g_{i}^{\mathrm{lin}}(\boldsymbol{x})\leq
0\\}\bigr{\\}}\geq 1-\epsilon,\quad i=1,\dots,m.$ (39)

Here, $\mathcal{P}$ is the set of possible realizations of the input
distribution (i.e., $\mathcal{P}$ here is either $\mathcal{P}_{\mathsf{N}}$ in
(32) or $\mathcal{P}$ in (34)). It is worth noting that in (39) the worst case
distributions are considered independently for each $i=1,\dots,m$. Constraint
(39) can be straightforwardly dealt with in the same manner as section 3.

In contrast, suppose that we consider a single (i.e., common) worst-case
distribution for all $i=1,\dots,m$. Then the distributionally-robust
reliability constraint is written as

$\displaystyle\sup_{p\in\mathcal{P}}\bigl{\\{}\mathrm{P}_{p}\\{g_{i}^{\mathrm{lin}}(\boldsymbol{x})\leq
0\ (i=1,\dots,m)\\}\bigr{\\}}\geq 1-\epsilon.$ (40)

Treatment of this constraint remains to be studied as future work. It is worth
noting that constraint (39) is conservative compared with constraint (40).

##  5 Numerical examples

In section 3.4 we have seen that an optimization problem of trusses under the
compliance constraint is reduced to problem (29). In this section we solve
this optimization problem numerically.

The algorithm presented in section 3.4 was implemented in Matlab ver.
9.8.0.666Source codes for solving the optimization problems presented in
section 5 are available on-line at https://github.com/ykanno22/moment_worst/.
The SDP problem in (31) was solved by CVX ver. 2.2 [13, 14] with SeDuMi ver.
1.3.4 [42, 40]. Computation was carried out on a 2.6 GHz Intel Core i7-9750H
processor with 32 GB RAM.

$1$$2$$1\,\mathrm{m}$$1\,\mathrm{m}$

Figure 1: Problem setting of example (I): 2-bar truss.  Table 1: Optimal solutions of example (I) with $\epsilon=0.01$, $\alpha=0.2$, and $\beta=0.01$ (probability distributions are assumed to be normal distributions).  |  $x_{1}$ ($\mathrm{mm}^{2}$) |  $x_{2}$ ($\mathrm{mm}^{2}$) | Obj. val. ($\mathrm{mm^{3}}$) |  $\pi(\tilde{\boldsymbol{x}})$ (J)
---|---|---|---|---
Nominal optim. | 1500.0 | 2121.3 | $4.5000\times 10^{6}$ | 100.000
$\ell_{\infty}$-norm unc. | 1558.0 | 2203.4 | $4.6741\times 10^{6}$ | 96.274
$\ell_{2}$-norm unc. | 1535.4 | 2171.4 | $4.6063\times 10^{6}$ | 97.692

(a)

(b)

(c)

Figure 2: Example of Monte Carlo simulation for a single sample of
$\boldsymbol{\mu}$ and $\varSigma$ (example (I) with the $\ell_{\infty}$-norm
uncertainty model; probability distributions are assumed to be normal
distributions). LABEL:sub@fig:2bar_ell1_x1 Samples of
$x_{1}=\tilde{x}_{1}+\zeta_{1}$; LABEL:sub@fig:2bar_ell1_x2 samples of
$x_{2}=\tilde{x}_{2}+\zeta_{2}$; and LABEL:sub@fig:2bar_ell1_sample_app_compl
computed values of $g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}$.

(a)

(b)

Figure 3: Results of double-loop Monte Carlo simulation (example (I) with the
$\ell_{\infty}$-norm uncertainty model; probability distributions are assumed
to be normal distributions). LABEL:sub@fig:2bar_ell1_prob_of_failure Failure
probability of linearly approximated constraint
$g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0$; and
LABEL:sub@fig:2bar_ell1_exact_failure_prob failure probability of constraint
$g(\tilde{\boldsymbol{x}}+\boldsymbol{\zeta})\leq 0$ without approximation.

(a)

(b)

Figure 4: Optimal value (example (I) with the $\ell_{\infty}$-norm uncertainty
model; probability distributions are assumed to be normal distributions)
versus LABEL:sub@fig:2bar_ell1_loop_reliability failure probability; and
LABEL:sub@fig:2bar_ell1_loop_uncertainty magnitude of uncertainty.  Figure 5:
Results of double-loop Monte Carlo simulation (example (I) with the
$\ell_{2}$-norm uncertainty model; probability distributions are assumed to be
normal distributions). Failure probability of linearly approximated constraint
$g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}\leq 0$.

(a)

(b)

Figure 6: Optimal value (example (I) with the $\ell_{2}$-norm uncertainty
model; probability distributions are assumed to be normal distributions)
versus LABEL:sub@fig:2bar_ell2_loop_reliability failure probability; and
LABEL:sub@fig:2bar_ell2_loop_uncertainty magnitude of uncertainty.

(a)

(b)

(c)

(d)

Figure 7: Optimal value (example (I); no restriction on distribution type is
assumed) versus LABEL:sub@fig:2bar_ell1_all_dist_reliability failure
probability (with the $\ell_{\infty}$-norm uncertainty model);
LABEL:sub@fig:2bar_ell1_all_dist_uncertainty magnitude of uncertainty
($\ell_{\infty}$-norm); LABEL:sub@fig:2bar_ell2_all_dist_reliability failure
probability (with the $\ell_{2}$-norm uncertainty model); and
LABEL:sub@fig:2bar_ell2_all_dist_uncertainty magnitude of uncertainty
($\ell_{2}$-norm).

###  5.1 Example (I): 2-bar truss

Consider a plane truss depicted in Figure 1. The truss has $n=2$ members and
$d=2$ degrees of freedom of the nodal displacements. The elastic modulus of
the members is $20\,\mathrm{GPa}$. A vertical external force of
$100\,\mathrm{kN}$ is applied at the free node. The upper bound for the
compliance is $\bar{\pi}=100\,\mathrm{J}$.

We first consider the uncertainty model with the $\ell_{\infty}$-norm, studied
in section 3.2. In the uncertainty model in (19) and (20), we put $A=B=I$ with
$m=k=n$, as considered in Example 3.1. The best estimates, or the nominal
values, of $\boldsymbol{\mu}$ and $\varSigma$ are set to

$\displaystyle\tilde{\boldsymbol{\mu}}=\boldsymbol{0},\quad\tilde{\varSigma}=\begin{bmatrix}0.07&0.02\\\
0.02&0.07\\\ \end{bmatrix}.$

The magnitude of uncertainty is $\alpha=0.2$ and $\beta=0.01$. The specified
upper bound for the failure probability is $\epsilon=0.01$. The optimal
solution obtained by the proposed method is listed in the row
“$\ell_{\infty}$-norm unc.” of Table 1, where “obj. val.” means the objective
value at the obtained solution. For comparison, the optimal solution of the
nominal optimization problem (i.e., the conventional structural volume
minimization under the compliance constraint without considering uncertainty)
is also listed.

The optimization result was verified as follows. We randomly generate
$\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}$ and $\varSigma\in U_{\varSigma}$,
and then generate $10^{6}$ samples drawn as
$\boldsymbol{\zeta}\sim\mathsf{N}(\boldsymbol{\mu},\varSigma)$. Figure 2a and
Figure 2b show the samples of
$\boldsymbol{x}=\tilde{\boldsymbol{x}}+\boldsymbol{\zeta}$ generated in this
manner. Figure 2c shows the values of the linearly approximated constraint
function,

$\displaystyle g(\tilde{\boldsymbol{x}})+\nabla
g(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}=\pi(\tilde{\boldsymbol{x}})+\nabla\pi(\tilde{\boldsymbol{x}})^{\top}\boldsymbol{\zeta}-\bar{\pi},$

for these samples. Therefore, the ratio of the number of samples of which
these function values are positive divided by the number of all samples (i.e.,
$10^{6}$) should be no greater than $\epsilon$ $(=0.01)$. We computed this
ratio for each of $10^{4}$ randomly generated samples of $\boldsymbol{\mu}\in
U_{\boldsymbol{\mu}}$ and $\varSigma\in U_{\varSigma}$, where the continuous
uniform distribution was used to generate samples of the components of
$\boldsymbol{\mu}$ and $\varSigma$. Figure 3a shows the histogram of the
values of this ratio computed in this manner, i.e., it shows distribution of
the failure probability estimated by double-loop Monte Carlo simulation. It is
observed in Figure 3a that, for every one of $10^{4}$ probability distribution
samples, the failure probability is no greater than $\epsilon$. Thus, it is
verified that the obtained solution satisfies the distributionally-robust
reliability constraint in (14). Indeed, among these samples of the failure
probability, the maximum value is $0.009054$ $(<\epsilon)$. For reference,
Figure 3b shows the histogram of failure probabilities computed for the
constraint function values without applying the linear approximation, i.e.,
$g(\boldsymbol{x})=\pi(\boldsymbol{x})-\bar{\pi}$. It is observed in Figure 3b
that, in only rare cases, the failure probability exceeds the target value
$\epsilon$ $(=0.01)$. Figure 4a shows the variation of the optimal value with
respect to the upper bound for the failure probability $\epsilon$, where
$\alpha=0.2$ and $\beta=0.01$ are fixed. As $\epsilon$ decreases, the optimal
value increases. In contrast, Figure 4b shows the variation of the optimal
value with respect to $\alpha$ and $\beta$, where $\epsilon=0.01$ is fixed.
Although in Figure 4b only values of $\alpha$ are shown, values of
$\beta\in[0,0.02]$ are also varied in a manner proportional to $\alpha$. As
the magnitude of uncertainty increases, the optimal value increases.

We next consider the uncertainty model with the $\ell_{2}$-norm, studied in
section 3.3. The uncertainty set is defined with $A$, $B$,
$\tilde{\boldsymbol{\mu}}$, $\tilde{\varSigma}$, $\alpha$, and $\beta$ used
above. The specified upper bound for the failure probability is
$\epsilon=0.01$. The obtained optimal solution is listed in the row
“$\ell_{2}$-norm unc.” of Table 1. It can be observed that the objective value
is small compared with the solution with the $\ell_{\infty}$-norm uncertainty
model. This is natural, because, with the common values of $\alpha$ and
$\beta$, the uncertainty set with the $\ell_{2}$-norm is included in the
uncertainty set with the $\ell_{\infty}$-norm. The optimization result is
verified in the same manner as above. Namely, Figure 5 shows $10^{4}$ samples
of the failure probability, each of which was computed with $10^{6}$ samples
of $\boldsymbol{\zeta}$. Among these samples, the maximum failure probability
is $0.009850$ $(<\epsilon)$, which verifies that the obtained solution
satisfies distributionally-robust reliability constraint (14). Figure 6a and
Figure 6b show the variations of the optimal value with respect to the failure
probability, $\epsilon$, and the magnitude of uncertainty, $\alpha$ and
$\beta$, respectively. These variations show trends similar to the ones with
the $\ell_{\infty}$-norm uncertainty model in Figure 4a and Figure 4b.

Finally, as discussed in section 4.1, we consider, not only the normal
distributions, but all the probability distributions with $\boldsymbol{\mu}$
and $\varSigma$ belonging to the uncertainty set. That is, the set of possible
realizations of probability distributions is given by (34). Figure 7 collects
the variations of the optimal value with respect to the failure probability,
$\epsilon$, and the magnitude of uncertainty, $\alpha$ and $\beta$ (in the
same manner as above, values of $\beta\in[0,0.02]$ are varied in a manner
proportional to $\alpha$). Compared with the results for normal distributions
in Figure 4 and Figure 6, the optimal value in Figure 7 is large, as expected.
Moreover, as $\epsilon$ decreases, the optimal value in Figure 7a and Figure
7c increases drastically compared with the cases in Figure 4a and Figure 6a.

###  5.2 Example (II): 29-bar truss

$1\,\mathrm{m}$$1\,\mathrm{m}$$1\,\mathrm{m}$$1\,\mathrm{m}$$1\,\mathrm{m}$

Figure 8: Problem setting of example (II): 29-bar truss.  Table 2: Optimal solutions of example (II) with $\epsilon=0.01$, $\alpha=0.2$, and $\beta=0.01$ (probability distributions are assumed to be normal distributions).  | Obj. val. ($\mathrm{mm^{3}}$) |  $\pi(\tilde{\boldsymbol{x}})$ (J)
---|---|---
Nominal optim. | $1.6616\times 10^{7}$ | 1000.00
$\ell_{\infty}$-norm unc. | $1.7918\times 10^{7}$ | 917.66
$\ell_{2}$-norm unc. | $1.7475\times 10^{7}$ | 944.21

(a)

(b)

Figure 9: Optimal value (example (II) with the $\ell_{2}$-norm uncertainty
model; probability distributions are assumed to be normal distributions)
versus LABEL:sub@fig:prog2_ell2_loop_reliability failure probability; and
LABEL:sub@fig:prog2_ell2_loop_uncertainty magnitude of uncertainty.

(a)

(b)

Figure 10: Optimal value (example (II) with the $\ell_{2}$-norm uncertainty
model; no restriction on distribution type is assumed) versus
LABEL:sub@fig:prog2ell2_all_dist_reliability failure probability; and
LABEL:sub@fig:prog2ell2_all_dist_uncertainty magnitude of uncertainty.

(a)

(b)

(c)

Figure 11: Obtained designs of example (II). Optimal solutions
LABEL:sub@fig:prog2_nominal without considering uncertainty;
LABEL:sub@fig:prog2_ell1_design with the $\ell_{\infty}$-norm uncertainty
model; and LABEL:sub@fig:prog2_ell2_design with the $\ell_{2}$-norm
uncertainty model. In LABEL:sub@fig:prog2_ell1_design and
LABEL:sub@fig:prog2_ell2_design, no restriction on distribution type is
assumed, $\epsilon=0.005$, $\alpha=0.2$, and $\beta=0.01$.

Consider a plane truss depicted in Figure 8, where $n=29$ and $d=20$. The
elastic modulus of the members is $20\,\mathrm{GPa}$. Vertical external forces
of $100\,\mathrm{kN}$ are applied at two nodes as shown in Figure 8. The upper
bound for the compliance is $\bar{\pi}=1000\,\mathrm{J}$. The lower bounds for
the member cross-sectional areas are $\bar{x}_{j}=200\,\mathrm{mm^{2}}$
$(j=1,\dots,n)$.

As for the uncertainty model, we consider both the model with the
$\ell_{2}$-norm, putting $A=B=I$ with $m=k=n$. The best estimates of
$\boldsymbol{\mu}$ and $\varSigma$ are

$\displaystyle\tilde{\boldsymbol{\mu}}=\boldsymbol{0},\quad\tilde{\varSigma}=0.05I+0.02(\boldsymbol{1}\boldsymbol{1}^{\top}),$

where $\boldsymbol{1}\in\mathbb{R}^{n}$ is an all-ones column vector. The
magnitude of uncertainty is $\alpha=0.2$ and $\beta=0.01$. The specified upper
bound for the failure probability is $\epsilon=0.01$.

The optimization results obtained by the proposed method are listed in Table
2. Figure 9a and Figure 9b show the variations of the optimal value with
respect to the failure probability and the magnitude of uncertainty,
respectively.

As done in section 4.1, we next require that the reliability constraint should
be satisfied for all the probability distributions satisfying
$\boldsymbol{\mu}\in U_{\boldsymbol{\mu}}$ and $\varSigma\in U_{\varSigma}$,
i.e., for any probability distribution belonging to $\mathcal{P}$ in (34). For
the $\ell_{2}$-norm uncertainty, Figure 10a and Figure 10b report the
variations of the optimal value with respect to the failure probability and
the magnitude of uncertainty, respectively. Figure 11 collects the optimal
solutions of the optimization problem without uncertainty, as well as the
distributionally-robust RBDO problems with the two uncertainty models. Here,
the width of each member in the figures are proportional to its cross-
sectional area.

##  6 Conclusions

This paper has dealt with reliability-based design optimization (RBDO) of
structures, in which knowledge of the input distribution that the design
variables follow is imprecise. Specifically, we only know that the expected
value vector and the variance-covariance matrix of the input distribution
belong to a specified convex set, and do not know their true values. Then we
attempt to optimize a structure, under the constraint that, even for the
worst-case input distribution, the failure probability of the structure is no
greater than the specified value. This constraint, called the
distributionally-robust reliability constraint, is equivalent to infinitely
many reliability constraints corresponding to all possible realizations of the
input distribution. Provided that change of a constraint function value is
well approximated as a linear function of uncertain perturbations of the
design variables, this paper has presented a tractable reformulation of the
distributionally-robust reliability constraint.

This paper has established the concept of distributionally-robust RBDO, and
developed fundamental results. Much remains to be studied. For instance, in
this paper we have considered uncertainty only in the design variables. Other
sources of uncertainty in structural optimization can be explored. Also, as
discussed in section 4.2, multiple performance requirement in the form of (40)
remains to be studied. Extension to topology optimization is of great
interest. Moreover, this paper relies on the assumption that quantity of
interest is approximated, with sufficient accuracy, as a linear function of
uncertainty perturbations of the design variables. Extension to nonlinear
cases can be attempted. Finally, development of a more efficient algorithm for
solving the optimization problem presented in this paper can be studied.

#### Acknowledgments

This work is supported by Research Grant from the Maeda Engineering Foundation
and JSPS KAKENHI (17K06633, 21K04351).

## References

  * Anjos and Lasserre [2012] M. F. Anjos, J. B. Lasserre (eds.): Handbook on Semidefinite, Conic and Polynomial Optimization. Springer, New York (2012).
  * Aoues and Chateauneuf [2010] Y. Aoues, A. Chateauneuf: Benchmark study of numerical methods for reliability-based design optimization. Structural and Multidisciplinary Optimization, 41, 277–294 (2010).
  * Ben-Haim [2006] Y. Ben-Haim: Information-gap Decision Theory: Decisions under Severe Uncertainty (2nd ed.). Academic Press, London (2006).
  * Ben-Tal et al. [2009] A. Ben-Tal, L. El Ghaoui, A. Nemirovski: Robust Optimization. Princeton University Press, Princeton (2009).
  * Ben-Tal and Nemirovski [1997] A. Ben-Tal, A. Nemirovski: Robust truss topology optimization via semidefinite programming. SIAM Journal on Optimization, 7, 991–1016 (1997).
  * Beyer and Sendhoff [2007] H.-G. Beyer, B. Sendhoff: Robust optimization—a comprehensive survey. Computer Methods in Applied Mechanics and Engineering, 196, 3190–3218 (2007).
  * Boyd and Vandenberghe [2004] S. Boyd, L. Vandenberghe: Convex Optimization. Cambridge University Press, Cambridge (2004).
  * Calafiore and El Ghaoui [2014] G. C. Calafiore, L. El Ghaoui: Optimization Models. Cambridge University Press, Cambridge (2014).
  * Cho et al. [2016] H. Cho, K. K. Choi, N. J. Gaul, I. Lee, D. Lamb, D. Gorsich: Conservative reliability-based design optimization method with insufficient input data. Structural and Multidisciplinary Optimization, 54, 1609–1630 (2016).
  * Choi et al. [2010] J. Choi, D. An, J. Won: Bayesian approach for structural reliability analysis and optimization using the Kriging dimension reduction method. Journal of Mechanical Design, 132, 051003 (2010).
  * Delage and Ye [2010] E. Delage, Y. Ye: Distributionally robust optimization under moment uncertainty with application to data-driven problems. Operations Research, 58, 595–612 (2010).
  * El Ghaoui et al. [2003] L. El Ghaoui, M. Oks, F. Oustry: Worst-case value-at-risk and robust portfolio optimization: a conic programming approach. Operations Research, 51, 543–556 (2003).
  * Grant and Boyd [2008] M. Grant, S. Boyd: Graph implementations for nonsmooth convex programs. In: V. Blondel, S. Boyd, H. Kimura (eds.), Recent Advances in Learning and Control (A Tribute to M. Vidyasagar), Springer, pp. 95–110 (2008).
  * Grant and Boyd [2021] M. Grant, S. Boyd: CVX: Matlab Software for Disciplined Convex Programming. http://cvxr.com/cvx/ (Accessed April 2021).
  * Goh and Sim [2010] J. Goh, M. Sim: Distributionally robust optimization and its tractable approximations. Operations Research, 58, 902–917 (2010).
  * Gunawan and Papalambros [2006] S. Gunawan, P. Y. Papalambros: A Bayesian approach to reliability-based optimization with incomplete information. Journal of Mechanical Design, 128, 909–918 (2006).
  * Guo et al. [2009] X. Guo, W. Bai, W. Zhang, X. Gao: Confidence structural robust design and optimization under stiffness and load uncertainties. Computer Methods in Applied Mechanics and Engineering, 198, 3378–3399 (2009).
  * Guo et al. [2011] X. Guo, J. Du, X. Gao: Confidence structural robust optimization by non-linear semidefinite programming-based single-level formulation. International Journal for Numerical Methods in Engineering, 86, 953–974 (2011).
  * Holmberg et al. [2015] E. Holmberg, C.-J. Thore, A. Klarbring: Worst-case topology optimization of self-weight loaded structures using semi-definite programming. Structural and Multidisciplinary Optimization, 52, 915–928 (2015).
  * Huan et al. [2019] Z. Huan, G. Zhenghong, X. Fang, Z. Yidian: Review of robust aerodynamic design optimization for air vehicles. Archives of Computational Methods in Engineering, 26, 685–732 (2019).
  * Ito et al. [2018] M. Ito, N. H. Kim, N. Kogiso: Conservative reliability index for epistemic uncertainty in reliability-based design optimization. Structural and Multidisciplinary Optimization, 57, 1919–1935 (2018).
  * Ito and Kogiso [2016] M. Ito, N. Kogiso: Information uncertainty evaluated by parameter estimation and its effect on reliability-based multiobjective optimization. Journal of Advanced Mechanical Design, Systems, and Manufacturing, 10, 16-00331 (2016).
  * Jekel and Haftka [2020] C. F. Jekel, R. T. Haftka: Risk allocation for design optimization with unidentified statistical distributions. AIAA Scitech 2020 Forum, Orlando (2020).
  * Jiang et al. [2013] Z. Jiang, W. Chen, Y. Fu, R.-J. Yang: Reliability-based design optimization with model bias and data uncertainty. SAE International Journal of Materials and Manufacturing, 6, 502–516 (2013).
  * Jung et al. [2019] Y. Jung, H. Cho, I. Lee: Reliability measure approach for confidence-based design optimization under insufficient input data. Structural and Multidisciplinary Optimization, 60, 1967–1982 (2019).
  * Kang and Zhang [2020] Z. Kang, W. Zhang: Construction and application of an ellipsoidal convex model using a semi-definite programming formulation from measured data. Computer Methods in Applied Mechanics and Engineering, 300, 461–489 (2016).
  * Kanno [2011] Y. Kanno: Nonsmooth Mechanics and Convex Optimization. CRC Press, Boca Raton (2011).
  * Kanno [2018] Y. Kanno: Robust truss topology optimization via semidefinite programming with complementarity constraints: a difference-of-convex programming approach. Computational Optimization and Applications, 71, 403–433 (2018).
  * Kanno [2019] Y. Kanno: A data-driven approach to non-parametric reliability-based design optimization of structures with uncertain load. Structural and Multidisciplinary Optimization, 60, 83–97 (2019).
  * Kanno [2020a] Y. Kanno: Dimensionality reduction enhances data-driven reliability-based design optimizer. Journal of Advanced Mechanical Design, Systems, and Manufacturing, 14, 19-00200 (2020a).
  * Kanno [2020b] Y. Kanno: On three concepts in robust design optimization: absolute robustness, relative robustness, and less variance. Structural and Multidisciplinary Optimization, 62, 979–1000 (2020b).
  * Kanno and Takewaki [2006] Y. Kanno, I. Takewaki: Sequential semidefinite program for robust truss optimization based on robustness functions associated with stress constraints. Journal of Optimization Theory and Applications, 130, 265–287 (2006).
  * Keshtegar and Lee [2016] B. Keshtegar, I. Lee: Relaxed performance measure approach for reliability-based design optimization. Structural and Multidisciplinary Optimization, 54, 1439–1454 (2016).
  * Lee et al. [2010] I. Lee, K. K. Choi, D. Gorsich: Sensitivity analyses of FORM-based and DRM-based performance measure approach (PMA) for reliability-based design optimization (RBDO). International Journal for Numerical Methods in Engineering, 82, 26–46 (2010).
  * Moon et al. [2018] M.-Y. Moon, H. Cho, K. K. Choi, N. Gaul, D. Lamb, D. Gorsich: Confidence-based reliability assessment considering limited numbers of both input and output test data. Structural and Multidisciplinary Optimization, 57, 2027–2043 (2018).
  * Moustapha and Sudret [2019] M. Moustapha, B. Sudret: Surrogate-assisted reliability-based design optimization: a survey and a unified modular framework. Structural and Multidisciplinary Optimization, 60, 2157–2176 (2019).
  * Noh et al. [2011a] Y. Noh, K. K. Choi, I. Lee, D. Gorsich, D. Lamb: Reliability-based design optimization with confidence level under input model uncertainty due to limited test data. Structural and Multidisciplinary Optimization, 43, 443–458 (2011a).
  * Noh et al. [2011b] Y. Noh, K. K. Choi, I. Lee, D. Gorsich, D. Lamb: Reliability-based design optimization with confidence level for non-Gaussian distributions using bootstrap method. Journal of Mechanical Design, 133, 091001 (2011b).
  * Oberkampf et al. [2004] W. L. Oberkampf, J. C. Helton, C. A. Joslyn, S. F. Wojtkiewicz, S. Ferson: Challenge problems: uncertainty in system response given uncertain parameters. Reliability Engineering and System Safety, 85, 11–19 (2004).
  * Pólik [2005] I. Pólik: Addendum to the SeDuMi User Guide: Version 1.1. Technical Report, Advanced Optimization Laboratory, McMaster University, Hamilton (2005). http://sedumi.ie.lehigh.edu/sedumi/ (Accessed April 2021).
  * Schöbi and Sudret [2017] R. Schöbi, B. Sudret: Structural reliability analysis for p-boxes using multi-level meta-models. Probabilistic Engineering Mechanics, 48, 27–38 (2017).
  * Sturm [1999] J. F. Sturm: Using SeDuMi 1.02, a MATLAB toolbox for optimization over symmetric cones. Optimization Methods and Software, 11–12, 625–653 (1999).
  * Takezawa et al. [2011] A. Takezawa, S. Nii, M. Kitamura, N. Kogiso, Topology optimization for worst load conditions based on the eigenvalue analysis of an aggregated linear system, Computer Methods in Applied Mechanics and Engineering, 200, 2268–2281 (2011).
  * Thore et al. [2017] C.-J. Thore, E. Holmberg, A. Klarbring: A general framework for robust topology optimization under load-uncertainty including stress constraints. Computer Methods in Applied Mechanics and Engineering, 319, 1–18 (2017).
  * Valdebenito and Schuëller [2010] M. A. Valdebenito, G. I. Schuëller: A survey on approaches for reliability-based optimization. Structural and Multidisciplinary Optimization, 42, 645–663 (2010).
  * Wang et al. [2020] Y. Wang, P. Hao, H. Yang, B. Wang, Q. Gao: A confidence-based reliability optimization with single loop strategy and second-order reliability method. Computer Methods in Applied Mechanics and Engineering, 372, 113436 (2020).
  * Wiesemann et al. [2014] W. Wiesemann, D. Kuhn, M. Sim: Distributionally robust convex optimization. Operations Research, 62, 1358–1376 (2014).
  * Yamashita and Yabe [2015] H. Yamashita, H. Yabe: A survey of numerical methods for nonlinear semidefinite programming. Journal of the Operations Research Society of Japan, 58, 24–60 (2015).
  * Yao et al. [2011] W. Yao, X. Chen, W. Luo, M. van Tooren, J. Guo: Review of uncertainty-based multidisciplinary design optimization methods for aerospace vehicles. Progress in Aerospace Sciences, 47, 450–479 (2011).
  * Youn and Wang [2008] B. D. Youn, P. Wang: Bayesian reliability-based design optimization using eigenvector dimension reduction (EDR) method. Structural and Multidisciplinary Optimization, 36, 107–123 (2008).
  * Zaman and Mahadevan [2017] K. Zaman, S. Mahadevan: Reliability-based design optimization of multidisciplinary system under aleatory and epistemic uncertainty. Structural and Multidisciplinary Optimization, 55, 681–699 (2017).
  * Zaman et al. [2011] K. Zaman, S. Rangavajhala, M. P. McDonald, S. Mahadevan: A probabilistic approach for representation of interval uncertainty. Reliability Engineering and System Safety, 96, 117–130 (2011).
  * Zhang et al. [2020] J. Zhang, L. Gao, M. Xiao, S. Lee, A. T. Eshghi: An active learning Kriging-assisted method for reliability-based design optimization under distributional probability-box model. Structural and Multidisciplinary Optimization, 62, 2341–2356 (2020).