Dimitri J. Mavriplis Department of Mechanical Engineering University of Wyoming Laramie, WY

1. Sensitivity Analysis Methods For Optimization, Error Estimation, and Uncertainty Quantification Dimitri J. Mavriplis Department of Mechanical Engineering University of Wyoming Laramie, WY

2. Motivation • Computational fluid dynamics analysis capabilities commonplace today • In addition to analysis capability, sensitivity capability is highly desirable • Design optimization • Error estimation • Parameter sensitivity • Uncertainty quantification • Sensitivities may be obtained by: • Perturb input, rerun analysis code (Finite difference) • Linearizing analysis code (tangent method) • Good for 1 input, many outputs • Adjoint method • Good for many inputs, one output

3. Overview • Adjoint methods well known for steady-state single disciplinary problems • Design optimization • Spatial discretization error/AMR • Extend to more complex problems • Unsteady flow problems • Multidisciplinary problems • Other error sources • Temporal error • Algebraic error • Coupling error • Uncertainty quantification • Investigate more sophisticated approaches: 2nd order • Hessian-based methods

4. Motivation • Continuous vs. Discrete Adjoint Approaches • Continuous: Linearize then discretize • Discrete: Discretize then Linearize • Continuous Approach: • More flexible adjoint discretizations • Framework for non-differentiable tasks (limiters) • Often invoked using flow solution as constraint using Lagrange multipliers • Discrete Approach: • Reproduces exact sensitivities of code • Verifiable through finite differences • Relatively simple implementation • Chain rule differentiation of analysis code • Transpose these derivates • (transpose and reverse order) • Includes boundary conditions

5. HIARMS Components

6. P0 P1 P2 PN PUNDIT: HI-ARMS Rotor-FSI: HI-ARMS MDM Mesh Deform 6DOF Mesh Motion Component Interfaces FSI Fluid-Structure Interface DCF Domain Connectivity Object Oriented Python Integration Framework Software Integration Framework (SIF) Distributed Memory processors communicating via MPI 20 wing twist shared data NBE Near-Body CFD OBE Off-Body CFD CSD Struct Dynamics NSU3D: Univ. of Wyoming SAMARC: LLNL and NASA-Ames RCAS: AFDD and ART HIARMS Components

7. Generalized Discrete Sensitivities • Consider a multi-phase analysis code: • L = Objective(s) • D = Design variable(s) • Sensitivity Analysis • Using chain rule:

8. Tangent Model • Special Case: • 1 Design variable D, many objectives L • Precompute all stuff depending on single D • Construct dL/dD elements as:

9. Adjoint Model • Special Case: • 1 Objective L, Many Design Variables D • Would like to precompute all left terms • Transpose entire equation:

10. Adjoint Model • Special Case: • 1 Objective L, Many Design Variables D • Would like to precompute all left terms • Transpose entire equation: precompute as:

11. Shape Optimization Problem • Multi-phase process:

12. Tangent Problem (forward linearization) • Examine Individual Terms: • : Design variable definition (CAD) • : Objective function definition

13. Tangent Problem (forward linearization) • Examine Individual Terms:

14. Sensitivity Analysis • Tangent Problem: • Adjoint Problem

15. Tangent Problem • 1: Surface mesh sensitivity: • 2: Interior mesh sensitivity: • 3: Residual sensitivity: • 4: Flow variable sensitivity: • 5: Final sensitivity

16. Adjoint Problem • 1: Objective flow sensitivity: • 2: Flow adjoint: • 3:Objective sens. wrt mesh: • 4: Mesh adjoint: • 5: Final sensitivity:

17. General Approach • Linearize each subroutine/process individually in analysis code • Possible use of automatic differentiation at this stage (only) • Check linearization by finite difference • Transpose, and check duality relation • Build up larger components • Check linearization, duality relation • Check entire process for FD and duality • Use single modular AMG solver for all phases

18. General Duality Relation • Analysis Routine: • Tangent Model: • Adjoint Model: • Duality Relation: • Necessary but not sufficient test • Check using series of arbitrary input vectors

19. Drag Minimization Problem • Total Optimization Time for 15 Design Cycles: 6 hours on 16 cpus of PC cluster • Flow Solver: 150 MG cycles • Flow Adjoint: 50 Defect-Correction cycles (x 4 MG) • Mesh Adjoint: 25 MG cycles • Mesh Motion: 25 MG cycles

20. Extension to Unsteady Problems:Governing Equations (ALE) • In ALE Form: • V(t) = control volume • Face integrated mesh velocity Formulated to obey GCL = f(xn,xn-1,xn-2,…) Mavriplis and Yang (AIAA 2005-5114)

21. Time Discretization • BDF1: • BDF2:

22. Unsteady Residual Form • BDF1: Due uniquely to ALE grid speed terms • BDF2: Similar expression depending on xn,xn-1,xn-2

23. Time-Dependent Adjoint • Adjoint equations now obtained using flow and mesh constraints at each time step: Objective at final time

24. Shape Optimization • Per design cycle (Steady Case) • One mesh deformation problem • One flow analysis • One flow adjoint solution • One mesh adjoint solution • Unsteady Shape Optimization • General functional dependence involves previous time step values • Chain rule results in forward time recurrence relation • When transposed (adjoint) results in backwards “integration” in time • Time history of solution must be stored for use by adjoint in reverse time integration • Write out to local cluster disks, read back in during adjoint phase

25. Validation of Sensitivities • Steady/Unsteady sensitivities compare well with finite difference values • Tangent/Adjoint values equivalent to machine precision • Duality principle

26. Extension to Unsteady Problems Pressure Contours for Pitching Airfoils Minf = 0.755, a0 = 0.016o, amax = 2.51o, w = 0.1628, t=0 to 54 27 time-steps with dt=2.0 Optimized Airfoil NACA0012 Baseline Airfoil

27. Time-Dependent Load Convergence/Comparison 90 optimization steps using LBFGS Mani and Mavriplis AIAA 2007-0060

28. Extension to Multidisciplinary Problems: Aeroelasticity • Formulation leads to : • Disciplinary adjoints: Fluids, Structures • Disciplinary adjoints are coupled at each time step • Coupled adjoint solver analogous (transpose) of coupled aeroelastic analysis solver • Used to demonstrate flutter suppression through shape optimization in 2D • Mani and Mavriplis AIAA 2008-6242

29. Adjoint-Based Error Estimation • Complex simulations have multiple error sources • Engineering simulations concerned with specific output objectives • Adjoint methods / Goal Oriented Approach • Use for a posteriori error estimation of specific objectives • Spatial error • Temporal error • Other error sources • Use to drive adaptive process

30. Adjoint-Based Spatial Error Estimation 30 • Formulation • Taylor series Li Wang, Dimitri Mavriplis (UW)

31. ADJOINT-BASED ERROR ESTIMATION • Approximated objective becomes error = adjoint . residual 31 • Formulation • Avoid solving the adjoint variables on the fine mesh, instead, • Solve on the coarse mesh • Reconstruct onto the fine mesh by using least squares method

32. Summary of Spatial Error Estimation and Refinement • Compute steady flow solution on coarse mesh H • Compute adjoint variables on coarse mesh H • Project adjoint variables, flow solution and mesh solution onto fine mesh h • Spatial error is then inner product of adjoint with corresponding non-zero residual on fine mesh h • Provides a prediction of objective value on fine mesh • Distribution of error in space is used to drive adaptation

33. Adjoint-Based Spatial Error Estimation + AMR Mach number contours Adjoint solution, Λ(2) • Adjoint Solution : Green’s Function for Objective (Lift) • Change in Lift for Point sources of Mass/Momentum • Error in objective ~ Adjoint . Residual (approx. solution) • Predicts objective value for new solution (on finer mesh) • Cell-wise indicator of error in objective (only)

34. h-refinement for target functional of lift Final h-adapted mesh (8387 elements) Close-up view of the final h-adapted mesh Fixed discretization order of p = 1

35. h-refinement for target functional of lift Error convergence history vs. degrees of freedom Functional Values and Corrected Values Comparison between h-refinement and uniform mesh refinement

36. Complex Geometry: Vehicle Stage Separation(CART3D/inviscid) Top View Initial Mesh Side View • Initial mesh contains only 13k cells • Final meshes contain between 8M to 20M cells

37. Pressure Contours M∞=4.5, α=0°

38. Minimal refinement of inter-stage region • Gap is highly refined • Overall, excellent convergence of functional and error estimate Cutaway view of inter-stage

39. Unsteady Multidisciplinary Problems Total error in solution Algebraic error Spatial error (discretization/resolution) Temporal error (discretization/resolution) Flow Mesh Other Flow Mesh Other Flow Mesh Other • Solution of time-dependent adjoint: backwards integration in time • Disciplinary adjoint inner product with disciplinary residual

40. Time-integrated functional • Interaction of isentropic vortex with slowly pitching NACA0012 • Mach number = 0.4225 • Reduced frequency = 0.001 • Center of pitch is quarter chord • Functional is 8,600 elements

41. Unsteady Adjoint Error Estimation Density contours of initial condition

42. Summary of Total Error Evaluation and Decomposition • Compute partially converged flow and mesh solution on coarse time domain • Compute adjoint variables on coarse domain using partially converged solution • Compute partial convergence error on coarse level time domain • Inner product of adjoint with partially converged (non-zero) residual • Project partially converged solution and adjoint variables onto fine time domain • Evaluate fine level error estimate as inner product between adjoint and residual on fine time domain • Combined temporal resolution and partial convergence error • Determine temporal resolution error by subtracting partial convergence error from total error estimate on fine time domain

43. Adaptation • Compute time-integrated averages of error component distributions • Adapt where error is greater than time-integrated average • Time resolution error: divide time step by 2 • Convergence error: tighten tolerance by predetermined factor

44. Temporal Error Adaptation Comparison of adapted temporal domain

45. Algebraic Error Adaptation Adapted Flow/Mesh convergence tolerances:

46. Adjoint-Based Refinement Results • Error in Lift versus CPU Time • Uniform cost is only finest solution cost • Adaptive cost is all solutions (+ adjoint cost) • Corrected value provides further improvement

47. Second-Order Methods • Adjoint is efficient approach for calculating first-order sensitivities (first derivatives) • Second-order (Hessian) information can be useful for enhanced capabilities: • Optimization • Hessian corresponds to Jacobian of optimization problem (Newton optimization) • Unsteady optimization seems to be hard to converge • Optimization for stability derivatives • Optimization under uncertainty • Uncertainty quantification • Method of moments (Mean of inputs = input of means) • Inexpensive Monte-Carlo (using quadratic extrapolation)

48. Forward-Reverse Hessian Construction • Hessian for N inputs is a NxN matrix • Complete Hessian matrix can be computed with: • One tangent/forward problem for each input • One adjoint problem • Inner products involving local second derivatives computed with automatic differentiation • Overall cost is N+1 solves for NxN Hessian matrix • Lower than double finite-difference: O(N2) • May be impractical for large number of inputs/design variables

49. Hessian Implementation • Implemented for steady and unsteady 2D airfoil problems • Validated against double finite difference for Hicks-Henne bump function design variables

50. Newton Optimization with Hessian • LBFGS is “best” gradient-based optimizer • Constructs approximate Hessian based on previous design iterations • KNITRO is Newton optimizer • Requires Hessian as input • Superior performance in terms of number of function calls • Added cost of Hessian recovered (2 to 6 design variables)