Shape and Topology Optimization of Composite Materials with the - - PowerPoint PPT Presentation

Inverse Problems, Control and Shape Optimization PICOF ’12

Shape and Topology Optimization

• f Composite Materials with the level-set method

Gabriel Delgado Gr´ egoire Allaire

EADS-Innovation Works Center of Applied Mathematics - ´ Ecole Polytechnique

Palaiseau, 3 April 2012

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 1 / 33

Outline

1

Motivation

2

Problem description

3

Main methods in Topology Optimization

4

The algorithm and its implementation

5

Results

6

Conclusion

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 2 / 33

Outline

1

Motivation

2

Problem description

3

Main methods in Topology Optimization

4

The algorithm and its implementation

5

Results

6

Conclusion

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 3 / 33

Motivation: Composite Materials in Aeronautics

Figure: Composite structures of an A380 and evolution of the use of composites in Airbus.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 4 / 33

Motivation: Laminated Composite Materials

Characteristics

Lamination of a sequence of unidirectionally reinforced plies. Each ply is typically a thin sheet

• f collimated fibers impregnated

with a polymer matrix material. Strong resistance against severe environmental conditions. Less expensive and lighter than metallic alloys. Greatest benefit and drawback: Espace of design possibilities.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 5 / 33

Motivation: Shape and Topology Optimization examples

Figure: Topology optimized wing structure of the A380 and a shape optimized airfoil.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 6 / 33

Outline

1

Motivation

2

Problem description

3

Main methods in Topology Optimization

4

The algorithm and its implementation

5

Results

6

Conclusion

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 7 / 33

2D Multi-layered composite description

Let be D ⊂ R2 fixed and Ω =

• Ωi

a laminated composite structure made from the superposition of N plies under plane stress state, of geometry Ωi ⊂ D, each one composed by two anisotropic elastic phases Ai and A0 (

• A0

<< 1). We consider the boundary of D made of two disjoint parts, ∂D = ΓD ∪ ΓN, |ΓD| = 0.

A0 A90 A−45 A45 Ω1 Ω2 Ω3 Ω4

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 8 / 33

Physical model description

On Ωi we denote by Σi the interface between both phases. The characteristic function of Ωi is denoted by χi(x), so the composite elastic tensor of Ω at x is described by Aχ =

N

• i=1

Aiχi +

N

• i=1

A0(1 − χi). The strain and stress tensors are related to the displacement field u as e(u) = 1

2(∇u + ∇T u) and σ(u) = Aχe(u). The external charges are

denoted by f (volume) and g (surface), and the displacement field uχ is the solution of the linearized elasticity system in D ⎧ ⎨ ⎩ −div(Aχe(u)) = f in D u = 0

• n ΓD

Ae(u)n = g

• n ΓN

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 9 / 33

Optimization problem

Find the best composite structure Ω according to a criteria J(Ω, uχ), where uχ is the state variable (displacement) and Ω is the control variable, among a set of admissible shapes of layers Uad (manufacture and geometric constraints), under a given state of external charges (f, g) by changing the geometry of each biphasic layer Ωi min

{Ωi}i=1..N∈ Uad, Ωi⊂D J(Ω, uχ)

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 10 / 33

Outline

1

Motivation

2

Problem description

3

Main methods in Topology Optimization

4

The algorithm and its implementation

5

Results

6

Conclusion

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 11 / 33

Main methods in Topology Optimization

Homogenization The idea is to find the optimal distribution of material inside a structure, by authorizing intermediates densities (0 ≤ θi ≤ 1) of each material. To do this, we define the set of homogenized tensors Gθ, corresponding to the set of composite materials A∗ made from the mixture of the phases Ai in θi proportion, where A∗ is the law tensor representing the homogenized micro-structure. Unfortunately the main problem of this method is to find an explicit description of the set Gθ.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 12 / 33

Main methods in Topology Optimization

Level set method Has as unique advantages to track topology changes and a clear and smooth boundary that can be easily managed. Problems that arise from a density approach like spurious eigenfrequencies and micro structure stress concentration are avoided. Coupled to the shape and topological derivative analysis, it makes the level set method a promising research direction in future applications on structures design. Let’s take a look of how it works.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 13 / 33

General explanation: Level set method

Ω0

φ < 0 φ > 0 φ = 0

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 14 / 33

General explanation: Shape sensitivity

θ Ωθ = (Id + θ)(Ω0) J(Ωθ) < J(Ω0)

Boundary movement

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 15 / 33

General explanation: Topological sensitivity

ω Aω

AΩ J(Ωω

θ ) < J(Ωθ)

Ωω

θ

Nucleation of an inhomogeneity

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 16 / 33

Outline

1

Motivation

2

Problem description

3

Main methods in Topology Optimization

4

The algorithm and its implementation

5

Results

6

Conclusion

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 17 / 33

Cost function and the adjoint problem

Let’s take a general cost function J(Ω, uχ(Ω)) as J(Ω, uχ(Ω)) =

• D

j(x, uχ(Ω))dx +

• ∂D

h(x, uχ(Ω))ds. This criteria can represent for example the compliance (j(x, u) = f · u, h(x, u) = g · u), the volume (j = 1), a least square target displacement (j(x, u) = |u − u0|2), a stress dependent function (j(x, u) = |σ(u)|2), etc. Furthermore, in order to avoid the calculation of the explicit variation of uχ w.r.t. the domain, we introduce the adjoint state problem ⎧ ⎨ ⎩ −div(Aχe(pχ)) = −j′(x, uχ) in D pχ = 0

• n ΓD

Aχe(pχ)n = −h′(x, uχ)

• n ΓN

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 18 / 33

The level-set method

We define the level set function φi(x) on D ⊂ R2 as ⎧ ⎨ ⎩ φi(x) = 0 ⇔ x ∈ ∂Ωi ∩ D, φi(x) < 0 ⇔ x ∈ Ωi, φi(x) > 0 ⇔ x ∈ (D \ ¯ Ωi). From this definition we can easily deduce Outward normal vector ni = ∇φi/|∇φi| Curvature κi = div ni This formulae have a meaning over all D and not only on Σi. If the domain Ωi(t) evolves in a pseudo-time t ∈ R+ according to a velocity field θ(x, t), then the level set transport equation is ∂φi ∂t + V|∇φi| = 0. where V(x, t) = θ · n is the normal component of the advection field.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 19 / 33

Level-set method: Reinitialization

In order to regularize the level-set function (which may become too flat or too steep), we reinitialize it periodically by solving

• − ∂φ

∂t + sign(φ0)(1 − |∇φ|) = 0

in D × R+ φ(t = 0, x) = φ0(x) in D which admits as a stationary solution the signed distance to the initial interface {φ0 = 0}. Reinitialize the level function is really important to

• btain a good approximation of the normal ni and the curvature κi of Σi.

Question: How we find V that minimizes J(Ω, uχ)? Hint: Shape derivative.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 20 / 33

Shape derivative: Variation of the cost function under boundary changes

We are interested on the variations of the functional J(Ω, uχ(Ω)) when the position of the interface Σi on each ply changes following a regular vector field θ. Let be ω0 an open smooth set such that ω0 ⊂ D. We denote by χω

0 the characteristic function of ω0, and we consider the

variations on the form χω

θ = χω 0 ◦ (Id + θ), i.e. χω θ = χω 0 ◦ (x + θ(x)).

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 21 / 33

Shape derivative

Definition

Let be the functional J(ω) : Uad → R. The shape derivative of J(ω) at ω0 is defined as the Fr´ echet derivative in W 1,∞(R2, R2) at 0 of the l’application θ → J(ω0 ◦ (Id + θ)), i.e. J(ω0 ◦ (Id + θ)) = J(ω0) + J′(ω0)(θ) + o(θ) with lim

θ→0

|o(θ)| θW 1,∞ = 0 where J′(ω) is a continuous linear form on W 1,∞(D; R2).

Theorem: Level-set adapted fixed mesh shape derivative

Let be Dh and Ωh =

• Ωi

h

• polygonal approximations of D and Ω. Given a

triangulation TΣ = {Kl}l of Dh, with Kl ∩ Σ = ∅, then J′(Ω)(θ) =

• Kl

Aχe(uχ) : e(pχ)(∇ · θ)dx −

• ∂Kl

Aχe(uχ) : e(pχ)θ · nds +

• Kl

Aχ{∇(e(uχ)) · θ} : e(pχ) + Aχ{∇(e(pχ)) · θ} : e(uχ)dx

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 22 / 33

How do we chose V using the shape sensitivity?

We take θ∗ smooth (θ∗W 1.∞ < 1) such that J′(Ω)(θ∗) < 0. Since the only values that matter are θ∗

|Σ (Hadamard), we solve the

regularization extension problem θ∗, vW (D) = −J′(Ω)(v), ∀v ∈ W(D), for some Hilbert space W(D) (H1 for instance). Then we chose as advection velocity for the level-set V = θ∗n. If we take a small enough advection step t << 1, we can assure that the objective function will decrease using the shape derivative definition: J(Ω0 ◦ (Id + tθ)) = J(Ω0) + tJ′(Ω0)(θ) + o(t). Drawback: Tendency to local minima, no change of topology (nucleation). Hint: Topological gradient.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 23 / 33

The topological derivative

Let ω a smooth open subset of R2. Let ρ > 0 be a small positive parameter which is intended to go to zero. For a point z ∈ D with d(z, ∂D) < ρ and d(z, Σi) < ρ, ∀i < N, we define the rescaled inclusion ωρ =

• x ∈ R2 : x − z

ρ ∈ ω

• ,

which, for a small enough ρ does not intersect any Σi. We are interested

• n the variation of the objective function when we embed an inclusion ωρ
• f Hooke law A∗ on the i-layer obtaining a composite Ωρ.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 24 / 33

Definition

If the objectif function J admits the following so-called topological asymptotic expansion for small ρ > 0: J(Ωρ) − J(Ω) − ρ2DJ(z) = o(ρ2), then the number DJ(z) is called the topological derivative of J at z for the inclusion shape ω.

Theorem

The topological derivative DJ(z) of the general cost function J, evaluated at z for an inclusion shape ω, has the following expression DJ(z) = −M(A, A∗, ω)e(uχ)(z) : e(pχ)(z) where M(A, A∗, ω) is the so-called elastic momentum tensor.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 25 / 33

The optimization algorithm

1

Initialization of each level set function φi

0 corresponding to an initial

guess χi

0 (most of time trivial guess χi 0 ≡ 1).

Iteration k

2

Solve the direct (uk

χ) and adjoint (pk χ) state problems posed on Ωk.

3

Shape/Topology sensitivities and level set advection For each layer i

◮ If k mod(Nopt) = 0 ◮ Calculate the topological derivative DJ and remove a small percentage

• f the area where it reaches the minimum.

◮ If k mod(Nopt) = 0 ◮ Computation of a regular velocity Vi

k such that the shape derivative

satisfies J′(χi

k)(Vi kn) ≤ 0.

◮ Deformation of the shape by solving the level-set equation with a time

step of ∆sk, chosen as J(χk+1) ≤ J(χk).

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 26 / 33

Outline

1

Motivation

2

Problem description

3

Main methods in Topology Optimization

4

The algorithm and its implementation

5

Results

6

Conclusion

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 27 / 33

Airbus simplified fuselage test case

Description

Find the stiffest (lightest) four layered (0o,90o,45o,-45o) composite structure, under a weight (compliance) constraint for the following cantilever-type problem

L=2 H=1 F1

This corresponds to an equivalent model of a flattened section of a fuselage.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 28 / 33

Algorithm parameters

Description

Unique structured isotropic rectangular mesh for the level-set and FE analysis. uχ, pχ = uχ and θ were taken as [P1(Ω), P1(Ω)] finite elements. Adapted augmented Lagrangian method used for the constrained

• ptimization.

Number of iterations 500, time 20 min (3200 elements, Intel(R) Core(TM) Duo CPU 2.8GHz, 2GB RAM).

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 29 / 33

Results

90o

• 45o

90o

• 45o

iteation 497; J= 2.79598; Constraint violation: -0.00645053

0o 45o 0o 45o

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 30 / 33

Outline

1

Motivation

2

Problem description

3

Main methods in Topology Optimization

4

The algorithm and its implementation

5

Results

6

Conclusion

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 31 / 33

Conclusions

1

A shape and topology optimization algorithm has been developed for composite materials, where the design variable is the shape and topology of each ply.

2

Difficulties related to the elastic low contrast and anisotropy of the materials were successfully overcome, defining the level-set adapted fixed mesh shape derivative and the elastic moment tensor for anisotropic materials.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 32 / 33

Perspectives

1

Introduce new pointwise cost functions dependent of the gradient of the state function (e.g. Von-mises,Reserve factor) and buckling load.

2

Enhance the numerical performance of the optimization algorithm using a feasible direction method (connexion to classical SQP

• ptimizers) instead of an Augmented Lagrangian method.

3

Explore two other design variables: the orientation of each ply and the stacking sequence.

4

Add more plies to the composite (8 for example) and parallelize the code.

Gabriel Delgado (EADS-X) PICOF ’12 Palaiseau, 3 April 2012 33 / 33