[PPT] - On the Importance of Adaptivity for Higher-order Discretizations in PowerPoint Presentation

SLIDE 1

On the Importance of Adaptivity for Higher-order Discretizations in Aerospace Applications

Joshua Krakos, Eric Liu, Huafei Sun, Masayuki Yano David Darmofal, Laslo Diosady, Bob Haimes Aerospace Computational Design Laboratory Massachusetts Institute of Technology August 6, 2012

SLIDE 2

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Motivation

CFD solutions for complex problems have been made possible by increases in computational and algorithmic power

Mach number distribution The Launch Abort Vehicle simulation above took 30 minutes to complete on 16 CPUs (Nemec et al, 2008, CART3D group)

SLIDE 3

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Mesh “convergence” comparison

(Chaffin, 2009, DPW4 Presentation)

Same CFD code (NSU3D) run on two “best practice” meshes of about 40 million nodes

Which solution is most realistic? How would this level of uncertainty be detected in practice? NASA mesh CD=277 counts CESSNA mesh CD=262 counts

SLIDE 4

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Error met or time over

Output-based adaptation

Objective: Increase reliability of CFD by estimating and autonomously controlling error in outputs (e.g. drag or lift)

Problem
Output
Max output error
Max time

Estimate error in

utputs

Calculate flow and

utputs

Adapt grid to control error

Flow solution
Outputs
Estimated errors

?

SLIDE 5

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Higher-order Discontinuous Galerkin Finite Element Method (DGFEM)

Approximations are degree p polynomials within

elements but discontinuous between elements DGFEM approximation: Find uh,p ∈ Vh,p such that Rh,p(uh,p, vh,p) = 0, ∀vh,p ∈ Vh,p

u(x,y) x y TH

uh,p(x, y)

SLIDE 6

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Why higher order for CFD?

Higher-order methods known to be more

efficient than lower-order for problems with smooth flows

Aerospace flows typically have limited

smoothness

Can higher-order methods be beneficial in

aerospace applications?

Adaptation key to realizing benefits of higher-
rder discretization on practical problems

SLIDE 7

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Outputs & adjoints

Jh,p(uh,p) is an output such as lift, drag, ef-

ficiency, etc.

Consider a (infinitessimal) perturbation to

the residual such that, Rh,p(uh,p + δuh,p, vh,p) + (δr, vh,p) = 0

The adjoint ψh,p ∈ Vh,p is the sensitivity of

the output to a residual perturbation, δJh,p ≡ (δr, ψh,p)

Interpretation: adjoint is transfer function

between δr and Jh,p.

In the infinite-dimensional case, the adjoint

satisfies a linear PDE.

Primal (Mach) Transonic RANS example Drag adjoint (mass)

SLIDE 8

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Outputs & adjoints

Jh,p(uh,p) is an output such as lift, drag, ef-

ficiency, etc.

Consider a (infinitessimal) perturbation to

the residual such that, Rh,p(uh,p + δuh,p, vh,p) + (δr, vh,p) = 0

The adjoint ψh,p ∈ Vh,p is the sensitivity of

the output to a residual perturbation, δJh,p ≡ (δr, ψh,p)

Interpretation: adjoint is transfer function

between δr and Jh,p.

In the infinite-dimensional case, the adjoint

satisfies a linear PDE.

Primal (Mach) Transonic RANS example Drag adjoint (mass)

Drag error indicator, ηκ

Transonic RANS example

SLIDE 9

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Continuous Optimization: Mesh-metric Duality

0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1

Metric field, M(x)

Mesh, Th Mesh generation Implied metric

(Intractable) discrete optimization problem

T ∗

h = arg inf Th

E(Th) s.t. C(Th) = Cost

Continuous relaxation (Loiselle, 2009)

M∗ = arg inf

M

E(M) s.t. C(M) = Cost

SLIDE 10

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Local sampling

       

For each configuration, solve local problems keeping states outside of κ0 fixed
Determine error estimate ηκi = Rh,p(uκi

h,p, ψh,p+1|κ0)

Produces a set of pairs, {Mκi, ηκi}

Yano & Darmofal, 2012

SLIDE 11

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

MOESS Algorithm

(Mesh Optimization via Error Sampling & Synthesis)

Solve flow and adjoint on current grid
Determine error-metric gradients via local sampling
Utilize steepest descents algorithm to improve metric
Remesh using improved metric

Yano & Darmofal, 2012

SLIDE 12

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Impact of Adaptation on Higher-order Efficiency

Subsonic Euler

When singularities are present, adaptive refinement critical to realize benefits of higher order

NACA 0012, M = 0.5, α = 2

2.5k 5k 10k 20k 10

−8

10

−7

10

−6

10

−5

10

−4

10

−3

degrees of freedom cd error estimate −1.44 −0.69 −3.28 p=1 (uniform) p=3 (uniform) p=1 (adapt) p=3 (adapt)

Adaptive refinement is per-

formed at 2, 500 and 5, 000 DOFs generating “optimal” meshes

Uniform refinement (each el-

ement subdivided into four) is performed

Uniform refinement compared

to adaptive refinement at 10, 000 and 20, 000 DOFs

SLIDE 13

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Impact of Adaptation on Higher-order Efficiency

Subsonic Euler

Adaptive mesh resolves trailing edge singularity more effectively

20K DOF mesh uniformally refined from 5K DOF mesh

Optimized 20K DOF mesh

NACA 0012, M = 0.5, α = 2 p = 3: Distribution of (log) elemental error (log10 ηκ)

SLIDE 14

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Adaptation, higher-order, and RANS

Subsonic RANS

Uniform refinement of 40K DOF to 160K DOF Error indicator for160K DOF, p=3

Yano, Modisette, Darmofal 2011

20k 40k 80k 160k 10

−3

10

−2

10

−1

10 10

1

degrees of freedom cd error estimate (counts) −0.98 −3.28 −1.71 p=1 (uniform) p=3 (uniform) p=1 (adapt) p=3 (adapt)

RAE2822, RANS, M = 0.3, Re = 6.5 × 106, α = 2.31

Adaptive 160K DOF

To see full benefit of higher-order approximations, solution irregularities must be controlled: adaptation is critical

SLIDE 15

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

MDA-3 RANS

Lift adaptation, boundary conforming

M = 0.2, Re = 9 × 106, 11 α values from 0 to 24.5

5 10 15 20 25 2 2.5 3 3.5 4 4.5 5 angle of attack cl 0.02 0.04 0.06 0.08 0.1 0.12 2 2.5 3 3.5 4 4.5 5 cd cl

α=0.0 α=8.1 α=16.21 α=23.28

p = 2 results

SLIDE 16

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

MDA-3 RANS

Lift adaptation, boundary conforming

M = 0.2, Re = 9 × 106, 11 α values from 0 to 24.5

5 10 15 20 25 10

−3

10

−2

10

−1

10 10

1

10

2

angle of attack cl error estimate fixed mesh adaptive 5 10 15 20 25 2 2.5 3 3.5 4 4.5 5 angle of attack cl fixed mesh adaptive

An appropriate mesh is critical (and subtle)
Consider α sweep using α = 8.1 optimized grid

cl cl error estimate

SLIDE 17

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Multi-element Airfoil RANS Case: Higher-order Workshop (2012)

cd convergence with WU

10

2

10

3

10

4

10

5

10

−6

10

−5

10

−4

10

−3

10

−2

10

−1

10

p = 1 |drag coefficient error| Work units

MIT UMich Bergano Wyoming

10

2

10

3

10

4

10

5

10

−6

10

−5

10

−4

10

−3

10

−2

10

−1

10

p = 2 |drag coefficient error| Work units

MIT UMich Bergano Wyoming

~1% error ~10% error

1000 work units

UMich and Bergamo results use

uniform mesh refinement

Wyoming is an hp adaptive result
MIT (h-adaptive) results 10-100 x

more efficient

MIT work includes adaptive time

from a coarse initial mesh

With adaptivity, p=2 error drops

rapidly with small additional work

SLIDE 18

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Cut-cell vs. boundary conforming

Transonic RANS

RAE2822, RANS, M = 0.729, Re = 6.5 × 106, α = 2.31

20k 40k 80k 160k 10

−7

10

−6

10

−5

10

−4

10

−3

10

−2

degrees of freedom cd error estimate −1.17 −2.08 −2.53 p=1 (cutcell) p=2 (cutcell) p=3 (cutcell) p=1 (conforming) p=2 (conforming) p=3 (conforming)

Final grid, p=3: 40K DOF Mach Initial grid

Yano, Modisette, Darmofal 2011

SLIDE 19

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Laminar delta wing example

Higher Order Workshop Test Case 2012:

Adaptation critical for achieving higher-order convergence (apriori meshes are uniform refinements)

M = 0.3, Re = 4000, α = 12.5

10

4

10

5

10

6

10

−5

10

−4

10

−3

10

−2

|CD − CD

ref|

dof 1% error 0.1% error HOW mesh (p=1) HOW mesh (p=2) L&H (p=1, hexa) MOESS (p=1) MOESS (p=2)

SLIDE 20

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Key take-away

Adaptation is critical to realize the performance benefits of higher-order discretizations on aerospace applications

SLIDE 21

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Adaptive higher-order methods: Current Status & Challenges

Higher-order 2D & 3D RANS (including shocks)

demonstrated on a priori meshes (Bassi; Darmofal; Fidkowski; Hartmann; Mavriplis; Peraire, etc.)

Robust anisotropic adaptation demonstrated for

2D steady RANS (Darmofal; Fidkowski; Hartmann)

Proof of concept demonstrations for adaptive 3D

RANS (Darmofal; Fidkowski; Hartmann)

Challenges: higher-order adaptive meshing;

robustness for under-resolved RANS; unsteadiness

SLIDE 22

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Questions?

SLIDE 23

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Affine-invariant metric framework

Employ affine-invariant description of a metric space (Pennec et al, 2006)

Sκ = log ⇣ M−1/2

κ0

MκM−1/2

κ0

⌘

Sκ (the step matrix) can be decomposed into Sκ = sκI + ˜

Sκ

sκ is isotropic and controls the area change
˜

Sκ controls orientation and stretching changes

First-order optimality conditions become

∂ηκ ∂sκ − λ∂ρκ ∂sκ = ∂ηκ ∂ ˜ Sκ =

SLIDE 24

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Error model synthesis

Define logarithmic error model fκi ≡ log(ηκi/ηκ0)

{Mκi, ηκi} → {Sκi, fκi}

Perform a least-squares fit to synthesis fκ(Sκ) = tr(RκSκ):

Rκ = arg min

Q∈Symd nconfig

X

i=1

(fκi − tr(QSκi))2

This gives ηκ(Sκ) = ηκ0 exp(rκsκd) exp

⇣ tr ⇣ ˜ Rκ ˜ Sκ ⌘⌘

For isotropic error and meshing this model reduces to,

ηiso

κ (h) = ηκ0

✓ h h0 ◆riso

κ

SLIDE 25

AEROSPACE COMPUTATIONAL DESIGN LABORATORY

Continuous Metric Optimization

Problem Statement: Seek optimal metric field

M∗ = arg min

M E(M)

s.t. C(M) = Cost

Choose E(M) = P

κ ηκ(Mκ)

Choose C(M) = P

κ ρκ(Mκ)

ρκ are the DOF in region κ.