I Hlaváček
email: hlavacek@math.cas.cz

A A Novotny
email: novotny@lncc.br

J Soko Lowski
email: jan.sokolowski@iecn.u-nancy.fr

A Żochowski

Energy change in elastic solids due to a spherical or circular cavity, considering uncertain input data

In the paper we consider topological derivative of shape functionals for elasticity, which is used to derive the worst and also the maximum range scenarios for behavior of elastic body in case of uncertain material parameters and loading. It turns out that both problems are connected, because the criteria describing this behavior have form of functionals depending on topological derivative of elastic energy. Therefore in the first part we describe the methodology of computing the topological derivative with some new additional conditions for shape functionals depending on stress. For the sake of fulness of presentation the explicit formulas for stress distribution around cavities are provided.

Introduction

In the paper we consider topological derivative of shape functionals for elasticity, which is used to derive the worst and also the maximum range scenarios for behavior of elastic body in case of uncertain material parameters.

It turns out that both problems are connected, because the criteria describing this behavior have form of functionals depending on topological derivative of elastic energy.

Therefore in the first part we describe the methodology of computing the topological derivative with some new additional conditions for shape functionals depending on stress.

For the sake of fulness of presentation the explicit formulas for stress distribution around cavities are provided.

Topological Derivative

The topological derivative T Ω of a shape functional J (Ω) is introduced in [START_REF] Soko Lowski | On topological derivative in shape optimization[END_REF] in order to characterize the infinitesimal variation of J (Ω) with respect to the infinitesimal variation of the topology of the domain Ω. The topological derivative allows us to derive the new optimality condition for the shape optimization problem:

J (Ω * ) = inf Ω J (Ω) .
The optimal domain Ω * is characterized by the first order condition [START_REF] Soko | Introduction to Shape Optimization[END_REF] defined on the boundary of the optimal domain Ω * , dJ(Ω * ; V ) ≥ 0 for all admissible vector fields V , and by the following optimality condition defined in the interior of the domain Ω * :

T Ω * (x) ≥ 0 in Ω * .
The other use of the topological derivative is connected with approximating the influence of the holes in the domain on the values of integral functionals of solutions, what allows us to solve a class of shape inverse problems.

In general terms the notion of the topological derivative (TD) has the following meaning. Assume that Ω ⊂ IR N is an open set and that there is given a shape functional J : Ω \ K → IR for any compact subset K ⊂ Ω. We denote by B ρ (x), x ∈ Ω, the ball of radius ρ > 0, B ρ (x) = {y ∈ IR N | y -x < ρ}, B ρ (x) is the closure of B ρ (x), and assume that there exists the following limit

T(x) = lim ρ↓0 J (Ω \ B ρ (x)) -J (Ω) |B ρ (x)|
which can be defined in an equivalent way by

T(x) = lim ρ↓0 J (Ω \ B ρ (x)) -J (Ω) ρ N
The function T(x), x ∈ Ω, is called the topological derivative of J (Ω), and provides the information on the infinitesimal variation of the shape functional J if a small hole is created at x ∈ Ω. This definition is suitable for Neumann-type boundary conditions on ∂B ρ . In several cases this characterization is constructive, i.e. TD can be evaluated for shape functionals depending on solutions of partial differential equations defined in the domain Ω.

For instance, TD may be computed for the 3D elliptic Laplace type equation, as well as for extremal values of cost functionals for a class of optimal control problems. All these examples have one common feature: the expression for TD may be calculated in the closed functional form.

As we shall see below, the 3D elasticity case is more difficult, since it requires evaluation of integrals on the unit sphere with the integrands which can be computed at any point, but the resulting functions have no explicit functional form. In the particular case of energy functional we obtain the closed formula. In section 5 we compare the results of the present paper with the formulae for 2D elasticity.

The main contribution of the present paper is the procedure for computations of the topological derivatives of shape functionals depending on the solutions of 3D elasticity systems. Therefore it constitutes an essential extension of the results given in [START_REF] Soko Lowski | On topological derivative in shape optimization[END_REF] for the 2D case.

Problem setting for elasticity systems

We introduce elasticity system in the form convenient for the evaluation of topological derivatives. Let us consider the elasticity equations in IR N , where N = 2 for 2D and

N = 3 for 3D,    div σ(u) = 0 in Ω u = g on Γ D σ(u)n = T on Γ N (1)
and the same system in the domain with the spherical cavity

B ρ (x 0 ) ⊂ Ω centered at x 0 ∈ Ω, Ω ρ = Ω \ B ρ (x 0 ),        div σ ρ (u ρ ) = 0 in Ω ρ u ρ = g on Γ D σ ρ (u ρ )n = T on Γ N σ ρ (u ρ )n = 0 on ∂B ρ (x 0 ) (2)
where n is the unit outward normal vector on ∂Ω ρ = ∂Ω ∪ ∂B ρ (x 0 ). Assuming that 0 ∈ Ω, we can consider the case x 0 = 0. Here u and u ρ denote the displacement vectors fields, g is a given displacement on the fixed part Γ D of the boundary, T is a traction prescribed on the loaded part Γ N of the boundary. In addition, σ is the Cauchy stress tensor given, for ξ = u (eq. 1) or ξ = u ρ (eq. 2), by

σ(ξ) = D∇ s ξ , (3)
where ∇ s (ξ) is the symmetric part of the gradient of vector field ξ, that is

∇ s (ξ) = 1 2 ∇ξ + ∇ξ T , (4)
and D is the elasticity tensor,

D = 2µII + λ (I ⊗ I) , (5)
with µ = E 2(1 + ν) , λ = νE (1 + ν)(1 -2ν) and λ = λ * = νE 1 -ν 2 (6)
being E the Young's modulus, ν the Poisson's ratio and λ * the particular case for plane stress.

In addition, I and II respectively are the second and fourth order identity tensors. Thus, the inverse of D is

D -1 = 1 2µ II - λ 2µ + N λ (I ⊗ I) ,
The first shape functional under consideration depends on the displacement field,

J u (ρ) = Ωρ F (u ρ ) dΩ , F (u ρ ) = (Hu ρ • u ρ ) p , ( 7
)
where F is a C 2 function, p ≥ 2 is an integer. It is also useful for further applications in the framework of elasticity to introduce the yield functional of the form

J σ (ρ) = Ωρ Sσ(u ρ ) • σ(u ρ ) dΩ , (8)
where S is an isotropic fourth-order tensor. Isotropicity means here, that S may be expressed as follows

S = 2mII + l (I ⊗ I) ,
where l, m are real constants. Their values may vary for particular yield criteria. The following assumption assures, that J u , J σ are well defined for solutions of the elasticity system.

(A) The domain Ω has piecewise smooth boundary, which may have reentrant corners with α < 2π created by the intersection of two planes. In addition, g, T must be compatible with u ∈ H 1 (Ω; IR N ).

The interior regularity of u in Ω is determined by the regularity of the right hand side of the elasticity system. For simplicity the following notation is used for functional spaces,

H 1 g (Ω ρ ) = {ψ ∈ [H 1 (Ω ρ )] N | ψ = g on Γ D }, H 1 Γ D (Ω ρ ) = {ψ ∈ [H 1 (Ω ρ )] N | ψ = 0 on Γ D }, H 1 Γ D (Ω) = {ψ ∈ [H 1 (Ω)] N | ψ = 0 on Γ D }
, here we use the convention that eg., in our notation H 1 g (Ω ρ ) stands for the Sobolev space of vector functions [H 1 g (Ω ρ )] N . The weak solutions to the elasticity systems are defined in the standard way.

Find

u ρ ∈ H 1 g (Ω ρ ) such that, for every φ ∈ H 1 Γ D (Ω ρ ), Ωρ D∇ s u ρ • ∇ s φ dΩ = Γ N T • φ dS (9)
We introduce the adjoint state equations in order to simplify the form of shape derivatives of functionals J u , J σ . For the functional J u they take on the form:

Find w ρ ∈ H 1 Γ D (Ω ρ ) such that, for every φ ∈ H 1 Γ D (Ω ρ ), Ωρ D∇ s w ρ • ∇ s φ dΩ = - Ωρ F u (u ρ ) • φ dΩ, (10)
whose Euler-Lagrange equation reads

       div σ ρ (w ρ ) = F u (u ρ ) in Ω ρ w ρ = 0 on Γ D σ ρ (w ρ )n = 0 on Γ N σ ρ (w ρ )n = 0 on ∂B ρ (x 0 ) , ( 11
) while v ρ ∈ H 1 Γ D (Ω ρ )
is the adjoint state for J σ and satisfies for all test functions φ ∈ H 1 Γ D (Ω) the following integral identity:

Ωρ D∇ s v ρ • ∇ s φ dΩ = -2 Ωρ DSσ(u ρ ) • ∇ s φ dΩ. ( 12
)
which associated Euler-Lagrange equation becomes

       div σ ρ (v ρ ) = -2div (DSσ ρ (u ρ )) in Ω ρ v ρ = 0 on Γ D σ ρ (v ρ )n = -2DSσ ρ (u ρ )n on Γ N σ ρ (v ρ )n = -2DSσ ρ (u ρ )n on S ρ (x 0 ) = ∂B ρ (x 0 ) . ( 13
)
Remark 1 We observe that DS can be written as

DS = 4µmII + γ (I ⊗ I) (14)
where

γ = λlN + 2 (λm + µl) (15)
Thus, when γ = 0, the boundary condition on ∂B ρ (x 0 ) in eq. ( 13) becomes homogeneous and the yield criteria must satisfy the constraint

m l = - µ λ + N 2 , ( 16
)
which is naturally satisfied for the energy shape functional, for instance. In fact, in this particular case, tensor S is given by

S = 1 2 D -1 ⇒ γ = 0 and 2m + l = 1 2E , ( 17
)
which implies that the adjoint solution associated to J σ can be explicitly obtained such that v ρ = -(u ρ -g).

Main result

We shall define the topological derivative of the functionals J u , J σ at the point x 0 as:

T J u (x 0 ) = lim ρ↓0 dJ u (ρ) d(|B ρ (x 0 )|) , (18)
T J σ (x 0 ) = lim ρ↓0 dJ σ (ρ) d(|B ρ (x 0 )|) . ( 19
)
Now we may formulate the following result, giving the constructive method for computing the topological derivatives:

Theorem 1 Assume that (A) is satisfied, then T J u (x 0 ) = - 1 2(N -1)π [ 2(N -1)πF (u) + Ψ(D -1 ; σ(u), σ(w))] x=x 0 , (20)
T J σ (x 0 ) = - 1 2(N -1)π [ Ψ(S; σ(u), σ(u)) + Ψ(D -1 ; σ(u), σ(v))] x=x 0 , (21)
where w, v ∈ H 1 Γ D (Ω) are adjoint variables satisfying the integral identities ( 10) and ( 12) for ρ = 0, i.e. in the whole domain Ω instead of Ω ρ , that is

Ω D∇ s w • ∇ s φ dΩ = - Ω F u (u) • φ dΩ. ( 22
) Ω D∇ s v • ∇ s φ dΩ = -2 Ω DSσ(u) • ∇ s φ dΩ. ( 23
)
for all test functions φ ∈ H 1 Γ D (Ω).

Some of the terms in (20), (21) require explanation. The function Ψ is defined as an integral over the unit sphere S 1 (0) = {x ∈ IR N | x = 1} of the following functions:

Ψ(S; σ(u(x 0 )), σ(u(x 0 ))) = S 1 (0) Sσ ∞ (u(x 0 ); x) • σ ∞ (u(x 0 ); x) dS (24) Ψ(D -1 ; σ(u(x 0 )), σ(v(x 0 ))) = S 1 (0) σ ∞ (u(x 0 ); x) • D -1 σ ∞ (v(x 0 ); x) dS (25) Ψ(D -1 ; σ(u(x 0 )), σ(w(x 0 ))) = S 1 (0) σ ∞ (u(x 0 ); x) • D -1 σ ∞ (w(x 0 ); x) dS (26)
The symbol σ ∞ (u(x 0 ); x) denotes the stresses for the solution of the elasticity system (2) in the infinite domain IR N \ B 1 (0) with the following boundary conditions:

• no tractions are applied on the surface of the ball, S 1 (0) = ∂B 1 (0);

• the stresses σ ∞ (u(x 0 ); x) tend to the constant value σ(u(x 0 )) as x → ∞.

In this notation σ ∞ (u(x 0 ); x) is a function of space variables depending on the functional parameter u(x 0 ), while σ(u(x 0 )) is a value of the stress tensor computed in the point x 0 for the solution u. The dependence between them results from the boundary condition at infinity listed above. The method for obtaining such solutions (and u ∞ ), based on [START_REF] Kachanov | Handbook of Elasticity Solutions[END_REF], is discussed in the next section.

In order to derive the above formulae (20), (21) we calculate the derivatives of the functional J u (ρ) with respect to the parameter ρ, which determines the size of the hole B ρ (x 0 ), by using the material derivative method [START_REF] Soko | Introduction to Shape Optimization[END_REF]. Then we pass to the limit ρ ↓ 0 using the asymptotic expansions of u ρ with respect to ρ. For the functional J u the shape derivative with respect to ρ is given by

J u (ρ) = Ωρ F u (u ρ ) • u ρ dΩ - Sρ(x 0 ) F (u ρ ) dS, (27)
and in the same way for the state equation:

Ωρ D∇ s u ρ • ∇ s φ dΩ - Sρ(x 0 ) D∇ s u ρ • ∇ s φ dS = 0, ( 28
)
where u ρ is the shape derivative, i.e. the derivative of u ρ with respect to ρ, [START_REF] Soko | Introduction to Shape Optimization[END_REF].

After substitution of the test functions φ = w ρ in the derivative of the state equation, φ = u ρ in the adjoint equation, we get

J u (ρ) = - Sρ(x 0 ) [F (u ρ ) + D∇ s u ρ • ∇ s w ρ ] dS = - Sρ(x 0 ) [F (u ρ ) + σ(u ρ ) • D -1 σ(w ρ )] dS, (29)
and similarly for J σ

J σ (ρ) = - Sρ(x 0 ) [Sσ(u ρ ) • σ(u ρ ) + D∇ s u ρ • ∇ s v ρ ] dS = - Sρ(x 0 ) [Sσ(u ρ ) • σ(u ρ ) + σ(u ρ ) • D -1 σ(v ρ )] dS. (30)
Observe, that both matrices D -1 and S are isotropic, and therefore the corresponding bilinear forms in terms of stresses are invariant with respect to the rotations of the coordinate system. Now we exploit the fact, that

dJ u (ρ) d(|B ρ (x 0 )|) = 1 2(N -1)πρ N -1 dJ u dρ ,
and use the existence of the asymptotic expansions for u ρ in the neighborhood of B ρ (x 0 ), namely

u ρ = u(x 0 ) + u ∞ + O(ρ 2 ). (31)
In addition, u ∞ is proportional to ρ, u ∞ IR N = O(ρ), on the surface S ρ (x 0 ) of the ball. The expansion of σ(u ρ ) corresponding to (31) has the form

σ(u ρ ) = σ ∞ (u(x 0 ); x) + O(ρ). (32)
It may be proved, that w ρ and v ρ have similar expansions.

Using the formulae (31),(32) we may justify the following passages to the limit:

lim ρ↓0 1 ρ N -1 Sρ(x 0 ) σ(u ρ ) • D -1 σ(v ρ ) dS = Ψ(D -1 ; σ(u(x 0 )), σ(v(x 0 ))), lim ρ↓0 1 ρ N -1 Sρ(x 0 ) σ(u ρ ) • D -1 σ(w ρ ) dS = Ψ(D -1 ; σ(u(x 0 )), σ(w(x 0 ))), lim ρ↓0 1 ρ N -1 Sρ(x 0 ) Sσ(u ρ ) • σ(u ρ ) dS = Ψ(S; σ(u(x 0 )), σ(u(x 0 ))), lim ρ↓0 1 ρ N -1 Sρ(x 0 ) F (u ρ ) dS = 2(N -1)πF (u(x 0 )).
This completes the proof of the theorem.

The main difficulty lies in the computation of the values of the functions denoted above as Ψ(S; σ(u(x 0 )), σ(u(x 0 ))), Ψ(D -1 ; σ(u(x 0 )), σ(w(x 0 ))) and Ψ(D -1 ; σ(u(x 0 )), σ(v(x 0 ))), which, in general, is difficult to obtain in the closed form, in contrast with the two dimensional case. Therefore we can approximate them using numerical quadrature. It is possible, because we may calculate the values of integrands at any point on the sphere. This makes the computations more involved, but does not increase the numerical complexity in comparison to evaluating single closed form expression.

Remark 2

The tensor S in the definition of J σ may, in fact, be arbitrary, not only isotropic. However, it is difficult to imagine such a need for the isotropic material. Anyway, in the general case, we would have to transform S according to the known rules for the fourth order tensor, connected with the rotation of the reference frame.

Topological derivatives in 3D elasticity

The shape functionals J u , J σ are defined in the same way as presented in section 2.2 with the exception, that J σ is now the energy stored in a 3D elastic body (see remark 1). The weak solutions to the elasticity system as well as adjoint equations are defined also analogously to the section 2.2. Then, considering the expansions presented in Appendix A.2, we may state the following result [START_REF] Novotny | Topological Sensitivity Analysis for three-dimensional linear elasticity problems[END_REF] (see also [START_REF] Garreau | The Topological Asymptotic for PDE Systems: The Elasticity Case[END_REF]):

Theorem 2 The expressions for the topological derivatives of the functionals J u , J σ have the form

T J u (x 0 ) = -F (u) + 3 2E 1 -ν 7 -5ν (10(1 + ν)σ(u) • σ(w) -(1 + 5ν)trσ(u)trσ(w)) x=x 0 , (33)
T J σ (x 0 ) = 3 4E 1 -ν 7 -5ν 10(1 + ν)σ(u) • σ(u) -(1 + 5ν)(trσ(u)) 2 x=x 0 . ( 34
)

Topological derivatives in 2D elasticity

For the convenience of the reader we recall here the results derived in [START_REF] Soko Lowski | On topological derivative in shape optimization[END_REF] for the 2D case. The principal stresses associated with the displacement field u are denoted by σ I (u), σ II (u), the trace of the stress tensor σ(u) is denoted by trσ(u) = σ I (u) + σ II (u). The shape functionals J u , J σ are defined in the same way as presented in section 2.2, with the tensor S isotropic (that is similar to D). The weak solutions to the elasticity system as well as adjoint equations are defined also analogously to the section 2.2. Then, from the expansions presented in Appendix A.1, we may formulate the following result [START_REF] Soko Lowski | On topological derivative in shape optimization[END_REF]:

Theorem 3
The expressions for the topological derivatives of the functionals J u ,J σ have the form

T J u (x 0 ) = -F (u) + 1 E (a u a w + 2b u b w cos 2δ) x=x 0 = -F (u) + 1 E (4σ(u) • σ(w) -trσ(u)trσ(w)) x=x 0 (35) T J σ (x 0 ) = -η(a 2 u + 2b 2 u ) + 1 E (a u a v + 2b u b v cos 2δ) x=x 0 = -η(4σ(u) • σ(u) -(trσ(u)) 2 ) + 1 E (4σ(u) • σ(v) -trσ(u)trσ(v)) x=x 0 (36)
Some of the terms in ( 35), (36) require explanation. According to eq. ( 15) for N = 2, constant η is given by

η = l + 2 m + γ ν E . (37)
Furthermore, we denote

a u = σ I (u) + σ II (u), b u = σ I (u) -σ II (u), a w = σ I (w) + σ II (w), b w = σ I (w) -σ II (w), a v = σ I (v) + σ II (v), b v = σ I (v) -σ II (v). (38)
Finally, the angle δ denotes the angle between principal stress directions for displacement fields u and w in (35), and for displacement fields u and v in (36).

Remark 3 For the energy stored in a 2D elastic body, tensor S is given by eq. (17), γ = 0 and η = 1/(2E). Thus, since v = -(u -g), we obtain the following well-known result

T J σ (x 0 ) = 1 2E 4σ(u) • σ(u) -(trσ(u)) 2 x=x 0 (39)
Compare these expressions to the 3D case. Their simplicity comes from the fact, that on the plane the rotation of one coordinate system with respect to the other is defined by the single value of the angle (here δ). This is a purely 2D phenomenon and it makes the explicit computations possible.

Uncertain input data

In reality, the values of input data(loading, material parameters) are guaranteed only in some given intervals. One of the simplest remedy is to apply the worst scenario or maximum range scenario method [START_REF] Hlaváček | Uncertain Input Data Problems and the Worst Scenario Method[END_REF]. In what follows, we present the methods for the traction problem (1) with ∂Ω = Γ N and the criterion corresponding to the topological derivatives (34) or (39), respectively.

Traction problem in 3D-elasticity

Let us consider a bounded domain Ω ⊂ R 3 with Lipschitz boundary ∂Ω ≡ Γ, occupied by a homogeneous and isotropic elastic body. Let the body be loaded by surface forces T ∈ [L ∞ (Γ)] 3 and the body forces be zero.

We introduce sets of admissible uncertain input data as follows :

(i) Lamé coefficients λ ∈ U λ ad = λ, λ , 0 ≤ λ < λ < ∞, µ ∈ U µ ad = µ, µ , 0 < µ < µ < ∞;
(ii) surface loading forces

T i ∈ U T i ad = τ ∈ L ∞ (Γ) : τ Γp ∈ C (0),1 Γ p , |τ | ≤ C 1 , |∂τ /∂s j | ≤ C 2 a.e. on Γ, j = 1, 2 , where Γ = P p=1 Γ p , Γ k ∩ Γ m = ∅ for k = m, i = 1, 2, 3, s j are local coordinates of the surface Γ p and C 1 , C 2 are given constants, T ≡ (T 1 , T 2 , T 3 ) ∈ U T ad = {T i ∈ U T i ad , i = 1, 2, 3 and Γ T dS = 0, Γ x × T dS = 0}.
Finally, we define U ad = U λ ad × U µ ad × U T ad and A ≡ {A, T }, A = {λ, µ}.

We will consider the following criterion-functional based on the topological derivative associated to the energy shape functional (34)

Φ(A, σ) = σ T H(A)σ
where σ ≡ σ(y) is the stress tensor of a full body at the center y ∈ Ω of a spherical cavity,

H(A) = 3(1 -ν) 4E (Λ 1 + 10(1 + ν) 7 -5ν Λ 2 ), (40)
Λ 1 = 1 3 I ⊗I, Λ 2 = II -Λ 1
, ν is the corresponding Poisson's constant and E the Young's modulus. Note that ν = λ 2(λ+µ) , E = µ(3λ+2µ) λ+µ .

Continuous dependence of the criterion on the input data

Our main result of the present section is given by the following theorem

Theorem 4 Let A n ∈ U ad , A n → A in R 2 × [L ∞ (Γ)] 3 as n → ∞. Then Φ(A n , σ(A n )) → Φ(A, σ(A)).
Proof. The estimate (47) follows from the Cauchy-Schwartz inequality and the boundedness of sets U λ ad , U µ ad . To justify (48), we write

a(A; u, u) ≥ 2µ Ω ε ij (u)ε ij (u)dx
and use the Korn's inequality

Ω ε ij (u)ε ij (u)dx ≥ c u 2 1,Ω ∀u ∈ V 0 (see e.g. [7]-Lemma 7.3.3). Lemma 2 Let λ n ∈ U λ ad , µ n ∈ U µ ad , λ n → λ and µ n → µ as n → ∞. Then ν n → ν and T * ij (ν n ) → T * ij (ν) in [L 2 (Γ)] 3 . ( 49
)
Proof. Since λ n + µ n ≥ λ + µ > 0, ν n = λ n 2(λ n + µ n ) → λ 2(λ + µ) = ν. We infer that s ij k (ν n ) → s ij k (ν) in L 2 (Γ), k = 1, 2, 3. (50)
Indeed, we have

λ n /κ n → λ/κ, µ n /κ n → µ/κ and u ij0 (ν n ) -u ij0 (ν) H 1 (Γ) ≤ C|ν n -ν| → 0,
so that (50) holds. Since the field w ij is independent of A, we arrive at (49).

Lemma 3 Let λ n ∈ U λ ad , µ n ∈ U µ ad , λ n → λ and µ n → µ as n → ∞ and u * ij (A n ) ∈ V 0 . Then u * ij (A n ) → u * ij (A) in [H 1 (Ω)] 3 .
Proof. For brevity, let us denote

T * n = T * ij (ν n ), T * = T * ij (ν), u * n = u * ij (A n ), u * = u * ij (A). By definition, we have a(A n ; u * n , v) = Γ T * n v dS (51) a(A; u * , v) = Γ T * v dS (52) for all v ∈ [H 1 (Ω)] 3 .
Let us consider also solutions ûn ∈ V 0 of the following problem

a(A; ûn , v) = Γ T * n v dS ∀v ∈ [H 1 (Ω)] 3 . (53)
From ( 53) and (52) we obtain

a(A; ûn -u * , v) = Γ (T * n -T * )v dS.
Proof. We may write

|σ kl (A n ) -σ kl (A)| = | Γ T n • (c klij (A n )G ij y (A n ))dS -Γ T • (c klij (A)G ij y (A))dS| ≤ | Γ T n • ((c klij (A n ) -c klij (A))G ij y (A n ))dS| +| Γ T n • (c klij (A)(G ij y (A n ) -G ij y (A)))dS| +| Γ (T n -T ) • (c klij (A)G ij y (A))dS| ≡ I 1 + I 2 + I 3 ,
where

I 1 ≤ Γ C A n -A 0,∞ |G ij y (A n )|dS → 0 and I 2 ≤ Γ C|G ij y (A n ) -G ij y (A)|dS → 0
due to Proposition 1 and the boundedness of T n in [L ∞ (Γ)] 3 . I 3 tend to zero by assumption.

Proof of Theorem 1. We have

|Φ(A n , σ(A n )) -Φ(A, σ(A))| ≤ |σ(A n ) T H(A n )(σ(A n ) -σ(A))| +|σ(A n ) T (H(A n ) -H(A))σ(A)| +|(σ(A n ) T -σ(A) T )H(A)σ(A)| = J 1 + J 2 + J 3 .
By Proposition 2 we infer that J 1 and J 3 tend to zero. We also use the continuity of the function A → H(A), which follows from Lemma 2 and the convergence

E n = 2µ n (1 + ν n ) → 2µ(1 + ν) = E ≥ 2µ > 0.
As a consequence, J 2 tends to zero, as well.

The worst scenario and the maximum range scenario

Suppose that we wish to be "on the safe side", taking uncertain input data A and T in consideration. Then we solve either the worst scenario problem

A 0 = arg max A∈U ad Φ(A, σ(A)) (59)
or the maximum range scenario problem: find (i) A 0 according to (59) and

(ii)

A 0 = arg min A∈U ad Φ(A, σ(A)). (60)
In other words, we seek exact upper and lower bounds of the criterion functional (see the monograph [START_REF] Hlaváček | Uncertain Input Data Problems and the Worst Scenario Method[END_REF] for applications of problem (60) within the frame of the fuzzy set theory).

Theorem 5 Problems (59) and ( 60) have at least one solution.

Proof. The set U ad is compact in R 2 × (

Traction problem in 2D-elasticity

Let us consider a plane elasticity, i.e., either the case of plane strain or that of plane stress. It is well-known, that both cases have the same stress-strain relations, where only the coefficient λ varies It is either λ or λ , see [START_REF] Novotny | Topological Sensitivity Analysis for three-dimensional linear elasticity problems[END_REF].

λ = Eν (1 + ν)(1 -2ν)
for plane strain, whereas

λ = λ * = Eν 1 -ν 2 for plane stress.
Let us consider a bounded domain Ω ⊂ R 2 with a Lipschitz boundary ∂Ω ≡ Γ, occupied by a homogeneous and isotropic elastic body, loaded only by surface loads T ∈ [L ∞ (Γ)] 2 . Assume that λ ∈ U λ ad , µ ∈ U µ ad and T i ∈ U T i ad , i = 1, 2, with U λ ad , U µ ad and U T i ad defined in section 1. Moreover, assume that the forces T are in equilibrium, i.e.

Γ T dS = 0, Γ (x 1 T 2 -x 2 T 1 )dS = 0. ( 61
)
We define

U T ad = {T ≡ (T 1 , T 2 ) : T i ∈ U T i ad , i = 1, 2, T satisfy (61)}, U ad = U λ ad × U µ ad × U T ad , A = {λ, µ}, A = {A, T }
and introduce the criterion-functional based on the topological derivative associated to the energy shape functional (39)

Φ(A, σ) = σ T H(A)σ, (62)
where σ ≡ σ(y) is the stress tensor of a full body at the center y ∈ Ω of a circular cavity, and

H(A) = (K + µ) 2Kµ (Λ 1 + 2Λ 2 ), (63)
where K = λ + µ is the bulk modulus.

Continuous dependence of the criterion on the input data

The main result of the present section will be represented by an analogue of Theorem 1 as follows.

Theorem 6 Let A n ∈ U ad , A n → A in R 2 × [L ∞ (Γ)] 2 as n → ∞. Then Φ(A n , σ(A n )) → Φ(A, σ(A)).
For the proof we shall employ the following integral representation formula, analogous to (41), namely

∂u i ∂y j (y) = Γ T • G ij y dS, i, j ∈ {1, 2}. (64)
By Reciprocity theorem, we obtain

Γ (T • u * ij -T * ij • u)dS = 0. (71)
Then ( 69) and (71) yield

∂u i ∂y j (y) = Γ T • (u * ij -u ij )dS + Γ u • (s ij -T * ij )dS.
The last integral vanishes by virtue of normalization conditions, since

s ij -T * ij = -w ij .
As a consequence, we arrive at the formula (64), where

G ij y = u * ij -u ij . (72)
Now we may go on in proving Theorem 3 as in the prooof of Theorem 1. We establish an analogue of Lemma 1, where the subspace V 0 is defined by

V 0 = v ∈ [H 1 (Ω)] 2 : Γ v dS = 0, Γ (x 1 v 2 -x 2 v 1 )dS = 0 .
For the Korn's inequality in V 0 , see e.g. Section 10.2.2 in [START_REF] Nečas | Mathematical Theory of Elastic and Elasto-plastic Bodies : An Introduction[END_REF].

As far as an analogue of Lemma 2 is concerned, we use the formula

ν = λ 2(λ + µ)
for plane strain and ν = λ * λ * + 2µ for plane stress.

It is readily seen that s ij ≡ s ij (ν), i.e., it does not depend on the modulus E. Then we can prove that λ * n → λ * , ν n → ν and

s ij (ν n ) → s ij (ν) in L 2 (Γ)] 2 as ν n → ν , since λ n (K n + 2µ n )/κ 0n → λ(K + 2µ)/κ 0 (73) and λ n K n /κ 0n → λK/κ 0 for K n = λ n + µ n , λ n ∈ U λ ad , µ n ∈ U µ ad , A n → A.
The field w ij is independent of A, so that we arrive at

T * ij (ν n ) → T * ij (ν) in L 2 (Γ)] 2 .
An analogue of Lemma 3 can be proved in the same way as Lemma 3. We infer that

u * ij (A n ) → u * ij (A) in [H 1 (Ω)] 2 . ( 74
)
Using again (73), we observe that

u ij (A n ) → u ij (A) in L 2 (Γ)] 2 . (75)
Combining (72) with (74), the Trace theorem and (75), we obtain

G ij y (A n ) → G ij y (A) in L 2 (Γ)] 2 . ( 76
)
Theorem 3 follows in a way parallel to the proof of Theorem 1, from (76), the uniform convergence of surface loads on Γ and the continuity of the function A → H(A).

The worst scenario and the maximum range scenario

Both the worst scenario problem (59) and the maximum range scenario problem (60) have at least one solution. This assertion is a consequence of Theorem 3 and the compactness of the set

U ad in R 2 × 2 i=1 P p=1
C(Γ p ).

• for i = 1

σ rr 1 (ξ ρ ) = σ I (ξ) 14 -10ν ρ 3 r 3 - ρ 5 r 5 + 14 -10ν -10(5 -ν) ρ 3 r 3 + 36 ρ 5 r 5 sin 2 θ sin 2 ϕ , (84)
σ rθ 1 (ξ ρ ) = σ I (ξ) 14 -10ν -5ν + 5(1 + ν) ρ 3 r 3 -12 ρ 5 r 5 sin 2θ sin 2 ϕ , (85)
σ rϕ 1 (ξ ρ ) = σ I (ξ) 14 -10ν -5ν + 5(1 + ν) ρ 3 r 3 -12 ρ 5 r 5 sin θ sin 2ϕ , (86)
• for i = 2 Remark 4 It is important to mention that the stress distribution for i = 1, 2 was obtained from a rotation of the stress distribution for i = 3. In addition, the derivation of this last result (for i = 3) can be found in [START_REF] Kachanov | Handbook of Elasticity Solutions[END_REF], for instance.

  p )), so that the assertion follows from Theorem 1.

r 5 -r 5 r 5

 555 14 -10ν + 25(1 -2ν) cos 2θ sin 2 ϕ , (87)σ θϕ 1 (ξ ρ ) = σ I (ξ) 14 -10ν -5ν + 5(1 -2ν) ρ 3 r 3 + 3 ρ 5 r 5 cos θ sin 2ϕ ,(88)σ ϕϕ 1 (ξ ρ ) = σ 1 (ξ) 56 -40ν -20ν + (11 -10ν) ρ 3 r 3 + 9 ρ 5 r 5 + 28 -20ν + 5(1 -2ν) cos 2θ sin 2 ϕ ,

r 5

 5 sin 2 θ , (101) where σ I (ξ), σ III (ξ) and σ III (ξ) are the principal stress values of tensor σ(ξ), for ξ = u, ξ = w or ξ = v associated to the original domain without hole Ω.

1,Ω ∀u ∈ V 0 ,(48)whereV 0 = {v ∈ [H 1 (Ω)] 3 : Γ v dS = 0, Γ v × x dS = 0}.

Γ (s ij • u -T • u ij )dS, (69)which follows by differentiating the so-colled Somigliana's identityu i (y) = Γ (T • u i y -u • s i y )dS,(70)where ∂s i y /∂y j = -s ij .

Acknowledgments 4 Conclusions

We have seen that the worst case and maximal range scenario problems are solvable with criterions of energy-based topological derivative. The same methodology, considering topological derivatives of different shape functional may be applied to derive similar analysis for criteria dependent for example on displacement (kinematic constraints) and yield constraints.

Proof is based on the formulas ([7]-Theorem 10.1.1)

and

where

and r = x -y. Since (κ(A n )) -1 -→ (κ(A)) -1 and the components u ij0 k are bounded on Γ,

The vector field u * ij (A) is the displacement solving the first boundary value problem with zero body forces and the equilibriated surface loading

where

and

represents a rigid body displacement such that

(e i denote unit vectors in the directions of Cartesian coordinates). The field w ij is uniquely determined by the conditions shown. Inserting (43) in (45), we observe that

since u ij0 , λ κ and µ κ are independent of the modulus E. Lemma 1 Let us define

Inserting v := ûn -u * and using Lemma 1, we infer that

so that ûn -u * 1,Ω → 0 follows from Lemma 2. We can show that

Indeed, (51) and Lemma 1 yield that

so that (55) follows from Lemma 2. We can use ( 51), ( 53) and Lemma 1 to obtain

Then ( 55) and ( 56) yield

The convergence u * n → u * in [H 1 (Ω)] 3 follows from the triangle inequality, (54) and (57).

as n → ∞.

Proof. Since by (42) we have

the assertion follows from Lemma 3, the Trace theorem and (44).

Proposition 2 Let the stress components at the point y be

We can construct the vector function G ij y in a way parallel to that of the proof of Theorem 10.1.1 in [START_REF] Nečas | Mathematical Theory of Elastic and Elasto-plastic Bodies : An Introduction[END_REF]. First, we consider the well-known Kelvin's solution

where κ 0 = 4πµ(K + µ), r = x -y and define u ij = -∂u i y /∂y j . The corresponding surface forces on Γ are then

We can find that

Let us construct the rigid body translation

Note that the field w ij is uniquely determined by conditions (67). If we define

the forces T * ij are in equilibrium, i.e., they satisfy conditions (61). There exists a unique displacement field u * ij , which solves the first boundary value problem of elasticity with zero body forces and surface loads T * ij and satisfies the normalization conditions

Next, we assume that the field u fullfils conditions (68) as well and consider the so-called Love's formula ∂u i ∂y j (y) =

A Stress distribution around cavities

We present in this appendix the analytical solution for the stress distribution around a circular (N = 2) and spherical (N = 3) cavities respectively for two and three-dimensional linear elastic bodies.

A.1 Circular cavity

Considering a polar coordinate system (r, θ), we have the following expansion for the stress distribution σ(ξ ρ ) around a free boundary circular cavity (σ rr (ξ ρ ) = σ rθ (ξ ρ ) = 0 on ∂B ρ (x 0 )), with

where the angle θ u = θ and θ w = θ + δ, with δ denoting the angle between principal stress directions for displacement fields u and w in (35). In addition, the following expansion for σ(v ρ ) satisfying the boundary condition on ∂B ρ (x 0 ) given by σ rθ (v ρ ) = 0 and σ rr (v ρ ) = -2γσ θθ (u ρ ), holds

where the angle θ v = θ + δ, with δ denoting the angle between principal stress directions for displacement fields u and v in (36). Finally,

where σ I (ξ) and σ II (ξ) are the principal stress values of tensor σ(ξ), for ξ = u, ξ = w or ξ = v associated to the original domain without hole Ω.

A.2 Spherical cavity

Let us introduce a spherical coordinate system (r, θ, ϕ). Then, the stress distribution around the spherical cavity B ρ (x 0 ) is given by σ rr (ξ ρ ) = σ rr 1 (ξ ρ ) + σ rr 2 (ξ ρ ) + σ rr 3 (ξ ρ ) + O(ρ) , σ rθ (ξ ρ ) = σ rθ 1 (ξ ρ ) + σ rθ 2 (ξ ρ ) + σ rθ 3 (ξ ρ ) + O(ρ) , σ rϕ (ξ ρ ) = σ rϕ 1 (ξ ρ ) + σ rϕ 2 (ξ ρ ) + σ rϕ 3 (ξ ρ ) + O(ρ) , σ θθ (ξ ρ ) = σ θθ 1 (ξ ρ ) + σ θθ 2 (ξ ρ ) + σ θθ 3 (ξ ρ ) + O(ρ) , σ θϕ (ξ ρ ) = σ θϕ 1 (ξ ρ ) + σ θϕ 2 (ξ ρ ) + σ θϕ 3 (ξ ρ ) + O(ρ) , σ ϕϕ (ξ ρ ) = σ ϕϕ 1 (ξ ρ ) + σ ϕϕ 2 (ξ ρ ) + σ ϕϕ 3 (ξ ρ ) + O(ρ) ,

where ξ ρ = u ρ , ξ ρ = w ρ or ξ ρ = v ρ ; σ rr i (ξ ρ ), σ rθ i (ξ ρ ), σ rϕ i (ξ ρ ), σ θθ i (ξ ρ ), σ θϕ i (ξ ρ ) and σ ϕϕ i (ξ ρ ), for i = 1, 2, 3, are written, as: