Compressible Navier-Stokes (Euler) Solver based on Deal.II Library - - PowerPoint PPT Presentation

Compressible Navier-Stokes (Euler) Solver based on Deal.II Library

Lei Qiao

Northwestern Polytechnical University Xi’an, China T exas A&M University College Station, T exas

1

Fifth deal.II Users and Developers Workshop T exas A&M University, College Station, TX, USA August 5, 2015

Outline

Motivation and Goal Current work Conclusion and TODOs

2

Motivation

This is what I typically need to compute

3

The pain point

Mesh generation: Time consuming Boring Experience depending Solution: Mesh adaptation — let solver tell what is good mesh

4

Deal.II the library

• Mesh adaptation with tree data structure
• FEM discretization
• Parallelization and excellent scalability

5

Goal

A N-S solver based on deal.II that:

• Starts from initial value on coarse mesh
• Converges to solution on a reasonable fine

mesh adaptively

• With High-order capability

6

Current work

Governing Equation Solving technique Solver implementation Solver verification

7

Governing Equation

Compressible Navier-Stokes Equation ∂Q(w) ∂t + r · F(w) = S(w) Use primitive variables as working var. w = (uj, ρ, p)T

8

Governing Equation

Fc(w) =   ρu ⊗ u + Ip ρu (E + p)u   =   ρuiuj + δijp ρui (E + p)ui  

9

The Flux F = Fc − Fv

Governing Equation

The Flux F = Fc − Fv With T = p γRρ τij = µ( ∂ui ∂xj + ∂uj ∂xi ) + λδij ∂uk ∂xk

10

Fv(w) = 1 Re∗

∞

  τij τijui + κ ∂T

∂xj

 

Governing Equation

Weak form: w = (uj, ρ, p)T Z

Ω

vl ∂Ql(w) ∂t + Z

Ω

vl ∂Fl,i(w) ∂xi = Z

Ω

vlSl(w) ∀v = vl ∈ T , Find w ∈ S that satisfy On domain Ω

11

Governing Equation

Weak form: Integrate by part to get boundary flux Z

Ω

vl ∂Ql(w) ∂t + Z

∂Ω

vlFl,i(w)ni − Z

Ω

∂vl ∂xi Fl,i(w) = Z

Ω

vlSl(w)

12

Governing Equation

T reat boundary conditions and discontinuities on hanging node with this boundary flux Compute numerical flux with Roe scheme[1]

13

Governing Equation

Boundary conditions:

Far field: Riemann invariant

14

Slip wall: Non-penetration

u · n = 0 ρout = ρin pout = pin

Solving technique

BDF-1 (Implicit Euler) time integration Newton linearization Direct linear solver from T rilinos*

* Iterative solver doesn’t stable enough, only used for scalability test

For steady run, time step size is determined according to norm of Newton update[2]

15

Solver implementation

• Starts from step-33 of deal.II tutorial
• Change working vars. from conservative ones to primitive ones
• viscous flux
• Riemann boundary condition
• Roe flux[2] to replace the too viscous Lax flux
• Parallelize: learned from step-40
• Version control with Git
• CMake project
• Regression test suit

16

Solver verification

• Scalability
• Manufactured solution[3]
• Subsonic
• Supersonic
• Adaptive simulation over circle
• Flow over foil NACA2412

17

Scalability

22,288 cells 91,312 DoFs

18

GMRES solver from T rilinos with DD+ILUT precond.

Solver verification

• Manufactured solution[3]

∂Q(w) ∂t + r · F(w) = S(w) w = J(x, y) =        u(x, y)= u0+ux sin(auxπx)+uy cos(auyπy)+uxy cos(auxyπxy) v(x, y)= v0 +vx cos(avxπx)+vy sin(avyπy) +vxy cos(avxyπxy) ρ(x, y)= ρ0+ρx sin(aρxπx)+ρy cos(aρyπy)+ρxy cos(aρxyπxy) p(x, y)= p0+px cos(apxπx)+py sin(apyπy) +pxy sin(apxyπxy)

S(x, y) = ∂Q(J) ∂t + r · F(J)

19

Solver

Jh(x, y)

Solver verification

• Manufactured solution[3]

w = J(x, y) =        u(x, y)= u0+ux sin(auxπx)+uy cos(auyπy)+uxy cos(auxyπxy) v(x, y)= v0 +vx cos(avxπx)+vy sin(avyπy) +vxy cos(avxyπxy) ρ(x, y)= ρ0+ρx sin(aρxπx)+ρy cos(aρyπy)+ρxy cos(aρxyπxy) p(x, y)= p0+px cos(apxπx)+py sin(apyπy) +pxy sin(apxyπxy)

Choosing different Constants to get subsonic

• r

supersonic case

20

Solver verification

• Manufactured solution
• Subsonic
• convergence 

• rder

21

Solver verification

• Manufactured solution
• Supersonic
• convergence 

• rder

22

Solver verification

• Adaptive simulation over circle
• C1 mapping on wall boundary
• Subsonic: Mach = 0.1
• Almost inviscid: Cv=1e-6 just for stabilization
• Converged Cd=0.00029 which ideal value is zero

23

Solver verification

• Adaptive simulation over cylinder

Converged density contour

24

Solver verification

• Adaptive simulation over cylinder

Converge history of density contour and mesh adaptation

25

Solver verification

• Flow over NACA2412 foil
• Subsonic: Mach=0.3

Converged pressure contour

26

Solver verification

• Flow over NACA2412 foil
• Subsonic: Mach=0.3

Oscillation may be caused by geometry non-smoothness

27

Solver verification

• Almost inviscid: Cv=3e-6 just for stabilization
• Converged Cd=0.00025 which ideal value is zero
• Converged Cl=0.2576, Cm=0.05235
• As a reference, xFoil[4] gives
• Cd=-0.00076
• Cl=0.2704
• Cm=0.0585

28

• Flow over NACA2412 foil
• Subsonic: Mach=0.3

Conclusion

• A prototype NS solver based on deal.II is constructed
• The solver could run in parallel but the solver doesn’t

scale perfectly

• High order convergence is confirmed by MMS
• Solving process could start from very coarsen mesh

and go on with mesh adaptation

• The solver can give out reasonable aerodynamics data

29

T

• Dos
• Stable and efficient preconditioner for iterative solver

(Now my hope is on MDF ordering of ILU)

• Describe wall boundary with NURBS geometry
• Non-slip boundary condition
• Anisotropic adaptation for boundary layer

30

Reference

1.

• J. Blazek. Computational Fluid Dynamics: Principles and

Applications (Second Edition). Elsevier Science, Oxford, second edition edition, 2005. ISBN 978-0-08- 044506-9. 2.

• J. Gatsis. Preconditioning T

echniques for a Newton–Krylov Algorithm for the Compressible Navier–Stokes Equations. PhD thesis, University of T

• ronto, 2013.

3.

• C. J. Roy, C. C. Nelson, T

. M. Smith, and C. C. Ober. Verification of Euler/Navier- stokes codes using the method of manufactured

• solutions. International Journal for Numerical Methods in Fluids,

44(6):599--620, 2004. 4.

• M. Drela, H. Youngren. http://web.mit.edu/drela/Public/web/

xfoil/

31

Thanks

32