# Risk-averse estimates of effective properties in heterogeneous elasticity

Jeremy Bleyer

<figure><img src="image_1.png" /></figure>

## To cite this version:

Jeremy Bleyer. Risk-averse estimates of effective properties in heterogeneous elasticity. Comptes Rendus. Mécanique, 2023, 351 (G1), pp.29-42. 10.5802/crmeca.171. hal-03964172

HAL Id: hal-03964172

https://enpc.hal.science/hal-03964172

Submitted on 30 Jan 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Comptes Rendus
## Mécanique

Jeremy Bleyer

Risk-averse estimates of effective properties in heterogeneous elasticity

Volume 351 (2023), p. 29-42

Published online: 24 January 2023

https://doi.org/10.5802/crmeca.171

<figure><img src="image_3.png" /></figure>

This article is licensed under the
CREATIVE COMMONS ATTRIBUTION 4.0 INTERNATIONAL LICENSE.
http://creativecommons.org/licenses/by/4.0/

<figure><img src="image_4.png" /></figure>

Les Comptes Rendus. Mécanique sont membres du
Centre Mersenne pour l'édition scientifique ouverte
www.centre-mersenne.org

e-ISSN : 1873-7234
# 1. Introduction

When considering the behavior of real materials, various factors can influence the resulting mechanical properties (e.g. defects, inclusions, mechanical properties of the constitutive phases, geometry of the microstructure, etc.). Many of such factors cannot be predicted exactly and are often considered as random. Even in the elastic regime, the mechanical properties of the material therefore exhibit stochastic fluctuations. A typical example is concrete which is an inherently multiscale material in which uncertainties on constitutive properties occur at various scales. For instance, at the cement paste level, uncertainties on the hydrate mechanical properties can be magnified during the upscaling process towards the macroscopic elastic properties [1]. At the mortar scale, the Interfacial Transition Zone (ITZ) is also known to be a major parameter impacting the effective properties whereas its constitutive properties are very difficult to determine precisely [2]. In practice, one often attempts at measuring experimentally or estimating numerically the average behavior of such stochastic properties. Developments in uncertainty quantification techniques now aim at also including additional statistical information such as variance or higher-order moments [3,4]. However, including such additional information in a global structural analysis (e.g. via stochastic finite-element methods [5]) remains challenging. Moreover, the generation of prior stochastic models of random elastic tensors is not a trivial task, especially when considering anisotropic symmetry classes [6, 7]. Finally, even if the variance of quantities of interest is valuable, engineering design often requires to estimate accurately worst-case scenarios with low probabilities of occurrence. In this context, information about the distribution tail of quantities of interest is required which is even more challenging to obtain and include in a global analysis. Various reliability analysis approaches (e.g. FORM/SORM [8]) have then been proposed to assess the risk of failure of a structure. Without being exhaustive, other illustrations of uncertainty propagation techniques in the context of heterogeneous materials can be found for instance in [9–13].

In the present work, we propose a novel approach to estimate pessimistic (and also optimistic) values of stochastic elastic material properties by exploiting the concept of *risk measures* which is widely used in financial mathematics, see [14] for a broad overview. More precisely, attention has been turned in recent years towards convex coherent risk measures which benefit from interesting mathematical properties, see [15]. A risk measure assigns to a given random variable a single value representing an estimate of its typical value. For instance, one can think about the mean value, the median value, the minimum value, etc. Some risk measures are obviously more useful than others. For instance, the mean value is not particularly useful if one looks for a safe design in engineering applications. Designing with respect to the minimum value of a mechanical property (such as elastic stiffness or strength for instance) would, on the contrary, be overly conservative. Additionally, the mean value $± k$ standard deviations is often used as a measure of data variability. However, it cannot be directly used as an accurate confidence interval measure unless the underlying distribution is precisely known. It is therefore not easy to use for describing low probability events without any additional information. A very popular risk measure which succeeds in describing the typical value in the distribution tail is the so-called Conditional Value-at-Risk (CVaR) which can be defined as the expected value of the $α$-quantile of the random variable, see Section 3.1. For $α$ going from 0 to 1, the CVaR evolves from the mean-value to the supremum value of the random variable. If high values of the latter are undesirable, we can then say that the CVaR evolves from a *risk-neutral* estimate (for $α = 0$) to a *risk-averse* estimate (for $α$ close to 1, typically 0.95, 0.99, etc.). The risk-averse estimate therefore corresponds to the expected value of the largest values of the random variable which occur with a probability $1 - α$. In our opinion, such a definition corresponds to a reasonable worst-case estimate without being overly conservative for engineering applications.
Our main contribution is to make use of the CVaR on the elastic free and complementary energies to produce optimistic and pessimistic estimates of elastic properties. Obviously, the latter seem more relevant from the engineering point of view than the latter. To do so, we will rely on a key aspect which is that the CVaR is a convex risk measure. This will ensure thermodynamic stability of the resulting risk-averse elastic behavior. To make our developments more explicit, we will formulate them in the context of homogenization by considering the effective macroscopic behavior of an elastic heterogeneous RVE depending on random microscopic properties (e.g. the elastic properties of its constituents). Nevertheless, the proposed formulations can be readily applied to other situations involving an elastic behavior.

The manuscript will be organized as follows: Section 2 recalls notations and generic homogenization results in random heterogeneous elasticity; Section 3 introduces the CVaR optimistic and pessimistic estimates; Section 4 then applies the proposed estimates to the case of isotropic composites; Section 5 discusses a linear elastic approximation of the proposed risk-averse effective behaviors; finally, Section 6 draws some conclusions and research perspectives.

## 2. Random homogenization setting

In the remainder, we will consider the homogenization setting of a heterogeneous elastic material. We will assume that a RVE $\Omega$ exists and that effective elastic properties $C^{\text{hom}}$ can be computed for given microscopic material properties, assuming in particular separation of scales. We will consider that the latter are stochastic and we are interested in the corresponding stochastic estimate of the macroscopic elastic properties. In particular, we do not investigate the limit of infinitely large RVE domain for which effective properties become deterministic but rather focus on a fixed RVE size. The RVE is subjected to a macroscopic strain $E$ and denote by $\boldsymbol{\Sigma}$ the corresponding macroscopic stress. The effective elastic properties $C^{\text{hom}}$ related to the elastic total strain energy $\Psi(E)$ can be obtained from the solution of a microscopic elasticity problem defined on the RVE since:

$$ \Psi(E) = \frac{1}{2} E : C^{\text{hom}} : E = \min_{\epsilon \in \mathcal{K} \Lambda(E)} \frac{1}{|\Omega|} \int_{\Omega} \psi(\epsilon; y) \, d\Omega \quad (1) $$

where $\psi(\epsilon; y)$ is the local elastic strain energy at a given point $y \in \Omega$ and $\epsilon$ denotes any kinematically admissible (KA) strain field with the macroscopic strain $E$.

Similarly, we have the corresponding characterization of the effective compliance $\mathbb{S}^{\text{hom}} = (C^{\text{hom}})^{-1}$ using the complementary elastic energy:

$$ \Psi^*(\boldsymbol{\Sigma}) = \frac{1}{2} \boldsymbol{\Sigma} : \mathbb{S}^{\text{hom}} : \boldsymbol{\Sigma} = \min_{\boldsymbol{\sigma} \in \mathcal{S}(\boldsymbol{\Sigma})} \frac{1}{|\Omega|} \int_{\Omega} \psi^*(\boldsymbol{\sigma}; y) \, d\Omega \quad (2) $$

where $\psi^*(\boldsymbol{\sigma}; y)$ is the local elastic stress energy for any statically admissible (SA) stress field $\boldsymbol{\sigma}$ with the macroscopic stress $\boldsymbol{\Sigma}$.

For a given macroscopic strain $E$, the corresponding macroscopic stress can therefore be obtained as follows:

$$ \boldsymbol{\Sigma} = \frac{\partial \Psi}{\partial E}(E) = \arg \min_{\hat{\boldsymbol{\Sigma}}} \Psi^*(\hat{\boldsymbol{\Sigma}}) - \hat{\boldsymbol{\Sigma}} : E \quad (3) $$

Let us now assume that the constitutive materials possess random properties characterized by a random variable $\zeta$ with a known probability distribution. The macroscopic strain and stress energies, and thus the corresponding macroscopic stiffness and compliance, are now random variables i.e.:

$$ \Psi(E; \zeta) = \frac{1}{2} E : C^{\text{hom}}(\zeta) : E \quad (4) $$

$$ \Psi^*(\boldsymbol{\Sigma}; \zeta) = \frac{1}{2} \boldsymbol{\Sigma} : \mathbb{S}^{\text{hom}}(\zeta) : \boldsymbol{\Sigma} \quad (5) $$
<figure><img src="image_5.png" /><figcaption>**Figure 1.** Illustration of the CVaR definition.</figcaption></figure>

It is then customary to consider the ensemble average of the above relations to compute effective macroscopic properties:

$$ \Psi_{\text{eff}}(E) = \mathbb{E}[\Psi(E; \zeta)] = \frac{1}{2} E : \mathbb{E}[C^{\text{hom}}(\zeta)] : E \quad (6) $$

$$ \Psi_{\text{eff}}^{*}(E) = \mathbb{E}[\Psi^{*}(\Sigma; \zeta)] = \frac{1}{2} \Sigma : \mathbb{E}[S^{\text{hom}}(\zeta)] : \Sigma \quad (7) $$

Note that one classically defines the corresponding effective elastic stiffness as follows:

$$ C_{\text{eff}}^{\text{hom}} = \mathbb{E}[C^{\text{hom}}(\zeta)] \quad (8) $$

so that $\Psi_{\text{eff}}(E) = (1/2)E : C_{\text{eff}}^{\text{hom}} : E$.

Our aim is now to provide some risk-averse estimates of the macroscopic stiffness and compliance.

# 3. CVaR estimates

## 3.1. The conditional value-at-risk

To introduce the CVaR, we first define the Value-at-Risk (VaR$_{α}$) of level $α \in [0; 1)$ (typically close to 1) of a random variable X as follows:

$$ \text{VaR}_{\alpha}(X) = \inf\{Z \text{ s.t. } F_{X}(Z) \geq \alpha\} \quad (9) $$

where $F_X(Z) = \mathbb{P}(Z \leq X)$ is the cumulative distribution function of X.

For a continuous distribution, the CVaR is defined as the expected value of X above VaR$_{α}$ (see Figure 1)

$$ \text{CVaR}_{\alpha}(X) = \mathbb{E}[X \text{ s.t. } X \geq \text{VaR}_{\alpha}(X)] \quad (10) $$

The above definition is slightly more technical when considering discrete distributions, we refer the reader to [16] for a rigorous definition.

A key result due to [17] is that CVaR benefits from the following convex optimization characterization:

$$ \text{CVaR}_{\alpha}(X) = \inf_{\lambda} \lambda + \frac{1}{1-\alpha} \mathbb{E}[\langle X - \lambda \rangle_{+}] \quad (11) $$

where $\langle \star \rangle_{+} = \max\{\star, 0\}$ denotes the positive part. From this definition, it is clear that $\text{CVaR}_0(X) = \mathbb{E}[X]$. Let us also point out that, if the minimum of (11) is unique, then the optimal value is exactly $\lambda = \text{VaR}_\alpha(X)$.

The CVaR$_{α}$ of a vector x in $\mathbb{R}^N$ can also be interpreted as a norm parameterized by $α$. Indeed, for $α = 0$, it reduces to the $L_1$-norm scaled by a factor $1/N$, whereas for $α = 1$ it reduces to the
$L_\infty$-norm. For intermediate values, it can be seen as an average of the $k$ largest values of $|x|$ where $k/N = 1 - \alpha$, see [18] for a more precise definition.

## 3.2. Optimistic estimate

We now introduce the notation $\Psi_{\bar{\alpha}}(E)$ to represent, for a given $E$, the CVaR of the stochastic macroscopic strain energy $\Psi(E; \zeta)$ for some confidence level $\alpha$:

$$ \Psi_{\bar{\alpha}}(E) = \text{CVaR}_{\alpha}(\Psi(E; \zeta)) \quad (12) $$

Owing to the convex representation formula (11), $\Psi_{\bar{\alpha}}$ is still a convex function of $E$. However, it is no longer quadratic but only piecewise-quadratic. It remains however quadratic for radial paths $E(\eta) = \eta E^0$ since $\Psi_{\bar{\alpha}}(\eta E^0) = \eta^2 \Psi_{\bar{\alpha}}(E^0)$. Note that:

$$ \Psi_{\bar{0}}(E) = \mathbb{E}[\Psi(E; \zeta)] \quad (13) $$

$$ \Psi_{\bar{\alpha}}(E) \le \Psi_{\bar{\beta}}(E) \quad \forall \alpha \le \beta \quad (14) $$

As a result, for a confidence level $\alpha$ close to 1, $\Psi_{\bar{\alpha}}(E)$ will provide an optimistic estimate for the random macroscopic strain energy $\Psi(E; \zeta)$, resulting in an optimistic estimate of the macroscopic stiffness.

Since $\Psi_{\bar{\alpha}}(E)$ is convex, we can define the corresponding risk-averse macroscopic stress for a given strain by generalizing (3) as follows:

$$ \Sigma_{\bar{\alpha}} = \frac{\partial \Psi_{\bar{\alpha}}}{\partial E}(E) = \argmin_{\hat{\Sigma}} (\Psi_{\bar{\alpha}})^*(\hat{\Sigma}) - \hat{\Sigma} : E \quad (15) $$

Again, owing to the piecewise-quadratic nature of $\Psi_{\bar{\alpha}}$, the resulting macroscopic stress-strain relationship is piecewise-linear elastic, except for radial loading paths on which it remains linear elastic.

## 3.3. Pessimistic estimate

Conversely, we can define the CVaR $(\Psi^*)_{\bar{\alpha}}$ of the stochastic macroscopic complementary energy $\Psi^*(\Sigma; \zeta)$:

$$ (\Psi^*)_{\bar{\alpha}}(\Sigma) = \text{CVaR}_{\alpha}(\Psi^*(\Sigma; \zeta)) = \mathbb{E}[\Psi^*|\Psi^* \ge \text{VaR}_{\alpha}(\Psi^*)] \quad (16) $$

Similarly to $\Psi_{\bar{\alpha}}(E)$, $(\Psi^*)_{\bar{\alpha}}(\Sigma)$ is a convex function of $\Sigma$ and provides a risk-averse estimate of the random complementary energy $\Psi^*(\Sigma; \zeta)$ for $\alpha$ close to 1. $(\Psi^*)_{\bar{\alpha}}(\Sigma)$ is also piecewise-quadratic and quadratic along radial paths and the resulting macroscopic strain-stress relationship is piecewise-linear elastic.

Considering a Legendre-Fenchel transformation, we can compute its corresponding conjugate function $((\Psi^*)_{\bar{\alpha}})^*$. Let us therefore introduce:

$$ \Psi_{\underline{\alpha}} = ((\Psi^*)_{\bar{\alpha}})^* \quad (17) $$

Using (47) for $\Psi^*$, one obtains the following convex optimization formulation for $\Psi_{\underline{\alpha}}$:

$$ \begin{array}{l} \Psi_{\underline{\alpha}}(E) = ((\Psi^*)_{\bar{\alpha}})^*(E) = \inf_{\hat{E}, z} \mathbb{E} \left[ \frac{1}{z} \Psi(\hat{E}) \right] \\ \text{s.t. } \mathbb{E}[\hat{E}] = E \\ 0 \le z \le \frac{1}{1-\alpha} \\ \mathbb{E}[z] = 1 \end{array} \quad (18) $$
which corresponds to some weighted average of $\Psi$. We also have:

$$ \Psi_0(E) = (\mathbb{E}[\Psi^*])^*(E) \quad (19) $$

$$ \Psi_{\underline{\alpha}}(E) \geq \Psi_{\underline{\beta}}(E) \quad \forall \alpha \leq \beta \quad (20) $$

As a result, for a confidence level $\alpha$ close to 1, $\Psi_{\underline{\alpha}}(E)$ provides a pessimistic estimate of the random strain energy, resulting in a pessimistic estimate of the macroscopic stiffness.

# 4. Application to isotropic composites

## 4.1. Bulk and shear moduli risk-averse estimates

Let us consider the case of an isotropic composite material described by apparent bulk and shear moduli $\kappa(\zeta)$ and $\mu(\zeta)$. One has:

$$ \Psi(E; \zeta) = \frac{1}{2}(\kappa(\zeta)(\text{tr}E)^2 + 2\mu(\zeta)E^d : E^d) \quad (21) $$

$$ \Psi^*(\Sigma; \zeta) = \frac{1}{2} \left( \frac{1}{9\kappa(\zeta)} (\text{tr}\Sigma)^2 + \frac{1}{2\mu(\zeta)} \Sigma^d : \Sigma^d \right) \quad (22) $$

where $E^d = \text{dev}(E)$ and $\Sigma^d = \text{dev}(\Sigma)$.

For purely deviatoric strain states i.e. $\text{tr}(E) = 0$, one has $\Psi(E; \zeta) = \mu(\zeta)E^d : E^d$ so that:

$$ \begin{aligned} \Psi_{\bar{\alpha}}(E) &= \inf_{\lambda} \lambda + \frac{1}{1-\alpha} \langle \mu(\zeta) E^d : E^d - \lambda \rangle_+ \\ &= \left( \inf_{\hat{\lambda}} \lambda + \frac{1}{1-\alpha} \langle \mu(\zeta) - \hat{\lambda} \rangle_+ \right) E^d : E^d \\ &= \text{CVaR}_{\alpha}(\mu) E^d : E^d \end{aligned} \quad (23) $$

As a consequence, $\Psi_{\bar{\alpha}}$ is still a quadratic form on the space of purely deviatoric strains and corresponds to a linear elastic material of shear modulus $\mu_{\bar{\alpha}} = \text{CVaR}_{\alpha}(\mu)$.

A similar reasoning can be made for purely spherical strains. $\Psi_{\bar{\alpha}}$ is still a quadratic form on such a space and corresponds to a linear elastic material of bulk modulus $\kappa_{\bar{\alpha}} = \text{CVaR}_{\alpha}(\kappa)$.

Similarly, considering the case of purely deviatoric stress states i.e. $\text{tr}(\Sigma) = 0$, then one has $\Psi^*(\Sigma; \zeta) = (1/4\mu(\zeta))\Sigma^d : \Sigma^d$ so that:

$$ \begin{aligned} \Psi_{\bar{\alpha}}^*(\Sigma) &= \inf_{\lambda} \lambda + \frac{1}{1-\alpha} \left\langle \frac{1}{4\mu(\zeta)} \Sigma^d : \Sigma^d - \lambda \right\rangle_+ \\ &= \frac{1}{4} \text{CVaR}_{\alpha} \left( \frac{1}{\mu} \right) \Sigma^d : \Sigma^d \end{aligned} \quad (24) $$

As a consequence, $\Psi_{\bar{\alpha}}^*$ is quadratic on the space of purely deviatoric stresses and corresponds to a linear elastic material of shear modulus $\mu_{\bar{\alpha}} = \text{CVaR}_{\alpha}(\mu^{-1})^{-1}$. Similarly, for purely spherical stress states, one obtains a linear elastic material of bulk modulus $\kappa_{\bar{\alpha}} = \text{CVaR}_{\alpha}(\kappa^{-1})^{-1}$.

Finally, let us recall that the optimal value of $\lambda$ corresponds to the VaR. Then, $1/\mu^* = \text{VaR}_{\alpha}(\mu^{-1})$ is such that $F_{\mu^{-1}}(1/\mu^*) = \alpha$. But one also has $F_{1/\chi}(x) = 1 - F_{\chi}(1/x)$ so that:

$$ F_{\mu^{-1}}(1/\mu^*) = 1 - F_{\mu}(\mu^*) = \alpha \quad \Rightarrow \quad F_{\mu}(\mu^*) = 1 - \alpha \quad (25) $$

We conclude that $\mu^* = 1/\text{VaR}_{\alpha}(\mu^{-1})$ corresponds to the $(1-\alpha)$-quantile of $\mu(\zeta)$. However, $\mu_{\underline{\alpha}}$ does not correspond to the expectation below this $(1-\alpha)$-quantile but rather the “harmonic expectation” $\mathbb{E}[\mu^{-1}]^{-1}$ below this quantile.
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Empirical PDF and CDF (light gray) and computed CVaR and VaR optimistic and pessimistic estimates for $f = 0.25$.</figcaption></figure>

<table><caption>**Table 1.** Numerical values of risk-averse and risk-neutral homogenized shear modulus estimates</caption><thead><tr><th>CVaR<sub>α</sub>(μ<sup>-1</sup>)<sup>-1</sup></th><th>VaR<sub>α</sub>(μ<sup>-1</sup>)<sup>-1</sup></th><th>E[μ]</th><th>CVaR<sub>α</sub>(μ)</th><th>VaR<sub>α</sub>(μ)</th></tr></thead><tbody><tr><td>12.856</td><td>13.290</td><td>15.174</td><td>17.211</td><td>17.810</td></tr></tbody></table>

## 4.2. *Numerical example*

As an illustrative example, let us consider an isotropic biphasic material consisting of two incompressible phases of shear moduli $\mu_1$ and $\mu_2$ and where $f$ is the volume fraction of phase 2. We estimate the effective modulus by resorting to a self-consistent scheme [19]. Both phases shear moduli are considered as random variables following a lognormal distribution of mean value $\mu_1^0 = 10$ and $\mu_2^0 = 50$ and a standard deviation of 10% for both phases. Note that we choose lognormal distributions for the shear modulus for illustration purposes only. More accurate choices on probability distribution of elastic behaviors can be found in the work of [20] for instance. We consider a set of $N = 10,000$ Monte-Carlo realizations to estimate the different expectations and CVaR. We choose a confidence level of $\alpha = 0.95$ unless stated otherwise.

For $f = 0.25$, Figures 2a and b respectively display the empirical PDF and CDF of the homogenized shear moduli estimate. The computed expected value and the various risk-averse estimates are reported in Table 1. One can see that both pessimistic and optimistic VaR levels indeed correspond to the values for which the empirical CDF reaches a value of $1 - \alpha$ and $\alpha$ respectively (thin horizontal dashed lines in Figure 2b) and that the CVaR corresponds to an average value of the remaining part (below or above the VaR) of the distribution. Figure 3 also reports the evolution of both pessimistic $\mu_\alpha$ and optimistic $\mu_\bar{\alpha}$ shear modulus estimates, normalized by the effective modulus $\mu^{\text{eff}} = E[\mu]$ as a function of the chosen confidence level $\alpha$. One can see that the risk-averse estimates deviate slowly from the effective modulus for low values of $\alpha$ and more rapidly when $\alpha$ approaches 1, capturing the effect of the flatter tails of the distribution.

Finally, Figure 4 also reports the evolution of the various estimates as a function of the volume fraction $f$. For each volume fraction, the dotted symbols represent a discrete value of the empirical distribution, the size and the color being proportional to the corresponding probability of occurrence. We can clearly see that the proposed risk-averse estimates provide a
<figure><img src="image_7.png" /><figcaption>**Figure 3.** Risk-averse shear modulus estimates as a function of the confidence level $\alpha$ for $f = 0.25$.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 4.** Risk-averse shear modulus estimates as a function of the volume fraction.</figcaption></figure>

nice bracketing of the random modulus around the central part of the distribution for the whole range of volume fractions.

## 5. A linear elastic approximation

As demonstrated previously, the CVaR estimates are an efficient risk-averse measure of the elastic properties of a random material. However, considering the CVaR on the total free or complementary energy does not yield a linear elastic material anymore but a piecewise-linear elastic which may be more cumbersome to use in practice. Nonetheless, we have seen that the corresponding behavior is still linear elastic for purely deviatoric and spherical states. It therefore appears natural to consider, as an optimistic (resp. pessimistic) estimate of the random material, an isotropic linear elastic material characterized by the bulk and shear moduli $\kappa_{\bar{\alpha}}, \mu_{\bar{\alpha}}$ (resp. $\kappa_{\underline{\alpha}}, \mu_{\underline{\alpha}}$).

Let us now give a more formal definition of such a proposition. We consider the spherical and deviatoric decomposition of an isotropic linear elastic material free energy:

$$ \Psi(E; \zeta) = \Psi^{\text{sph}}(E; \zeta) + \Psi^{\text{dev}}(E; \zeta) = \frac{1}{2} \kappa(\zeta) (\text{tr} E)^2 + \mu(\zeta) E^d : E^d \quad (26) $$

$$ \Psi^*(\Sigma; \zeta) = (\Psi^{\text{sph}})^*(\Sigma; \zeta) + (\Psi^{\text{dev}})^*(\Sigma; \zeta) = \frac{1}{18\kappa(\zeta)} (\text{tr}\Sigma)^2 + \frac{1}{4\mu(\zeta)} \Sigma^d : \Sigma^d \quad (27) $$
Owing to the sub-additivity property of the CVaR, we have that:

$$ \Psi_{\bar{\alpha}}(E) \le \Psi_{\bar{\alpha}}^{\text{sph}}(E) + \Psi_{\bar{\alpha}}^{\text{dev}}(E) \quad \forall E \quad (28) $$

$$ (\Psi^*)_{\bar{\alpha}}(\Sigma) \le ((\Psi^{\text{sph}})^*)_{\bar{\alpha}}(\Sigma) + ((\Psi^{\text{dev}})^*)_{\bar{\alpha}}(\Sigma) \quad \forall \Sigma \quad (29) $$

where the right-hand sides are given by:

$$ \Psi_{\bar{\alpha}}^{\text{sph}}(E) = \frac{1}{2} \kappa_{\bar{\alpha}} (\text{tr} E)^2 \quad (30) $$

$$ \Psi_{\bar{\alpha}}^{\text{dev}}(E) = \mu_{\bar{\alpha}} E^d : E^d \quad (31) $$

$$ ((\Psi^{\text{sph}})^*)_{\bar{\alpha}}(\Sigma) = \frac{1}{18}(\kappa^{-1})_{\bar{\alpha}}(\text{tr}\,\Sigma)^2 = \frac{1}{18\kappa_{\bar{\alpha}}} (\text{tr}\,\Sigma)^2 \quad (32) $$

$$ ((\Psi^{\text{dev}})^*)_{\bar{\alpha}}(\Sigma) = \frac{1}{4}(\mu^{-1})_{\bar{\alpha}}\Sigma^d : \Sigma^d = \frac{1}{4\mu_{\bar{\alpha}}}\Sigma^d : \Sigma^d \quad (33) $$

Since the Legendre–Fenchel transform is order-reversing and since $(\Psi^{\text{sph}})^*$ (resp. $(\Psi^{\text{dev}})^*$) depends only on tr $\Sigma$ (resp. $\Sigma^d$), (27) is equivalent to:

$$ \Psi_{\bar{\alpha}} = ((\Psi^*)_{\bar{\alpha}})^*(E) \ge ((((\Psi^{\text{sph}})^*)_{\bar{\alpha}})^*(E) + (((\Psi^{\text{dev}})^*)_{\bar{\alpha}})^*(E)) \quad (34) $$

with:

$$ (((\Psi^{\text{sph}})^*)_{\bar{\alpha}})^*(E) = \Psi_{\bar{\alpha}}^{\text{sph}}(E) = \frac{1}{2} \kappa_{\bar{\alpha}} (\text{tr} E)^2 \quad (35) $$

$$ (((\Psi^{\text{dev}})^*)_{\bar{\alpha}})^*(E) = \Psi_{\bar{\alpha}}^{\text{dev}}(E) = \mu_{\bar{\alpha}} E^d : E^d \quad (36) $$

As a result, we have $\forall E$:

$$ \Psi_{\bar{\alpha}}^{\text{sph}}(E) + \Psi_{\bar{\alpha}}^{\text{dev}}(E) \le \Psi_{\bar{\alpha}}(E) \le \Psi^{\text{eff}}(E) \le \Psi_{\bar{\alpha}}(E) \le \Psi_{\bar{\alpha}}^{\text{sph}}(E) + \Psi_{\bar{\alpha}}^{\text{dev}}(E) \quad (37) $$

where the left-most estimate $\Psi_{\bar{\alpha}}^{\text{sph}}(E) + \Psi_{\bar{\alpha}}^{\text{dev}}(E)$ is a quadratic form corresponding indeed to an isotropic elastic material of moduli $(\kappa_{\bar{\alpha}}, \mu_{\bar{\alpha}})$ and the right-most estimate $\Psi_{\bar{\alpha}}^{\text{sph}}(E) + \Psi_{\bar{\alpha}}^{\text{dev}}(E)$ is a quadratic form corresponding indeed to an isotropic elastic material of moduli $(\kappa_{\bar{\alpha}}, \mu_{\bar{\alpha}})$.

## 5.1. Generalization to other symmetry classes

The previous approximate linear elastic formulation can also be easily generalized to other material symmetry classes by considering the Kelvin decomposition of the elastic moduli [21]. Let $C(\zeta) = \{c^{(i)}(\zeta)\}$ (resp. $\mathbb{S}(\zeta) = \{s^{(i)}(\zeta)\}$) be defined by its stochastic eigen-moduli $c^{(i)}(\zeta)$ (resp. $s^{(i)}(\zeta)$). The stochastic elastic free and complementary energies therefore write as:

$$ \Psi(E; \zeta) = \frac{1}{2} \sum_{i=1}^{K} c^{(i)}(\zeta) E^{(i)} : E^{(i)} \quad (38) $$

$$ \Psi^*(\Sigma; \zeta) = \frac{1}{2} \sum_{i=1}^{K} s^{(i)}(\zeta) \Sigma^{(i)} : \Sigma^{(i)} \quad (39) $$

where $E^{(i)} = \mathbb{P}^{(i)} : E$ and $\Sigma^{(i)} = \mathbb{P}^{(i)} : \Sigma$ where $\mathbb{P}^{(i)}$ denotes the projector onto mode number (i) in the Kelvin decomposition.

Note that in (38)–(39), we assume that the uncertainty $\zeta$ does not affect the symmetry class<fnref n="1" /> and the corresponding projectors $\mathbb{P}^{(i)}$. The optimistic elastic material can then be defined as $C_{\bar{\alpha}} = \{c_{\bar{\alpha}}^{(i)}\}$. Similarly, the pessimistic elastic material can be defined from $\mathbb{S}_{\bar{\alpha}} = \{s_{\bar{\alpha}}^{(i)}\}$ and thus $C_{\alpha} = (\mathbb{S}_{\bar{\alpha}})^{-1} = \{c_{\alpha}^{(i)}\}$, thereby generalizing the previously discussed approximation.

<footnote n="1">¹In the case when this is no longer true, one can still use the generic CVaR estimates of Section 3.</footnote>
<table><caption>**Table 2.** Numerical parameters of Case 1</caption><thead><tr><th>Property</th><th>Distribution</th><th>Numerical values</th></tr></thead><tbody><tr><td>E<sub>1</sub></td><td>Lognormal</td><td>Mean 10, std 20%</td></tr><tr><td>E<sub>2</sub></td><td>Lognormal</td><td>Mean 50, std 20%</td></tr><tr><td>ν<sub>1</sub></td><td>Uniform</td><td>[0.15, 0.25]</td></tr><tr><td>ν<sub>2</sub></td><td>Uniform</td><td>[0.25, 0.35]</td></tr></tbody></table>

<table><caption>**Table 3.** Numerical parameters of Case 2</caption><thead><tr><th>Property</th><th>Distribution</th><th>Numerical values</th></tr></thead><tbody><tr><td>κ<sub>1</sub></td><td>Lognormal</td><td>Mean 5.56, std 20%</td></tr><tr><td>κ<sub>2</sub></td><td>Lognormal</td><td>Mean 41.67, std 20%</td></tr><tr><td>μ<sub>1</sub></td><td>Lognormal</td><td>Mean 4.17, std 20%</td></tr><tr><td>μ<sub>2</sub></td><td>Lognormal</td><td>Mean 19.23, std 20%</td></tr></tbody></table>

<table><caption>**Table 4.** Numerical parameters of Case 3</caption><thead><tr><th>Property</th><th>Distribution</th><th>Numerical values</th></tr></thead><tbody><tr><td>E<sub>1</sub></td><td>Lognormal</td><td>Mean 10, std 20%</td></tr><tr><td>E<sub>2</sub></td><td>Lognormal</td><td>Mean 50, std 20%</td></tr><tr><td>ν<sub>1</sub></td><td>Uniform</td><td>[-0.45, 0.45]</td></tr><tr><td>ν<sub>2</sub></td><td>Uniform</td><td>[-0.45, 0.45]</td></tr></tbody></table>

## 5.2. *Numerical example*

We now consider an isotropic biphasic composite consisting of two linear elastic materials of Young modulus and Poisson ratio $E_i$, $v_i$ for $i = 1, 2$. $f$ still denotes the volume fraction of phase 2 and is taken as $f = 0.25$ in the following. We again estimate the effective moduli using a self-consistent scheme for the various cases of Tables 2–4.

We first consider Case 1 of Table 2 where both phases Young modulus to be random variables following a lognormal distribution of mean value $E_1^0 = 10$ and $E_2^0 = 50$ and a standard deviation of 20% for both phases. The Poisson ratios are assumed to be random variables with a uniform distribution centered around $v_1^0 = 0.2$ and $v_2^0 = 0.3$ with a $\pm 0.05$ maximum deviation. In such a context, the empirical PDF of the effective bulk and shear moduli of the composite are represented in Figure 5a. The corresponding pessimistic and optimistic CVaR estimates for a 0.95-confidence level have also been reported on this figure. In order to assess the accuracy of the linear elastic approximation which has just been proposed, we represent, for a given free energy density $\Psi$, the level-set $2\Psi^*(\Sigma) = 1$. Results are represented in the normalized stress space $\bar{\Sigma}_m = (1/3 \text{tr} \Sigma) / \sqrt{\kappa^{\text{eff}}}$ and $\bar{\Sigma}_d = \|\text{dev}(\Sigma)\| / \sqrt{2\mu^{\text{eff}}}$ so that the level set corresponds to a unit circle for $\Psi^{\text{eff}} = E[\Psi]$. Linear elastic materials possess level sets which form an ellipse in such a space. Clearly from Figure 5b, the proposed linear elastic approximation $\Psi_{\bar{\alpha}}^{\text{sph}} + \Psi_{\bar{\alpha}}^{\text{dev}}$ is extremely close to the CVaR optimistic estimate $\Psi_{\bar{\alpha}}$ in the present case. Similarly, the linear elastic approximation $\Psi_{\bar{\alpha}}^{\text{sph}} + \Psi_{\bar{\alpha}}^{\text{dev}}$ is also extremely close to the CVaR pessimistic estimate $\Psi_{\bar{\alpha}}$ for this example.

To assess the quality of the approximation in other situations, we also considered Case 2 of Table 3 where $\kappa_i$ and $\mu_i$ are assumed to be lognormal random variable with a 20% standard deviation and a mean value corresponding to the previous values of $E_i^0$ and $v_i^0$. Results are reported in Figure 6 where we can first observe that the empirical PDF are slightly different, especially regarding the bulk modulus distribution. As regards the linear elastic approximations,
<figure><img src="image_9.png" /><figcaption>(b) Level set $2\Psi^*(\Sigma) = 1$ for the different free-energy estimates.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>(b) Level set $2\Psi^*(\Sigma) = 1$ for the different free-energy estimates.</figcaption></figure>

<caption>**Figure 6.** Results for Case 2.</caption>

the quality is still very good, exhibiting only a small overestimation (resp. underestimation) compared to the original CVaR estimates for mixed spheric–deviatoric stress states. As a result, we can conclude that the proposed approximations will furnish very good risk-averse linear elastic estimates of the stochastic composite material.

Finally, to further assess under which conditions the proposed approximation becomes too crude, we investigated the very artificial and non-realistic setting of Case 3 of Table 4 in which the Young modulus is again lognormal with the same mean and 20% standard deviation as before but in which both phases Poisson ratios are assumed to be uniformly distributed in the interval [-0.45; 0.45]. The corresponding results are reported in Figure 7. In this artificial case, the corresponding PDF of the bulk modulus is extremely heavy-tailed and the proposed linear elastic approximations to the CVaR estimates deviate significantly from the latter for mixed spheric/deviatoric stress states. One can clearly see in this case that the CVaR estimates $\Psi_{\bar{\alpha}}$
<figure><img src="image_11.png" /><figcaption>**Figure 7.** Results for Case 3.</figcaption></figure>

and $\Psi_{\alpha}$ are not quadratic functions since their level set is far from being an ellipse, whereas the proposed elastic approximation are the tightest ellipse which coincides for purely spherical and deviatoric states.

**6. Conclusions and perspectives**

In this work, we proposed to use the conditional value-at-risk (CVaR) as a risk measure to replace a stochastic elastic model with either an optimistic or a pessimistic model depending on a given confidence level $\alpha$. The main interest of the CVaR is that it benefits from a definition based on convex optimization. More precisely, when applied to the stochastic convex free energy of a stochastic elastic material, the CVaR free energy is still convex (but not necessarily quadratic) and defines an optimistic elastic material. The pessimistic counterpart is obtained when taking the CVaR of the complementary energy. Our approach enables to define in a systematic and consistent way risk-averse estimates of a stochastic convex potential. Finally, since the resulting risk-averse behaviors are elastic but not necessarily linear anymore, we also proposed lower and upper approximations based on the Kelvin decomposition of the considered symmetry class.

These developments obviously pave the way to their extension towards nonlinear constitutive models. The convex property of the CVaR risk measure makes it natural to consider the case of generalized standard materials (GSM) characterized by stochastic convex free energy and dissipation pseudo-potential. Upon proper definition of the CVaR-based risk averse estimates of such free energy and dissipation potentials, one could therefore aim at proposing risk-averse estimates of various stochastic nonlinear behavior entering the GSM framework. Finally, the computation of the CVaR might require a large number of samples when using a Monte-Carlo approximation. Computing such samples might be prohibitive in the case of full-field simulations for instance. In this case, various approaches attempting at reducing this cost could be used such as variance reduction, importance sampling or surrogate modeling.

**Conflicts of interest**

The author has no conflict of interest to declare.
## Supplementary data

Supporting information for this article is available on the journal's website under https://doi.org/10.5802/crmeca.171 or from the author.

## Appendix A. Dual characterization of CVaR

Starting from the convex representation formula (11), $\Psi_{\bar{\alpha}}$ also benefit from the following convex dual formulation (see for instance [22]):

$$ \begin{array}{l} \Psi_{\bar{\alpha}}(\mathbf{E}) = \sup_z \mathbb{E}[z\Psi(\mathbf{E})] \\ \text{s.t.} \quad 0 \le z \le \frac{1}{1-\alpha} \\ \mathbb{E}[z] = 1 \end{array} \qquad (40) $$

Applying a Legendre-Fenchel transform, one has:

$$ \begin{array}{l} (\Psi_{\bar{\alpha}})^*(\boldsymbol{\Sigma}) = \sup_{\mathbf{E}} \boldsymbol{\Sigma} : \mathbf{E} - \Psi_{\bar{\alpha}}(\mathbf{E}) = \sup_{\mathbf{E}} \inf_{z} \boldsymbol{\Sigma} : \mathbf{E} - \mathbb{E}[z\Psi(\mathbf{E})] \\ \text{s.t.} \quad 0 \le z \le \frac{1}{1-\alpha} \\ \mathbb{E}[z] = 1 \end{array} \qquad (41) $$

Restricting to a discrete setting with $N$ random variables for simplicity and exchanging the inf/sup order, the objective function is given by:

$$ \inf_{z} \sup_{\mathbf{E}} \boldsymbol{\Sigma} : \mathbf{E} - \frac{1}{N} \sum_{i=1}^{N} z_i \Psi(\mathbf{E}; \boldsymbol{\zeta}_i) \qquad (42) $$

$$ = \inf_{z} \sup_{\mathbf{E}} \boldsymbol{\Sigma} : \mathbf{E} - \frac{1}{N} \sum_{i=1}^{N} g_i(\mathbf{E}) \qquad (43) $$

$$ = \infinf_{z} \inf_{\boldsymbol{\sigma}_i} \frac{1}{N} \sum_{i=1}^{N} g_i^*(\hat{\boldsymbol{\Sigma}}_i) \qquad (44) $$
$$ \text{s.t.} \quad \frac{1}{N} \sum_{i=1}^{N} \hat{\boldsymbol{\Sigma}}_i = \boldsymbol{\Sigma} $$

where $g_i(\mathbf{E}) = z_i \Psi_i(\mathbf{E})$ with $\Psi_i(\mathbf{E})$ denoting the $i$-th realization of $\Psi(\mathbf{E}; \boldsymbol{\zeta})$. Moreover, one has:

$$ g_i^*(\boldsymbol{\Sigma}) = z_i \Psi_i^* \left( \frac{\boldsymbol{\Sigma}}{z_i} \right) \qquad (45) $$

Accounting for the fact that $\Psi_i^*$ is quadratic, one finally has:

$$ \begin{array}{l} (\Psi_{\bar{\alpha}})^*(\boldsymbol{\Sigma}) = \inf_{\boldsymbol{\Sigma}_i, z_i} \frac{1}{N} \sum_{i=1}^{N} \frac{1}{z_i} \Psi_i^*(\hat{\boldsymbol{\Sigma}}_i) \\ \text{s.t.} \quad \frac{1}{N} \sum_{i=1}^{N} \hat{\boldsymbol{\Sigma}}_i = \boldsymbol{\Sigma} \\ \qquad 0 \le z_i \le \frac{1}{1-\alpha} \\ \qquad \frac{1}{N} \sum_{i=1}^{N} z_i = 1 \end{array} \qquad (46) $$
which can be generalized to a continuous setting as follows:

$$ (\Psi_{\bar{\alpha}})^*(\Sigma) = \inf_{\hat{\Sigma}, z} \mathbb{E} \left[ \frac{1}{z} \Psi^*(\hat{\Sigma}) \right] \\ \text{s.t.} \quad \mathbb{E}[\hat{\Sigma}] = \Sigma \\ 0 \le z \le \frac{1}{1-\alpha} \\ \mathbb{E}[z] = 1 \qquad (47) $$

## References

[1] L. Göbel, T. Lahmer, A. Osburg, "Uncertainty analysis in multiscale modeling of concrete based on continuum micromechanics", *Eur. J. Mech. A Solids* **65** (2017), p. 14-29.

[2] A. Socié, Y. Monerie, F. Péralès, "Effects of the microstructural uncertainties on the poroelastic and the diffusive properties of mortar", *J. Theor. Comput. Appl. Mech.* (2022), p. 1-21.

[3] S.-i. Sakata, F. Ashida, T. Kojima, M. Zako, "Three-dimensional stochastic analysis using a perturbation-based homogenization method for elastic properties of composite material considering microscopic uncertainty", *Int. J. Solids Struct.* **45** (2008), no. 3-4, p. 894-907.

[4] S. Sakata, F. Ashida, M. Zako, "Kriging-based approximate stochastic homogenization analysis for composite materials", *Comput. Methods Appl. Mech. Eng.* **197** (2008), no. 21-24, p. 1953-1964.

[5] G. Stefanou, "The stochastic finite element method: past, present and future", *Comput. Methods Appl. Mech. Eng.* **198** (2009), no. 9-12, p. 1031-1051.

[6] J. Guilleminot, A. Noshadravan, C. Soize, R. G. Ghanem, "A probabilistic model for bounded elasticity tensor random fields with application to polycrystalline microstructures", *Comput. Methods Appl. Mech. Eng.* **200** (2011), no. 17-20, p. 1637-1648.

[7] J. Guilleminot, C. Soize, "Stochastic modeling of anisotropy in multiscale analysis of heterogeneous materials: A comprehensive overview on random matrix approaches", *Mech. Mater.* **44** (2012), p. 35-46.

[8] R. Rackwitz, "Reliability analysis—a review and some perspectives", *Struct. Saf.* **23** (2001), no. 4, p. 365-395.

[9] M. Kamiński, M. Kleiber, "Perturbation based stochastic finite element method for homogenization of two-phase elastic composites", *Comput. Struct.* **78** (2000), no. 6, p. 811-826.

[10] M. Tootkaboni, L. Graham-Brady, "A multi-scale spectral stochastic method for homogenization of multi-phase periodic composites with random material properties", *Int. J. Numer. Methods Eng.* **83** (2010), no. 1, p. 59-90.

[11] B. Hiriyur, H. Waisman, G. Deodatis, "Uncertainty quantification in homogenization of heterogeneous microstructures modeled by XFEM", *Int. J. Numer. Methods Eng.* **88** (2011), no. 3, p. 257-278.

[12] M. Thapa, S. B. Mulani, R. W. Walters, "Stochastic multi-scale modeling of carbon fiber reinforced composites with polynomial chaos", *Compos. Struct.* **213** (2019), p. 82-97.

[13] V. Dubey, A. Noshadravan, "A probabilistic upscaling of microstructural randomness in modeling mesoscale elastic properties of concrete", *Comput. Struct.* **237** (2020), article no. 106272.

[14] R. T. Rockafellar, "Coherent approaches to risk in optimization under uncertainty", in *OR Tools and Applications: Glimpses of Future Technologies*, Informs, Catonsville, MD, 2007, p. 38-61.

[15] P. Artzner, F. Delbaen, J.-M. Eber, D. Heath, "Coherent measures of risk", *Math. Financ.* **9** (1999), no. 3, p. 203-228.

[16] R. T. Rockafellar, S. Uryasev, "Conditional value-at-risk for general loss distributions", *J. Bank. Financ.* **26** (2002), no. 7, p. 1443-1471.

[17] R. T. Rockafellar, S. Uryasev et al., "Optimization of conditional value-at-risk", *J. Risk* **2** (2000), p. 21-42.

[18] K. Pavlikov, S. P. Uryasev, "CVaR norm and applications in optimization", *Optim. Lett.* **8** (2014), p. 1999-2020.

[19] B. Budiansky, "On the elastic moduli of some heterogeneous materials", *J. Mech. Phys. Solids* **13** (1965), no. 4, p. 223-227.

[20] J. Guilleminot, C. Soize, "On the statistical dependence for the components of random elasticity tensors exhibiting material symmetry properties", *J. Elast.* **111** (2013), no. 2, p. 109-130.

[21] M. M. Mehrabadi, S. C. Cowin, "Eigentensors of linear anisotropic elastic materials", *Q. J. Mech. Appl. Math.* **43** (1990), no. 1, p. 15-41.

[22] M. Ang, J. Sun, Q. Yao, "On the dual representation of coherent risk measures", *Ann. Oper. Res.* **262** (2018), no. 1, p. 29-46.