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 S IMO & H UGHES , 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 ǫ in = γ ( σ/σ ref ) ǫ in ) Maxwell creep model ˙ σ = E (˙ ˙ ǫ − ˙ σ + ( Eγ/σ ref ) σ = E ˙ ˙ ǫ y + λy = f, y (0) = y 0 ˙ λ = Eγ/σ ref ≥ 0 , f = ηλσ ref , η = ˙ ǫ/γ
AMPLIFICATION FACTOR ( f = 0 and λ constant) 1 A bE = A dG(0) = 1 + λ ∆ t 1 − 1 3 λ ∆ t A dG(1) = A IRKR2A − 3 = 1 + 2 3 λ ∆ t + 1 6 ( λ ∆ t ) 2 1 − 2 5 λ ∆ t + 1 20 ( λ ∆ t ) 2 A dG(2) = A IRKR2A − 5 = 20 ( λ ∆ t ) 2 + 1 1 + 3 5 λ ∆ t + 3 60 ( λ ∆ t ) 3 1 A IRKL3C − 2 = 1 + λ ∆ t + 1 2 ( λ ∆ t ) 2
1 0.8 0.6 A 0.4 dG(0) = BE 0.2 IRKL3C ☛ dG(2) ❂ ❘ 0 dG(1) ✰ 0 2 4 6 8 10 12 14 16 18 20 ∆ tλ 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 � t n +1 � t n +1 y + λy ) w dt + [ y n ] w + ( ˙ n = fw dt t n t n [ y n ] = y + n − y − where n
PETROV-GALERKIN y = N 1 y + w = N w 1 ω 1 + N w n + N 2 y − and λ = N 1 λ n + N 2 λ n +1 n +1 , 2 ω 2 � ω 1 � T �� A 11 � 1 � f 1 y + � � � � �� A 12 + ( y + n n − y − n ) = 0 − y − 0 ω 2 A 21 A 22 f 2 n +1 tn +1 tn +1 � � N w ˙ λN w A ij = m ij + k ij , m ij = k ij = N j dt, i N j dt i tn tn − A 21 (1 + A 11 ) f 2 − A 21 f 1 y − n +1 = y − n + (1 + A 11 ) A 22 − A 12 A 21 (1 + A 11 ) A 22 − A 12 A 21 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) = y 0 ˙ ( 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) β = 0 . 5 λ 0 , φ = − 1 . 0 , η = 2 . 0 1.5 1.4 1.3 A exact IRKL3C 1.2 bE 1.1 1 0 2 4 6 8 10 12 14 16 18 20 λ 0 ∆ t
Vanishing diffusivity (strain hardening) β = 0 . 5 λ 0 , φ = 1 . 0 , η = 2 . 0 40 35 30 25 A 20 15 exact IRKL3C 10 bE 5 0 2 4 6 8 10 12 14 16 18 20 λ 0 ∆ t
MATERIAL MODEL ǫ in = 3 2 γ n , where n = s / ¯ σ ˙ The scalar ¯ σ is the equivalent stress � ¯ � − Q � σ � γ = f ∗ exp sinh m Rθ σ y
Bar in uniaxial tension (strain rate 10 − 5 1/s) θ ( t ) = θ 0 ± ∆ θ ( t/t max ) , where t max = ǫ max / ˙ ǫ θ 0 = 293 K , ∆ θ = 40 K Binary near eutectic Sn40Pb solder: E = 33 GPa Q = 12 kcal/mol ν = 0 . 3 2 · 10 − 3 kcal/mol · K R = σ y 20 MPa = 10 5 s − 1 f = m = 3 . 5
ǫ = 10 − 5 1/s) Thermally softening case ( ˙ 0.5 × 0.4 + � + × � + × 0.3 σ � + × � + × σ y � 0.2 exact dG(0) � 0.1 + bE × IRKL3C + × � 0 0 0.002 0.004 0.006 0.008 0.01 ǫ
ǫ = 10 − 5 1/s) Thermally hardening case ( ˙ � + × 1 � + × 0.8 � + × � + σ × 0.6 � × + σ y 0.4 exact dG(0) � + 0.2 bE × IRKL3C + × � 0 0 0.002 0.004 0.006 0.008 0.01 ǫ
ǫ = 10 − 5 1/s) Relative error ( ˙ 10 − 1 � � 10 − 2 � � 10 − 3 + × + × 10 − 4 + × dG(0) � + + bE 10 − 5 × IRKL3C × 10 − 6 10 − 2 10 − 1 ∆ t/t max
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
Recommend
More recommend