Jeremy Bleyer
email: jeremy.bleyer@enpc.fr

Stochastic generalized standard materials and risk-averse effective behavior

Keywords: uncertainty, dissipative materials, convexity, risk measure, stochastic method, variational principle, conditional value-at-risk

In this work, we develop a theoretical formulation for describing dissipative material behaviors in a stochastic setting, using the framework of Generalized Standard Materials (GSM). Our goal is to capture the variability inherent in the material model while ensuring thermodynamic consistency, by employing the mathematical framework of stochastic programming. We first show how average behaviors can be computed using the expected value of the free energy and dissipation pseudo-potentials. We then introduce the concept of a risk-averse effective measure, which provides both an optimistic and a pessimistic estimate of the uncertain material behavior. To this end, we utilize the Conditional Value-at-Risk, a widely used risk measure in mathematical finance. We also demonstrate how these concepts can be extended to variational problems at the structure scale, allowing us to compute the effective response of a structure composed of a stochastic material.

Introduction

When modeling and simulating real materials, it is important to account for uncertainties in material properties that are inevitable. In traditional deterministic constitutive modeling, material models are typically calibrated using a set of experimental data. However, even though these experimental data sets generally exhibit statistical distribution, the models are often calibrated based solely on the mean of the data, disregarding any information about uncertainties.

Advanced techniques of uncertainty quantification (UQ) aim to propagate the uncertainties in material properties from the material scale to the structural scale. These techniques are particularly important for designing and optimizing structures under uncertain operating conditions or environmental loads. They can be classified into two different categories: intrusive and non-intrusive ones, see [START_REF] Ghanem | Handbook of uncertainty quantification[END_REF][START_REF] Sudret | Uncertainty propagation and sensitivity analysis in mechanical models-contributions to structural reliability and stochastic spectral methods[END_REF][START_REF] Sullivan | Introduction to uncertainty quantification[END_REF] for a general review. Non-intrusive approaches are a class of techniques used to propagate uncertainties through computational models without modifying the models themselves. The basic idea behind non-intrusive UQ methods is to use the model as a black box and construct a surrogate model, which is a computationally efficient approximation of the original model. The surrogate model is then used to perform the uncertainty quantification. The main advantage of non-intrusive UQ is that it can be applied to any existing computational model without the need for additional code development. Typical techniques involve Monte-Carlo sampling, polynomial chaos expansion, surface response methods or machine-learning techniques.

Conversely, intrusive UQ methods involve modifying the original computational model by including additional terms to propagate uncertainty through the system. They are called intrusive because they require direct modification of the original simulation code, which may be computationally expensive and time-consuming. However, they might be more accurate than non-intrusive approaches. A typical example is stochastic finite element method (SFEM) [START_REF] Ghanem | Stochastic finite elements: a spectral approach[END_REF], which involves introducing probabilistic models for material properties, loads, and boundary conditions into the finite element model of the structure. The probabilistic models are then propagated through the finite element analysis to obtain probabilistic estimates of the structural response. This approach is however difficult to generalize to highly non-linear material behavior such as elastoplasticity. This is notably due to the important coupling between elasticity, plastic yield surfaces and hardening which can all be stochastic but also due to the history-dependent character of elastoplastic behavior. Some works have attempted to include a stochastic modeling of elastoplastic behaviors. In [START_REF] Einav | A thermomechanical framework of plasticity based on probabilistic micromechanics[END_REF], a thermomechanical elasto-plastic framework for random heterogeneous materials is proposed . In [START_REF] Jeremić | Probabilistic elasto-plasticity: formulation in 1d[END_REF],the probabilistic yielding of 1D elastoplasticity is considered and models the elastoplastic evolution using nonlocal Fokker-Planck-Kolmogorovn equations. This methodology has then been extended to other plastic behaviors in 3D Jeremić and Sett [START_REF] Jeremić | On probabilistic yielding of materials[END_REF], Karapiperis et al. [START_REF] Karapiperis | Fokker-planck linearization for non-gaussian stochastic elastoplastic finite elements[END_REF]. In Arnst and Ghanem [START_REF] Arnst | A variational-inequality approach to stochastic boundary value problems with inequality constraints and its application to contact and elastoplasticity[END_REF], stochastic variational inequalities were considered with application to unilateral contact and elastoplasticity. The latter have been approximated using a polynomial chaos expansion, collocation of the inequality constraints and subsequently solved using discrete optimization solvers.

In an engineering context, it is often desirable to characterize the tail of the probability distribution of some quantity of interest. In this case, the problem is often reformulated as computing the probability that the quantity of interest exceeds a certain threshold, which is characterized by the "failure" of the system. Such structural reliability methods [START_REF] Ditlevsen | Structural reliability methods[END_REF][START_REF] Lemaire | Fiabilité des structures: Couplage mécano-fiabiliste statique[END_REF] may rely on Monte-Carlo sampling but they quickly become prohibitive when computing events of low probability. Other techniques may involve importance sampling or approximations based on first-order/second-order reliability method (FORM/SORM). In the present work, we will also aim at characterizing the tail of the material probability distribution to provide engineers with a pessimistic estimate of the response which can be used for safe design.

The proposed contribution falls into the class of intrusive UQ and aims at developing a stochastic variational formulation for the response of a structure consisting of non-linear and non-smooth material behavior, possibly rate and history-dependent. In particular, our work will be positioned in the context of standard generalized materials which characterize complex material non-linear behavior by a set of internal state variables, a free energy and dissipation pseudo-potential. By considering that both of them are convex potentials, the fulfillment of thermodynamic consistency is guaranteed by construction. This level of generality permits the description of a broad range of non-linear material constitutive relations including viscoplasticity, elastoplasticity, damage or a combination thereof. Despite the theoretical importance of this framework, it is surprising that only little attention has been paid to the study of GSM in a stochastic framework where material properties characterizing the material behavior are uncertain. Our first goal is therefore to provide a thermodynamically consistent formulation of stochastic GSM. Let us point out that a large body of work has been devoted to the modeling and generation of uncertain material parameters, e.g. in elasticity [START_REF] Guilleminot | A probabilistic model for bounded elasticity tensor random fields with application to polycrystalline microstructures[END_REF][START_REF] Guilleminot | Stochastic modeling of anisotropy in multiscale analysis of heterogeneous materials: A comprehensive overview on random matrix approaches[END_REF], hyperelasticity [START_REF] Caylak | Stochastic hyperelastic modeling considering dependency of material parameters[END_REF][START_REF] Mihai | Stochastic isotropic hyperelastic materials: constitutive calibration and model selection[END_REF], etc. This aspect is not covered in this work and assumes that we can describe a stochastic material behavior with certain parameters of known probability distribution.

Second, we will use the framework of stochastic programming [START_REF] Birge | Introduction to stochastic programming[END_REF][START_REF] Shapiro | A tutorial on stochastic programming[END_REF] to formulate an incremental stochastic variational principle describing the evolution of such materials. To our knowledge, this is the first contribution which formulates the evolution equations of stochastic dissipative materials using the concepts of stochastic programming. Convex stochastic programming is a well studied field which considers convex optimization problem formulated on cost functions and constraints which depend on random parameters of known probability distributions. In our case, the use of stochastic programming offers several advantages. First, it will enable to preserve a key aspect of the GSM framework which is the convexity of the thermodynamic potential. We can therefore hope to formulate an effective behavior based on convex effective potentials which will still belong to the GSM framewok, ensuring thermodynamic consistency. Second, it will allow us to distinguish more clearly the dependence of various state variables with respect to the uncertainty.

Finally, stochastic programming often relies on the use of the expected value of the cost function over all realizations. However, a structure is not expected to operate well on average only, but with respect to all realizations up to a given confidence level. The general aim is that its probability of failure should be sufficiently low i.e. that the risk of failure is acceptable. The expected value therefore ignores the occurrence of undesirable events in the tail distribution. Similarly, in finance, a portfolio must be chosen in order to safeguard against large losses. It then becomes more relevant to consider a risk-averse measure of the objective function which will also minimize the risk of undesirable events in the design process rather than the average behavior. In the following, we will discuss the choice of a risk-averse measure to be used in a stochastic programming approach. Our objective is to obtain a pessimistic estimate of the structural performance for some given confidence level. Typically, the resulting pessimistic estimate should reflect a material with weaker stiffness and strength properties than the average behavior, see Figure 1.

To do so, we will rely on developments in convex risk measures used in finance such as the Conditional Value-at-Risk (CVaR), which we recently introduced in the context of stochastic elastic materials [START_REF] Bleyer | Risk-averse estimates of effective properties in heterogeneous elasticity[END_REF]. To apply this concept to the underlying GSM framework, we will also need to define a dual CVaR (dCVaR) which, to our knowledge, was never investigated in the previous mathematical literature.

To summarize, the proposed manuscript contributions are threefold.

• We propose and analyze stochastic programming formulations for deriving the average behavior of stochastic GSM models.

• We introduce CVaR and dCVaR risk measures to formulate risk-averse (either optimistic or pessimistic) effective behaviors of a stochastic GSM.

• We extend the derived risk-averse formulations to global variational principles for a structure consisting of a stochastic dissipative material.

Section 2 provides a brief overview of the GSM framework in a deterministic setting and discusses primal and dual incremental principle which can be used to solve the corresponding evolution equations. Section 3 introduces the stochastic programming formalism of random GSM and discusses different formulations of the effective behavior using the expected value operator. Section 4.1 introduces the concept of coherent risk measures and studies the mathematical properties of the CVaR. The novel dual CVaR risk measure is presented and analyzed. Both risk measures are then used in Section 5 to formulated optimistic and pessimistic effective behaviors of a stochastic GSM. Finally, Section 6.2 extends the same concepts to the computation of a risk-averse effective response of a structure. Section 7 closes the manuscript with some conclusions and perspectives for future research.

Notations: Bold-face symbols indicate tensors (stress, strain) or vectors. E [X] denotes the expected value of a random variable X. E [P ] [Y ] denotes the probability of an event Y . The convex conjugate f * (y) of a convex function f (x) is defined as:

f * (y) = sup x x • y -f (x) (1)
For other concepts of convex analysis and convex optimization, we refer the reader to the seminal books of Rockafellar [START_REF] Rockafellar | Convex analysis[END_REF] and Boyd and Vandenberghe [START_REF] Boyd | Convex optimization[END_REF].

Deterministic generalized standard materials and variational principles

Building upon the works of Ziegler [START_REF] Ziegler | Some extremum principles in irreversible thermodynamics, with application to continuum mechanics[END_REF] and Moreau [START_REF] Moreau | Sur les lois de frottement, de plasticité et de viscosité[END_REF] who introduced the notion of a pseudo-potential of dissipation, Halphen and Nguyen [START_REF] Halphen | Sur les matériaux standard généralisés[END_REF] formalized such concepts more generally with the frameworks of Generalized Standard Materials (GSM). Such a setting enables in particular to easily construct constitutive laws which automatically satisfy the principles of thermodynamics.

Generalized Standard Materials

Let us consider a material described by the total strain ε and a set of internal state variables α which can describe various irreversible phenomena such as plasticity, viscoplasticity or damage. The GSM framework postulates the existence of:

• a Helmholtz free energy density ψ(ε, α) depending on the state variables (ε, α)

• a dissipation pseudo-potential φ( ε, α) depending on the rate of state variables ( ε, α).

Note that φ might also depend on the state itself but we do not discuss this case for simplicity.

Moreover, both potentials ψ and φ are assumed to be convex function of their respective variables and are supposed to be non-negative and to vanish at the origin:

ψ(ε, α) ≥ 0 ∀(ε, α) and ψ(0, 0) = 0 (2) φ( ε, α) ≥ 0 ∀( ε, α) and φ(0, 0) = 0 (3)
In duality to internal state variables, there exist associated thermodynamic forces. The constitutive equations (state laws) relating forces to state variables are obtained from the corresponding potentials:

σ = σ nd + σ d (4a) 0 = Y nd + Y d (4b) (σ nd , Y nd ) ∈ ∂ (ε,α) ψ(ε, α) (4c) (σ d , Y d ) ∈ ∂ ( ε, α) φ( ε, α) (4d)
where the subscript nd (resp. d ) refers to a non-dissipative (resp. dissipative) quantity and where ∂ generally denotes the subdifferential. Note that, in most cases, the free energy is a smooth function so that the subdifferential often reduces to a classical differentiation whereas the pseudo-potential is often non-smooth in many applications (e.g. plasticity). Moreover, the Young-Fenchel property also provides the inverse constitutive relations using the corresponding conjugate potentials:

(ε, α) ∈ ∂ (σ nd ,Y nd ) ψ * (σ nd , Y nd ) (5a) ( ε, α) ∈ ∂ (σ d ,Y d ) φ * (σ d , Y d ) (5b)
Finally, if equations (4) describe the constitutive relations for a generic GSM materials, we must mention two important specific cases. The first one is the case where the total strain ε is not considered as a dissipative variable e.g. (visco)elasticity, elasto(visco)plasticity, damage, etc. Then, σ d = 0 and σ = σ nd . The second important case is that of rate-independent behaviors such as elastoplasticity, as opposed to elastoviscoplasticity for instance. Such cases correspond to a homogeneous dissipation pseudo-potential. This specific situation results in (4d) being an inclusion in a convex set and (5b) being an inclusion in the corresponding normal cone.

Primal incremental variational principle

The quasi-static evolution over [0; T ] of a system Ω made of a GSM material can be expressed as the following time-continuous variational principle:

inf u(t),α(t) T 0 Ω d dt ψ(ε, α) + φ( ε, α; α) dΩ dt - T 0 P ext ( u) dt (6)
under appropriate initial conditions α(0) = α 0 and kinematic boundary conditions. The minimizers (u(t), α(t)) provide the evolution of the system over the whole time interval [0; T ]. In practice, solving the above evolution is achieved using a time discretization scheme. Dividing the total interval into finite time intervals [t n ; t n+1 ], an approximate solution (u n+1 , α n+1 ) is computed based on the knowledge of the solution at the previous time step (u n , α n ). There may exist different strategies for choosing a time discretization scheme. The most common strategy is to discretize in time the evolution equations (4), using a fully implicit Euler scheme or a θ-scheme. The nonlinear resolution of the constitutive equations is then embedded into a global Newton-Raphson method solving for the system equilibrium. In some instances, the resulting time-discrete system may loose some structural properties of the original evolution equation, leading for instance to non-symmetric tangent operators.

Another approach relies on a discretization of the original time-continuous principle [START_REF] Bleyer | Applications of conic programming in non-smooth mechanics[END_REF]. Following [START_REF] Maso | Quasistatic evolution problems for linearly elastic-perfectly plastic materials[END_REF][START_REF] Miehe | Strain-driven homogenization of inelastic microstructures and composites based on an incremental variational formulation[END_REF][START_REF] Mielke | Evolution of rate-independent systems[END_REF][START_REF] Ortiz | The variational formulation of viscoplastic constitutive updates[END_REF], we replace in [START_REF] Bleyer | Applications of conic programming in non-smooth mechanics[END_REF] written over [t n , t n+1 ], the rates with finite increments i.e.: ε ≈ ∆ε ∆t , α ≈ ∆α ∆t [START_REF] Bleyer | Risk-averse estimates of effective properties in heterogeneous elasticity[END_REF] where ∆ = n+1n . Moreover, we approximate the power of external forces over [t n ; t n+1 ] by the corresponding work at time t n+1 :

t n+1 tn P ext ( u) dt ≈ W ext,n+1 (u) = Ω f n+1 • u dΩ + ∂Ω N T n+1 • u dS (8)
where f and T are applied body forces and surface tractions respectively.

The solution (u n+1 , α n+1 ) is computed from the resolution of the following incremental variational problem:

(u n+1 , α n+1 ) = arg inf u∈U ad ,α Ω ψ(ε, α) + ∆tφ ∆ε ∆t , ∆α ∆t dΩ -W ext,n+1 (u) (9)
where U ad is the space of kinematically admissible displacements satisfying, for instance, zero Dirichlet boundary conditions u = 0 on ∂Ω D = ∂Ω \ ∂Ω N .

The incremental variational problem (9) therefore amounts to solving for displacement and internal state variables increments a potential which is parameterized by the previous state. We see that the resulting problem is convex and can also be formally written in a more compact form as:

u n+1 = arg inf u∈U ad J(ε; ε n , α n ) -W ext,n+1 (u) (10)
where

J(ε; ε n , α n ) = Ω j(ε; ε n , α n ) dΩ (11) j(ε; ε n , α n ) = inf α ψ(ε, α) + ∆tφ ε -ε n ∆t , α -α n ∆t ( 12
)
where J can be seen as an incremental global potential and j is the associated volume density. The latter is obtained through an implicit minimization over the state variables α as a function of the observable state variable ε. Both of them are parameterized by the values (ε n , α n ) of the state variables at the previous time step. The optimality condition:

σ ∈ ∂j(ε; ε n , α n ) ( 13
)
represents the final stress/strain relation where internal state variables are implicitly solved for in the partial minimization. Akin to variational time integrators, using such variational formulations of constitutive relation updates [START_REF] Heuzé | A variational formulation of thermomechanical constitutive update for hyperbolic conservation laws[END_REF][START_REF] Ortiz | The variational formulation of viscoplastic constitutive updates[END_REF] enables to preserve properties of the underlying variational structure, such as symmetries of the tangent operator for instance.

Dual incremental variational principle

Using classical results from convex analysis, there exist a dual principle to (9) which involves stresses σ and thermodynamic forces Y rather than strain and internal variables as primary unknowns:

(σ n+1 , Y n+1 ) = arg inf σ,Y ,σ nd ,σ d ,Y nd ,Y d Ω ψ * (σ nd , Y nd ) + ∆tφ * (σ d , Y d ) + σ d : ε n + Y d • α n dΩ s.t. div σ + f n+1 = 0 in Ω σn = T n+1 on ∂Ω N σ = σ nd + σ d 0 = Y nd + Y d (14)
In the case where ε is not considered to be a dissipative variable, it reduces to:

(σ n+1 , Y n+1 ) = arg inf σ,Y Ω (ψ * (σ, -Y ) + ∆tφ * (Y ) + Y • α n ) dΩ s.t. div σ + f n+1 = 0 in Ω σn = T n+1 on ∂Ω N (15)
where we introduced Y = Y d = -Y nd . This expression is consistent with the one derived by De Angelis and Cancellara [START_REF] De Angelis | Multifield variational principles and computational aspects in rate plasticity[END_REF] in the viscoplastic case.

In the case of rate-independent materials, φ is a homogeneous function so that (9) simplifies as:

(u n+1 , α n+1 ) = arg inf u∈U ad ,α Ω (ψ(ε, α) + φ (α -α n )) dΩ -W ext,n+1 (u) (16)
This results in φ * being the indicator of some convex set G and (15) reduces to:

(σ n+1 , Y n+1 ) = arg inf σ,Y Ω (ψ * (σ, -Y ) + Y • α n ) dΩ s.t. div σ + f n+1 = 0 in Ω σn = T n+1 on ∂Ω N Y ∈ G in Ω (17)

Stochastic generalized standard materials

We will now revisit the above GSM formulation of material behavior in the stochastic setting. As we have seen before, dissipative materials can be described by the incremental potential j which is obtained through an implicit minimization over the state variables α as a function of the observable state variable ε. The incremental potential j is parameterized by the values (ε n , α n ) of the state variables at the previous time step. For the sake of simplicity, we will drop the explicit mention of such a dependency in the following. In this section, we work at the material point level and thus do not explicitly mention the dependency of the density j upon the position vector x ∈ Ω.

Formulation of an effective material behavior

We now assume that the incremental potential depends upon uncertain parameters ζ i.e. j(ε; ζ) is now a stochastic convex function. In practice, uncertainty can affect material parameters describing the behavior. In this case, such parameters become random variables. Alternatively, material properties can also vary in space so that uncertain material parameters become random fields. In either case, uncertainty is encoded via the vector of random parameters ζ which are supposed to have a known probability distribution. Now, let us discuss our main objective. Many UQ techniques aim at describing the probability distribution of some quantities of interest based on the input probability distribution of uncertain input parameters. Here, our goal is different in the sense that we would like to be able to describe the effective behavior of the material with respect to the uncertainty ζ. In other words, we would like to formulate a deterministic behavior which captures the essential features of the underlying uncertainties in some way, which remains to be specified. Formally, we introduce the following effective potential:

j eff (ε) = R [j(ε; ζ)] (18)
where R refers to some effective measure of j over all realizations ζ. Since j is defined as a minimum problem involving an uncertain objective function, the computation of j eff therefore falls into the scope of stochastic programming. Working with a stochastic programming framework will indeed prove to be extremely fruitful for defining, in a consistent manner, the effective state and evolutions associated with the effective material behavior.

Obviously, the choice of the effective measure R heavily relies on our wanted definition of the corresponding behavior. The most classical choice is the expected value R[j] = E [j] which considers the average behavior of the material. In this case, variance and low probability events are discarded when using the corresponding behavior for a computation at the structure scale. Borrowing from the finance vocabulary, the expected value is a risk-neutral effective measure.

However, engineering applications often require to consider scenarios of low probability which might be more detrimental to the structural safety. In such a context, the choice of the effective measure could be aimed at reproducing some kind of pessimistic estimate of the underlying material behavior (e.g. smaller stiffness and strength than the average behavior etc.). In this case, we look for a risk-averse effective measure. Conversely, one could also be interested in computing an optimistic estimate of the corresponding uncertain behavior using a risk-seeking measure. If optimistic estimates may be less useful from a safety analysis point of view, they could serve as defining some kind of confidence interval when used in conjunction with a pessimistic estimate. One major difficulty in choosing such measures is that we would like to preserve fundamental mathematical properties such as convexity and positivity of the material behavior potentials. Fortunately, we will be able to do so by introducing later in section 4.1 the notion of convex risk measures.

Before that, let us first discuss the more classical case where R[j] = E [j]. In this case, most observations will bear some similarities with the setting of homogenization of random media. In particular, we show that there exist two possible ways of defining an effective behavior, based either on the primal strain-based variational principle [START_REF] Caylak | Stochastic hyperelastic modeling considering dependency of material parameters[END_REF] or on the corresponding dual stress-based principle [START_REF] Georghiou | Generalized decision rule approximations for stochastic programming via liftings[END_REF]. Moreover, we show that a particular attention must be paid to the definition of internal state variables when using both formulations, especially when yield conditions are involved.

Internal state variables in the stochastic setting

In the stochastic setting, the first question which arises indeed concerns the dependency of the different state variables with respect to the uncertainty. Since we are interested in computing an effective behavior of the material, we consider the total strain ε to be imposed by the external observer. In a global structural analysis, ε will for instance be obtained by minimizing the effective potential associated with j eff , see section 6.2. In the vocable of stochastic programming, ε is referred to as a first-stage variable since its value is chosen independently from the exact realization of the uncertainty ζ. On the contrary, the implicit state variables α are referred to as second-stage variables i.e. their value directly depends on the uncertainty realization. To make this dependency more explicit, we will now use the subscripted notation α ζ .

Another important step in the GSM framework is to define the list of internal state variables. In the elastoplastic deterministic case, one can equivalently introduce the elastic ε el and plastic strain ε p as internal state variables or one could also just introduce ε p and replace ε el with εε p . In the stochastic setting, this choice is less innocent than it appears. Assuming that ε p ζ is the uncertainty-dependent plastic strain, a first possibility is to define the elastic strain as:

ε el ζ = ε -ε p ζ ( 19
)
so that:

ε = ε el ζ + ε p ζ almost surely (20)
i.e. the sum of the random elastic and plastic strains is equal to the deterministic total strain ε with probability 1. Note that a similar approach is considered in Einav and Collins [START_REF] Einav | A thermomechanical framework of plasticity based on probabilistic micromechanics[END_REF].

Another approach would consist in enforcing the last relation on average only. More precisely, we assume that there exists a random elastic strain ε el ζ such that:

ε = E ε el ζ + ε p ζ ( 21
)
Both approaches are obviously not equivalent. In [START_REF] Jeremić | On probabilistic yielding of materials[END_REF], we have explictly introduced an additional internal state variable ε el ζ which is related to the total strain in average only.

Both approaches can also be reformulated as follows. Considering the observable strain ε, we also introduce the uncertain total strain deviation from average e ζ with the property E [e ζ ] = 0. This strain deviation is such that the real uncertain total strain is ε ζ = ε + e ζ and such that E [ε ζ ] = ε. We therefore consider e ζ as an "internal" state variable which is not observable from the perspective of the effective structural behavior. The elastic strain is then defined as:

ε el ζ = ε ζ -ε p ζ = ε + e ζ -ε p ζ ( 22
)
We then see that the first approach ( 19) corresponds to the situation where we impose that the deviation is zero e ζ = 0 whereas the second approach ( 21) corresponds to the assumption that the deviation is zero on average only.

In fact, this point of view involving e ζ as an additional state variable is not specific to elastoplasticity but generalizes to any material behavior in the sense that:

j eff (ε) = inf e ζ ,α ζ E [ψ(ε + e ζ , α ζ )] + E [φ (ε + e ζ -ε n -e n,ζ , α ζ -α n,ζ )] s.t. e ζ = 0 (23)
for the first approach or:

j eff (ε) = inf e ζ ,α ζ E [ψ(ε + e ζ , α ζ )] + E [φ (ε + e ζ -ε n -e n,ζ , α ζ -α n,ζ )] s.t. E [e ζ ] = 0 (24)
for the second approach.

For instance, in the case of linear elasticity characterized by an uncertain elasticity tensor C ζ , the first approach with e ζ = 0 yields:

j eff (ε) = E 1 2 ε : C ζ : ε = 1 2 ε : E [C ζ ] : ε (25)
whereas the second approach yields:

j eff (ε) = inf e ζ E 1 2 (ε + e ζ ) : C ζ : (ε + e ζ ) = 1 2 ε : E C -1 ζ -1 : ε s.t. E [e ζ ] = 0 (26)
We therefore see that the first approach will produce an elastic behavior characterized by the average elastic moduli whereas the second expression is associated with the harmonic average of the elastic moduli. Note that both quantities are not equal in general, as we only have

E [C ζ ] E C -1 ζ -1
where inequality is understood in the sense of quadratic forms. Note that a result similar to (26) can also be found in [START_REF] Junker | A relaxation approach to modeling the stochastic behavior of elastic materials[END_REF].

In the general case, the effective dual potential (j eff ) * (σ) where σ is the observable stress is given by:

(j eff ) * (σ) = sup ε,e ζ ,α ζ σ : ε -E [ψ(ε + e ζ , α ζ )] -E [φ (ε + e ζ -ε n -e n,ζ , α ζ -α n,ζ )] s.t. e ζ = 0 or E [e ζ ] = 0 ( 27
)
For the first approach (23), one obtains:

(j eff ) * (σ) = inf σ nd ζ ,σ d ζ ,Y ζ E ψ * (σ nd ζ , -Y ζ ) + E φ * (σ d ζ , Y ζ ) + E σ d ζ : ε n + Y ζ • α n,ζ s.t. σ = E σ nd ζ + σ d ζ ( 28
) whereas for the second approach [START_REF] Karapiperis | Fokker-planck linearization for non-gaussian stochastic elastoplastic finite elements[END_REF], one obtains:

(j eff ) * (σ) = inf σ nd ζ ,σ d ζ ,Y ζ E ψ * (σ nd ζ , -Y ζ ) + E φ * (σ d ζ , Y ζ ) + E σ d ζ : (ε n + e n,ζ ) + Y ζ • α n,ζ s.t. σ = σ nd ζ + σ d ζ a.s. ( 29
)
We therefore see that when we do not consider the strain deviation e ζ as an internal state variable, the thermodynamic stress σ nd ζ + σ d ζ should be considered as a second-stage variable which is equal to the observable stress σ only on average. On the contrary, when we do consider the strain deviation e ζ as an internal state variable, the thermodynamic stress becomes deterministic. For instance, in the case where ε is not a dissipative variable,

σ d ζ = 0 and σ ζ = σ nd ζ . Thus, one has either E [σ ζ ] = σ or σ ζ = σ respectively.
To conclude, this discussion bears striking similarities with Voigt and Reuss bounds in periodic homogenization. Indeed, if dependence upon the uncertainty ζ is replaced with a spatial dependence and expectation is replaced with a spatial average, ( 23) is akin to compute the average potentials with a uniform strain ε in space. This indeed leads to a "stiff" Voigt estimate of the elastic effective modulus as in [START_REF] Lemaire | Fiabilité des structures: Couplage mécano-fiabiliste statique[END_REF]. On the contrary, ( 29) is akin to imposing a uniform stress field σ in space, leading to a "soft" Reuss estimate as in [START_REF] Mafusalov | Conditional value-at-risk (cvar) norm: Stochastic case[END_REF] in the elastic case.

Illustrative application in 1D

We close this section by illustrating the different stochastic formulations on the example of a one-dimensional elasto-plastic behavior with isotropic power-law hardening. The internal state variables are α = (ε p , p) where p is the cumulated plastic strain. The corresponding 1D potentials are given by:

ψ(ε, ε p , p) = ψ el (ε, ε p ) + ψ h (p) (30)
ψ el (ε, ε p ) = 1 2 E(ζ)(ε -ε p ) 2 (31)
ψ h (p) = 1 1 + 1/m H(ζ)p 1+1/m (32) φ( εp , ṗ) = σ 0 (ζ) ṗ for | εp | ≤ ṗ +∞ otherwise ( 33
)
We consider that the Young modulus E(ζ), hardening modulus H(ζ) and yield stress σ 0 (ζ) are independent random variables following a lognormal distribution. We assume the powerlaw exponent to be deterministic (m = 3 in the following numerical simulations).

The stochastic programming formulation [START_REF] Junker | A relaxation approach to modeling the stochastic behavior of elastic materials[END_REF] where we consider the total strain to be fully deterministic reads:

j eff (ε) = inf ε p ζ ,p ζ 1 2 E E(ζ)(ε -ε p ζ ) 2 + 1 1 + 1/m E H(ζ)p 1+1/m ζ + E [σ 0 (ζ)(p ζ -p n,ζ )] s.t. |ε p ζ -ε p n,ζ | ≤ p ζ -p n,ζ (34)
The corresponding dual potential is given by:

(j eff ) * (σ) = inf σ ζ ,R ζ 1 2 E 1 E(ζ) σ 2 ζ + 1 1 + m E 1 H(ζ) m R 1+m ζ + E σ ζ ε p n,ζ -R ζ p n,ζ s.t. |σ ζ | ≤ σ 0 (ζ) ∀ζ σ = E [σ ζ ] (35)
The second stochastic programming formulation [START_REF] Karapiperis | Fokker-planck linearization for non-gaussian stochastic elastoplastic finite elements[END_REF] where we consider the additional total strain deviations reads:

j eff (ε) = inf e ζ ,ε p ζ ,p ζ 1 2 E E(ζ)(ε + e ζ -ε p ζ ) 2 + 1 1 + 1/m E H(ζ)p 1+1/m ζ + E [σ 0 (ζ)(p ζ -p n,ζ )] s.t. |ε p ζ -ε p n,ζ | ≤ p ζ -p n,ζ (36)
The corresponding dual potential is given by:

(j eff ) * (σ) = inf R ζ 1 2 E 1 E(ζ) σ 2 + 1 1 + m E 1 H(ζ) m R 1+m ζ + E σε p n,ζ -R ζ p n,ζ s.t. |σ| ≤ σ 0 (ζ) ∀ζ (37)
Figure 2a represents the corresponding effective behavior for both formulations in the case where E [E] = 20, E [H] = 1 and E [σ 0 ] = 1 and a 20% standard deviation for both E and H and 10% for σ 0 . The effective behavior is obtained by discretizing the imposed total strain in 30 increments. Each increment is obtained by solving either formulation (34) (in blue) or [START_REF] Rockafellar | Convex analysis[END_REF] (in red). In both cases, expectations have been computed using a Monte-Carlo sampling approximation with a sample size N = 500. Minimization over the N different state variables is performed using the cvxpy package [START_REF] Agrawal | A rewriting system for convex optimization problems[END_REF][START_REF] Diamond | CVXPY: A Python-embedded modeling language for convex optimization[END_REF].

First, it can be seen that both approaches lead to a very similar behavior, the second formulation [START_REF] Rockafellar | Convex analysis[END_REF] being slightly softer than the first one due to the presence of additional total strain variables. Interestingly, the obtained effective behavior is close to the nominal deterministic behavior except that it exhibits an earlier and progressive onset of plasticity since it accounts for a distribution of yield stress with different values. This is further confirmed in Figure 2b which corresponds to an almost perfectly-plastic behavior with E [H] = 0.01. We can observe that the effective behavior corresponding to (34) admits an initial hardening phase before reaching a perfectly plastic plateau. A similar behavior is obtained from the approach proposed in [START_REF] Einav | A thermomechanical framework of plasticity based on probabilistic micromechanics[END_REF]. Moreover, in the ultimate perfectly plastic regime, the elastic term vanishes so that ε p ζ = ε. One can then see that [START_REF] Ortiz | The variational formulation of viscoplastic constitutive updates[END_REF] gives approximately plastic case. The final plastic plateau is therefore given by the average yield stress

j eff (ε) ≈ E [σ 0 ] |ε -ε n |
E [σ 0 ] = 1.
Conversely, formulation [START_REF] Rockafellar | Convex analysis[END_REF] clearly fails at reproducing a reasonable effective behavior. This can indeed be well understood when inspecting the dual formulation [START_REF] Rockafellar | Coherent approaches to risk in optimization under uncertainty[END_REF] in which the last constraint will enforce the observable stress to satisfy the plasticity criterion for any possible realization of the yield stress σ 0 (ζ). As a result, this constraint can in fact be reformulated as |σ| ≤ inf ζ σ 0 (ζ) which is overly conservative. Note that we obtain here a non-zero yield stress due to the fact that we have a finite set of samples but, in theory for a lognormal distribution, this minimum should be zero.

A soft-constrained formulation

We might therefore be tempted to conclude that the second formulation (36) allowing such non-zero total strain deviations e ζ is not relevant. However, when inspecting again the dual formulation [START_REF] Rockafellar | Coherent approaches to risk in optimization under uncertainty[END_REF], we can see that a minor modification could lead to a relevant effective estimate, yet different from ( 34)- [START_REF] Pavlikov | Cvar norm and applications in optimization[END_REF]. Indeed, the last constraint can also be formulated as g(σ; ζ) ≤ 1 ∀ζ where g(σ; ζ) = |σ|/σ 0 (ζ) is the gauge function of the corresponding yield criterion. Since the second formulation involves a deterministic stress, we then see that satisfying the yield condition almost surely is too restrictive. Another classical approach in stochastic programming is to replace such a hard constraint by a soft counterpart such as imposing that the constraint is satisfied on average only, hence E [g(σ; ζ)] ≤ 1. In the present case, this amounts to imposing |σ| ≤

1 E 1 σ 0 (ζ)
. The new dual formulation therefore reads:

(j eff ) * (σ) = inf R ζ 1 2 E 1 E(ζ) σ 2 + 1 1 + m E 1 H(ζ) m R 1+m ζ + E σε p n,ζ -R ζ p n,ζ s.t. E [|σ|/σ 0 (ζ)] ≤ 1 ∀ζ (38)
Going back to the corresponding primal formulation, [START_REF] Rockafellar | Conditional value-at-risk for general loss distributions[END_REF] is the dual to the following primal formulation:

j eff (ε) = inf e ζ ,ε p ζ ,p ζ 1 2 E E(ζ)(ε + e ζ -ε p ζ ) 2 + 1 1 + 1/m E H(ζ)p 1+1/m ζ + max ζ [σ 0 (ζ)(p ζ -p n,ζ )] s.t. |ε p ζ -ε p n,ζ | ≤ p ζ -p n,ζ (39)
where the only difference with [START_REF] Rockafellar | Convex analysis[END_REF] is that the effective measure of the dissipation potential is no longer the expectation but the worst-case value.

Figure 3 shows the obtained effective estimate using this third formulation involving soft yield constraints. Clearly, we can see that this formulation removes the deficiencies of (36) which we identified earlier. Besides, it yields an effective behavior which is close to formulation [START_REF] Ortiz | The variational formulation of viscoplastic constitutive updates[END_REF], the main difference being that there is no progressive onset of plasticity in this case. This is due to the fact that the yield constraint is determinisitic. Note that for a lognormal variable X we have:

E X -1 -1 = E [X] 1 + var[X] E [X] 2 -1 ≈ E [X] (40)
when the variance is small as in the present case. This explains for instance in Figure 3b why the effective yield stress are very close for both formulations (roughly only 1% difference since we have a 10% standard deviation on σ 0 ). However, let us point out that we might see more differences for other types of distributions than a lognormal one.

To conclude, we have seen that there exist different possible formulations of an effective behavior for a stochastic dissipative material. Depending on the chosen hypothesis regarding whether total strain is considered as a first-stage or second-stage variable, the associated thermodynamic stress should be respectively considered as a second-stage or a first-stage variable. Depending on the chosen hypothesis, we have also seen that the retained effective measure plays an important role, especially concerning the dissipation potential. These various formulations seem to yield very close effective behavior which can question the usefulness of having such a discussion. However, we considered here the expectation as an effective measure of the behavior. In the following, we will now discuss other choices, in particular risk measures which can take into account low probability events rather than average behavior. In such a case, the distinction between these different formulations will be crucial.

Risk-averse measures

As discussed before, we now aim at departing from a risk-neutral effective behavior and rather aim at representing some pessimistic (or risk-averse) effective behavior of the underlying stochastic material. To do so, we first introduce the concept of coherent risk measures.

Coherent risk measures

In financial mathematics, a risk measure is a functional which amounts to quantify the level of risk of a given random variable X where X can be, for instance, some financial loss or cost, see [START_REF] Rockafellar | Coherent approaches to risk in optimization under uncertainty[END_REF] for a broad overview. In the subsequent presentation, we follow the tradition in optimization where costs should be minimized. Therefore, we consider that large positive values of X are disliked. The choice of a proper risk measure will aim at guarding against large values of X. Many possibilities of risk measures can be considered, depending on how the decision-maker evaluates risks. Typical examples are for instance:

• the expected value: R[X] = E [X] (41)
• the safety margin using k standard deviations:

R[X] = E [X] + k std[X], for k > 0 ( 42
)
• the worst-case value:

R[X] = sup X (43)
• the Value-at-Risk (VaR) for a level β ∈ [0; 1] (or the β-quantile):

R[X] = VaR β (X) = inf{Z s.t. P [Z ≥ X] ≥ β} (44)
• and many more...

Recently, the finance community has devoted its attention to the use of so-called coherent risk measures which are a specific class of risk measures benefiting from desirable properties proposed by Artzner et al. [START_REF] Artzner | Coherent measures of risk[END_REF] such as:

• R[C] = C for all constants C • convexity: R[(1 -λ)X + λY ] ≤ (1 -λ)R[X] + λR[Y ] for all X, Y and λ ∈ [0; 1] • monotonicity: R[X] ≤ R[Y ] if X ≤ Y almost surely • positive homogeneity: R[λX] = λR[X] for λ > 0
Based on the first two properties, a coherent measure of risk will also be translationally invariant i.e. R[X + C] = R[X] + C for any constant C. Finally, one key property in the present work is that convexity is preserved under coherent risk measures i.e.:

if f (x; ζ) is convex wrt x, then R[f ](x) is convex in x ( 45
)
As a first example, the expected value ( 41) is a coherent risk measure but it is not very useful in quantifying the amount of risk since it does not take into account events in the distribution tail. This measure is therefore termed risk-neutral. Risk measures satisfying R[X] > E [X] are instead called risk-averse measures.

The safety margin measure ( 42) is for instance risk-averse but is not coherent. Conversely, the worst-case value (43) fulfills both coherency and risk-aversion but is not useful in practice due to being severely conservative. The Value-at-Risk VaR β (X) might be an interesting risk-averse measure as it gives an estimate of X which is exceeded only with probability 1β. It is however not convex and hence does not fall into the class of coherent risk measures. As a result, stochastic programs in which an uncertain convex function j(ε, ζ) is replaced with VaR β (j)(ε) would not be convex anymore which is extremely undesirable from the computational and theoretical perspective. A very popular risk-averse measure closely linked to the VaR and which overcomes this drawback is the Conditional Value-at-Risk (CVaR).

Conditional Value-at-Risk and risk-averse estimates 4.2.1. Definition

For a continuous distribution, the CVaR is defined as the expected value of X above VaR β (see Figure 4):

CVaR β (X) = E [X s.t. X ≥ VaR β (X)] (46)
The above definition is slightly more technical when considering discrete distributions, we refer the reader to [START_REF] Rockafellar | Conditional value-at-risk for general loss distributions[END_REF] for a rigorous definition. The CVaR is then a coherent risk measure and a key result due to [START_REF] Rockafellar | Optimization of conditional value-at-risk[END_REF] is that CVaR benefits from the following convex optimization characterization:

CVaR β (X) = inf λ λ + 1 1 -β E X -λ + (47)
where + = max{ , 0} denotes the positive part. From this definition, it is clear that CVaR 0 (X) = E [X]. Let us also point out that, if the minimum of (47) is unique, then the optimal value is exactly λ = VaR β (X).

The CVaR β of |x| for a vector x ∈ 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 ∞ -norm. For intermediate values, it can be seen as an average of the k largest values of |x| where k/N = 1β, see [START_REF] Mafusalov | Conditional value-at-risk (cvar) norm: Stochastic case[END_REF][START_REF] Pavlikov | Cvar norm and applications in optimization[END_REF] for a more precise definition.

Application to an elastic potential

Applying such concepts from financial optimization to mechanics requires to decide what would be the "loss" function against which we want to be immunized. This not necessarily an obvious question to answer. Let us first investigate what would be obtained when considering the CVaR as an effective measure R in [START_REF] Guilleminot | Stochastic modeling of anisotropy in multiscale analysis of heterogeneous materials: A comprehensive overview on random matrix approaches[END_REF] in the simple case where the behavior is linear elastic, hence here j(ε) = ψ(ε) = 1 2 ε : C : ε. We assume that uncertainty affects only the material Young modulus so that E(ζ) is a random variable and C(ζ) = E(ζ)C 0 for some reference C 0 material with unit Young modulus. Then, one easily sees that:

ψ eff (ε) = CVaR β (ψ) (ε) = 1 2 ε : CVaR β (E) C 0 : ε (48)
i.e. one obtains a linear elastic material with an effective Young modulus E eff = CVaR β (E).

In the following, we use the more compact notation E β = CVaR β (E) to denote the CVaR effective value. The CVaR estimate enjoys the following properties:

E 0 = E [E] (49) E 1 = sup E (50) E β ≥ VaR β (E) (51) E β ≤ E β ∀β ≤ β (52)
For β ∈ [0; 1], the effective modulus therefore varies from its mean value to its supremum (which can be infinite). In practice, β is often chosen close to 1 e.g. β = 0.95 or 0.99 so that E β provides a best-case estimate of the material Young modulus. It might therefore appear strange to consider this stiff estimate as being a useful risk measure from the engineering point of view, since stiff materials are generally looked for. As it will be clear later, considering a CVaR in a primal displacement-based principle will yield an overestimation of the mechanical performance compared to its mean value. Structural performance underestimation will instead be obtained through duality. The confidence level β will enable us to reach tail behaviors (both left and right) when β is close to 1. We will therefore use in the following the term risk or risk-averse measure to denote a measure of the deviation from the mean value, irrespective of the fact that such a value can be beneficial or detrimental from the engineering point of view. We will also refer to β as the risk-aversion level.

Generic effective stress-strain relationship

Let us now study the stress-strain relationship which would be obtained when considering an effective free energy CVaR β (ψ) (ε) for a generic uncertain free energy ψ(ε; ζ). To avoid technicalities, let us consider that uncertainty is represented by a finite distribution of N scenarios ζ i of probabilities p i . The effective free-energy therefore reads:

ψ eff (ε) = CVaR β (ψ) (ε) = inf λ λ + 1 1 -β N i=1 p i ψ(ε; ζ i ) -λ + (53)
Using classical convex analysis arguments, the associated stress σ ∈ ∂ψ eff (ε) can be obtained as the solution to the following maximization problem:

sup ε σ : ε -ψ eff (ε) = -inf ε,λ λ + 1 1 -β N i=1 p i ψ(ε; ζ i ) -λ + -σ : ε (54)
which is characterized by the following optimality conditions:

1 = 1 1 -β N i=1 p i Heav (ψ(ε; ζ i ) -λ) (55) σ = 1 1 -β N i=1 p i ∂ψ ∂ε (ε; ζ i ) Heav (ψ(ε; ζ i ) -λ) (56)
where Heav (x) = 1 if x > 0 and 0 otherwise. As a result, we can see from (55) that λ is obtained, for a given ε, as the optimal level set of the uncertain free energies such that the set A of active scenarios, i.e. scenarios which have ψ(ε; ζ i ) ≥ λ occur with probability 1β. Knowing such a partition A of active scenarios, condition (56) can also be rewritten as:

σ = 1 1 -β i∈A p i σ i (57)
That is, the effective stress σ is obtained as the weighted sum of the stress states σ i = ∂ψ(ε; ζ i )/∂ε corresponding to each active scenario. Note that if ψ is non-smooth the latter relations should be understood in the sense of set membership.

Dual Conditional Value-at-Risk (dCVaR)

We see from characterization (57), that the obtained CVaR constitutive relation is indeed optimistic in the sense that it retains the scenarios with the highest free energy levels, resulting in a stiff estimate in general. In this section, we would like to derive a corresponding pessimistic estimate which would retain only scenarios with low free energy levels.

Quadratic elastic potentials

In Bleyer [START_REF] Bleyer | Risk-averse estimates of effective properties in heterogeneous elasticity[END_REF], we proposed to obtain the pessimistic estimate of an elastic potential ψ by taking the conjugate function of the CVaR of the conjugate i.e. (ψ * ) * β . In the elastic case, we easily see that:

(ψ * ) β (σ) = 1 2 σ : CVaR β E -1 S 0 : σ (58) (ψ * ) * β (ε) = 1 2 ε : (CVaR β E -1 ) -1 C 0 : ε (59)
where S 0 = (C 0 ) -1 . This defines an elastic material with an effective Young modulus given by:

E eff = E β = CVaR β E -1 -1 = ((E -1 ) β ) -1 (60)
where we introduced the dual Conditional Value-at-Risk (dCVaR) notation E β which enjoys the following properties:

E 0 = E E -1 -1 (61) E 1 = inf E (62) E β ≤ VaR β (E -1 ) -1 (63) E β ≤ E β ∀β ≤ β (64)
In the generic case of a convex function f (x), one has the following convex formulation:

(f * ) * β (x) = inf x,z E [zf ( x/z)] s.t. E [ x] = x 0 ≤ z ≤ 1 1 -β E [z] = 1 (65)
Unfortunately, we see that this definition is not suitable for all convex functions. Indeed, for positive homogeneous functions for instance, the objective simplifies to E [f ( x)] and the auxiliary variable z, and thus β, play no role. For instance, in the case where f (x) = k(ζ) x p , f (x) is the support function of a L q -ball of radius k(ζ) with 1 p + 1 q = 1, one can easily see that (65) yields the ball of smallest radius i.e. (f * ) * β (x) = {inf k(ζ)} x p . This is obviously not what we are looking for as this worst-case estimate is way too conservative and does not depend on β. We can remark that we recover an issue similar to that encountered previously with formulation (37) of the effective behavior.

Polar-based definition

One important property of the CVaR is that if f ≤ g, then CVaR β (f ) ≤ CVaR β (g). The use of f * in the definition of the dCVaR was motivated by the fact f * ≥ g * , so that:

(CVaR β (f * )) * ≤ (CVaR β (g * )) * (66)
Moreover, we then have:

(CVaR β (f * )) * ≤ (CVaR β (f * )) * ∀β ≤ β (67)
which inverts the ordering of the standard CVaR with respect to β. Hence, R[f ] = (CVaR β (f * )) * would define another convex risk measure. Unfortunately, it is not coherent because it lacks the homogeneity property.

To obtain a coherent convex risk-measure, we look for another operation on convex functions which would invert ordering. It turns out that for the specific class of non-negative convex function which vanish at the origin, as our thermodynamical potentials, there exists another transform than the convex conjugate which enjoys this property [START_REF] Milman | Hidden structures in the class of convex functions and a new duality transform[END_REF]. This transform is the polar f • of a positive convex function f which is defined as [START_REF] Rockafellar | Convex analysis[END_REF]Th. 15.4]:

f • (x) = inf{µ ≥ 0 s.t. x • y -µf (y) ≤ 1 ∀y} (68) = inf{µ ≥ 0 s.t. µf * (x/µ) ≤ 1} (69) which also satisfies f • ≥ g • for f ≤ g. Moreover f •• (x) = f (x).
As a result, we define the dual CVaR as:

dCVaR β (f ) (x) = (CVaR β (f • )) • (x) (70)
and use f β as a short-hand notation for dCVaR β (f ).

Properties

In the case of quadratic convex functions, f • = f * so that the polar definition (70) coincides with that given in Bleyer [START_REF] Bleyer | Risk-averse estimates of effective properties in heterogeneous elasticity[END_REF] in the context of linear elasticity.

Moreover, in the case where the uncertain function f

(x; ζ) = k(ζ)g(x)
where k(ζ) > 0 and g is a non-negative convex function, we have that:

f • (x; ζ) = (k(ζ)g(x)) • = 1 k(ζ) g • (x) (71)
so that:

dCVaR β (f ) (x) = CVaR β 1 k(ζ) g • (x) • (72) = 1 CVaR β (k -1 (ζ)) g •• (x) (73) = 1 CVaR β (k -1 (ζ)) g(x) (74)
In particular, for the previous example involving a homogeneous function

f (x) = k(ζ) x p , the new dCVaR estimates results in dCVaR β (f ) (x) = CVaR β (k -1 (ζ)) -1
x p which is much less pessimistic than inf ζ k(ζ). Moreover, we therefore recover the same definition as for quadratic elastic potentials.

As a conclusion, the dCVaR can hence be defined for any constant function (random variable k(ζ)) as follows:

k β := dCVaR β (k(ζ)) = (CVaR β k -1 (ζ) ) -1 = (k -1 β ) -1 (75)
Note that one has:

f β ≤ f β ≤ f 0 = E [f • ] • ≤ f 0 = E [f ] ≤ f β ≤ f β ∀ 0 ≤ β ≤ β < 1 (76)
Let us also recall that the optimal value of λ in the CVaR definition (47) corresponds to the VaR. Then,

g • = VaR β (f • ) is such that F f • (g • ) = β. But one also has F f • (g • ) = P [f • ≤ g • ] = P [f ≥ g] so that: F f • (g • ) = 1 -F f (g) = β ⇒ F f (g) = 1 -β (77)
We conclude that g = VaR β (f • ) • corresponds to the (1β)-quantile of f (x; ζ). However, f β does not correspond to the expectation below this (1β)-quantile but rather the "polar expectation"

E [f • ] • below this quantile.
Finally, we have the following convex representation of dCVaR β (f ) (x):

dCVaR β (f ) (x) = inf v≥0, x max E [vf ( x/v; ζ)] ; (1 -β) sup ζ {vf ( x/v; ζ)} s.t. E [ x] = x E [v] = 1 (78)
Besides, the Lengendre-Fenchel conjugate of the dCVaR is given by:

dCVaR β (f ) * (y) = inf z≥0 sup ζ {zf * (y/z; ζ)} s.t. CVaR β (z) ≤ 1 (79)
Note that when f is a homogeneous function, we can ignore the auxiliary variables v in (78). Besides, in this case f * is the indicator of some convex set G which is represented via its gauge function g = f • as follows: G = {y s.t. f • (y) ≤ 1}. In this case, (79) now becomes:

dCVaR β (f ) * (y) = 0 if CVaR β (f • ) (y) ≤ 1 +∞ otherwise ( 80
)
[ToDo: proofs]

Risk-averse behavior of stochastic GSM

In this section, we discuss the application of the CVaR and dCVaR to derive a risk-averse effective behavior of a stochastic GSM. To illustrate this objective on the example of the elastoplastic 1D behavior of section 3.3, we would like to formulate:

• optimistic effective free energy and dissipation potentials which would result in an effective 1D behavior corresponding to an upper enveloppe of the bundle of stochastic material responses shown in Fig. 2 for instance. The resulting effective behavior is expected to exhibit a stiffer elastic response, a stronger hardening and a larger yield stress than the average or the nominal behavior.

• pessimistic effective free energy and dissipation potentials which would result in an effective 1D behavior corresponding to a lower enveloppe of the same bundle of stochastic material responses. The resulting effective behavior is expected to exhibit a more compliant elastic response, a weaker hardening and a smaller yield stress than the average or the nominal behavior.

Optimistic estimate

A natural extension of the CVaR concepts introduced in section 4.1 is to replace both stochastic free energy potential and stochastic dissipation pseudo-potential with their CVaR risk measure:

ψ eff (ε, α ζ ) = CVaR β (ψ(ε, α ζ ; ζ)) (81) φ eff ( ε, αζ ) = CVaR β (φ( ε, αζ ; ζ)) (82)
Note however that the above notation is slightly abusive in the sense that internal state variables α ζ are second-stage variables which depend on the uncertainty and cannot be replaced a priori with a finite-set of observable internal state variables which do not depend on the uncertainty anymore. The resulting effective potentials are not practically effective since we still have to track the uncertainty-dependent state variables. This issue is very similar to the homogenization of a heterogenous GSM. Indeed, the mathematical structure of GSM is preserved by up-scaling [START_REF] Suquet | Elements of homogenization for inelastic solid mechanics[END_REF][START_REF] Suquet | Local and global aspects in the mathematical theory of plasticity[END_REF] but this has to be done at the expense of introducing an infinite number of internal variables which represent the local state variables microscopic fields. Here again, the dependence on the local spatial variable is similar to the dependence on the uncertainty.

Nevertheless, taking into account these optimistic effective potentials in lieu of the expected value in the incremental potential [START_REF] Junker | A relaxation approach to modeling the stochastic behavior of elastic materials[END_REF] results in the following definition of the optimistic incremental potential j eff (ε):

j eff (ε) = inf α ζ CVaR β (ψ(ε, α ζ )) + CVaR β (φ (ε -ε n , α ζ -α n,ζ )) ( 83
)
Similarly to section 4.2.3, we consider here the case of finite scenarios i = 1, . . . , N with probabilities p i . The set of internal state variables is therefore a collection {α i } of N state variables. In this case, the effective potential can be written as:

j eff (ε) = min λ,η,α i λ + 1 1 -β N i=1 ψ(ε, α i ; ζ i ) -λ + + η + 1 1 -β N j=1 φ ε -ε n , α j -α n,j ; ζ j -η + (84)
We remark in particular that the definition of this optimistic effective potential introduces only two additional scalar auxiliary variables λ and η arising from the CVaR convex representation (47).

The optimality condition with respect to λ and η results in the characterization of two sets of active scenarios A ψ and A φ associated with both potentials:

1 = 1 1 -β N i=1 p i Heav (ψ(ε, α i ; ζ i ) -λ) = 1 1 -β i∈A ψ p i (85) 1 = 1 1 -β N j=1 p j Heav φ(ε -ε n , α j -α n,j ; ζ j ) -η = 1 1 -β j∈A φ p j (86)
Note that the sets of active scenarios A ψ and A φ are not necessarily the same since some realizations can represent materials associated with a large free energy and a low dissipation potential, and conversely. One realization can therefore be active with respect to the free energy but inactive with respect to the dissipation potential.

The optimality condition with respect to the total stress σ results in:

σ = 1 1 -β i∈A ψ p i ∂ψ ∂ε (ε, α i ; ζ i ) + 1 1 -β j∈A φ p j ∂φ ∂ ε (ε -ε n , α j -α j,n ; ζ j ) (87) = 1 1 -β i∈A ψ p i σ nd i + 1 1 -β j∈A φ p j σ d j (88)
where we see that the total stress is still the sum of a non-dissipative stress σ nd and a dissipative stress σ d . The non-dissipative stress is obtained as a weighted sum of nondissipative stresses σ nd i associated with the active scenarios i ∈ A ψ whereas the dissipative stress is obtained as a weighted sum of dissipative stresses σ d j associated with the active scenarios j ∈ A φ . Finally, the complementary laws characterizing the thermodynamic forces Y nd k , Y d k are obtained from the optimality conditions with respect to the α k as follows:

0 = i∈A ψ δ ik ∂ψ ∂α k (ε, α i ; ζ i ) + j∈A φ δ jk ∂φ ∂ αk (ε -ε n , α j -α j,n ; ζ j ) (89) =          Y nd k + Y d k if k ∈ A ψ , k ∈ A φ Y nd k if k ∈ A ψ , k / ∈ A φ Y d k if k / ∈ A ψ , k ∈ A φ 0 if k / ∈ A ψ , k / ∈ A φ

Pessimistic estimate

The dual pessimistic estimate is obtained by replacing both stochastic free energy potential and stochastic dissipation pseudo-potential with their dCVaR risk measure:

ψ eff (ε, α ζ ) = dCVaR β (ψ(ε, α ζ ; ζ)) (90) φ eff ( ε, αζ ) = dCVaR β (φ( ε, αζ ; ζ)) (91)
Again taking into account these pessimistic effective potentials in lieu of the expected value in the incremental potential [START_REF] Junker | A relaxation approach to modeling the stochastic behavior of elastic materials[END_REF] results in the following definition of the optimistic incremental potential j eff (ε):

j eff (ε) = inf α ζ dCVaR β (ψ(ε, α ζ )) + dCVaR β (φ (ε -ε n , α ζ -α n,ζ )) (92)
In order to have a better feeling of how this dual formulation works, let us now consider the case where ε is non-dissipative. In the case of discrete scenarios and using the convex characterization (78), we obtain:

j eff (ε) = min v i ≥0,r i ,s i , ε i α i r β + s β s.t. N i=1 p i v i = 1 N i=1 p i ε i = ε v i ψ( ε i /v i , α i /v i ; ζ i ) ≤ r i ∀i = 1, . . . , N φ(α i -α n,i ; ζ i ) ≤ s i ∀i = 1, . . . , N (93)
where x β is the composite norm of a vector x ∈ R N defined as:

x β = max N i=1 p i x i , (1 -β) x ∞ (94)
Variables r i and s i have been introduced to express the problem in convex form.

Interestingly, the use of dCVaR naturally introduces an auxiliary field ε which must be equal to ε on average. We therefore recover the formulation the same ideas as in formulation [START_REF] Karapiperis | Fokker-planck linearization for non-gaussian stochastic elastoplastic finite elements[END_REF] where the total strain was considered to be adjustable. Here, we can consider that ε i = ε + e i where e are the total strain deviations introduced earlier.

Moreover, in the case β = 0, the composite norm degenerates to the ∞-norm: x 0 = x ∞ . We therefore recover that the term related to the dissipation pseudo-potential is in fact dCVaR 0 (φ) = sup ζ φ(α ζα n,ζ ) as in the soft-constrained formulation of section 3.4.

Illustration on the 1D elastoplastic behavior

We reconsider the elastoplastic behavior with power-law hardening introduced in section 3.3. Clearly, applying formulation (84) to this case with β = 0 is equivalent to formulation [START_REF] Ortiz | The variational formulation of viscoplastic constitutive updates[END_REF] since the CVaR degenerates to the expected value. However, the above dCVaR formulation (92) for β = 0 is slightly different than formulation [START_REF] Rockafellar | Optimization of conditional value-at-risk[END_REF] in the case β = 0 beacause dCVaR 0 (ψ) = E [ψ] only when ψ is quadratic, otherwise the expression is slightly different. Since hardening is not quadratic in this case, ψ is not quadratic and the expression is slightly different, although very close.

Figure 5 represents the corresponding optimistic and pessimistic effective behavior for various values of the risk-aversion level β. First, we can notice, as expected, that the risk-neutral case yields effective responses very similar to those of Figure 3. Second, we observe that, for β = 0.5 and β = 0.95, we indeed obtain risk-averse effective behaviors either on the optimistic side with the CVaR or on the pessimistic side with the dCVaR. The resulting behavior exhibits a similar elastoplastic behavior characterized by a first elastic stage followed by a hardening stage. With a moderate risk-aversion level β = 0.5, we obtain responses which are in-between extremal and nominal behaviors. On the contrary, for a stronger risk-aversion level β = 0.95, the resulting behavior is much closer to extremal behaviors. Note however that we plotted only 250 realizations over 500 for better readability of the figures. Both estimates therefore capture the effective behavior of both tails of the stress distribution at fixed strain level but without being associated with a single best or worst case realization. This is further confirmed by results of Figure 6 which correspond to the nearly perfectly plastic case. Again, the risk-aversion level enables to yield effective behavior which progressively evolve from being close to the nominal behavior for β = 0 to the extreme upper and lower envelopes of the material responses when β progressively approaches 1. As before, it is worth noting that the optimistic behaviors result in a pseudo-hardening even if the underlying stochastic behavior is always nearly perfectly plastic. This is attributed to the progressive yielding of different realizations associated with different yield stresses when combining the total stress using the partial weighted sum of (57). However, as regards the pessimistic estimate, we have seen that the stress is deterministic and that we have yielding when reaching the yield stress (σ 0 ) β . The resulting behavior is therefore also perfectly plastic. Finally, let us point out that this 1D example is quite simple since each uncertain potential is of the form f (

x; ζ) = k(ζ)g(x)
where k is an uncertain material parameter and g a reference convex potential. Moreover, the three material parameters (E, H and σ 0 ) are independent random variables. In such a case, a simpler way of formulating an effective behavior could consist in replacing each material parameter with its corresponding CVaR (resp. dCVaR) estimate k β (resp. k β ) resulting in an equivalent deterministic behavior with optimistic/pessimistic material parameters. Figure 7 illustrates the difference between such a behavior based on equivalent parameters with the results from Figure 6 based on CVaR/dCVaR effective potentials. We can notice that the pessimistic behaviors are quite close to each other, especially in the perfectly plastic case where they match exactly. As regards optimistic estimates, the differences are larger, certainly owing to the fact that the coupling between internal state variables is absent in the equivalent deterministic behavior. This results in a perfectly plastic behavior in absence of hardening whereas the CVaR variational formulation exhibits a stochastic hardening due to the coupling between state variables.

Risk-averse stochastic programming of the structural response

To ease notations, we will consider in the subsequent developments only the rate-independent case when ε is not a dissipative variable such as in ( 16) and [START_REF] Guilleminot | A probabilistic model for bounded elasticity tensor random fields with application to polycrystalline microstructures[END_REF].

In this section, we investigate how the previous conceps can be applied to compute a risk-averse global response of a structure consisting of a stochastic GSM.

A truss example

We will illustrate the proposed formulations on the 2D truss structure of Figure 8 made of M = 30 members. Each member is of identical cross-section S = 1 and assumed to follow the hardening elastoplastic behavior of ( 30)- [START_REF] Mosek Aps | The MOSEK Optimizer API for Python[END_REF]. As before, the Young modulus E, the hardening modulus H and the yield stress σ 0 for each member are assumed to be independent random variables following a lognormal distribution of mean E [E] = 20, E [H] = 1 and E [σ 0 ] = 1. Unless stated otherwise, the corresponding standard deviations are again taken to be of 20% for the elastic moduli E and H and 10% for the yield stress σ 0 . The reference loading F represented in Figure 8 consists of vertical forces applied to the upper face. We will use a displacement-controlled path-following strategy by driving the associated work-conjugate displacement U = F T u and report the corresponding load factor.

To fix ideas, Figure 9 represents the bundle of structural responses for each of N = 500 realizations. It also reports the resolution of the corresponding deterministic problem when each bar material property is affected an optimistic value (E β , H β , (σ 0 ) β ) or a pessimistic (E β , H β , (σ 0 ) β ). Clearly, we see that simply taking such values is way too optimistic (resp. pessimistic) since it assumes that each bar will take vary large (or very small) material properties simultaneously, which is highly unlikely.

Structural effective response

To extend the previous concepts to the structural effective response, we consider the incremental variational principle [START_REF] De Angelis | Multifield variational principles and computational aspects in rate plasticity[END_REF] where the global potential J is replaced by a certain risk measure R[J].

When using the risk-neutral measure R = E, we have:

R[J](ε) = E Ω j(ε; ζ) dΩ = Ω E [j] (ε) dΩ = Ω R[j](ε) dΩ (95)
In this case, the effective global potential J eff = R[J] is equal to the global potential associated with the effective density j eff = R[j] = E [j]. However, this is no longer true when considering more complex nonlinear risk measures such as CVaR which are not additive. In the more general case, we therefore need to consider the following stochastic program:

u n+1 = arg inf u∈U ad R [J(ε; ζ)] -W ext,n+1 (u) ( 96
)
where J is defined in [START_REF] Diamond | CVXPY: A Python-embedded modeling language for convex optimization[END_REF] and is stochastic due to the stochastic free energy and dissipation potentials. For later use, we also introduce the global potentials:

Ψ(ε, α; ζ) = Ω ψ(ε, α; ζ) dΩ (97) Φ(ε, α; ζ) = Ω φ(ε, α; ζ) dΩ (98)
We assume for simplicity that the loading and boundary conditions, and thus the set of kinematically admissible displacements U ad and W ext,n+1 , are deterministic. The displacement field u should, in theory, change its value according to each realization of the uncertainty ζ. We are however interested in finding an effective response of the structure accounting for the stochastic nature of its constitutive material. For this reason, the displacement field u and the corresponding strain ε = ∇ s u are considered to be firststage variables in (96). As before, the implicit state variables α are considered as a secondstage variables.

Primal risk-averse formulation

Similarly to what has been proposed in section 5, we propose to define a risk-averse effective structural response by considering the CVaR of the global potentials Ψ and Φ. This results in the following stochastic program:

u n+1 = arg inf u∈U ad inf α ζ CVaR β (Ψ) (ε, α ζ ) + CVaR β (Φ) (α ζ -α n,ζ ) -W ext,n+1 (u) (99)
Introducing the convex formulation (47) in lieu of the CVaR operator will result in (99) being a convex stochastic program where α ζ are second-stage variables and u, ε and the two λ parameters associated with each CVaR terms are first-stage variables.

This risk-averse stochastic program can be numerically solved by adopting a Monte-Carlo sampling approximation based on N realizations ζ i of the uncertainty. Introducing auxiliary variables to reformulate the CVaR operator yields the following equivalent deterministic problem:

inf u∈U ad ,α i ,Ψ (i) + ,Φ (i) + ,λ Ψ ,λ Φ λ Ψ + 1 1 -β 1 N N i=1 Ψ (i) + + λ Φ + 1 1 -β 1 N N k=1 Φ (i) + -W ext,n+1 (u) s.t. Ω ψ(ε, α i ; ζ i ) dΩ -λ Ψ ≤ Ψ (i) + ∀i = 1, . . . , N Ω φ(α i -α n,i ; ζ i ) dΩ -λ Φ ≤ Φ (i) + 0 ≤ Ψ (i) + 0 ≤ Φ (i) + (100)
In the case when β = 0, we recover a classical equivalent deterministic problem based on the expected value:

inf u∈U ad ,α i ,Ψ (i) ,Φ (i) 1 N N i=1 Ψ (i) + 1 N N k=1 Φ (i) -W ext,n+1 (u) s.t. Ω ψ(ε, α i ; ζ i ) dΩ ≤ Ψ (i) ∀i = 1, . . . , N Ω φ(α i -α n,i ; ζ i ) dΩ ≤ Φ (i) (101)
In all cases, spatial discretization of the displacement field is performed using the finiteelement method. Moreover, most GSM behaviors, such as the considered elastoplastic behavior ( 30)- [START_REF] Mosek Aps | The MOSEK Optimizer API for Python[END_REF], can be expressed using conic-representable functions, see [START_REF] Bleyer | Applications of conic programming in non-smooth mechanics[END_REF]. The resulting problem therefore belongs to the class of convex conic programs which can be easily formulated using tools such as cvxpy and solved using dedicated interior-point algorithms such as MOSEK [START_REF] Mosek Aps | The MOSEK Optimizer API for Python[END_REF]. Note that the introduction of the CVaR measure only adds two additional variables scalar variables λ Ψ and λ Φ compared to a formulation based on the expectation. It introduces some additional non-linearity in the problem due to the presence of the positive part in (47). The latter can however be considered as very minor compared to the material non-linearity. Obviously, solving the above effective problem is computationally intensive since a large number of realizations must be considered in order for the Monte-Carlo sampling approximation to converge and the problem size scales with N . A large body of works in the literature is therefore devoted to improving the efficiency of computing the corresponding effective behavior using various techniques such as variance reduction, adaptive sampling, scenario reductions, etc. We will however leave this important aspect to future research.

Dual risk-averse formulation

Similarly and using the developments of section 4.3, the dual risk-averse formulation which is aimed to yield pessimistic estimates of the structural response is obtained by considering the dCVaR measure of the global potentials Ψ and Φ as follows:

u n+1 = arg inf u∈U ad inf α ζ dCVaR β (Ψ) (ε, α ζ ) + dCVaR β (Φ) (α ζ -α n,ζ ) -W ext,n+1 (u) (102)
In particular, since Φ is a homogeneous potential, we have:

dCVaR β (Φ) (α ζ -α n,ζ ) = max E [Φ(α ζ -α n,ζ )] , (1 -β) sup ζ {Φ(α ζ -α n,ζ )} (103)
Besides, when computing the corresponding dual problem, we obtain:

σ n+1 = arg inf σ,Y ζ dCVaR β (Ψ) * (σ, -Y ζ ) + E [Y ζ • α ζ,n ] s.t. div σ + f n+1 = 0 in Ω σn = T n+1 on ∂Ω N CVaR β (Φ • ) (Y ζ ) ≤ 1 (104)
in which we used result (80).

In the present case, we have

Φ * (Y ζ ) = Ω φ * (Y ζ ) dΩ (105)
where each φ * (Y ζ ) defines a convex yield domain in the space of thermodynamic forces Y ζ which we can represent by its gauge function g = φ • as g(Y ζ ) ≤ 1. We therefore have:

Φ * (Y ζ ) = 0 if g(Y ζ ) ≤ 1 ∀x ∈ Ω +∞ otherwise (106) Φ • (Y ζ ) = sup Ω g(Y ζ ) (107)
We see in particular that the corresponding risk-averse yield criterion will involve the CVaR of the maximum value of the gauge function g attained on Ω. In the risk-neutral case β = 0, this reduces to E [g(Y )] ≤ 1 which coincides with the yield constraint proposed in [START_REF] Rockafellar | Conditional value-at-risk for general loss distributions[END_REF].

Finally, using (78), we have:

σ n+1 = arg inf σ,Y ζ ,z ζ sup ζ Ω z ζ ψ * (σ/z ζ , -Y ζ /z ζ ; ζ) dΩ + E [Y ζ • α ζ,n ] s.t. div σ + f n+1 = 0 in Ω σn = T n+1 on ∂Ω N CVaR β sup Ω g(Y ζ ) ≤ 1 CVaR β (z ζ ) ≤ 1 z ζ ≥ 0 (108)
As before, both problems (102) and (108) are convex stochastic problems which can be approximated using a Monte-Carlo sampling approximation and discretized in space using a FE formulation. In practice, we use the same FE approximation and Monte-Carlo sampling to approximate (99) and (102). We simply change the CVaR risk measure to the dCVaR measure between both problems. As (100), the equivalent deterministic problem resulting from the discretized version of ( 102) is a conic optimization problem which is formulated and solved in the same fashion.

Illustrative application

We first analyze the effective response of the truss structure in the risk-neutral β = 0. Figure 10 shows the corresponding response obtained when solving the risk-neutral stochastic programs (99) and (102). As in Figure 5a, a small duality gap between both formulations is obtained since CVaR and dCVaR yield slightly different effective measures. However, both formulations yield a good approximation of the average structural response when comparing with the empirical distribution of load-displacement curves (colored symbols).

We then consider the resolution of both problems in the risk-averse case with β = 0.95 with a sample size of N = 100. Figure 11 shows the obtained optimistic and pessimistic responses when considering different types of uncertainty.

On Figure 11a, we assume that uncertainty affects only the Young modulus. In this case, the bundle of responses is quite thin and the behavior has much less variance in the final hardening stage of the response. Both formulations are however able to correctly account for a best-case and a pessimistic estimate of the response. On Figure 11b, uncertainty affects only the hardening modulus which results in a much wider spread of the responses after a first deterministic elastic stage. Again, both risk-averse estimates are of very good quality compared to the empirical distribution. They coincide in the first deterministic elastic stage and then exhibit a hardening behavior with a different slope, following that of the response distribution. A similar observation can also be made concerning Figure 11c which considers only a yield strength uncertainty. Here, the difference is that both risk-averse estimates show a similar hardening response but the onset of plasticity occurs at different load levels, closely matching that of the empirical distribution. Finally, Figure 11d is the most interesting since it considers a simultaneous uncertainty on all three mechanical parameters. Again, the agreement is also excellent in this case. As a conclusion, these results show that the characteristics of the risk-averse estimates which have been observed on the material point response in section 5 transpose also to the global structural behavior.

As regards this last case of uncertainty on all mechanical parameters, Figure 12 displays the deformed configurations of the truss structure for two different load levels corresponding to an imposed displacement U = 2 and U = 4. The effective responses obtained from the resolution of both risk-averse formulations (99) and (102) have been compared against 25 of the uncertain realizations. First, we can see that the spread of nodal displacements can be quite large for some nodes and such a spread increases with the load level. Let us point out that, since we have a displacement control, all deformed configuration have the same equivalent displacement U which corresponds to a weighted-average of the vertical displacement of the top surface. However, each deformed configuration will be associated with a very different state of internal forces. This is particularly highlighted by Fig. 12b which shows that the pessimistic effective response using dCVaR (in red) is subject to large straining, see the vertical members on both supports. On the contrary, the optimistic estimate (in blue) seems to be much stiffer for the same imposed displacement.

Figure 13 investigates the influence on the risk-aversion level β on the obtained riskaverse estimates. As expected, going from 0 to 1 yields a stiffer and stronger response for the CVaR optimistic estimate. Conversely, a softer and weaker response is obtained with the dCVaR pessimistic estimate. One can see that going from β = 0.95 to β = 0.995 does not have too much of an influence here. However, we must moderate this observation since we used only a small sample size in order to limit the computational cost. To further assess this aspect, Figure 13b investigates the influence of the sample size N for β = 0.95. It seems that reasonable estimates are already obtained with N = 50. In particular, it does not seem that one formulation is more sensitive than the other to the chosen sample size.

Finally, we further assess the efficiency of the proposed formulations on the much more challenging problem of a cyclic loading. In Figure 14, the displacement amplitude is varied as U : 0 → 3 → 0 → 5 → 0. We can first observe that the cyclic behavior of the truss is correctly reproduced. Moreover, in the plastic evolution phases, both formulation indeed produce the correct pessimistic and best-case estimates, even after load direction has flipped signs. However, we can notice that in the elastic unloading/reloading phases, there is a crossover between the pessimistic and best-case estimates. To explain this, let us consider the first unloading from U = 3 to 0. Even if the best-case estimate is stiffer than the pessimistic one, at U = 3, the best-case response is associated with smaller plastic strain levels than the pessimistic. Elastic unloading therefore occurs at different levels of plastic strain and this strain difference is then later compensated by stiffer/stronger response in the unloading stage.

Conclusions and perspectives

In this work, we have used the tools of stochastic programming to propose a new approach to the treatment of stochastic behaviors in the framework of Generalized Standard Materials. The proposed methodology tackles the case of convex free-energy and dissipation pseudopotential which both depend on uncertain parameters in a very general manner.

In first part, we have proposed a formulation describing the effective material behavior by considering the expected value for both potentials. We have also proposed an alternate dual formulation which must consider the total strain to be an adjustable variable with a fixed expected value associated with the observable total strain of the effective behavior. Moreover, we have shown, for such a dual formulation, that the dissipation pseudo-potential should be replaced with its supremum over all realizations. When inspecting the resulting formulation for an elastoplastic behavior, we have seen that this formulation is associated with satisfying the yield condition in an average sense.

In a second part, we have introduced the notion of coherent risk measures which have been successfully introduced in the financial mathematics community. The Conditional Value-at-Risk (CVaR) is one example of a coherent risk measure which is able to account for the tail behavior of a stochastic potential while preserving convexity. Applying this concept to both potentials results in an optimistic estimate of the effective material behavior.

To extend this concept in order to obtain a pessimistic estimate, we had to introduce a new risk measure, the dual Conditional Value-at-Risk (dCVaR), which is defined using the polar transform of positive convex functions. To the best of our knowledge, we did not find any reference of such a definition in the mathematical literature. The use of the dCVaR is consistent with the dual formulation which we studied in the first part and is indeed able to produce a pessimistic estimate of the material response.

One of the key point in the success of such formulations is that convexity of both potentials is preserved. By construction, our optimistic and pessimistic effective behaviors both satisfy the definition of a GSM which ensures consistent thermodynamical properties. Moreover, the use of an incremental variational principle accounting for the history of internal state variables, as opposed to a deformation theory of plasticity for instance, enables to treat cyclic loadings or non-proportional loadings without any additional difficulty.

Finally, even if our contribution does not address the modeling of constitutive laws since we assume given uncertain parameters with given probability distribution, our framework provides a way to perform up-scaling and propagate uncertainty from one spatial level to an other while preserving the GSM structure.

The treatment of stochastic material behavior within the tools of stochastic programming paves the way to numerous developments which could benefit from the developments of this mathematical community. One key challenge is however the numerical cost associated with an extensive number of internal state variables when using Monte-Carlo sampling. Although the equilibrium equations do not depend on the sample size since we are computing an effective structural response, the nonlinear constitutive law must be solved simultaneously for N samples, resulting in significant computational expense. To alleviate this issue, one potential strategy is to identify the set of active scenarios, for instance during global Newton iterations, which would significantly reduce the cost of evaluating the constitutive law. Indeed, when using the CVaR with a large risk-aversion level, we should expect only a small fractions of scenarios to be active at a given increment. Additionally, using linear or non-linear decision rules [START_REF] Georghiou | Generalized decision rule approximations for stochastic programming via liftings[END_REF] for second-stage variables instead of sampling could be a viable alternative to reduce this cost. Finally, we should also mention current works developing Domain Specific Languages which aim at simplifying the formulation of stochastic programs [START_REF] Ali | Disciplined convex stochastic programming: A new framework for stochastic optimization[END_REF].

Figure 1 :

 1 Figure1: A stochastic material or structural behavior and its average, optimistic and pessimistic effective behavior. Dashed blue/red lines represent slightly less optimistic/pessimistic estimates than the solid ones.

  Figure2arepresents the corresponding effective behavior for both formulations in the case where E [E] = 20, E [H] = 1 and E [σ 0 ] = 1 and a 20% standard deviation for both E and H and 10% for σ 0 . The effective behavior is obtained by discretizing the imposed total strain in 30 increments. Each increment is obtained by solving either formulation (34) (in blue) or[START_REF] Rockafellar | Convex analysis[END_REF] (in red). In both cases, expectations have been computed using a Monte-Carlo sampling approximation with a sample size N = 500. Minimization over the N different state variables is performed using the cvxpy package[START_REF] Agrawal | A rewriting system for convex optimization problems[END_REF][START_REF] Diamond | CVXPY: A Python-embedded modeling language for convex optimization[END_REF].First, it can be seen that both approaches lead to a very similar behavior, the second formulation[START_REF] Rockafellar | Convex analysis[END_REF] being slightly softer than the first one due to the presence of additional total strain variables. Interestingly, the obtained effective behavior is close to the nominal deterministic behavior except that it exhibits an earlier and progressive onset of plasticity since it accounts for a distribution of yield stress with different values. This is further confirmed in Figure2bwhich corresponds to an almost perfectly-plastic behavior with E [H] = 0.01. We can observe that the effective behavior corresponding to (34) admits an initial hardening phase before reaching a perfectly plastic plateau. A similar behavior is obtained from the approach proposed in[START_REF] Einav | A thermomechanical framework of plasticity based on probabilistic micromechanics[END_REF]. Moreover, in the ultimate perfectly plastic regime, the elastic term vanishes so that ε p ζ = ε. One can then see that[START_REF] Ortiz | The variational formulation of viscoplastic constitutive updates[END_REF] gives approximately j eff (ε) ≈ E [σ 0 ] |εε n | in the asymptotic regime of the perfectly

Figure 2 :

 2 Figure2: Effective behavior of a stochastic hardening elastoplastic material. The black dashed line corresponds to the nominal deterministic case where all parameters take their average value. Thin black lines correspond to 100 independent realizations of the stochastic behavior.

Figure 3 :

 3 Figure 3: Effective behavior of a stochastic hardening elastoplastic material using the soft constraint formulations (39).

Figure 4 :

 4 Figure 4: Illustration of the CVaR definition

  β = 0.95 (strong risk aversion)

Figure 5 :

 5 Figure 5: Optimistic (CVaR) and pessimistic (dCVaR) effective behaviors of a stochastic elastoplastic behavior (hardening case with E [H] = 1) for various values of the risk aversion level β.

Figure 6 :

 6 Figure 6: Optimistic (CVaR) and pessimistic (dCVaR) effective behaviors of a stochastic elastoplastic behavior (nearly perfectly plastic case with E [H] = 0.01) for various values of the risk aversion level β.

  (a) hardening case (E [H] = 1) (b) nearly perfectly plastic case (E [H] = 0.01)

Figure 7 :

 7 Figure 7: Comparison with elastoplastic behaviors associated with optimistic/pessimistic values of the material parameters (β = 0.95).

Figure 8 :

 8 Figure 8: Truss structure with reference loading consisting of downwards vertical forces of intensity 1 on the middle nodes and 0.5 on the end nodes.

Figure 9 :

 9 Figure 9: Stochastic structural response for 500 realizations and deterministic responses using optimistic and pessimistic values for each material property.

Figure 10 :

 10 Figure 10: Risk-neutral β = 0 structural responses using the primal CVaR measure (99) and the dual dCVaR measure (102). Colored symbols indicate the empirical distribution of the structural response for fixed values of the imposed displacement. Color and size of the marker denote the frequency of the observed response.

Figure 11 :

 11 Figure 11: Risk-averse estimates for various types of uncertainty

4 Figure 12 :

 412 Figure 12: Effective deformed configurations: optimistic CVaR (blue), pessimistic dCVaR (red) and 25 random realizations

Figure 13 :

 13 Figure 13: Influence of risk-aversion level and number of sampling realizations on risk-averse estimates

Figure 14 :

 14 Figure 14: Risk-averse structural responses for a cyclic loading.

Acknowledgements

The author would like to thank Vincent Leclère, Djimédo Kondo and R. Tyrell Rockafellar for discussing some aspects of this work.