Scalable methods for optimal control of systems governed by PDEs - - PowerPoint PPT Presentation

scalable methods for optimal control of systems governed
SMART_READER_LITE
LIVE PREVIEW

Scalable methods for optimal control of systems governed by PDEs - - PowerPoint PPT Presentation

Scalable methods for optimal control of systems governed by PDEs with random coefficient fields Alen Alexanderian, 1 Peng Chen, 2 Omar Ghattas , 2 emi Petra 3 Georg Stadler, 4 Umberto Villa 2 No 1 Department of Mathematics North Carolina State


slide-1
SLIDE 1

Scalable methods for optimal control of systems governed by PDEs with random coefficient fields

Alen Alexanderian,1 Peng Chen,2 Omar Ghattas,2 No´ emi Petra3 Georg Stadler,4 Umberto Villa2

1Department of Mathematics

North Carolina State University

2Institute for Computational Engineering and Sciences

The University of Texas at Austin

3School of Natural Sciences

University of California, Merced

3Courant Institute of Mathematical Sciences

New York University

April 20, 2017 ICES Babuˇ ska Forum The University of Texas at Austin

slide-2
SLIDE 2

From data to decisions under uncertainty

Uncertain parameter m Mathematical Model A(m, u) = f Quantity of Interest (QoI) q(m)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 2 / 41

slide-3
SLIDE 3

From data to decisions under uncertainty

Uncertain parameter m Bayesian inversion πpost(m|y) ∝ πlike(y|m)π0(m) Experimental data y Mathematical Model A(m, u) = f Quantity of Interest (QoI) q(m) posterior prior

πlike(y|m)=πnoise(Bu − y)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 2 / 41

slide-4
SLIDE 4

From data to decisions under uncertainty

Uncertain parameter m Bayesian inversion πpost(m|y) ∝ πlike(y|m)π0(m) Experimental data y Mathematical Model A(m, u) = f Design of experiments min

ξ

  • Ψ
  • πpost(m|y; ξ)
  • π(y)dy

ξ : experimental design Ψ: design objective/criterion Quantity of Interest (QoI) q(m) posterior prior experiments

πlike(y|m)=πnoise(Bu − y)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 2 / 41

slide-5
SLIDE 5

From data to decisions under uncertainty

Uncertain parameter m Bayesian inversion πpost(m|y) ∝ πlike(y|m)π0(m) Experimental data y Mathematical Model A(m, u) = f(·; z) z := control function Design of experiments min

ξ

  • Ψ
  • πpost(m|y; ξ)
  • π(y)dy

ξ : experimental design Ψ: design objective/criterion Quantity of Interest (QoI) q(z, m) Optimal control/design under uncertainty e.g. min

z

  • q(z, m)µ(dm)

posterior prior experiments

πlike(y|m)=πnoise(Bu − y)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 2 / 41

slide-6
SLIDE 6

Example: Groundwater contaminant remediation

Source: Reed Maxwell, CSM

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 3 / 41

slide-7
SLIDE 7

Example: Groundwater contaminant remediation

Inverse problem

Infer (uncertain) soil permeability from (uncertain) measurements of pressure head at wells and from a (uncertain) model of subsurface flow and transport

Prediction (or forward) problem

Predict (uncertain) evolution of contaminant concentration at municipal wells from (uncertain) permeability and (uncertain) subsurface flow/transport model

Optimal experimental design problem

Where should new observation wells be placed so that permeability is inferred with the least uncertainty?

Optimal design problem

Where should new remediation wells be placed so that (uncertain) contaminant concentrations at municipal wells are minimized?

Optimal control problem

What should the rates of extraction/injection at remediation wells be so that (uncertain) contaminant concentrations at municipal wells are minimized?

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 4 / 41

slide-8
SLIDE 8

Applications of inverse problems in CCGO

Antarctic ice sheet flow (+ ocean dynamics) Joint with Patrick Heimbach, Tom Hughes, Tobin Isaac (Georgia Tech), Tom O’Leary-Roseberry, Noemi Petra (UC-Merced), Georg Stadler (NYU), Umberto Villa, Alice Zhu Global and regional seismic inversion, joint seismic–EM inversion, inverse scattering Joint with Hossein Aghakhani, Nick Alger, Tan Bui, Ben Crestel, David Keyes (KAUST), George Turkiyyah (KAUST), Georg Stadler (NYU), Umberto Villa Global mantle convection Joint with Mike Gurnis (Caltech), Johann Rudi, Georg Stadler (NYU) Poroelastic subsurface flow inversion and management of induced seismicity Joint with Amal Alghamdi, Marc Hesse, Georg Stadler (NYU), Umberto Villa, Karen Willcox (MIT) Turbulent combustion: inference and control Joint with George Biros, Peng Chen, Matthias Heinkenschloss (Rice), Myoungkyu Lee, Bob Moser, Todd Oliver, Chris Simmons, David Sondak, Andrew Stuart (Caltech), Umberto Villa, Karen Willcox (MIT) Reservoir inversion Joint with George Biros, Tan Bui, Clint Dawson, Sam Estes, John Lee, Umberto Villa Soft tissue biomechanical inversion Joshua Chen, Michael Sacks, Umberto Villa

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 5 / 41

slide-9
SLIDE 9

Forward and inverse global mantle convection modeling

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 6 / 41

slide-10
SLIDE 10

Scalable solver (2015 Gordon Bell Prize)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 7 / 41

slide-11
SLIDE 11

Bayesian inversion for basal friction field in Antarctica

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 8 / 41

slide-12
SLIDE 12

Bayesian global seismic inversion

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 9 / 41

slide-13
SLIDE 13

Bayesian poroelastic inversion

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 10 / 41

slide-14
SLIDE 14

Joint seismic-electromagnetic inversion

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 11 / 41

slide-15
SLIDE 15

Optimal control of systems governed by PDEs with uncertain parameter fields

PDE-constrained control objective: q = q(u(z, m), z) where u depends on z and m through: A(u, m) = f(z) q: control objective A: forward PDE operator u: state variable m: uncertain parameter field z: control function Control of injection wells in porous medium flow (SPE10 permeability data) Problem: given uncertainty model for m, find z that “optimizes” q(u(z, m), z)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 12 / 41

slide-16
SLIDE 16

Optimization under uncertainty (OUU)

H : parameter space, infinite-dimensional separable Hilbert space q(z, m): control objective functional m ∈ H : uncertain model parameter field, z: control function Optimization under uncertainty (OUU): min

z

q(z, m)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 13 / 41

slide-17
SLIDE 17

Optimization under uncertainty (OUU)

H : parameter space, infinite-dimensional separable Hilbert space q(z, m): control objective functional m ∈ H : uncertain model parameter field, z: control function Risk-neutral optimization under uncertainty (OUU): min

z

Em{q(z, m)} Em{q(z, m)} =

  • H

q(z, m) µ(dm)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 13 / 41

slide-18
SLIDE 18

Optimization under uncertainty (OUU)

H : parameter space, infinite-dimensional separable Hilbert space q(z, m): control objective functional m ∈ H : uncertain model parameter field, z: control function Risk-averse (Mean-Var) optimization under uncertainty (OUU): min

z

Em{q(z, m)} + β varm{q(z, m)} Em{q(z, m)} =

  • H

q(z, m) µ(dm) varm{q(z, m)} = Em{q(z, m)2} − Em{q(z, m)}2

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 13 / 41

slide-19
SLIDE 19

Optimization under uncertainty (OUU)

H : parameter space, infinite-dimensional separable Hilbert space q(z, m): control objective functional m ∈ H : uncertain model parameter field, z: control function Risk-averse (Mean-Var) optimization under uncertainty (OUU): min

z

Em{q(z, m)} + β varm{q(z, m)} Em{q(z, m)} =

  • H

q(z, m) µ(dm) varm{q(z, m)} = Em{q(z, m)2} − Em{q(z, m)}2 Main challenges:

Integration over infinite/high-dimensional parameter space Evaluation of q requires PDE solves

Standard Monte Carlo approach (Sample Average Approximation) is prohibitive

Numerous (nmc) samples required, each requires PDE solve Resulting PDE-constrained optimization problem has nmc PDE constraints

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 13 / 41

slide-20
SLIDE 20

Some existing approaches for PDE-constrained OUU

Methods based on stochastic collocation, sparse/adaptive sampling, POD, ...

Schulz & Schillings, Problem formulations and treatment of uncertainties in aerodynamic design, AIAA J, 2009. Borz` ı & von Winckel, Multigrid methods and sparse-grid collocation techniques for parabolic optimal control problems with random coefficients, SISC, 2009. Borz` ı, Schillings, & von Winckel, On the treatment of distributed uncertainties in PDE-constrained optimization, GAMM-Mitt. 2010. Borz` ı & von Winckel, A POD framework to determine robust controls in PDE optimization, Computing and Visualization in Science, 2011. Gunzburger & Ming, Optimal control of stochastic flow over a backward-facing step using reduced-order modeling, SISC 2011. Hou, Lee, & Manouzi, Finite element approximations of stochastic optimal control problems constrained by stochastic elliptic PDEs, J Math Anal Appl, 2011. Gunzburger, Lee, & Lee, Error estimates of stochastic optimal Neumann boundary control problems, SINUM, 2011. Rosseel & Wells, Optimal control with stochastic PDE constraints and uncertain controls, CMAME, 2012. Tiesler, Kirby, Xiu, & Preusser, Stochastic collocation for optimal control problems with stochastic PDE constraints, SICON, 2012. Kouri, Heinkenschloss, Ridzal, & Van Bloemen Waanders, A trust-region algorithm with adaptive stochastic collocation for PDE optimization under uncertainty, SISC, 2013. Chen, Quarteroni, & Rozza, Stochastic optimal Robin boundary control problems of advection-dominated elliptic equations, SINUM, 2013. Kunoth & Schwab, Analytic regularity and gPC approximation for control problems constrained by linear parametric elliptic and parabolic PDEs, SICON, 2013. Kouri, A multilevel stochastic collocation algorithm for optimization of PDEs with uncertain coefficients, JUQ, 2014. Chen & Quarteroni, Weighted reduced basis method for stochastic optimal control problems with elliptic PDE constraint, JUQ, 2014. Ng & Willcox, Multifidelity approaches for optimization under uncertainty, IJNME, 2014. Kouri, Heinkenschloss, Ridzal, & van Bloemen Waanders, Inexact objective function evaluations in a trust-region algorithm for PDE-constrained

  • ptimization under uncertainty, SISC, 2014.

Chen, Quarteroni, & Rozza, Multilevel and weighted reduced basis method for stochastic optimal control problems constrained by Stokes equations, Num. Math. 2015. Ng & Willcox, Monte Carlo information-reuse approach to aircraft conceptual design optimization under uncertainty, J Aircraft, 2015. Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 14 / 41

slide-21
SLIDE 21

Example: Control of injection wells in porous medium

1 2 3 4

¯ m = mean of log permeability field ¯ q = target pressure at production wells

state PDE: single phase flow in a porous medium −∇ · (em∇u) =

nc

  • i=1

zifi(x) with Dirichlet lateral & Neumann top/bottom BCs uncertain parameter: log permeability field m control variables: zi, mass source at injection wells; fi, mollified Dirac deltas control objective: q(z, m) := 1

2Qu(z, m) − ¯

q2, ¯ q: target pressure dimensions: ns = nm = 3242, nc = 20, nq = 12

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 15 / 41

slide-22
SLIDE 22

Porous medium with random permeability field

Distribution law of m: µ = N( ¯ m, C) (Gaussian measure on Hilbert space H ) Take covariance operator as square of inverse of Poisson-like operator: C = (−κ∆ + αI)−2 κ, α > 0 C is positive, self-adjoint, of trace-class; µ well-defined on H (Stuart ’10)

κ α ∝ correlation length; the larger α, the smaller the variance

Random draws for κ = 2 × 10−2, α = 4

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 16 / 41

slide-23
SLIDE 23

OUU with linearized parameter-to-objective map

Mean-var risk-averse optimal control problem (including cost of controls): min

z

Em{q(z, m)} + β varm{q(z, m)} + γ z2 Linear approximation to parameter-to-objective map about ¯ m: qlin(z, m) = q(z, ¯ m) + gm(z, ¯ m), m − ¯ m gm(z, ·) := dq(z, ·) dm is the gradient with respect to m The moments of the linearized objective: Em{qlin(z, ·)} = q(z, ¯ m) varm{qlin(z, ·)} = gm(z, ¯ m), C[gm(z, ¯ m)] qlin(z, ·) ∼ N

  • q(z, ¯

m), gm(z, ¯ m), C[gm(z, ¯ m)]

  • Omar Ghattas (ICES, UT Austin)

Optimal control under uncertainty Mar 24, 2017 17 / 41

slide-24
SLIDE 24

Mean-var risk-averse optimal control problem with linearized parameter-to-objective map

State-and-adjoint-PDE constrained optimization problem (quartic in z):

min

z∈Z J (z) := 1

2Qu − ¯ q2 + β 2 gm( ¯ m), C[gm( ¯ m)] + γ 2 z2 with gm( ¯ m) =e ¯

m∇u · ∇p,

where −∇ · (e ¯

m∇u) = nc

  • i=1

zifi state equation −∇ · (e ¯

m∇p) = −Q∗(Qu − ¯

q) adjoint equation

Lagrangian of the risk-averse optimal control problem with qlin:

L (z, u, p, u⋆, p⋆) = 1 2Qu − ¯ q2 + β 2

  • e ¯

m∇u · ∇p, C[e ¯ m∇u · ∇p]

  • + γ

2 z2 +

  • e ¯

m∇u, ∇u⋆

nc

  • i=1

zifi, u⋆ +

  • e ¯

m∇p, ∇p⋆

+ Q∗(Qu − ¯ q), p⋆

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 18 / 41

slide-25
SLIDE 25

Gradient computation for risk averse optimal control

“State problem” for risk-averse optimal control problem with qlin:

  • e ¯

m∇u, ∇˜

u

  • =

nc

  • i=1

zifi, ˜ u

  • e ¯

m∇p, ∇˜

p

  • = −Q∗(Qu − ¯

q), ˜ p for all test functions ˜ p and ˜ u

“Adjoint problem” for risk-averse optimal control problem with qlin:

  • e ¯

m∇p⋆, ∇˜

p

  • = −β
  • e ¯

m∇u · ∇˜

p, C[e ¯

m∇u · ∇p]

  • e ¯

m∇u⋆, ∇˜

u

  • = −Q∗(Qu − ¯

q), ˜ u − β

  • e ¯

m∇˜

u · ∇p, C[e ¯

m∇u · ∇p]

  • − Q∗Qp⋆, ˜

u for all test functions ˜ p and ˜ u

Gradient: ∂J ∂zj = γzj − fj, u⋆, j = 1, . . . , nc Cost of objective = 2 PDE solves; cost of gradient = 2 PDE solves

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 19 / 41

slide-26
SLIDE 26

Risk-averse optimal control with linearized objective

−20 20 40 60 80 100 120 140 0.01 0.02 0.03 0.04 distribution Θ(z0, ·) Θlin(z0, ·) Θquad(z0, ·)

initial (suboptimal) control z0

  • distrib. of exact & approx objectives at z0

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 20 / 41

slide-27
SLIDE 27

Risk-averse optimal control with linearized objective

−20 20 40 60 80 100 120 140 0.01 0.02 0.03 0.04 distribution Θ(z0, ·) Θlin(z0, ·) Θquad(z0, ·)

initial (suboptimal) control z0

  • distrib. of exact & approx objectives at z0

2 4 6 8 10 12 14 16 1 2 3 4 5 distribution Θ(zopt

lin , ·)

Θlin(zopt

lin , ·)

Θquad(zopt

lin , ·)

  • ptimal control zopt

lin based on qlin

  • distrib. of exact & approx objectives at zopt

lin Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 20 / 41

slide-28
SLIDE 28

Quadratic approximation to parameter-to-objective map

Quadratic approximation to the parameter-to-control-objective map: qquad(z, m) = q(z, ¯ m) + gm(z, ¯ m), m − ¯ m + 1 2Hm(z, ¯ m)(m − ¯ m), m − ¯ m gm: gradient of parameter-to-objective map Hm: Hessian of parameter-to-objective map Observations: Quadratic approximation does not lead to a Gaussian control objective However, can derive analytic formulas for the moments of qquad in the infinite-dimensional Hilbert space setting

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 21 / 41

slide-29
SLIDE 29

Analytic expressions for mean and variance with qquad

Mean: Em{qquad(z, ·)} = q(z, ¯ m)+1 2tr ˜ Hm(z, ¯ m)

  • Variance:

varm{qquad(z, ·)} = gm(z, ¯ m), C[gm(z, ¯ m)]+1 2tr ˜ Hm(z, ¯ m)2 where ˜ Hm = C1/2HmC1/2 is the covariance-preconditioned Hessian Risk averse optimal control objective with qquad: J(z) = q(z, ¯ m) + 1 2tr ˜ Hm(z, ¯ m)

  • + β

2

  • gm(z, ¯

m), C[gm(z, ¯ m)] + 1 2tr ˜ Hm(z, ¯ m)2

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 22 / 41

slide-30
SLIDE 30

Randomized trace estimator

Randomized trace estimation: tr( ˜ Hm) ≈ 1 ntr

ntr

  • j=1
  • ˜

Hmξj, ξj

  • = 1

ntr

ntr

  • j=1

Hmζj, ζj tr( ˜ H2

m) ≈ 1

ntr

ntr

  • j=1

Hmζj, C[Hmζj] where ξj are Gaussian random fields and ζj = C1/2ξj In computations, we use draws ζj ∼ N(0, C) =: ν Straightforward to show:

  • H

Hmζ, ζ ν(dζ) = tr( ˜ Hm),

  • H

Hmζ, C[Hmζ] ν(dζ) = tr( ˜ H2

m)

Finite dimensional algorithm: H. Avron and S. Toledo, Randomized algorithms for estimating the trace of an implicit symmetric positive semi-definite matrix, Journal of the ACM, 2011.

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 23 / 41

slide-31
SLIDE 31

Risk-averse optimal control with quadraticized objective

min

z∈Z

1 2Qu − ¯ q2

2+ 1

2ntr

ntr

  • j=1

ζj, ηj+ β 2

  • gm( ¯

m), C[gm( ¯ m)]+ 1 2ntr

ntr

  • j=1
  • C1/2ηj
  • 2

with gm( ¯ m) =e ¯

m∇u · ∇p

ηj = e ¯

m(ζj∇u · ∇p + ∇υj · ∇p + ∇u · ∇ρj)

  • Hmζj

j ∈ {1, . . . , ntr} where −∇ · (e ¯

m∇u) = N i=1 zifi

−∇ · (e ¯

m∇p) = −Q∗(Qu − ¯

q) −∇ · (e ¯

m∇υj) = ∇ · −(ζje ¯ m∇u)

−∇ · (e ¯

m∇ρj) = −Q∗Qυj + ∇ · (ζje ¯ m∇p)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 24 / 41

slide-32
SLIDE 32

Risk-averse optimal control with quadraticized objective

min

z∈Z

1 2Qu − ¯ q2

2+ 1

2ntr

ntr

  • j=1

ζj, ηj+ β 2

  • gm( ¯

m), C[gm( ¯ m)]+ 1 2ntr

ntr

  • j=1
  • C1/2ηj
  • 2

with gm( ¯ m) =e ¯

m∇u · ∇p

ηj = e ¯

m(ζj∇u · ∇p + ∇υj · ∇p + ∇u · ∇ρj)

  • Hmζj

j ∈ {1, . . . , ntr} where −∇ · (e ¯

m∇u) = N i=1 zifi

−∇ · (e ¯

m∇p) = −Q∗(Qu − ¯

q) −∇ · (e ¯

m∇υj) = ∇ · −(ζje ¯ m∇u)

−∇ · (e ¯

m∇ρj) = −Q∗Qυj + ∇ · (ζje ¯ m∇p)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 24 / 41

slide-33
SLIDE 33

Risk-averse optimal control with quadraticized objective

min

z∈Z

1 2Qu − ¯ q2

2+ 1

2ntr

ntr

  • j=1

ζj, ηj+ β 2

  • gm( ¯

m), C[gm( ¯ m)]+ 1 2ntr

ntr

  • j=1
  • C1/2ηj
  • 2

with gm( ¯ m) =e ¯

m∇u · ∇p

ηj = e ¯

m(ζj∇u · ∇p + ∇υj · ∇p + ∇u · ∇ρj)

  • Hmζj

j ∈ {1, . . . , ntr} where −∇ · (e ¯

m∇u) = N i=1 zifi

−∇ · (e ¯

m∇p) = −Q∗(Qu − ¯

q) −∇ · (e ¯

m∇υj) = ∇ · −(ζje ¯ m∇u)

−∇ · (e ¯

m∇ρj) = −Q∗Qυj + ∇ · (ζje ¯ m∇p)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 24 / 41

slide-34
SLIDE 34

Risk-averse optimal control with quadraticized objective

min

z∈Z

1 2Qu − ¯ q2

2+ 1

2ntr

ntr

  • j=1

ζj, ηj+ β 2

  • gm( ¯

m), C[gm( ¯ m)]+ 1 2ntr

ntr

  • j=1
  • C1/2ηj
  • 2

with gm( ¯ m) =e ¯

m∇u · ∇p

ηj = e ¯

m(ζj∇u · ∇p + ∇υj · ∇p + ∇u · ∇ρj)

  • Hmζj

j ∈ {1, . . . , ntr} where −∇ · (e ¯

m∇u) = N i=1 zifi

−∇ · (e ¯

m∇p) = −Q∗(Qu − ¯

q) −∇ · (e ¯

m∇υj) = ∇ · −(ζje ¯ m∇u)

−∇ · (e ¯

m∇ρj) = −Q∗Qυj + ∇ · (ζje ¯ m∇p)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 24 / 41

slide-35
SLIDE 35

Risk-averse optimal control with quadraticized objective

min

z∈Z

1 2Qu − ¯ q2

2+ 1

2ntr

ntr

  • j=1

ζj, ηj+ β 2

  • gm( ¯

m), C[gm( ¯ m)]+ 1 2ntr

ntr

  • j=1
  • C1/2ηj
  • 2

with gm( ¯ m) =e ¯

m∇u · ∇p

ηj = e ¯

m(ζj∇u · ∇p + ∇υj · ∇p + ∇u · ∇ρj)

  • Hmζj

j ∈ {1, . . . , ntr} where −∇ · (e ¯

m∇u) = N i=1 zifi

−∇ · (e ¯

m∇p) = −Q∗(Qu − ¯

q) −∇ · (e ¯

m∇υj) = ∇ · −(ζje ¯ m∇u)

−∇ · (e ¯

m∇ρj) = −Q∗Qυj + ∇ · (ζje ¯ m∇p)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 24 / 41

slide-36
SLIDE 36

Lagrangian for risk-averse optimal control with qquad

L (z, u, p,{υj}ntr

j=1, {ρj}ntr j=1, u⋆, p⋆, {υ⋆ j}ntr j=1, {ρ⋆ j}ntr j=1)

= 1 2 Qu − ¯ q2

2

+ 1 2ntr

ntr

  • j=1
  • ζj,
  • e ¯

m(ζj∇u · ∇p + ∇υj · ∇p + ∇u · ∇ρj)

  • + β

2

  • e ¯

m∇u · ∇p, C[e ¯ m∇u · ∇p]

  • +

β 4ntr

ntr

  • j=1
  • C1/2

e ¯

m(ζj∇u · ∇p + ∇υj · ∇p + ∇u · ∇ρj)

  • 2

+

  • e ¯

m∇u, ∇u⋆

− N

i=1 zifi, u⋆

+

  • e ¯

m∇p, ∇p⋆

+ Q∗(Qu − ¯ q), p⋆ +

ntr

  • j=1
  • e ¯

m∇υj, ∇υ⋆ j

  • +
  • ζje ¯

m∇u, ∇υ⋆ j

  • +

ntr

  • j=1
  • e ¯

m∇ρj, ∇ρ⋆ j

  • +
  • Q∗Qυj, ρ⋆

j

  • +
  • ζje ¯

m∇p, ∇ρ⋆ j

  • Omar Ghattas (ICES, UT Austin)

Optimal control under uncertainty Mar 24, 2017 25 / 41

slide-37
SLIDE 37

Adjoint & gradient for risk-averse optimal control w/qquad

Adjoint problem for qquad approximation

− ∇ · (e ¯

m∇ρ⋆ j) = b(j) 1

j ∈ {1, . . . , ntr} − ∇ · (e ¯

m∇υ⋆ j) + Q∗Qρ⋆ j = b(j) 2

j ∈ {1, . . . , ntr} − ∇ · (e ¯

m∇p⋆) − ntr

  • j=1

∇ · (ζje ¯

m∇ρ⋆ j) = b3

− ∇ · (e ¯

m∇u⋆) + Q∗Qp⋆ − ntr

  • j=1

∇ · (ζje ¯

m∇υ⋆ j) = b4

Gradient for qquad approximation

∂L ∂zj = γzj − fj, u⋆, j = 1, . . . , nc Cost of objective = 2 + 2ntr PDE solves; cost of gradient = 2 + 2ntr PDE solves

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 26 / 41

slide-38
SLIDE 38

Risk-averse optimal control with quadraticized objective

−20 20 40 60 80 100 120 140 0.01 0.02 0.03 0.04 distribution Θ(z0, ·) Θlin(z0, ·) Θquad(z0, ·)

initial (suboptimal) control z0

  • distrib. of exact & approx objectives at z0

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 27 / 41

slide-39
SLIDE 39

Risk-averse optimal control with quadraticized objective

−20 20 40 60 80 100 120 140 0.01 0.02 0.03 0.04 distribution Θ(z0, ·) Θlin(z0, ·) Θquad(z0, ·)

initial (suboptimal) control z0

  • distrib. of exact & approx objectives at z0

2 4 6 8 10 12 14 16 0.2 0.4 0.6 distribution Θ(zopt

quad, ·)

Θquad(zopt

quad, ·)

  • ptimal control zopt

quad based on qquad

  • distrib. of exact & approx objectives at zopt

quad Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 27 / 41

slide-40
SLIDE 40

Effect of number of trace estimator vectors on distribution

  • f control objective evaluated at optimal controls

1 2 3 4 5 6 7 8 0.2 0.4 value of control objective distribution ntr = 5 ntr = 20 ntr = 40 ntr = 60

Optimal controls zopt

quad computed for each value of trace estimator using quadratic

approximation of control objective, qquad Each curve based on 10,000 samples of distribution of q(zopt

quad, m)

(control objective evaluated at optimal control)

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 28 / 41

slide-41
SLIDE 41

Comparison of distribution of control objective for optimal controls based on linearized and quadraticized objective

1 2 3 4 5 6 7 8 0.2 0.4 0.6 value of control objective distribution Θ(zopt

quad, ·)

Θ(zopt

lin , ·)

1000

5000 10000 1 2 3 4 5

sample size sample mean

1000 5000 10000 5 10 15

sample size sample variance

Comparison of the distributions of q(zopt

lin , m) with q(zopt quad, m)

β = 1, γ = 10−5 and ntr = 40 trace estimation vectors KDE results are based on 10,000 samples Inserts show Monte Carlo sample convergence for mean and variance

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 29 / 41

slide-42
SLIDE 42

Optimal control of stochastic turbulent combustion

Turbulent non-premixed Hydrogen-Air slot-jet flame, with co-flow Sources of uncertainty

Inputs: inlet flow conditions Model parameters: turbulence models, chemical kinetics Model inadequacy: chemical kinetics, flamelet model, turbulent transport

Quantities of interest

Heat released over a specified distance Net rate of NOx production

Design Problem

Objective: Maximize heat release with constraint on NOx production Control: Inlet profiles and optionally volumetric heat sink

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 30 / 41

slide-43
SLIDE 43

UQ challenges

Forward solves are relatively expensive

Naive sampling methods are cost prohibitive

Model inadequacies are significant, complex, and deeply embedded in physical model

Requires physics-based approaches to inadequacy Leads to high/infinite dimensional stochastic models

Multiple QoIs, some related to rare events

Tail probabilities challenging to compute Realistic risk measures often result in non-smooth

  • bjective functions

These features are common in multiscale and multiphysics problems across scientific disciplines They combine to make UQ for inference, prediction, and

  • ptimization extremely challenging

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 31 / 41

slide-44
SLIDE 44

Extensions (Peng Chen, Umberto Villa)

Optimal control of turbulent jet governed by a turbulence model with embedded stochastic inadequacy model Randomized SVD to compute trace of ˜ H, combined with eigenvalue sensitivities to computing gradient Monte Carlo corrections to quadratic approximation: Use of Taylor approximation as a control variate

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 32 / 41

slide-45
SLIDE 45

Dirichlet boundary control for turbulent jet flow: problem

Control is horizontal velocity profile at inlet boundary ΓI Objective is to maximize jet width at Γ0 Constraint on inlet momentum:

  • ΓI (u · n)2ds = MI

Random input is an inadequacy field for turbulent viscosity (5151 dimensions)

Left: sketch of the physical domain of the turbulence jet flow, with inlet boundary ΓI, outlet boundary ΓO, top and bottom wall ΓW , the center axis ΓC, and the cross-section Γ0. The computational domain D is the top part of the physical domain. Right: two random samples drawn from the Gaussian measure N(0, C) with C = (−α1△ + α2I)−2 where α1 = α2 = 0.5.

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 33 / 41

slide-46
SLIDE 46

Dirichlet boundary control for turbulent jet flow: model

−∇ ·

  • (ν + γνt,0)
  • ∇u + ∇u⊤

+ (u · ∇) u + ∇p = 0, in D, ∇ · u = 0, in D, −∇ · ((ν + (γ + em) νt,0) ∇γ) + u · ∇γ − 1 2 u · e1 x1 + bγ = 0, in D, σn(u) · τ = 0, u · n + χW z = 0,

  • n ΓI,

σn(u) · n = 0, u · τ = 0,

  • n ΓO ∪ ΓW ,

σn(u) · τ = 0, u · n = 0,

  • n ΓC,

γ − γ0 = 0,

  • n ΓI ∪ ΓW ,

σγ

n(γ) · n = 0,

  • n ΓO ∪ ΓC.

σn(u) = (ν + γνt,0)

  • ∇u + ∇u⊤

· n σγ

n(γ) = (ν + (γ + em)νt,0)∇γ · n

νt,0 = C √ M(x1 + aW)1/2 with M =

  • ΓI

||udns||2ds γ0 = 0.5 − 0.5 tanh

  • 5

30 − x1 30

  • (x2 − 1 − 0.5x1)
  • Omar Ghattas (ICES, UT Austin)

Optimal control under uncertainty Mar 24, 2017 34 / 41

slide-47
SLIDE 47

Dirichlet boundary control for turbulent jet flow: results

The velocity field of the turbulent jet flow extracted from the DNS data (top) and

  • btained from the optimal control with quadratic approximation with variance reduction

(bottom).

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 35 / 41

slide-48
SLIDE 48

Computing the trace of ˜ H via randomized SVD

Large number of randomized trace vectors might be required to compute tr( ˜ H) Randomized SVD estimates trace at cost

  • f 2r products of ˜

H with random vectors (r is numerical rank of ˜ H) Resulting cost is 2r incremental forward/adjoint solves and 4r Poisson solves (Steps 2 and 4) Covariance operator is compact and Hessian is often low-rank (QoI is sensitive to limited number of modes) so composition is low-rank

Randomized SVD (double pass algorithm)

1

Generate i.i.d. Gaussian matrix R ∈ Rn×r with r = numerical rank of ˜ H (r ≪ n)

2

Form Y = ˜ HR

3

Compute Q = orthonormal basis for Y

4

Define B ∈ Rr×r := QT ˜ HQ

5

Decompose B = ZΛZT

6

Low-rank approximation: ˜ H ≈ V ΛV T , where V ∈ Rn×r := QZ

7

Trace estimation: tr( ˜ H) ≈ tr(B)

Thus r ≪ n, independent of parameter dimension n, and with high probability |tr( ˜ H) − tr(B)| ≤

  • r<i≤n

|λi( ˜ H)| Quadratic-based approximations of E[q] and var[q] require 4r linearized forward solves (small multiple of nonlinear forward solve, for highly nonlinear forward problems) Computing gradient of tr( ˜ H) wrt controls via eigenvalue sensitivity

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 36 / 41

slide-49
SLIDE 49

Control of turbulent jet: Computing tr( ˜ H) via randomized SVD vs. randomized trace estimator

20 40 60 80 100

N

10-8 10-7 10-6 10-5 10-4 10-3

λ + λ− error1 error2

20 40 60 80 100

N

10-8 10-7 10-6 10-5 10-4 10-3

λ + λ− error1 error2

The decay of the generalized eigenvalues (λ+ and λ−) of ˜ H and the randomized trace estimator error (error1) and the randomized SVD error (error2) with the number terms

  • N. Left: at the initial control. Right: at the optimal control.

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 37 / 41

slide-50
SLIDE 50

Computation of the derivative of ˆ T( ˜ Hm) wrt control z

Recall the trace estimation ˆ T( ˜ Hm) =

N

  • j=1

λj( ˜ Hm) and ˆ T( ˜ H2

m) = N

  • j=1

λ2

j( ˜

Hm) Eigenvalue sensitivity equation: find (λ′

j, φ′ j) s.t. ∀φ ∈ M,

φ, H′

m(m0)ψj + φ, Hm(m0)ψ′ j = λ′ jφ, C−1ψj + λjφ, C−1ψ′ j

where ′ represents the derivative w.r.t. z Setting φ = ψj, and noting ψj, C−1ψj = 1 as well as symmetry of Hm(m0) and C−1, we obtain the expression for the eigenvalue sensitivity in terms of the Hessian derivative, λ′

j = ψj, H′ m(m0)ψj,

which allows us to compute the gradient of tr( ˜ Hm) at any optimization iteration at a cost of N + p pairs of linearized forward/adjoint PDE solves Similarly, the gradient of tr( ˜ H2

m) can be computed with

(λ2

j)′ = 2λjλ′ j = 2ψj, Hm(m0) ψj ψj, H′ m(m0) ψj.

where we have used the fact λj = ψj, Hm(m0) ψj.

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 38 / 41

slide-51
SLIDE 51

Taylor approximation as a variance reduction

Statistics computed by Taylor approximations may be biased Use Monte Carlo quadrature to correct Taylor approximation, e.g., E[Q] = E[Qlin] + E[Q − Qlin

  • Y

] ≈ Q(m0) + ˆ Y

  • MC estimator

Similar MC correction for quadratic approximation Qquad Mean squared error (MSE) of MC estimate of E[Q] and E[Y ] MSE(Q) ≍ 1 N Var[Q] vs. MSE(Y ) ≍ 1 N Var[Y ] A much smaller number of MC samples is required for E[Y ] as Var[Y ] ≪ Var[Q] provided Qlin is a good approximation of (i.e., highly correlated to) Q see Multifidelity Monte Carlo work of Karen Willcox et al.

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 39 / 41

slide-52
SLIDE 52

Turbulent jet model with embedded inadequacy

Estimation of E[q] using Monte Carlo and variance reduction Monte Carlo where E[qlin] = q( ¯ m)= 1.04854; E[qquad] = q( ¯ m) + 1

2tr

˜ Hm

  • ≈ 1.04825.
  • q

E[qlin] + Ylin E[qquad] + Yquad E[q] 1.04800 1.04817 1.04818 MSE 1.3168E-06 5.1655E-08 1.1868E-08 N 100 100 100

Sample variance and number N of MC samples to obtain an accuracy τ = 10−4

q Ylin Yquad variance 1.3168E-04 5.1655E-06 1.1868E-06 N 13168 516 119 Variance reduction with low-rank-Hessian-based quadratic approximation of parameter-to-QoI map leads to 100X reduction in # of Monte Carlo samples for 5151-dimensional uncertain inputs.

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 40 / 41

slide-53
SLIDE 53

Conclusions

Scalable computational framework based on Taylor approximation and variance reduction for PDE-constrained optimal mean-var control Cost of objective/gradient evaluation measured in PDE solves is independent

  • f parameter dimension when covariance preconditioned Hessian is a compact
  • perator

Randomized SVD is more efficient than randomized trace estimator; gradient computed via eigenvalue sensitivity at no additional cost in PDE solves Taylor approximation used as control variate leads to 100X reduction in MC samples and unbiased estimate

  • A. Alexanderian, N. Petra, G. Stadler, O. Ghattas, Mean-variance risk-averse optimal

control of systems governed by PDEs with random coefficient fields using local

  • approximations. arXiv:1602.07592, 2016.
  • P. Chen, U. Villa, O. Ghattas, Taylor approximation and variance reduction for

PDE-constrained optimal control under uncertainty. Preprint, 2017. Optimization under uncertainty research supported by DARPA EQUiPS (turbulent combustion), DOE MMICCs (subsurface flow), and DOE Extreme Scale UQ programs

Omar Ghattas (ICES, UT Austin) Optimal control under uncertainty Mar 24, 2017 41 / 41