RADGUM: The Recovery-Assisted DG code of the - PowerPoint Presentation

Slide1

RADGUM:The Recovery-Assisted DG code of theUniversity of Michigan(WS1 case only)

January 6

th

, 2017

5

th

International Workshop on High-Order CFD Methods

Kissimmee, Florida

Philip E. Johnson & Eric Johnsen

Scientific Computing and Flow Physics Laboratory

Mechanical Engineering Department

University of Michigan, Ann Arbor

Slide2

Code Overview

Basic Features:

Spatial Discretization:

Discontinuous Galerkin, nodal basisTime Integration: Explicit Runge-Kutta (4th order and 8th order available)Riemann solver: Roe, SLAU2†Quadrature: One quadrature point per basis functionNon-Standard Features:ICB reconstruction: compact technique, adjusts Riemann solver argumentsCompact Gradient Recovery (CGR): Mixes Recovery with traditional mixed formulation for viscous termsShock Capturing: PDE-based artificial dissipation inspired by C-method†† of Reisner et al.Discontinuity Sensor: Detects shock/contact discontinuities, tags “troubled” elements

†Kitamura & Shima, JCP 2013††Reisner et al., JCP 2013

1

Slide3

Ω

A

Exact Distribution U

DG solution: ,

 

Recovered solution:

 

 

Schematic from [Johnson & Johnsen, APS DFD 2015]

Recovery Concept†

†Van Leer & Nomura, AIAA Conf. 2005

Ω

A

Ω

B

Ω

B

Ω

B

Ω

A

2

Slide4

Recovery Demonstration:

 

3

 

 

Slide5

Recovery Demonstration:

 

Recovered solution

(degree

polynomial) more accurate at interface 

3

 

 

Slide6

Recovery Demonstration:

 

ICB reconstructions (degree

equal at closest quadrature points

 

3

 

 

Slide7

For diffusive fluxes: CGR maintains compact stencil†

, offers advantages over BR2

Larger allowable explicit timestep size

Improved wavenumber resolutionFor advection problems: DG weak form: Must calculate flux along interfacesConventional approach (upwind DG): plug in left/right values of DG solutionConventional approach:Our approach: ICB reconstruction scheme††

Replace left/right solution values with ICB reconstruction:Our Approach vs. Conventional DG

†† Khieu & Johnsen, AIAA Aviation 2014

† Johnson & Johnsen, AIAA Aviation 2017

4

Slide8

Taylor-Green Test (WS1)

Code setup

: p2 elements, uniform hex mesh (27 DOF/element), RK4 time integration

Reference result taken from HiOCFD3 workshopOur approach allows larger stable time step

5ICB+CGR: 2.5 CPU-hoursConventional: 9.2 CPU-Hours

ICB+CGR: 75 CPU-hours

Conventional: 304 CPU-Hours

Slide9

Energy Spectrum Computation

Populate velocity

on evenly-spaced 3D grid

Build discrete

;

For each

:

average over entire grid (all

) for velocity correlation

Open

Matlab

 

6

Slide10

Energy Spectrum Computation

Build 3D Fourier transform of each correlation:

,

,

Calculate energy spectrum:

Normalize: scale

to achieve

 

 

7

Slide11

Conclusions

Were the verification cases helpful and which ones were used?

TGV: First 3D simulation, demonstrates value of ICB+CGR for nonlinear problem

What improvements are needed to the test case?TGV: Standardize energy spectrum calculation and make reference data more easily accessibleDid the test case prompt you to improve your methods/solverYes: added 3D capabilityWhat worked well with your method/solver?Feature resolution on Cartesian meshes (ICB very helpful)What improvements are necessary to your method/solver?ICB/CGR robustness on non-Cartesian elements

8

Slide12

9

SciTech Talk

Title:

A Compact Discontinuous Galerkin Method for Advection-Diffusion ProblemsSession: FD-33, High-Order CFD Methods 1Setting: Sun 2, January 10, 9:30 AMAcknowledgementsComputing resources were provided by the NSF via grant 1531752 MRI: Acquisition of Conflux, A Novel Platform for Data-Driven Computational Physics (Tech. Monitor: Ed Walker).

Slide13

References

Kitamura, K. &

Shima

, E., “Towards shock-stable and accurate hypersonic heating computations: A new pressure flux for AUSM-family schemes,” Journal of Computational Physics, Vol. 245, 2013.Reisner, J., Serensca, J., Shkoller, S., “A space-time smooth artificial viscosity method for nonlinear conservation laws,” Journal of Computational Physics, Vol. 235, 2013.Johnson, P.E. & Johnsen, E., “A New Family of Discontinuous Galerkin Schemes for Diffusion Problems,” 23rd AIAA Computational Fluid Dynamics Conference, 2017.Khieu, L.H. & Johnsen, E., “Analysis of Improved Advection Schemes for Discontinuous Galerkin Methods,” 7th AIAA Theoretical Fluid Dynamics Conference

, 2011.Cash, J.R. & Karp, A.H., “A Variable Order Runge-Kutta Method for Initial Value Problems with Rapidly Varying Right-Hand Sides,” ACM Transactions on Mathematical Software, Vol. 16, No. 3, 1990.

Slide14

Spare Slides

Slide15

Vortex Transport Case (VI1)

6

Setup 1

:

,

RK4, SLAU Riemann solver

Setup 2

:

, RK8

†

(13 stages), SLAU Riemann solverICB usage: Apply ICB on Cartesian meshes, conventional DG otherwise

 

 

 

 

 

 

 

 

 

EQ:

Global

error of

:

Convergence:

order

on Cartesian mesh, order

on perturbed quad mesh

 

† Cash & Karp, ACMTMS 1990

Slide16

Shock-Vortex Interaction (CI2)

7

Configurations:

Cartesian (), Cartesian (, Irregular Simplex

Setup: RK4 time integration, SLAU (Cartesian) and Roe (Simplex) Riemann solversShock Capturing: PDE-based artificial dissipationICB usage: Only on Cartesian grids

 

Quad

 

Quad

 

Simplex

 

Slide17

CGR = Mixed Formulation + Recovery

Must choose interface

approximation from available data

BR2: Take average of left/right solutions at the interfaceCompact Gradient Recovery (CGR): = recovered solutionInterface gradient: CGR formulated to maintain compact stencil 

Gradient approximation in

:

 

Weak equivalence with

:

 

Integrate by parts for

weak form

:

 

5

Slide18

Recovery: reconstruction technique introduced by Van Leer and Nomura

†

in 2005

Recovered solution () and DG solution () are equal in the weak senseGeneralizes to 3D hex elements via tensor product basis The Recovery Concept

Representations of

 

 

 

 

Recovered Solution for :

constraints for

:

 

Interface Solution along :

†Van Leer & Nomura, AIAA Conf. 2005

Slide19

Recovery Demonstration: All Solutions

Slide20

Each interface gets a pair of ICB reconstructions, one for each element:

coefficients per element:

Constraints for

(Similar for

Choice of

affects behavior of ICB scheme

Illustration uses

 

The

ICB reconstruction

 

 

 

Example:

(2 DOF/element)

 

 

Slide21

The

Function: ICB-Modal vs. ICB-Nodal

 

ICB-Modal (original): is lowest mode in each element’s solutionICB-Nodal (new approach): is degree

Lagrange interpolantUse Gauss-Legendre quadrature nodes as interpolation pointsTake nonzero at closest quadrature point

 

Sample

choice for

:

Each

is unity at quadrature point nearest interface

 

Slide22

The

Function: ICB-Modal vs. ICB-Nodal

 

ICB-Modal: Each matches the average of in neighboring cell 

ICB-Nodal: Each matches

at near quadrature point

 

Slide23

Fourier analysis performed on 2 configurations:Conventional: Upwind DG + BR2

New: ICB-Nodal + CGR

Linear advection-diffusion, 1D:

Define element Peclet number:Set Initial condition:Cast numerical scheme in matrix-vector form: Fourier Analysis

Scheme

uDG

+ BR2

ICB

+ CGR

Scheme

uDG

+ BR2

ICB

+ CGR

Analysis Procedure

†

:

† Watkins et al., Computers & Fluids 2016

Slide24

Eigenvalue corresponding to exact solution:

Fourier Analysis

Diagonalize

the update matrix: Calculate initial expansion weights, : Watkins et al. derived estimate for initial error growth: = eigenvalue of  

Eigenvalue Example:

ICB+CGR,

,

,

 

† Watkins et al., Computers & Fluids 2016

Slide25

Wavenumber Resolution

To calculate wavenumber resolution:

Define some error tolerance(

and Peclet number ()Identify cutoff wavenumber, according to:Calculate resolving efficiency:

 

† Watkins et al., Computers & Fluids 2016

Slide26

Scheme Comparison:

 

P

Conventional

ICB

+ CGR

1

0.0296

0.1103

2

0.0531

0.0776

3

0.0844

0.1113

4

0.1022

0.1225

5

0.1196

0.1304

P

Conventional

ICB

+ CGR

1

0.0940

0.2389

2

0.1200

0.1793

3

0.1451

0.1755

4

0.1677

0.2628

5

0.1743

0.1874

Fourier analysis, Linear advection-diffusion

Resolving efficiency measures effectiveness of update scheme’s consistent eigenvalue

Slide27

Compact Gradient Recovery (CGR) Approach

Similar to BR2: Manage flow of information by altering gradient reconstruction

1D Case shown for simplicity: Let

, be gradient reconstructions in

Perform Recovery over

,

for

on the shared interface

 

 

 

Slide28

Example with

elements:

Representations of

 

Ω

A

Ω

B

Process Description:

Start with the DG polynomials

in

and

in

.

 

The ICB Approach (Specifically,

ICBp

[0])

Recovery is applicable ONLY for viscous terms; unstable for advection terms.

Interface-Centered Binary (ICB) reconstruction scheme modifies Recovery approach for hyperbolic PDE.

Slide29

Process Description:

Start with the DG polynomials

in

and in .

Obtain reconstructed solution in

, containing

DOF.

 

Example with

elements:

Representations of

 

Ω

A

Ω

B

The ICB Approach (Specifically,

ICBp

[0])

Slide30

Process Description:

Start with the DG polynomials

in

and in

.Obtain reconstructed solution

in

, containing

DOF.

Perform similar operation for

Use ICB solutions as inputs to

 

Example with

elements:

Representations of

 

Ω

A

Ω

B

The ICB Approach (Specifically,

ICBp

[0])

ICB Method achieves

order of accuracy

Generalizes to 2D via tensor-product basis

 

Slide31

Discontinuity Sensor

Approach:

Check cell averages for severe density/pressure jumps across element interfaces

Calculate =cell average for each elementAt each interface, use sensor of Lombardini to check for shock wave:If Lax entropy condition satisfied (hat denotes Roe average at interface):Check pressure jump:If , tag both elements as “troubled”At each interface, check for contact discontinuityCalculate wave strength propagating the density jump:Check relative strength:If

, tag both elements as “troubled”