A time discontinuous Petrov-Galerkin method for the integration of - - PowerPoint PPT Presentation

A time discontinuous Petrov-Galerkin method for the integration of inelastic constitutive equations

Reijo Kouhia Laboratory of Structural Mechanics Helsinki University of Technology

OUTLINE

• Motivation
• Integration algorithms

– Amplification factors of a model problem – Discontinuous Galerkin approach

• Example
• Concluding remarks

MOTIVATION

SIMO & HUGHES, Computational Inelasticity, Remark 3.3.2.2 on page 125: “The overall superiority of the radial return method relative to other return schemes is conclusively established in Krieg and Krieg [1977]; Schreyer, Kulak and Kramer [1979] and Yoder and Whirley [1984].”

WHY ??

SCALAR MODEL PROBLEM Maxwell creep model ˙ ǫin = γ(σ/σref) ˙ σ = E(˙ ǫ − ˙ ǫin) ˙ σ + (Eγ/σref)σ = E ˙ ǫ ˙ y + λy = f, y(0) = y0 λ = Eγ/σref ≥ 0, f = ηλσref, η = ˙ ǫ/γ

AMPLIFICATION FACTOR (f = 0 and λ constant) AbE = AdG(0) = 1 1 + λ∆t AdG(1) = AIRKR2A−3 = 1 − 1

3λ∆t

1 + 2

3λ∆t + 1 6(λ∆t)2

AdG(2) = AIRKR2A−5 = 1 − 2

5λ∆t + 1 20(λ∆t)2

1 + 3

5λ∆t + 3 20(λ∆t)2 + 1 60(λ∆t)3

AIRKL3C−2 = 1 1 + λ∆t + 1

2(λ∆t)2

0.2 0.4 0.6 0.8 1 2 4 6 8 10 12 14 16 18 20

A ∆tλ

dG(0) = BE IRKL3C dG(1) dG(2)

☛ ❂ ✰ ❘

IRKL3C method is the most accurate when ∆t is large enough !!!! BE also good when λ∆t > 15

SOME REQUIREMENTS

Ideal integrator for inelastic constitutive models should be:

• 1. L-stable
• 2. For ˙

y + λy = 0 (λ constant) the amplification factor should be (a) strictly positive (b) monotonous (convex) Pad´ e (0, q)-approximations of exp(−λt) are positive and monotonous. IRKL3C-2 = Pad´ e-(0,2) for ˙ y + λy = 0

QUESTION Can we design a discontinuous Galerkin method with these properties ?? ANSWER yes, if using discontinuous Petrov-Galerkin approach or underintegration tn+1

tn

( ˙ y + λy)w dt + [yn] w+

n =

tn+1

tn

fw dt where [yn] = y+

n − y− n

PETROV-GALERKIN

y = N1y+

n + N2y− n+1,

w = N w

1 ω1 + N w 2 ω2

and λ = N1λn + N2λn+1 ω1 ω2 T A11 A12 A21 A22 y+

n

y−

n+1

• + (y+

n − y− n )

1

• −

f1 f2

• = 0

Aij = mij + kij, mij =

• tn+1

tn

N w

i

˙ Nj dt, kij =

• tn+1

tn

λN w

i Nj dt

y−

n+1 =

−A21 (1 + A11)A22 − A12A21 y−

n +

(1 + A11)f2 − A21f1 (1 + A11)A22 − A12A21

UNDERINTEGRATION Using 2-point Lobatto integration for the dG(1)-scheme we get the two stage IRKL3C-method

MODEL PROBLEM

˙ y + λy = f, y(0) = y0 (f = ηλσref, η = ˙ ǫ/γ, λ = Eγ/σref) λ(t) = λ0 [1 − φ + φ exp(−βt)] Two special cases increasing diffusivity φ = −1 vanishing diffusivity φ = 1, then λ → 0 when t → ∞

Increasing diffusivity (strain softening)

1 1.1 1.2 1.3 1.4 1.5 2 4 6 8 10 12 14 16 18 20

A λ0∆t β = 0.5λ0, φ = −1.0, η = 2.0

exact IRKL3C bE

Vanishing diffusivity (strain hardening)

5 10 15 20 25 30 35 40 2 4 6 8 10 12 14 16 18 20

A λ0∆t β = 0.5λ0, φ = 1.0, η = 2.0

exact IRKL3C bE

MATERIAL MODEL

˙ ǫin = 3

2γn,

where n = s/¯ σ The scalar ¯ σ is the equivalent stress γ = f ∗ exp −Q Rθ

• sinhm

¯ σ σy

Bar in uniaxial tension (strain rate 10−5 1/s)

θ(t) = θ0 ± ∆θ(t/tmax), where tmax = ǫmax/˙ ǫ θ0 = 293K, ∆θ = 40K Binary near eutectic Sn40Pb solder: E = 33 GPa Q = 12 kcal/mol ν = 0.3 R = 2 · 10−3 kcal/mol·K σy = 20 MPa f = 105 s−1 m = 3.5

Thermally softening case (˙ ǫ = 10−5 1/s)

0.1 0.2 0.3 0.4 0.5 0.002 0.004 0.006 0.008 0.01

σ σy ǫ

exact dG(0)

• bE

+ + + + + + +

IRKL3C

× × × × × × ×

Thermally hardening case (˙ ǫ = 10−5 1/s)

0.2 0.4 0.6 0.8 1 0.002 0.004 0.006 0.008 0.01

σ σy ǫ

exact dG(0)

• bE

+ + + + + + +

IRKL3C

× × × × × × ×

Relative error (˙ ǫ = 10−5 1/s)

10−6 10−5 10−4 10−3 10−2 10−1 10−2 10−1 ∆t/tmax

dG(0)

• bE

+ + + + +

IRKL3C

× × × × ×

CONCLUDING REMARKS

• Two-stage IRKL3C method seems to be an accurate integrator

also for large time steps

• Discontinuous Petrov-Galerkin approach can produce a method

similar to IRKL3C

• Underintegrated dG(1) method produces a method similar to

IRKL3C

• Asymptotically third order accuracy can be obtained by

switching full integration in the dG(1) scheme if the time step is small